sam文件是短序列比对生成的文件，是二代测序中最核心的文件。在RNAseq，变异检测等分析中，都需要首先生成sam文件格式。bam文件是sam格式的二进制格式，转换为二进制之后，可以减小文件的存储。掌握sam/bam文件的操作是处理二代测序数据的非常重要的内容，例如sam与bam的转换，排序，建立索引，reads计数等等操作。
21 mpileup
mpileup是将reads堆叠起来，作用是展示每一个位点比对的细节，然后可以通过这个细节来推测比对情况，pileup是之前maq工具查找SNP重要的步骤，目前已经逐渐将mpileup功能转移到bcftools中了。 samtools mpileup --reference ref/ref.fna A1.sorted.bam
22 MarkDuplication
Dupliacation reads会对变异检测造成干扰，得到一些假阳性的结果，因此，需要将这些reads去除掉，只留一份即可。可以在比对之前去除，但是比较消耗内存。也可以在比对之后进行标记。这一步骤只是在每一行比对结尾出添加一些CIGAR标志，并不过滤数据。samtools可以标记Duplication，也可以去除掉reads，GATK也可以进行标记。 gatk MarkDuplicates -I A1.sorted.bam -M Sample1.markdup_metrics.txt -O A1.sorted.markdup.bam
samtools index A1.sorted.markdup.bam
23 BQSR
BQSR是Base Quality Score Recalibratio的简称，是利用GATK做人全基因组分析的必要步骤。首先利用BaseRecalibrator进行机器学习，然后利用ApplyBQSR应用于比对结果。但是需要注意，BQSR至少要达到100万个突变位点，因此，外显子和目标区域捕获一般达不到这个要求，也就做不了这步奏。 #进行学习
gatk BaseRecalibrator \
-R Homo_sapiens_assembly38.fasta \
-I A1.sorted.markdup.bam \
--known-sites 1000G_phase1.snps.high_confidence.hg38.vcf.gz \
--known-sites Mills_and_1000G_gold_standard.indels.hg38.vcf.gz \
--known-sites dbsnp_146.hg38.vcf.gz \
-O A1.sorted.markdup.recal_data.table >bqsr.log
#进行BQSR
gatk ApplyBQSR \
--bqsr-recal-file A1.sorted.markdup.recal_data.table \
-R Homo_sapiens_assembly38.fasta \
-I A1.sorted.markdup.bam \
-O A1.sorted.markdup.BQSR.bam
#对新的bam创建索引
samtools index A1.sorted.markdup.BQSR.bam
24 利用bcftools进行SNP检测
可以使用bcftools直接来筛选SNP，输入排序并建立索引的bam即可，如果能做Mark Duplication则更好了。 #利用bcftools进行变异检测
bcftools mpileup -f ref.fna A1.sorted.bam | bcftools call -c -v -o A1.bcftool.vcf
25 利用freebayes进行snp检测
与bcftools类似，使用freebayes来筛选SNP也非常容易，输入排序并建立索引的bam。 freebayes -f ref.fna -C 5 A1.sorted.bam -v A1.freebayes.vcf
---------- END ----------
|
