代码分享|| 同一套代码和思路已经发了三篇，提供完整代码和数据的泛癌分析文献分享！

代码分享|| 同一套代码和思路已经发了三篇，提供完整代码和数据的泛癌分析文献分享！ by 新叶生物

最近大家做代码资源分享的公众号非常多，感觉这个赛道已经白热化，你的收藏夹里面是不是已经收藏了N多代码，但是都在吃灰？有个朋友告诉我他已经收藏了100多套了！！！后面我打算做视频讲解复现分析好了，虽然很多文献有代码有数据，但是对大部分没有生信基础的人来说，复现起来应该还是有很大的门槛！

今天分享的文献有一套代码，使用这套代码类似的思路差不多发了三篇文献啦！来看看吧！

提供代码的文献：

其他两篇文献如下：

H.-R. Feng, X.-N. Shen, X.-M. Zhu, W.-T. Zhong, D.-X. Zhu, J. Zhao, Y.-J. Chen, F. Shen, K. Liu, L. Liang, Unveiling Major Histocompatibility Complex-mediated Pan-Cancer Immune Features by Integrated Single-Cell and Bulk RNA Sequencing, Cancer Letters, https://doi.org/10.1016/j.canlet.2024.217062

Wu, Zhengquan et al. “A transcriptomic pan-cancer signature for survival prognostication and prediction of immunotherapy response based on endothelial senescence.” Journal of biomedical science vol. 30,1 21. 28 Mar. 2023, doi:10.1186/s12929-023-00915-5

代码和数据

https://figshare.com/articles/online_resource/Integrated_analysis_of_single-cell_and_bulk_RNA_sequencing_data_reveals_a_pan-cancer_stemness_signature_predicting_immunotherapy_response/17654633

每一幅都有对应的代码：

文献主要内容

看看思路，你是不是也能做一篇出来~

文章做的癌细胞干性泛癌分析，基本假设：免疫检查点抑制剂（ICI）是癌症治疗的一个突破，但只有有限的一部分患者从中受益。癌症干细胞性可能是ICI耐药性的潜在罪魁祸首，但目前还缺乏直接的临床证据。作者收集了来自ICI治疗的患者的公开scRNA-Seq数据集，以阐明癌症干性和ICI反应之间的关系。

所以本文还为大家提供了ICI治疗的单细胞数据队列，可以利用起来，为自己的科学研究课题所用！

作者使用大规模泛癌数据：包括34个scRNA-Seq数据集、TCGA数据集、10 ICI转录组队列，开发和验证了一个新的 stemness signature (Stem.Sig) 。此外，Stem.Sig中的基因（EMC3, BECN1, VPS35, PCBP2, VPS29, PSMF1, GCLC, KXD1, SPRR1B）使用 17个筛选潜在免疫治疗靶点的CRISPR数据集进行了治疗方面的探索。

ICI治疗的单细胞数据队列

用于开发Stem.Sig的scRNA数据集

免疫治疗队列

CRISPR数据集

免疫治疗的预测基因表达特征列表

图1：识别和验证癌细胞干细胞性和ICI结果之间的负相关

图1的代码：

###==Figure 1===###

rm (list=ls())

gc()

#====loading data===####
library(CytoTRACE)
library(Seurat)
library(dplyr)
library(ggpubr)
source('Code/plotCytoTRACE.v2.0.R') #modified plotCytoTRACE function implemented in CytoTRACE
source('Code/R_rainclouds.R')


####====Fig. 1A and 1B===####
#GSE115978
cohort <- "SKCM_GSE115978_aPD1" 
load(paste0("data/",cohort,'/scRNA.Rdata'))
load(paste0('data/',cohort,'/phenotype.Rdata'))

cell <- phenotype[phenotype$cellType == 'Malignant cells' ,'Cell'] 
response <- phenotype[phenotype$cellType == 'Malignant cells' ,'Response']
names(response) <- phenotype[phenotype$cellType == 'Malignant cells' ,'Cell']

scRNA <- scRNA[,colnames(scRNA) %in% cell]
scRNA <- as.data.frame(scRNA)
results <- CytoTRACE(scRNA)
  
Fig.1A <- plotCytoTRACE.v2.0(results,response)



