CN113838528A - 基于单细胞免疫组库数据的单细胞水平耦合可视化方法 - Google Patents

Abstract

本发明公开了一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法。方法包括：筛选收集已有的单细胞免疫组库原始数据集，获得多个单细胞免疫组库原始数据集，对各个单细胞免疫组库原始数据集进行特征提取后，获得各个单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息，再搭建数据库和终端浏览器界面，实现单细胞免疫组库原始数据集的所有信息在终端浏览器界面中进行可交互的可视化。本发明收集大量的单细胞免疫组库原始数据集，允许短时间且低计算资源耗费的情况下，允许以无技术偏差的方式比较不同样本、不同疾病模型、不同文献来源的单细胞免疫组库样本数据；同时实现了对交互式可视化的TCR/BCR‑转录组耦合分析。

Description

基于单细胞免疫组库数据的单细胞水平耦合可视化方法

技术领域

本发明涉及了单细胞免疫组库数据的分析方法，尤其是涉及了一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法。

背景技术

T细胞和B细胞是体内负责免疫的两类细胞。T细胞受体(TCR)和B细胞受体(BCR)分别是表达在T细胞膜表面和B细胞膜表面的一种受体(蛋白质)。TCR/BCR的物理结构决定了他们所在细胞所能识别的抗原或病原体。理论上，每个T细胞和B细胞所表达的TCR/BCR都是独一无二的(一种TCR/BCR只在一个细胞上表达)。因此，每个人体体内总共10⁶～10⁸种T细胞和B细胞能够应对大量的外来病原体。在大量的TCR/BCR种类中，如果能够找到专门消灭某种病原体的TCR/BCR，将可以促进新疾病治疗方案的研究。

单细胞免疫组库技术能够将疾病条件下的TCR/BCR库复现。已有方法可以解析单细胞免疫组库技术所得到的样本TCR/BCR库中克隆扩增的TCR/BCR和clonotype，提供样本clonotype的部分特征。

但是，TCR/BCR所在的T或B细胞各基因的表达情况尚不能被分析，TCR/BCR在样本条件下的作用目前只能通过其富集程度推测，而TCR/BCR在样本条件下所发挥的作用极大程度上需要依赖TCR/BCR所在的T或B细胞的转录组特征(例如所在细胞的类型，特定基因的表达等)来进行推断。

同时，已有的方法只能提取并展示样本TCR/BCR库或样本富集的clonotype内TCR/BCR共有特征。需要解释的是，clonotype是多个有共同特征的TCR/BCR形成的TCR/BCR群，目前的方法是提取和分析样本每个clonotype所代表的TCR/BCR群的共同特征，如共同CDR3特征、共同V(D)JC基因使用特征等。目前已有方法的精度尚不能达到对单个细胞、单个TCR/BCR的特征提取和分析，如单个TCR/BCR的CDR3编码序列、V(D)JC基因使用信息，单个TCR/BCR全长编码序列。精度不够高可能导致错失有生物学功能的TCR/BCR信息。同时，由于不提供单个TCR/BCR V(D)J全长序列的信息，无法通过湿实验的手段重构TCR/BCR进行功能研究。

此外，目前的方法只能实现对单个样本的单细胞TCR/BCR测序数据的解析，在实现跨样本和跨文献的数据的比较时，需要对每个样本从头处理，处理结束后才能进行比较分析。这种方法需要耗费大量时间和计算资源。

发明内容

为了解决背景技术中存在的问题和需求，本发明提供了一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，是一种针对单细胞免疫组库数据的，大规模、高精度、耦合细胞转录组特征的分析方法。本发明将解决以下技术问题：

1.耦合分析单细胞TCR/BCR测序数据和单细胞转录组测序数据，全方位提取、分析、展示样本内的TCR/BCR的特征；

2.提取、分析、展示样本内全部TCR/BCR的全长编码序列，为TCR/BCR的功能分析提供条件；

3.大规模整合不同样本、不同文献来源、不同疾病模型的单细胞免疫组库数据，实现短时间及低资源耗费的免疫组库数据的跨样本、跨文献比较。

本发明采用的技术方案是：

本发明包括以下步骤：

步骤S1：筛选收集已有的单细胞免疫组库原始数据集，获得多个单细胞免疫组库原始数据集，其中每个单细胞免疫组库原始数据集均由单细胞TCR/BCR测序数据与其相关联的单细胞转录组测序数据组成；

步骤S2：选择一个单细胞免疫组库原始数据集，对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据进行标签以及特征提取，获得各个TCR/BCR的细胞标签和TCR/BCR特征，TCR/BCR特征包括各个TCR/BCR的全长编码序列、TCR/BCR所使用的V(D)JC基因、抗原互补决定区3的核苷酸序列和氨基酸序列以及TCR/BCR双链信息；

