Bio.codonalign.codonalignment module

Code for dealing with Codon Alignment.

CodonAlignment class is inherited from MultipleSeqAlignment class. This is the core class to deal with codon alignment in biopython.

class Bio.codonalign.codonalignment.CodonAlignment(records='', name=None, alphabet=CodonAlphabet(Standard))

Bases: Bio.Align.MultipleSeqAlignment

Codon Alignment class that inherits from MultipleSeqAlignment.

>>> from Bio.Alphabet import generic_dna
>>> from Bio.SeqRecord import SeqRecord
>>> from Bio.Alphabet import IUPAC, Gapped
>>> a = SeqRecord(CodonSeq("AAAACGTCG", alphabet=default_codon_alphabet), id="Alpha")
>>> b = SeqRecord(CodonSeq("AAA---TCG", alphabet=default_codon_alphabet), id="Beta")
>>> c = SeqRecord(CodonSeq("AAAAGGTGG", alphabet=default_codon_alphabet), id="Gamma")
>>> print(CodonAlignment([a, b, c]))
CodonAlphabet(Standard) CodonAlignment with 3 rows and 9 columns (3 codons)
AAAACGTCG Alpha
AAA---TCG Beta
AAAAGGTGG Gamma
__init__(self, records='', name=None, alphabet=CodonAlphabet(Standard))

Initialize the class.

__str__(self)

Return a multi-line string summary of the alignment.

This output is indicated to be readable, but large alignment is shown truncated. A maximum of 20 rows (sequences) and 60 columns (20 codons) are shown, with the record identifiers. This should fit nicely on a single screen. e.g.

__getitem__(self, index, alphabet=None)

Return a CodonAlignment object for single indexing.

__add__(self, other)

Combine two codonalignments with the same number of rows by adding them.

The method also allows to combine a CodonAlignment object with a MultipleSeqAlignment object. The following rules apply:

  • CodonAlignment + CodonAlignment -> CodonAlignment

  • CodonAlignment + MultipleSeqAlignment -> MultipleSeqAlignment

get_aln_length(self)

Get aligment length.

toMultipleSeqAlignment(self)

Convert the CodonAlignment to a MultipleSeqAlignment.

Return a MultipleSeqAlignment containing all the SeqRecord in the CodonAlignment using Seq to store sequences

get_dn_ds_matrix(self, method='NG86', codon_table=NCBICodonTable(id=1, names=['Standard', 'SGC0'], ...))

Available methods include NG86, LWL85, YN00 and ML.

Argument:
  • method - Available methods include NG86, LWL85, YN00 and ML.

  • codon_table - Codon table to use for forward translation.

get_dn_ds_tree(self, dn_ds_method='NG86', tree_method='UPGMA', codon_table=NCBICodonTable(id=1, names=['Standard', 'SGC0'], ...))

Cnstruct dn tree and ds tree.

Argument:
  • dn_ds_method - Available methods include NG86, LWL85, YN00 and ML.

  • tree_method - Available methods include UPGMA and NJ.

classmethod from_msa(align, alphabet=CodonAlphabet(Standard))

Convert a MultipleSeqAlignment to CodonAlignment.

Function to convert a MultipleSeqAlignment to CodonAlignment. It is the user’s responsibility to ensure all the requirement needed by CodonAlignment is met.

Bio.codonalign.codonalignment.mktest(codon_alns, codon_table=NCBICodonTable(id=1, names=['Standard', 'SGC0'], ...), alpha=0.05)

McDonald-Kreitman test for neutrality.

Implement the McDonald-Kreitman test for neutrality (PMID: 1904993) This method counts changes rather than sites (http://mkt.uab.es/mkt/help_mkt.asp).

Arguments:
  • codon_alns - list of CodonAlignment to compare (each CodonAlignment object corresponds to gene sampled from a species)

Return the p-value of test result.