I need help to read fastq file on my server virtual machine im using the following python code
def readFastq(filename):
sequences = []
qualities = []
with open(filename) as fh:
while True:
fh.readline()
seq = fh.readline().rstrip()
fh.readline()
qual = fh.readline().rstrip()
if len(seq) == 0:
break
sequences.append(seq)
qualities.append(qual)
return sequences, qualities
seqs, quals = readFastq('sample-seqs.fastq')
print (seqs)
giving me the following
syntax error: invalid syntax (pointing to the line before the last line)
thank you
•
written
by
sarmademad •
I fully agree with the others, that unless there is a very good reason to write your own fastq parser, you should use Biopython etc.
However, for completeness, here's some working code. I've changed it quite a lot since there really isn't a good reason to use readline() since files are inherently iterable.
# Usage: python readfastq.py myseqs.fastq
import sys
from itertools import izip
def readFastq(filename):
headers = []
sequences = []
qualities = []
with open(filename, 'r') as ifh:
for header, sequence, sep, quality in izip(ifh, ifh, ifh, ifh):
headers.append(header.rstrip())
sequences.append(sequence.rstrip())
qualities.append(quality.rstrip())
return headers, sequences, qualities
h, s, q = readFastq(sys.argv[1])
# Do whatever with h, s and q....
There may be more efficient ways to go about it, but use of izip at least means this won't eat up all your memory if you ran it on large files.
This is elegant for parsing clean fastq file. And in python3, zip can be used to replace izip.
This parser assumes the user needs the entire file (h, s, q) in memory, while it's (I guess) more likely iterating over the fastq is more desirable.
Your code has a lot of issues associated with data structures. Please, see below code,
def readFastq(filename):
fastq_open = open(filename, "rU")
for line in fastq_open:
header1 = line.rstrip()
sequences = next(fastq_open).rstrip()
header2 = next(fastq_open).rstrip()
qualities = next(fastq_open).rstrip()
if len(seq) == 0:
break
return sequences, qualities
fastq_open.close()
seqs, quals = readFastq('sample-seqs.fastq')
print (seqs)
with is the preferred method for working with files, rather than "manually" opening and closing.
if you ever tried any test files, you should know this code is not working.
Do not try to write a code to parse the fastq file by yourself (unless you just want to practice), its a tricky one. Just use biopython like skbrimer mentioned above.
Log in to answer this question.
If it's pointing to the line before the last line, is your input filename/path correct?
You code runs for me, though the output only returns
[+], [+], [+]- so some of your.readline()logic is off too.