2

I am handling fastq files from next generation sequencing data that are organized as follow:

  • Row 1 information about the sequencing machine, the lane, the "tiles" and the number of the read
  • Row 2 information about the sequence
  • Row 3 is a symbol + used as a separator
  • Row 4 information about the quality of the reads

I would like to remove all the reads that are in the interval of specific numbers in the 5th position of the first row (in black).

Here an example of reads to remove between 1101 to 1103. Input:

@ST-E00204:114:HHKTJALXX:4:1101:22962:1538_1:N:0:1/1
NGTGTTTTTAATTATTAAGTTATTTTTTTAGTTTTTTAAGGATTTTTATAGTAGTAATAGAAATTTAATTAAGATAGAAAATTTTAAGTGTGGTTAGGATTGTAGTTTTGTTGGTATTATGTTGATTTAGTATAAGTAAAGTTTTGATTTT
+
AAAAJJJJFJJAJJJJAJAJJJJJJJJJAJJ-FJJJJF--FJJJJFJJJFFJJJFFJ-JJJJFFFFJ-AJ7AJJJJJJJJJJFJJJJFJFAJFFJJF-AAFAJFJJ7AJAJJFJFJJJ7FFFFFJFJJ-7F-77AJF--7FJJ


@ST-E00204:114:HHKTJALXX:4:1102:7101:2012 2:N:0:1
NATTTAAAAATACCCACTATAAAACATAAAATATAACAAAAAAACTAAAAATCATAAAAAATAAAAAAAATCCACTTCACGTCTTTTAACAATTTCGTCATTTTTAACATCCTCAAATAAATTATTCTCATTTTCCATAACTTCCAATTTT
+
!A-A-FJJJJJJ-FJAJFJJJJJFJJAJJJJJJ-F-AJJJJJJ-F-FJJFJJFJFFFFF<-F<FJJJF-<7<JF<-7AAFFJ--A<A77--7FAAF-A----7FF-7-7<F-J<A-7--<F7---77<----7-<<FA---7<<---7---

@ST-E00204:114:HHKTJALXX:4:1103:7141:2012 2:N:0:1
NAAAACATAAAATATAACAAACAAACTAAAAATCATAAAAAATAAAAAACATCCACTTAACAACTTAAAAAATAACAAAATCACTAATTATAATAAAAAATAAAAAATACACACTCTAACACCTAAAACACAACCAAAAAACTAAAACTCC
+
!AAFFFFJJJJA-F--AFFJJ-F<JJF<AJFJ<JF-7<JJAA7-J-FFFJ7JJJFJ-F<AJJJJFFJ-A-F-AJ<FF-JFFF-77<JJ---777<7----7-A<J-A-7<<FFF<--7--7-FFFF-<---7---7A-<A7FA------7-

@ST-E00204:114:HHKTJALXX:4:1104:7101:2012 2:N:0:1
NATTTAAAAATACCCACTATAAAACATAAAATATAACAAAAAAACTAAAAATCATAAAAAATAAAAAAAATCCACTTCACGTCTTTTAACAATTTCGTCATTTTTAACATCCTCAAATAAATTATTCTCATTTTCCATAACTTCCAATTTT
+
!A-A-FJJJJJJ-FJAJFJJJJJFJJAJJJJJJ-F-AJJJJJJ-F-FJJFJJFJFFFFF<-F<FJJJF-<7<JF<-7AAFFJ--A<A77--7FAAF-A----7FF-7-7<F-J<A-7--<F7---77<----7-<<FA---7<<---7---

Desired output:

@ST-E00204:114:HHKTJALXX:4:1104:7101:2012 2:N:0:1
NATTTAAAAATACCCACTATAAAACATAAAATATAACAAAAAAACTAAAAATCATAAAAAATAAAAAAAATCCACTTCACGTCTTTTAACAATTTCGTCATTTTTAACATCCTCAAATAAATTATTCTCATTTTCCATAACTTCCAATTTT
+
!A-A-FJJJJJJ-FJAJFJJJJJFJJAJJJJJJ-F-AJJJJJJ-F-FJJFJJFJFFFFF<-F<FJJJF-<7<JF<-7AAFFJ--A<A77--7FAAF-A----7FF-7-7<F-J<A-7--<F7---77<----7-<<FA---7<<---7---

