empty vcf from vg call
2
0
Entering edit mode
3 months ago
philgev2 ▴ 10

Hi

I am getting an empty vcf file when I run below commands

vg construct -t 10 -r /data/public-data/hg38/hg38_v2_Homo_sapiens_assembly38.fa -v /data/public-data/1000genomes/ALL.wgs.shapeit2_integrated_snvindels_v2a.GRCh38.27022019.sites_chr_1:Y.vcf >1000g.vg

vg prune -k 12 -e 2 -t 16 -p 1000g.vg >1000g.pruned.vg

#### store the graph in the xg/gcsa index pair
vg index -t 10 -p -b /data/users/philge/tmp/ -x 1000g.xg -g 1000g.gcsa -k 16 1000g.pruned.vg

####align a read to the indexed version of the graph
vg map -t 10 -b /data/users/philge/genomegraph/vg/dbsnp/sorted_A019.bam -x 1000g.xg -g 1000g.gcsa > ./alignment/sorted_A019_1000g.gam

####Compute the read support from the gam (with ignoring qualitiy < 5)
vg pack -t 10 -x 1000g.xg -g ./alignment/sorted_A019_1000g.gam -Q 5 -o ./alignment/sorted_A019_1000g.pack

####Generate a VCF from the support
vg call -t 10 -r ./1000g.snarls ./1000g.xg -k ./alignment/sorted_A019_1000g.pack > ./alignment/sorted_A019_1000g_no_aug.vcf

All my input files for vg call contain data. But output vcf is empty. Any suggestions on what could be the reason?

Thanks
Philge

vg • 226 views
ADD COMMENT
0
Entering edit mode

any suggestions, please? I am stuck with the analysis.

ADD REPLY
1
Entering edit mode
3 months ago
glenn.hickey ▴ 230

The commands look okay (though vg snarls seems to be missing). Assuming your snarls file is fine, the other reason for an empty VCF could be not enough mappings found in the pack file. You can use vg depth -k -b 10000 to verify the coverage in your pack file (and with -g to verify the coverage in your GAM).

ADD COMMENT
0
Entering edit mode
11 weeks ago
kingcohn ▴ 20

I have a similar issue, did you resolve this?

ADD COMMENT

Login before adding your answer.

Traffic: 1106 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6