一千个读者心中有一千个哈姆雷特

分享是一种态度

序言

前面我们分享过这篇文章，见月经周期过程人类子宫内膜的单细胞转录图谱 。但正如一千个读者心中有一千个哈姆雷特，文献阅读也是同样的道理，不同的人关注的角度可能不一样，那么趁热乎让我们换个角度再了解下这篇文章吧！

文章信息

文章题目：Single-cell transcriptomic atlas of the human endometrium during the menstrual cycle
发表时间：2020.09.14
期刊：nature medicine
影响因子：36.13
作者：Stephen R. Quake
单位：斯坦福大学
原文链接：https://www.nature.com/articles/s41591-020-1040-z
DOI ：10.1038/s41591-020-1040-z

摘要

子宫内膜（endometrium）在月经周期（menstrual cycle）中发生重构、脱落和再生，这一系列过程伴随着内膜中不同细胞的基因表达改变。本文首次在单细胞分辨率上描绘了人类子宫内膜的细胞异质性，以及在月经周期中的转录动态。作者依据转录组数据构建了子宫内膜的4个主要转化阶段，并发现了一种未曾报道过的纤毛细胞（ciliated epithelia）。作者还发现，人类胚胎移植的窗口期起始于上皮细胞（epithelia）的突发性转录激活，并伴有基质成纤维细胞（stromal fibroblasts）的广泛蜕膜化。总之，该研究提供了一份高精度的细胞分子图谱，有助于理解人类子宫内膜在月经周期中的稳态维持和转化，以及生殖和再生医学的发展。

测序数据

• 组织：月经周期 Day4 - 27 的人类子宫内膜

• 建库方法：Fluidigm C1 system（Day4 - 27）和 10x Genomics Single cell 3′ reagent kit v.3.1（Day16 - 26）。作者文章的重心主要放在 C1 的数据。

基本数据分析

Fluidigm C1 数据

• 参考基因组：Ensembl human reference genome GRCh38.87 (dna.primary_assembly)

• 上游分析流程：fastqx (v.0.11.7) 质控 + STAR (v.2.5) 比对 + picard (2.9) MarkDuplicates 去除重复 + HTSeq (v.0.7.0) 定量

• 质控：在 R 中利用函数 ecdf() 对 External RNA Controls Consortium（ERCC）的比例（）拟合经验累积分布（empirical cumulative distribution），保留2,148个 ecdf_null () < 0.05 的细胞。

• 归一化：（Seurat 中 NormalizeData() 的 scale.factor 设置成 100000）

• 高可变基因：对给定 ，根据 拟合 Overdispersion，选择离散度 top 1,000 的基因

• 聚类：根据 1 – Pearson’s correlation 确定细胞在高维空间的距离，利用 tSNE 降维后，使用 density-based spatial clustering of applications with noise (DBSCAN) 的方法聚类细胞

• 差异表达基因：Benjamini–Hochberg’s procedure (BH) adjusted P values with Wilcoxon’s rank-sum test （作者直接用 wilcox.test 包完成的，用 Seurat 的 FindMarkers() 应该也能实现）

• 细胞周期分析：

1. 根据 43 个 G1/S genes 和 55 G2/M genes 的平均表达量计算细胞周期评分
2. 对推定的 cycling cells 做差异表达（Wilcoxon’s rank-sum test），鉴定出 cell-cycle-associated transcriptome signatures
3. 利用这些 signatures 重新计算细胞周期评分

10x 数据

• 参考基因组：10x Genomics 提供的 GRCh38-3.0.0 也就是 Ensembl human reference genome GRCh38.93 (dna.primary_assembly)

• 上游分析流程：Cell Ranger (v.3.1.0) with --expect-cells=10000

• 质控：最终得到 71,032 个细胞。质控过程相当复杂：

1. 对 Cell Ranger 输出的全部细胞进行降维、聚类，根据 nUMI、nGene 和线粒体基因百分比，给每个聚类的质量评分（原文未给出计算方式）

2. 对每个聚类做差异表达分析，鉴定 uniquely expressed genes

3. 移除没有鉴定到 uniquely expressed genes，并且质量评分低（未给出阈值）的聚类

