TCGA计划的4个找somatic mutation的软件使用体验

体细胞突变（somatic mutation）是指患者某些组织或者器官后天性地发生了体细胞变异，虽然它不会遗传给后代个体，却可以通过细胞分裂，遗传给子代细胞。体细胞突变对肿瘤的发生发展有关键性的作用，并且它也是制定肿瘤癌症靶向治疗措施的关键所在。

NGS使体细胞变异的检测更加全面，成本更低，在检测多种体细胞变异上具有很大的优势，但在使用过程中还存在着挑战：如样品降解、覆盖度不足、遗传异质性和组织污染（杂质)等问题。 为应对以上挑战，降低错误率，科学家采取了不同的算法和统计模型用于检测体细胞突变。目前最受欢迎的有*Varscan、SomaticSniper、 Strelka 和MuTect2 *。

这些软件大都是直接对肿瘤-正常样本的每个位点进行比较，对肿瘤样本中明显高于正常样本的次等位基因进行标记，作为体细胞变异，同时排除种系突变和杂合性丢失（LOH）情况。虽然这些软件具有较高的引用率，并在不断地更新，但仍存在不足：

a 缺乏完整可靠的实验来评估检测结果；

b 缺乏金标准，不能保证检测到的灵敏度和特异性最高；

c 在实际应用中，各软件的相对优缺点在很大程度上是未知的。

下面是TCGA计划采取的软件：

  • MuSE

  • varscan

  • MuTect

  • SomaticSniper

当然，能找somatic mutation的软件还有很多，比如Strelka等，就不一一讲解啦。其实最基础的原理都是应该是 除去在normal样本里面出现过的germline变异位点，可以很简单的GATK UnifiedGenotyper followed by simple subtraction即可。

首先用mutect2

这个软件已经被整合到GATK里面啦，所以下载GATK即可使用它。

java软件，下载即可使用，GATK软件下载以前需要自行注册，目前是开放下载了，使用代码是：

  1. GENOME=/home/jianmingzeng/reference/genome/human_g1k_v37/human_g1k_v37.fasta

  2. GATK=/home/jianmingzeng/biosoft/GATK/GenomeAnalysisTK.jar

  3. DBSNP=/home/jianmingzeng/annotation/variation/human/dbSNP/All_20160601.vcf.gz

  4. normal_bam=NPC10F-N_recal.bam

  5. tumor_bam=NPC10F-T_recal.bam

  6. sample=NPC10F

  7. java -Djava.io.tmpdir=$TMPDIR   -Xmx40g -jar $GATK  -T MuTect2 \

  8. -R $GENOME -I:tumor $tumor_bam-I:normal $normal_bam \

  9. --dbsnp  $DBSNP   -o ${sample}-mutect2.vcf

出来的结果里面有些比较陌生的tags，需要仔细理解，这样才能看懂vcf文件，并进行进一步的过滤。

  1. ##FILTER=<ID=alt_allele_in_normal,Description="Evidence seen in the normal sample">

  2. ##FILTER=<ID=clustered_events,Description="Clustered events observed in the tumor">

  3. ##FILTER=<ID=germline_risk,Description="Evidence indicates this site is germline, not somatic">

  4. ##FILTER=<ID=homologous_mapping_event,Description="More than three events were observed in the tumor">

  5. ##FILTER=<ID=multi_event_alt_allele_in_normal,Description="Multiple events observed in tumor and normal">

  6. ##FILTER=<ID=panel_of_normals,Description="Seen in at least 2 samples in the panel of normals">

  7. ##FILTER=<ID=str_contraction,Description="Site filtered due to contraction of short tandem repeat region">

  8. ##FILTER=<ID=t_lod_fstar,Description="Tumor does not meet likelihood threshold">

  9. ##FILTER=<ID=triallelic_site,Description="Site filtered because more than two alt alleles pass tumor LOD">

旧版的mutect只是对一个位点REJECT或者PASS，但是新版增加了多种情况来解释为什么REJECT，就是上面的那些tag的组合。

即使我把N-T反过来用mutect2来call somatic mutation，仍然会有125个位点PASS，只需要是在tumor里面纯粹的野生型，在normal里面是AF非常低的杂合即可。

其次是varscan

java软件都是下载即可使用，官网可以下载，我把它下载到了 ~/biosoft/VarScan/VarScan.v2.3.9.jar目录，使用代码是：

  1. TMPDIR=/home/jianmingzeng/tmp/software

  2. normal_bam=NPC10F-N_recal.bam

  3. tumor_bam=NPC10F-T_recal.bam

  4. sample=NPC10F

  5. normal_pileup="samtools mpileup -q 1 -f $reference $normal_bam";

  6. tumor_pileup="samtools mpileup -q 1 -f $reference $tumor_bam";

  7. # Next, issue a system call that pipes input from these commands into VarScan :

  8. java -Djava.io.tmpdir=$TMPDIR   -Xmx40g  -jar ~/biosoft/VarScan/VarScan.v2.3.9.jar \

  9. somatic <($normal_pileup) <($tumor_pileup) $sample

  10. java -jar ~/biosoft/VarScan/VarScan.v2.3.9.jar processSomatic $sample.snp

