40. 免疫受体分析#
关键要点
VDJ 测序测量 AIR 的序列信息，从而让我们得以洞察 B 细胞和 T 细胞的功能。
为了检测双细胞，应当识别或过滤掉具有多个 AIR 的细胞。
根据具体分析的不同，我们还会额外过滤掉 AIR 信息不完整的细胞。
环境设置
安装 conda:
在创建环境之前，请确保 conda 已安装在您的系统中。
保存 yml 内容:
将 yml 选项卡中的内容复制到名为
environment.yml的文件中。
创建环境:
打开终端或命令提示符。
运行以下命令：
conda env create -f environment.yml
激活环境:
创建好环境后，使用以下命令激活它：
conda activate <environment_name>
替换
<environment_name>，名称就是在environment.yml文件中指定的那个。在 yml 文件里，它看起来像这样：name: <environment_name>
验证安装:
通过运行以下命令，检查环境是否创建成功：
conda env list
40.1. IR 在细胞中的作用#
免疫受体（Immune Receptor，IR）是识别入侵机体的潜在有害抗原和毒素的关键。这些受体通常位于特化细胞的膜上，这些细胞能识别有害物质并启动保护机制。针对识别外来因子的特定结构，存在不同种类的 IR：
模式识别受体（PRRs）：病原体相关分子模式（PAMPS）
杀伤激活/抑制受体（KARs/KIRs）：宿主细胞异常
补体受体：补体蛋白
Fc 受体：表位-抗体复合物
细胞因子受体：细胞因子
B 细胞受体（BCR）：表位
T 细胞受体（TCR）：与主要组织相容性复合体（MHC）结合的线性表位
40.2. 适应性免疫系统的IRs#
在适应性免疫系统中，有两大类带有 IR 的淋巴细胞，它们分别在胸腺（Thymus）和骨髓（Bone marrow）中产生，称为 T 细胞和 B 细胞 [Cooper and Alder, 2006]。这两类细胞的适应性免疫受体（AIR）能检测来自不同病原体和肿瘤细胞的抗原。不过，它们与抗原相互作用的方式不同：B 细胞受体（BCR）直接识别可溶性或膜结合的表位，而 T 细胞受体则与结合在表面蛋白 MHC 上的线性肽（pMHC）相互作用——后者通过把细胞内容物呈递到细胞表面，使免疫系统得以“窥视”细胞内部。被抗原激活后，B 细胞和 T 细胞会执行多种功能，如对抗病原体、调节免疫应答，或通过增殖形成记忆。
图40.1 两种 AIR 类型的示意图。创建工具： BioRender.com.#
40.3. 通过V(D)J重组形成适应性免疫受体#
AIR 是由两条链组成的蛋白复合体。根据 T 细胞类型，这两条链分别称为 α 链和 β 链，或 γ 链和 δ 链。类似地，BCR 的链称为重链和轻链（亚型：κ 和 λ）。其中 α、γ、κ、λ 链由各自的 V 基因和 J 基因形成，而其余的链还额外包含 D 基因。为简单起见，本书中我们统一称之为 VJ 链和 VDJ 链。
适应性 IR 谱系高度多态，以覆盖庞大的抗原空间——换言之，TCR 与 BCR 序列极为多样。据估计，不同 TCR 的总数可达 1020 ；而一个人体内携带 107 个不同的序列 [Zarnitsyna et al., 2013]。同样，自然会出现多达 1018 个不同的 BCR [Briney et al., 2019]。一个个体的全部受体集合构成其 AIR 库（repertoire）。
序列多样性由 V(D)J 重组机制引入：
组合多样性 由不同的 V、(D、) 和 J 基因的组合，以及不同 VJ 链与 VDJ 链的组合而引入。
连接多样性 通过在基因连接处插入核苷酸来实现。
图40.2 以 αβ-TCR 为例的 VDJ 重组。创建工具： BioRender.com.#
此外，在细胞快速增殖期间，BCR 的结合位点会发生突变，这一过程称为体细胞超突变（somatic hypermutation），用于亲和力成熟。
在最终形成的 AIR 中，每条链上都检测到三个高变区——互补决定区（Complementary Determining Region，CDR）1–3，AIR 正是在这里与其靶标相互作用。CDR1 和 CDR2 由 V 基因编码，而 CDR3 横跨 V、(D、) 和 J 基因的交界处，因此是 AIR 中最具多样性的部分。因此，人们通常认为，一个 AIR 的特异性主要由其更具多样性的 VDJ 链的 CDR3 区域决定。
40.4. VDJ-sequencing#
40.4.1. 细胞分离#
在对单细胞数据进行定量之前，应当有针对性地分选并分离出感兴趣的细胞。常见的选择包括：
荧光激活细胞分选（FACS） 是一种基于流式细胞术的方法，能够依据荧光探针，从原始细胞悬液中标记出感兴趣的细胞。细胞悬液由快速流动的液流承载，这股细胞流通过振动机制被分解成单个液滴。就在液流断裂成液滴之前，它会流经一个荧光检测站，在那里测量每个细胞的荧光信号。随后还可以给液滴充电，以便进一步分选。
磁激活细胞分选（MACS） 可以使用连接在磁珠上的抗体、酶、凝集素（lectin）或链霉亲和素来标记目标细胞。一旦细胞从原始悬液中被标记，就施加磁场来吸引磁珠，并把悬液中其余的细胞弃去；关闭磁场后即可收集到目标细胞。这种方法的一个优点是，也能收集没有特异标记可供标记的目标细胞：这种情况下，用一组（鸡尾酒式的）标记去标记非目标细胞，待磁场把连接了磁珠的非目标细胞捕获后，再把感兴趣的细胞洗脱出来收集。
激光捕获显微切割（LCM） 能够从显微镜标本中提取细胞群或单细胞，而不损伤周围组织。进行 LCM 所需的部件包括：一台倒置显微镜、一个激光控制单元、一个用于稳定载物台的显微操纵杆、一台 CCD 相机和一台彩色显示器。LCM 的基本思路是：通过目视检测目标细胞的形态特征来标定细胞，固定载物台，再用激光脉冲熔化薄薄的热塑膜，从而在不损伤周围组织的情况下去除感兴趣的细胞。
微流控 是一种用途广泛的方法，即使是纳升量级的少量原始悬液也能处理。微流控方法有多种，包括基于细胞亲和层析的微流控、基于细胞物理特性的微流控、基于免疫磁珠的微流控，以及基于某些细胞类型介电特性进行分离的微流控。最常用的微流控方法是层析分离：用经过修饰、包含必要抗体的芯片作为固定相，以在流动相中捕获目标细胞；当缓冲液从芯片流走后，再用一种溶液把结合在抗体上的细胞分离下来，收集起来供进一步分析 [Hu et al., 2016].
40.4.2. 免疫受体测序#
从单细胞分离物中辨别 V(D)J 链的一种常见方法，是基于全长单细胞 RNA 测序，对不同链的序列进行计算重建；其中基于 5′ 端 RNA 模板的 Smart-seq2 是应用最广泛的方案之一。在计算方法方面，TRAPeS、TraCeR 和 VDJPuzzle 通常用于基于 scRNA-seq 数据重建 TCR 序列，而 BALDR [Upadhyay et al., 2018], BASIC [Canzar et al., 2016] 和 BraCer [Lindeman et al., 2018] 则被证明能够稳健地恢复 BCR 序列。不过，这些方法往往会忽略 V(D)J 区域中重组产物和可变剪接产物的全貌。为解决这一问题，出现了一些替代方法，例如 RAGE-seq：它基于为免疫受体测序设计的 PCR 模板来捕获特定的 TCR 和 BCR 片段，并使用长读长的 Oxford Nanopore 捕获完整序列，而其余的 cDNA 则用例如 Illumina 提供的短读长方案来处理。 [Singh et al., 2019].
40.5. AIR 库分析#
VDJ 测序为我们提供核苷酸序列、进而提供配对两条链的 AIR 蛋白序列，由此可确定 V、(D、) J、C 基因以及 CDR3 序列。总体而言，AIR 序列决定了单个 B 细胞和 T 细胞的特异性。因此，VDJ 测序所获得的信息为我们提供了细胞功能的指标，而功能又与 AIR 的靶抗原直接关联。这使我们能以三种主要方式利用 AIR 信息：
表型分析：我们可以通过识别具有相同或相似 AIR（从而具有相同特异性）的细胞来对免疫细胞分组。有了这些分组，就可以观察特定疾病的细胞在不同条件下如何反应（例如受刺激时的转录变化）、免疫细胞是否发生了增殖，或免疫库的多样性在免疫应答后如何变化。
序列分析：在识别出若干 AIR 分组（例如在其他模态中检测到的反应性聚类）之后，我们可以提取其序列属性，如 V、D、J 基因的使用情况或富集的序列基序（motif），它们与特定疾病或疗法相关。
特异性推断（Specificity Inference）：最后，我们可以利用序列，通过数据库查询、序列距离或预测器，把 AIR 与其靶抗原匹配起来。这能直接识别出对特定传染病、肿瘤或自身抗原产生反应的细胞。
40.6. 数据集#
为了展示对 IR 进行预处理和分析的不同方法，我们将使用 Haniffa 实验室公布的数据集 [Stephenson et al., 2021]，它包含来自 130 名患者、在严重急性呼吸综合征冠状病毒 2（SARS-CoV-2）感染背景下的超过 750,000 个外周血单个核细胞（PBMC）的转录组测序数据。
该数据集的患者来自三个来源（纽卡斯尔、剑桥和伦敦），疾病严重程度各不相同（无症状、轻度、中度、重度和危重），同时还包含阴性对照（健康者、患其他严重呼吸道疾病者，以及静脉注射脂多糖以模拟全身炎症反应的健康者）。此外，还提供了年龄、性别、吸烟状况等患者层面的信息。
选择该数据集，是因为它是一个包含 VDJ 测序的大规模单细胞数据集的范例，且为社区所熟知。它包含超过 150,000 个 B 细胞和 200,000 个带 IR 注释的 T 细胞，我们将在本分析教程中加以考察。
要下载数据，请运行以下命令：
path_data = "data/"
path_bcr_input = f"{path_data}/BCR_00_read_aligned.csv"
path_tcr_input = f"{path_data}/TCR_00_read_aligned.tsv"
! wget -O $path_bcr_input -nc https://figshare.com/ndownloader/files/35574338
! wget -O $path_tcr_input -nc https://figshare.com/ndownloader/files/35574539
--2022-06-07 18:30:52-- https://figshare.com/ndownloader/files/35574338
Resolving figshare.com (figshare.com)... 54.72.163.193, 52.50.42.102, 2a05:d018:1f4:d000:7421:2135:71b2:6a10, ...
Connecting to figshare.com (figshare.com)|54.72.163.193|:443... connected.
HTTP request sent, awaiting response... 302 Found
Location: https://s3-eu-west-1.amazonaws.com/pfigshare-u-files/35574338/bcr_cellranger.csv?X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Credential=AKIAIYCQYOYV5JSSROOA/20220607/eu-west-1/s3/aws4_request&X-Amz-Date=20220607T163052Z&X-Amz-Expires=10&X-Amz-SignedHeaders=host&X-Amz-Signature=b18cb63ca03820e07cb8b3d9375466ca6405e655822269f5ae36d7449a28b810 [following]
--2022-06-07 18:30:52-- https://s3-eu-west-1.amazonaws.com/pfigshare-u-files/35574338/bcr_cellranger.csv?X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Credential=AKIAIYCQYOYV5JSSROOA/20220607/eu-west-1/s3/aws4_request&X-Amz-Date=20220607T163052Z&X-Amz-Expires=10&X-Amz-SignedHeaders=host&X-Amz-Signature=b18cb63ca03820e07cb8b3d9375466ca6405e655822269f5ae36d7449a28b810
Resolving s3-eu-west-1.amazonaws.com (s3-eu-west-1.amazonaws.com)... 52.218.62.171
Connecting to s3-eu-west-1.amazonaws.com (s3-eu-west-1.amazonaws.com)|52.218.62.171|:443... connected.
HTTP request sent, awaiting response... 200 OK
Length: 84451806 (81M) [text/csv]
Saving to: ‘data//BCR_00_read_aligned.csv’
100%[======================================>] 84,451,806 39.1MB/s in 2.1s
2022-06-07 18:30:54 (39.1 MB/s) - ‘data//BCR_00_read_aligned.csv’ saved [84451806/84451806]
--2022-06-07 18:30:55-- https://figshare.com/ndownloader/files/35574539
Resolving figshare.com (figshare.com)... 54.72.163.193, 52.50.42.102, 2a05:d018:1f4:d000:7421:2135:71b2:6a10, ...
Connecting to figshare.com (figshare.com)|54.72.163.193|:443... connected.
HTTP request sent, awaiting response... 302 Found
Location: https://s3-eu-west-1.amazonaws.com/pfigshare-u-files/35574539/TCR_mergedUpdated.tsv?X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Credential=AKIAIYCQYOYV5JSSROOA/20220607/eu-west-1/s3/aws4_request&X-Amz-Date=20220607T163055Z&X-Amz-Expires=10&X-Amz-SignedHeaders=host&X-Amz-Signature=141ab989ae794394fb6c441d50c0ea6771f96a8d048a8a4c400c10dba267b97c [following]
--2022-06-07 18:30:55-- https://s3-eu-west-1.amazonaws.com/pfigshare-u-files/35574539/TCR_mergedUpdated.tsv?X-Amz-Algorithm=AWS4-HMAC-SHA256&X-Amz-Credential=AKIAIYCQYOYV5JSSROOA/20220607/eu-west-1/s3/aws4_request&X-Amz-Date=20220607T163055Z&X-Amz-Expires=10&X-Amz-SignedHeaders=host&X-Amz-Signature=141ab989ae794394fb6c441d50c0ea6771f96a8d048a8a4c400c10dba267b97c
Resolving s3-eu-west-1.amazonaws.com (s3-eu-west-1.amazonaws.com)... 52.218.62.171
Connecting to s3-eu-west-1.amazonaws.com (s3-eu-west-1.amazonaws.com)|52.218.62.171|:443... connected.
HTTP request sent, awaiting response... 200 OK
Length: 247085830 (236M) [text/tab-separated-values]
Saving to: ‘data//TCR_00_read_aligned.tsv’
100%[======================================>] 247,085,830 42.4MB/s in 5.8s
2022-06-07 18:31:01 (40.9 MB/s) - ‘data//TCR_00_read_aligned.tsv’ saved [247085830/247085830]
40.7. 加载数据#
在本教程中，我们主要使用两个 Python 包来进行加载、细胞级排序和可视化：
Scanpy:用于单细胞分析的通用软件包(theislab/scanpy, [Wolf et al., 2018])
Scirpy：用于免疫受体分析的 scanpy 扩展（scverse/scirpy, [Sturm et al., 2020])
这里，我们仅用 Scirpy 来演示 IR 分析。不过，也存在若干功能类似的工具，例如 immunarch（R， [ImmunoMind Team, 2019]), scRepertoire (R, [Borcherding et al., 2020])、dandelion（R， [Stephenson et al., 2021])、Platypus（R， [Yermanos et al., 2021]）综述见 [Valkiers et al., 2022].
警告
Scirpy 更改了 其数据结构的格式 （在 v0.13 中）。虽然整体分析流程没有改变，但本章中展示的一些输出可能不再准确。
See scirpy 发布说明 关于这个变化的更多细节。在更新本章之前，请参考 scirpy 官方文档.
import warnings
warnings.filterwarnings(
"ignore",
".*IProgress not found*",
)
warnings.simplefilter(action="ignore", category=FutureWarning)
import pandas as pd
import scanpy as sc
import scirpy as ir
warnings.simplefilter(action="ignore", category=pd.errors.DtypeWarning)
让我们设定数据输入和输出路径。
path_bcr_out = f"{path_data}/BCR_01_preprocessed.h5ad"
path_tcr_csv = f"{path_data}/TCR_00_read_aligned.csv"
path_tcr_out = f"{path_data}/TCR_01_preprocessed.h5ad"
40.7.1. 原始数据#
我们先查看 Cell Ranger 流程的原始输出，以便更好地理解所处理的数据。我们将载入 filtered_contig_annotations.csv" 以查看其内容。每行将代表一个序列的一个测量。
df_bcr_raw = pd.read_csv(path_bcr_input, index_col=0)
# The column 'productive' contains mixed data types which are not compatible with downstream tools.
# We correct them by casting them to strings.
df_bcr_raw["productive"] = df_bcr_raw["productive"].astype(str)
print(f"Total measurements: {len(df_bcr_raw)}")
df_bcr_raw.head(5)
Total measurements: 373670
| barcode | is_cell | contig_id | high_confidence | length | chain | v_gene | d_gene | j_gene | c_gene | full_length | productive | cdr3 | cdr3_nt | reads | umis | raw_clonotype_id | raw_consensus_id | patient_id | newBarcode | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | AAACCTGAGCGATATA-MH9143270 | True | AAACCTGAGCGATATA-MH9143270_contig_1 | True | 646 | IGL | IGLV1-44 | None | IGLJ1 | IGLC2 | True | True | CAAWDDSLSASYVF | TGTGCAGCATGGGATGACAGCCTGAGTGCCTCTTATGTCTTC | 154 | 14 | clonotype9 | clonotype9_consensus_2 | COVID-002 | NaN |
| 1 | AAACCTGAGCGATATA-MH9143270 | True | AAACCTGAGCGATATA-MH9143270_contig_2 | True | 579 | IGH | IGHV3-33 | IGHD3-3 | IGHJ4 | IGHM | True | True | CARDNRVYYDFWSGYPDYW | TGTGCGAGAGACAACCGAGTCTATTACGATTTTTGGAGTGGTTATC... | 144 | 13 | clonotype9 | clonotype9_consensus_1 | COVID-002 | NaN |
| 2 | AAACCTGTCTACCAGA-MH9143270 | True | AAACCTGTCTACCAGA-MH9143270_contig_1 | True | 574 | IGH | IGHV1-69D | IGHD2-15 | IGHJ5 | IGHM | True | True | CARDRRSAYCSGGSCWGGDWFDPW | TGTGCGAGAGATCGTCGTTCAGCTTATTGTAGTGGTGGTAGCTGCT... | 564 | 42 | clonotype10 | clonotype10_consensus_1 | COVID-002 | NaN |
| 3 | AAACCTGTCTACCAGA-MH9143270 | True | AAACCTGTCTACCAGA-MH9143270_contig_2 | True | 744 | IGL | IGLV1-51 | None | IGLJ3 | IGLC3 | True | False | CGTWDSSLSAWVF | TGCGGAACATGGGATAGCAGCCTGAGTGCTTGGGTGTTC | 188 | 21 | clonotype10 | None | COVID-002 | NaN |
| 4 | AAACGGGCATTTGCTT-MH9143270 | True | AAACGGGCATTTGCTT-MH9143270_contig_1 | True | 565 | IGH | IGHV3-33 | None | IGHJ4 | IGHM | True | True | CANPTRPYSSSWWYFDYW | TGTGCGAACCCCACCCGTCCATATAGCAGCAGCTGGTGGTACTTTG... | 234 | 20 | clonotype11 | clonotype11_consensus_1 | COVID-002 | NaN |
该表载有与我们预处理或分析有关的下列条目。有关 Cell Ranger 具体条目的更详细解释，可在此查阅（https://support.10xgenomics.com/single-cell-vdj/software/pipelines/latest/output/annotation）。在其他数据格式中，条目（如 AIRR）的名称会有所不同，但底层信息基本相似。
barcode：测得该 contig 的细胞标签
is_cell：表示该 Barcode 是否与某个细胞相关
high_confidence：该测量结果为 IR 的置信度
chain：IR 的链（例如 TRA：T 细胞受体 α 链；IGH：免疫球蛋白重链）
{v,d,j,c}_gene：用于形成 IR 特定片段的基因
full_length：是否捕获了完整的 IR（见下文）
productive：IR 是否为有效（productive）（见下文）
cdr3{_nt}：IR 链 CDR3 的{核苷酸}序列
patient_id：标识不同患者的 ID
40.7.2. 有效的 AIR#
尽管我们能检测到 AIR 序列，但它未必是“有效的（productive）”，也就是说未必能形成一个有效的 AIR。那些不能产生功能性 AIR 的序列会被标记为“无效（non-productive）”。在用 Scirpy 等工具加载数据时，这些序列通常会被忽略，不用于任何下游分析。有效（productive）免疫受体由 10x Genomics 定义为 此处 ：
从 V 基因跨越到 J 基因的序列
含有起始 密码子 （位于前导区）
在起始密码子的阅读框内包含一个 CDR3。
在 V–J 跨度内不包含终止密码子
示例 1：序列未跨越完整的 IR（V 到 J 基因），如 full_length 列所示。可以看到此处存在基因注释缺失，尤其是 V 基因和 J 基因。
columns = [
"barcode",
"v_gene",
"d_gene",
"c_gene",
"j_gene",
"productive",
"full_length",
]
df_bcr_raw[~df_bcr_raw["full_length"]][columns].head()
| barcode | v_gene | d_gene | c_gene | j_gene | productive | full_length | |
|---|---|---|---|---|---|---|---|
| 34 | AACTCCCCAAGTACCT-MH9143270 | None | None | IGLC3 | IGLJ5 | False | False |
| 38 | AACTCTTCAGATCGGA-MH9143270 | IGKV1-39 | None | IGKC | None | False | False |
| 43 | AACTGGTCAGGGAGAG-MH9143270 | IGLV4-3 | None | IGLC2 | None | False | False |
| 65 | ACACCAAGTCGATTGT-MH9143270 | None | IGHD5-12 | IGHD | IGHJ4 | False | False |
| 76 | ACATCAGCATACTACG-MH9143270 | None | IGHD3-3 | IGHD | IGHJ6 | False | False |
示例 2：Contig 表达全长，但没有可识别的 CDR3。
columns += ["cdr3", "cdr3_nt"]
df_bcr_raw[(df_bcr_raw["productive"] == "False") & (df_bcr_raw["full_length"])][
columns
].head(5)
| barcode | v_gene | d_gene | c_gene | j_gene | productive | full_length | cdr3 | cdr3_nt | |
|---|---|---|---|---|---|---|---|---|---|
| 3 | AAACCTGTCTACCAGA-MH9143270 | IGLV1-51 | None | IGLC3 | IGLJ3 | False | True | CGTWDSSLSAWVF | TGCGGAACATGGGATAGCAGCCTGAGTGCTTGGGTGTTC |
| 13 | AAAGCAAGTCGAAAGC-MH9143270 | IGHV5-51 | IGHD2-15 | IGHD | IGHJ4 | False | True | None | None |
| 23 | AACGTTGAGCTGCGAA-MH9143270 | IGHV3-30 | None | IGHM | IGHJ4 | False | True | CESGRATITMIVVTDLDYW | TGCGAAAGCGGGCGGGCTACGATTACTATGATAGTAGTGACTGACC... |
| 27 | AACTCAGTCAAAGTAG-MH9143270 | IGHV3-23 | IGHD3-16 | IGHD | IGHJ4 | False | True | CAKAWIGLDSEIYYDYIWGSYRIDFDYW | TGTGCGAAAGCTTGGATCGGACTCGACTCCGAGATTTATTATGATT... |
| 30 | AACTCCCAGGTGCACA-MH9143270 | IGKV6-21 | None | IGKC | IGKJ3 | False | True | None | None |
我们现在将加载 TCR 数据并将其从标签分隔格式转换为逗号分隔格式，供以后使用。由于数据集内多个样本的 Barcode 存在重叠，我们进一步复制了唯一 CellID 作为 barcode 注释。除了前面讨论过的 contig 注释外，这些数据还包含患者层面（如年龄、结局……）和细胞层面（细胞类型）的注释。
df_tcr_raw = pd.read_csv(path_tcr_input, sep="\t")
df_tcr_raw["barcode"] = df_tcr_raw.pop("CellID")
df_tcr_raw.to_csv(path_tcr_csv)
print(f"Total measurements: {len(df_tcr_raw)}")
df_tcr_raw.head()
Total measurements: 547630
| Centre | barcode | is_cell | contig_id | high_confidence | length | chain | v_gene | d_gene | j_gene | ... | Days_from_onset | time_after_LPS | Worst_Clinical_Status | Outcome | leiden | consensus | initial_clustering | study_id | AgeRange | Age | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | Newcastle | AAACCTGAGAAACCTA-MH9179824 | True | AAACCTGAGAAACCTA-1_contig_2 | True | 512 | TRB | TRBV5-1 | None | TRBJ2-3 | ... | 8 | NaN | Moderate | Home | 26 | CD8 | CD8 | COVID-014 | (50, 59] | 54 |
| 1 | Newcastle | AAACCTGAGAAACCTA-MH9179824 | True | AAACCTGAGAAACCTA-1_contig_1 | True | 573 | TRA | TRAV8-4 | None | TRAJ31 | ... | 8 | NaN | Moderate | Home | 26 | CD8 | CD8 | COVID-014 | (50, 59] | 54 |
| 2 | Newcastle | AAACCTGAGAACAATC-MH9143324 | True | AAACCTGAGAACAATC-1_contig_1 | True | 635 | TRB | TRBV19 | None | TRBJ1-5 | ... | 16 | NaN | Moderate | Home | 25 | gdT | gdT | COVID-029 | (20, 29] | 26 |
| 3 | Newcastle | AAACCTGAGAACTCGG-MH8919329 | True | AAACCTGAGAACTCGG-1_contig_1 | True | 499 | TRA | TRAV34 | None | TRAJ30 | ... | 5 | NaN | Death | Death | 23 | Treg | Treg | COVID-045 | (70, 79] | 77 |
| 4 | Newcastle | AAACCTGAGAACTCGG-MH8919329 | True | AAACCTGAGAACTCGG-1_contig_2 | True | 496 | TRB | TRBV5-1 | None | TRBJ2-5 | ... | 5 | NaN | Death | Death | 23 | Treg | Treg | COVID-045 | (70, 79] | 77 |
5 rows × 47 columns
我们会发现这份数据已经处理过：所有 contig 都是全长的，且代表一个有效（productive）的 IR。
df_tcr_raw["full_length"].value_counts()
True 547630
Name: full_length, dtype: int64
df_tcr_raw["productive"].value_counts()
True 547630
Name: productive, dtype: int64
40.8. 细胞对齐数据#
为了在细胞层面进行分析，接下来我们需要通过 Cell Barcode，把先前按 contig 排序的 IR 信息对齐成基于细胞的格式。注意，通常会有许多 contig 映射到同一个细胞，因为每个 IR 由 VJ 链和 VDJ 链组成。此外，已有研究表明一个细胞可以表达多个 IR [Schuldt and Binstadt, 2019].
我们将使用 Python 包 Scirpy，它在从各种格式读取 TCR 信息时会自动执行细胞对齐：
adata_tcr = ir.io.read_10x_vdj(path_tcr_csv)
print(f"Amount cells: {len(adata_tcr)}")
WARNING: Non-standard locus name ignored: Multi
Amount cells: 280045
IR 信息将以排序后的形式存储在 adata_tcr.obs 中，存为 Pandas DataFrame，其每个索引代表一个细胞。该 DataFrame 包含以下列：
has_ir：表示是否检测到 IR
multi_chain：表示是否检测到超过 2 个 IR
extra_chains：如果 multi_chain 为 true，则列出额外的链
我们在原始数据中观测到的 contig 数据也会被一并存储。它以两条链（VJ 和 VDJ）的单独列形式存储，每个细胞每条链最多 2 个 contig：
IR_V{D}J_{1,2}_locus
IR_V{D}J_{1,2}_productive
IR_V{D}J_{1,2}_{v,d,j,c}_call
IR_V{D}J_{1,2}_junction{_aa}
adata_tcr.obs.head(5)
| is_cell | high_confidence | multi_chain | extra_chains | IR_VJ_1_c_call | IR_VJ_2_c_call | IR_VDJ_1_c_call | IR_VDJ_2_c_call | IR_VJ_1_consensus_count | IR_VJ_2_consensus_count | ... | IR_VDJ_2_locus | IR_VJ_1_productive | IR_VJ_2_productive | IR_VDJ_1_productive | IR_VDJ_2_productive | IR_VJ_1_v_call | IR_VJ_2_v_call | IR_VDJ_1_v_call | IR_VDJ_2_v_call | has_ir | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| cell_id | |||||||||||||||||||||
| AAACCTGAGAAACCTA-MH9179824 | True | True | False | [] | TRAC | NaN | TRBC2 | NaN | 960.0 | NaN | ... | NaN | True | None | True | None | TRAV8-4 | NaN | TRBV5-1 | NaN | True |
| AAACCTGAGAACAATC-MH9143324 | True | True | False | [] | NaN | NaN | TRBC1 | NaN | NaN | NaN | ... | NaN | None | None | True | None | NaN | NaN | TRBV19 | NaN | True |
| AAACCTGAGAACTCGG-MH8919329 | True | True | False | [] | TRAC | NaN | TRBC2 | NaN | 5122.0 | NaN | ... | NaN | True | None | True | None | TRAV34 | NaN | TRBV5-1 | NaN | True |
| AAACCTGAGAAGCCCA-MH9143373 | True | True | False | [] | TRAC | NaN | NaN | NaN | 25914.0 | NaN | ... | NaN | True | None | None | None | TRAV16 | NaN | NaN | NaN | True |
| AAACCTGAGAAGGACA-MH9143274 | True | True | False | [] | TRAC | NaN | TRBC2 | NaN | 3212.0 | NaN | ... | NaN | True | None | True | None | TRAV8-1 | NaN | TRBV7-6 | NaN | True |
5 rows × 45 columns
请注意，这里不会自动添加患者层面的信息。让我们通过加载原始数据、在细胞级别上对齐、并按 Barcode 索引来补充这些信息。
patient_information = [
"barcode",
"Centre",
"Sample",
"patient_id",
"Collection_Day",
"Sex",
"Swab_result",
"Status",
"Smoker",
"Status_on_day_collection",
"Status_on_day_collection_summary",
"Days_from_onset",
"time_after_LPS",
"Worst_Clinical_Status",
"Outcome",
"initial_clustering",
"study_id",
"AgeRange",
"Age",
]
df_patient = df_tcr_raw[patient_information].copy()
df_patient["Days_from_onset"] = df_patient["Days_from_onset"].astype(
str
) # mixed type (str, int)
df_patient = df_patient.drop_duplicates().reset_index(drop=True)
# Assigning barcode as index
df_patient.index = df_patient.pop("barcode")
df_patient.index.name = None
df_patient.head()
| Centre | Sample | patient_id | Collection_Day | Sex | Swab_result | Status | Smoker | Status_on_day_collection | Status_on_day_collection_summary | Days_from_onset | time_after_LPS | Worst_Clinical_Status | Outcome | initial_clustering | study_id | AgeRange | Age | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| AAACCTGAGAAACCTA-MH9179824 | Newcastle | MH9179824 | COVID-014 | D0 | Female | Positive | Covid | Non-smoker | Ward_O2 | Moderate | 8 | NaN | Moderate | Home | CD8 | COVID-014 | (50, 59] | 54 |
| AAACCTGAGAACAATC-MH9143324 | Newcastle | MH9143324 | COVID-029 | D0 | Male | Positive | Covid | Non-smoker | Ward_O2 | Moderate | 16 | NaN | Moderate | Home | gdT | COVID-029 | (20, 29] | 26 |
| AAACCTGAGAACTCGG-MH8919329 | Newcastle | MH8919329 | COVID-045 | D0 | Female | Positive | Covid | Non-smoker | ITU_intubated | Critical | 5 | NaN | Death | Death | Treg | COVID-045 | (70, 79] | 77 |
| AAACCTGAGAAGCCCA-MH9143373 | Newcastle | MH9143373 | COVID-060 | D0 | Male | Negative | Non_covid | Smoker | Ward_O2 | Non_covid | 3 | NaN | Non-covid | Home | CD4 | COVID-060 | (50, 59] | 54 |
| AAACCTGAGAAGGACA-MH9143274 | Newcastle | MH9143274 | COVID-037 | D0 | Female | Positive | Covid | Non-smoker | ITU_O2 | Severe | 21 | NaN | Severe | Home | CD4 | COVID-037 | (50, 59] | 59 |
我们现在可以把这些信息添加到 adata 对象中，并检查该注释是否已作为观测项添加。
adata_tcr.obs[df_patient.columns] = df_patient
adata_tcr.obs.head()
| is_cell | high_confidence | multi_chain | extra_chains | IR_VJ_1_c_call | IR_VJ_2_c_call | IR_VDJ_1_c_call | IR_VDJ_2_c_call | IR_VJ_1_consensus_count | IR_VJ_2_consensus_count | ... | Status_on_day_collection | Status_on_day_collection_summary | Days_from_onset | time_after_LPS | Worst_Clinical_Status | Outcome | initial_clustering | study_id | AgeRange | Age | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| cell_id | |||||||||||||||||||||
| AAACCTGAGAAACCTA-MH9179824 | True | True | False | [] | TRAC | NaN | TRBC2 | NaN | 960.0 | NaN | ... | Ward_O2 | Moderate | 8 | NaN | Moderate | Home | CD8 | COVID-014 | (50, 59] | 54 |
| AAACCTGAGAACAATC-MH9143324 | True | True | False | [] | NaN | NaN | TRBC1 | NaN | NaN | NaN | ... | Ward_O2 | Moderate | 16 | NaN | Moderate | Home | gdT | COVID-029 | (20, 29] | 26 |
| AAACCTGAGAACTCGG-MH8919329 | True | True | False | [] | TRAC | NaN | TRBC2 | NaN | 5122.0 | NaN | ... | ITU_intubated | Critical | 5 | NaN | Death | Death | Treg | COVID-045 | (70, 79] | 77 |
| AAACCTGAGAAGCCCA-MH9143373 | True | True | False | [] | TRAC | NaN | NaN | NaN | 25914.0 | NaN | ... | Ward_O2 | Non_covid | 3 | NaN | Non-covid | Home | CD4 | COVID-060 | (50, 59] | 54 |
| AAACCTGAGAAGGACA-MH9143274 | True | True | False | [] | TRAC | NaN | TRBC2 | NaN | 3212.0 | NaN | ... | ITU_O2 | Severe | 21 | NaN | Severe | Home | CD4 | COVID-037 | (50, 59] | 59 |
5 rows × 63 columns
以同样的方式，我们可以加载 BCR 数据并检查 IR 注释。
adata_bcr = ir.io.read_10x_vdj(path_bcr_input)
adata_bcr.obs.head(5)
WARNING: Non-standard locus name ignored: Multi
WARNING: Non-standard locus name ignored: None
| is_cell | high_confidence | multi_chain | extra_chains | IR_VJ_1_c_call | IR_VJ_2_c_call | IR_VDJ_1_c_call | IR_VDJ_2_c_call | IR_VJ_1_consensus_count | IR_VJ_2_consensus_count | ... | IR_VDJ_2_locus | IR_VJ_1_productive | IR_VJ_2_productive | IR_VDJ_1_productive | IR_VDJ_2_productive | IR_VJ_1_v_call | IR_VJ_2_v_call | IR_VDJ_1_v_call | IR_VDJ_2_v_call | has_ir | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| cell_id | |||||||||||||||||||||
| AAACCTGAGAACTCGG-MH9143420 | True | True | False | [] | IGKC | NaN | IGHM | NaN | 6982.0 | NaN | ... | NaN | True | None | True | None | IGKV3-11 | NaN | IGHV3-23 | NaN | True |
| AAACCTGAGAATAGGG-MH8919227 | True | True | False | [] | IGKC | NaN | IGHA1 | NaN | 626.0 | NaN | ... | NaN | True | None | True | None | IGKV3-15 | NaN | IGHV1-24 | NaN | True |
| AAACCTGAGAATGTGT-MH8919326 | True | True | False | [{"c_call": "IGLC1", "consensus_count": 748, "... | IGKC | NaN | IGHM | NaN | 1028.0 | NaN | ... | NaN | True | None | True | None | IGKV1-33 | NaN | IGHV1-2 | NaN | True |
| AAACCTGAGAATGTGT-MH9143424 | True | True | False | [] | IGLC2 | NaN | NaN | NaN | 736.0 | NaN | ... | NaN | True | None | None | None | IGLV3-21 | NaN | NaN | NaN | True |
| AAACCTGAGACACTAA-MH9143423 | True | True | False | [] | IGKC | NaN | IGHG1 | NaN | 66.0 | NaN | ... | NaN | True | None | True | None | IGKV3-20 | NaN | IGHV3-23 | NaN | True |
5 rows × 45 columns
与之前一样，我们可以提取为每个 contig 注释的额外信息。
patient_information = ["barcode", "patient_id"]
df_patient = df_bcr_raw[patient_information].copy()
df_patient = df_patient.drop_duplicates().reset_index(drop=True)
# Assigning barcode as index
df_patient.index = df_patient.pop("barcode")
df_patient.index.name = None
df_patient.head()
| patient_id | |
|---|---|
| AAACCTGAGCGATATA-MH9143270 | COVID-002 |
| AAACCTGTCTACCAGA-MH9143270 | COVID-002 |
| AAACGGGCATTTGCTT-MH9143270 | COVID-002 |
| AAACGGGGTAGCACGA-MH9143270 | COVID-002 |
| AAACGGGTCGTACCGG-MH9143270 | COVID-002 |
并把它们加入到细胞注释中：
adata_bcr.obs[df_patient.columns] = df_patient
adata_bcr.obs.head()
| is_cell | high_confidence | multi_chain | extra_chains | IR_VJ_1_c_call | IR_VJ_2_c_call | IR_VDJ_1_c_call | IR_VDJ_2_c_call | IR_VJ_1_consensus_count | IR_VJ_2_consensus_count | ... | IR_VJ_1_productive | IR_VJ_2_productive | IR_VDJ_1_productive | IR_VDJ_2_productive | IR_VJ_1_v_call | IR_VJ_2_v_call | IR_VDJ_1_v_call | IR_VDJ_2_v_call | has_ir | patient_id | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| cell_id | |||||||||||||||||||||
| AAACCTGAGAACTCGG-MH9143420 | True | True | False | [] | IGKC | NaN | IGHM | NaN | 6982.0 | NaN | ... | True | None | True | None | IGKV3-11 | NaN | IGHV3-23 | NaN | True | COVID-011 |
| AAACCTGAGAATAGGG-MH8919227 | True | True | False | [] | IGKC | NaN | IGHA1 | NaN | 626.0 | NaN | ... | True | None | True | None | IGKV3-15 | NaN | IGHV1-24 | NaN | True | COVID-075 |
| AAACCTGAGAATGTGT-MH8919326 | True | True | False | [{"c_call": "IGLC1", "consensus_count": 748, "... | IGKC | NaN | IGHM | NaN | 1028.0 | NaN | ... | True | None | True | None | IGKV1-33 | NaN | IGHV1-2 | NaN | True | COVID-003 |
| AAACCTGAGAATGTGT-MH9143424 | True | True | False | [] | IGLC2 | NaN | NaN | NaN | 736.0 | NaN | ... | True | None | None | None | IGLV3-21 | NaN | NaN | NaN | True | COVID-039 |
| AAACCTGAGACACTAA-MH9143423 | True | True | False | [] | IGKC | NaN | IGHG1 | NaN | 66.0 | NaN | ... | True | None | True | None | IGKV3-20 | NaN | IGHV3-23 | NaN | True | COVID-036 |
5 rows × 46 columns
40.9. 质量控制#
为了分析，我们依赖高质量的输入数据。因此，识别出 AIR 信息有误或不完整的细胞至关重要：
不完整的 AIR：某个细胞只被分配到 VJ 链或 VDJ 链中的一条，因为另一条链在测序过程中丢失了。这些细胞虽然仍然有效，但在需要完整 AIR 序列信息时，我们无法将其用于下游分析。
多个 AIR：一个细胞也可能被分配到多个 AIR。虽然已观察到 T/B 细胞可以表达双 AIR [Schuldt and Binstadt, 2019]，但带有两条以上 IR 的细胞往往提示是双细胞，不应用于下游分析。
通过汇总各 AIR 链的数据，我们可以很容易地为每个细胞标定这种 AIR 状态。这样，之后我们就能过滤掉那些 AIR 无法用于分析的细胞。具体需要哪些 AIR 信息，在很大程度上取决于你的下游分析。
40.9.1. 标定 AIR 状态#
ir.tl.chain_qc(adata_tcr)
ir.tl.chain_qc(adata_bcr)
下面给出生成的列，并说明其可能的取值：
chain_pairing
orphan {VJ}/{VDJ}（孤链）：只有 VJ 或 VDJ 链其中之一的细胞
single pair（单对）：拥有完整配对（VJ 与 VDJ）AIR 的细胞
extra {VJ}/{VDJ}（额外链）：在完整配对 AIR 之外、还多出一条 VJ 或 VDJ 链的细胞
two full chains（两对完整 AIR）：带有两对配对 AIR 的细胞
multichain（多链）：具有两条以上 VJ 或 VDJ 链的细胞（很可能是双细胞）
receptor_type：BCR、TCR、无 IR，或模糊（ambiguous，TCR 与 BCR 的混合链）
receptor_subtype：α/β、ɣ/δ、IG-λ 和 IG-κ 等链的组合配置
40.9.2. 可视化#
根据数据的不同，链配对的质量可能有差异。你应当跨不同条件（如样本、数据来源等）可视化结果，因为这些条件可能给数据集带来质量差异。下面我们按不同的数据采集中心，可视化所收集数据的质量。
_ = ir.pl.group_abundance(adata_tcr, groupby="Centre", target_col="chain_pairing")
三个中心的数据主要都是带单对 TCR 的细胞，这是数据质量良好的标志。通常，当我们查看每种条件下细胞的占比、而非细胞的绝对数量时，能更好地比较不同条件：
_ = ir.pl.group_abundance(
adata_tcr, groupby="Centre", target_col="chain_pairing", normalize=True
)
可以看到，所有中心提供的数据总体质量相当。大多数 T 细胞表达单个 TCR，其次是孤立（orphan）的 VDJ 链和 VJ 链，这在 TCR 测序中很典型。质量良好的样本，其典型值大致为：单 AIR 大于 60%、孤链 10–20%、额外链约 10%。注意，在联合测序多种模态时，没有 AIR 的细胞所占比例自然取决于研究中非 AIR 细胞的数量。
_ = ir.pl.group_abundance(
adata_tcr, groupby="Centre", target_col="receptor_type", normalize=True
)
_ = ir.pl.group_abundance(
adata_tcr, groupby="Centre", target_col="receptor_subtype", normalize=True
)
可以看到，该数据集只包含 TCR（第一张图），这并不意外，因为我们已经把 T 细胞和 B 细胞的 DataFrame 分开了。此外，该数据集只包含 α/β-T 细胞。由于图中无法看出没有 IR 的细胞数量，我们查看绝对数值，发现只有 22 个细胞没有 AIR 信息。在联合测序多种模态时，这一比例在大多数数据集中会更大，具体取决于 AIR 细胞的数量。
adata_tcr.obs["chain_pairing"].value_counts()
single pair 196957
orphan VDJ 45266
extra VJ 19034
orphan VJ 7937
extra VDJ 7473
two full chains 3356
no IR 22
Name: chain_pairing, dtype: int64
对于 BCR，我们将从患者层面来看。为便于可视化，我们下采样到两名患者：
adata_bcr_tmp = adata_bcr[
adata_bcr.obs["patient_id"].isin(["COVID-003", "AP11"])
].copy()
_ = ir.pl.group_abundance(
adata_bcr_tmp, groupby="patient_id", target_col="chain_pairing"
)
_ = ir.pl.group_abundance(
adata_bcr_tmp, groupby="patient_id", target_col="chain_pairing", normalize=True
)
虽然我们看到两名患者之间 B 细胞的数量有相当大的差异，但链配对的质量是相当的。
_ = ir.pl.group_abundance(
adata_bcr_tmp, groupby="patient_id", target_col="receptor_type", normalize=True
)
_ = ir.pl.group_abundance(
adata_bcr_tmp, groupby="patient_id", target_col="receptor_subtype", normalize=True
)
同样，数据集中只有一种类型的 AIR（BCR），但这次轻链是 κ 链和 λ 链的混合。
40.9.3. 过滤#
你可能想根据不同的 AIR 状态来过滤细胞。在决定如何过滤时，要根据你的研究在“可用数据量”与“数据质量”之间权衡。通常，所有细胞都会保留在数据集中，并标记其 AIR 状态；对于不同的下游分析，再只使用信息正确的细胞。例如，如果基于 CDR3-VDJ 链做数据库查询（见第 X 章），可以使用孤立的 VDJ 链，而在搜索完整 AIR 时则不行。下面我们将演示如何进行过滤。至于是在此处过滤，还是在各分析环节之前过滤，由读者自行决定。
过滤可以采用不同程度的严格度和不同的组合，例如:
拥有 IR：只保留表达了 IR 的细胞，因为我们无法对其余细胞做序列层面的分析
无多链：过滤掉带有超过 2 个 AIR 的细胞，以排除双细胞
只有一个 AIR：过滤掉带有两个 AIR 或额外链的细胞，因为此时无法把特异性归因于某一条特定序列
无单链：过滤掉只有 VJ 或 VDJ 链的细胞，由于缺少信息，无法用于所有下游分析
由于本书后续会使用完整数据，我们会针对上述每个步骤，把数据过滤到一个临时的 AnnData 对象中。从绝对细胞数可以看出每一步会损失多少细胞。
adata_bcr_tmp.obs["chain_pairing"].value_counts()
single pair 1489
orphan VJ 312
extra VJ 46
ambiguous 31
orphan VDJ 25
two full chains 6
extra VDJ 4
no IR 1
Name: chain_pairing, dtype: int64
print(f"Amount of all B cells:\t\t\t\t{len(adata_bcr)}")
adata_bcr_tmp = adata_bcr[adata_bcr.obs["chain_pairing"] != "no IR"]
print(f"Amount of B cells with AIR:\t\t\t{len(adata_bcr_tmp)}")
adata_bcr_tmp = adata_bcr_tmp[adata_bcr_tmp.obs["chain_pairing"] != "multi_chain"]
print(f"Amount of B cells without doublets:\t\t{len(adata_bcr_tmp)}")
adata_bcr_tmp = adata_bcr_tmp[
~adata_bcr_tmp.obs["chain_pairing"].isin(
["two full chains", "extra VJ", "extra VDJ"]
)
]
print(f"Amount of B cells with unique AIR per cell:\t{len(adata_bcr_tmp)}")
adata_bcr_tmp = adata_bcr_tmp[adata_bcr_tmp.obs["chain_pairing"] == "single pair"]
print(f"Amount of B cells with sinlge complete AIR:\t{len(adata_bcr_tmp)}")
Amount of all B cells: 159446
Amount of B cells with AIR: 159185
Amount of B cells without doublets: 159185
Amount of B cells with unique AIR per cell: 153936
Amount of B cells with single complete AIR: 108395
最终，我们得到一个只含单对 BCR 的 AnnData 对象。不过，这样数据量可能会大幅减少，从而限制下游分析。
adata_bcr_tmp.obs["chain_pairing"].value_counts()
single pair 108395
Name: chain_pairing, dtype: int64
同样，我们可以为TCR数据进行这种过滤。
adata_tcr.obs["chain_pairing"].value_counts()
single pair 196957
orphan VDJ 45266
extra VJ 19034
orphan VJ 7937
extra VDJ 7473
two full chains 3356
no IR 22
Name: chain_pairing, dtype: int64
print(f"Amount of all T cells:\t\t\t\t{len(adata_tcr)}")
adata_tcr_tmp = adata_tcr[adata_tcr.obs["chain_pairing"] != "no IR"]
print(f"Amount of T cells with AIR:\t\t\t{len(adata_tcr_tmp)}")
adata_tcr_tmp = adata_tcr_tmp[adata_tcr_tmp.obs["chain_pairing"] != "multi_chain"]
print(f"Amount of T cells without doublets:\t\t{len(adata_tcr_tmp)}")
adata_tcr_tmp = adata_tcr_tmp[
~adata_tcr_tmp.obs["chain_pairing"].isin(
["two full chains", "extra VJ", "extra VDJ"]
)
]
print(f"Amount of T cells with unique AIR per cell:\t{len(adata_tcr_tmp)}")
adata_tcr_tmp = adata_tcr_tmp[adata_tcr_tmp.obs["chain_pairing"] == "single pair"]
print(f"Amount of T cells with sinlge complete AIR:\t{len(adata_tcr_tmp)}")
Amount of all T cells: 280045
Amount of T cells with AIR: 280023
Amount of T cells without doublets: 280023
Amount of T cells with unique AIR per cell: 250160
Amount of T cells with sinlge complete AIR: 196957
adata_tcr_tmp.obs["chain_pairing"].value_counts()
single pair 196957
Name: chain_pairing, dtype: int64
最后，我们保存好已注释的数据，供后续分析使用。
sc.write(adata=adata_tcr, filename=path_tcr_out)
sc.write(adata=adata_bcr, filename=path_bcr_out)
40.10. Quiz#
40.11. 参考文献#
Nicholas Borcherding, Nicholas L Bormann, and Gloria Kraus. Screpertoire: an r-based toolkit for single-cell immune receptor analysis. F1000Research, 2020.
Bryan Briney, Anne Inderbitzin, Collin Joyce, and Dennis R Burton. Commonality despite exceptional diversity in the baseline human antibody repertoire. Nature, 566(7744):393–397, 2019.
Stefan Canzar, Karlynn E Neu, Qingming Tang, Patrick C Wilson, and Aly A Khan. BASIC: BCR assembly from single cells. Bioinformatics, 33(3):425–427, 10 2016. URL: https://doi.org/10.1093/bioinformatics/btw631, arXiv:https://academic.oup.com/bioinformatics/article-pdf/33/3/425/25145999/btw631.pdf, doi:10.1093/bioinformatics/btw631.
Max D Cooper and Matthew N Alder. The evolution of adaptive immune systems. Cell, 124(4):815–822, 2006.
Ping Hu, Wenhua Zhang, Hongbo Xin, and Glenn Deng. Single cell isolation and analysis. Frontiers in cell and developmental biology, 4:116, 2016.
Ida Lindeman, Guy Emerton, Lira Mamanova, Omri Snir, Krzysztof Polanski, Shuo-Wang Qiao, Ludvig M. Sollid, Sarah A. Teichmann, and Michael J. T. Stubbington. Bracer: b-cell-receptor reconstruction and clonality inference from single-cell rna-seq. Nature Methods, 15(8):563–565, Aug 2018. URL: https://doi.org/10.1038/s41592-018-0082-3, doi:10.1038/s41592-018-0082-3.
Nathaniel J Schuldt and Bryce A Binstadt. Dual tcr t cells: identity crisis or multitaskers? The Journal of Immunology, 202(3):637–644, 2019.
Mandeep Singh, Ghamdan Al-Eryani, Shaun Carswell, James M Ferguson, James Blackburn, Kirston Barton, Daniel Roden, Fabio Luciani, Tri Giang Phan, Simon Junankar, and others. High-throughput targeted long-read single cell sequencing reveals the clonal and transcriptional landscape of lymphocytes. Nature communications, 10(1):1–13, 2019.
Emily Stephenson, Gary Reynolds, Rachel A Botting, Fernando J Calero-Nieto, Michael D Morgan, Zewen Kelvin Tuong, Karsten Bach, Waradon Sungnak, Kaylee B Worlock, Masahiro Yoshida, and others. Single-cell multi-omics analysis of the immune response in covid-19. Nature medicine, 27(5):904–916, 2021.
Gregor Sturm, Tamas Szabo, Georgios Fotakis, Marlene Haider, Dietmar Rieder, Zlatko Trajanoski, and Francesca Finotello. Scirpy: a scanpy extension for analyzing single-cell t-cell receptor-sequencing data. Bioinformatics, 36(18):4817–4818, 2020.
Amit A. Upadhyay, Robert C. Kauffman, Amber N. Wolabaugh, Alice Cho, Nirav B. Patel, Samantha M. Reiss, Colin Havenar-Daughton, Reem A. Dawoud, Gregory K. Tharp, Iñaki Sanz, Bali Pulendran, Shane Crotty, F. Eun-Hyung Lee, Jens Wrammert, and Steven E. Bosinger. Baldr: a computational pipeline for paired heavy and light chain immunoglobulin reconstruction in single-cell rna-seq data. Genome Medicine, 10(1):20, Mar 2018. URL: https://doi.org/10.1186/s13073-018-0528-3, doi:10.1186/s13073-018-0528-3.
Sebastiaan Valkiers, Nicky de Vrij, Sofie Gielis, Sara Verbandt, Benson Ogunjimi, Kris Laukens, and Pieter Meysman. Recent advances in t-cell receptor repertoire analysis: bridging the gap with multimodal single-cell rna sequencing. ImmunoInformatics, pages 100009, 2022.
F Alexander Wolf, Philipp Angerer, and Fabian J Theis. Scanpy: large-scale single-cell gene expression data analysis. Genome biology, 19(1):1–5, 2018.
Alexander Yermanos, Andreas Agrafiotis, Raphael Kuhn, Damiano Robbiani, Josephine Yates, Chrysa Papadopoulou, Jiami Han, Ioana Sandu, Cédric Weber, Florian Bieberich, and others. Platypus: an open-access software for integrating lymphocyte single-cell immune repertoires with transcriptomes. NAR genomics and bioinformatics, 3(2):lqab023, 2021.
Veronika Zarnitsyna, Brian Evavold, Louie Schoettle, Joseph Blattman, and Rustom Antia. Estimating the diversity, completeness, and cross-reactivity of the t cell repertoire. Frontiers in immunology, 4:485, 2013.
ImmunoMind Team. immunarch: An R Package for Painless Bioinformatics Analysis of T-Cell and B-Cell Immune Repertoires. August 2019. URL: https://doi.org/10.5281/zenodo.3367200, doi:10.5281/zenodo.3367200.