cmcpy Package

cmcpy Package

__main__ Module

amino_acid_spaces Module

class amino_acid_spaces.RegionAminoAcidSpace(num_aas=None, coords=None, num_dims=1, seed=42, labels=None)[source]

Bases: amino_acid_spaces._AminoAcidSpace

Region amino acid spaces model amino acid (dis)similarities in bounded regions of a finite number of dimensions.

>>> aa = RegionAminoAcidSpace(num_aas = 5,num_dims = 2)
>>> map(lambda x:x.round(2),aa.coords)
[array([[ 0.16,  0.71]]), array([[ 0.37,  0.16]]), array([[ 0.6,  0.6]]), array([[ 0.73,  0.87]]), array([[ 0.95,  0.06]])]
>>> dm =  aa.get_distance_matrix()
>>> dm.round(3)
array([[ 0.   ,  0.594,  0.455,  0.597,  1.027],
       [ 0.594,  0.   ,  0.498,  0.795,  0.584],
       [ 0.455,  0.498,  0.   ,  0.297,  0.647],
       [ 0.597,  0.795,  0.297,  0.   ,  0.837],
       [ 1.027,  0.584,  0.647,  0.837,  0.   ]])
class amino_acid_spaces.RingAminoAcidSpace(num_aas=None, seed=42, coords=None, labels=None)[source]

Bases: amino_acid_spaces._AminoAcidSpace

Ring amino acid spaces model amino acid (dis)similarities in a one-dimensional circular physicochemical amino acid space

>>> aa = RingAminoAcidSpace(num_aas = 5)
>>> map(lambda x: x.round(3),aa.coords)
[array([[ 0.156]]), array([[ 0.375]]), array([[ 0.599]]), array([[ 0.732]]), array([[ 0.951]])]
>>> dm = aa.get_distance_matrix()
>>> dm.round(3)
array([[ 0.   ,  0.219,  0.443,  0.424,  0.205],
       [ 0.219,  0.   ,  0.224,  0.357,  0.424],
       [ 0.443,  0.224,  0.   ,  0.133,  0.352],
       [ 0.424,  0.357,  0.133,  0.   ,  0.219],
       [ 0.205,  0.424,  0.352,  0.219,  0.   ]])

codon_spaces Module

class codon_spaces.RingCodonSpace(num_codons, mu)[source]

Bases: codon_spaces._CodonSpace

Ring codon spaces are wrapped linear mutation spaces where codons mutate only to their two immediate neighbors. For ring codon models, mu defines the probability of change to one of two codon neighbors. The probability of no change is [1 - (2*mu)]

>>> codons = RingCodonSpace(num_codons = 5,mu = 0.1)
>>> mm = codons.get_mutation_matrix()
>>> mm.round(3)
array([[ 0.9 ,  0.05,  0.  ,  0.  ,  0.05],
       [ 0.05,  0.9 ,  0.05,  0.  ,  0.  ],
       [ 0.  ,  0.05,  0.9 ,  0.05,  0.  ],
       [ 0.  ,  0.  ,  0.05,  0.9 ,  0.05],
       [ 0.05,  0.  ,  0.  ,  0.05,  0.9 ]])
get_derivative_matrix()[source]
post_process_perturbative_solution(lpert, vpert)[source]
class codon_spaces.WordCodonSpace(num_bases, num_positions, mu, kappa=1.0)[source]

Bases: codon_spaces._CodonSpace

Word codon spaces model natural codons with a finite number of bases and a finite word-length.

For word codon models, mu defines the total probability of change of a base to any neighbor. The probability of no change of a single base is defined as (1 - mu).

If kappa is not equal to 1.0, then num_bases must be even.

>>> codons = WordCodonSpace(num_bases = 2,num_positions = 2, mu = 0.1)
>>> mm = codons.get_mutation_matrix()
>>> mm.round(3)
array([[ 0.81 ,  0.045,  0.045,  0.003],
       [ 0.045,  0.81 ,  0.003,  0.045],
       [ 0.045,  0.003,  0.81 ,  0.045],
       [ 0.003,  0.045,  0.045,  0.81 ]])