步骤S3：对各个TCR/BCR的抗原互补决定区3的核苷酸序列进行克隆型定义，获得各个TCR/BCR的克隆型；

步骤S4：对当前单细胞免疫组库原始数据集中的单细胞转录组测序数据进行标签提取，获得各个细胞的转录组信息的细胞标签，计算各个细胞的细胞转录组特征，获得单细胞转录组测序数据的基因表达计数矩阵，以及预测单细胞转录组测序数据中各个细胞的类型；

步骤S5：根据当前单细胞免疫组库原始数据集的各个TCR/BCR的细胞标签和各个细胞的转录组信息的细胞标签，实现当前单细胞免疫组库原始数据集的TCR/BCR—转录组耦合，获得当前单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息；

步骤S6：对剩余的单细胞免疫组库原始数据集同样重复步骤S2-S5，获得各个单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息；

步骤S7：数据库搭建；对各个单细胞免疫组库原始数据集进行标注后每个单细胞免疫组库原始数据集对应有唯一的识别码，将各个单细胞免疫组库原始数据集的基本信息、各个细胞的TCR/BCR—转录组信息以及各个单细胞免疫组库原始数据集的来源论文均存储到数据库中，通过识别码在数据库中进行数据集的索引；

步骤S8：数据库的可视化；搭建终端浏览器界面，终端浏览器界面与数据库相结合，在终端浏览器界面中选择单细胞免疫组库原始数据集，根据所选择的单细胞免疫组库原始数据集的识别码调用数据库中所选择的单细胞免疫组库原始数据集的所有信息，所选择的单细胞免疫组库原始数据集的所有信息在终端浏览器界面中进行可交互的可视化。

所述步骤S2具体为：

S2.1：选择一个单细胞免疫组库原始数据集，其中单细胞TCR/BCR测序数据主要由多个TCR/BCR的序列信息组成；根据TCR/BCR参考基因组，对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据中各个TCR/BCR的序列信息进行比对，获得各个TCR/BCR的单链所使用的V(D)JC基因以及各个TCR/BCR的单链的完整性；

S2.2：从各个TCR/BCR的序列信息中提取各个TCR/BCR的单链的抗原互补决定区3的核苷酸序列和氨基酸序列；

S2.3：从各个TCR/BCR的序列信息中提取各个TCR/BCR的细胞标签；对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据进行筛选过滤，获得筛选后的多个TCR/BCR的序列信息；

S2.4：根据各个TCR/BCR的单链所使用的V(D)JC基因和各个TCR/BCR的细胞标签，从筛选后的各个TCR/BCR的序列信息中提取各个TCR/BCR的单链的从5’非编码区到C基因区的全长编码序列；

S2.5：根据各个TCR/BCR的细胞标签，将TCR的α和β链进行匹配以及BCR的轻链和重链进行匹配，获得TCR/BCR双链信息。

所述步骤S4具体为：

S4.1：当前单细胞免疫组库原始数据集中的单细胞转录组测序数据主要由多个细胞的转录组信息组成，提取单细胞转录组测序数据中各个细胞的转录组信息的细胞标签；

S4.2：根据人类蛋白表达参考基因组，对每个细胞的转录组信息进行比对，计算每个细胞中各基因的表达量并作为当前细胞的细胞转录组特征，由各个细胞的细胞转录组特征组成基因表达计数矩阵；

S4.3：根据各个细胞的细胞转录组特征，利用有监督法或者无监督法对各个细胞的类型进行预测，获得各个细胞的类型。

所述S4.3具体为：

A利用有监督法预测细胞类型

根据各个细胞的细胞转录组特征，计算各个细胞的细胞转录组特征和已知各类细胞的参考转录组特征的相似度，根据相似度高低判断各个细胞的类型；

B利用无监督法预测细胞类型

根据各个细胞的细胞转录组特征，计算各个细胞的细胞转录组特征之间的相似度，根据细胞与细胞的细胞转录组特征的相似度划分细胞群，获得当前单细胞转录组测序数据的各个细胞群，将各个细胞群和当前单细胞转录组测序数据的全部细胞的细胞转录组特征做对比，获得各个细胞群的特异性表达基因。

所述步骤S5具体为：

根据当前单细胞免疫组库原始数据集的各个TCR/BCR的细胞标签和各个细胞的转录组信息的细胞标签，将各个细胞的TCR/BCR特征、克隆型、细胞转录组特征和类型进行耦合，过滤仅有TCR/BCR信息或仅有转录组信息的细胞，过滤双链无法匹配的TCR/BCR，获得当前单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息。

