0

I have to retrieve a specific information with known length from my fastq file and add it in another position.
For example given as input the following fastq file:

@SRR5394526.1 1 length=150  
CGATGTTAAATCAACGATAACTACACCG    
+SRR5394526.1 1 length=150  
AA<AFJFJJJJJJJJJJAJJJJJJJJJF    

I would like as output:

@SRR5394526.1.CGATGT 1 length=150    
TAAATCAACGATAACTACACCG    
+SRR5394526.1.CGATGT 1 length=150    
FJJJJJJJJJJAJJJJJJJJJF  

As you can notice the first 6 nucleotides were removed from either the sequence in the second line but also the sequence in the 4th line and appended after the first number 1 in the 1st and the 3rd lines. I have millions of chunk of this size (4 rows) within the file and this is just an example.

I have already found how to add/append info in a file: sed 's/myinfo/&,/4' and how to remove info in a file sed -e '423s!//!!; 424s!printf!//&!' but it's not enough. Any idea is much appreciated.

3
  • 1
    Just checking to see if you are aware of the StackExchange Bioinformatics site at bioinformatics.stackexchange.com ?
    – Kusalananda
    Apr 2, 2019 at 22:51
  • oh, I was not. Should I move my question it there? Apr 2, 2019 at 22:53
  • 1
    No need for this particular one I think.
    – Kusalananda
    Apr 2, 2019 at 22:55

2 Answers 2

3

Using awk:

awk '(FNR-1) % 2 == 0 { name=$1; chr=$2; len=$3; next }
     (FNR-2) % 4 == 0 { seq=substr($0,1,6) }
                      { print name "." seq, chr, len
                        print substr($0,7) }' file.fastq >newfile.fastq

This awk program has three blocks.

  1. The first block triggers on every second row (the sequence and the quality data header lines), starting on the first one. It saves the three bits of information on that row into three variables. It then skips to the next line of input immediately.

  2. The second block extracts the 6 first characters from the sequence line into seq, but only for ever fourth row starting at row 2 (the sequence lines only).

  3. The last block runs only on rows not handled by the first block (every sequence or quality data line) and constructs the output.

To use this on gzip-compressed files (or bgzip-compressed files, commonly used in bioinformatics projects), use

zcat file.fastq.gz | awk '...' | bgzip -c >newfile.gz

To use a variable as the value used for cutting, consider

awk -v n=6 '(FNR-1) % 2 == 0 { name=$1; chr=$2; len=$3; next }
            (FNR-2) % 4 == 0 { seq=substr($0,1,n) }
                             { print name "." seq, chr, len
                               print substr($0,n+1) }'

Where the -v n=6 controls the length of the cut.

You may also put the actual awk code (everything inside the single quotes) into its own script file and use it as

awk -v n=6 -f script.awk file.fastq
0
0

data in fastq file, use gnu sed on 4 by 4 lines,

$ sed -nE ' N;N;N;s/(.+\.1)(\s.+\n)(.{6})(\w+)\s*(\n.+\.1)(.+\n).{6}(\w+)/\1.\3\2\4\5.\3\6\7/p' fastq

@SRR5394526.1.CGATGT 1 length=150
TAAATCAACGATAACTACACCG
+SRR5394526.1.CGATGT 1 length=150
FJJJJJJJJJJAJJJJJJJJJF

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .