該腳本適用于featureCounts數(shù)counts后的cpm和tpm計(jì)算

setwd("~/Desktop")#設(shè)置工作目錄

countdata<-read.table("counts.txt",skip = 1,sep="\t",header = T,row.names = 1)

metadata <- countdata[,1:5]#提取基因信息count數(shù)據(jù)前的幾列

countdata <- countdata[,6:ncol(countdata)]#提取counts數(shù)，counts數(shù)據(jù)主題部分

prefix<-"couts"#設(shè)置輸出文件前綴名

cpm <- t(t(countdata)/colSums(countdata) * 1000000)#參考cpm定義

avg_cpm <- data.frame(avg_cpm=rowMeans(cpm))

#-----TPM Calculation------

kb <- metadata$Length / 1000

rpk <- countdata / kb

tpm <- t(t(rpk)/colSums(rpk) * 1000000)

avg_tpm <- data.frame(avg_tpm=rowMeans(tpm))

write.csv(avg_tpm,paste0(prefix,"_avg_tpm.csv"))

write.csv(avg_cpm,paste0(prefix,"_avg_cpm.csv"))

write.csv(tpm,paste0(prefix,"_tpm.csv"))

write.csv(cpm,paste0(prefix,"_cpm.csv"))