所述步骤S8中，使用t-SNE方法对所选择的单细胞免疫组库原始数据集的基因表达计数矩阵进行降维处理后，获得所选择的单细胞免疫组库原始数据集的全部细胞的降维后的细胞转录组特征，以二维散点图的方式展示全部细胞的降维后的细胞转录组特征；使用小提琴图显示各细胞类型的各个基因的表达量；使用直方图显示各个细胞类型的数量；使用网络图显示各个TCR/BCR的克隆型；使用直方图显示各个克隆型的细胞数量；使用复合直方图显示各个克隆型内部的各细胞类型的细胞数量；获得多种可视化图；通过细胞标签对多种可视化图进行可视化交互，实现了有相同克隆型的细胞的细胞转录组特征的提取和比较以及单细胞的TCR/BCR特征-细胞转录组特征的耦合。

所述S2.4具体为：

对单细胞TCR/BCR测序数据进行处理后，获得包含所有TCR/BCR全长序列文件，将各个TCR/BCR的细胞标签以及各个TCR/BCR的单链的完整性作为索引依据，在包含所有TCR/BCR全长序列文件中找到对应TCR/BCR的全长编码序列的位置并提取，获得TCR/BCR双链未匹配的全长编码序列；

再根据各个TCR/BCR的单链所使用的V(D)JC基因中的C基因信息判断TCR/BCR双链未匹配的全长编码序列所来自的具体肽链，最后根据各个TCR/BCR的细胞标签实现TCR/BCR双链未匹配的全长编码序列的双链匹配，获得各个TCR/BCR的单链的从5’非编码区到C基因区的全长编码序列。

本发明的有益效果是：

1.耦合分析。本发明对单细胞免疫组库数据进行交互式可视化的TCR/BCR-转录组耦合分析，提取、分析了前10大克隆型的转录组特征(特定基因表达、细胞类型组成)，帮助定位疾病样本中每个克隆型乃至每个TCR/BCR在样本条件下的生物学功能，发现疾病特异性、样本特异性的TCR/BCR；

2.单细胞精度。本发明提供每个样本中全部TCR/BCR的特征情况，所提取、分析、展示的单TCR/BCR全长核苷酸编码序列，帮助对任何有潜在生物学的TCR/BCR进行后续的湿实验功能研究。

3.大规模。本发明提供了对多种疾病模型样本快速分析单细胞免疫组库的方法，允许在短时间且低计算资源耗费的情况下，探索多种疾病条件下的单细胞免疫组库数据，允许以无技术偏差的方式比较不同样本、不同疾病模型、不同文献来源的单细胞免疫组库样本数据。

附图说明

图1为本发明的整体流程图。

图2为本发明对经过免疫点治疗罹患肠炎病人的样本分析得到的单细胞免疫组库数据中单细胞TCR/BCR测序数据集的可视化结果示意图，提示了病人体内引发肠炎的特定细胞和TCR/BCR克隆型。

图3为本发明对经过免疫点治疗罹患肠炎病人的样本分析得到的单细胞免疫组库数据中单细胞RNA测序数据集的可视化结果示意图，提示了病理情况下的细胞转录组特征。

图4为本发明对经过免疫点治疗罹患肠炎病人的样本分析得到的单细胞免疫组库数据中单细胞TCR/BCR测序数据集耦合单细胞转录组测序数据集的可视化结果示意图，提示了病人体内特定TCR/BCR所在细胞的转录组特征和这些细胞在生理病理过程中的潜在生物学作用。

图5为本发明的数据集分析的示意图。

图6为本发明对经过免疫点治疗但未罹患肠炎病人的样本分析得到的单细胞免疫组库数据中单细胞TCR/BCR测序数据集耦合单细胞转录组测序数据集的可视化结果示意图，作为对照组，进一步提示了免疫检查点疗法造成肠炎的机制。

具体实施方式

下面结合附图及具体实施例对本发明作进一步详细说明。

如图1所示，本发明包括以下步骤：

步骤S1：筛选收集已有的单细胞免疫组库原始数据集，获得多个单细胞免疫组库原始数据集，其中每个单细胞免疫组库原始数据集均由单细胞TCR/BCR测序数据与其相关联的单细胞转录组测序数据组成；其中单细胞TCR/BCR测序数据与其相关联的单细胞转录组测序数据为同一样本所采集的数据，单细胞TCR/BCR测序数据主要由多个TCR/BCR的序列信息组成，单细胞转录组测序数据主要由多个细胞的转录组信息组成。

具体实施中，步骤S1具体为：

筛选收集多篇文献以及各篇文献中多个疾病模型下的单细胞免疫组库原始数据集，将各个单细胞免疫组库原始数据集的数据格式进行统一，检索单细胞免疫组库原始数据集所在数据库，如果来源于SRA或GEO数据库，则该单细胞免疫组库原始数据集的数据格式为SRA数据格式，使用NCBI SRA工具包将该单细胞免疫组库原始数据集的数据格式转换为FASTQ数据格式，获得多个单细胞免疫组库原始数据集；

