3
$\begingroup$

I have multiplexed pair-end fastq reads with dual barcodes. The issue is that one barcode is present in the header and one is present at the beginning of the read. I need a method to demultiplex this data, but in order to assign a read to an individual, both barcodes are required, as there is overlap between the barcodes. It seems there are packages available to demultiplex using header ID or in-line barcodes to demultiplex, but not both.

example reads:

@700819F:525:HT235BCXX:2:1101:1139:2144 1:N:0:ATCACGAT
CGAATTGCAGATTTTTTCTGAATAAAGCAGTGCAATAAAATTCCCCGCAAAAACACTTNANNNGNNNNNNNNNNNNNNNNNNNNANNNNNNGNTAATAAA
+
GGGGGIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII#<###.####################.######<#.<<GGGI
@700819F:525:HT235BCXX:2:1101:1212:2172 1:N:0:ATCACGAT
AAGGATGCAGGGCATCTCCCTCAGGCTGCGCTCTATCGAAGTCATCCCAGAATTAGATTCCGACCACAGACCAGTCTTAGTCAAACTAGGACCCGAGTGT
+
GGGGGIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIGIIIIIIIIIIIIIIIIIIIIGIIIIIIGGIIIGIIIIIIIIIIIIIGGIGIGIIIGIGI
@700819F:525:HT235BCXX:2:1101:1110:2173 1:N:0:ATCACGAT
GGTTGTGCAGAAAGAGTTGCTGATAAACTTAGCCATGCAGAACAGAATTATGAGTTAGAAGTATGTATATATATACCAATCACTATATCAACCCATTACC
+
<G.G<GGIIIG.GA.GAGAAGG<AA<A<.<<<GA<.<<.G<.G<<A<GGAA....<G.G..<<.GA..<A.<GG<<<.<..<<.GG.A..G..<<<.<<G

One of the barcodes is the sequence found at the end of the header line (ATCACGAT). The second barcode is found at the beginning of the sequence (CGAAT in read 1). Both are required for assignment to an individual.

$\endgroup$
3
  • $\begingroup$ Can you define where in the read the barcode should be found? $\endgroup$ – Bioathlete Jan 10 '18 at 20:38
  • $\begingroup$ One of the barcodes is the sequence found at the end of the header line (ATCACGAT). The second barcode is found at the beginning of the sequence (CGAAT in read 1). Both are required for assignment to an individual. This is ddRADSeq data, if that is important. $\endgroup$ – Caleb Benson Jan 10 '18 at 20:51
  • $\begingroup$ Ok that make sense. I didn't know the expected length of the barcode in the read. $\endgroup$ – Bioathlete Jan 10 '18 at 21:03
2
$\begingroup$

You could try something like that:

#!/bin/sh

fq_in="${1}"

bioawk -c fastx '{print substr($comment, 7, 15)}' ${fq_in} | uniq \
    > "barcode1.txt"
bioawk -c fastx '{print substr($seq, 1, 5)}' ${fq_in} | uniq \
    > "barcode2.txt"

for bc1 in $(cat "barcode1.txt")
do
    for bc2 in $(cat "barcode2.txt")
    do
        fq_out="${bc1}_${bc2}.fq"
        bioawk -c fastx -v bc1=${bc1} -v bc2=${bc2} \
            'substr($comment, 7, 15) == bc1 && substr($seq, 1, 5) == bc2 {print "@"$name" "$comment"\n"$seq"\n+\n"$qual}' \
            ${fq_in} \
            > ${fq_out}
    done
done

This might be slow. You can a little bit of time if you directly provide lists of barcodes for the for loops.

It could be more efficient to write a custom program that can write to multiple files. It happens that I wrote something that does this in python just yesterday, but only reading the inline barcode. Modifying it to also use the in-header barcode is likely simple.

$\endgroup$
1
$\begingroup$

Since you've found separate programs to demultiplex by Illumina barcode and in-line barcode (I presume this is actually just the sequence from the cut site), just run one and then the other on its output. Typically the Illumina barcode would have been used to demultiplex the files before you even received them (if not, your sequencing facility is incompetent), in which case you just need to demultiplex by the in-line barcode

$\endgroup$
1
$\begingroup$

The program extract from UMI-Tools will extract the barcode from the read and add it to the ID. It will also handle pair-end data so that even if the barcode is only on one of the reads, it will be added to the ID of both reads.

$\endgroup$

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.