One Idea was to use:

split -l 4 myfile.fq

and then remove each file according the numbers in the 5th position like:

grep -v ":1104"
grep -v ":1105"

etc. but the problem is that the file is very big. Also I have to remove big intervals like from 1000 to 2000 and each number correspond to a lot of reads.

1
  • no they are not separated by the blank lines but thanks for the suggestion. – fusion.slope Jun 16 '17 at 9:56
2

awk solution:

awk -F':' -v RS="@" 'NR>1 && ($5<1101 || $5>1103){ print RS$0 }' myfile.fq

The output:

@ST-E00204:114:HHKTJALXX:4:1104:7101:2012 2:N:0:1
NATTTAAAAATACCCACTATAAAACATAAAATATAACAAAAAAACTAAAAATCATAAAAAATAAAAAAAATCCACTTCACGTCTTTTAACAATTTCGTCATTTTTAACATCCTCAAATAAATTATTCTCATTTTCCATAACTTCCAATTTT
+
!A-A-FJJJJJJ-FJAJFJJJJJFJJAJJJJJJ-F-AJJJJJJ-F-FJJFJJFJFFFFF<-F

Details:

  • -F':' - field separator :

  • -v RS="@" - considering @ as record separator

  • ($5<1101 || $5>1103) - check if the needed field fits the condition "remove between 1101 to 1103"

1

Here's an awk solution that assumes that every block is exactly 4 lines long.

awk -F ':' '
    NR % 4 == 1 { keep = ($5 < 1101 || $5 > 1103) }
    keep
'

Explanations:

  • -F : sets the field separator to :. $5 is the fifth field on the current line (numbered from 1).
  • NR % 4 == 1 { … } executes the code in the braces if the condition NR % 4 == 1 is true, i.e. if the current line number is a multiple of 4 plus 1.
  • keep = (…) sets the variable keep to the specified condition.
  • The line with just keep says to print the current line if the variable keep is true. Since keep is only set every 4 lines, it applies identically to all the lines within a 4-line block.
2
  • o nice also this solution with "keep". Also if i would like to select multiple intervals of numbers, should i do: { keep = ($5 < 1101 || $5 > 1103) ($5 < 1600 || $5 > 1500)($5 < 2100 || $5 > 2000)} ? – fusion.slope Jun 17 '17 at 13:15
  • 1
    @fusion.slope No, that doesn't work. It's only syntactically correct by coincidence, because (a) (b) is the string a concatenated with b and booleans are automatically converted to a string "0" or "1". To keep everything except for the intervals 1103–1105, 1500–1600 and 2000–2100, use keep = ! ((1103 <= $5 && $5 <= 1105) || (1500 <= $5 && $5 <= 1600) || (2000 <= $5 && $5 <= 2100)) – Gilles 'SO- stop being evil' Jun 17 '17 at 15:41
0
sed -e '
   /@/s/:/\n/4; /@/s/:/\n/4
   /\n110[1-3]\n/,/^$/d
   y/\n/:/
' fastq.data

  1. Mark the serial number on lines that begin with @.
  2. Delete the chunk which holds the serial number 1101 or 1102 or 1103.
  3. On all other chunks take away the markers and let it go to stdout.

Result

@ST-E00204:114:HHKTJALXX:4:1104:7101:2012 2:N:0:1
NATTTAAAAATACCCACTATAAAACATAAAATATAACAAAAAAACTAAAAATCATAAAAAATAAAAAAAATCCACTTCACGTCTTTTAACAATTTCGTCATTTTTAACATCCTCAAATAAATTATTCTCATTTTCCATAACTTCCAATTTT
+
!A-A-FJJJJJJ-FJAJFJJJJJFJJAJJJJJJ-F-AJJJJJJ-F-FJJFJJFJFFFFF<-F

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.