Single-cell RNA-seq dissects the intratumoral heterogeneity of triple-negative breast cancer based on gene regulatory networks

基于基因调控网络的单细胞RNA-seq剖析三阴性乳腺癌的瘤内异质性

发表期刊：Mol Ther Nucleic Acids

发表日期：2021 Jan 1

影响因子：7.032

DOI:  10.1016/j.omtn.2020.12.018

一、研究背景

乳腺癌总是伴随着高度的腔内异质性，是女性最常诊断的癌肿之一。传统上，乳腺癌患者可根据ER（雌激素受体）、PR（孕激素受体）和HER2状态的免疫组化特征，分为腔内A（LumA）、腔内B（LumB）、人表皮生长因子受体-2+（HER2+）和基底样亚型。 大多数基底样乳腺癌表现为三阴性表型， 三阴性乳腺癌（TNBC）总是比其他亚型具有更高的转移概率和更差的预后。

随着单细胞RNA测序(scRNA-seq)技术的发展，我们可以深入了解单细胞分辨率下的遗传异质性，为癌症研究提供参考。有研究证明了乳腺癌的恶性细胞，来源于上皮细胞。此外，恶性细胞可与免疫细胞和基质细胞相互作用，构成了复杂的肿瘤微环境。应用全局调控模型构建了基于scRNA-seq数据的基因调控网络（GRNs），这可以帮助我们认识疾病中的关键调控因子。

二、材料与方法 

1 数据来源

1）单细胞数据：GSE118390中六名TNBC患者(PT039、PT058、PT081、PT084、PT089和PT126)的scRNA-seq数据，其中包括1,189个高质量的细胞。乳腺癌的恶性细胞来源于上皮细胞，因此，只保留了868个上皮细胞用于识别恶性细胞。从以前的研究中下载了240个正常乳腺上皮细胞的基因表达谱，作为参考。

3）TCGA数据库40中获取105个正常和115个TNBC样本的基因表达谱进行差异化基因表达。

4）METABRIC数据库收集了186例TNBC患者的临床资料。

 2 分析流程

流程图

三、结果展示

01 - TNBC患者的内在分子亚型情况

作者获得了868个上皮细胞，并利用基因表达谱和基因组信息推断其CNVs。利用推断的CNVs将上皮细胞分为3个不同的亚组，集群1和集群3有明显的高变量的CNVs。因此推测这两个集群可能是恶性细胞，其中分别包括134和411个细胞。簇2中的细胞有较低的CNVs变量，这可能是非恶性细胞，由323个细胞组成（图1A和1B）。接下来，采用tSNE进行降维。结果显示，上皮细胞被很好地分割成四个不同的亚群（图1C），这证实了层次聚类的结果。

为了鉴别每个恶性上皮细胞的分子亚型，作者采用PAM50模型将恶性上皮细胞划分为5个分子亚型，分别命名为正常样、基底样、Her2+、LumA和LumB（图1D），发现6名患者的亚型组成呈现出极大的多样性（图1E）。其中PT039和PT081由相对较高比例的基底样亚型细胞组成，通过查阅6例患者的临床资料，发现PT039和PT081的肿瘤体积相对较大（图1E），这表明基底样亚型的高构成可能与乳腺癌患者的进展有关。

图1 TNBC患者的细胞组成

为了研究TNBC肿瘤微环境中的细胞-细胞相互作用网络，采用python CellPhoneDB软件包分析不同细胞类型的配体-受体相互作用。p<0.05的配体-受体对被认为在两种细胞类型中存在显著的相互作用。作者发现，巨噬细胞在细胞-细胞通讯网络中起主导作用（图2A），揭示了巨噬细胞在肿瘤微环境中的重要作用。

然后，对五个分子亚型进行了基因组变异分析（GSVA），以确定细胞在50个标志基因组中的活性。发现五个分子亚型表现出极大的多样性（图2B），这揭示了不同分子亚型的瘤内异质性。例如，基底样亚群在G2/M检查点显著富集（图2C）。

图2 细胞与细胞之间的相互作用和标志基因组活性分析

02 - 分子亚型特异性GRNs