其中，实施例中共收集了16种组织，包括肺泡灌洗液、结肠上皮驻留淋巴细胞、外周血单个核细胞、多种肿瘤组织；来源于7种疾病状态，包括新冠肺炎、乳腺癌、结肠癌、自发性溃疡性结肠炎、炎症性肠病、免疫检查点治疗引发的结肠炎、肾细胞癌、基底细胞癌等。

步骤S2：选择一个单细胞免疫组库原始数据集，对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据进行标签以及特征提取，获得各个TCR/BCR的细胞标签和TCR/BCR特征，TCR/BCR特征包括各个TCR/BCR的全长编码序列、TCR/BCR的单链所使用的V(D)JC基因、抗原互补决定区3(CDR3)的核苷酸序列和氨基酸序列以及TCR/BCR双链信息；

步骤S2具体为：

S2.1：以Adrienne M Luoma等人在2020年发表的单细胞免疫组库数据挖掘经免疫检查点治疗病人罹患肠炎的新机制中的结肠上皮驻留淋巴细胞的单细胞免疫组库原始数据集为例，样本名为CPIc_C2。单细胞免疫组库原始数据集中的细胞TCR/BCR测序数据主要由多个TCR/BCR的序列信息组成；根据TCR/BCR参考基因组，对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据中各个TCR/BCR的序列信息进行比对，即对各个TCR/BCR的V(D)JC基因的使用做注释，获得各个TCR/BCR的单链所使用的V(D)JC基因以及各个TCR/BCR的单链的完整性；

S2.2：利用Cell Ranger方法从各个TCR/BCR的序列信息中提取各个TCR/BCR的单链的抗原互补决定区3(CDR3)的核苷酸序列和氨基酸序列；

S2.3：利用Cell Ranger方法从各个TCR/BCR的序列信息中提取各个TCR/BCR的细胞标签(Barcode)；对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据中多个TCR/BCR的序列信息进行筛选，过滤不完整的和细胞外的TCR/BCR序列信息，获得筛选后的多个TCR/BCR的序列信息；

S2.4：根据各个TCR/BCR的单链所使用的V(D)JC基因和各个TCR/BCR的细胞标签，从筛选后的各个TCR/BCR的序列信息中提取各个TCR/BCR的单链的从5’非编码区到C基因区的全长编码序列(包含了TCR/BCR全长可变区的序列)；

S2.4具体为：

对单细胞TCR/BCR测序数据进行处理后，获得包含所有TCR/BCR全长序列文件(all_contig_annotations.json)，将各个TCR/BCR的细胞标签以及各个TCR/BCR的单链的完整性作为索引依据，在包含所有TCR/BCR全长序列文件中找到对应TCR/BCR的全长编码序列的位置并提取，获得TCR/BCR双链未匹配的全长编码序列；

再根据各个TCR/BCR的单链所使用的V(D)JC基因中的C基因信息判断TCR/BCR双链未匹配的全长编码序列所来自的具体肽链，最后根据各个TCR/BCR的细胞标签实现TCR/BCR双链未匹配的全长编码序列的双链匹配，获得各个TCR/BCR的单链的从5’非编码区到C基因区的全长编码序列；

S2.5：根据各个TCR/BCR的细胞标签，利用Cell Ranger方法将TCR的α和β链进行匹配以及BCR的轻链和重链进行匹配，获得TCR/BCR双链信息，其中，TCR的单链为α或β链，BCR的单链为轻链或重链。

步骤S3：对各个TCR/BCR的抗原互补决定区3(CDR3)的核苷酸序列进行克隆型定义，即拥有相同CDR3的核苷酸序列的TCR/BCR为同一克隆型，获得各个TCR/BCR的克隆型；

如图2所示，图2的A展示了高富集TCR的克隆型及克隆型ID，克隆型15(clonotype15)为该数据集来源样本中最高富集TCR的克隆型，提示了属于克隆型15的TCR/BCR以及这些TCR/BCR所在细胞可能在患者炎症过程中起到作用。图2的B以网络图的形式展示了通过本发明分析得到的各TCR/BCR和克隆型(clonotype)的富集特征，其中每个点代表一个TCR/BCR，一个小网络代表内部TCR/BCR同属于一个克隆型。图2的C展示了本发明所提取的高富集TCR/BCR的高精度特征信息，例如对于单个TCR/BCR(图2的B中的单个结点)双链V(D)JC基因使用情况和全长编码序列。高精度的信息对可能发挥生理病理作用的特定细胞和TCR/BCR进行研究。图2的D展示了TCR/BCR的高精度信息，对于每个TCR/BCR的V(D)JC基因使用情况，CDR3序列信息，以及全长编码序列，都展示在该表格中，以便提取特定信息，以及允许不同TCR/BCR之间信息的对比。

