免疫浸润分析大盘点,解决小鼠和人的免疫浸润分析难点
Hi,大家好,我是晨曦
这期的晨曦碎碎念,我们来展开一个新话题——免疫浸润分析
这期的话题其实也是小伙伴很早以前就跟晨曦说想要一个系统教程的话题,正好晨曦最近对于免疫浸润分析也做了一些总结,所以本期推文也可以算是晨曦的一个学习总结
欢迎各位小伙伴在评论区进行交流哦~
那么,我们就开始吧

晨曦碎碎念系列传送门
免疫浸润分析的相关背景知识
首先,让我们来给免疫浸润分析下一个宏观的定义,来回答一下究竟什么是免疫浸润分析呢?
其实我们可以很简单的理解,免疫浸润分析其实就是为了明确免疫细胞在人体微环境内的组成情况,进而明确究竟是哪种免疫细胞在疾病的发生发展中产生了重要作用,这就是免疫浸润分析
好了,了解了概念,我们下面再来提出一个问题:我们通过什么算法来明确免疫细胞的组成呢?
大体的算法思想我们可以归结为两类:
1.marker gene;
2.反卷积思想
算法本身的逻辑不需要掌握,只需要知道两种算法得到的免疫浸润结果是不同的
反卷积算法得到的是每个样本免疫细胞的比例
marker gene算法得到的是每一个样本某个细胞的富集分数
也就是说尽管免疫浸润的方法有很多,但是基本上都是基于以上两种算法核心开发出来的,就好比汽车引擎的整体构造基本上是相似的,但是却延申出来了很多汽车品牌以及汽车类型
然后我们继续往下探索,那么目前来说免疫浸润分析都有哪些方法呢?
目前来说,我们用的最多的就是:CIBERSORT、ssGSEA、Xcell、MCPcounter、immundeconv等等
这里我们主要介绍的就是CIBERSORT、ssGSEA、以及MCPcounter,相信掌握了以上免疫浸润分析的方法,会让你的文章有一个更广阔的天空
CIBERSORT分析
首先,我们从实用的角度出发,算法的逻辑本质不做掌握,我们只需要知道一句话:作者为了让我们更方便的使用,把代码整理成了一个代码文件并写成了一个函数,我们只需要提供所需要的文件即可以完成分析
CIBERSORT基于线性支持向量回归(linear support vector regression)的原理进行反卷积分析
我们的输入数据需要两个文件:
1.表达矩阵文件
(观测为Gene symbol;变量为样本ID)
输入数据的要求:
1.不可以有缺失值和负数;
2.不可以取log(CIBERSORT会对于数据是否取过log有一个判断,简单来说就是如果取了log也会给你去掉,因为输入数据接受的是没有log的形式),这个在CIBERSORT函数的源代码中讲解的很清楚,至于为什么要取log,其实我们可以去看我们的参考数据集本身就是没有log的,我们是要通过表达量彼此映照来分析的,所以不需要log
#anti-log if max < 50 in mixture fileif(max(Y) < 50) {Y <- 2^Y}
3.芯片数据对应不同平台有不同的要求,Affymetrix芯片建议进行RMA标准化,illumina芯片建议用limma包进行标准化处理,Agilent芯片建议用limma包进行标准化处理
4.对于测序数据,不支持原始数据,而是建议使用TPM或者FPKM都可以
2.参考数据集文件
这里其实就是CIBERSORT分析比较自由的地方,也就是说,我们能否进行各种物种的免疫浸润分析,其实是依赖于参考数据集文件,我们可以把算法进行一个简单的理解,我们通过参考数据集然后比对我们的表达矩阵文件,进而得到每一个样本的免疫细胞比例,说到这里各位应该都理解了,我的参考数据集不同,就可以进行不同物种的免疫浸润分析
那么目前来说我们都有哪些参考数据集可以选择呢?
LM22——这个数据集来自Robust enumeration of cell subsets from tissue expression profiles,这篇文献对于这个数据集有一个比较清晰的解答:we designed and validated a leukocyte gene signature matrix, termed LM22. It contains 547 genes that distinguish 22 human hematopoietic cell phenotypes, including seven T cell types, naïve and memory B cells, plasma cells, NK cells, and myeloid subsets ,简单的理解,这个数据集的构建并没有划归恶性还是良性,所以肿瘤和非肿瘤,物种为人的时候我们都可以使用CIBERSORT算法进行免疫浸润分析(LM22作为参考数据集)
mice——这个数据集来自Inference of immune cell composition on the expression profiles of mouse tissue,这篇文献中的参考数据集是一个小鼠的免疫细胞数据集,可以替换LM22来进行肿瘤和非肿瘤,物种为小鼠的免疫浸润分析,当然这篇文献也提出了一种方法ImmunCC,有时间也会给大家进行介绍
掌握了上面这两个参考数据集以后,我们就可以正式进入代码部分的讲解,其实代码部分就很简单了,因为作者已经把复杂的地方都已经准备妥当,我们实际上只需要提供两个文件即可
#代码讲解——针对物种为人#### 1.准备工作####rm(list = ls())source("CIBERSORT.R") sig_matrix <- "LM22.txt"#参考数据集文件mixture_file = 'example.txt'#表达矩阵文件####2.免疫浸润分析####res_cibersort <- CIBERSORT(sig_matrix, mixture_file, perm=10, QN=TRUE)#perm参数为置换系数,理论上来说越大、运行越慢、结果越准确、适度即可(10-1000)#QN参数为数据类型,芯片为TRUE,测序为FALSE####3.额外话题——输入数据如何准备####example <- read.table("example.txt", header=T, sep="\t", check.names=F,row.names = 1)#上述示例数据包含四个样本的一个物种为人的芯片数据#导出函数write.table(cbind(rownames(example), example),"example2.txt",quote = F, sep = "\t", row.names=FALSE)#本质上就是对write.table函数的使用#导出文件之前需要把行名变成一列,不然后面就会有报错?write.table
上面的参考数据集我们选择是基于物种为人的情况,我们下面再来看一下基于物种为小鼠的时候应该如何进行操作
示例表达矩阵文件来自GSE7762,选择了其中六个样本
注意:下面表达矩阵文件名字虽然一样,但是物种是不一样的
#代码讲解——针对物种为小鼠#### 1.准备工作####rm(list = ls())source("CIBERSORT.R") sig_matrix <- "mice.txt"#参考数据集文件mixture_file = 'example.txt'#表达矩阵文件####2.免疫浸润分析####res_cibersort <- CIBERSORT(sig_matrix, mixture_file, perm=10, QN=TRUE)#perm参数为置换系数,理论上来说越大、运行越慢、结果越准确、适度即可(10-1000)#QN参数为数据类型,芯片为TRUE,测序为FALSE####3.额外话题——参考数据集如何准备####example <- read.csv("mice.csv", header=T,check.names=F,row.names = 1)#导出函数write.table(cbind(rownames(example), example),"mice.txt",quote = F, sep = "\t", row.names=FALSE)#本质上就是对write.table函数的使用#导出文件之前需要把行名变成一列,不然后面就会有报错(本质是因为CIBERSORT函数对于文件读取比较单一)?write.table
这样的话,通过替换参考数据集,我们就可以完成针对肿瘤/非肿瘤不同物种(人或者小鼠)的免疫浸润分析,至于如何进行相关的可视化,网络上都有很多相关的教程,完全可以运用我们的扒图技能来进行更加丰富的可视化展示
那么,下面我们再来扩展一下其它类型的免疫浸润分析,并且如果各位小伙伴对于三大免疫评分感兴趣,在评论区留言就可以后续安排更新哦
Ps:疯狂暗示QAQ
ssGSEA分析
ssGSEA分析即单样本基因富集分析,是GSEA方法的扩展,通过定义富集分数来提供输入数据内每一个样本中基因集的富集程度
简单理解,通过ssGSEA富集分析得到的结果是每一个样本在每一种免疫细胞内的富集分数,这个富集分数为作者定义的产物,至于什么是富集分数不需要理解的很透彻,只需要明白一个最浅显的道理,富集分数为1的大于富集分数为0.1的即可
要进行ssGSEA分析自然也就需要两个文件:
1.表达矩阵文件
(观测为Gene symbol;变量为样本ID)
这里对于输入数据的要求可以参考CIBERSROT的要求,基本上来说是没有问题的
2.参考数据集文件
通过与上面CIBERSORT进行对比,我们可以发现,ssGSEA的参考数据集仅仅只是告诉你哪些基因在哪些免疫细胞上,这就是marker gene算法的思想,因为没有表达量,所以就定义了一个富集分数,所以,既然我们已经知道了参考数据集的本质(格式不一定和上面的类似,只需要有对应关系的信息即可),那么自然就可以延申出来不同参考基因集所计算出来的ssGSEA
那么目前来说我们都有哪些参考数据集可以选择呢?
Pan-cancer:这个数据集来自Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade,为我们提供了泛癌的基因和免疫之间的关系,所以我们在进行肿瘤相关的、物种为人的ssGSEA分析的时候可以使用这个数据集
TISIDB数据库:这个包含一系列数据集来自TISIDB (hku.hk)数据库,我们通过前面的知识可以知道,ssGSEA的参考基因集数据只需要得到一组对应关系即可,所以只要我们定义一种对应关系即可以进行ssGSEA免疫浸润分析(其实就是单样本富集分析)
掌握了上面这两个参考数据集以后,我们就可以正式进入代码部分的讲解,其实代码部分就很简单了,因为R包已经把复杂的地方都已经准备妥当,我们实际上只需要提供两个文件即可
####1.准备工作#####针对TISIDB数据库library(data.table)library(GSVA)#这个包可以完成ssGSEA分析rm(list = ls())### 2.准备参考基因文件####cellMarker <- data.table::fread("cellmarker_TISIDB.txt",data.table = F)#其实只需要两列信息:Gene和免疫细胞的对应关系colnames(cellMarker)[3] <- "celltype"cellMarker <- cellMarker[,-c(1,4)]cellMarker$celltype <- gsub('[-]', ' ', cellMarker$celltype)list<- split(as.matrix(cellMarker)[,1], cellMarker[,2])#ssGSEA参考数据集要求的格式就是一个包含对应信息的list#save(cellMarker,file = "cellMarker_ssGSEA.Rdata")#### 3.准备表达量矩阵####### 行是基因,列是样本(形式参考CIBERSORT)expr <- data.table::fread("example.txt",data.table = F)rownames(expr) <- expr[,1]expr <- expr[,-1]expr <- as.matrix(expr)### 4.使用ssGSEA量化免疫浸润####gsva_data <- gsva(expr,list, method = "ssgsea")上面选择的是TISIDB数据库中的对T细胞介导敏感/不敏感的数据集,当然其中还有其它数据集,只需要获得相关的对应关系即可进行上述重复操作####1.准备工作#####针对Pancancer数据集library(data.table)library(GSVA)#这个包可以完成ssGSEA分析rm(list = ls())### 2.准备参考基因文件####cellMarker <- data.table::fread("cellMarker_pancancer.csv",data.table = F)#其实只需要两列信息:Gene和免疫细胞的对应关系colnames(cellMarker)[2] <- "celltype"list<- split(as.matrix(cellMarker)[,1], cellMarker[,2])#ssGSEA参考数据集要求的格式就是一个包含对应信息的list#save(cellMarker,file = "cellMarker_ssGSEA.Rdata")#### 3.准备表达量矩阵####### 行是基因,列是样本(形式参考CIBERSORT)expr <- data.table::fread("example.txt",data.table = F)rownames(expr) <- expr[,1]expr <- expr[,-1]expr <- as.matrix(expr)### 4.使用ssGSEA量化免疫浸润####gsva_data <- gsva(expr,list, method = "ssgsea")
至此,ssGSEA就给大家介绍到了这里,其实我们通过学习CIBERSORT和ssGSEA可以发现,我们大部分的时间只需要准备两个文件,而决定我们的分析是否可以应用于肿瘤或者非肿瘤、人或者小鼠的前提其实就是参考基因集的类型
前面介绍的CIBSERSORT提供了肿瘤和非肿瘤,物种为人或者小鼠的免疫浸润分析
ssGSEA因为我们的参考基因集本身就是肿瘤的且物种为人,所以我们ssGSEA在这两个参考基因集下只适用于肿瘤且物种为人的情况
下面,我们就继续介绍最后一种免疫浸润分析——MCPcounter
MCPcounter分析
通过前面的学习,我们可以知道,不管是我们要进行CIBERSORT还是ssGSEA,能否适配于我们研究的物种或者是说我们研究的疾病,关键在于参考数据集的获取,所以以后,我们在想要完成我们个性化的免疫浸润分析的时候,可以先问一下自己,有相关的参考数据集吗?这,其实就是免疫浸润方法的本质
那么下面,我们继续介绍最后一个免疫浸润方法——MCPcounter
我们依旧首先来明确我们的输入数据,MCPcounter只需要我们提供一个文件即可完成分析(即表达矩阵)
1.表达矩阵
(观测是基因,变量是样本)
我们通过对方法学文献Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression的阅读,对该免疫浸润的输入数据进行定义
1.数据来源为转录组(芯片、测序均可);
2.MCP-counter is the first validated computational method that enables the robust quantification of the abundance of multiple immune and non-immune stromal populations in the transcriptome of cellularly heterogeneous tissues such as normal or malignant tissues(即MCPcounter可以针对肿瘤/非肿瘤组织,物种为人)
3.MCPcounter估计是“单一样本”的分数,因为它们是在每个样本上独立计算的。这些分数可以用来直接比较队列中不同样本中相应细胞类型的丰度(即该方法可以比较同一个细胞类型在不同sample的丰度,但是不能用于不同细胞类型之间的比较)
4.从方法学中我们可以明确(Already-normalized “level 3” data were downloaded separately for each project)即针对测序数据需要经过标准化的处理(TPM即可),芯片数据我们依旧可以选择RMA方法标准化后的数据
5.那么这里就会存在一个争议,是否需要对数据进行log,别急,我们后面会实际对比一下情况,然后给出建议
掌握了上面的背景知识,我们就可以正式进入代码部分的讲解,其实代码部分就很简单了,因为R包已经把复杂的地方都已经准备妥当,我们实际上只需要提供一个文件即可(另外两个文件在github上下载)
#代码讲解####1.准备工作####rm(list = ls())#install.packages(c('devtools','curl')) library(curl)library(devtools)install_github('ebecht/MCPcounter',ref='master', subdir='Source',force = TRUE)####2.验证安装是否成功####library(MCPcounter)?MCPcounter.estimate#通过帮助文档进行学习####3.免疫浸润分析#####需要连接github在线读取两个文件,所以如果网速卡,可能会影响分析#已经从网站下载好两个文件,直接读取即可#https://github.com/ebecht/MCPcountergenes <- data.table::fread("genes.txt",data.table = F)probesets <- data.table::fread("probesets.txt",data.table = F,header = F)example <- read.table("example.txt", header=T, sep="\t", check.names=F,row.names = 1)##示例数据包含四个样本的一个物种为人的芯片数据,经过RMA方法,没有log#dat <- log(example+1)results<- MCPcounter.estimate(example, featuresType= "HUGO_symbols", probesets=probesets, genes=genes)
可以很清楚的看到,我们只需要提供一个表达矩阵文件就可以完成免疫浸润分析,那么,我们下面就来探讨一下log与否对我们的结果有什么影响?
没有log的结果如下:
经过log的结果如下
首先我们需要明确,CIBERSORT得到的结果为同一个样本不同免疫细胞类型所占的比例;ssGSEA得到的则是同一个样本在不同免疫细胞上的富集分数;那么MCPcounter得到的则是同一个细胞类型在不同样本中的丰度,也就是说MCPcounter关注的点在同一个细胞类型上,不同细胞类型是无法比较的,这就是MCPcounter和CIBERSORT以及ssGSEA的区别
那么我们就知道了,毕竟我们得到的数据大部分是用在可视化上,也就是说只要趋势没有改变,那么log其实就是可有可无的,我们这回来看一下上面两个结果,同一个细胞类型在不同样本上的丰度大小的趋势并没有发生改变,但是如果我们进行热图可视化等方式的时候,因为彼此差距过于悬殊,很有可能会导致可视化的颜色出现偏移,所以这里我们采取折中说,你可以log,然后看可视化,如果颜色配比悬殊,也可以不log,这个并没有影响到结果
好,到这里,我们又产生了一个问题,那如果针对小鼠呢?有没有向前面那样更换参考数据集就可以完成的分析?
其实,因为我们这个分析需要用到两个从github上下载的数据,所以整理数据,莫不如直接去找一个计算细胞丰度的小鼠免疫浸润的方法,很幸运,我们发现了——mMCP-counter
作者亲切的在github上对这个工具进行了概括:小鼠版本的MCPcounter,一种从转录组数据中估计异质组织的免疫和基质组成的工具
而且如果大家阅读了参考文献The murine Microenvironment Cell Population counter method to estimate abundance of tissue-infiltrating immune and stromal cell populations in murine samples using gene expression可以发现,这个工具跟MCPcounter本质上来说其实没有太多的区别,只不过针对的物种不同而已
那么,我们依旧首先对输入数据进行讨论,与MCPcounter类似仅仅只需要一个文件就可以完成分析
1.表达矩阵
(观测是基因,变量是样本)
针对输入数据的形式,在帮助文档中详细阐明,需要对数据进行标准化并且进行对数转换
这里的对数转换可以理解为缩短数值之间彼此的差距而已,我们差异分析的时候取的是log2,这里我们直接取log就可以
至于其它的数据要求参考MCPcounter即可
掌握了上面的背景知识,我们就可以正式进入代码部分的讲解,其实代码部分就很简单了,因为R包已经把复杂的地方都已经准备妥当,我们实际上只需要提供一个文件即可
#代码实战####1.准备工作####install.packages("devtools")library(devtools)devtools::install_github("cit-bioinfo/mMCP-counter")library(mMCPcounter)library(data.table)library(tidyverse)?mMCPcounter::mMCPcounter.estimate####2.免疫浸润####expressionData <- fread("example.txt",data.table = FALSE)rownames(expressionData) <- expressionData[,1]expressionData <- expressionData[,-1]expressionData <- log(expressionData+1)class(expressionData)data <- mMCPcounter.estimate(expressionData, features = c("Gene.Symbol","ENSEMBL.ID","Probes")[1])
然后我们就得到了小鼠的免疫浸润丰度的结果,在同一个细胞类型上查看不同样本的丰度计数即可
好,到这里,我们的推文本来应该就结束了,但是因为过年,所以给还在学习的小伙伴们加加餐,我们再来深入研究一下,我们如果对小鼠的基因进行同源转换到人上后然后进行免疫浸润分析是否可行呢?
先说一下个人观点吧,是可行的,但是如果大家阅读过晨曦写的同源转换的推文大家就应该知道,转换并不是100%的,会有部分基因无法转换,这就会导致这部分基因没有办法参与到免疫浸润分析当中,但是,这个也算是一个思路,感兴趣的小伙伴,可以结合这篇推文以及晨曦写的下面这篇推文一起学习哦~
不同物种也能合并做生信?给你支个妙招,让数据起死回生!
然后如果大家对于免疫浸润的方法还想要了解更多,推荐阅读下面这篇文献
Quantifying tumor-infiltrating immune cells from transcriptomics data
好啦,到这里,本期推文就全部结束了,感谢大家这一年来对挑圈联靠的支持,也希望新的一年可以和大家一起学习更多的知识,也欢迎大家在评论区对推文内容进行讨论和补充,同时也祝大家新的一年、心想事成、身体健康、万事如意、学业事业双丰收~
我是晨曦,我们下期再见
Ps:回复“晨曦免疫浸润”就可以获得本期推文全部代码+示例数据+参考文献
晨曦单细胞笔记系列传送门
晨曦从零开始学画图系列传送门
晨曦单细胞数据库系列传送门

END

撰文丨晨   曦
排版丨四金兄
主编丨小雪球
欢迎大家关注解螺旋生信频道-挑圈联靠公号~
继续阅读
阅读原文