文章主要讨论了在处理转录组测序的count矩阵时，手动计算基因表达量的变化倍数与金标准算法（如DESeq2，edgeR，limma-voom）之间的差异。文章强调了金标准算法在处理这种情况时的优势，特别是在处理含有离群点基因的样品时。
文章指出在处理转录组测序的count矩阵时，手动计算基因表达量的变化倍数与金标准算法存在差异，这是正常现象。
文章描述了一种手动计算变化倍数的流程，包括加载count矩阵、数据预处理、差异表达基因的识别等步骤。
文章强调了金标准算法（如DESeq2，edgeR，limma-voom）在处理转录组测序数据时的优势，特别是在处理含有离群点基因的样品时。
文章提到了在某些情况下手动计算和算法计算的变化倍数会出现大面积冲突，主要是由于样品中可能存在离群点基因的影响。
一个小伙伴表示他在处理两分组的转录组测序后的count矩阵的时候，发现手动计算的变化倍数跟金标准算法（DESeq2，edgeR，limma-voom）计算的不一样！
这个超级正常的，我们首先让人工智能大模型解释一下手动计算两分组的count矩阵中基因表达量的变化倍数（fold change）的过程：
变化倍数（fold change）我自己写过一个代码：
load(file = './symbol_matrix.Rdata')
table(group_list)
group_list
symbol_matrix[1:4,1:4] ## 基因名字的样品，矩阵
dat = log2(edgeR::cpm(symbol_matrix)+1)
dat[1:4,1:4]
# 转换为长格式
library(reshape2)
melted_data varnames = c("Gene", "Sample"),
value.name = "ExpressionLevel")
library(ggplot2)
ggplot(melted_data, aes(x = ExpressionLevel, fill = Sample)) +
geom_density(alpha = 0.5) +
labs(title = "Density Plot with Genes of Interest",
x = "Expression Level", y = "Density") +
theme_minimal()
mylogFC=rowMeans(dat[,4:6])-rowMeans(dat[,1:3])
然后我们很容易去拿到金标准算法（DESeq2，edgeR，limma-voom）的差异分析结果进行对比：
load( file = 'deg-raw/DEG_deseq2.Rdata' )
head(DEG_deseq2)
df=data.frame(
deseq2=DEG_deseq2$log2FoldChange,
my=mylogFC[match(rownames(DEG_deseq2),names(mylogFC))]
)
head(df)
plot(df)
boxplot(dat['MKNK1',] ~group_list)
DEG_deseq2['MKNK1',]
df['MKNK1',]
plot(abs(DEG_deseq2$log2FoldChange),log(DEG_deseq2$baseMean+1))
就可以看到，两者其实很容易冲突，因为金标准算法（DESeq2，edgeR，limma-voom）考虑到的事情更多，毕竟是引用近10万的算法啊！
引用近10万的算法如果加上另外的两个包，轻轻松松破十万啊 ：
- 作者：MD Robinson · 2010 · 被引用次数：39135
- 作者：ME Ritchie · 2015 · 被引用次数：31689
如果大家简简单单的随随便便的就能计算变化倍数， 那要这些引用近10万的算法干啥子？
其实大部分基因并不会冲突
人类的gtf文件里面记录的五六万基因里面，只有不到两万是蛋白质编码基因，去除5000低表达量的，剩下的一万多基因其实使用金标准算法（DESeq2，edgeR，limma-voom）拿到的变化倍数跟自己简单单的随随便便的就能计算变化倍数，并不会冲突很多。
但是有一些时候确实是会出现大面积冲突，因为可能会某个组里面的某一些样品有离群点基因。正常情况下每个转录组测序的样品会两千万个reads，平均到两万个基因应该是每个基因1000左右的reads，但是等于核糖体和线粒体它们每个基因可以有几百万个reads，非常的浪费你的测序。这个时候，如果有一个样品的两千万个reads里面的一千万都是核糖体和线粒体，它就影响了你手动计算变化倍数，但是并不会影响金标准算法（DESeq2，edgeR，limma-voom）。因为金标准算法（DESeq2，edgeR，limma-voom）会考虑到这样的离群点基因。
