单细胞测序技术作为发表高分文章的利器,吸引了很多科研工作者的关注,但在开展单细胞项目的过程中,科研人员对于测序后的分析挖掘工作或多或少都有些担忧,毕竟单细胞测序数据如此丰富,现有的分析方法又五花八门。当前有哪些好用的高级分析?单细胞测序文献是如何利用这些高级分析的?利用高级分析进行数据挖掘的方向有哪些?

……

工欲善其事,必先利其器,看完这一篇推文,解答以上所有疑惑。

单细胞测序数据分析之1、细胞轨迹分析–揭示细胞发育分化动态轨迹

细胞轨迹分析可以在单细胞分辨率验证已知的细胞分化关系,推断未知的细胞分化路径,挖掘一些稀少的中间状态细胞,解析细胞分化过程中的起调控作用的关键基因,在发育生物学中细胞分化、谱系发育研究方向、肿瘤/疾病微环境中免疫细胞的动态变化研究中均有广泛应用。目前进行细胞轨迹分析的方法和软件非常之多,大致可以概括为两种方法,一种是以monocle[1]软件为代表的拟时序分析(pseudotime analysis),另一种则是以velocyto /scVelo为主的RNA速度分析(RNA velocity)[2]。详细内容可以看这篇推文:单细胞数据分析没有思路?试试细胞轨迹分析~(内附代码)
2021年Nature Communication发表的Single-cell sequencing of immune cells from anticitrullinated peptide antibody positive and negative rheumatoid arthritis[3],研究首先通过比较两种亚型类风湿性关节炎(RA)患者与健康对照外周血免疫细胞图谱差异,关注到了B细胞类型丰度变化,发现ACPA-RA组IGHG4+ plasma B细胞明显升高,HLA-DRB5+ plasma B细胞明显降低(图1,左图);随后对B细胞进行拟时序分析探究B细胞亚型分化关系,发现在健康对照和RA患者之间B细胞亚型存在不同的分化路径:在健康对照中naïve B细胞有两条分化路径,1条是分化为HLA-DRB5+plasma B细胞再分化为plasma B细胞,另1条分化为Mem-unsw B细胞到Mem-sw B细胞;而RA患者中,naïve B细胞直接分化为IGHG4+ plasma B细胞(图1,右图)。
图1 利用Monocle推断类风湿性关节炎PBMCs B细胞的分化轨迹[3]

图1 利用Monocle推断类风湿性关节炎PBMCs B细胞的分化轨迹[3]

2020年Nature Communication发表的Single-cell RNA-seq reveals that glioblastoma recapitulates a normal neurodevelopmental hierarchy[4],利用RNA速度表征了胶质母细胞瘤的分化过程,证实了胶质母细胞瘤层次结构顶点是肿瘤祖细胞,星形细胞、间充质、少突胶质细胞和神经元癌细胞比肿瘤祖细胞的分化程度更高,并且在每个患者中均可见特定细胞分化方向(图2)。

图2 基于RNA velocity分析胶质母细胞瘤细胞分化路径[4]

图2 基于RNA velocity分析胶质母细胞瘤细胞分化路径[4]

单细胞测序数据分析之2、细胞通讯分析–解析细胞间的信号通讯关系

多细胞生物是由很多不同类型细胞组成的开放而复杂体系,配体受体复合物介导的细胞间通讯对协调发育、分化和炎症等多种生物学过程至关重要。细胞通讯分析,又称细胞受体-配体互作分析,是以细胞亚群的基因表达量数据为研究对象,通过获得细胞中配体及受体基因的表达信息,比较细胞类型之间的配体与受体基因表达差异,分析得到细胞亚群间的信号通讯关系,在阐明生物学过程中细胞间通讯的复杂性、多样性和动态性方面有重要意义。

