# 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)
                        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(
                        1, N_CODON, N_CODON)
                expDxx = numpy.array([expDyy[r].transpose() for r in range(1)])
                V = (expDxx - expDyy) / Dxx_Dyy
            with numpy.errstate(under='ignore'):  # OK if some values 0
                numpy.copyto(V, self.mu * t * expDxx,
                             where=Dxx_Dyy_lt_ALMOST_ZERO)
            self._cached[('V', t)] = V
        V = self._cached[('V', t)]

        with numpy.errstate(under='ignore'):  # don't worry if some values 0
            dM_param = broadcastMatrixMultiply(self.A, broadcastMatrixMultiply(self.B[param] * V, self.Ainv))  # noqa: E501
            if tips is None:
                return numpy.tile(dM_param, (self.nsites, 1, 1))
            else:
                newdM_param = numpy.zeros((len(tips), N_CODON))
                for i in range(len(tips)):
                    newdM_param[i] = (dM_param[0][:, tips[i]])
                if gaps is not None:
                    newdM_param[gaps] = numpy.zeros(N_CODON, dtype='float')
                return newdM_param

    def _update_Pxy(self):
        """Update `Pxy` using current `omega`, `kappa`, and `Phi_x`."""
        numpy.copyto(self.Pxy_no_omega, self.Phi_x.transpose(),
                     where=CODON_SINGLEMUT)

# 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)

                        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(
                        1, N_CODON, N_CODON)
                expDxx = numpy.array([expDyy[r].transpose() for r in range(1)])
                V = (expDxx - expDyy) / Dxx_Dyy
            with numpy.errstate(under='ignore'):  # OK if some values 0
                numpy.copyto(V, self.mu * t * expDxx,
                             where=Dxx_Dyy_lt_ALMOST_ZERO)
            self._cached[('V', t)] = V
        V = self._cached[('V', t)]

        with numpy.errstate(under='ignore'):  # don't worry if some values 0
            dM_param = broadcastMatrixMultiply(self.A, broadcastMatrixMultiply(self.B[param] * V, self.Ainv))  # noqa: E501
            if tips is None:
                return numpy.tile(dM_param, (self.nsites, 1, 1))
            else:
                newdM_param = numpy.zeros((len(tips), N_CODON))
                for i in range(len(tips)):
                    newdM_param[i] = (dM_param[0][:, tips[i]])
                if gaps is not None:
                    newdM_param[gaps] = numpy.zeros(N_CODON, dtype='float')
                return newdM_param

    def _update_Pxy(self):
        """Update `Pxy` using current `omega`, `kappa`, and `Phi_x`."""
        numpy.copyto(self.Pxy_no_omega, self.Phi_x.transpose(),
                     where=CODON_SINGLEMUT)

        self.Pxy_no_omega[0][CODON_TRANSITION] *= self.kappa