ixxmu / mp_duty

抓取网络文章到github issues保存
https://archives.duty-machine.now.sh/
106 stars 30 forks source link

泛癌的免疫浸润可视化(附小豆花讲课视频) #2922

Closed ixxmu closed 1 year ago

ixxmu commented 1 year ago

https://mp.weixin.qq.com/s/NYnjfWPvjsQVlGnZlFxryg

ixxmu commented 1 year ago

泛癌的免疫浸润可视化(附小豆花讲课视频) by 生信星球

 今天是生信星球陪你的第880天

   大神一句话,菜鸟跑半年。我不是大神,但我可以缩短你走弯路的半年~

   就像歌儿唱的那样,如果你不知道该往哪儿走,就留在这学点生信好不好~

   这里有豆豆和花花的学习历程,从新手到进阶,生信路上有你有我!

好久不见!我的产假还剩一个多月,小娃娃两个半月了,天使宝宝,省爹省妈,所以我也有些时间了,陆续开始学习和备课,整理资料咯。

从Timer2.0的主页发现了一个优秀的资源:

http://timer.comp-genomics.org/infiltration_estimation_for_tcga.csv.gz

这里面是各种方法计算出的TCGA免疫浸润结果。有现成的结果,做图就完全是R语言技巧的输出咯。

可以做接下来的几件事:

1.将某一方法计算出来的免疫浸润结果拆分成单独的数据框。

2.绘制箱线图,展示某种免疫细胞在所有癌症中,tumor和normal的丰度比较。

3.绘制箱线图,展示某种免疫细胞在所有癌症中,tumor样本的丰度,并按照中位数从小到大排好顺序。

你可以先自己试试,再看我写的答案。

1.拆分

rm(list = ls())
a = read.csv("infiltration_estimation_for_tcga.csv.gz",
             row.names = 1,
             check.names = F)
head(colnames(a))
## [1] "B cell_TIMER"                 "T cell CD4+_TIMER"           
## [3] "T cell CD8+_TIMER"            "Neutrophil_TIMER"            
## [5] "Macrophage_TIMER"             "Myeloid dendritic cell_TIMER"
library(stringr)
b = str_split(colnames(a),"_",simplify = T) %>% as.data.frame()
colnames(b) = c("cell","method")
table(b$method)
#
##     CIBERSORT CIBERSORT-ABS          EPIC    MCPCOUNTER     QUANTISEQ 
##            22            22             8            11            11 
##         TIMER         XCELL 
##             6            39
k = b$method=="CIBERSORT";table(k)
## k
## FALSE  TRUE 
##    97    22
ciber = a[,k]
colnames(ciber) = str_remove(colnames(ciber),"_CIBERSORT")
ciber[1:4,1:4]
##                 B cell naive B cell memory B cell plasma T cell CD8+
## TCGA-OR-A5J1-01  0.002937309   0.002282572    0.00000000  0.11294112
## TCGA-OR-A5J2-01  0.046380291   0.000000000    0.15149520  0.07351904
## TCGA-OR-A5J3-01  0.061034875   0.000000000    0.19053762  0.01649453
## TCGA-OR-A5J5-01  0.253704390   0.000000000    0.03964469  0.03734717

2.泛癌的tumor和normal某种细胞丰度比较

需要泛癌临床信息,才能知道每个样本属于哪种癌症。可以从xena下载

https://tcga-pancan-atlas-hub.s3.us-east-1.amazonaws.com/download/Survival_SupplementalTable_S1_20171025_xena_sp

#整理临床信息,使其与免疫细胞丰度矩阵行名一致
m = read.delim("Survival_SupplementalTable_S1_20171025_xena_sp",row.names = 1)
y = intersect(rownames(ciber) , rownames(m));length(y)
## [1] 11031
ciber = ciber[y,]
m = m[y,]
identical(rownames(ciber),rownames(m))
## [1] TRUE
#做图
library(ggplot2)
ciber$type = m$cancer.type.abbreviation
ciber$Group =  ifelse(str_sub(rownames(ciber),14,15)<10,"tumor","normal")
ggplot(ciber,aes(x = type,y = `B cell naive`))+
  geom_boxplot(aes(fill = Group),alpha = 0.7)+
  scale_fill_manual(values = c("navy","firebrick3"))+
  theme_bw()+
  theme(axis.text.x = element_text(angle = 90,vjust = 0.5))

3.某种免疫细胞的丰度在所有tumor中的比较

随便拿一个细胞B cell naive来画。

library(dplyr)
dat = ciber[ciber$Group=="tumor",]
#利用因子,调整横坐标顺序
dat2 = group_by(dat,type) %>% 
  summarise(median = median(`B cell naive`)) %>% 
  arrange(median)
dat$type=factor(dat$type,levels = dat2$type)
library(ggplot2)
ggplot(dat,aes(type,`B cell naive`))+
  geom_boxplot(outlier.size = 0.5)+
  theme_bw()+
  theme(axis.text.x = element_text(angle = 90,vjust = 0.5))
具体的代码让小豆花老师讲讲吧。
如果因为代码看不懂,而跟不上正文的节奏,可以来找我,系统学习。以下课程都是循环开课。下一期的时间,点进去咨询微信咯
生信零基础入门学习小组
生信入门班(四周线上直播课)
数据挖掘班(医生/医学生首选)