论文阅读之方法： Single-cell transcriptomics of 20 mouse organs creates a Tabula Muris

The Tabula Muris Consortium., Overall coordination., Logistical coordination. et al. Single-cell transcriptomics of 20 mouse organs creates a Tabula Muris. Nature 562, 367–372 (2018).

论文地址：https://doi.org/10.1038/s41586-018-0590-4

代码地址：http://tabula-muris.ds.czbiohub.org

方法

小鼠和器官采集
四只10至15周龄的雄性和四只处女雌性C57BL/6JN小鼠从查尔斯·里弗的国家衰老研究所种群运输至VA帕洛阿尔托的兽医医学单元（VMU）。小鼠在这两个地点均在12小时光/暗周期下饲养，并提供自由进食和饮水。查尔斯·里弗的饮食为NIH-31，而VA VMU则使用Teklad 2918饮食。未记录或追踪小鼠的同窝关系，小鼠在VA VMU饲养不超过2周，之后进行安乐死。在组织收集前，小鼠于早上8点被放入无菌采集室15分钟，以收集新鲜的粪便颗粒。麻醉使用2.5% v/v Avertin，小鼠体重被称量，毛发被剃除，并通过心脏穿刺采血，之后进行经心脏灌注，使用20ml PBS。然后立即收集肠系膜脂肪组织，避免暴露于肝脏和胰腺灌注液中，这会影响细胞分选。由于无法从同一只小鼠中分离肝脏和胰腺的可活细胞，因此每个器官使用两只雄性和两只雌性小鼠。然后按照以下顺序解剖整个器官：大肠、脾脏、胸腺、气管、舌头、大脑、心脏、肺、肾脏、生殖腺脂肪组织、膀胱、隔膜、肢体肌肉（胫前肌）、皮肤（背部）、皮下脂肪组织（腹股沟脂肪垫）、乳腺（脂肪垫2、3和4）、褐色脂肪组织（肩胛间脂肪垫）、主动脉和骨髓（脊柱和肢体骨）。器官收集于上午10点前完成。经过单细胞解离后，细胞悬液要么用于FACS分选至384孔板，要么用于微流控液滴文库的准备。所有动物护理和程序均按照VA帕洛阿尔托动物研究委员会批准的机构指南进行。

组织解离与样品准备
每种组织的具体协议在论文补充信息中描述。

样品大小、随机化和盲法
研究前未进行样品大小选择。未进行随机化和盲法处理：作者在整个研究过程中均知晓所有数据和与元数据相关的变量。

单细胞方法

裂解板准备
裂解板通过向384孔硬壳PCR板（Bio-Rad HSP3901）中分配0.4μl裂解缓冲液（0.5 U重组RNA酶抑制剂（Takara Bio, 2313B）、0.0625% Triton™ X-100（Sigma, 93443-100ML）、3.125 mM dNTP混合物（Thermo Fisher, R0193）、3.125μM Oligo-dT30VN（Integrated DNA Technologies, 5′AAGCAGTGGTATCAACGCAGAGTACT30VN-3′）和1:600,000 ERCC RNA内标混合物（Thermo Fisher, 4456740））进行制备，使用Tempest液体处理器（Formulatrix）。96孔裂解板也使用4μl裂解缓冲液制备。所有板被用AlumaSeal CS薄膜（Sigma-Aldrich Z722634）封口，并以3220g离心1分钟后快速冷冻在干冰上。板被存储在-80°C直到分选。

FACS分选
单细胞解离后，通过FACS将每个器官和组织中的单细胞分选到384孔或96孔板中。大多数器官使用SH800S（Sony）分选器分选到384孔板中。心脏和肝脏分选到96孔板，心肌细胞手工挑选到96孔板。肢体肌肉和隔膜则使用Aria III（Becton Dickinson）分选器分选到384孔板中。每个384孔板的最后两列故意留空。对于大多数器官，单细胞通过前向散射（forward scatter）选取，死细胞和常见细胞类型通过单一荧光通道排除。使用荧光抗体组合对大多数器官进行稀有细胞群体富集（详见补充信息），但也有一些仅用作可活细胞的染色。根据需要使用颜色补偿。在SH800上，除了最稀有的细胞类型外，其他所有细胞都使用“单细胞”设置，稀有细胞类型使用“超纯”设置。分选器每天使用FACS缓冲液进行校准，每分选8个板后再次校准。对于典型的分选，1-3 ml预染色细胞悬液经过过滤，温和震荡后加载到FACS机器上。首先以低压流动少量细胞进行检查，然后调整压力，暂停流动，解封第一目标板并加载，随后开始分选。如果细胞悬液浓度过高，则使用FACS缓冲液或1X PBS稀释。对于某些细胞类型，如肝细胞，使用96孔板分选，因为无法精确将单个细胞分选到384孔板中。分选后，立即封板并用预标记铝膜封口，离心后迅速冷冻在干冰上。每个384孔板的分选平均需要8分钟。