根据基因表达谱和TF结合座信息，构建了正常上皮细胞和5个分子亚型特异性GRNs。首先，计算了显著共表达的TF-基因对。其次，保留显著富集TF结合基团的靶基因进行进一步研究。最后，富集分数（ES）>1被认为是显著的TF-目标对的GRN建设。构建的GRN如图3A-3F所示。

图3 正常上皮细胞和五个分子亚型的GRNs

作者发现，正常上皮细胞和5种分子亚型的GRNs的度数分布大致符合幂律分布（图S1A-S1F），这表明6个网络可能存在高度连接的枢纽节点。

图S1 正常上皮和五种乳腺癌亚型特异性GRNs的程度分布情况

字云图说明了亚型特异性GRNs中TFs的程度（图S2A-S2F）。发现TFs的程度在每个亚型特定的GRN中表现出很大的差异。

图S2 表示正常上皮和5个乳腺癌亚型特异性GRN中TFs程度的字云图

此外，还发现在5个分子亚型特异性GRN中没有保守的调节器（TF-靶对）（图S3A），但在5个GRN中有10个保守调节器（TF）（图S3B）。这些保守的TFs可以在不同的分子亚型中发挥不同的作用，这可能与瘤内异质性有关。

图S3 韦恩图

然后，作者确定了6个GRNs中的关键基因。首先，采用中心性指标来评估基因的重要性（中心性指标包括程度，PageRank，之间，特征值，和封闭性）。然后，采用Q统计量对五个中心性指标进行整合。每个网络中由综合Q统计量排名前1%的基因被认为是关键基因。

接下来，为了评估关键基因是否富集于癌症基因，从癌症基因普查中获得已知的癌症基因。发现每个分子亚型特异性GRN中的关键基因都显著富集于已知的癌症基因。每个亚型中已知癌症基因所涉及的关键基因表现出极大的多样性（图S4）。例如，作者发现致癌基因MYC是基底样和LumB亚型特异性网络中的关键基因，具有较高的节点度。

图S4 每个亚型的癌症基因所涉及的关键基因（紫色格子）

03 - 关键基因的多种功能

为了研究5种分子亚型中关键基因的功能，作者重点研究了这些分子亚型中保守的关键基因。发现ETV6在五个分子亚型特异性GRN中被激活（图3B-3F）。推测ETV6可能在五个分子亚型中发挥重要作用。

然后，重点研究了ETV6的功能多样性。发现ETV6在不同亚型中调控着不同的基因。每个亚型中都有几个亚型特异性基因，这可能与乳腺癌的瘤内异质性有关（图4A）。此外，ETV6在基底样亚型特异性GRN中显示出最高的度（度=411），这可能是基底样亚型恶性特征的关键。

接下来，采用基因本体（GO）注释方案来评估每个亚型中ETV6调控基因的功能。发现ETV6 调控的基因在不同亚型中参与了不同的角色，这表现出高度的瘤内异质性（图 4B）。例如，基底样亚型在血管生成生物过程中被特异性注释，这是由于ETV6对基底样亚型中SGPL1的特异性调控。

然后，作者研究了ETV6在不同亚型中的表达水平。在单细胞分辨率下，发现ETV6在5个分子亚型中的表达明显高于正常上皮组（图4C）。此外，检查了ETV6在The Cancer Genome Atlas（TCGA）基底样患者和正常乳腺样本中的表达。发现与正常对照相比，ETV6在基底样患者中显著上调（图4D）。

此外，有研究表明，强核ETV6染色可能与子宫血清性肿瘤的不良预后有关，怀疑ETV6的激活可能与TNBC患者的生存有关。根据METABRIC数据集，发现ETV6表达较高的患者临床结局相对较差（图4E）。因此，失调的ETV6可能是肿瘤治疗的潜在靶点。

图4 ETV6的功能多样性

四、结论

在本研究中，作者提出了对乳腺癌亚型的GRNs进行综合分析，并对每个分子亚型的关键基因进行了剖析。利用PAM50模型对恶性细胞的亚型进行分类。保留共表达和富含转录结合基团的转录因子（TF）-基因对进行GRN构建。发现ETV6在乳腺癌的五个分子亚型中无处不在地被激活。失调的ETV6可以调控每个亚型中不同的基因，这可能与TNBC的瘤内异质性和进展有关。