Thank you once again! Your suggestions really helped!
Hello all, Where can I download all the chromosomal coordinates for the entire hg19 reference genome? Additionally, can anyone help me generate a python script that will split the whole genome into 50kbp windows? My draft python code is given below
step = 50000
with open('Chr.bed', 'w') as outfp:
for val in range(0, 300000000, step):
outfp.write('{0}\t{1}\t{2}\n'.format('chr', val, val+step))
I would appreciate any help and suggestions. Thank you
3 answers
The genome is available from e.g. UCSC and Ensembl and only a google search away plus there are plenty of questions here already on that matter, please use the search function. For the binning, take a look at bedtools makewindows. Not python but hey, why reinventing the wheel if a superfast option already exists.
You can use Kent tools and BEDOPS bedops --chop with standard Unix streams to do this quickly and efficiently:
$ fetchChromSizes hg19 | awk '{print $1"\t0\t"$2}' | sort-bed - | bedops --chop 50000 - > hg19.50k.bed
Note that the bounds of chromosomes will almost certainly not split evenly by 50kb units. In other words, chopping up a chromosome will leave some chunk at the end that is less than the desired length.
If you're doing calculations that involve aggregation, you may want to leave out that last < 50kb bin, by adding the -x operator to bedops --chop:
$ fetchChromSizes hg19 | awk '{print $1"\t0\t"$2}' | sort-bed - | bedops --chop 50000 -x - > hg19.50k.chunkless.bed
As a bonus, the output is sorted and ready to use for fast set operations with BEDOPS.
for hg19 chromosome sizes:
https://genome.ucsc.edu/goldenPath/help/hg19.chrom.sizes
for making windows of 10kb:
$ wget -np -r -nd https://genome.ucsc.edu/goldenPath/help/hg19.chrom.sizes
$ bedtools makewindows -g hg19.chrom.sizes -w 10000
try using pybedtools. https://daler.github.io/pybedtools/
Log in to answer this question.
Thank you for your suggestions! I will check out bedtools for making the genomic windows. Thank you once again
Thank you all for your expert insights! Bedtools and bedops are a better option than just writing my own python script. My purpose of generating these windows is to find the CNV frequency in each genomic intervals for 200,000 CNV coordinates in my dataset. Bedtools intersect function and -c flag seems a good option to find the overlap of CNV coordinates within those windows. The end goal is to see which interval or coordinates of genomes do not have any CNVs at all. And thank you once again guys for your expert insights. Really appreciate it.
If you want bins without CNVs:
That
awkbit can be adjusted to get bins with varied numbers of CNVs.Thank you so much for your help. Realy appreciate it. I will give it a try. Thanks once again