晨曦碎碎念系列之
一文揭示批次效应并提供解决思路
Hi,大家好,我是晨曦
今天这期推文是晨曦碎碎念的第二期,今天我们聚焦在单细胞文献的检索方面,来探索一下,如何更高效、更全面的检索单细胞相关文献,以及如何检索到更全面的单细胞数据挖掘类文献?
这次我们聚焦的问题在不管是初入生信领域的萌新,还是混迹在生信领域有一段时间的老萌新(Ps:晨曦可能就是老萌新吧QAQ)都不可避免会产生疑惑的问题——批次效应
不管是在分析芯片、bulk-seq还是在分析scRNA-seq的数据时,我们都会遇到这个问题,而且会在脑海中不断涌现各种疑问
1.什么是批次效应?
2.批次效应和生物学差异有什么区别呢?
3.我可以直接去除批次效应吗?
4.有哪些方法可以去除批次效应?
5.芯片、bulk-seq、scRNA-seq去除批次效应的原理和工具都是什么呢?
6.......
基于但不局限在以上问题,这篇推文应一些小伙伴的期望诞生,将为大家剖析一下批次效应的种种细节,希望可以让大家在遇到这类问题的时候可以从容应对~
那么,我们今天的学习就开始啦(文字数量较多,建议边听歌边阅读QAQ)

