3
// Created by: Laurent Gueguen
4
// Created on: July 2009
8
Copyright or © or Copr. CNRS, (November 16, 2004)
10
This software is a computer program whose purpose is to provide classes
11
for phylogenetic data analysis.
13
This software is governed by the CeCILL license under French law and
14
abiding by the rules of distribution of free software. You can use,
15
modify and/ or redistribute the software under the terms of the CeCILL
16
license as circulated by CEA, CNRS and INRIA at the following URL
17
"http://www.cecill.info".
19
As a counterpart to the access to the source code and rights to copy,
20
modify and redistribute granted by the license, users are provided only
21
with a limited warranty and the software's author, the holder of the
22
economic rights, and the successive licensors have only limited
25
In this respect, the user's attention is drawn to the risks associated
26
with loading, using, modifying and/or developing or reproducing the
27
software by the user in light of its specific status of free software,
28
that may mean that it is complicated to manipulate, and that also
29
therefore means that it is reserved for developers and experienced
30
professionals having in-depth computer knowledge. Users are therefore
31
encouraged to load and test the software's suitability as regards their
32
requirements in conditions enabling the security of their systems and/or
33
data to be ensured and, more generally, to use and operate it in the
34
same conditions as regards security.
36
The fact that you are presently reading this means that you have had
37
knowledge of the CeCILL license and that you accept its terms.
43
#include "CodonAsynonymousFrequenciesReversibleSubstitutionModel.h"
49
* @brief The Yang and Nielsen (1998) substitution model for codons.
50
* @author Laurent Guéguen
52
* This model has one rate of transitions and one rate of
53
* transversion. It also allows distinct equilibrium frequencies
54
* between codons. A multiplicative factor accounts for the selective
55
* restraints at the amino acid level, depending on the synonymy of
58
* For codons @f$i=i_1i_2i_3@f$ and @f$j=j_1j_2j_3@f$, the generator
59
* term @f$Q_{ij} (i \neq j)@f$ is:
61
* 0 if 2 or 3 of the pair @f$(i_1,j_1) (i_2,j_2) (i_3,j_3) @f$ are different.
63
* @f$\mu \pi_j \omega@f$ if exactly 1 of the pairs @f$(i_1,j_1)
64
* (i_2,j_2) (i_3,j_3) @f$ is different, that difference is a
65
* transversion and amino acids coded by i and j are different.
67
* @f$\mu \pi_j @f$ if exactly 1 of the pairs @f$(i_1,j_1) (i_2,j_2)
68
* (i_3,j_3) @f$ is different, that difference is a transversion and
69
* amino acids coded by i and j are the same.
71
* @f$\mu \kappa \pi_j \omega@f$ if exactly 1 of the pairs
72
* @f$(i_1,j_1) (i_2,j_2) (i_3,j_3) @f$ is different, that difference
73
* is a transition and amino acids coded by i and j are different.
75
* @f$\mu \kappa \pi_j @f$ if exactly 1 of the pairs @f$(i_1,j_1)
76
* (i_2,j_2) (i_3,j_3) @f$ is different, that difference is a
77
* transition and amino acids coded by @f$i@f$ and @f$j@f$ are the
80
* @f$\mu@f$ is a normalization factor.
82
* This model includes 2 parameters (@f$\kappa@f$ and @f$\omega@f$).
83
* The codon frequencies @f$\pi_j@f$ are either observed or infered.
86
* - Yang Z. and Nielsen R. (1998), _Journal of Molecular Evolution_ 46:409--418.
89
public AbstractReversibleSubstitutionModel
93
CodonAsynonymousFrequenciesReversibleSubstitutionModel pmodel_;
97
* @brief Tools to make the link between the Parameters of the
98
* object and those of pmodel_.
102
std::map<std::string,std::string> mapParNamesFromPmodel_;
104
ParameterList lParPmodel_;
107
YN98(const GeneticCode* gc, FrequenciesSet* codonFreqs);
109
YN98(const YN98& yn98) : AbstractReversibleSubstitutionModel(yn98),
110
pmodel_(yn98.pmodel_),
111
mapParNamesFromPmodel_(yn98.mapParNamesFromPmodel_),
112
lParPmodel_(yn98.lParPmodel_) {}
115
YN98* clone() const { return new YN98(*this); }
119
std::string getName() const { return "YN98"; }
121
const Vdouble& getFrequencies() const { return pmodel_.getFrequencies(); }
123
const Matrix<double>& getGenerator() const { return pmodel_.getGenerator(); }
125
const Vdouble& getEigenValues() const { return pmodel_.getEigenValues(); }
127
const Matrix<double>& getRowLeftEigenVectors() const { return pmodel_.getRowLeftEigenVectors(); }
129
const Matrix<double>& getColumnRightEigenVectors() const { return pmodel_.getColumnRightEigenVectors(); }
131
double freq(unsigned int i) const { return pmodel_.freq(i); }
133
double Qij(unsigned int i, unsigned int j) const { return pmodel_.Qij(i,j); }
135
double Pij_t (unsigned int i, unsigned int j, double t) const { return pmodel_.Pij_t(i, j, t); }
136
double dPij_dt (unsigned int i, unsigned int j, double t) const { return pmodel_.dPij_dt(i, j, t); }
137
double d2Pij_dt2(unsigned int i, unsigned int j, double t) const { return pmodel_.d2Pij_dt2(i, j, t); }
139
const Matrix<double>& getPij_t (double d) const { return pmodel_.getPij_t(d); }
140
const Matrix<double>& getdPij_dt (double d) const { return pmodel_.getdPij_dt(d); }
141
const Matrix<double>& getd2Pij_dt2(double d) const { return pmodel_.getd2Pij_dt2(d); }
143
void setFreq(std::map<int, double>& m);
145
unsigned int getNumberOfStates() const { return pmodel_.getNumberOfStates(); }
147
double getInitValue(unsigned int i, int state) const throw (BadIntException) { return pmodel_.getInitValue(i,state); }
149
void enableEigenDecomposition(bool yn) { eigenDecompose_ = 1; }
151
bool enableEigenDecomposition() { return pmodel_.enableEigenDecomposition(); }
153
const FrequenciesSet& getFreq() const { return pmodel_.getFreq(); }
155
double getRate() const { return pmodel_.getRate();}
157
void setRate(double rate) { pmodel_.setRate(rate);}
159
void addRateParameter();
162
void updateMatrices();
166
} //end of namespace bpp.