CytoTRACE推测细胞分化状态

python环境

CytoTRACE的iCytoTRACE函数需要调用python去除批次效应，因此需要先设置好python环境

mamba create -n SC &
    &
     mamba activate SC

mamba install -y -c conda-forge python=3.10 notebook ipywidgets pandas numpy seaborn matplotlib ipykernel openpyxl pyarrow scanpy python-igraph leidenalg pytables jaxlib leidenalg
pip install scanoramaCT -i https://pypi.tuna.tsinghua.edu.cn/simple
which python
# /opt/homebrew/Caskroom/mambaforge/base/envs/SC/bin/python

修改bug

找到intervaltree/intervaltree.py这个脚本修改bug

#25行的
import collections
#替换为
import collections.abc as collections
#我的intervaltree/intervaltree.py文件路径在
/opt/homebrew/Caskroom/mambaforge/base/envs/SC/lib/python3.10/site-packages/intervaltree/intervaltree.py

如果你不知道去哪里找intervaltree/intervaltree.py文件可以在terminal中，到入python库，报错的时候有文件路径，如果没报错就不需要修改了。

mamba activate SC
ipython
import scanoramaCT

安装CytoTRACE

using函数是我写在$HOME/.Rprofile中的函数，因此每次打开R就能使用。

using的功能是一次加载多个包，并且使用了suppressPackageStartupMessages函数，因此不会显示加载包过程中的信息。

wget https://cytotrace.stanford.edu/CytoTRACE_0.3.3.tar.gz
using(remotes)
remotes::install_local("CytoTRACE_0.3.3.tar.gz")

从anndata导出数据

adata是注释好细胞类型的数据，CellType是细胞类型，library_id是不同样本编号代表批次效应。

这里使用了Arrow格式作为R和Python的中间数据，可以参考。

adata=sc.read("adata.h5ad").raw.to_adata()
# AnnData object with n_obs × n_vars = 14268 × 32285
#     obs: 'CellType', 'library_id', 'n_genes', 'doublet_score', 'predicted_doublet', 'Cluster'
#     uns: 'Cluster_colors', 'CellType_colors', 'hvg', 'leiden', 'library_id_colors', 'log1p', 'neighbors', 'pca', 'scrublet', 'umap'
#     obsm: 'X_harmony', 'X_pca', 'X_umap'
#     obsp: 'connectivities', 'distances'
adata.to_df().reset_index().to_feather("matrix.arrow",compression='zstd', compression_level=1)
# Phenotype tables, NO headers
adata.obs.to_csv("pd.csv",index=True,header=False)
# batch info
adata.obs.loc[:,['library_id']].to_csv("batch.csv",index=False,header=True)

R中设置python路径

using(reticulate)
reticulate::use_python("/opt/homebrew/Caskroom/mambaforge/base/envs/SC/bin/python")
reticulate::py_exe()

整理input

• 读入

using(data.table, arrow, CytoTRACE, dplyr, tidyr, purrr)
df - arrow::read_ipc_file("matrix.arrow")
batch - data.table::fread("batch.csv")
pd - data.table::fread("pd.csv")

• 整理

datasets是没有名字的列表，里边是表达矩阵Matrix类型，phe是有名字的向量，名字是细胞barcode，值是细胞类型

a - split(df, batch$library_id)
datasets - purrr::map(a, ~ column_to_rownames(.x, "index") %>
    % t())
names(datasets) - NULL
phe - pd$V2
names(phe) - pd$V1
rm(df, batch, pd) # 删除不需要的变量

iCytoTRACE

results - CytoTRACE::iCytoTRACE(datasets, enableFast = TRUE, ncores = 8, subsamplesize = 1000)
CytoTRACE::plotCytoTRACE(results, phenotype = phe, gene = "Top2a", outputDir = './')
CytoTRACE::plotCytoGenes(results, numOfGenes = 10, outputDir = './')

Reference

https://mp.weixin.qq.com/s/S1-ClJEtR0ro0sYnIF6iaw
https://mp.weixin.qq.com/s/Al-FqOLMPBlrrT-JNchVhw
https://mp.weixin.qq.com/s/WiqU3nUFUysXf1M519l_Dg

声明：本文内容由网友自发贡献，本站不承担相应法律责任。对本内容有异议或投诉，请联系2913721942#qq.com核实处理，我们将尽快回复您，谢谢合作！