head(mycounts)
rownames(mycounts)<-mycounts[,1]
mycounts<-mycounts[,-1]
head(mycounts)
kb <- mycounts$Length / 1000
kb
countdata <- mycounts[,1:9]
rpk <- countdata / kb
rpk
tpm <- t(t(rpk)/colSums(rpk) * 1000000)
head(tpm)
write.table(tpm,file="2020武汉加油_tpm.xls",sep="\t",quote=F)
fpkm <- t(t(rpk)/colSums(countdata) * 10^6) （之前这里写成了10^9，多谢@不爱说话的生物狗 提醒，现在已经修改）
head(fpkm)
write.table(fpkm,file="2020武汉加油_fpkm.xls",sep="\t",quote=F)
fpkm_to_tpm = t(t(fpkm)/colSums(fpkm))*10^6
head(fpkm_to_tpm)
当然，已知所有基因的FPKM情况下，可以通过上述公式直接在excel里计算相应基因的TPM值。
以PBMC数据集进行探索
Seurat通过CreateSeuratObject函数创建对象后，将我们导入的UMI count原始稀疏矩阵储存在pbmc@assays[["RNA"]]@counts，此外Seurat自动计算每个细胞总的UMI count，即每一列数字之和，储存在[email protected][["nCount_RNA"]]；计算每个细胞总的基因数，每一列非0的行数，储存在[email protected][["nFeature_RNA"]]
我们也可以自己计算验证：
[email protected][["nCount_RNA"]]结果为向量，pbmc[["nCount_RNA"]]结果为矩阵
手动计算前六个细胞的基因Feature_RNA数
与Seurat自动计算的结果作对比
与Seurat自动计算的结果作对比
手动计算第一列（细胞）标准化后的值
Seurat自动计算
标准化类似于RNA-seq的FPKM，计算公式为：In( 1 + ( UMIA ÷ UMITotal ) × 10000 )，R语言中log表示自然对数In（以e为底的对数），此处的 UMITotal即为每个细胞对应的nCount_RNA。
pbmc@assays[["RNA"]]@data在未标准化之前，储存的是原始稀疏矩阵，标准化后，储存的为标准化的矩阵。
目前还不知道手动计算公式。
归一化后的矩阵储存在pbmc@assays[["RNA"]]@scale.data
