CytoTRACE2：单细胞转录组细胞分化潜能推断-拟时起点参考（R语言版及Python版）

###加载数据及安装包library(Seurat)DimPlot(sce1, label = T)sce_sub <- sce1[,sce1$cluster %in% c("YSMP","GMP","Myeloblast","Monocyte")]
devtools::install_github("digitalcytometry/cytotrace2", subdir = "cytotrace2_r")library(CytoTRACE2)

#data running-主要函数cytotrace2cytotrace2_sce <- cytotrace2(sce_sub, #seurat对象                             is_seurat = TRUE,                              slot_type = "counts", #counts和data都可以                             species = 'human')#物种要选择，默认是小鼠

class(cytotrace2_sce)# [1] "Seurat"# attr(,"package")# [1] "SeuratObject"


# cytotrace2_res <- cytotrace2(sce_sub@assays$RNA$data, #seurat对象#                              species = 'human')#物种要选择，默认是小鼠# # class(cytotrace2_res)# [1] "data.frame"

annotation <- data.frame(phenotype = sce_sub@meta.data$cluster) %>%   set_rownames(., colnames(sce_sub))
# plotting-一次性生成多个图，然后储存在一个list，用$查看即可plots <- plotData(cytotrace2_result = cytotrace2_sce,                   annotation = annotation,                   is_seurat = TRUE)

#如果这些图您需要放在文章中，需要修饰也是可以的#因为是基于ggplot的作图，所以修饰就很简单了#比如我们修饰一下主题library(ggplot2)for(i in 1:(length(plots)-1)) {

  plots[[i]] <- plots[[i]]+theme_bw()

}


#可以一个个查看图并保存# #p1# plots$CytoTRACE2_UMAP# #p2# plots$CytoTRACE2_Potency_UMAP# #p3# plots$CytoTRACE2_Relative_UMAP# #p4# plots$Phenotype_UMAP# #p5# plots$CytoTRACE2_Boxplot_byPheno
#我们这里为了方便展示，组合展示library(cowplot)plot_grid(plots[[1]],plots[[3]],plots[[4]],          plots[[5]],ncol=2)#ncol=4表示图片排为几列

cd data_analysis/cytotrace2_py/git clone https://github.com/digitalcytometry/cytotrace2cd cytotrace2/cytotrace2_pythonconda env create -f environment_py.ymlconda activate cytotrace2-pypip install .

gene_exp <- as.matrix(GetAssayData(sce_sub, layer = "counts"))write.table(gene_exp, file = "gene_exp.txt", sep = '\t',quote=F)
cell_anno <- data.frame(cellid = rownames(sce_sub@meta.data),                        celltype = sce_sub@meta.data$cluster)write.table(cell_anno, file = "cell_anno.txt", sep = '\t',quote=F, row.names = F)

getwd()setwd("/home/tq_ziv/data_analysis/cytotrace2_py/")

# sce_sub <- sce1[,sce1$cluster %in% c("YSMP","GMP","Myeloblast","Monocyte")]# save(sce_sub, file = "sce_sub.RData")
library(sceasy)library(reticulate)use_condaenv('sceasy')loompy <- reticulate::import('loompy')sceasy::convertFormat(sce_sub, from="seurat", to="anndata", outFile='sce_sub.h5ad')

import scanpy as scadata=sc.read_h5ad("./sce_sub.h5ad")

expression_matrix = pd.DataFrame(adata.to_df().T) #需要转置一下expression_matrix.head()expression_matrix.to_csv('expression_matrix.txt',sep="\t")

cell_annotations  = pd.DataFrame(data=adata.obs["cluster"])cell_annotationscell_annotations.to_csv('cell_annotations.txt',sep="\t")

#直接终端运行cytotrace2 --input-path gene_exp.txt --annotation-path cell_anno.txt --species human

#python中运行from cytotrace2_py.cytotrace2_py import *exp_path = "./expression_matrix.txt"annotation_path = "./cell_annotations.txt"species = "human"
results =  cytotrace2(exp_path,                     annotation_path=annotation_path,                     species=species)