[gvcf liftover] INFO/END does not update

[YingZhou001]

Hi Giulio,

Thanks for developing this liftover tool, it is easy to use and very efficient.

I was trying to liftover g38 to chm13 with command

bcftools +liftover -Ou ${gvcf_g38} --\
  --src-fasta-ref ${ref_g38} \
  --fasta-ref ${ref_chm13} \
  --chain ${chain} \
  --write-reject  \
  --lift-end \
  --reject-type z \

the input g38 data was

chr19   60963   .       G       <NON_REF>       .       .       END=60969       GT:DP:GQ:MIN_DP:PL      0/0:19:57:19:0,57,671
chr19   60970   .       A       <NON_REF>       .       .       END=60970       GT:DP:GQ:MIN_DP:PL      0/0:19:41:19:0,41,648
chr19   60971   .       T       <NON_REF>       .       .       END=60976       GT:DP:GQ:MIN_DP:PL      0/0:19:54:19:0,54,810

but the output didn't update the 'END' field

chr19   11946   .       G       <NON_REF>       .       .       END=60969       GT:DP:GQ:MIN_DP:PL      0/0:19:57:19:0,57,671
chr19   11953   .       A       <NON_REF>       .       .       END=60970       GT:DP:GQ:MIN_DP:PL      0/0:19:41:19:0,41,648
chr19   11954   .       T       <NON_REF>       .       .       END=60976       GT:DP:GQ:MIN_DP:PL      0/0:19:54:19:0,54,810

The variable 'end_info->v1.i' had the updated position value but didn't show in the final output.

I tested with both bcftools 1.20 and 1.21 and all shared the same results.

Could you help me check this?

Best,

Ying

[freeseek]

The problem shows up because you are using -Ou. I have added the line:

bcf_update_info_int32(args->out_hdr, rec, "END", &end_info->v1.i, 1);

See if the development version now works for you

[YingZhou001]

Thanks, it fixed the problem on the testing data, but for a larger data set, it stopped with

[W::bcf_update_info] INFO/END=0 is smaller than POS at chr1:26434592

According to #35, I think it could be reset as "MISSING", right?

[update]

the 'end_info->v1.i' value was -9223372036854775808, which is the same as bcf_int64_missing.

[freeseek]

Can you show me the variant that caused you this issue so that I can reproduce it?

[YingZhou001]

yes, the chain file is hg38ToHs1.over.chain.gz. Please let me know if you need other information.

chr1    26596641        .       C       <NON_REF>       .       .       END=26596641    GT:DP:GQ:MIN_DP:PL      0/0:1:3:1:0,3,42
chr1    26596642        .       A       <NON_REF>       .       .       END=42709142    GT:DP:GQ:MIN_DP:PL      0/0:0:0:0:0,0,0
chr1    42709143        .       T       <NON_REF>       .       .       END=42709173    GT:DP:GQ:MIN_DP:PL      0/0:1:3:1:0,3,15

[freeseek]

It should be fixed now

[YingZhou001]

Thanks Giulio but unfortunately it still failed. It said:

Writing to /tmp/bcftools.Tn6Xv1
INFO/MLEAC is handled by AGR rule
INFO/MLEAF is handled by AGR rule
FORMAT/AD is handled by AGR rule
FORMAT/PL is handled by AGR rule
Warning: input VCF includes symbolic alleles that might not properly lift over
[W::bcf_update_info] INFO/END=-128 is smaller than POS at chr1:26434592
[W::bgzf_read_block] EOF marker is absent. The input may be truncated
Merging 0 temporary files
Done
Cleaning

Here I attached the input data:
HG002.g38.g.vcf.gz

[freeseek]

It should be fixed now

[YingZhou001]

Thanks Giulio, the update works for the testing data I shared but failed again with a larger one. Still there are lots of END values that less than the POS values leading to failure in converting gvcf to vcf with PICARD. Let me investigate more about the failures and back to you later.
By the way, is there any document I can learn how the bcftools/liftover deals with gVCF?

[freeseek]

Okay, happy to debug if you can send me another example. It really does not deal with gVCF in any way different from a regular VCF. It is discouraged to liftover gVCFs though. See the main paper