I have downloaded a dataset (which documents nucleosome postioning in hESCs) that took almost 4 hours to transfer to my system. And I am using the following script to convert the data to a bed file:
curl ftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByExp/sra/SRX/SRX645/SRX645331/SRR1781834/SRR1781834.sra>SRR1781834.sra
fastq-dump ./SRR1781834.sra
# in the case of single-end reads, there are two input files for mapping:
bowtie -t -v 2 -p 8 -m 1 --solexa-quals mm9 SRR1781834.fastq SRR1781834.map
# in case you have used Boutie at the previous step, now you have to convert the file fom .map to .bed
perl -w /storage/projects/teif/scripts/bowtie2bed.pl SRR1781834.map SRR1781834.bed
I was told that this was a fairly logical way to convert sra files to bed files and yet the output in the final bed file is as follows:
chr1 16972827 16972927 U2 0 +
chr1 169261752 169261852 U2 0 -
chr10 14618525 14618625 U2 0 +
chr10 32927157 32927257 U2 0 -
chr10 36893563 36893663 U2 0 +
chr10 71759925 71760025 U2 0 -
chr10 71759929 71760029 U2 0 +
NOTE: the start and end position only documents the analysis of single nucleotide positions (this dataset should be documenting entire portions of chromosomes)....
Below is the documented descent in the size of the files as the script progresses (notice how the file size decreases with each new ouput):
-rw-r--r-- 1 15G Dec 9 16:49 SRR1781834.sra
-rw-r--r-- 1 73G Dec 9 17:03 SRR1781834.fastq
-rw-r--r-- 1 18K Dec 9 18:02 SRR1781834.map
-rw-r--r-- 1 1.9K Dec 11 14:19 SRR1781834.bed
-rw-r--r-- 1 1.9K Dec 11 14:25 SRR1781834_sorted.bed
My question is how is the data of such a large SRA file translating into such a small bed file?..... Am I doing this wrong and how can I improve my file conversion strategy to get more informative data?
0 answers
No answers yet.
Log in to answer this question.
Beside the issue with the .bed file, the .map file is also surprisingly small. What is the percentage of mapped reads ?
You are right:
reads processed: 238967267
reads with at least one reported alignment: 62 (0.00%)
reads that failed to align: 238966752 (100.00%)
reads with alignments suppressed due to -m: 453 (0.00%)
Reported 62 alignments to 1 output stream(s) Time searching: 00:56:39 Overall time: 00:56:39
Is there any way that you know of to get around this or is it just a faulty dataset?
Well first question : are u aligning to the correct genome ? I think mm9 stands for M.musculus (mouse) so that may explain why your alignment rate is so low. If that doesn't fix your problem, have a look at a few lines of your .fastq file. The file might be corrupted. You can also run fastqc (quality control) on your .fastq file to check the read quality, length, etc...
thank you will try this