boxdata <- data.frame(Cell = names(results$CytoTRACE),cytoTRACE=results$CytoTRACE)
boxdata <- left_join(boxdata,phenotype[,c('Cell','Response')],by='Cell')
boxdata$Response <- factor(boxdata$Response,levels = c("NR","TN"))


Fig.1B <- ggplot(boxdata, aes(x = Response, y = cytoTRACE, fill = Response)) +
  geom_flat_violin(aes(fill = Response),position = position_nudge(x = 0.1, y = 0),  trim = TRUE, alpha = .5, colour = NA)+
  geom_point(aes(x = .55, y = cytoTRACE, colour = Response),position = position_jitter(width = .05), size = 1, shape = 20)+
  geom_boxplot(aes(x = Response, y = cytoTRACE, fill = Response),outlier.shape = NA, alpha = .5, width = .1, colour = "black")+
  scale_colour_manual(values =  c("#E7B800","#2E9FDF"))+
  scale_fill_manual(values = c("#E7B800","#2E9FDF"))+
  theme_classic()+
  geom_signif(comparisons = list(c("NR","TN")),map_signif_level = TRUE,test = "wilcox.test")+
  theme(axis.title = element_text(size=15),
        axis.title.x = element_blank(),
        axis.text = element_text(size=12))+
  # ylim(0,1)+
  ylab('CytoTRACE scores')

Fig.1B

####====Fig. 1C and 1D===####

cohort <- "BCC_GSE123813_aPD1"
load(paste0("data/",cohort,'/scRNA.Rdata')) 
load(paste0("data/",cohort,'/phenotype.Rdata'))

cell <- phenotype[phenotype$cellType == 'Malignant cells' ,'Cell'] 
response <- phenotype[phenotype$cellType == 'Malignant cells' ,'Response']
names(response) <- phenotype[phenotype$cellType == 'Malignant cells' ,'Cell']

scRNA <- scRNA[,colnames(scRNA) %in% cell]
scRNA <- as.data.frame(scRNA)
results <- CytoTRACE(scRNA)


Fig.1C <- plotCytoTRACE.v2.0(results,response)

boxdata <- data.frame(Cell = names(results$CytoTRACE),cytoTRACE=results$CytoTRACE)
boxdata <- left_join(boxdata,phenotype[,c('Cell','Response')],by='Cell')
boxdata$Response <- factor(boxdata$Response,levels = c('NR','R'))


Fig.1D <- ggplot(boxdata, aes(x = Response, y = cytoTRACE, fill = Response)) +
  geom_flat_violin(aes(fill = Response),position = position_nudge(x = 0.1, y = 0),  trim = TRUE, alpha = .5, colour = NA)+
  geom_point(aes(x = .55, y = cytoTRACE, colour = Response),position = position_jitter(width = .05), size = 1, shape = 20)+
  geom_boxplot(aes(x = Response, y = cytoTRACE, fill = Response),outlier.shape = NA, alpha = .5, width = .1, colour = "black")+
  scale_colour_manual(values =  c("#E7B800","#2E9FDF"))+
  scale_fill_manual(values = c("#E7B800","#2E9FDF"))+
  theme_classic()+
  geom_signif(comparisons = list(c("NR","R")),map_signif_level = TRUE,test = "wilcox.test")+
  theme(axis.title = element_text(size=15),
        axis.title.x = element_blank(),
        axis.text = element_text(size=12))+
  ylab('CytoTRACE scores')

Fig.1D

图2：stemness signature的开发和描述

图2的代码：

####===Figure 2===####
rm (list=ls())
gc()


####===Fig. 2A===####
library(circlize)
library(RColorBrewer)
library(png)
library(graphics)
library(circlize)
library(RColorBrewer)
library(stringr)
library(ComplexHeatmap)
library(grImport2)
library(gridBase)
library(ggplot2)

dir <- data.table::fread('data/dir/dir.csv')
circos.clear()
Set3 <- brewer.pal(12,"Set3")
Set2 <- brewer.pal(8,"Set2")
col <- c(Set3,Set2)

sectors = dir$dataset
names(sectors) <- sectors

dataset <- data.frame(dataset = sectors)
rownames(dataset) <- sectors

col_cancer <- c(col[1], rep(col[2],2), col[3],col[4],rep(col[5],11),col[6],col[7],col[8],rep(col[9],2),col[10],col[11],rep(col[12],4),col[13],rep(col[14],2),
                rep(col[15],2),col[16],col[19])