晨曦碎碎念系列传送门
背景知识
首先,我们需要先理解一下概念即什么是批次效应?
在2010年发表的文献Tackling the widespread and critical impact of batch effects in high-throughput data为我们详细解答了这个问题,部分原文如下:
Batch effects are subgroups of measure-ments that have qualitatively different behaviour across conditions and are unre-lated to the biological or scientific variables in a study. For example, batch effects may occur if a subset of experiments was run on Monday and another set on T uesday, if two technicians were responsible for different subsets of the experiments or if two different lots of reagents, chips or instruments were used. These effects are not exclusive to high-throughput biology and genomics research1, and batch effects also affect low-dimensional molecular measurements, such as northern blots and quantitative PCR. Although batch effects are difficult or impossible to detect in low-dimensional assays, high-throughput technologies provide enough data to detect and even remove them. However, if not properly dealt with, these effects can have a particularly strong and pervasive impact. Specific examples have been documented in published studies2,3 in which the biological variables were extremely correlated with technical variables, which subsequently led to serious concerns about the validity of the biological conclusions4
针对上述内容我们简单概括
批次效应的概念为批次效应是实验结果中出现并导致数据偏离正确性的一种广泛存在的性质
提问:研究变量的选择和批次效应的产生或者消失有关系吗?
回答:没有半毛钱关系~
既然了解的批次效应的概念,那么接下里我们就需要明确批次效应都可能出现在哪些情况?
1.一个实验在不同时间进行,比如,实验的前半部分在周一完成,实验的后半部分在周二完成
2.一个实验由不同的操作者进行,比如,实验的前半部分由A同学完成,实验的后半部分由B同学进行
3.实验的设备不同,举例:测序的平台、试剂、技术、方法等
4.自己的测序数据与公共数据库的数据天生存在着批次效应
针对第四点晨曦给大家解释一下
首先,上传到公共数据库的数据时间上肯定是有一个时间差的,同时很难保证公共数据库的数据和自己的测序数据所采用的外在环境是一致的,所以其天生存在着批次效应
既然我们了解完了批次效应的概念和批次效应所出现的情况,接下来,晨曦想提出一个问题
提问:批次效应的去除究竟是去除了真正的批次效应还是顺手也把部分生物学差异也给去除了?
回答:这个问题的回答其实是比较矛盾的,但是这里晨曦依旧是拿出来文献来给大家进行介绍,我们依旧是选择我们前面拿出来的文献,然后我们往下看,就可以看到下面这句话
这里晨曦简单概括一下,在传统的生物学分析标准流程中,我们习惯用归一化来去除批次效应,所谓的归一化就是调整单个样本测量的全局属性以便更适当地进行比较,但是归一化并不能够完全消除批次效应,因为批次效应是影响特定的基因集,并可能以不同方式影响不同基因,算法显然并不会考虑的如此“个性化”
好啦,到这里,小伙伴们应该有疑问了
提问:曦曦,什么是归一化?标准化?中心化?
回答:要说到这三个概念,晨曦真的是上网查了很多资料,也是踩了很多坑,最后总结出来一个还算可靠的答案,来跟小伙伴们做一个分享
首先,我们说一个让小伙伴们很无语的情况,归一化(Normalization)和标准化(Standardization,这两个概念的英文单词就让晨曦很无语,标准化可以理解,为什么归一化是Normalization啊?
为了解答上面这个问题,下面这段英文就说的很好
This process of converting a raw score into a standard score is called standardizing or normalizing.                                                  
                                                                                                             ---from Wikipedia
好了,从维基百科中我们知道,原来标准化有两个英文名字,一个是StandardizationNormalization
而我们说的归一化却只有一个英文名字,也就是Normalzation
所以从概念上推衍,我们所说的标准化其实是包括归一化的,这就是为什么我们常常说数据标准化,这是一个大概念
然后我们这里探讨一下归一化和标准化的概念,至于为什么不说中心化,是因为这个和前面两个概念很好区分,晨曦放在后面简单带一笔即可
首先既然标准化包含归一化,势必这两个概念有着某些相似的地方即这两个概念其实本质上来说都是对数据的按照某种特定比例进行缩放,但是既然两者又不是等同,势必会有着各自的特点,我们通过概念就可以一窥其中奥妙
归一化(Normalization):将一列数据的值缩放到[0,1]之间,当然这是我们最常见的情况,其实从概念上来讲你缩放到[0,10000]也没人管你,但是一般来说我们都缩放到[0,1]之间,归一化也是不改变数据分布的
Ps:归一化是否改变数据分布在网络上有很多不一样的答案,晨曦搜索到的教程也是各种说法都有,但是从概念上出发归一化和标准化都属于线性变换,既然标准化不改变数据分布,同属于线性变换的归一化也是不改变数据分布的
标准化(Standardization):将一列数据变换为均值为0,标准差为1的状态,且不改变数据的分布,之所以我们标准化后的数据常常符合正态分布,是因为我们的数据在标准化之前往往都是正态分布的形式
而且标准化本身就是一个数据处理的方法,所以不管你的数据是正态还是非正态分布都是可以运用标准化的
了解完两者的概念之后,我们还需要提出两个问题来让我们的理解更进一步
提问归一化和标准化是对整体数据进行缩放还是对数据的某一特征进行缩放?
回答:这个问题看似复杂,其实我们举个例子就可以很好的解决,假设我们有一个数据框,变量为身高、性别、体重、年龄、血压等等个体信息,观测我们则是A同学、B同学、C同学等等,我们归一化或者标准化一定是对数据的某个特征进行缩放,如果我对于整个数据进行缩放,那么毫无疑问试图把身高、性别、体重、年龄等信息融合到一起本身就是个错误,因为它们本身就是不同类型的信息
提问:线性代数中的标准化和我们这里说的标准化是一个意思吗?
回答:其实在线性代数中,将一个向量除以向量的长度也被称为标准化,不过这里的标准化含义为将向量变为长度为1的单位向量,所以叫法是一样的,但是概念确实不同的
概念理解到这里其实就足够了,因为在生物信息中,我们往往只是选择方法和调动R包,细致的算法不是我们学习的重点,但是了解它们却可以让我们充满安全感~
然后我们继续,了解完概念后,我们就需要往下推进,谈一下,两者的联系和差异
联系部分我们在前面都谈过了,本质上来说两者都是对数据进行线性变换(不改变原始数据排列),我们常说的min-max归一化、mean归一化以及Z-score标准化都是隶属于概念下的方法
重点我们来谈一下差异
1. 归一化严格限定了数据范围,但是标准化只是规定了均值为0,标准差为1,可以形象理解为归一化是把数据强行挤在了一块而不考虑数据之间的距离
这里我们可以扩展一个现实情况,通常在机器学习中,标准化是更常用的手段,而相比较来说归一化的使用场景就不是很大,是因为标准化由于只限定了均值和标准差,所以其对样本之间的间距就可以很好的进行划分,当样本存在异常的时候,标准化也不会因为个别样本的异常而导致正常样本出现大范围位置上的变动
但是归一化不一样,因为是把数据挤在了一起,所以当个别样本出现问题的时候就会挤压正常样本,缩短正常样本之间的距离,如果正常样本是一个整体还好,如果正常样本存在亚型那么就可能导致再后续的模型训练中,模型需要更长时间的收敛,因为需要将样本分开,这就是为什么我们经常选择标准化而不是归一化
2.归一化对数据的缩放仅仅取决于极值,而标准化不光取决于极值还取决于每一个数值
这部分在机器学习中可能会用到,与我们常规分析无关,简单理解即可
那么,上面我们就说完了标准化和归一化的概念、联系还有差异,至于其中更细微的算法晨曦就不在这里过多介绍了,因为往往我们只需要选择就可以
既然都已经讲了这么多了,我们不如把这方面给展开进行讲解,接下来我们讲一下标准化、归一化的用途
其实在我们跑标准流程,诸如GEO、TCGA甚至于scRNA-seq的标准流程过程中,各位小伙伴可能会注意到,我们用到上述思想的机会很少,甚至于已经定型的流程压根就不给我们思考的时间,往往一键运行即可
但是,这样做只会把我们推向流水线,会让我们的思想越来越单一,其实标准化和归一化在统计模型、机器学习中运用的还是十分广泛的,我们平时很少注意到,是因为R包往往封装了这类算法,再加上我们运行标码所以也不会注意到了
1.比如在回归模型中,如果自变量x的量纲的不同就会导致回归系数无法直接解读或者报错,所以我们需要把 咨询变量x统一到可比的状态下,这时候我们需要标准化
2.再比如我们常规的主成分分析(PCA)中,第一步就需要把数据整体平移到[0,0]然后去找,因为这样做我 们可以减少一个维度的信息,好达到降维的目的,但是如果数据的量纲不同,则会导致距离计算的时候量纲 大的会影响量纲小的造成数据的偏倚,所以我们需要标准化
3.再比如线性模型的构建如果数据进行过标准化或者归一化,可以加快模型的构建速度,提高分析速率
以上这些都是我们对数据进行处理的好处,既然用途我们也提到了,那么我们该如何选择呢?
这一块其实就涉及到一些线性代数和机器学习方面的知识了,所以这方面晨曦也是看了网络上一些相关的教程最后得到了一个较为准确的答案(欢迎在评论区里讨论)
1.如果对于数据有着严格范围的要求,选择归一化
2.如果不知道选择什么,直接选择标准化
3.如果数据分布不均匀,不要用归一化
4.只有数据彼此量纲不同的时候我们才需要进行标准化,如果量纲相同则不需要,不是我们拿到一个数据就需要去标准化,比如决策树就不同标准化,究其原因我们进行标准化的目的其实在组内就是为了平衡组内的差异,在组间就是平衡组间的差异,仅此而已~
好啦,唠唠叨叨说了这么多,其实能看到这里的晨曦觉得应该不多吧QAQ
我们把标准化还有归一化都给大家介绍的很透彻了,至于一个小尾巴中心化,晨曦就直接说概念了
中心化:也被称为零均值处理,就是将每个原始数据减去这些数据的均值即可
至此,我们背景知识方面就给大家介绍到了这里,相信各位小伙伴再看完这篇推文后,对于自己后续应该选择什么方法应该会有一个大概的方向而不是随机选择
那么下面我们再花点时间讲解一下如何实现吧
毕竟光说不练假把式嘛
实战
我们涉及到批次效应无外乎就是多数据集的合并,下面晨曦将通过代码为大家解析一下我们常规的多数据联合的过程(相关的流程代码依旧来源于最初的哪篇文献)
我们选择消除批次效应的R包为SVA包,函数则是选择ComBat函数
首先我们看一下官方文档对于ComBat函数的示例
# ComBat(dat, batch, mod=NULL, par.prior = TRUE, prior.plots = FALSE)# dat: 基因组测量矩阵(探针维度 X 样本数),探针维度例如marker数、基因数.....,例如表达谱矩阵# batch: 批次协变量,只能传入一个批次协变量!# mod: 这是一个模式矩阵,包含的变量是为了告诉函数我这个变量存在的是生物学差异而不是批次差异,不要动,一般来说就是我们合并前数据的分组信息(Turmor OR Normal)# par.prior: 基于参数/非参数,默认为基于参数然后我们这里拿内置数据集进行一下操作##加载数据BiocManager::install('bladderbatch')library(bladderbatch) #查看内置数据集data()#加载数据集data(bladderdata) #这是一项膀胱癌研究中相关的57个样本的基因表达数据,这些数据已经使用RMA标准化,并且其形式类似于ExpressionSet所以可以通过类似的函数去调用表达矩阵和临床信息(注释信息)pheno <- pData(bladderEset) # 使用pData函数加载注释信息edata <- exprs(bladderEset) # 使用exprs函数加载表达数据
然后我们依次看一下我们最后获得这两个数据是什么样子的
首先是注释信息,如下:
其中的batch就是批次信息
注意:有的数据集是不提供这个信息的,那我们就把数据集看成一个整体,比如说两个数据集合并,那么A数据集的batch就全是1,B数据集的batch就全是2,当然如果各个样本测序有时间,那么我们就需要提取时间并赋值自定义的批次信息,其目的就是把数据划分一下
其次我们看一下表达信息,如下:
然后接下来我们只要几步函数就可以完成批次矫正
#假设这是两个数据集,我们在确保样本行数一致的时候,直接合并即可#exp = cbind(数据集1,数据集2)#去除批次效应pheno$hasCancer <- pheno$cancer == "Cancer"model <- model.matrix(~hasCancer, data=pheno)combat_edata <- ComBat(dat = edata, batch = pheno$batch, mod = model)#评估批量画图dist_mat_combat <- dist(t(combat_edata))clustering_combat <- hclust(dist_mat_combat, method = "complete")plot(clustering, labels = pheno$batch)plot(clustering, labels = pheno$cancer))
至此,批次效应的由来、定义、用途以及如何解决就给各位小伙伴介绍到了这里
本推文算的上是晨曦学习的一个笔记,也欢迎各位小伙伴如果有不同的观点欢迎在评论区进行交流
同时各位小伙伴如果有感兴趣的问题也欢迎在评论区发言哦!
最后,附一下本推文的参考资料:
RNA-seq workflow: gene-level exploratory analysis and differential expression (bioconductor.org)
Tackling the widespread and critical impact of batch effects in high-throughput data | Nature Reviews Genetics
那么这期推文到这里就结束啦,我是晨曦,我们下期再见~
晨曦单细胞笔记系列传送门
晨曦从零开始学画图系列传送门
晨曦单细胞数据库系列传送门

END

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