Thank you very much! One other question. Before trim, I had sequence length = 100. After trimming, I had a sequence length of 3-100. Can you make any sense of this? I think I am having a hard time figuring out what length it is measuring.
Hello! I am new to RNA seq. I quality trimmed my fastq sequences via fastX toolkit using a phred score of 30. I would like to figure out the read length after the quality trim (phred score < 30 was removed). Any help would be much appreciated.
-Nikelle
2 answers
How about this:
sed -n '1~4p' filename.fastq | perl -ne 'chomp;print length($_) . "\n"' | sort -n | uniq -c >length.dist
Hi, Thanks Chris. If I were to input this, what would this be generating?
Let's break it down:
sed -n '1~4p' filename.fastq
Gives you every 4th line of the file (the sequence line)
perl -ne 'chomp;print length($_) . "\n"'
outputs the length of that line
sort -n | uniq -c
condenses it into a table of counts like this:
3 98
123 99
22 100
Thanks very much, that was helpful. What are the two different columns in the table?
Count and read length (see man uniq)
Chris, the fourth line of FASTQ is the quality score, not sequence (but it should be trimmed to the same length as the sequence string, so results should be the same).
Log in to answer this question.
This kind of trimming is very severe and likely unnecessary. I assume you are aligning to a reference genome. You are probably throwing away lot of good data if you did lose a lot if bases after trimming.
See this thread for a paper referenced in there about trimming using quality: Which Phred value to use in trimming
You might want to consult www.rnaseq.wiki for a really nice description of the steps involved in processing RNAseq. IIRC, they cover read trimming in the appropriate section