>>> codons = WordCodonSpace(num_bases = 4,num_positions = 2, mu = 0.2,kappa = 2)
>>> mm = codons.get_mutation_matrix()
>>> mm.round(3)
get_derivative_matrix()[source]

This function exists to serve the perturbative solution in evolvers.py

post_process_perturbative_solution(lpert, vpert)[source]

evolvers Module

evolvers – cmcpy module for abstract base class of Ardell Sella Evolvers

class evolvers.ArdellSellaEvolverAbstractBase(initial_code, site_types, delta, epsilon, observables)[source]

Bases: object

Abstract Base Class for ArdellSellaEvolvers for Code-Message Coevolution corresponding to models published in Ardell and Sella (2001, 2002) and Sella and Ardell (2002, 2006). Concrete Implementations subclass from this for different implementations of Eigenvalue solutions.

compute_code_fitness_given_messages(equilibrated_messages, effective_code_matrix)[source]

Implement eg eqns. 2-7 from Sella and Ardell(2006)

compute_max_fitness_code_mutation()[source]
equilibrate_messages()[source]

Compute eigensystems in site-types for an established genetic code.

This finds eigensystems (codon frequencies and growth rates) for different site-types given the genetic code.

Abstract method: subclasses must:

1) store their results by setting self.eigenvalues and self.eigenmatrix

  1. set self.equilibrated to True

3) call super(<<SubClass>>, self).equilibrate_messages() to print observables at end of subclass method where <<SubClass>> is the subclass name to print observables.

evolve_code_unless_frozen()[source]

Mutate genetic code.

Unless genetic code is frozen, mutate it according to the Ardell and Sella models, and update code to most fit mutant if it exists. If no more fit mutant code exists, set the “frozen” attribute to True.

evolve_one_step()[source]

Mutate genetic code and equilibrate messages.

Unless genetic code is frozen, mutate it according to the Ardell and Sella models, update code to most fit mutant if it exists, and equilibrate messages to the new mutant genetic code. If no more fit mutant code exists, set the “frozen” attribute to True.

evolve_until_frozen()[source]

Iteratively evlove genetic code and messages.

Until genetic code is frozen, mutate it according to the Ardell and Sella models, update code to most fit mutant if it exists, and equilibrate messages to the new mutant genetic code. Once no more fit mutant code exists, set the “frozen” attribute to True.

get_eigenvalue(msm, eigenvec)[source]
get_mutation_selection_matrix(alpha)[source]
get_selection_matrix(alpha)[source]
growth_rate()[source]
growth_rate_from_lambda()[source]
initial_equilibrate_messages()[source]

Compute eigensystems in site-types for an established genetic code.

This finds eigensystems (codon frequencies and growth rates) for different site-types given the genetic code.

messages()[source]
print_initial_observables()[source]
print_observables()[source]
print_observables_header()[source]
class evolvers.ArdellSellaEvolverNumpy(initial_code, site_types, num_processes, observables, delta=1e-32, epsilon=1e-12)[source]

Bases: evolvers.ArdellSellaEvolverAbstractBase

compute_eigensystem(alpha, max_time=60, numpy_type=<type 'numpy.float64'>)[source]
equilibrate_messages()[source]
class evolvers.ArdellSellaEvolverNumpyMulticore(initial_code, site_types, num_processes, observables, delta=1e-32, epsilon=1e-12)[source]

Bases: evolvers.ArdellSellaEvolverAbstractBase

equilibrate_messages()[source]
class evolvers.ArdellSellaEvolverNumpyProcessChild(in_queue, out_queue)[source]

Bases: multiprocessing.process.Process

Finds the eigensystem (message equilibrium and growth rate) for an established genetic code

run()[source]
class evolvers.ArdellSellaEvolverPerturbative(initial_code, site_types, num_processes, observables, delta=1e-32, epsilon=1e-12)[source]

Bases: evolvers.ArdellSellaEvolverAbstractBase

choose(n, k)[source]

A fast way to calculate binomial coefficients by Andrew Dalke (contrib).