2020年Nature Communications发表的Single-cell RNA sequencing highlights the role of inflammatory cancer-associated fibroblasts in bladder urothelial carcinoma[5],对8例膀胱癌肿瘤样本和3例配对癌旁粘膜样本进行了单细胞转录组测序分析,利用CellphoneDB[6]软件进行细胞通讯分析,分析结果发现,炎性癌症相关成纤维细胞iCAFs与其他细胞类型的相互作用最多,且与内皮细胞ECs的相互作用最强(图3a);iCAF表达更高水平的CXCL12,其受体包括 DPP4、CXCR3、CXCR4 和 ACKR3 (CXCR7),由于CXCR4和CXCR3 在免疫细胞上广泛表达,猜测iCAFs分泌 CXCL12是膀胱癌免疫浸润状态的原因(图3b);iCAFs产生VEGF(包括VEGFA和VEGFB),与ECs的VEGF受体(FLT1、KDR、MET和FLT4)结合,促进血管生成,与肿瘤的恶性进展有关(图3c)。

图3 利用CellPhoneDB分析膀胱癌不同细胞类型间配体-受体相互作用[6]

图3 利用CellPhoneDB分析膀胱癌不同细胞类型间配体-受体相互作用[6]

单细胞测序数据分析之3、GSEA/GSVA分析–基于功能基因集的富集分析策略

GSEA( Gene Set Enrichment Analysis),是2005年由Broad Institute研究开发的一种基于基因集的富集分析方法[7],用来评估一个预先定义的基因集的基因在与表型相关度排序的基因表中的分布趋势,从而判断其对表型的贡献。GSEA是从所有基因的表达丰度出发,分析在不同的通路中的基因的整体表达影响,这也是区别于GO/KEGG富集分析的地方,GSEA不需要设定差异阈值筛选目标基因集,理论上更容易囊括细微但协调性的变化对生物通路的影响。

上述膀胱癌的单细胞文献中[5],利用GSEA对iCAFs和肌癌相关成纤维细胞mCAF两种细胞类型的基因进行了GSEA富集分析,结果显示,iCAFs与胞外基质降解基因集显著相关,暗示iCAFs可能在肿瘤转移中发挥作用,同时iCAFs与细胞因子-细胞因子受体相互作用基因集也显著相关(图4,左图);mCAFs则与肌肉收缩、PGC1A通路显著相关,这一结果与之前的体外研究结果一致(图5,右图)。

图4 利用GSEA分析膀胱癌肿瘤相关成纤维细胞显著富集的通路[5]

图4 利用GSEA分析膀胱癌肿瘤相关成纤维细胞显著富集的通路[5]

GSVA分析(Gene Set Variation Analysis)[8]与GSEA(Gene Set Enrichment Analysis)都是一种基于基因集的富集分析方法,能够有效弥补传统富集分析对微效基因的有效信息挖据不足等问题,更为全面地对某一功能单位的调节作用进行解释,Broad研究所在GSEA发布8年之后,开发了GSVA算法来拓展基因集分析的应用。GSEA分析主要用于两两组间比较的方案设计中,对于分组比较复杂的方案设计则比较适合GSVA分析,GSVA不需要预先进行样本之间的差异分析,依据表达矩阵就可以计算每个样本中特定基因集的变异分数。

2020年Cell Research发表的Single-cell transcriptomics reveals regulators underlying immune cell diversity and immune subtypes associated with prognosis in nasopharyngeal carcinoma[9],利用GSVA分析了鼻咽癌NPC不同细胞亚型间的通路活性差异,结果显示,参与干扰素 (IFN)-α和IFN-γ通路在巨噬细胞中相对上调,而血管生成、NF-κB 介导的 TNF-α 信号传导和缺氧相关信号通路在单核细胞中上调(图5,左图);生发中心(GC)B细胞表现出丙酮酸代谢和MYC途径的激活,而浆细胞中的糖酵解途径被上调,在浆细胞中也检测到巨噬细胞趋化途径的高活性,表明这些途径在促进巨噬细胞募集到 NPC 中的 TME 中的潜在作用(图5,右图)。

图5 利用GSVA分析鼻咽癌不同髓系细胞亚型间的通路活性差异[9]