步骤S4：对当前单细胞免疫组库原始数据集中的单细胞转录组测序数据进行标签提取，获得各个细胞的转录组信息的细胞标签，计算各个细胞的细胞转录组特征，获得单细胞转录组测序数据的基因表达计数矩阵，以及预测单细胞转录组测序数据中各个细胞的类型；

步骤S4具体为：

S4.1：当前单细胞免疫组库原始数据集中的单细胞转录组测序数据主要由多个细胞的转录组信息组成，利用Cell Ranger方法提取单细胞转录组测序数据中各个细胞的转录组信息的细胞标签Barcode；

S4.2：根据人类蛋白表达参考基因组，利用Cell Ranger方法对每个细胞的转录组信息进行比对，计算每个细胞中各基因的表达量并作为当前细胞的细胞转录组特征，由各个细胞的细胞转录组特征组成基因表达计数矩阵(count matrix)；

S4.3：根据各个细胞的细胞转录组特征，利用有监督法或者无监督法对各个细胞的类型进行预测，获得各个细胞的类型。

S4.3具体为：

A利用有监督法预测细胞类型

根据各个细胞的细胞转录组特征，通过SingleR算法计算各个细胞的细胞转录组特征和已知各类细胞的参考转录组特征的相似度(斯皮尔曼相关度)，根据相似度高低，具体实施中，将与当前细胞的细胞转录组特征最相似的参考转录组特征所在细胞的类型定义为当前细胞的细胞类型，从而判断各个细胞的类型；

其中，各类细胞的参考转录组特征为对已知的各类T细胞和B细胞的测序所得到的基因表达信息。

B利用无监督法预测细胞类型

根据各个细胞的细胞转录组特征，计算各个细胞的细胞转录组特征之间的相似度，根据细胞与细胞的细胞转录组特征的相似度划分细胞群，获得当前单细胞转录组测序数据的各个细胞群，利用scanpy方法将各个细胞群和当前单细胞转录组测序数据的全部细胞的细胞转录组特征做对比，获得各个细胞群的特异性表达基因，根据各个细胞群的特异性表达基因来判断各个细胞群所代表的细胞类型，从而获得当前单细胞转录组测序数据的各个细胞的类型。

当某细胞群内细胞表达某基因的均值显著高于或显著低于当前单细胞转录组测序数据的全部细胞表达该基因的均值，该基因被列为该细胞群的特异性表达基因。根据各个细胞群的特异性表达基因来判断各个细胞群所代表的细胞类型，从而获得当前单细胞转录组测序数据的各个细胞的类型。

细胞转录组特征上，图3的A展示了该测试数据集内细胞t-SNE分布，每个点代表数据集中每个细胞，颜色由预测得细胞类型决定。细胞与细胞之间的距离代表了细胞与细胞之间细胞转录组特征的相似度。图3的B展示了各类型细胞在本数据集中的数量。如图3的B所示，效应记忆CD8+T细胞(effector memory CD8 T cell)在本数据集中大量富集，其他类型的T细胞占比较少，这提示了效应记忆CD8+T细胞在免疫点治疗引发肠炎中的重要作用。

同时，通过本发明可以得到数据集单个基因表达水平。图3的C展示了本数据集中效应记忆CD8+T细胞在t-SNE图中的分布。图3的D展示了本数据集中，高表达IFNG基因的细胞。图3的C和图3的D的细胞所在t-SNE图的位置一致，提示了本数据集来源的样本中，效应记忆CD8+T细胞中高表达IFNG基因。图3的E小提琴图同样展示了本数据集来源的样本中，效应记忆CD8+T细胞中高表达IFNG基因。既往研究表明，IFNG基因的蛋白产物IFNγ可促进各类T细胞在多种组织和疾病中发挥炎症作用。该数据集来源的细胞高水平表达IFNG基因提示了患者结肠组织正在经历炎症反应，同时，IFNG在效应记忆CD8+T细胞的表达水平进一步证实了效应记忆CD8+T细胞在这一病理过程中的作用。

此外，图3的F显示PDCD1基因在样本中整体表达水平较低，而主要集中在各类CD4+T细胞(Th1，Follicular helper T cells)中表达。既往研究表明，PDCD1基因产物PD1是免疫点抑制剂，这可能提示了PD1主要由CD4+T细胞在肠炎中发挥抑制的抗炎作用，而以CD8+T细胞为主的特征可能提示了免疫抑制作用被抑制的情况。

步骤S5：根据当前单细胞免疫组库原始数据集的各个TCR/BCR的细胞标签和各个细胞的转录组信息的细胞标签，实现当前单细胞免疫组库原始数据集的TCR/BCR—转录组耦合，获得当前单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息；

