在RNA-Seq分析中，为了获取基因表达量差异，于是产生了RPKM、FPKM、TPM定量方法，去除测序深度和基因长度的影响。RPKM用于单端测序；FPKM用于双端测序。而TPM计算方法类似于RPKM，但更具有优势，是目前使用的较多的定量方法。详细的解释及计算公式：
https://www.jianshu.com/p/1940c5954c81
（建议学习一下，至少了解清楚TPM计算过程，否则得到比对结果后自己不知道怎么算TPM）
简单来说，TPM定量分为三步：
1. 基因长度标准化。测序reeads数除以基因长度length，得到RPK
2. 进行百万转换。因为计算的单位为M
3. 测序深度标准化。reads数除以总reads数即总RPK
宏基因组的定量用TPM是个不错的选择。
进行宏基因组的TPM定量，需要文件预测的CDS序列文件和对应的fastq（cleandata）。
软件安装：BWA和Samtools
#傻瓜式安装conda
#samtools
conda install -c bioconda samtools#bwa
conda install -c bioconda bwa
#如果你会自己编译环境可以不用conda安装，反正能装上并且能使用就行以bwa构建CDS序列索引进行比对
#构建索引
bwa index Unigene.fasta#例如我有5个样本A,B,C,D,E
for i in {A,B,C,D,E} \
do bwa mem -k 19 – t 24 Unigene.fasta \
${i}_clean_R1.fastq ${i}_clean_R2.fastq > ./${i}.sam \ #压缩的fastq文件也可以，即后缀为gz
donesamtools提取比对到clean上的reads数
#Samtools转换sam文件为bam文件 1.3版本前 单个样本举例，多样本for循环
samtools view -bS sample.sam > sample.bam
samtools sort sample.bam > sample_sort.bam
samtools index sample_sort.bam#1.3版本后（sort将sam转化为bam与排序同时进行）
samtools sort sample.sam > sample_sort.bam
samtools index sample_sort.bam#获取每个ORF比对的read数
for i in {A,B,C,D,E} \
do samtools idxstats ${i}_sort.bam > ${i}_mapped.txt \
done
#注意，这一步之前需要经过sort和index获得的txt结果有四列，从左至右分别是：基因名、基因长度、mapped_read、unmapped_read
因为输出的文件不含表头，我们手动添加一个
for i in {A,B,C,D,E} \
do \
sed -i "1 i GeneID\t\length\tmapped_read\tunmapped_read " ${i}_mapped.txt \
done
因为TPM计算是不需要unmapped_read的，所以只保留前面三列
for i in {A,B,C,D,E} \
do \
cut -f 1-3 ${i}_read.txt > ${i}_read_cut.txt \
done这时候的文件就是我们需要的了，只需根据公式来计算即可以获得TPM ，这里我贴上用R语言和python计算的代码，比较简单。
R的脚本：
#读取文件
df <- read.delim("BJS_read.txt",header = T,row.names = 1)
#截取需要的部分
df_cut <- df[,1:2]
#计算RPK
df_cut$RPK <- df_cut$mapped_read*1000/df_cut$length
n <- nrow(df_cut)
result <- df_cut[-n,]#文件最后一行存在一个*行，含有NULL，需要去除掉才能计算，否则结果都为NULL
TotalRPK <- sum(result$RPK)
result$TPM <- (result$RPK*10e6)/TotalRPK
write.csv(result,file = "BJS_TPM.csv")
python：
import pandas as pd
count = pd.read_table("BJS_read_test.txt",index_col=0)
count_part = count.iloc[:,[0,1]]
#n = count_part.shape[0]#获取数据框行列数
count_cut=count_part.drop(index="*")
#创建一个新的数据框
new_data = pd.DataFrame(columns=("RPK","TPM"))
df = pd.concat([count_part,new_data],axis=1)
df["RPK"] = df.apply(lambda x: (x["mapped_read"]*1000)/x["length"] ,axis=1)
TotalRPK = sum(df["RPK"])
df["TPM"] = df.apply(lambda x: (x["mapped_read"]*10e6)/TotalRPK,axis=1)
result = df.iloc[:,[0,3]]
result.to_csv("BJStmp.csv",index=True,sep="\t")两个的结果都一样，需要注意的是由sam文件到最后的txt文件的最后一行会含有一个*的行，代表的是未匹配到任何CDS的reads，需要在计算中将其去除。
最终会得到每个样本的TPM文件，可以使用python或者R中的merge函数将所有的样本TPM拼接到一起，得到一个最终的表格
