1
$\begingroup$

I have paired end reads and want to extract only those which have this sequence within them:

TGTATGTAAACTTCCGACTTCAACTGTA

I tried using grep -B1 -A2 "TGTATGTAAACTTCCGACTTCAACTGTA" file_1.fq | grep -v "^--$" and the same from the mate pair. After this I tried aligning them , but they dont align because the reads headers are different. Can someone help in extracting the mate pair of one read from the other file (2.fq), so that they can be aligned ? regards.

Improve this question
$\endgroup$
  • $\begingroup$ Can you sit down and write out exactly what you want to do? I think it's pretty clear why what you did is not returning matched pairs of reads, but it's not clear to me if you see why. $\endgroup$ – swbarnes2 Jan 14 at 0:55
  • $\begingroup$ As I wrote, that I have paired end reads. I want to extract reads from both pairs which have in inserted sequence TGTATGTAAACTTCCGACTTCAACTGTA in them. It belongs to a transposon along with guideDNA. The fastq files are stored as forward and reverse reads. When I executed the script it brings out the reads which have insertion within them. However, they dont align with Bowtie2 as the reads are not from the same pairs. $\endgroup$ – user1738234 Jan 14 at 8:33
  • $\begingroup$ What I'm trying to get at is, do you want reads which lack the sequence if their mate has it? Because if the answer is "no", why are you trying to align paired reads when you are throwing away a member of a pair? And if the answer is "yes", you are clearly throwing those reads away. $\endgroup$ – swbarnes2 Jan 14 at 21:37
2
$\begingroup$

Your question can be rewritten as, "how do I get the union of reads in two files?" Since you have already subset the mates, the steps are:

  1. Process each file, storing the list of reads names in each.
  2. Take the union of these.
  3. Iterate through the original/non-subset files and write only the reads with names in the union to the output.

The files will be both subset correctly and synced properly such that bowtie2 or any other aligner can correctly process them.

Improve this answer
$\endgroup$
  • $\begingroup$ Hi Devon....is there a tool which can do this ? Kindly let me know. $\endgroup$ – user1738234 Jan 14 at 12:47
  • $\begingroup$ Probably, but in the time it would take to search for it it'd be faster to write the 10-20 lines of python (or whatever you prefer) to just do it. $\endgroup$ – Devon Ryan♦ Jan 14 at 13:40
  • 1
    $\begingroup$ If you want to stick with shell scripting: comm -12 file1ids file2ds should give the union, followed by a similar grep like you used but with -f union.txt to give the patterns in a file. The undocumented --no-group-separator argument to grep may help you too if you're using GNU grep. $\endgroup$ – Jesse Jan 14 at 17:21
  • $\begingroup$ Alternately, you could go through both files simultaneously, and output both reads if either of them matches. $\endgroup$ – swbarnes2 Jan 14 at 21:38
  • $\begingroup$ Based on my reading of the question, I'd guess they're interested in the intersection of reads, rather than the union $\endgroup$ – gringer 2 days ago
0
$\begingroup$

It's a little bit unclear what you're asking, but if you want to repair paired-end files that have become disordered (or the linking has been lost), then the repair function of bbtools should work. This will pair up any coincident read IDs, and create a separate file for any reads that have no pairs in the other file:

https://jgi.doe.gov/data-and-tools/bbtools/bb-tools-user-guide/repair-guide/

repair.sh in1=broken1.fq in2=broken2 out1=fixed1.fq out2=fixed2.fq outs=singletons.fq repair

If grep is used on one file, and repair used on that file in combination with the original file from the matching read, then the relevant read pairs should appear in the output:

repair.sh in1=filtered_1.fq in2=file_2 out1=fixed_1.fq out2=fixed_2.fq outs=unpaired_2.fq repair
Improve this answer
$\endgroup$
  • $\begingroup$ It's not clear to me that both halves of any pairs will be caught by that grep, so there won't be anything repairable. $\endgroup$ – swbarnes2 Jan 15 at 5:07
  • $\begingroup$ That's the "unclear what you're asking" bit. $\endgroup$ – gringer 2 days ago

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.