步骤S5具体为：

根据当前单细胞免疫组库原始数据集的各个TCR/BCR的细胞标签和各个细胞的转录组信息的细胞标签，将各个细胞的TCR/BCR特征、克隆型、细胞转录组特征和类型进行耦合，过滤仅有TCR/BCR信息或仅有转录组信息的细胞，过滤双链无法匹配的TCR/BCR，获得当前单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息。

图4的A、图4的B展示了本发明分析各个克隆型内的细胞的细胞转录组特征的结果。如图4的A所示，在克隆型15中，效应记忆CD8+T细胞是主要细胞类型，也就是说克隆型15中的TCR主要表达在效应记忆CD8+T细胞中。同时如图4的B所示，在克隆型15中，IFNG基因呈现高表达水平的状态。转录组信息所提供的克隆型15中的细胞类型组成和特定基因表达信息进一步提示了克隆型15中TCR通过接受抗原刺激后触发某种IFNG表达通路的致炎作用机制。

图4的C展示了本发明分析的高富集TCR(属于前10大克隆型的TCR/BCR)所在的效应记忆CD8+T细胞特定类型的细胞转录组特征结果。高富集TCR所在的效应记忆CD8+T细胞表现出IFNG高表达的情况，提示了效应记忆CD8+T细胞在免疫点治疗后发挥致炎的病理作用机制，即能够靶向识别致病抗原的TCR所在的CD8+T细胞将特异扩增为效应记忆T细胞，后通过胞内特定分子通路激活IFNG的表达，最终引起肠炎。

在每个克隆型内部，本发明提供了高精度的单个TCR/BCR所在细胞的细胞转录组信息。如图4的D所示，单个细胞的TCR/BCR信息和单个细胞的细胞类型和特定基因表达量(在此处是IFNG的表达量)被展示。根据该细胞的TCR/BCR特征以及其耦合的细胞转录组特征，该细胞在生理病理中的功能能够被更准确的被推断。根据本发明所提供的高精度全长序列信息，该TCR的功能可被继续研究。

步骤S6：对剩余的单细胞免疫组库原始数据集同样重复步骤S2-S5，获得各个单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息；

步骤S7：数据库搭建；对各个单细胞免疫组库原始数据集进行标注后每个单细胞免疫组库原始数据集对应有唯一的识别码，将各个单细胞免疫组库原始数据集的基本信息(包括数据发表者、发表时间、样本疾病模型等)、各个细胞的TCR/BCR—转录组信息以及各个单细胞免疫组库原始数据集的来源论文均存储到数据库中，通过识别码在数据库中进行数据集的索引；其中，各个细胞的TCR/BCR—转录组信息包括各个细胞的细胞标签、TCR/BCR特征(即全长编码序列、TCR/BCR所使用的V(D)JC基因、抗原互补决定区3(CDR3)的核苷酸序列和氨基酸序列)、TCR/BCR双链信息、克隆型、细胞转录组特征、细胞类型。

具体实施中，各个单细胞免疫组库原始数据集的基因表达计数矩阵存储为HDF5文件，各个单细胞免疫组库原始数据集中的各个细胞的细胞标签、TCR/BCR特征、克隆型和细胞类型存储为二进制格式，存储在云服务器PostgreSQL数据库引擎提供的数据库服务上。识别码作为各个数据集的索引，可以用来快速检索和提取数据库中各个数据集的TCR/BCR特征、克隆型、细胞转录组特征。

图5展示了数据库网页主页和在数据库中检索该数据集的方法。通过选择研究特定样本(本例中所选择的样本为CPIc_C2)，由本发明得到的各项结果信息将被立刻展示。图5A展示了本发明索引特定样本数据集的方法，图5B展示了所索引样本来源文章的基本信息，图5C展示了所索引样本的基本信息。

步骤S8：数据库的可视化；利用D3.js框架搭建终端浏览器界面，终端浏览器界面与数据库相结合，在终端浏览器界面中选择单细胞免疫组库原始数据集，根据所选择的单细胞免疫组库原始数据集的识别码调用数据库中所选择的单细胞免疫组库原始数据集的所有信息，所选择的单细胞免疫组库原始数据集的所有信息在终端浏览器界面中进行可交互的可视化。

步骤S8中，使用t-SNE方法对所选择的单细胞免疫组库原始数据集的基因表达计数矩阵进行降维处理后，获得所选择的单细胞免疫组库原始数据集的全部细胞的降维后的细胞转录组特征，以二维散点图的方式展示全部细胞的降维后的细胞转录组特征；使用小提琴图显示各细胞类型的各个基因的表达量；使用直方图显示各个细胞类型的数量；使用网络图显示各个TCR/BCR的克隆型；使用直方图显示各个克隆型的细胞数量；使用复合直方图显示各个克隆型内部的各细胞类型的细胞数量；获得多种可视化图；通过细胞标签对多种可视化图进行可视化交互，实现了有相同克隆型的细胞的细胞转录组特征的提取和比较以及单细胞的TCR/BCR特征-细胞转录组特征的耦合。

