1. RNA-seq数据归一化为什么Raw Counts不够用刚拿到RNA-seq的Raw Counts矩阵时很多新手会直接用它做差异分析结果发现样本间比较总是出现偏差。我刚开始做生物信息分析时也踩过这个坑——直到发现同一个基因在不同样本中的count值差异可能完全来自技术因素。Raw Counts本质上只是比对到每个基因的reads数它隐藏着三个关键干扰项第一是基因长度陷阱。就像用桶接雨水基因越长就像开口越大的桶自然能接到更多reads雨滴。比如一个10kb的基因和1kb的基因即使真实表达量相同前者的count值可能是后者的10倍。去年我分析小鼠转录组时最长的基因Ttn约100kb的count值总是异常高但这并不代表它真的高表达。第二是测序深度差异。假设样本A测了50M reads样本B只测了10M reads那么A中所有基因的count值天然会是B的5倍左右。这就像比较两场降雨量时一个用大桶接1小时另一个用小桶接10分钟直接比较水量毫无意义。第三是片段化偏差。特别是双端测序中文库制备时的随机打断会导致某些基因区域更容易被捕获。有次做植物RNA-seq发现某些基因的5端counts异常高后来发现是GC含量影响了片段化效率。提示覆盖率Coverage ratio和测序深度Coverage depth常被混淆。前者是测到碱基占基因组的比例后者是每个碱基平均被测到的次数。例如10x深度理论上覆盖率约99.9%但实际由于不均匀分布某些区域可能未被覆盖。2. 主流归一化方法原理与实战2.1 RPKM长度与深度双重校正RPKMReads Per Kilobase per Million是最早提出的归一化方法计算公式为RPKM (10^6 * 10^3 * raw count) / (total reads * gene length)其中10^6对应百万级测序深度标准化10^3对应千碱基基因长度标准化。在R中实现RPKM计算calculate_rpkm - function(counts, lengths) { # counts: 基因表达矩阵行基因列样本 # lengths: 基因长度向量 million_reads - colSums(counts) / 1e6 kb_length - lengths / 1e3 rpkm - sweep(counts, 2, million_reads, /) rpkm - sweep(rpkm, 1, kb_length, /) return(rpkm) }但RPKM有个致命缺陷它先除以总reads数再做长度校正。这会导致不同样本的RPKM总和差异巨大使得样本间比较失真。就像先按人数平分披萨再根据食量调整——最后瘦子可能比胖子分得更多。2.2 TPM更合理的校正顺序TPMTranscripts Per Million改进了计算顺序先对每个基因做长度归一化count / gene_length再对样本做深度归一化上一步结果 / sum(normalized_counts) * 10^6Python实现示例import numpy as np def calculate_tpm(counts, lengths): # 长度归一化 rate counts / lengths[:, None] # 深度归一化 tpm 1e6 * rate / np.sum(rate, axis0) return tpm实测对比用同一套人类肝细胞数据RPKM各样本总和差异达3倍而TPM总和均为百万。这正是TPM成为当前金标准的原因——它使得不同样本间具有可比性就像把所有人的收入先按工作时间标准化再比较时薪。2.3 FPKM双端测序的特化版本FPKMFragments Per Kilobase per Million本质上是RPKM的双端测序版本区别在于单端测序FPKMRPKM双端测序成对read计为1个fragment在Linux中用featureCounts统计fragmentsfeatureCounts -p -a annotation.gtf -o counts.txt aligned.bam参数-p表示统计fragments而非reads。注意某些工具如HTSeq默认统计reads会导致双端数据计算不一致。2.4 RPM/CPM小RNA分析的利器RPMReads Per Million/CPMCounts Per Million是最简单的归一化方法RPM (raw count / total reads) * 10^6适用于miRNA等长度差异小的场景。我在肺癌sRNA分析中发现用TPM反而会过度校正短miRNA而RPM能更好保留真实差异。3. 方法对比与选择指南3.1 四种方法数学关系对比方法校正顺序适用场景总和特性RPKM深度→长度早期单端测序各样本不一致TPM长度→深度现代差异分析总和1,000,000FPKM深度→长度双端测序各样本不一致RPM仅深度miRNA/sRNA分析总和1,000,0003.2 选择决策树根据我的项目经验推荐以下选择路径是否为小RNA是 → 选择RPM/CPM否 → 进入下一步测序类型单端 → TPM或RPKM双端 → TPM或FPKM是否需要跨样本比较是 → 优先TPM否 → RPKM/FPKM也可接受特别注意某些工具如DESeq2要求输入Raw Counts而非归一化值。我曾犯过将TPM矩阵输入DESeq2的错误导致差异分析结果完全失真。4. 实战中的常见陷阱4.1 基因长度定义争议不同工具对gene_length的定义不同RSEM使用有效长度effective length考虑片段化偏差HTSeq直接使用外显子长度和StringTie考虑转录本异构体建议保持上下游工具一致。有次我用HTSeq的counts但用RSEM的长度计算TPM导致高表达基因被异常压低。4.2 零长度基因处理某些lncRNA在注释文件中长度为0会导致除以零错误。我的处理方案# 将零长度替换为全基因组中位数 lengths[lengths 0] np.median(lengths[lengths 0])4.3 跨物种比较的坑比较人和小鼠数据时直接使用TPM仍可能存在偏差。因为两物种的平均转录本长度不同小鼠基因通常更短GC含量分布不同 这时需要更复杂的校正方法如TMM或RLE。在最近的项目中我采用edgeR的calcNormFactors先做样本间标准化再用TPM做基因间标准化效果显著优于单一方法。具体代码片段library(edgeR) dge - DGEList(countsraw_counts) dge - calcNormFactors(dge, methodTMM) tpm - calculate_tpm(dge$counts, gene_lengths)
【RNA-seq实战】从Raw Counts到归一化矩阵：RPKM、TPM、FPKM与RPM/CPM的选择与计算全解
发布时间：2026/5/28 16:23:09
1. RNA-seq数据归一化为什么Raw Counts不够用刚拿到RNA-seq的Raw Counts矩阵时很多新手会直接用它做差异分析结果发现样本间比较总是出现偏差。我刚开始做生物信息分析时也踩过这个坑——直到发现同一个基因在不同样本中的count值差异可能完全来自技术因素。Raw Counts本质上只是比对到每个基因的reads数它隐藏着三个关键干扰项第一是基因长度陷阱。就像用桶接雨水基因越长就像开口越大的桶自然能接到更多reads雨滴。比如一个10kb的基因和1kb的基因即使真实表达量相同前者的count值可能是后者的10倍。去年我分析小鼠转录组时最长的基因Ttn约100kb的count值总是异常高但这并不代表它真的高表达。第二是测序深度差异。假设样本A测了50M reads样本B只测了10M reads那么A中所有基因的count值天然会是B的5倍左右。这就像比较两场降雨量时一个用大桶接1小时另一个用小桶接10分钟直接比较水量毫无意义。第三是片段化偏差。特别是双端测序中文库制备时的随机打断会导致某些基因区域更容易被捕获。有次做植物RNA-seq发现某些基因的5端counts异常高后来发现是GC含量影响了片段化效率。提示覆盖率Coverage ratio和测序深度Coverage depth常被混淆。前者是测到碱基占基因组的比例后者是每个碱基平均被测到的次数。例如10x深度理论上覆盖率约99.9%但实际由于不均匀分布某些区域可能未被覆盖。2. 主流归一化方法原理与实战2.1 RPKM长度与深度双重校正RPKMReads Per Kilobase per Million是最早提出的归一化方法计算公式为RPKM (10^6 * 10^3 * raw count) / (total reads * gene length)其中10^6对应百万级测序深度标准化10^3对应千碱基基因长度标准化。在R中实现RPKM计算calculate_rpkm - function(counts, lengths) { # counts: 基因表达矩阵行基因列样本 # lengths: 基因长度向量 million_reads - colSums(counts) / 1e6 kb_length - lengths / 1e3 rpkm - sweep(counts, 2, million_reads, /) rpkm - sweep(rpkm, 1, kb_length, /) return(rpkm) }但RPKM有个致命缺陷它先除以总reads数再做长度校正。这会导致不同样本的RPKM总和差异巨大使得样本间比较失真。就像先按人数平分披萨再根据食量调整——最后瘦子可能比胖子分得更多。2.2 TPM更合理的校正顺序TPMTranscripts Per Million改进了计算顺序先对每个基因做长度归一化count / gene_length再对样本做深度归一化上一步结果 / sum(normalized_counts) * 10^6Python实现示例import numpy as np def calculate_tpm(counts, lengths): # 长度归一化 rate counts / lengths[:, None] # 深度归一化 tpm 1e6 * rate / np.sum(rate, axis0) return tpm实测对比用同一套人类肝细胞数据RPKM各样本总和差异达3倍而TPM总和均为百万。这正是TPM成为当前金标准的原因——它使得不同样本间具有可比性就像把所有人的收入先按工作时间标准化再比较时薪。2.3 FPKM双端测序的特化版本FPKMFragments Per Kilobase per Million本质上是RPKM的双端测序版本区别在于单端测序FPKMRPKM双端测序成对read计为1个fragment在Linux中用featureCounts统计fragmentsfeatureCounts -p -a annotation.gtf -o counts.txt aligned.bam参数-p表示统计fragments而非reads。注意某些工具如HTSeq默认统计reads会导致双端数据计算不一致。2.4 RPM/CPM小RNA分析的利器RPMReads Per Million/CPMCounts Per Million是最简单的归一化方法RPM (raw count / total reads) * 10^6适用于miRNA等长度差异小的场景。我在肺癌sRNA分析中发现用TPM反而会过度校正短miRNA而RPM能更好保留真实差异。3. 方法对比与选择指南3.1 四种方法数学关系对比方法校正顺序适用场景总和特性RPKM深度→长度早期单端测序各样本不一致TPM长度→深度现代差异分析总和1,000,000FPKM深度→长度双端测序各样本不一致RPM仅深度miRNA/sRNA分析总和1,000,0003.2 选择决策树根据我的项目经验推荐以下选择路径是否为小RNA是 → 选择RPM/CPM否 → 进入下一步测序类型单端 → TPM或RPKM双端 → TPM或FPKM是否需要跨样本比较是 → 优先TPM否 → RPKM/FPKM也可接受特别注意某些工具如DESeq2要求输入Raw Counts而非归一化值。我曾犯过将TPM矩阵输入DESeq2的错误导致差异分析结果完全失真。4. 实战中的常见陷阱4.1 基因长度定义争议不同工具对gene_length的定义不同RSEM使用有效长度effective length考虑片段化偏差HTSeq直接使用外显子长度和StringTie考虑转录本异构体建议保持上下游工具一致。有次我用HTSeq的counts但用RSEM的长度计算TPM导致高表达基因被异常压低。4.2 零长度基因处理某些lncRNA在注释文件中长度为0会导致除以零错误。我的处理方案# 将零长度替换为全基因组中位数 lengths[lengths 0] np.median(lengths[lengths 0])4.3 跨物种比较的坑比较人和小鼠数据时直接使用TPM仍可能存在偏差。因为两物种的平均转录本长度不同小鼠基因通常更短GC含量分布不同 这时需要更复杂的校正方法如TMM或RLE。在最近的项目中我采用edgeR的calcNormFactors先做样本间标准化再用TPM做基因间标准化效果显著优于单一方法。具体代码片段library(edgeR) dge - DGEList(countsraw_counts) dge - calcNormFactors(dge, methodTMM) tpm - calculate_tpm(dge$counts, gene_lengths)
