Speed up extract_tracts.py step for large vcf

[tulsi92]

I am currently running the extract_tracts.py step just on chr19 broken down into 4cM chunks. There are >30k individuals in the VCF and still ~500k variants in each chunk so it is taking a long time to run. It has been running for a few days already and is not even halfway through the variants. Just wondering if it should be taking this long and if there is there any way to speed this step up?

Example of code I am running below in case any of the output params are slowing it down:

python3 extract_tracts.py \
    --vcf hg38_chr19.chunk_chr19_37780539-41942739.phased.vcf.gz \
    --msp hg38_chr19.msp \
    --num-ancs 5 \
    --output-dir output/19 \
    --output-vcf \
    --compress-output

[shahnirav1005]

Hi @tulsi92,

The extract_tracts.py step can be time-consuming because it cross-references each variant with ancestry data from the ancestral windows in the MSP file and generates haplotype count (hapcount) and dosage files. This process is explained in the Tractor tutorial here. In short, the script performs these operations for each individual in your dataset, which requires significant amount of time.

Given the size of your dataset, the process is likely taking longer. The good news is that you only need to run extract_tracts.py once. Afterward, you can use the output files (hapcount and dosage) with any number of phenotypes.

Here are a few suggestions to speed it up:

1. Consider chunking your files further by variants, then merging them. You could also chunk by individuals, but this requires creating temp. MSP and VCF files, as extract_tracts.py expects a 1-to-1 match between these files. Once processed, you can merge them using paste and cat.
2. The --compress-output flag compresses files during saving, which can slow things down. If you have enough storage, try removing this option, saving files uncompressed, and then manually compressing them later.
3. While not a recommended solution (use caution), if some ancestries account for a small percentage of your dataset, you could reduce the number of ancestries used in the local ancestry step. However, this may slightly increase error rates in your LAI calls, so proceed with caution. Reducing the number of ancestries should also improve speed, as it's currently generating 5 hapcount and 5 dosage files.

I hope this helps! Let us know if you have any other questions. You can also email me at nirav.shah@bcm.edu.
Nirav

[pettyalex]

My team has had some incredibly slow extract_tracts.py runtimes, so I took a look at it and found at least one major speed improvement: #45