本发明提供直接、快速的跨样本对比。在主页索引相对应的对照组数据(样本：CPIc_NC1,样本来源：接受免疫点治疗但未表现肠炎症状的病人结肠组织)并提交，可以快速查看可视化结果，如图6所示。图6的A显示，在对照组样本中高富集的克隆型所包含的TCR数量为5-30，该数量远小于肠炎病人结肠组织内的高富集克隆型TCR数量，提示了肠炎病人的结肠组织内存在的严重炎症状态。图6的B显示，高富集的TCR主要表达在终末记忆CD8+T细胞，而不是效应记忆CD8+T细胞，和肠炎病人形成显著对比。既往研究表明，终末记忆T细胞可能代表了历史炎症而非近发炎症，因此进一步证明了肠炎病人结肠组织内的严重炎症状态，以及效应记忆CD8+T细胞在免疫点治疗引发肠炎病理过程中的重要机制。

另外，如图6的C所示，对照组数据集中，最高程度富集TCR的克隆型47的细胞IFNG基因表达水平明显低于肠炎病人克隆型15中的细胞IFNG的表达，提示了肠炎病人结肠组织中属于克隆型15的细胞和TCR的致炎作用。然而，如图6的D所示，该数据集中虽然效应记忆CD8+T细胞数量不高，但为数不多的高富集TCR所在的效应记忆CD8+T细胞表达了较高水平的IFNG。考虑到既往研究发现的免疫点疗法对结肠组织的副作用，这一现象可能提示了该对照组样本来源的病人正在肠炎前期未发病阶段，可能的原因是引发炎症反应的克隆型和TCR，以及其所在的效应记忆CD8+T细胞还未大量扩增。在这一数据集中，克隆型135的细胞类型组成(效应记忆CD8+T细胞占主要类型，提示近发炎症)以及其所包含的TCR数量(有扩增但扩增量还不大)提示了本猜想的正确性。

Claims (7)

1.一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，包括以下步骤：

步骤S1：筛选收集已有的单细胞免疫组库原始数据集，获得多个单细胞免疫组库原始数据集，其中每个单细胞免疫组库原始数据集均由单细胞TCR/BCR测序数据与其相关联的单细胞转录组测序数据组成；

步骤S2：选择一个单细胞免疫组库原始数据集，对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据进行标签以及特征提取，获得各个TCR/BCR的细胞标签和TCR/BCR特征，TCR/BCR特征包括各个TCR/BCR的全长编码序列、TCR/BCR所使用的V(D)JC基因、抗原互补决定区3的核苷酸序列和氨基酸序列以及TCR/BCR双链信息；

步骤S3：对各个TCR/BCR的抗原互补决定区3的核苷酸序列进行克隆型定义，获得各个TCR/BCR的克隆型；

步骤S4：对当前单细胞免疫组库原始数据集中的单细胞转录组测序数据进行标签提取，获得各个细胞的转录组信息的细胞标签，计算各个细胞的细胞转录组特征，获得单细胞转录组测序数据的基因表达计数矩阵，以及预测单细胞转录组测序数据中各个细胞的类型；

步骤S5：根据当前单细胞免疫组库原始数据集的各个TCR/BCR的细胞标签和各个细胞的转录组信息的细胞标签，实现当前单细胞免疫组库原始数据集的TCR/BCR—转录组耦合，获得当前单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息；

步骤S6：对剩余的单细胞免疫组库原始数据集同样重复步骤S2-S5，获得各个单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息；

步骤S7：数据库搭建；对各个单细胞免疫组库原始数据集进行标注后每个单细胞免疫组库原始数据集对应有唯一的识别码，将各个单细胞免疫组库原始数据集的基本信息、各个细胞的TCR/BCR—转录组信息以及各个单细胞免疫组库原始数据集的来源论文均存储到数据库中，通过识别码在数据库中进行数据集的索引；

步骤S8：数据库的可视化；搭建终端浏览器界面，终端浏览器界面与数据库相结合，在终端浏览器界面中选择单细胞免疫组库原始数据集，根据所选择的单细胞免疫组库原始数据集的识别码调用数据库中所选择的单细胞免疫组库原始数据集的所有信息，所选择的单细胞免疫组库原始数据集的所有信息在终端浏览器界面中进行可交互的可视化。

2.根据权利要求1所述的一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，所述步骤S2具体为：

