Trees | Indices | Help |
---|
|
object --+ | SeqRecord
A SeqRecord object holds a sequence and information about it.
Main attributes:
Additional attributes:
You will typically use Bio.SeqIO to read in sequences from files as SeqRecord objects. However, you may want to create your own SeqRecord objects directly (see the __init__ method for further details):
>>> from Bio.Seq import Seq >>> from Bio.SeqRecord import SeqRecord >>> from Bio.Alphabet import IUPAC >>> record = SeqRecord(Seq("MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF", ... IUPAC.protein), ... id="YP_025292.1", name="HokC", ... description="toxic membrane protein") >>> print record ID: YP_025292.1 Name: HokC Description: toxic membrane protein Number of features: 0 Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein())
If you want to save SeqRecord objects to a sequence file, use Bio.SeqIO for this. For the special case where you want the SeqRecord turned into a string in a particular file format there is a format method which uses Bio.SeqIO internally:
>>> print record.format("fasta") >YP_025292.1 toxic membrane protein MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF <BLANKLINE>
You can also do things like slicing a SeqRecord, checking its length, etc
>>> len(record) 44 >>> edited = record[:10] + record[11:] >>> print edited.seq MKQHKAMIVAIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF >>> print record.seq MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
|
|||
Inherited from |
|
|||
letter_annotations Dictionary of per-letter-annotation for the sequence. |
|||
seq The sequence itself, as a Seq or MutableSeq object. |
|||
Inherited from |
|
Add another sequence or string to this sequence. The other sequence can be a SeqRecord object, a Seq object (or similar, e.g. a MutableSeq) or a plain Python string. If you add a plain string or a Seq (like) object, the new SeqRecord will simply have this appended to the existing data. However, any per letter annotation will be lost: >>> from Bio import SeqIO >>> handle = open("Quality/solexa_faked.fastq", "rU") >>> record = SeqIO.read(handle, "fastq-solexa") >>> handle.close() >>> print record.id, record.seq slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN >>> print record.letter_annotations.keys() ['solexa_quality'] >>> new = record + "ACT" >>> print new.id, new.seq slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNNACT >>> print new.letter_annotations.keys() [] The new record will attempt to combine the annotation, but for any ambiguities (e.g. different names) it defaults to omitting that annotation. >>> from Bio import SeqIO >>> handle = open("GenBank/pBAD30.gb") >>> plasmid = SeqIO.read(handle, "gb") >>> handle.close() >>> print plasmid.id, len(plasmid) pBAD30 4923 Now let's cut the plasmid into two pieces, and join them back up the other way round (i.e. shift the starting point on this plasmid, have a look at the annotated features in the original file to see why this particular split point might make sense): >>> left = plasmid[:3765] >>> right = plasmid[3765:] >>> new = right + left >>> print new.id, len(new) pBAD30 4923 >>> str(new.seq) == str(right.seq + left.seq) True >>> len(new.features) == len(left.features) + len(right.features) True When we add the left and right SeqRecord objects, their annotation is all consistent, so it is all conserved in the new SeqRecord: >>> new.id == left.id == right.id == plasmid.id True >>> new.name == left.name == right.name == plasmid.name True >>> new.description == plasmid.description True >>> new.annotations == left.annotations == right.annotations True >>> new.letter_annotations == plasmid.letter_annotations True >>> new.dbxrefs == left.dbxrefs == right.dbxrefs True However, we should point out that when we sliced the SeqRecord, any annotations dictionary or dbxrefs list entries were lost. You can explicitly copy them like this: >>> new.annotations = plasmid.annotations.copy() >>> new.dbxrefs = plasmid.dbxrefs[:] |
Implements the 'in' keyword, searches the sequence. e.g. >>> from Bio import SeqIO >>> record = SeqIO.read(open("Fasta/sweetpea.nu"), "fasta") >>> "GAATTC" in record False >>> "AAA" in record True This essentially acts as a proxy for using "in" on the sequence: >>> "GAATTC" in record.seq False >>> "AAA" in record.seq True Note that you can also use Seq objects as the query, >>> from Bio.Seq import Seq >>> from Bio.Alphabet import generic_dna >>> Seq("AAA") in record True >>> Seq("AAA", generic_dna) in record True See also the Seq object's __contains__ method. |
Returns the record as a string in the specified file format. This method supports the python format() function added in Python 2.6/3.0. The format_spec should be a lower case string supported by Bio.SeqIO as an output file format. See also the SeqRecord's format() method.
|
Returns a sub-sequence or an individual letter. Slicing, e.g. my_record[5:10], returns a new SeqRecord for that sub-sequence with approriate annotation preserved. The name, id and description are kept. Any per-letter-annotations are sliced to match the requested sub-sequence. Unless a stride is used, all those features which fall fully within the subsequence are included (with their locations adjusted accordingly). However, the annotations dictionary and the dbxrefs list are not used for the new SeqRecord, as in general they may not apply to the subsequence. If you want to preserve them, you must explictly copy them to the new SeqRecord yourself. Using an integer index, e.g. my_record[5] is shorthand for extracting that letter from the sequence, my_record.seq[5]. For example, consider this short protein and its secondary structure as encoded by the PDB (e.g. H for alpha helices), plus a simple feature for its histidine self phosphorylation site: >>> from Bio.Seq import Seq >>> from Bio.SeqRecord import SeqRecord >>> from Bio.SeqFeature import SeqFeature, FeatureLocation >>> from Bio.Alphabet import IUPAC >>> rec = SeqRecord(Seq("MAAGVKQLADDRTLLMAGVSHDLRTPLTRIRLAT" ... "EMMSEQDGYLAESINKDIEECNAIIEQFIDYLR", ... IUPAC.protein), ... id="1JOY", name="EnvZ", ... description="Homodimeric domain of EnvZ from E. coli") >>> rec.letter_annotations["secondary_structure"] = " S SSSSSSHHHHHTTTHHHHHHHHHHHHHHHHHHHHHHTHHHHHHHHHHHHHHHHHHHHHTT " >>> rec.features.append(SeqFeature(FeatureLocation(20,21), ... type = "Site")) Now let's have a quick look at the full record, >>> print rec ID: 1JOY Name: EnvZ Description: Homodimeric domain of EnvZ from E. coli Number of features: 1 Per letter annotation for: secondary_structure Seq('MAAGVKQLADDRTLLMAGVSHDLRTPLTRIRLATEMMSEQDGYLAESINKDIEE...YLR', IUPACProtein()) >>> print rec.letter_annotations["secondary_structure"] S SSSSSSHHHHHTTTHHHHHHHHHHHHHHHHHHHHHHTHHHHHHHHHHHHHHHHHHHHHTT >>> print rec.features[0].location [20:21] Now let's take a sub sequence, here chosen as the first (fractured) alpha helix which includes the histidine phosphorylation site: >>> sub = rec[11:41] >>> print sub ID: 1JOY Name: EnvZ Description: Homodimeric domain of EnvZ from E. coli Number of features: 1 Per letter annotation for: secondary_structure Seq('RTLLMAGVSHDLRTPLTRIRLATEMMSEQD', IUPACProtein()) >>> print sub.letter_annotations["secondary_structure"] HHHHHTTTHHHHHHHHHHHHHHHHHHHHHH >>> print sub.features[0].location [9:10] You can also of course omit the start or end values, for example to get the first ten letters only: >>> print rec[:10] ID: 1JOY Name: EnvZ Description: Homodimeric domain of EnvZ from E. coli Number of features: 0 Per letter annotation for: secondary_structure Seq('MAAGVKQLAD', IUPACProtein()) Or for the last ten letters: >>> print rec[-10:] ID: 1JOY Name: EnvZ Description: Homodimeric domain of EnvZ from E. coli Number of features: 0 Per letter annotation for: secondary_structure Seq('IIEQFIDYLR', IUPACProtein()) If you omit both, then you get a copy of the original record (although lacking the annotations and dbxrefs): >>> print rec[:] ID: 1JOY Name: EnvZ Description: Homodimeric domain of EnvZ from E. coli Number of features: 1 Per letter annotation for: secondary_structure Seq('MAAGVKQLADDRTLLMAGVSHDLRTPLTRIRLATEMMSEQDGYLAESINKDIEE...YLR', IUPACProtein()) Finally, indexing with a simple integer is shorthand for pulling out that letter from the sequence directly: >>> rec[5] 'K' >>> rec.seq[5] 'K' |
Create a SeqRecord. Arguments:
You will typically use Bio.SeqIO to read in sequences from files as SeqRecord objects. However, you may want to create your own SeqRecord objects directly. Note that while an id is optional, we strongly recommend you supply a unique id string for each record. This is especially important if you wish to write your sequences to a file. If you don't have the actual sequence, but you do know its length, then using the UnknownSeq object from Bio.Seq is appropriate. You can create a 'blank' SeqRecord object, and then populate the attributes later.
|
Iterate over the letters in the sequence. For example, using Bio.SeqIO to read in a protein FASTA file: >>> from Bio import SeqIO >>> record = SeqIO.read(open("Fasta/loveliesbleeding.pro"),"fasta") >>> for amino in record: ... print amino ... if amino == "L" : break X A G L >>> print record.seq[3] L This is just a shortcut for iterating over the sequence directly: >>> for amino in record.seq: ... print amino ... if amino == "L" : break X A G L >>> print record.seq[3] L Note that this does not facilitate iteration together with any per-letter-annotation. However, you can achieve that using the python zip function on the record (or its sequence) and the relevant per-letter-annotation: >>> from Bio import SeqIO >>> rec = SeqIO.read(open("Quality/solexa_faked.fastq", "rU"), ... "fastq-solexa") >>> print rec.id, rec.seq slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN >>> print rec.letter_annotations.keys() ['solexa_quality'] >>> for nuc, qual in zip(rec,rec.letter_annotations["solexa_quality"]): ... if qual > 35: ... print nuc, qual A 40 C 39 G 38 T 37 A 36 You may agree that using zip(rec.seq, ...) is more explicit than using zip(rec, ...) as shown above. |
Returns the length of the sequence. For example, using Bio.SeqIO to read in a FASTA nucleotide file: >>> from Bio import SeqIO >>> record = SeqIO.read(open("Fasta/sweetpea.nu"),"fasta") >>> len(record) 309 >>> len(record.seq) 309 |
Returns True regardless of the length of the sequence. This behaviour is for backwards compatibility, since until the __len__ method was added, a SeqRecord always evaluated as True. Note that in comparison, a Seq object will evaluate to False if it has a zero length sequence. WARNING: The SeqRecord may in future evaluate to False when its sequence is of zero length (in order to better match the Seq object behaviour)! |
Add another sequence or string to this sequence (from the left). This method handles adding a Seq object (or similar, e.g. MutableSeq) or a plain Python string (on the left) to a SeqRecord (on the right). See the __add__ method for more details, but for example: >>> from Bio import SeqIO >>> handle = open("Quality/solexa_faked.fastq", "rU") >>> record = SeqIO.read(handle, "fastq-solexa") >>> handle.close() >>> print record.id, record.seq slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN >>> print record.letter_annotations.keys() ['solexa_quality'] >>> new = "ACT" + record >>> print new.id, new.seq slxa_0001_1_0001_01 ACTACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN >>> print new.letter_annotations.keys() [] |
A concise summary of the record for debugging (string). The python built in function repr works by calling the object's ___repr__ method. e.g. >>> from Bio.Seq import Seq >>> from Bio.SeqRecord import SeqRecord >>> from Bio.Alphabet import generic_protein >>> rec = SeqRecord(Seq("MASRGVNKVILVGNLGQDPEVRYMPNGGAVANITLATSESWRDKAT" ... +"GEMKEQTEWHRVVLFGKLAEVASEYLRKGSQVYIEGQLRTRKWTDQ" ... +"SGQDRYTTEVVVNVGGTMQMLGGRQGGGAPAGGNIGGGQPQGGWGQ" ... +"PQQPQGGNQFSGGAQSRPQQSAPAAPSNEPPMDFDDDIPF", ... generic_protein), ... id="NP_418483.1", name="b4059", ... description="ssDNA-binding protein", ... dbxrefs=["ASAP:13298", "GI:16131885", "GeneID:948570"]) >>> print repr(rec) SeqRecord(seq=Seq('MASRGVNKVILVGNLGQDPEVRYMPNGGAVANITLATSESWRDKATGEMKEQTE...IPF', ProteinAlphabet()), id='NP_418483.1', name='b4059', description='ssDNA-binding protein', dbxrefs=['ASAP:13298', 'GI:16131885', 'GeneID:948570']) At the python prompt you can also use this shorthand: >>> rec SeqRecord(seq=Seq('MASRGVNKVILVGNLGQDPEVRYMPNGGAVANITLATSESWRDKATGEMKEQTE...IPF', ProteinAlphabet()), id='NP_418483.1', name='b4059', description='ssDNA-binding protein', dbxrefs=['ASAP:13298', 'GI:16131885', 'GeneID:948570']) Note that long sequences are shown truncated. Also note that any annotations, letter_annotations and features are not shown (as they would lead to a very long string).
|
A human readable summary of the record and its annotation (string). The python built in function str works by calling the object's ___str__ method. e.g. >>> from Bio.Seq import Seq >>> from Bio.SeqRecord import SeqRecord >>> from Bio.Alphabet import IUPAC >>> record = SeqRecord(Seq("MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF", ... IUPAC.protein), ... id="YP_025292.1", name="HokC", ... description="toxic membrane protein, small") >>> print str(record) ID: YP_025292.1 Name: HokC Description: toxic membrane protein, small Number of features: 0 Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein()) In this example you don't actually need to call str explicity, as the print command does this automatically: >>> print record ID: YP_025292.1 Name: HokC Description: toxic membrane protein, small Number of features: 0 Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein()) Note that long sequences are shown truncated.
|
Returns the record as a string in the specified file format. The format should be a lower case string supported as an output format by Bio.SeqIO, which is used to turn the SeqRecord into a string. e.g. >>> from Bio.Seq import Seq >>> from Bio.SeqRecord import SeqRecord >>> from Bio.Alphabet import IUPAC >>> record = SeqRecord(Seq("MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF", ... IUPAC.protein), ... id="YP_025292.1", name="HokC", ... description="toxic membrane protein") >>> record.format("fasta") '>YP_025292.1 toxic membrane protein\nMKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\n' >>> print record.format("fasta") >YP_025292.1 toxic membrane protein MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF <BLANKLINE> The python print command automatically appends a new line, meaning in this example a blank line is shown. If you look at the string representation you can see there is a trailing new line (shown as slash n) which is important when writing to a file or if concatenating mutliple sequence strings together. Note that this method will NOT work on every possible file format supported by Bio.SeqIO (e.g. some are for multiple sequences only). |
|
letter_annotationsDictionary of per-letter-annotation for the sequence. For example, this can hold quality scores used in FASTQ or QUAL files. Consider this example using Bio.SeqIO to read in an example Solexa variant FASTQ file as a SeqRecord: >>> from Bio import SeqIO >>> handle = open("Quality/solexa_faked.fastq", "rU") >>> record = SeqIO.read(handle, "fastq-solexa") >>> handle.close() >>> print record.id, record.seq slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN >>> print record.letter_annotations.keys() ['solexa_quality'] >>> print record.letter_annotations["solexa_quality"] [40, 39, 38, 37, 36, 35, 34, 33, 32, 31, 30, 29, 28, 27, 26, 25, 24, 23, 22, 21, 20, 19, 18, 17, 16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5] The letter_annotations get sliced automatically if you slice the parent SeqRecord, for example taking the last ten bases: >>> sub_record = record[-10:] >>> print sub_record.id, sub_record.seq slxa_0001_1_0001_01 ACGTNNNNNN >>> print sub_record.letter_annotations["solexa_quality"] [4, 3, 2, 1, 0, -1, -2, -3, -4, -5] Any python sequence (i.e. list, tuple or string) can be recorded in the SeqRecord's letter_annotations dictionary as long as the length matches that of the SeqRecord's sequence. e.g. >>> len(sub_record.letter_annotations) 1 >>> sub_record.letter_annotations["dummy"] = "abcdefghij" >>> len(sub_record.letter_annotations) 2 You can delete entries from the letter_annotations dictionary as usual: >>> del sub_record.letter_annotations["solexa_quality"] >>> sub_record.letter_annotations {'dummy': 'abcdefghij'} You can completely clear the dictionary easily as follows: >>> sub_record.letter_annotations = {} >>> sub_record.letter_annotations {}
|
seqThe sequence itself, as a Seq or MutableSeq object.
|
Trees | Indices | Help |
---|
Generated by Epydoc 3.0.1 on Wed Jul 14 01:37:23 2010 | http://epydoc.sourceforge.net |