摘要
按照正常情况，送去测序的样品最好是同一个批次上机测序，避免外部干扰。最近接到一个项目，拿到手的数据就是分了四批。组长提醒我研究一下批量效应的处理方式。因此，这里总结一下批量处理的分析流程。
环境配置
R版本：3.6.1
依赖R包：limma
使用代码：
library(limma) #调用limma包，线性分析主要包
data <- read.table("all_count.txt",header = T, sep = "\t", row.names = 1) #输入定量后的count值
data <- as.matrix(data) # 转化为矩阵
group <- read.table("group.list",header = T, sep = "\t", row.names = 1) #输入样品对应的批次
batch = group$batch #读取批次列
modcombat = model.matrix(~1,data = group) #生成批次设计矩阵，不懂可以百度model.matrix()函数
all_count_cor <- removeBatchEffect(data, batch = batch ,design = modcombat) #消除批量效应
all_count_cor_zero_int <-round(pmax(all_count_cor,0)) #对负值取0（pmax()函数），对小数点四舍五入取整（round()函数）。
head(all_count_cor_zero_int) #检验结果
write.table(all_count_cor_zero_int,"all_count_cor2.txt" , quote = FALSE,sep = "\t") #输出校正后结果
结果展示：
样品名称前缀为生物学重复，中间数字为批次（1-4）。
处理前，通过纵列可以看到，2，3，4批次的样品聚类都是比较接近的，说明可能存在的批量效应
处理后，尽管4批次聚类还是比较近，但其他批次已经更多地依照生物学重复聚类。
总结
批量效应值得警惕，然而，矫枉过正导致差异基因数量地减少也会存在问题。如何平衡是之后仍需研究的地方，比如最开始我使用的sva包的combat函数，同样可以处理批量效应。但是我这边不知道什么问题，使用报错。似乎批量效应的消除主要还是通过线性分析，这方面的基础知识可以再补一补。