S2.1：选择一个单细胞免疫组库原始数据集，其中单细胞TCR/BCR测序数据主要由多个TCR/BCR的序列信息组成；根据TCR/BCR参考基因组，对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据中各个TCR/BCR的序列信息进行比对，获得各个TCR/BCR的单链所使用的V(D)JC基因以及各个TCR/BCR的单链的完整性；

S2.2：从各个TCR/BCR的序列信息中提取各个TCR/BCR的单链的抗原互补决定区3的核苷酸序列和氨基酸序列；

S2.3：从各个TCR/BCR的序列信息中提取各个TCR/BCR的细胞标签；对当前单细胞免疫组库原始数据集的单细胞TCR/BCR测序数据进行筛选过滤，获得筛选后的多个TCR/BCR的序列信息；

S2.4：根据各个TCR/BCR的单链所使用的V(D)JC基因和各个TCR/BCR的细胞标签，从筛选后的各个TCR/BCR的序列信息中提取各个TCR/BCR的单链的从5’非编码区到C基因区的全长编码序列；

S2.5：根据各个TCR/BCR的细胞标签，将TCR的α和β链进行匹配以及BCR的轻链和重链进行匹配，获得TCR/BCR双链信息。

3.根据权利要求1所述的一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，所述步骤S4具体为：

S4.1：当前单细胞免疫组库原始数据集中的单细胞转录组测序数据主要由多个细胞的转录组信息组成，提取单细胞转录组测序数据中各个细胞的转录组信息的细胞标签；

S4.2：根据人类蛋白表达参考基因组，对每个细胞的转录组信息进行比对，计算每个细胞中各基因的表达量并作为当前细胞的细胞转录组特征，由各个细胞的细胞转录组特征组成基因表达计数矩阵；

S4.3：根据各个细胞的细胞转录组特征，利用有监督法或者无监督法对各个细胞的类型进行预测，获得各个细胞的类型。

4.根据权利要求3所述的一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，所述S4.3具体为：

A利用有监督法预测细胞类型

根据各个细胞的细胞转录组特征，计算各个细胞的细胞转录组特征和已知各类细胞的参考转录组特征的相似度，根据相似度高低判断各个细胞的类型；

B利用无监督法预测细胞类型

根据各个细胞的细胞转录组特征，计算各个细胞的细胞转录组特征之间的相似度，根据细胞与细胞的细胞转录组特征的相似度划分细胞群，获得当前单细胞转录组测序数据的各个细胞群，将各个细胞群和当前单细胞转录组测序数据的全部细胞的细胞转录组特征做对比，获得各个细胞群的特异性表达基因。

5.根据权利要求1所述的一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，所述步骤S5具体为：

根据当前单细胞免疫组库原始数据集的各个TCR/BCR的细胞标签和各个细胞的转录组信息的细胞标签，将各个细胞的TCR/BCR特征、克隆型、细胞转录组特征和类型进行耦合，过滤仅有TCR/BCR信息或仅有转录组信息的细胞，过滤双链无法匹配的TCR/BCR，获得当前单细胞免疫组库原始数据集中的各个细胞的TCR/BCR—转录组信息。

6.根据权利要求1所述的一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，所述步骤S8中，使用t-SNE方法对所选择的单细胞免疫组库原始数据集的基因表达计数矩阵进行降维处理后，获得所选择的单细胞免疫组库原始数据集的全部细胞的降维后的细胞转录组特征，以二维散点图的方式展示全部细胞的降维后的细胞转录组特征；使用小提琴图显示各细胞类型的各个基因的表达量；使用直方图显示各个细胞类型的数量；使用网络图显示各个TCR/BCR的克隆型；使用直方图显示各个克隆型的细胞数量；使用复合直方图显示各个克隆型内部的各细胞类型的细胞数量；获得多种可视化图；通过细胞标签对多种可视化图进行可视化交互，实现了有相同克隆型的细胞的细胞转录组特征的提取和比较以及单细胞的TCR/BCR特征-细胞转录组特征的耦合。

7.根据权利要求2所述的一种基于单细胞免疫组库数据的单细胞水平耦合可视化方法，其特征在于，所述S2.4具体为：

对单细胞TCR/BCR测序数据进行处理后，获得包含所有TCR/BCR全长序列文件，将各个TCR/BCR的细胞标签以及各个TCR/BCR的单链的完整性作为索引依据，在包含所有TCR/BCR全长序列文件中找到对应TCR/BCR的全长编码序列的位置并提取，获得TCR/BCR双链未匹配的全长编码序列；

再根据各个TCR/BCR的单链所使用的V(D)JC基因中的C基因信息判断TCR/BCR双链未匹配的全长编码序列所来自的具体肽链，最后根据各个TCR/BCR的细胞标签实现TCR/BCR双链未匹配的全长编码序列的双链匹配，获得各个TCR/BCR的单链的从5’非编码区到C基因区的全长编码序列。

Images