Entering edit mode
3.9 years 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
any suggestions, please? I am stuck with the analysis.