Is variant calling done on a per-position basis? I've read recommendations to split the BAM file by chromosomes and parallel call the chromosomes for speed
Could I further segment each chromosome into chunks and do calling on each chunk? For example if I: 1) Split a chromosome into 1MB segments. 2) Parallel variant call each 1MB segment. 3) Concatenate the VCF.
Would the resulting concatenated file be correct? Would I be missing any information that might be shared among sites on the same chromosome that variant callers use?
•
written
by
Damian Kao
2 answers
Brad Chapman has described this procedure previously on his blog here (search for "Parallelism by genomic regions"). This at least used to be part of the bcbio-nextgen workflow, though I don't know if it's still included. That procedure is slightly more involved since it finds appropriate segments, rather then using 1MB chunks, but the principle is the same.
Hello Damian,
Is variant calling done on a per-position basis?
This depends on your variant caller. GATKS UnifiedGenotyper and samtools mpileup do so as far as I know. GATKs HaplotypeCaller and freebayes are doing local denovo assembly. So these variant caller need information around a suspected variant.
Instead of splitting the bam file, I would use the possibilty of most variant caller to do the calling within a given genomic region. So every process have than every alignment information it needs.
fin swimmer
Log in to answer this question.
-
GATK - Output differences in Calls when Split by chromosome intervals and without splitting
written by Maverick •
-
Further splitting single-chromosomes VCF files to decrease their size to annotate with VEP
written by mmorad •
-
Merge overlapping regions of a single chromosome which were phased separately
written by Georgi •
-
Split bam into chunks but keep reads mapped to the same region together
written by M9Q6Ixc31Fc8 •
-
VCF merge or concatenate?
written by rokragna295 •
-
Merge contigs in fasta file
written by Rubal
-
How does samtools mpileup handle multiple .bam inputs?
written by traviata •
-
Abyss genome assembly input
written by Damian Kao
-
Split Bam Files By Region For Parallel Variant Calling
written by William
-
Extract Fastq From Large Bam File
written by michealsmith
If you plan on calling SNPs only then I don't see a problem. However, if you are looking for structural variants, there would be missing data on the edges of the chunks, especially for SVs spanning multiple chunks. Additionally, how do you plan to keep track of the size or exact position, additional liftUp files?
Thanks for the reply. I see what you mean with the SVs and possibly even indels. I am not interested in SVs for now, but do want to preserve indel information if I can.
The BAM files I am working with are low coverage. I guess I'll to write a script to chunk the BAM file based on coverage "islands" where each island should be at least 100bp apart.