4.1.6: final chunk fails with "failed loading target haplotypes" depending on target BGZF block layout

[mattbarton]

Summary

On minimac4 4.1.6, imputation of the final chunk of a chromosome fails with

Imputing 20:40000001-58114749 ...
Loading target haplotypes ...
Error: failed loading target haplotypes

while every earlier chunk of the same run loads and imputes normally.

The failure appears to depend on the BGZF block layout of the target VCF rather than on its contents. Rewriting the same target through bcftools view --output-type z and re-indexing makes the identical command on the identical region succeed, even though the two files contain byte-identical records.

This looks distinct from #76, which was diagnosed as inconsistent sample ploidy between adjacent variants. Here the record bodies of the working and failing targets are identical, so no content difference can explain it, and the error arrives with no accompanying sample/position detail.

Environment

• minimac4 4.1.6 (conda-forge, linux-64)
• Linux x86_64, 8 vCPU / 32 GiB, no memory pressure (29 GiB free, no OOM events)
• Target VCF written by Eagle 2.3.5 with --vcfOutFormat=z
• Reference in .msav form; chromosome length 58,114,749 bp; ~3.9k target samples
• Command:

minimac4 --all-typed-sites --format GT,DS,GP,HDS --threads 7 \
    -o out.vcf.gz reference_panel_chr20.msav target_chr20.vcf.gz

Default --chunk 20000000 and --overlap 3000000, so the chunks are 20:1-20000000, 20:20000001-40000000, 20:40000001-58114749. The first two impute all samples; the third fails as above.

The distinguishing observation

Rewriting the target with no intended change to its contents fixes it:

bcftools view target_chr20.vcf.gz --output target_rewritten.vcf.gz --output-type z
tabix -f -p vcf target_rewritten.vcf.gz
# same minimac4 command against target_rewritten.vcf.gz -> all three chunks impute

Comparing the two files:

	Eagle output	bcftools rewrite
bcftools view --no-header md5	identical	identical
record count	6154	6154
BGZF blocks	1493	1543
max uncompressed block	65280 B	65280 B
blocks > 64 KiB	0	0
28-byte BGZF EOF marker	present	present

Both files are valid BGZF and bgzip -t passes on both. The only difference I can find is where the block boundaries fall (and the ##bcftools_viewCommand provenance line bcftools adds to the header).

Region bisection

Using -r against the original Eagle-written target:

region	result
20:40000001-50000000	OK (2033 typed sites)
20:56000000-58114749	OK (749 typed sites)
20:40000001-58000000	fails
20:40000001-58114748	fails
20:40000001-58114749 with -w 0	OK

So it is not a single malformed record, not a specific position, and not region size — the two 20 Mb chunks earlier in the same file are larger than the 18 Mb region that fails.

The -w 0 result is the interesting one: with the default 3 Mb overlap the final chunk's read window extends to ~20:61114749, i.e. about 3 Mb past the last variant in the file, whereas -w 0 stops at it. That is the only structural difference I can see between the final chunk and the earlier ones.

Ruled out by test

• Contig header length. The target declared ##contig=<ID=20,length=58090568> while holding variants out to 58,114,749. Reheadering to 58,114,750 and again to 62,000,000 (comfortably past chunk end + overlap) both still failed.
• Memory — 29 GiB free, no OOM events; the window holds only 2782 sites x ~3.9k samples.
• Missing BGZF EOF marker — present and correct in both files.
• Stale or bad index — a fresh tabix -f over the original blocks still failed; only rewriting the blocks helped.
• Genotype content — missing (./.) and unphased calls appear in all three chunks, including the two that impute fine.

Question

Is the target reader sensitive to compressed-block boundaries when the requested read window extends beyond the last variant in the file? If so, a final chunk whose overlap runs past the end of the data would fail or succeed depending only on how the writer happened to lay out its BGZF blocks, which matches everything observed here.

Happy to run further diagnostics. I can't attach the original files as the genotypes are not mine to share, but I can test patches or produce additional measurements on request.

[jonathonl]

Thanks for providing this thorough report. This seems like a bug in processing some corner case of the tabix index. Can you rerun the original VCF after generating a CSI index (bcftools index target_chr20.vcf.gz)?