图5 利用GSVA分析鼻咽癌不同髓系细胞亚型间的通路活性差异[9]

单细胞测序数据分析之4、肿瘤拷贝数变异分析–揭示恶性细胞表型

拷贝数变异(Copy number variation, CNV)是由基因组发生重排而导致的,基因组大片段的拷贝数增加或者减少,基因组结构变异(Structural variation, SV) 的重要组成部分,也是人类疾病的重要致病因素之一。与正常细胞相比,肿瘤基因组部分区域呈现过表达或低表达状态,通过与一组参考的“正常”细胞相比,比较不同样本间或不同细胞类型之间的CNV基因表达差异,探索肿瘤基因组位置上基因的表达强度,最终反映基因大片段区域的CNV事件,鉴定体细胞整个染色体或大片段染色体的扩增或缺失。

2021年Cell death discovery发表的Single-cell RNA transcriptome reveals the intra-tumoral heterogeneity and regulators underlying tumor progression in metastatic pancreatic ductal adenocarcinoma[10],利用inferCNV[11]对胰腺导管腺癌转移灶的5种导管细胞亚型进行inferCNV分析,鉴别恶性细胞,分析结果发现,五种导管细胞亚型都经历了显著的CNV事件,均为转移灶中的恶性细胞(图6E);计算CNV scores并绘制小提琴,发现1型导管细胞的CNV水平显著高于其他类型的导管细胞,提示1型导管细胞恶性程度更高(图6F)。

图6 利用inferCNV揭示胰腺导管腺癌恶性细胞类型[11]

图6 利用inferCNV揭示胰腺导管腺癌恶性细胞类型[11]

单细胞测序数据分析之5、转录因子活性分析–挖掘关键调控转录因子

 单细胞研究通常会涉及到一个核心关键问题:细胞的异质性以及这种异质性是如何发展和维持的。这种细胞异质性很大程度上是由潜在的基因调控网络决定的,特定转录因子(transcription factor,TF)集合的协同表达驱动各自靶标基因的表达,从而建立特定的基因表达谱。因此,单细胞的基因调控网络对于深入挖掘细胞异质性背后的生物学意义是至关重要的。SCENIC Suite[12]是一套基于SCENIC(Single-Cell Regulatory Network Inference and Clustering)来研究和破译基因调控的工具,能从单细胞转录组数据中推断TF、基因调控网络和细胞类型,基本原理是基于共表达和DNA调控保守序列(motif)分析推断基因调控网络,然后在每个细胞中分析网络活性以鉴定细胞状态。

上述鼻咽癌的单细胞文献[9],利用SCENIC分析髓系细胞的转录因子活性,发现BNR1H3 和 TFEC 转录因子活性在巨噬细胞中显著上调(图7e);进一步利用生存分析揭示了NR1H3和TFEC的高表达水平与 NPC 患者预后改善之间的关联(图k, n),暗示NR1H3和TFEC在鼻咽癌发展中的重要调控作用。

图7  基于SCENIC分析鼻咽癌关键调控转录因子[9]

图7  基于SCENIC分析鼻咽癌关键调控转录因子[9]

单细胞测序数据分析之6、加权基因共表达网络分析–筛选表型相关核心调控网络

加权基因共表达网络分析(Weighted Gene Co-expression Network Analysis, WGCNA)[13]是构建基因共表达网络的常用方法,可以探索模块与特定表型或疾病的关联关系,最终达到鉴定基因网络的目的;单细胞测序技术可以揭示特定肿瘤组织中的细胞特异性,对细胞进行分类,并且识别特定的标志物,但其检测的细胞数量和病例来源都是有限的。利用WGCNA分析单细胞转录组测序数据,可以提供一套有别于高変基因、差异分析的方法,不依赖于数据库直接用表达量的相关性值预测调控关系,筛选某些细胞亚群中有关联作用的基因集(称为模块),可以从成千上万的基因中挑选出高度相关的基因的模块,并将模块与表型进行关联,寻找marker gene或治疗靶点。