4. 移除包含两个不同谱系 uniquely expressed genes 的聚类（doublet clusters）

5. 利用DoubletFinder 去除 Homotypic doublets（由转录组相似的细胞类型构成）

6. 对于包含非纤毛上皮（unciliated epithelia）和基质成纤维细胞（stromal fibroblasts）的聚类，利用 R 包 mixtools 对 nGene 拟合 Gaussian mixture model

7. 对于每种细胞类型，鉴定出 Gaussian distribution with the lowest mean，并对 计算一个阈值 ，去除基因数小于 的细胞

归一化：（Seurat 中 NormalizeData() 的 scale.factor 设置成 100000）

去除批次效应：Seurat v3.1.2 函数 FindIntegrationAnchors() 和 IntegrateData()

去除细胞周期效应：ScaleData()

聚类：FindNeighbors() 和 FindClusters()

可以看到，作者的数据分析夹带了不少私货。对于推断子宫内膜转录动态的部分，我将在主要结果中叙述。

主要结果

解析人类子宫内膜在月经周期中的细胞异质性

作者通过 Fluidigm C1 数据鉴定出 6 种细胞类型：纤毛和非纤毛上皮细胞（unciliated/ciliated epithelium）、基质成纤维细胞（stromal fibroblast）、内皮细胞（endothelium）、巨噬细胞（macrophage）和淋巴细胞（lymphocyte）。在鉴定两群上皮细胞时，作者先是注意到它们均表达上皮细胞的标记基因，进一步的差异表达分析和 Gene Ontology 分析显示，其中一群显著富集了纤毛相关的 cellular component 或 biological process 字段，故命名为纤毛上皮。

作者同时通过补充 10x 数据鉴定出了平滑肌细胞（smooth muscle cell），表达 PDGFRB、MCAM 和 SUSD2，提示其具有间充质干细胞的特征。

根据月经周期中的转录动态划分人类子宫内膜的不同时期

目前通常使用自上一次月经出血开始的天数来区分内膜的不同状态，但由于月经周期长度在同一个体内和不同个体间的差异，以及同一个体内膜本身的细胞异质性，传统方法可能存在很大偏倚。作者分别对非纤毛上皮和基质成纤维细胞利用全转录组做 tSNE 降维，大致将月经周期的转录动态划分为 4 个阶段（Phase 1 - 4）。作者发现，有两名女性的分期和其取样时间点相违背，且同一时间点取样的内膜竟然被划分到不同的阶段，提示了仅利用天数区分内膜状态的局限性。

作者认为子宫内膜周期属于不完全连续的生物学过程，于是进一步使用基于 mutual information (MI) 的方法对上皮细胞和基质细胞构建拟时间轨迹。

1. 计算每个 在 的表达值和时间变量 （月经时间点）之间的 MI：
2. 基于经验累积分布（empirical cumulative distribution）计算 的统计显著性，对非纤毛上皮细胞和基质成纤维细胞分别鉴定出 3,198 和 1,156 个时间相关基因（time-associated genes）
3. 利用时间相关基因重新进行 tSNE 降维，拟合主曲线（principal curve）代表两种细胞在月经周期的拟时间轨迹
4. 以拟时间值作为时间变量 ，重复步骤 1 和 2，对非纤毛上皮细胞和基质成纤维细胞分别鉴定出 1,382 和 527 个基因的转录动态

作者发现，利用时间相关基因进行 tSNE 降维，既明显地区分了 4 个分期的特征，又在分期内或分期间保留了部分连续的特征，这与子宫内膜周期的天然特征更相符。拟时间轨迹和取样时间点之间高度相关，且同一个体的两种细胞亚群的拟时间值也高度相关，提示这种轨迹成功地削弱了个体间和个体内差异造成的偏倚。

作者观察到非纤毛上皮细胞在 Phase 4 和其他 Phase 之间有明显的不连续性，而在基质成纤维细胞中未观察到这一特征。通过分析拟时间轴上的基因表达动态，作者发现以 PAEP、GPX3 和 CXCL14 为代表的基因模块在进入 Phase 4 时出现突然的转录激活，提示 Phase 4 可能代表子宫内膜进入胚胎种植窗口期（the window of implantation, WOI）或分泌中期。在同一时期，基质成纤维细胞除了表达经典的 WOI 相关基因 DKK1 和 CRYAB，还表达蜕膜化相关的转录因子 FOXO1 和标记基因 IL15，并且这些基因的上调早在 Phase 3 就已经开始，到 Phase 4 达到高峰，提示基质成纤维细胞以一种更渐进的方式进入 WOI。

