Hi,
I am trying to search for the presence of couple sequences (around 400) each with a size of 23 bps,in different fastq files, while allowing 1-2 mismatches at maximum. I am not sure if turning the fastq to a genome(transcriptome) would be a nice approach? I have tried making the fastq -> fasta -> building blast database -> running blastn, however it did not run as my query is not only one sequence.
Example part of my query.file :
ATTTTTCTGAAAAACCCCCTACGA
AACAGGAAGTCAAAAAAAGCCAA
AGGATTTTTTTTTTTCTGGGGACA
The output I am aiming to have is, for each read in my query.file, which of these sequences are having 100% (or having 1-2 mismatches) match in fastq file, and possibly where in the fastq file.
I would appreciate your suggestions! Thank you!
1 answer
You can grab the fastq sequences containing these 23-mers with BBDuk like this:
bbduk.sh in=file.fastq outm=matched.fastq ref=23mers.fa k=23 hdist=2
"hdist=2" allows 2 mismatches; you can alteratively set that to 1 or 0. This does not tell you where the match is, but you can do that like this:
bbduk.sh in=matched.fastq out=masked.fastq ref=23mers.fa k=23 hdist=2 kmask=lc
That will convert the matched regions to lowercase.
Log in to answer this question.
-
Aligning very short sequences using bowtie2, I need fine control over how many mismatches are allow…
written by Mauro •
-
Multiple alignments with same score in Smith-Waterman
written by amansinha.sw •
-
Blast 2.11 not blasting all sequences in a file?
written by DNAngel
-
How to build 16S rRNA matrix?
written by buce_buper •
-
Mapping short reads (23 b) to large catalog of short (20 b) sequences
written by NaskoMak •
-
Searching for the common sequences in multiple (>2) fastq files
written by BPors •
-
What happens during seed stitching when STAR initially gets the wrong MMP?
written by eric.kern13
-
Read fragments receive maximum MAPQ in bowtie2
written by valerie
-
Error while performing MACS program
written by Mohit •
-
Bowtie Issues
written by Namibia •
Hi,
May be you can try to ta align with bwa aln your 23 bps seq against your fastq files as ref after you transformed it as fasta ?
Best
Thank you for your suggestion. Would this work if my reads are in text format?
You could use bowtie instead of blast. Make a fasta from the fastq, build a bowtie index from it, then align the query. Bowtie has an option that controls how many mismatches are allowed in the seed (-n). As the seed (28bp) is longer than your queries, setting the max seed mismatches to 1 or 2 should be sufficient for your goal.
Thank you for your answer. I would like to try, but I have these reads in just text format, therefore I cannot turn it to fastq. I think in Bowtie I have use reads in fastq format
No, several formats are accepted:
Thank you! I have eventually used BBDUK but I will give bowtie a try soon with these options. ( -r).
I was not aware of that these is a function in BB. This BB stuff is really a jack-of-all-trades.