4
$\begingroup$

Okay so I'm using grep to try and get a preview of some trimming operations that are not going as expected.. Lets say that my sequence in the FastQ file is: ATNGCNATCG

What I want to do is..

grep "ATCGCTATCG" my.fastq

..and match the sequence given above

Surely there is some way or some existing tool that I can use besides doing: grep "[A|N][T|N][C|N][G|N][C|N]...etc."

$\endgroup$
3
  • 2
    $\begingroup$ Your fix is going to match NNNNN, ANNNC and etc. - are you sure this is what you want? $\endgroup$
    – Greg
    Mar 24, 2021 at 11:37
  • $\begingroup$ Technically, no it's not what I want, but as this is just trying to get a picture for what's going on under the hood I think it wont be a problem :) $\endgroup$
    – RPINerd
    Mar 24, 2021 at 14:06
  • $\begingroup$ Note that this will fail if the target sequence is split across multiple fastq lines. The fastq format does not require only one line of sequence per entry, you can have multiple ones. $\endgroup$
    – terdon
    Mar 30, 2021 at 17:20

1 Answer 1

5
$\begingroup$

If you want to stick to grep, use a scripting language such as Perl to generate the regex programmatically. For example:

perl -le 'print join "", map "[${_}N]", split //, $ARGV[0];' ATCGCTATCG 

Prints:

[AN][TN][CN][GN][CN][TN][AN][TN][CN][GN]

You can use it in grep like so:

grep '[AN][TN][CN][GN][CN][TN][AN][TN][CN][GN]' <<< ATNGCNATCG

Prints:

ATNGCNATCG

If that works for you, you could make it into a little bash function that also runs the grep. Add these lines to your ~/.bashrc:

grepN(){
    seq="$1"
    file="$2"
    pattern=$(perl -le 'print join "", map "[${_}N]", split //, $ARGV[0];' "$seq")
    grep "$pattern" "$file"
}

You can now run

grepN ATCGCTATCG my.fastq

Of course, this is not a good idea since the sequence might be in different lines, but that's what you were doing originally.

$\endgroup$
4
  • $\begingroup$ Awesome, I'll try this out. It will at least save me the time typing out all those N's hah. $\endgroup$
    – RPINerd
    Mar 24, 2021 at 14:16
  • 1
    $\begingroup$ It does involve one extra step, but saves me tons of time and gets me the data I need! Thanks! $\endgroup$
    – RPINerd
    Mar 24, 2021 at 14:42
  • 1
    $\begingroup$ @RPINerd I edited Timur's very clever solution to add a shell function that you can use to do the whole thing in one step. I hope you don't mind, Timur, please let me know if you'd rather I posted a separate answer. It just seemed like a trivial extension of yours and you did all the work so... $\endgroup$
    – terdon
    Mar 30, 2021 at 17:19
  • $\begingroup$ @terdon Your edits and suggestions are always welcome! Thank you. $\endgroup$ Mar 30, 2021 at 18:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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