cDNA合成与文库准备
cDNA合成使用Smart-seq2协议进行。简而言之，含单细胞裂解液的384孔板在冰上解冻，随后进行第一链合成。每孔加入0.6μl反应混合液（16.7 U μl-1 SMARTScribe逆转录酶（Takara Bio, 639538）、1.67 U μl-1重组RNA酶抑制剂（Takara Bio, 2313B）、1.67X第一链缓冲液（Takara Bio, 639538）、1.67μM TSO（Exiqon, 5′-AAGCAGTGGTATCAACGCAGAGTGAATrGrGrG-3′）、8.33 mM二硫苏糖醇（Bioworld, 40420001-1）、1.67 M甜菜碱（Sigma, B0300-5VL）和10 mM MgCl2（Sigma, M1028-10X1ML））使用Tempest液体处理器加入每个孔。逆转录反应在ProFlex 2 × 384热循环仪（Thermo Fisher）中以42°C孵育90分钟，并通过在70°C加热5分钟停止。
随后，使用Mantis液体处理器（Formulatrix）向每孔加入1.5μl PCR混合液（1.67X KAPA HiFi HotStart ReadyMix（Kapa Biosystems, KK2602）、0.17μM IS PCR引物（IDT, 5′-AAGCAGTGGTATCAACGCAGAGT-3′）和0.038 U μl-1 Lambda外切酶（NEB, M0262L）），并在ProFlex 2x384热循环仪上进行第二链合成，使用以下程序：1）37°C孵育30分钟，2）95°C孵育3分钟，3）98°C孵育20秒，67°C孵育15秒，72°C孵育4分钟，共进行23个循环，4）72°C孵育5分钟。
扩增产物用1:10的比例稀释为cDNA与10mM Tris-HCl（Thermo Fisher, 15568025），并通过荧光染料测定法（Quant-iT dsDNA High Sensitivity kit；Thermo Fisher, Q33120）在SpectraMax i3x微孔板读取器（Molecular Devices）上测量浓度。选定具有较高浓度的样本板进行下游处理，如果空白样本（含ERCC的非细胞孔）均值浓度大于0 ng μl-1，并且通过Quant-iT dsDNA标准曲线的线性回归得出R2值大于0.98，则选择这些样本孔。随后将这些孔重新格式化至新的384孔板中，浓度为0.3 ng μl-1，最终体积为0.4μl，使用Echo 550声波液体分配器（Labcyte）进行。

Illumina 测序文库的制备
如之前所述，Illumina 测序文库的制备步骤如下：
简而言之，使用 Nextera XT Library Sample Preparation kit（Illumina，FC-131-1096）对双链 cDNA 进行 tagmentation。每个孔内加入 0.8 μl Nextera tagmentation DNA 缓冲液（Illumina）和 0.4 μl Tn5 酶（Illumina），然后在 55°C 下孵育 10 分钟。通过加入 0.4 μl Neutralize Tagment Buffer（Illumina）并在室温下以 3,220g 离心 5 分钟停止反应。
接着进行索引 PCR 反应，加入 0.4 μl 5 μM i5 索引引物、0.4 μl 5 μM i7 索引引物和 1.2 μl Nextera NPM 混合液（Illumina）。在 ProFlex 2x384 热循环仪上进行 PCR 扩增，使用以下程序：

1. 72°C 保温 3 分钟，
2. 95°C 30 秒，
3. 12 个循环：95°C 10 秒，55°C 30 秒，72°C 1 分钟，
4. 72°C 保温 5 分钟。

文库合并、质量控制与测序
文库制备完成后，使用 Mosquito 液体处理器（TTP Labtech）将每个文库板的孔进行合并。合并后，通过使用 0.7x AMPure 珠（Fisher，A63881）进行两次纯化。文库质量通过 Fragment Analyzer（AATI）进行毛细管电泳评估，并使用 qPCR（Kapa Biosystems，KK4923）在 CFX96 Touch 实时 PCR 检测系统（Biorad）上进行定量。将合并的板池文库标准化至 2 nM，然后从 10 或 20 个板中混合等体积的文库，生成测序样本池。在测序前，加入 0.2% 的 PhiX 控制文库。

384孔和96孔板的测序文库
文库在 NovaSeq 6000 测序系统（Illumina）上进行测序，使用 2 × 100-bp 的双端读取和 2 × 8-bp 或 2 × 12-bp 的索引读取，使用 200 或 300 循环的试剂盒（Illumina，20012861 或 20012860）。

