一文了解TooManyCells

分享是一种态度

1. TooManyCells 简介

• TooManyCells 是由宾夕法尼亚大学的 Gregory W. Schwartz 等人开发的一种聚类算法，开发的本意是用于 scRNA-seq 数据的分析，当然也可以用于任意的以 observations 为列，以 features 为行的数据

• 相关文献：TooManyCells identifies and visualizes relationships of single-cell clades(https://www.nature.com/articles/s41592-020-0748-5?tdsourcetag=s_pctim_aiomsg) 发表在 Nature Methods

• 据我了解，该算法有 2 种包装形式，一种是 too-many-cells 软件，一种是 TooManyCellsR R 包，由于 R 包报错问题尚未解决，所以本教程只涉及 too-many-cells 软件

• 参考 github 文档

• https://github.com/GregorySchwartz/tooManyCellsR?tdsourcetag=s_pctim_aiomsg

• https://github.com/GregorySchwartz/too-many-cells

2. 软件安装

下面仅介绍一种安装方式，更多方式请见 github 文档

1. 可以选择用 conda 安装 curl

2. 用 curl 安装 stack

## 安装 stack
curl -sSL https://get.haskellstack.org/ | sh #这个 curl 可以用 conda 安装
stack setup

1. 用 stack 安装 too-many-cells

## 安装 too-many-cells
git clone https://github.com/GregorySchwartz/too-many-cells.git
cd too-many-cells
stack install

1. 完成后打开 nohup.out 看到最后提示将 /blabla/.local/bin 路径加入 PATH 变量，因为这是软件的安装路径

3. 文件架构

个人习惯，project 下建三个文件夹：

• scpt，用于存放脚本；

• input，用于存放输入信息；

• out，用于存放输出信息

4. 关于常规表达矩阵转 cellranger 结果

• 一开始由于我的疏忽，以为这个软件只能读取 cellranger 结果中的 3 个文件这样格式的输入文件，但我手头只有普通表达矩阵，在向曾老师求助之后他很快发了推文：表达矩阵逆转为10X的标准输出3个文件

• 从这段非常秀的代码中我学到这么几点：

• R 中的“常规写入”

file="matrix.mtx"
  sink(file)
  cat("%%MatrixMarket matrix coordinate integer general\n")
  cat("%\n")
  sink()

• 操作循环的输出结果（如每次循环返回一个小数据框，最后rbind到一起）

tmp=do.call(rbind,lapply(1:ncol(ct),function(i){
return(data.frame(row=1:nrow(ct),
                  col=i,
                  exp=ct[,i]))
  }))

• write.csv中的append参数，避免覆盖原有内容
  

5.默认参数

• 但是，在我又读了一遍 github 文档之后，发现输入既可以是一个文件夹（里面放 cellranger 的 3 个文件），也可以是一个 csv 格式的普通表达矩阵…于是还是采用后者读取数据

• 除了表达矩阵之外，还需要一个输入文件 labels.csv，大致长下面这个样子：

  如果已知细胞有不同的来源，或者数据分析之后对细胞有注释需求都可以通过这个输入文件实现

• 代码

too-many-cells make-tree \
        --matrix-path ../input/expr_count.csv \
        --labels-file ../input/labels.csv \
        --draw-collection "PieRing" \
        --output ../out \
        > ../out/clusters.csv

• 对于我的需求来说最后的输出只有两个文件有用，一个是 clusters.csv，记录聚类结果，一个是 dendrogram.svg，可视化聚类结果

• 可视化效果：
  
  

6. “修剪”树枝

• 显然，默认参数下的分支太细了，我们可以通过两种方式来调整

• 直接设置 --min-size 参数为一个值，如 100，以规定最小分支细胞数

• 设置 --smart-cutoff 参数为一个值，如 4，以规定最小分支细胞数为 4*median absolute deviation，这个我还不太明白，但是我猜测这是根据每个分支的情况决定分支大小，可能会更合理

• 代码

too-many-cells make-tree \
    --prior ../out \
    --labels-file ../input/labels.csv \
    --smart-cutoff 1 \ #经调试，我的数据最合适的值是1
    --min-size 1 \
    --draw-collection "PieChart" \
    --output ../out_pruned \
    > ../out_pruned/clusters_pruned.csv

• 可视化结果：
  
  

• 只有一种颜色是因为我的 labels 只标了一种，如果 labels 有多种，那么就会呈现这样的效果：
  
  

7. 取细胞子集

• 如果有取出一部分有一定特点的细胞进一步分析，就需要用到 clusters_pruned.csv 了

• 首先查看这个文件的结构

$ head clusters_pruned.csv
cell,cluster,path
AAACGGGAGGTGTTAA.1,9,9/8/7/6/5/4/3/2/1/0
AACACGTTCGGCGGTT.1,9,9/8/7/6/5/4/3/2/1/0
AACCGCGGTATATGAG.1,9,9/8/7/6/5/4/3/2/1/0
ACACCCTTCTGGTTCC.1,9,9/8/7/6/5/4/3/2/1/0
ACCTTTAAGGTGTTAA.1,9,9/8/7/6/5/4/3/2/1/0
ACGAGGACACGTTGGC.1,9,9/8/7/6/5/4/3/2/1/0
AGGGAGTCAGGCTCAC.1,9,9/8/7/6/5/4/3/2/1/0
AGGGATGAGCGATAGC.1,9,9/8/7/6/5/4/3/2/1/0
AGTGGGAAGATGTAAC.1,9,9/8/7/6/5/4/3/2/1/0

• 看来是通过数字的形式记录 cluster 信息的，并且记录了从小到大的每一个分支，但是这个数字和图怎么对应呢？可以把数字标出来

too-many-cells make-tree \
    --prior ../out \
    --labels-file ../input/labels.csv \
    --smart-cutoff 1 \
    --min-size 1 \
    --draw-collection "PieChart" \
    --draw-node-number \ #只需多加这个参数
    --output ../out_pruned \
    > ../out_pruned/clusters_pruned.csv

• 结果：
  
  

• 那么接下来就可以取任意分支的所有细胞的 barcodes 并由此去除该细胞子集的表达信息了

8. 其他可视化选项

• 由于我没有这些需求，所以这里仅仅是搬运 github 文档

8.1. 分支末端画成饼图

• 代码

too-many-cells make-tree \
    --prior out \
    --labels-file labels.csv \
    --smart-cutoff 4 \
    --min-size 1 \
    --draw-collection "PieChart" \
    --output out_pruned \
    > clusters_pruned.csv

• 可视化结果：
  
  

8.2. 调整树枝宽度

• 代码

too-many-cells make-tree \
    --prior out \
    --labels-file labels.csv \
    --smart-cutoff 4 \
    --min-size 1 \
    --draw-collection "PieChart" \
    --draw-max-node-size 40 \
    --output out_pruned \
    > clusters_pruned.csv

• 可视化结果：
  
  

8.3. 不按分支大小缩放

• 代码

too-many-cells make-tree \
    --prior out \
    --labels-file labels.csv \
    --smart-cutoff 4 \
    --min-size 1 \
    --draw-collection "PieChart" \
    --draw-max-node-size 40 \
    --draw-no-scale-nodes \
    --output out_pruned \
    > clusters_pruned.csv

• 可视化结果：
  
  

9. 体现特定基因表达量

• 代码

too-many-cells make-tree \
    --prior ../out \
    --matrix-path ../input/expr_count.csv \
    --labels-file ../input/labels.csv \
    --smart-cutoff 1 \
    --min-size 1 \
    --feature-column 2 \
    --draw-leaf "DrawItem (DrawThresholdContinuous [(\"gene1\", 0), (\"gene2\", 0)])" \
    --draw-colors "[\"#e41a1c\", \"#377eb8\", \"#4daf4a\", \"#eaeaea\"]" \
    --draw-scale-saturation 10 \ #如果不加这个参数，很可能表达量普遍较低以至于整张图没有颜色，至于这个值多少比较合适我还没有试过
    --output ../out_gene_expression \
    > ../out_gene_expression/clusters_pruned.csv

• 可视化结果：
  
  

10. 任意分支之间的差异分析

• 比如我想比较上图中两根蓝色（也就是 gene1 高表达，gene2 低表达）的分支，一看他们的 cluster 号码分别是 110 和 148（当然可以直接空缺一个数字，表示比较一个分支和其他所有细胞）

• 代码

too-many-cells differential \
    --matrix-path ../input/expr_count.csv \
    -n "([110], [148])" \
    +RTS -N24
    > ../out/differential.csv

11. diversity

• 比较两个细胞群的多样性（需要先跑过 make-tree 得到结果）

• 代码

too-many-cells diversity\
    --priors ../out1 \
    --priors ../out2 \
    -o ../out_diversity_stats

12. 伪时序分析

• 代码

too-many-cells paths\
    --prior ../out \
    --labels-file ../input/labels.csv \
    --bandwidth 3 \
    -o ../out_paths

• 结果：
  
  

• 个人认为还是用 slingshot 吧…

往期回顾

---

---

---

如果你对单细胞转录组研究感兴趣，但又不知道如何入门，也许你可以关注一下下面的课程

看完记得顺手点个“在看”哦！

生物 | 单细胞 | 转录组丨资料每天都精彩

长按扫码可关注