# Path: phydmslib/constants.py
# N_CODON = len(CODON_TO_INDEX)
# 
# CODON_NONSYN = numpy.full((N_CODON, N_CODON), False, dtype='bool')
# 
# CODON_TRANSITION = numpy.full((N_CODON, N_CODON), False, dtype='bool')
# 
# CODON_SINGLEMUT = numpy.full((N_CODON, N_CODON), False, dtype='bool')
# 
# ALMOST_ZERO = 1.0e-6
# 
# CODON_NT_INDEX = numpy.full((3, N_CODON), -1, dtype='int')
# 
# STOP_CODON_TO_NT_INDICES = []
# 
# N_STOP = 0
# 
# N_NT = len(INDEX_TO_NT)
# 
# N_AA = len(INDEX_TO_AA)
# 
# INDEX_TO_NT = dict(enumerate('ACGT'))
# 
# CODON_NT_MUT = numpy.full((N_NT, N_CODON, N_CODON), False, dtype='bool')
# 
# CODON_NT_COUNT = numpy.zeros((N_NT, N_CODON), dtype='int')
# 
# CODON_TO_AA = []
# 
# INDEX_TO_AA = dict(enumerate('ACDEFGHIKLMNPQRSTVWY'))
# 
# AA_TO_INDEX = {aa: i for (i, aa) in INDEX_TO_AA.items()}
# 
# CODON_NT = numpy.full((3, N_NT, N_CODON), False, dtype='bool')


import copy
import six
import abc
import warnings
import numpy
import scipy.misc
import scipy.optimize
import scipy.linalg
import scipy.stats
import scipy.special
    import doctest
from phydmslib.numutils import (broadcastMatrixMultiply,
                                broadcastMatrixVectorMultiply,
                                broadcastGetCols)
from phydmslib.constants import (N_CODON, CODON_NONSYN, CODON_TRANSITION,
                                 CODON_SINGLEMUT, ALMOST_ZERO, CODON_NT_INDEX,
                                 STOP_CODON_TO_NT_INDICES, N_STOP, N_NT, N_AA,
                                 INDEX_TO_NT, CODON_NT_MUT, CODON_NT_COUNT,
                                 CODON_TO_AA, INDEX_TO_AA, AA_TO_INDEX,
                                 CODON_NT)
                if gaps is not None:
                    dM_param[gaps] = numpy.zeros(N_CODON, dtype='float')
            return dM_param

        paramval = getattr(self, param)
        if isinstance(paramval, float):
            paramisvec = False
        else:
            assert isinstance(paramval, numpy.ndarray) and paramval.ndim == 1
            paramisvec = True
            paramlength = paramval.shape[0]

        if ('expD', t) not in self._cached:
            self._cached[('expD', t)] = numpy.exp(self.D * self.mu * t)
        expD = self._cached[('expD', t)]

        if ('V', t) not in self._cached:
            if 'Dxx_Dyy' not in self._cached:
                Dyy = numpy.tile(self.D, (1, N_CODON)).reshape(
                        self.nsites, N_CODON, N_CODON)
                Dxx = numpy.array([Dyy[r].transpose() for r in
                                   range(self.nsites)])
                self._cached['Dxx_Dyy'] = Dxx - Dyy
            Dxx_Dyy = self._cached['Dxx_Dyy']
            if 'Dxx_Dyy_lt_ALMOST_ZERO' not in self._cached:
                self._cached['Dxx_Dyy_lt_ALMOST_ZERO'] = numpy.fabs(
                        Dxx_Dyy) < ALMOST_ZERO
            Dxx_Dyy_lt_ALMOST_ZERO = self._cached['Dxx_Dyy_lt_ALMOST_ZERO']
            with numpy.errstate(divide='raise', under='ignore',
                                over='raise', invalid='ignore'):

# Path: phydmslib/constants.py
# N_CODON = len(CODON_TO_INDEX)
# 
# CODON_NONSYN = numpy.full((N_CODON, N_CODON), False, dtype='bool')
# 
# CODON_TRANSITION = numpy.full((N_CODON, N_CODON), False, dtype='bool')
# 
# CODON_SINGLEMUT = numpy.full((N_CODON, N_CODON), False, dtype='bool')
# 
# ALMOST_ZERO = 1.0e-6
# 
# CODON_NT_INDEX = numpy.full((3, N_CODON), -1, dtype='int')
# 
# STOP_CODON_TO_NT_INDICES = []
# 
# N_STOP = 0
# 
# N_NT = len(INDEX_TO_NT)
# 
# N_AA = len(INDEX_TO_AA)
# 
# INDEX_TO_NT = dict(enumerate('ACGT'))
# 
# CODON_NT_MUT = numpy.full((N_NT, N_CODON, N_CODON), False, dtype='bool')
# 
# CODON_NT_COUNT = numpy.zeros((N_NT, N_CODON), dtype='int')
# 
# CODON_TO_AA = []
# 
# INDEX_TO_AA = dict(enumerate('ACDEFGHIKLMNPQRSTVWY'))
# 
# AA_TO_INDEX = {aa: i for (i, aa) in INDEX_TO_AA.items()}
# 
# CODON_NT = numpy.full((3, N_NT, N_CODON), False, dtype='bool')


# Path: phydmslib/models.py
import copy
import six
import abc
import warnings
import numpy
import scipy.misc
import scipy.optimize
import scipy.linalg
import scipy.stats
import scipy.special
    import doctest
from phydmslib.numutils import (broadcastMatrixMultiply,
                                broadcastMatrixVectorMultiply,
                                broadcastGetCols)
from phydmslib.constants import (N_CODON, CODON_NONSYN, CODON_TRANSITION,
                                 CODON_SINGLEMUT, ALMOST_ZERO, CODON_NT_INDEX,
                                 STOP_CODON_TO_NT_INDICES, N_STOP, N_NT, N_AA,
                                 INDEX_TO_NT, CODON_NT_MUT, CODON_NT_COUNT,
                                 CODON_TO_AA, INDEX_TO_AA, AA_TO_INDEX,
                                 CODON_NT)

                if gaps is not None:
                    dM_param[gaps] = numpy.zeros(N_CODON, dtype='float')
            return dM_param

        paramval = getattr(self, param)
        if isinstance(paramval, float):
            paramisvec = False
        else:
            assert isinstance(paramval, numpy.ndarray) and paramval.ndim == 1
            paramisvec = True
            paramlength = paramval.shape[0]

        if ('expD', t) not in self._cached:
            self._cached[('expD', t)] = numpy.exp(self.D * self.mu * t)
        expD = self._cached[('expD', t)]

        if ('V', t) not in self._cached:
            if 'Dxx_Dyy' not in self._cached:
                Dyy = numpy.tile(self.D, (1, N_CODON)).reshape(
                        self.nsites, N_CODON, N_CODON)
                Dxx = numpy.array([Dyy[r].transpose() for r in
                                   range(self.nsites)])
                self._cached['Dxx_Dyy'] = Dxx - Dyy
            Dxx_Dyy = self._cached['Dxx_Dyy']
            if 'Dxx_Dyy_lt_ALMOST_ZERO' not in self._cached:
                self._cached['Dxx_Dyy_lt_ALMOST_ZERO'] = numpy.fabs(
                        Dxx_Dyy) < ALMOST_ZERO
            Dxx_Dyy_lt_ALMOST_ZERO = self._cached['Dxx_Dyy_lt_ALMOST_ZERO']
            with numpy.errstate(divide='raise', under='ignore',
                                over='raise', invalid='ignore'):

                expDyy = numpy.tile(expD, (1, N_CODON)).reshape(