转录组Count格式数据转化为FPKM/TPM格式
很多时候我们得到的转录组格式为Count,例如在TCGA数据库下载的数据,如果我们想使用FPKM格式或者TPM,那么就需要转换,不过TCGA数据库也提供了FPKM的格式,貌似miRNA数据只有Count格式数据,如果想使用FPKM,就需要进行格式转换。
1、计算基因的长度,可以计算基因在染色体的起始和结束之差,也可以计算每个基因的最长转录本或所有外显子之和:
a、计算基因在染色体的起始和结束之差需要用到R包biomaRt,安装方法如下:
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("biomaRt")
b、计算每个基因的最长转录本或所有外显子之和需要用到R包GenomicFeatures,安装方法如下:
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("GenomicFeatures")
2、得到基因长度的文件之后,使用R读取基因长度文件与被转换格式文件,按照下面的公式进行转换
a、Count2FPKM
count2fpkm <- function(countnum, genelen){
total <- sum(countnum)
exp( log(countnum) + log(1e9) - log(genelen) - log(total) )
}
fpkm<- apply(count_data, 2, count2fpkm, genelen = gene_l$gene_l)
fpkm_df<-data.frame(fpkm)
colnames(fpkm_df)<-colnames(count_data)
write.table(fpkm, "fpkm.txt", sep="\t", quote=F, row.names=T)
b、Count2TPM
n <- count_data / gene_l$gene_l
tpm <- t( t(n) / colSums(n) ) * 1e6
tpm<-data.frame(tpm)
write.table(tpm, "tpm.txt", sep="\t", quote=F, row.names=T)
c、FPKM2TPM
FPKM_data<-read.table("fpkm.txt",header = T,row.names = 1)
fpkm2tpm <- function(fpkm)
{
exp(log(fpkm) - log(sum(fpkm)) + log(1e6))
}
tpm2 <- apply(FPKM_data,2,fpkm2tpm)
write.table(tpm2,"tpm2.txt",sep="\t", quote=F, row.names=T)
赞 (0)