Hi - I have been doing some analysis of published cancer sequencing data, and I have generated lists of somatic (cancer-specific) variants, filtered based on various criteria, for which I would like to query gnomAD. (Essentially to ask, how many of them have never/almost never been seen in germline sequencing of healthy individuals, either due to strong functional effects that have been selected against, or because they arise via mutagenic processes that are not operating in the germline?) I put "uncharacterized" in scare quotes because of course many of these variants have been characterized, but not by me and not in terms of their frequency in the general population.
I can format my lists in whatever format is appropriate, e.g. as VCF or MAF files. Right now I'm doing the analysis in R, and the tables are represented as GRanges objects - so I have the genome position, control allele, mutant allele, etc., but no SNP IDs. I'd like to be able to input that file and get back a report of the gnomAD allele frequency (and potentially other annotations) for each one. Is there a simple way to do this? (Bonus points if it can be done within R, but I'm fine with command-line approaches if necessary.)
Read the gnomAD VCF file in R and do a lookup? Heavy-handed but will work.
Thanks - this is probably what I will have to do, but my understanding is that the gnomAD VCF file is really enormous (almost 500 GB, assuming I’m looking at the right thing), so this will present some major practical difficulties.
Save just the CHR, POS, REF, ALT, ID and required frequency fields in a tab-delimited file, then use
data.table::freadto read it in. Do it on a cluster node with ~64GB RAM if possible. If not, switch to command line (bcftools/vt/vep) but in any case you'll need a bit of compute power to get this done.