微流控液滴单细胞分析
单细胞通过 GemCode 单细胞仪（10x Genomics）在液滴乳液中捕获，并按照 10x Genomics 协议使用 GemCode 单细胞 3′ Gel Bead 和 Library V2 Kit 构建 scRNA-seq 文库。简而言之，单细胞悬液在倒置显微镜下检查，若样本质量合格，则将样本在含 2% FBS 的 PBS 中稀释至 1000 个细胞每 μl。如果细胞悬液中含有细胞聚集体或碎片，则在 4°C 下以 300g 离心 5 分钟，使用含 2% FBS 的 PBS 进行两次额外的洗涤。使用 Moxi GO II（Orflo Technologies）或血球计数板测量细胞浓度。将细胞加载到每个通道中，目标输出为每个样本 5000 个细胞。所有反应均在 Biorad C1000 Touch 热循环仪上使用 96 深孔反应模块进行。使用 12 个循环进行 cDNA 扩增和样本索引 PCR。扩增的 cDNA 和最终文库通过 Fragment Analyzer 使用高灵敏度 NGS 分析试剂盒（Advanced Analytical）进行评估。10x cDNA 文库的平均片段长度通过 Fragment Analyzer（AATI）和 Kapa Library Quantification 试剂盒（Illumina）通过 qPCR 进行定量。每个文库被稀释至 2 nM，16 个文库的等体积混合，准备每次 NovaSeq 测序运行。使用 100 循环试剂盒进行测序，Read 1 为 26 个碱基，Index 1 为 8 个碱基，Read 2 为 90 个碱基（Illumina 20012862）。在测序前加入 0.2% 到 1% 的 PhiX 控制文库。文库在 NovaSeq 6000 测序系统（Illumina）上进行测序。

数据处理
NovaSeq 测序数据使用 bcl-2fastq 版本 2.19.0.316 进行去重。使用 STAR 版本 2.5.2b，参数 TK，将读取数据对齐到 mm10plus 基因组。使用 HTSEQ 版本 0.6.1p1（默认参数）生成基因计数，除非将‘stranded’设置为‘false’，并将‘mode’设置为‘intersection-nonempty’。
微流控液滴平台的数据使用 CellRanger 版本 2.0.1（来自 10x Genomics）进行去重和对齐，使用默认参数。

聚类
使用 Seurat 包版本 2.2.1 执行标准的过滤、变量基因选择、降维和聚类过程。详细的工作示例，包括每个操作的数学公式，可以在 Organ Annotation Vignette 中找到。根据组织的不同，每个组织的聚类过程（分辨率和主成分数量（PCs））在 GitHub 上的组织特定 Rmd 文件中可以查看。对于每种组织和每种测序方法（FACS 和微流控液滴），执行了以下步骤：

1. 细胞排序：按细胞 ID 进行字典排序，以确保可重复性。
2. 细胞筛选：排除了检测到基因数少于 500 的细胞。（一个基因被认为是检测到的，如果至少有一个读数映射到该基因）。排除了读数少于 50,000（FACS）或 1,000 UMI（微流控液滴）的细胞。
3. 计数归一化：对每个细胞的计数进行对数归一化，使用自然对数计算 1 + 每百万计数（FACS）或 1 + 每万计数（微流控液滴）。
4. 选择变量基因：使用标准化对数离散度（阈值为 0.5）选择变量基因，其中标准化是根据对数均值表达的分箱值分别进行的。
5. 主成分分析（PCA）：将变量基因投影到低维子空间。主成分的数量是根据解释的方差图来选择的。
6. 共享最近邻图：基于低维子空间内由主要主成分构成的欧几里得距离构建共享最近邻图。使用 Louvain 方法的变体进行细胞聚类，该方法在模块化函数中包含分辨率参数。
7. 细胞可视化：使用二维 t 分布随机邻域嵌入（t-SNE）对 PC 投影数据进行可视化。
8. 细胞类型分配：使用已知标记基因的丰度将细胞类型分配给每个聚类。每个组织的标记基因表达图展示在扩展数据中。
9. 细胞类型混合聚类的精细化：当聚类似乎包含多种细胞类型时，通过增加聚类的分辨率参数或对数据进行子集化并重新运行步骤 3 到 7 进行精细化。

对所有 FACS 处理细胞和所有微流控液滴处理细胞进行类似的全局分析，以产生无偏的聚类结果。

异质性评分
设 C 为一个簇，分解为注释的细胞类型 C=T1∪T2∪...∪Tk。对于每一对细胞类型 Ti 和 Tj，我们计算它们成员之间的平均距离：

其中，dij为细胞类型 Ti和 Tj之间的距离。异质性评分 C 为所有至少含有五个细胞的细胞类型 TT 之间距离的最大值。对于 FACS 数据，细胞的向量 x 是上文步骤 5 中的主成分（PC）投影。扩展数据图 9 包含了选定簇的细胞类型距离矩阵 dij的热图，以及包含多个细胞类型的所有簇的异质性评分条形图。