结果文件如下：

  1.  101353 NPC10F.snp

  2.   95145 NPC10F.snp.Germline

  3.    4831 NPC10F.snp.Germline.hc

  4.    3503 NPC10F.snp.LOH

  5.    3485 NPC10F.snp.LOH.hc

  6.    2343 NPC10F.snp.Somatic

  7.     557 NPC10F.snp.Somatic.hc

其中值得关注的就是 NPC10F.snp.Somatic.hc文件里面的位点，打开浏览可以看出其实就是normal样本里面的野生型变化到tumor里面的杂合型，就是allel frequency的增加。

接着是muse

首先需要下载安装

  1. cd ~/biosoft

  2. mkdir muse &&  cd muse

  3. wget http://bioinformatics.mdanderson.org/Software/MuSE/MuSEv1.0rc_submission_b391201 muse

  4. mv MuSEv1.0rc_submission_b391201 muse

  5. chmod 777 muse

该软件需要的input跟mutect2一样，都是：

  • (1) the indexed reference genome FASTA file,

  • (2) the binary sequence alignment/map formatted (BAM) sequence data from the pair of tumor and normal DNA samples, and

  • (3) the dbSNP variant call format (VCF) file that should be bgzip compressed, tabix indexed and based on the same reference genome as (1).

它给的测试代码是：

  1. ./MuSE call –O Output.Prefix –f Reference.Genome Tumor.bam Matched.Normal.bam

  2. ./MuSE sump -I Output.Prefix.MuSE.txt -G –O Output.Prefix.vcf –D dbsnp.vcf.gz

需要分成两步来运行，修改一下就可以处理自己的样本啦。

  1. reference=/home/jianmingzeng/reference/genome/human_g1k_v37/human_g1k_v37.fasta

  2. DBSNP=/home/jianmingzeng/annotation/variation/human/dbSNP/All_20160601.vcf.gz

  3. normal_bam=NPC10F-N_recal.bam

  4. tumor_bam=NPC10F-T_recal.bam

  5. sample=NPC10F

  6. ~/biosoft/muse/muse call -O $sample -f $reference $tumor_bam $normal_bam

  7. # ~/biosoft/muse/muse sump -I ${sample}.MuSE.txt -E –O ${sample}_muse.vcf –D $DBSNP

看看结果文件，就是典型的VCF格式，而且tags不多值得注意的只有 ID=SS,Number=1,Type=Integer,Description="Variant status relative to non-adjacent Normal,0=wildtype,1=germline,2=somatic,3=LOH,4=post-transcriptional modification,5=unknown

前5列 CHROM POS ID REF ALT很正常，第6列QUAL全部是点，第7列FILTER 把位点分级了，统计如下：

  1.     81 Tier1

  2.     26 Tier2

  3.     12 Tier3

  4.     26 Tier4

  5.     39 Tier5

第8列是 INFO 信息，全部是SOMATIC

第9,10,11列是 GT:DP:AD:BQ:SS 格式的tumor和normal，可以看到normal都是野生型0/0, tumor全部是杂合突变1/0，只是allel frequency不同而已，介于0~1之间。

最后是SomaticSniper

第一次用这个软件，也是需要先安装

  1. cd ~/biosoft

  2. mkdir SomaticSniper &&  cd SomaticSniper

  3. wget https://github.com/genome/somatic-sniper/archive/v1.0.5.0.tar.gz

  4. tar zxvf v1.0.5.0.tar.gz

  5. cd somatic-sniper-1.0.5.0

  6. mkdir build && cd build

  7. cmake ../

  8. make deps

  9. make -j

  10. make test

服务器的make好像有点问题，没有安装成功，先略过吧。

还看看strelka吧

这款软件专门设计给临床使用的，说明书非常详细 https://github.com/Illumina/strelka/tree/master/docs/userGuide 在github上面可以找到最新版地址，安装如下：

  1. cd ~/biosoft

  2. mkdir strelka &&  cd strelka

  3. wget https://github.com/Illumina/strelka/releases/download/v2.8.2/strelka-2.8.2.centos5_x86_64.tar.bz2

  4. tar xvfj strelka-2.8.2.centos5_x86_64.tar.bz2

其实就是下载了python脚本，如下：

也是需要tumor和normal的bam文件，还有参考基因组文件，就可以进行找somatic mutation啦。

  1. reference=/home/jianmingzeng/reference/genome/human_g1k_v37/human_g1k_v37.fasta

  2. DBSNP=/home/jianmingzeng/annotation/variation/human/dbSNP/All_20160601.vcf.gz

  3. normal_bam=NPC10F-N_recal.bam

  4. tumor_bam=NPC10F-T_recal.bam

  5. sample=NPC10F

  6. cd ~/data/public/npc/

  7. mkdir -p somatic/strelka

  8. ~/biosoft/strelka/strelka-2.8.2.centos5_x86_64/bin/configureStrelkaSomaticWorkflow.py \

  9. --normalBam $normal_bam \

  10. --tumorBam $tumor_bam \

  11. --referenceFasta $reference \

  12. --runDir  somatic/strelka

上面的代码只是生成了一个python脚本，还需要自行提交该脚本才能得到结果。

Successfully created workflow run script. To execute the workflow, run the following script and set appropriate options:

需要多读软件的FAQ：https://sites.google.com/site/strelkasomaticvariantcaller/home/faq

赞 (0)

相关推荐