names(col_cancer) <- sectors



image = 'data/png/venn.png'
image = as.raster(readPNG(image))

circos.clear()
circos.par(gap.degree = 2, cell.padding = c(0, 0, 0, 0),
           track.margin = c(0.01, 0.01))
circos.heatmap(dataset, split = dataset$dataset, col = col_cancer, track.height = 0.02,rownames.side = 'outside')

circos.track(ylim = c(0, 1),track.height=0.15,bg.border = "#EEEEEE", panel.fun = function(x, y) {
  circos.raster(image, CELL_META$xcenter, CELL_META$ycenter, 
                width = "0.7cm",
                facing = "inside")
  
})

circos.track(ylim = c(0, 1), sectors = sectors,
             bg.col = "#8190A5", bg.border = "#EEEEEE" , track.height = 0.25)

circos.trackText(x = rep(0.5, 34), y = rep(0.5, 34),
                 labels = paste0(rep('G',34), 1:34),
                 cex = 0.5, sectors = sectors, col = "white", font = 2, facing = "clockwise",
                 niceFacing=T)

draw.sector(center = c(0, 0), start.degree = 0, end.degree = 360,
            rou1 = 0.25, col = "#4D6381", border = "#EEEEEE")

text(0,0,"Stem.sig",col = 'white',cex = 2,font = 2)


cancer <- unique(col_cancer)
names(cancer) <- unique(dir$cancer)
lgd_cancer = Legend(title = "Cancer", at = names(cancer),
                 legend_gp = gpar(fill = cancer))
draw(lgd_cancer, x = unit(1.2, "snpc"), just = "left")



####===Fig. 2B===####
library(clusterProfiler)
library(org.Hs.eg.db)
library(ReactomePA)

load('data/sig/Stem.Sig.Rdata')

sig <- Stem.Sig

sig_ID <- bitr(sig,fromType = 'SYMBOL',toType = 'ENTREZID',OrgDb = "org.Hs.eg.db")

eReac <- enrichPathway(gene = sig_ID$ENTREZID,
                       organism = 'human',
                       pvalueCutoff = 0.05)

bardata <- eReac@result[1:20,]
bardata <- bardata[order(bardata$qvalue),]
bardata <- bardata[order(bardata$Description),]
bardata$Description <- factor(bardata$Description,levels=rev(bardata$Description))
barplot <- ggplot(bardata,aes(y=Description,x=Count))+geom_bar(stat = "identity",aes(fill=qvalue))+scale_fill_gradientn(colours =  c("#4D6381", "#BAC2CC"))+theme_bw()

eReacx <- setReadable(eReac, 'org.Hs.eg.db', 'ENTREZID')

cnetplot <-  cnetplot(eReacx, circular = TRUE, colorEdge = TRUE,categorySize="pvalue",showCategory = 20,layout = 'kk') 

barplot 
cnetplot

图3：使用泛癌TCGA队列分析Stem.Sig与免疫耐药之间的潜在联系

图4：使用Stem.Sig预测ICI结果

图5：Stem.Sig的AUC与其他预测基因特征的比较

图6：利用CRISPR筛选数据探索来自Stem.Sig的潜在治疗靶点

Cite：

Zhang, Zhen et al. “Integrated analysis of single-cell and bulk RNA sequencing data reveals a pan-cancer stemness signature predicting immunotherapy response.” Genome medicine vol. 14,1 45. 29 Apr. 2022, doi:10.1186/s13073-022-01050-w

往期精彩

代码分享|| IF10+的 python版本的单细胞数据分析，图的颜值超高！

代码分享|| 学习一个单细胞与空转数据库 ssREAD（阿尔兹海默症）网站的构建：提供前端+后端代码，文章图表分析代码

代码分享|| 单细胞亚群注释软件你还有不知道的吗？快来看看这个发表在Nature Communications的工具！

代码分享|| Nature肝脏疾病进展单细胞分析超详细代码：上游分析+下游，详细到fq文件重命名，数据预处理等！快上车！

代码分享|| 单细胞上车尤未晚！囊括了单细胞分析方方面面包括大类注释、亚群分析、细胞通讯、轨迹分析、生存、多因素cox分析等