compute_eigensystem(alpha, max_time=60, numpy_type=<type 'numpy.float64'>)[source]
equilibrate_messages()[source]
class evolvers.ArdellSellaEvolverPowerCUDA(initial_code, site_types, num_processes, observables, delta=1e-32, epsilon=1e-12)[source]

Bases: evolvers.ArdellSellaEvolverAbstractBase

equilibrate_messages()[source]
class evolvers.ArdellSellaEvolverPowerMethod(initial_code, site_types, num_processes, observables, delta=1e-32, epsilon=1e-12)[source]

Bases: evolvers.ArdellSellaEvolverAbstractBase

compute_eigensystem(alpha, max_time=60, numpy_type=<type 'numpy.float64'>)[source]
equilibrate_messages()[source]
class evolvers.ArdellSellaEvolverPowerMethodProcessChild(in_queue, out_queue, num_codons, delta, max_time=60)[source]

Bases: multiprocessing.process.Process

Finds the eigensystem (message equilibrium and growth rate) for an established genetic code

run()[source]
class evolvers.ArdellSellaEvolverPowerMulticore(initial_code, site_types, num_processes, observables, delta=1e-32, epsilon=1e-12)[source]

Bases: evolvers.ArdellSellaEvolverAbstractBase

equilibrate_messages()[source]
in_queue = <multiprocessing.queues.JoinableQueue object at 0x104eb6290>
out_queue = <multiprocessing.queues.JoinableQueue object at 0x104e42510>
class evolvers.eigensystem_CUDA_implementation(parent, max_time=60, delta=1e-32)[source]
calculate()[source]
done()[source]
error_check()[source]
get_eigenmatrix()[source]
get_eigenvalue(alpha)[source]
get_eigenvalues()[source]

genetic_codes Module

class genetic_codes.GeneticCodeMutation(code, codon, aa)[source]
get_effective_code_matrix()[source]
class genetic_codes.InitiallyAmbiguousGeneticCode(codons, amino_acids, misreading=None)[source]

Bases: genetic_codes._GeneticCode

class genetic_codes.UserInitializedGeneticCode(codons, amino_acids, code_matrix=None, code_dict=None, misreading=None)[source]

Bases: genetic_codes._GeneticCode

User-Initialized Genetic Codes are initialized with a numpy.ndarray code matrix or a dict of codons mapping to indices (not labels) of amino acids

>>> codons = codon_spaces.WordCodonSpace(num_bases = 4,num_positions = 2, mu = 0.2,kappa = 2)
>>> aas    = amino_acid_spaces.RegionAminoAcidSpace(num_aas = 20, seed = 40)
>>> cm = numpy.eye(16)
>>> cm = numpy.hstack((cm,numpy.zeros((16,4))))
>>> cm.shape
(16, 20)
>>> cm[0][1] = 1
>>> cm /= cm.sum(axis = 1).reshape(16,1)
>>> gc = UserInitializedGeneticCode(codons,aas,code_matrix = cm)
>>> gc.num_codons
16
>>> gc.num_amino_acids
20
>>> gc.ambiguous_codons()
set([0])
>>> gc.encoded_aas
set([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15])
>>> print gc
|* b c d|
|e f g h|
|i j k l|
|m n o p|
>>> gc.as_labelled_dict()
{0: '*', 1: 'b', 2: 'c', 3: 'd', 4: 'e', 5: 'f', 6: 'g', 7: 'h', 8: 'i', 9: 'j', 10: 'k', 11: 'l', 12: 'm', 13: 'n', 14: 'o', 15: 'p'}
>>> gc.as_dict()
{0: '*', 1: 1, 2: 2, 3: 3, 4: 4, 5: 5, 6: 6, 7: 7, 8: 8, 9: 9, 10: 10, 11: 11, 12: 12, 13: 13, 14: 14, 15: 15}
>>> gc2 = UserInitializedGeneticCode(codons,aas,code_dict = {0: '*', 1: 1, 2: 2, 3: 3, 4: 4, 5: 5, 6: 6, 7: 7, 8: 8, 9: 9, 10: 10, 11: 11, 12: 12, 13: 13, 14: 14, 15: 15})
>>> print gc2
|* b c d|
|e f g h|
|i j k l|
|m n o p|
>>> print gc2.code_matrix[0]
[ 0.05  0.05  0.05  0.05  0.05  0.05  0.05  0.05  0.05  0.05  0.05  0.05
  0.05  0.05  0.05  0.05  0.05  0.05  0.05  0.05]

