How Do I Create A Seqrecord In Biopython?

Question

[Jeremy Leipzig]

	id1 ="HWI-EAS380:8:1:16:830/1"
	seq1 ="AGGGCGTTCAGCAGCCAGCTTGCGGCAAAACTGCGTAACCGTCTTCTCGTTCTCT
	AAAAACCATTTTTCGTCCCCTTCGGGGCGGTGGTCTATAGTGTTATTAATATCAA
	GTTGGGGGAGCACATTGTAGCATTG"
	qual1="abbbbbaab`abaabbaabaaaab^E^``^aaabaa_\_abaaaaaaaa`aaaa`
	Z^`^^aaaaaaa`aa^aaa``_aa_aaaaaaaaaaa`aaaa`aaaaaabaaabba
	aaaaaaaaaaa_baaaabbbbbaba"

view raw biostars-966.py hosted with ❤ by GitHub

assume these are fastq-illumina quality scores and the sequence is unambiguous dna

I want to create a SeqRecord object from this. Thanks.

Answer 1

[Peter]

This answer is a bit of a cheat, but currently Biopython doesn't expose the FASTQ quality string decoding functionality in the public API:

from Bio import SeqIO
from StringIO import StringIO
fastq_string = "@%s\n%s\n+\n%s\n" % (id1, seq1, qual1)
record = SeqIO.read(StringIO(fastq_string), "fastq-illumina")

Why are you trying to do this?

[Jeremy Leipzig]

i am merging overlapping paired end sequences

Answer 2

[Yu]

>>> id1  ="HWI-EAS380:8:1:16:830/1"
>>> seq1 ="AGGGCGTTCAGCAGCCAGCTTGCGGCAAAACTGCGTAACCGTCTTCTCGTT"
>>> from Bio.SeqRecord import SeqRecord
>>> from Bio.Seq import Seq
>>> sr = SeqRecord(Seq(seq1), id1, '', '')
>>> sr
SeqRecord(seq=Seq('AGGGCGTTCAGCAGCCAGCTTGCGGCAAAACTGCGTAACCGTCTTCTCGTT', Alphabet()), id='HWI-EAS380:8:1:16:830/1', name='', description='', dbxrefs=[])

[Jeremy Leipzig]

i need the quality in the record the whole per-letter-annotation thing is what is giving me trouble