I have RNA-seq data for a specific cell line. I also have a specific sequence of interest. I map RNA-seq reads to the sequence of interest using bowtie and generate a .sam file:
First I build an index for my sequence:
bowtie-build campy.fa campy
I then map the reads:
bowtie -s -p 2 campy campy-pre-1m.fastq campy-pre-1m.sam
After generating a .sam file I sort it and then try to use cufflinks to quantity the expression level.
sort -k 3,3 -k 4,4n campy-pre-1m.sam > campy-pre-1m.sam.sorted
cufflinks -o campy-pre-1m.cufflinks campy-pre-1m.sam.sorted
Im getting the following error when I use cufflinks : Cufflinks requires that if your file has SQ records in the SAM header that they appear in the same order as the chromosomes names in the alignments. If there are no SQ records in the header, or if the header is missing, the alignments must be sorted lexicographically by chromsome name and by position.
Since Im not mapping my reads to a genome but only to a 4000 nucleotide sequence which is located on chromosome 12, is there a way to use cufflinks with out sorting since I don't really need to sort anything. Even if i tried to sort first it still doen't work. any help would be greatly appreciated
1 answer
How about if you try sorting your SAM file with samtools? This requires conversion of SAM to BAM, but conversion and sorting can be accomplished in one step via piping:
samtools view -uS campy-pre-1m.sam | samtools sort - campy-pre-1m
The command above converts SAM to BAM in the first step, taking SAM as input, sending uncompressed BAM to stdout. The second step takes the input on stdin (-), sorts it, and places it in a BAM file with the given name. Then you can run cufflinks on the bam file:
cufflinks -o campy-pre-1m.cufflinks campy-pre-1m.bam
With samtools, there are also other ways to view aspects of your alingment quickly without running cufflinks.
Log in to answer this question.