samtools help

[guang918]

Hi,

I had two Bam (1.bam, 2.bam) files from Bowtie2 and tried to call variants by Samtools. I have tried two procedures and they gave me different answers.

Procedure #1:
samtools mpileup -uf reference.fa bam1 bam2 | bcftools view -bvcg - > 12.bcf
bcftools view 12.bcf | vcfutils.pl varFilter -D100 12.flt.vcf

Procedure #2:
samtools merge 12_merged.bam 1.bam 2.bam
samtools mpileup -uf reference.fa 12_merged.bam | bcftools view -bvcg - > 12_merged.bcf
bcftools view 12_merged.bcf | vcfutils.pl varFilter -D100 12_merged.vcf

I can't figure out why. Please give some suggestions.

Thanks very much. Happy Thanks giving day.

[dpryan]

Well, in one case you're calling SNPs on multiple samples and in the other on a single sample with higher depth. I wouldn't expect them to give the same results.

[guang918]

Thanks very much for your reply.

I have one biological sample. And the Illumina sequencing will give me several FASTQ files, 1.fq, 2.fq... If I want to identify SNPs to reference genome, should I align them to genome separately, or combine the sequence first and align them as a single file?

Thanks a million.

[dpryan]

It won't much matter if you align them separately and then merge or concatenate the files and then align the results (if you have paired-end reads, some aligners re-estimate the insert-size distribution throughout the alignment, so that could change things a bit). Whatever you do, don't create multiple BAM files from the same sample and then treat them as multiple biological samples.