免疫组库分析软件MiXCR的安装和使用
安装
docker镜像安装
docker pull milaboratory/mixcr
镜像使用方法
docker run -it --rm -m 4g -v /path/to/data:/work milaboratory/mixcr:latest mixcr align test_R1.fastq.gz test_R2.fastq.gz test.vdjca
普通安装
在Mac OS X上使用Homebrew或者在Linux系统上使用linuxbrew安装
brew install milaboratory/all/mixcr
将已经安装的MiXCR更新到最新版本
brew update
brew upgrade mixcr
手动安装
- 在release page下载最新版稳定编译的MiXCR
- 解开压缩包
- 将解压缩生成的文件夹添加到环境变量中
- 添加mixcr脚本的符号链接到bin文件夹
- 或者在执行脚本中指定MiXCR的全路径来直接使用
快速使用
综述
典型的MiXCR工作流程主要由三个部分构成：
- align：将测序结果比对到T细胞或B细胞受体的V、D、J、C基因参考序列上
- assemble：利用前一步骤获得的比对结果拼接clonotypes（为了提取特定的基因区域信息，比如CDR3）
- export：输出比对结果（exportAlignments模块）或者clones信息（exportClones模块），生成可读文件
MiXCR的assemble模块有几种不同的拼接方法可以选择：
- assembleContigs：拼接完整的TCR或者IG受体clonotype序列
- 对于RNA-Seq or non-targeted DNA data, 工作流程可能包括以下两部分：
- assemblePartial：将有重叠区域的序列片段拼接成相对较长的包含CDR3区域的contigs
- extend: 估算测序和比对质量较好但长度较短的TCR比对序列的germline序列
为了简化输入命令，MiXCR提供了analyze命令模块，打包了整个分析流程
MiXCR支持一下若干种数据类型：fasta，fastq，fastq.gz，paired-end fastq和fastq.gz。作为每一步骤的输出结果，MiXCR生成包含各种信息的二进制压缩文件（比对生成alignments，拼接生成clones）。利用exportAlignments和exportClones命令模块，每一个二进制文件都可以转化成tab分割的可读文本文件。
实例
默认流程 / multiplex-PCR
利用analyze amplicon命令分析multiplex-PCR扩增的TCR/BCR基因DNA片段
mixcr analyze amplicon --species hs --starting-material dna --5-end v-primers --3-end j-primers --adapters adapters-present --receptor-type IGH input_R1.fastq input_R2.fastq analysis
只有一个参数修改为非默认值（--receptor-type IGH），这个参数的改变可以让MiXCR 调用针对B细胞优化的比对模块并且只输出IGH序列。其实这个参数是可以缺省的，缺省状态下MiXCR会调用默认的比对模块并输出样本中所有的TCR/BCR序列。
生成的文件（analysis.clonotypes.IGH.txt）是一个tab分隔的表格，包含CDR3序列拼接的所有clonotypes（克隆丰度，CDR3序列， VDJ基因等）。
详细流程
利用analyze amplicon模块与执行下面的命令是等价的
> mixcr align -s hs -p kAligner2 input_R1.fastq input_R2.fastq alignments.vdjca
... Building alignments
> mixcr assemble alignments.vdjca clones.clns
... Assembling clones
> mixcr exportClones --chains IGH clones.clns clones.txt
... Exporting clones to tab-delimited file
基于5’RACE扩增实验的数据分析
考虑基于5’RACE（一个read覆盖CDR3区域和临近序列，另一个read覆盖V基因的5’UTR和下游序列）实验准备的IGH基因cDNA文库双端测序的数据处理流程，全部分析流程可以通过analyze amplicon命令实现
> mixcr analyze amplicon --species hs --starting-material rna --5-end v-primers --3-end j-primers --adapters adapters-present input_R1.fastq input_R2.fastq analysis
结果文件（analysis.clonotypes.<chains>.txt）包含详细的clonotypes信息。
详细流程
利用analyze amplicon模块将执行下面的流程：
- Align 把原始序列比对到IGH基因的VDJ基因序列片段上
> mixcr align -s hs -OvParameters.geneFeatureToAlign=VTranscript --report analysis.report input_R1.fastq input_R2.fastq analysis.vdjca
用来比对V基因的非默认基因特征（-OvParameters.geneFeatureToAlign=VTranscript）同时利用了两个reads的信息，为了让MiXCR利用CDR3反向read比对V基因的5’UTRS和部分5’端编码区域。MiXCR还会生成report文件（通过可选参数--report指定），其中包含的具体运行统计信息如下
Analysis Date: Mon Aug 25 15:22:39 MSK 2014
Input file(s): input_r1.fastq,input_r2.fastq
Output file: alignments.vdjca
Command line arguments: align --report alignmentReport.log input_r1.fastq input_r2.fastq alignments.vdjca
Total sequencing reads: 323248
Successfully aligned reads: 210360
Successfully aligned, percent: 65.08%
Alignment failed because of absence of V hits: 4.26%
Alignment failed because of absence of J hits: 30.19%
Alignment failed because of low total score: 0.48%
可以利用exportAlignments命令将比对生成的二进制结果（analysis.vdjca）转化为可读的文本文件。
- Assemble 拼接clonotypes
> mixcr assemble --report analysis.report analysis.vdjca -a analysis.clna
这一步骤会校正PCR和测序错误并建立clonotypes，默认情况下clonotypes会拼接CDR3序列；可以通过设置assemble模块的参数来制定其他的基因区域（参考assemble documentation），可选的report文件analysis.report包含各种调试信息
Analysis Date: Mon Aug 25 15:29:51 MSK 2014
Input file(s): alignments.vdjca
Output file: clones.clns
Command line arguments: assemble --report assembleReport.log alignments.vdjca clones.clns
Final clonotype count: 11195
Total reads used in clonotypes: 171029
Reads used, percent of total: 52.89%
Reads used as core, percent of used: 92.04%
Mapped low quality reads, percent of used: 7.96%
Reads clustered in PCR error correction, percent of used: 0.04%
Clonotypes eliminated by PCR error correction: 72
Percent of reads dropped due to the lack of clonal sequence: 2.34%
Percent of reads dropped due to low quality: 3.96%
Percent of reads dropped due to failed mapping: 5.87%
- Export 将包含clones列表的二进制文件（analysis.clna）导出为可读的文本文件
> mixcr exportClones --chains TRA analysis.clna analysis.clonotypes.TRA.txt
> mixcr exportClones --chains TRB analysis.clna analysis.clonotypes.TRB.txt
> ...
导出的clones信息如下表所示
| Clone count | Clone fraction | … | V hits | J hits | seq. CDR3 | AA. seq. CDR3 | … |
|---|---|---|---|---|---|---|---|
| 4369 | 2.9E-3 | … | IGHV4-39*00(1388) | IGHJ6 *00(131) | TGTGTGAG… | CVRHKPM… | … |
| 3477 | 2.5E-3 | … | IGHV4-34*00(1944) | IGHJ4 *00(153) | TGTGCGAT… | CAIWDVGL… | … |
| … | … | … | … | … | … | … | … |
导出的各种选项详见export文档，上述的所有步骤都可以根据特定研究的分析流程进行个性化设置。
高质量全长IG免疫组库分析
对于基于cDNA全长的IG免疫组库分析，我们一般推荐UMI标签制备文库并使用非对称双端测序350 bp + 100 bp Illumina MiSeq测序方法（详情参考Nature Protocols paper）。这种方法可以获得长片段高质量测序结果，而且可以利用MiGEC software有效去除PCR和测序错误。获得的高质量数据可以进一步利用MiXCR处理，以提取全长IGH或IGL组库。
全部分析流程可以使用analyze amplicon命令
> mixcr analyze amplicon --species hs --starting-material rna --5-end v-primers --3-end j-primers --adapters adapters-present --receptor-type BCR --region-of-interest VDJRegion --only-productive --align "-OreadsLayout=Collinear" --assemble "-OseparateByC=true" --assemble "-OqualityAggregationType=Average" --assemble "-OclusteringFilter.specificMutationProbability=1E-5" --assemble "-OmaxBadPointsPercent=0" input_R1.fastq input_R2.fastq analysis
这一步骤会生成以下结果文件（analysis.clonotypes.IGH.txt，analysis.clonotypes.IGK.txt和analysis.clonotypes.IGL.txt），其中包括详细的clonotypes信息。这里我们要强调几个可选参数：
- --receptor-type BCR 需要MiXCR调用B细胞优化的比对模块（等同于对align模块使用-p kAligner2参数）并且只输出IG序列。
- region-of-interest VDJRegion 对assemble模块使用-OassemblingFeatures=VDJRegion参数
- --only-production 在export输出的clonotypes中过滤掉out-of-frame和stop codon
- --align <option> 在align过程中设置其他的参数
- --assemble <option> 在assemble过程中设置其他参数
详细流程
使用analyze amplicon命令等同于执行下面的MiXCR步骤
> mixcr align -p kaligner2 -s hs -r alignmentReport.txt -OreadsLayout=Collinear -OvParameters.geneFeatureToAlign=VTranscript read_R1.fastq.gz read_R2.fastq.gz alignments.vdjca
选项-s用来指定物种（e.g. homo sapiens - hsa, mus musculus - mmu），参数-OreadsLayout 用来设定reads方向（Collinear, Opposite, Unknown）。这里需要注意的是，经过MiGEC分析的双端reads方向是Collinear。
除了KAligner2，也可以使用默认的MiXCR比对模块，只是也许会忽略一些亚变异类型，这些变异类型是由V基因片段的若干核苷酸插入形成的。
- Assemble拼接clones
> mixcr assemble -r assembleReport.txt -OassemblingFeatures=VDJRegion -OseparateByC=true -OqualityAggregationType=Average -OclusteringFilter.specificMutationProbability=1E-5 -OmaxBadPointsPercent=0 alignments.vdjca clones.clns
-OseparateByC=true 把clones按照不同的抗体亚型分类
-OcloneClusteringParameters=null 关闭基于频率的PCR错误校正
根据数据质量，可以通过设置-ObadQualityThreshold参数来调节输入数据的阈值来优化clonotypes的提取。
- Export输出clones结果
> mixcr exportClones -c IGH -o -t clones.clns clones.txt
选项-o和-t用于过滤包含out-of-frame和stop codon的clonotypes，-c指定哪条链的数据应该被提取（e.g. IGH, IGL）。
RNA-Seq数据分析
MiXCR可以用于提取RNA-Seq数据中TCR和BCR的CDR3组库，提取效率取决于样本红T/B细胞的丰度和测序长度。推荐2x150bp或者2x100bp的双端测序方法。不过在双端2x50bp的RNA-Seq数据（比如肿瘤样本中），主要clonotypes信息也可以被提取。
单一analyze shotgun命令可以完成分析
> mixcr analyze shotgun --species hs --starting-material rna --only-productive input_R1.fastq input_R2.fastq analysis
生成的结果文件（analysis.clonotypes.TRA.txt, analysis.clonotypes.IGH.txt等）包含clonotypes的详细信息。
详细流程
- Align 比对reads
> mixcr align -s hs -p rna-seq -OallowPartialAlignments=true data_R1.fastq.gz data_R2.fastq.gz alignments.vdjca
所有mixcr align的参数都可以在这里使用（比如-s来指定物种）：
-OallowPartialAlignments=true选项保留部分比对结果用于后续的assemblePartial模块
- Assemble parial reads拼接部分reads
> mixcr assemblePartial alignments.vdjca alignmentsRescued.vdjca
为了获得包含CDR3全长序列的拼接reads，建议使用迭代mixcr的assemblePartial模块多次迭代拼接结果。多次迭代需要-p参数，根据我们的经验，两次迭代后结果最优
> mixcr assemblePartial alignments.vdjca alignmentsRescued_1.vdjca
> mixcr assemblePartial alignmentsRescued_1.vdjca alignmentsRescued_2.vdjca
- 利用已有V和J基因延长TCR比对结果，基于germline序列不全覆盖度不完全的CDR3s
> mixcr extendAlignments alignmentsRescued_2.vdjca alignmentsRescued_2_extended.vdjca
- Assemble拼接clones
> mixcr assemble alignmentsRescued_2_extended.vdjca clones.clns
所有mixcr assemble的参数都可以在这里使用：
- 对于低质量数据，建议降低输入质量阈值（e.g. -ObadQualityThreshold=15）
- 为了克隆丰度与错误校正算法相结合，增加下面的选项：-OaddReadsCountOnClustering=true
- Exporting导出clones
> mixcr exportClones -c TRA -o -t clones.clns clones.txt
可以指定导出感兴趣的免疫受体链（-c TRA 或者 -c TRB等），也可以去除包含out-of-frame（选项-o）和stop codon的突变体（选项-t）。
参数解读
| 模块名称 | 模块功能 |
|---|---|
| analyze | 对指定输入文件执行MiXCR整套分析流程 |
| align | 对输入测序reads生成V/D/J/C基因比对序列 |
| assemble | 拼接clones |
| assembleContigs | 拼接全长序列 |
| assemblePartial | 拼接部分比对reads 生成更长的序列 |
| extend | 用germline序列预测比对序列或clones |
| exportAlignments | 将V/D/J/C比对结果导出为tab分隔文件 |
| exportAlignmentsPretty | 导出比对结果的详细信息 |
| exportClones | 将拼接的clones导出为tab分隔文件 |
| exportClonesPretty | 导出clones的详细信息 |
| exportReadsForClones | 从clones&比对结果（*.clna）中导出特定clone的reads，如果没有指定clone，所有对应的reads都会被导出 |
| exportAlignmentsForClones | 从clones&比对结果（*.clna）中导出特定clone的比对结果 |
| exportReads | 从vdjca文件导出原始reads |
| mergeAlignments | 将若干*.vdjca文件合并为一个比对文件 |
| filterAlignments | 过滤比对结果 |
| sortAlignments | 根据read ID排序vdjca文件中的比对结果 |
| alignmentsDiff | 计算两个vdjca文件的差异 |
| clonesDiff | 计算两个clns文件的差异 |
| slice | 分割clna文件 |
结果解读
输出文件包含内容
| 输出表头 | 注释内容 |
|---|---|
| cloneId | clone识别号码 |
| cloneCount | clone数量 |
| cloneFraction | clone比例 |
| targetSequences | 目标序列 |
| targetQualities | 目标质量 |
| allVHitsWithScore | 所有V基因命中和评分 |
| allDHitsWithScore | 所有D基因命中和评分 |
| allJHitsWithScore | 所有J基因命中和评分 |
| allCHitsWithScore | 所有C基因命中和评分 |
| allVAlignments | 所有V基因比对结果 |
| allDAlignments | 所有D基因比对结果 |
| allJAlignments | 所有J基因比对结果 |
| allCAlignments | 所有C基因比对结果 |
| nSeqFR1 | FR1核苷酸序列 |
| minQualFR1 | FR1最小质量 |
| nSeqCDR1 | CDR1核苷酸序列 |
| minQualCDR1 | CDR1最小质量 |
| nSeqFR2 | FR2核苷酸序列 |
| minQualFR2 | FR2最小质量 |
| nSeqCDR2 | CDR2核苷酸序列 |
| minQualCDR2 | CDR2最小质量 |
| nSeqFR3 | FR3核苷酸序列 |
| minQualFR3 | FR3最小质量 |
| nSeqCDR3 | CDR3核苷酸序列 |
| minQualCDR3 | CDR3最小质量 |
| nSeqFR4 | FR4核苷酸序列 |
| minQualFR4 | FR4最小质量 |
| aaSeqFR1 | FR1氨基酸序列 |
| aaSeqCDR1 | CDR1氨基酸序列 |
| aaSeqFR2 | FR2氨基酸序列 |
| aaSeqCDR2 | CDR2氨基酸序列 |
| aaSeqFR3 | FR3氨基酸序列 |
| aaSeqCDR3 | CDR3氨基酸序列 |
| aaSeqFR4 | FR4氨基酸序列 |
| refPoints | 参考点 |