misreading Module

class misreading.PositionalMisreading(codons, misreading)[source]

Bases: misreading._Misreading

Positional misreading models misreading on word codon spaces which model natural codons with a finite number of bases and a finite word-length.

For positional misreading, the misreading parameter is a list of positional misreading parameters mr_i which define the total misreading probability of a base at position i to any neighbor. The probability of no misreading of a single base at position i is defined as (1 - mr_i).

>>> codons = codon_spaces.WordCodonSpace(num_bases = 2,num_positions = 2, mu = 0.1)
>>> misreading = PositionalMisreading(codons,[0.1,0.01])
>>> mr = misreading.get_misreading_matrix()
>>> mr.round(3)
array([[ 0.891,  0.099,  0.009,  0.001],
       [ 0.099,  0.891,  0.001,  0.009],
       [ 0.009,  0.001,  0.891,  0.099],
       [ 0.001,  0.009,  0.099,  0.891]])    
class misreading.RingMisreading(codons, misreading)[source]

Bases: misreading._Misreading

Ring misreading is one-dimensional misreading uniform over all other codons. For ring misreading, (mr/(nc - 1)) is the probability of misreading as a specific codon. The probability of no misreading is (1 - (mr)).

The misreading parameter is a list with one element, mr.

>>> codons = codon_spaces.RingCodonSpace(num_codons = 5,mu = 0.1)
>>> misreading = RingMisreading(codons,[0.1])
>>> mr = misreading.get_misreading_matrix()
>>> mr.round(3)
array([[ 0.9  ,  0.025,  0.025,  0.025,  0.025],
       [ 0.025,  0.9  ,  0.025,  0.025,  0.025],
       [ 0.025,  0.025,  0.9  ,  0.025,  0.025],
       [ 0.025,  0.025,  0.025,  0.9  ,  0.025],
       [ 0.025,  0.025,  0.025,  0.025,  0.9  ]])

observables Module

Control and select output from CMCpy simulations

class observables.Observables(show_codes=True, show_messages=False, show_initial_parameters=True, show_matrix_parameters=False, show_fitness_statistics=False, show_code_evolution_statistics=False, show_frozen_results_only=False, print_precision=6, show_all=False)[source]

site_type_spaces Module

Site-type fitness matrices are intended as site-types over rows and amino acids over columns

class site_type_spaces.MirroringSiteTypeSpace(amino_acids, phi, weights=None)[source]

This class models site-types in one-to-one correspondence with amino acids as according to the published models of Ardell and Sella.

>>> aa = amino_acid_spaces.RingAminoAcidSpace(num_aas = 5)
>>> dm = aa.get_distance_matrix()
>>> dm.round(3)
array([[ 0.   ,  0.219,  0.443,  0.424,  0.205],
       [ 0.219,  0.   ,  0.224,  0.357,  0.424],
       [ 0.443,  0.224,  0.   ,  0.133,  0.352],
       [ 0.424,  0.357,  0.133,  0.   ,  0.219],
       [ 0.205,  0.424,  0.352,  0.219,  0.   ]])
>>> st = MirroringSiteTypeSpace(aa,phi = 0.96)
>>> fm = st.get_fitness_matrix()
>>> fm.round(3)
array([[ 1.   ,  0.991,  0.982,  0.983,  0.992],
       [ 0.991,  1.   ,  0.991,  0.986,  0.983],
       [ 0.982,  0.991,  1.   ,  0.995,  0.986],
       [ 0.983,  0.986,  0.995,  1.   ,  0.991],
       [ 0.992,  0.983,  0.986,  0.991,  1.   ]])
get_fitness_matrix()[source]
get_site_type_weights()[source]