差异表达重叠分析
对于 FACS 和微流控液滴数据，每个器官的差异表达分析通过 Wilcoxon 秩和检验进行，该检验由 Seurat 包中的 FindAllMarkers 函数实现。差异表达分析在细胞本体学组之间进行，结果生成了每个细胞本体学组与同一器官中所有其他本体学组之间差异表达的基因列表（ln(FoldChange)>0.25）。对于 microwell-seq 数据，我们使用了每种细胞类型和每个器官对应的已发布基因列表。然后，我们评估了三种方法之间这些列表的重叠情况。由于命名法不完全相同，分析是在三种方法之间可以匹配的细胞类型之间进行的（见补充表 2）。

批量基因表达谱的相关性
对于 FACS 和微流控液滴共享的 33 个细胞群体，计算每个群体的平均基因表达谱。此类批量基因表达谱的质量取决于检测到的分子总数。FACS 每个细胞检测到更多分子，但细胞数较少；微流控液滴每个细胞检测到的分子较少，但细胞数较多。为了评估在注释的细胞类型上不同方法之间的一致性，使用 Pearson 相关系数对每个共享细胞群体的对数表达谱进行计算。（只有在至少一个数据集中每百万计数大于等于 1 的基因才会被考虑。取对数之前加上伪计数 1 每百万。）

解离评分的计算
对于每个器官，对 140 个与解离相关的基因进行主成分分析（PCA）。第一主成分作为“解离评分”，因为它对应于这些基因的方差。

定义细胞类型富集的转录因子
转录因子被定义为 1,140 个基因，这些基因被注释为 Gene Ontology 术语“DNA 结合转录因子活性”，这些数据从 Mouse Genome Informatics 数据库（http://www.informatics.jax.org/mgihome/GO/project.shtml，2017年11月10日访问）下载。细胞类型被定义为细胞本体学和器官注释的唯一组合（例如，Lung_Endothelial_cell）。所有分析都在完整的数据集上进行，唯一例外的是相关图，这些数据通过随机选择每种细胞类型的 60 个细胞进行子抽样。富集转录因子通过 Seurat FindMarkers 函数进行定义，使用 Wilcoxon 显著性检验，比较目标细胞类型与所有其他细胞类型的组合。筛选标准为 p 值 < 10^-3、平均差异 > 0.2、pct.1 - pct.2 > 0.1（检测差异百分比 > 0.1），以及 pct.1 > 0.3（在 > 30% 的目标细胞中被检测到）。

使用细胞本体学类比较不同方法的细胞类型
我们使用 OntologyX R 包族版本 2.4（包括 ontologyIndex、ontologyPlot 和 ontologySimilarity 库）绘制了代表性的细胞本体学树状图（函数 onto_plot）。为了计算纠缠图（tanglegram），我们使用从所有表达基因构建的树状图作为参考，与通过特定基因本体学细胞功能（如转录因子、细胞表面标记、RNA 剪接因子）生成的树状图进行比较。纠缠得分使用 step2side 方法计算（untangle 函数来自 dendextend R 包）。纠缠是衡量两棵树状图之间对齐程度的指标。纠缠得分的范围是 0（完全对齐）到 1（无对齐）29。

通过随机森林定义转录因子网络
我们使用随机森林（将多个单一决策树结合起来的分类器）来计算每个基因在定义细胞类型中的重要性30。使用 varSelRF R 包版本 0.7-8，该包使用袋外误差作为最小化标准，并通过逐步消除最不重要的变量（即从随机森林分析中返回的变量重要性）来进行变量消除。该算法迭代地拟合随机森林，每次迭代都丢弃变量（基因）中最不重要的部分；选定的基因集是能够产生最小袋外误差的基因集合。这样可以选择出一组非冗余的较小变量集。

报告摘要
关于研究设计的更多信息，请参阅本文附带的 Nature Research 报告摘要。

代码可用性
所有用于分析的代码可以在 GitHub 上找到（https://github.com/czbiohub/tabula-muris）。

数据可用性
所有来自 Tabula Muris 的数据、协议和分析脚本作为公共资源共享（http://tabula-muris.ds.czbiohub.org/）。FACS（https://doi.org/10.6084/m9.figshare.5829687.v7）和微流控液滴（https://doi.org/10.6084/m9.figshare.5968960.v2）所有单细胞的基因计数和元数据，以及所有生成的 R 对象（https://doi.org/10.6084/m9.figshare.5821263.v1），以及 FACS 索引数据（https://doi.org/10.6084/m9.figshare.5975392）可在 Figshare 上访问（https://figshare.com/projects/Tabula_Muris_Transcriptomic_characterization_of_20_organs_and_tissues_from_Mus_musculus_at_single_cell_resolution/27733），原始数据可从基因表达资源库（Gene Expression Omnibus，GSE109774）获得。

这应该是科研工作者的梦中情刊吧。。。学习