6
$\begingroup$

I'm currently trying to run bwa-mem on Influenza substrains using the following command:

~/bwa mem h5n1_1_cons.fa h5n1_1_read1.fq h5n1_1_read2.fq

h5n1_1_cons.fa is the consensus sequence for substrain h5n1_1, and the fq files are paired end reads 1 and 2. Running the program gives me this error:

[mem_sam_pe] paired reads have different names: "JQ714243-920-1", "JQ714243-920-2"

Visual inspection shows that the read names are identical until the last character. I have found a package (bbmap) with a utility, bbrename.sh, that renames the sequences, producing:

./bbrename.sh in=h5n1_1_read1.fq out=fixed prefix=h5n1_1

which yields the following output

-bash-4.2$ head fixh5n1_1_read1.fq
@h5n1_0
gaaaatgtctcaggagttgatgcatggtgtttaatctttggtttgcgctg
+
@4>C5?7FA?9@@@;1,-4B+E??4E9,DA>6CB)*2.)?7/.'65+75/
@h5n1_1
catcagttggattttgcctaaggatgtccaccatccttatcccaccaatt
+
7(9?)=C3AEF+@@F=C-E>EE57CD:;'?.*?7;>:3>2.;29..)=2'
@h5n1_2
catgcatttgttggaataatgttctcacaaatcaactgtattgaccgctg

This fixes the problem and allows BWA-MEM to work. However, header information is discarded. Is there another way to ensure the header names are the same while retaining all header information?

$\endgroup$
2
  • 1
    $\begingroup$ Please edit your question and show us a few examples of the read names in your original input files. The output of awk 'NR % 4 == 1{print $1}' h5n1_1_read1.fq | head and awk 'NR % 4 == 1{print $1}' h5n1_1_read2.fq | head for example. $\endgroup$
    – terdon
    Nov 29 '17 at 11:08
  • $\begingroup$ Since you mention bbmap already, why not use that for mapping your reads altogether? It is after all called bbmap $\endgroup$
    – holmrenser
    Nov 30 '17 at 7:47
7
$\begingroup$

Assuming the read name lines only contain things like @JQ714243-920-1 and nothing else, then sed -i '1~4 s/-[12]$//' file.fastq (N.B., example updated a bit for safety, thanks to @terdon) will remove the -1 and -2 suffixes.

However, you should try to figure out why the read names were munged into that form to begin with (at least Illumina machines won't produce that, I don't know about others).

$\endgroup$
1
  • $\begingroup$ You could also replace the the second - with a space. Most mapping tools only regard the part of the id before the first space as the read name. $\endgroup$ Nov 30 '17 at 10:57

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.