/vcf_to_gds_overview

Description of VCF to GDS conversion on the UK Biobank RAP

Primary LanguageShell

VCF to GDS Conversion

Description of VCF to GDS conversion on the UK Biobank RAP for STAAR processing

Due to the number of vcfs and their file sizes on the DNAnexus platform, it may be necessary to reduce the amount of data contained within VCFs if wanting to merge all VCFs associated with a chromosome on a workstation for subsequent STAAR annotations and processing. Some applets are provided here that may help facilitate this process. The examples provided below focus on WGS data. More examples will be provided soon to handle WES data.

Step 1: Trimming down data in the VCFs

See vcf_trimmer for an applet that removes fields and performs required filtering through bcftools. This applet takes in a list of files as an input to process. This applet could be used across multiple jobs submitted on the DNAnexus platform which would require unique lists of VCFs to be split. Note, left-alignment and normalisation not presently set.

Example: Processing WGS (200K) VCFs for chromosome 17 over 100 jobs:

Generating the input file lists:

If you wanted to process the 200K WGS data release for chromosome 17 then vcf_trimmer would require a list of VCFs associated with that chromosome. Furthermore, you may also want to submit 100 jobs whereby the list of VCFs associated with chromosome 17 is split into 100 unique and equally sized VCF lists. Running the get_vcf_file_list_by_chr_and_split.sh bash script in this repo will produce all the relevant input files for chromosome 17 (chr17_vcf_list_[1:100]) as well as for all other chromosomes for the 60,648 VCFs currently available.

bash get_vcf_file_list_by_chr_and_split.sh

This script has been provided as a starting point that can be amended to suit your requirements. Note, the file lists generated by this bash script will need to be subsequently uploaded to a project folder on the DNAnexus platform.

Submitting the jobs:

Once you have a set of 100 files listing the VCFs uploaded to a project folder on the RAP, you can run the vcf_trimmer with each of the files generated above. For example, if we wanted to run vcf_trimmer on chromosome 17 to remove FORMAT fields except GT and keep only variants with AAscore > 0.5, then the command would look something like this:

for i in {1..100}
do
  dx run /path/to/vcf_trimmer \
    -ivcf_file_list=/path/to/chr17_vcf_list_${i} \
    -ifile_label=trimmed \
    -ioutput_dir=/path/to/output/dir \
    -iqc_thresholds="INFO/AAScore>0.5" \
    -ifields_to_remove="FORMAT/FT,FORMAT/AD,FORMAT/MD,FORMAT/DP,FORMAT/RA,FORMAT/PP,FORMAT/GQ,FORMAT/PL" \
    -y
done

Step 2: Merging VCFs

See vcf_merger. Another applet wrapper for bcftools. This applet also takes in a list of VCF files. Continuing the example above, to obtain a list of files in block order of a directory on the RAP for chromosome 17 after running vcf_trimmer, type:

dx ls "/path/to/output/dir/ukb24304_c17_b*_v1_trimmed.vcf.gz" | sort -t"b" -k3.1 -n | awk -v d="/path/to/output/dir" '{print d $0}' > chr17_vcfs_to_merge

The output file should be subsequently uploaded to the RAP for processing with vcf_merger.

Step 3: Converting VCF to GDS for subsequent annotation/analysis in the STAARpipeline

See vcf2gds. This applet comes with an R library that will be unpacked during runtime and used for data conversion.