Vg call in a specific region
1
0
Entering edit mode
17 months ago

Hi, I am trying to call variants in a narrow region using vg chunk and vg call. It worked with an older version of vg tools but does not seem to work with the newer one.

In a pervious version, this worked for me

vg chunk -x ref_alt.xg -a SRR5647740_paired_wasp_mapped_sorted.gam -g -p 3:3716672-3718990 -c 10 > chunk.vg

for file in chunk_0_3_3716352_3719310.gam; do vg augment -a pileup -t 20 -Z ${file/.gam/_aug.trans} -S${file/.gam/_aug.support} chunk.vg $file >${file/.gam/_aug.vg} 2>err_vg_aug ; done

for file in chunk_0_3_3716352_3719310_aug.vg; do vg call -z ${file/_aug.vg/_aug.trans} -s${file/_aug.vg/_aug.support} $file >${file/_aug.vg/_calls.vcf}; done

I tried to adapt it to the current release of vg tools and this is what I have tried

vg chunk -x ref_alt.xg -a SRR5647740_paired_wasp_mapped_sorted.gam -g -p 3:3716672-3718990 -c 10 > chunk.vg

vg augment chunk.vg chunk_0_3_3716352_3719310.gam -A aug.gam > aug.vg

vg index aug.vg -x aug.xg

vg pack -x aug.xg -g aug.gam -Q 5 -o aln_aug.pack

vg call aug.xg -k aln_aug.pack > SRR5647740_calls.vcf

Is this correct? Its producing a vcf file but with no variants. The old version for the same input file did. Thanks!

Vg call vg tools vg chunk vg • 733 views
1
Entering edit mode
17 months ago
glenn.hickey ▴ 250

It's less necessary to break the graph into tiny chunks with vg chunk -p. Here's an example on chunking by chromosome:

https://github.com/vgteam/vg/wiki/Whole-genome-calling-and-genoyping#calling-the-graph-by-first-splitting-into-components

That said, workflows like that should still work. You can use vg depth to verify the coverage in you GAM and pack files. Lack of coverage is the only reason I could think of why call would produce no output. Perhaps the -Q 5 filter in the new workflow is filtering everything out?

0
Entering edit mode

Thanks so much! I found out the error had nothing to do with vg tools. It was some odd error encountered when running the program on a HPC. Got it resolved. Thanks!