非纤毛上皮细胞中，一部分 WOI 基因在分泌晚期即迅速下调，而另一部分的激活一直持续到新一轮周期的开始（Phase 1），提示 WOI 的关闭是一个连续渐进的过程。同时，基质成纤维细胞的蜕膜化相关基因逐渐下调，结束分泌期的蜕膜化状态，最终造成子宫内膜脱落出血。

在成功构建了非纤毛上皮细胞和基质成纤维细胞的月经周期轨迹之后，作者进一步鉴定了月经周期调控相关的转录因子，并且比较了经典的和新鉴定的标记基因。

进一步区分腺上皮和管腔上皮的转录动态

作者注意到，在非纤毛上皮细胞的 tSNE 图中，细胞在与拟时间轨迹垂直的两个方向上呈现另一种分布。进一步分析发现两个亚群分别表达腺上皮（glandular）和管腔上皮（luminal）的特征基因。作者于是进一步阐述了两个具有一定空间意义的细胞分群在月经周期的转录动态，并通过功能富集分析揭示了管腔上皮在增殖期表现出解剖结构改变、小管生成和 Wnt 通路激活等特征。

蜕膜化过程中淋巴细胞和基质成纤维细胞的相互作用

在蜕膜化（Phase 4）的内膜中，淋巴细胞显著表达 NK 细胞的标记基因，并伴有多种 NK 细胞受体的激活或抑制。作者进一步将淋巴细胞分为 CD3-（NK only） 和CD3+ （NK and T）的亚群。在蜕膜化内膜中，与CD3+亚群相比，CD3- 亚群的细胞毒性相关基因、趋化因子基因上调更加明显，提示其功能激活以及与基质成纤维细胞的相互作用更加活跃（图 a 和图 b）。作者鉴定了蜕膜化内皮中基质成纤维细胞和淋巴细胞各自高表达的受-配体基因，并通过免疫荧光验证了蜕膜化过程中二者的相互作用。

总结

该研究提供了一个相对完整的多样本、多细胞类型和跨数据平台的人类子宫内膜单细胞图谱，报道了新的纤毛上皮细胞亚群；通过推断上皮细胞和基质细胞的转录轨迹，作者建立了一个相对无偏倚区分月经周期的模型，并阐述了胚胎种植窗口期的启动、关闭及其基因表达的时空调控；最后，本文揭示了淋巴细胞和基质成纤维细胞的相互作用对内膜的蜕膜化可能具有重要的意义。本文的数据能为子宫内膜疾病的研究提供正常水平的基准，新的周期计算模型将有助于日后在月经周期某些特定阶段对内膜生理学和生殖医学相关的课题进行集中探索，并结合其他生理指标更精确地研究分子和行为水平的改变对女性月经周期的影响。

在测序方面，本文的大量篇幅都用于阐述 Fluidigm C1 的数据，而商业化平台 10x Genomics 的数据仅在适当的时候用作补充。在数据分析方面，作者没有选择无脑走当下流行的 Seurat 或 Monocle 等流程，而是充分考虑了子宫内膜周期的自然特征，并选用相应的数学模型构建拟时间轨迹，可谓私货满满。我上一次读到类似的文章还是 2018 年发表在 Nature，来自 Seurat 开发团队 Satija lab 的 《Developmental diversification of cortical inhibitory interneurons》，尤其推断细胞轨迹时采用的思路（鉴定时间相关基因 - 降维 - 拟合主曲线）简直如出一辙。本文的关键代码已经上传到了 github（https://github.com/wanxinw/endometrium），原始数据也上传到了 GEO（GSE111976）。值得一提的是，作者上传到 SRA 的 10x 数据和某个学徒作业的文章是类似的，同样需要 fastq-dump 从 .sra 文件提取出 .fastq 才能用于 Cell Ranger 流程。

往期精彩

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

单细胞天地欢迎你