2 answers
Assuming you have the positions of your genes in a sorted BED file, you can skip the first step of downloading genes and filtering them for those of interest.
Download a list of genes:
$ wget -qO- ftp://ftp.sanger.ac.uk/pub/gencode/Gencode_human/release_26/gencode.v26.basic.annotation.gff3.gz \
| gunzip -c \
| convert2bed --input=gff - \
> gencode.v26.bed
Then grep for your genes of interest (names should match Gencode gene names, ENSG* etc.):
$ LC_ALL=C && grep -F -f genes-of-interest.txt gencode.v26.bed > gencode.v26.goi.bed
Once you have genes of interest and their positions in a sorted BED file, you can grab the SNPs and put those into BED format:
$ wget -qO- http://hgdownload.cse.ucsc.edu/goldenpath/hg38/database/snp147.txt.gz \
| gunzip -c \
| awk -v OFS="\t" '{ print $2,$3,($3+1),$5 }' \
| sort-bed - \
> hg38.snp147.bed
Now you have genes, SNPs, and their positions. You can map SNPs to genes within a 1M window:
$ bedmap --echo --echo-map-id --range 1000000 --delim '\t' gencode.v26.goi.bed hg38.snp147.bed > answer.bed
The file answer.bed will contain the gene records on each line. The last column of each line will be a tab-delimited list of SNP IDs, which overlap the gene within 1Mb.
BEDOPS and GNU tools: http://bedops.readthedocs.io/en/latest/
Those look like HGNC names. If you're working with hg38:
$ wget -qO- http://hgdownload.cse.ucsc.edu/goldenpath/hg38/database/refGene.txt.gz \
| gunzip -c \
| awk -v OFS="\t" '($9>1){ print $3,$5,$6,$13,$9,$4 }' \
| sort-bed - \
> hg38.hgnc.bed
Then grep as before:
$ LC_ALL=C && grep -F -f genes-of-interest.txt hg38.hgnc.bed > hg38.hgnc.goi.bed
The rest is the same.
I don't understand your question. You already have a text file containing the genes you are interested in, correct? Your original question says you have a list of genes.
I'm afraid I don't understand. Check that you have the files listed in these commands, perhaps?
Hi mms140130,
If an answer was helpful you should upvote it, if the answer resolved your question you should mark it as accepted.
After you have made a final BED file, you could convert it to CSV with something like:
$ cat in.bed | tr '\t' ',' > in.csv
Otherwise, it is already a (tab-delimited) text file. If you need to rename it from some-file.bed to some-file.txt, use the mv command.
Hi, You can use BiomaRt package to grab this information
library(biomaRt)
listMarts()
snps<-useMart("ENSEMBL_MART_SNP")
listDatasets(snps)
Then you choose whatever data-set you are interested in, lets say "hsapiens_snp" which includes Human Short Variants (SNPs and indels excluding flagged variants) (GRCh38.p10).
snps<-useDataset("hsapiens_snp",mart = snps)
listFilters(snps)
Then you can choose "chromosaomal_region" filter, to get the all information you need about SNPs in that region.
P.S. you have to get the coordinates for your genes first, then put this coordinates in an object (name it for example coordinates)
snps_info<-getBM(filters="chromosomal_region",attributes=c("refsnp_id","chr_name","allele","minor_allele","minor_allele_freq","clinical_significance","polyphen_prediction","polyphen_score","sift_prediction","sift_score"),values=coordinates,mart=snps)
If you have the coordinates (chromosomal region) of your genes, you can just add 1Mb to the start and the end coordinates.
Log in to answer this question.
-
Investigating differentially expressed genes
written by vanessaward876 •
-
Is there a GRCh37 equivalent of EnsDb.Hsapiens.v86 or some other method for getting a complete list…
written by SpacemanSpiffo •
-
what is the difference between snps files in UCSC download
written by mms140130
-
Cis- eQTL Question
written by mms140130
-
Finding Low Quality SNPs
written by mms140130
-
Extract genes from gene list
written by Bioblazer
-
how I can find the exact position of genes encode polymerases in genomes?
written by zizigolu
-
Tissue Specific Gene-Gene Similarity (Co-Expression)
written by GenomicEnthusiast •
-
How To Use Plink To Compute The Ld Between Two Sets Of Snps?
written by Fayue1015 •
-
Calculating Gene Densities From Rnaseq Data
written by Ahmetz •
actually I applied your code it gave me one vector of snps while I need for each gene the snps that is within 1mb
Please use
ADD COMMENT/ADD REPLYto properly place this comment in the correct answer thread below. This helps keep threads logically organized.Hello mms140130!
Questions similar to yours can already be found at:
We have closed your question to allow us to keep similar content in the same thread.
If you disagree with this please tell us why in a reply below. We'll be happy to talk about it.
Cheers!
I don't think it is the same thing