Hello, I have 3000 WGS CRAM files and I want to split them into 1Mbp chunks. I want to split with exact genomic coordinate locations, e.g. starting from 1 to 1000000bp, 1000001bp to 2000000bp, 2000001bp to 3000000 etc. for all chromosomes. Therefore, each chunks have the similar corresponding region in each sample. Is there any way that I can do this?
this is untested, but a possible quick shell script approach
# make one megabase windows and convert to 'locstrings' that samtools can use
bedtools makewindows -g hg19.chrom.sizes -w 1000000 | awk -F '\t' '{printf("%s:%d-%s\n",$1,int($2)+1,$3);}' > out.txt
for i in *.cram; do
mkdir $i_folder # make subfolders, because otherwise you may exceed the amount of files allowed in a single directory on linux
while read p; do
samtools view -T ref.fa $i $p -o $i_folder/$i_$p.cram
done < out.txt;
done;
i would question whether this is truly necessary, but it's your call
a really minor point : bed is 0 based but samtools use 1-based coordinates. So it's should be something like:
bedtools makewindows -g hg19.chrom.sizes -w 1000000 | awk -F '\t' '{printf("%s:%d-%s\n",$1,int($2)+1,$3);}'
Log in to answer this question.
in addition to the good answers provided, might be useful/helpful to look at tabix. Basically, for a one time cost, tabix will index a file, making subsequent searches for a range of data ultra fast.