2021年Nature communications发表的Cellular and molecular landscape of mammalian sinoatrial node revealed by single-cell RNA sequencing[14],利用WGCNA分析成年小鼠窦房结SAN细胞基因表达网络,确定了十个共表达模块(图8a);功能富集分析结果显示有六个模板与离子通道的调节有关,其中红色模块基因主要集中在K+和Ca2+通道上(图8c)。

图8 利用WGCNA分析小鼠窦房结细胞基因共表达调控网络[14]

图8 利用WGCNA分析小鼠窦房结细胞基因共表达调控网络[14]

 看完上述内容,单细胞测序分析思路与方法get起来,属于你的单细胞故事将由你自己讲述~

百迈客拥有专业且强大的单细胞生信分析团队,在项目实战中积累了大量的数据挖掘经验,能提供但不限于以上内容的高级分析,欢迎点击下方按钮联系我们。

图9  百迈客部分分析实例展示

图9  百迈客部分分析实例展示

 

 

参考文献:

[1] http://cole-trapnell-lab.github.io/monocle-release/docs/

[2] La Manno, Gioele et al. “RNA velocity of single cells.” Nature vol. 560,7719 (2018): 494-498. doi:10.1038/s41586-018-0414-6

[3] Wu, Xunyao et al. “Single-cell sequencing of immune cells from anticitrullinated peptide antibody positive and negative rheumatoid arthritis.” Nature communications vol. 12,1 4977. 17 Aug. 2021, doi:10.1038/s41467-021-25246-7

[4] Bergen, Volker et al. “Generalizing RNA velocity to transient cell states through dynamical modeling.” Nature biotechnology vol. 38,12 (2020): 1408-1414. doi:10.1038/s41587-020-0591-3

[5] Chen, Zhaohui et al. “Single-cell RNA sequencing highlights the role of inflammatory cancer-associated fibroblasts in bladder urothelial carcinoma.” Nature communications vol. 11,1 5077. 8 Oct. 2020, doi:10.1038/s41467-020-18916-5

[6] Subramanian, Aravind et al. “Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles.” Proceedings of the National Academy of Sciences of the United States of America vol. 102,43 (2005): 15545-50. doi:10.1073/pnas.0506580102

[7] Chen, Zhaohui et al. “Single-cell RNA sequencing highlights the role of inflammatory cancer-associated fibroblasts in bladder urothelial carcinoma.” Nature communications vol. 11,1 5077. 8 Oct. 2020, doi:10.1038/s41467-020-18916-5

[8] Hänzelmann, Sonja et al. “GSVA: gene set variation analysis for microarray and RNA-seq data.” BMC bioinformatics vol. 14 7. 16 Jan. 2013, doi:10.1186/1471-2105-14-7

[9] Chen, Yu-Pei et al. “Single-cell transcriptomics reveals regulators underlying immune cell diversity and immune subtypes associated with prognosis in nasopharyngeal carcinoma.” Cell research vol. 30,11 (2020): 1024-1042. doi:10.1038/s41422-020-0374-x

[10] Xu, Qianhui et al. “Single-cell RNA transcriptome reveals the intra-tumoral heterogeneity and regulators underlying tumor progression in metastatic pancreatic ductal adenocarcinoma.” Cell death discovery vol. 7,1 331. 3 Nov. 2021, doi:10.1038/s41420-021-00663-1

[11] https://github.com/broadinstitute/inferCNV

[12] https://scenic.aertslab.org/

[13] Langfelder, Peter, and Steve Horvath. “WGCNA: an R package for weighted correlation network analysis.” BMC bioinformatics vol. 9 559. 29 Dec. 2008, doi:10.1186/1471-2105-9-559

[14] Liang, Dandan et al. “Cellular and molecular landscape of mammalian sinoatrial node revealed by single-cell RNA sequencing.” Nature communications vol. 12,1 287. 12 Jan. 2021, doi:10.1038/s41467-020-20448-x

 

最近文章