Mapping SNP rsIds/ positions to genes coding for proteins

Question

[Harsh]

I have started a new project and am a beginner in bioinformatics. I have a summary statistics file from a GWAS study. The file has information like the Chromosome number, position of the SNP. The rsids are not included. Now I want to map these SNPs to Genes that code for proteins. With the condition: SNP should 50kb upstream or downstream of the gene boundary, or inside the gene boundary.

Thank You

[Pierre Lindenbaum]

extract genes and functions +/-50kb around a SNP ; How To Map A Snp To A Gene Around +/- 60Kb ?

Answer 1

[Alex Reynolds]

Here's one way to do it:

To get SNPs, assuming you are working with hg38 human GWAS data:

wget -qO- https://ftp.ncbi.nlm.nih.gov/snp/organisms/human_9606_b151_GRCh38p7/VCF/All_20180418.vcf.gz \
    | gunzip -c \
    | vcf2bed --sort-tmpdir=/tmp --max-mem=2G - \
    | awk -v FS="\t" -v OFS="\t" '{ print "chr"$0 }' \
    > hg38.dbSNP155.bed

Take the difference of padded and unpadded gene regions to get the genomic space covering upstream and downstream of the gene boundaries, but not including the gene (if I am understanding your problem's condition correctly):

bedops --difference <(bedops --range 50000 --everything genes.bed | bedops --merge -) genes.bed > regions.bed

Then map the SNPs to the regions:

bedmap --echo --echo-map-id --delim '\t' regions.bed SNPs.bed > answer.bed

If your condition is that the SNP can be within the gene or 50kb upstream or downstream of the gene boundary, then you can skip a lot of the complication and just map directly:

bedmap --echo --echo-map-id --delim '\t' <(bedops --range 50000 --everything genes.bed) SNPs.bed > answer.bed

If you want to remove the 50kb of padding in the answer:

bedmap --echo --echo-map-id --delim '\t' <(bedops --range 50000 --everything genes.bed) SNPs.bed | bedops --range -50000 - > answer.bed

[Harsh]

Actually the condition is, SNP should be within the gene boundary or 50kb upstream or 50kb downstream of the gene boundary.

[Alex Reynolds]

See edited answer.

[Harsh]

thank you!