bpp-phyl3 3.0.0
HKY85.cpp
Go to the documentation of this file.
1// SPDX-FileCopyrightText: The Bio++ Development Group
2//
3// SPDX-License-Identifier: CECILL-2.1
4
5#include <Bpp/Numeric/Matrix/MatrixTools.h>
6
7#include "../FrequencySet/NucleotideFrequencySet.h"
8#include "HKY85.h"
9
10using namespace bpp;
11
12// From the STL:
13#include <cmath>
14
15using namespace std;
16
17/******************************************************************************/
18
19HKY85::HKY85(
20 shared_ptr<const NucleicAlphabet> alpha,
21 double kappa,
22 double piA,
23 double piC,
24 double piG,
25 double piT) :
26 AbstractParameterAliasable("HKY85."),
27 AbstractReversibleNucleotideSubstitutionModel(alpha, make_shared<CanonicalStateMap>(alpha, false), "HKY85."),
28 kappa_(kappa), k1_(), k2_(), r_(),
29 piA_(piA), piC_(piC), piG_(piG), piT_(piT), piY_(), piR_(),
30 theta_(piG + piC), theta1_(piA / (1. - theta_)), theta2_(piG / theta_),
31 exp1_(), exp21_(), exp22_(), l_(), p_(size_, size_)
32{
33 addParameter_(new Parameter("HKY85.kappa", kappa, Parameter::R_PLUS_STAR));
34 addParameter_(new Parameter("HKY85.theta", theta_, FrequencySetInterface::FREQUENCE_CONSTRAINT_SMALL));
35 addParameter_(new Parameter("HKY85.theta1", theta1_, FrequencySetInterface::FREQUENCE_CONSTRAINT_SMALL));
36 addParameter_(new Parameter("HKY85.theta2", theta2_, FrequencySetInterface::FREQUENCE_CONSTRAINT_SMALL));
37 updateMatrices_();
38}
39
40/******************************************************************************/
41
42void HKY85::updateMatrices_()
43{
44 kappa_ = getParameterValue("kappa");
45 theta_ = getParameterValue("theta");
46 theta1_ = getParameterValue("theta1");
47 theta2_ = getParameterValue("theta2");
48 piA_ = theta1_ * (1. - theta_);
49 piC_ = (1. - theta2_) * theta_;
50 piG_ = theta2_ * theta_;
51 piT_ = (1. - theta1_) * (1. - theta_);
52 piR_ = piA_ + piG_;
53 piY_ = piT_ + piC_;
54 k1_ = kappa_ * piY_ + piR_;
55 k2_ = kappa_ * piR_ + piY_;
56
57 freq_[0] = piA_;
58 freq_[1] = piC_;
59 freq_[2] = piG_;
60 freq_[3] = piT_;
61
62 generator_(0, 0) = -( piC_ + kappa_ * piG_ + piT_);
63 generator_(1, 1) = -( piA_ + piG_ + kappa_ * piT_);
64 generator_(2, 2) = -(kappa_ * piA_ + piC_ + piT_);
65 generator_(3, 3) = -( piA_ + kappa_ * piC_ + piG_ );
66
67 generator_(1, 0) = piA_;
68 generator_(3, 0) = piA_;
69 generator_(0, 1) = piC_;
70 generator_(2, 1) = piC_;
71 generator_(1, 2) = piG_;
72 generator_(3, 2) = piG_;
73 generator_(0, 3) = piT_;
74 generator_(2, 3) = piT_;
75
76 generator_(2, 0) = kappa_ * piA_;
77 generator_(3, 1) = kappa_ * piC_;
78 generator_(0, 2) = kappa_ * piG_;
79 generator_(1, 3) = kappa_ * piT_;
80
81 // Normalization:
82 r_ = isScalable() ? 1. / (2. * (piA_ * piC_ + piC_ * piG_ + piA_ * piT_ + piG_ * piT_ + kappa_ * (piC_ * piT_ + piA_ * piG_))) : 1;
83
84 setScale(r_);
85
86 // Exchangeability:
87 exchangeability_(0, 0) = generator_(0, 0) / piA_;
88 exchangeability_(0, 1) = generator_(0, 1) / piC_;
89 exchangeability_(0, 2) = generator_(0, 2) / piG_;
90 exchangeability_(0, 3) = generator_(0, 3) / piT_;
91
92 exchangeability_(1, 0) = generator_(1, 0) / piA_;
93 exchangeability_(1, 1) = generator_(1, 1) / piC_;
94 exchangeability_(1, 2) = generator_(1, 2) / piG_;
95 exchangeability_(1, 3) = generator_(1, 3) / piT_;
96
97 exchangeability_(2, 0) = generator_(2, 0) / piA_;
98 exchangeability_(2, 1) = generator_(2, 1) / piC_;
99 exchangeability_(2, 2) = generator_(2, 2) / piG_;
100 exchangeability_(2, 3) = generator_(2, 3) / piT_;
101
102 exchangeability_(3, 0) = generator_(3, 0) / piA_;
103 exchangeability_(3, 1) = generator_(3, 1) / piC_;
104 exchangeability_(3, 2) = generator_(3, 2) / piG_;
105 exchangeability_(3, 3) = generator_(3, 3) / piT_;
106
107 // Eigen values:
108 eigenValues_[0] = 0;
109 eigenValues_[1] = -r_ * (kappa_ * piY_ + piR_);
110 eigenValues_[2] = -r_ * (kappa_ * piR_ + piY_);
111 eigenValues_[3] = -r_;
112
113 // Eigen vectors:
114 leftEigenVectors_(0, 0) = piA_;
115 leftEigenVectors_(0, 1) = piC_;
116 leftEigenVectors_(0, 2) = piG_;
117 leftEigenVectors_(0, 3) = piT_;
118
119 leftEigenVectors_(1, 0) = 0.;
120 leftEigenVectors_(1, 1) = piT_ / piY_;
121 leftEigenVectors_(1, 2) = 0.;
122 leftEigenVectors_(1, 3) = -piT_ / piY_;
123
124 leftEigenVectors_(2, 0) = piG_ / piR_;
125 leftEigenVectors_(2, 1) = 0.;
126 leftEigenVectors_(2, 2) = -piG_ / piR_;
127 leftEigenVectors_(2, 3) = 0.;
128
129 leftEigenVectors_(3, 0) = piA_ * piY_ / piR_;
130 leftEigenVectors_(3, 1) = -piC_;
131 leftEigenVectors_(3, 2) = piG_ * piY_ / piR_;
132 leftEigenVectors_(3, 3) = -piT_;
133
134 rightEigenVectors_(0, 0) = 1.;
135 rightEigenVectors_(0, 1) = 0.;
136 rightEigenVectors_(0, 2) = 1.;
137 rightEigenVectors_(0, 3) = 1.;
138
139 rightEigenVectors_(1, 0) = 1.;
140 rightEigenVectors_(1, 1) = 1.;
141 rightEigenVectors_(1, 2) = 0.;
142 rightEigenVectors_(1, 3) = -piR_ / piY_;
143
144 rightEigenVectors_(2, 0) = 1.;
145 rightEigenVectors_(2, 1) = 0.;
146 rightEigenVectors_(2, 2) = -piA_ / piG_;
147 rightEigenVectors_(2, 3) = 1.;
148
149 rightEigenVectors_(3, 0) = 1.;
150 rightEigenVectors_(3, 1) = -piC_ / piT_;
151 rightEigenVectors_(3, 2) = 0.;
152 rightEigenVectors_(3, 3) = -piR_ / piY_;
153}
154
155/******************************************************************************/
156
157double HKY85::Pij_t(size_t i, size_t j, double d) const
158{
159 l_ = rate_ * r_ * d;
160 exp1_ = exp(-l_);
161 exp22_ = exp(-k2_ * l_);
162 exp21_ = exp(-k1_ * l_);
163
164 switch (i)
165 {
166 case 0: // A
167 switch (j)
168 {
169 case 0: return piA_ * (1. + (piY_ / piR_) * exp1_) + (piG_ / piR_) * exp22_; // A
170 case 1: return piC_ * (1. - exp1_); // C
171 case 2: return piG_ * (1. + (piY_ / piR_) * exp1_) - (piG_ / piR_) * exp22_; // G
172 case 3: return piT_ * (1. - exp1_); // T, U
173 default: return 0;
174 }
175 case 1: // C
176 switch (j)
177 {
178 case 0: return piA_ * (1. - exp1_); // A
179 case 1: return piC_ * (1. + (piR_ / piY_) * exp1_) + (piT_ / piY_) * exp21_; // C
180 case 2: return piG_ * (1. - exp1_); // G
181 case 3: return piT_ * (1. + (piR_ / piY_) * exp1_) - (piT_ / piY_) * exp21_; // T, U
182 default: return 0;
183 }
184 case 2: // G
185 switch (j)
186 {
187 case 0: return piA_ * (1. + (piY_ / piR_) * exp1_) - (piA_ / piR_) * exp22_; // A
188 case 1: return piC_ * (1. - exp1_); // C
189 case 2: return piG_ * (1. + (piY_ / piR_) * exp1_) + (piA_ / piR_) * exp22_; // G
190 case 3: return piT_ * (1. - exp1_); // T, U
191 default: return 0;
192 }
193 case 3: // T, U
194 switch (j)
195 {
196 case 0: return piA_ * (1. - exp1_); // A
197 case 1: return piC_ * (1. + (piR_ / piY_) * exp1_) - (piC_ / piY_) * exp21_; // C
198 case 2: return piG_ * (1. - exp1_); // G
199 case 3: return piT_ * (1. + (piR_ / piY_) * exp1_) + (piC_ / piY_) * exp21_; // T, U
200 default: return 0;
201 }
202 default: return 0;
203 }
204}
205
206/******************************************************************************/
207
208double HKY85::dPij_dt(size_t i, size_t j, double d) const
209{
210 l_ = rate_ * r_ * d;
211 exp1_ = exp(-l_);
212 exp22_ = exp(-k2_ * l_);
213 exp21_ = exp(-k1_ * l_);
214
215 switch (i)
216 {
217 case 0: // A
218 switch (j)
219 {
220 case 0: return rate_ * r_ * (piA_ * -(piY_ / piR_) * exp1_ - (piG_ / piR_) * k2_ * exp22_); // A
221 case 1: return rate_ * r_ * (piC_ * exp1_); // C
222 case 2: return rate_ * r_ * (piG_ * -(piY_ / piR_) * exp1_ + (piG_ / piR_) * k2_ * exp22_); // G
223 case 3: return rate_ * r_ * (piT_ * exp1_); // T, U
224 default: return 0;
225 }
226 case 1: // C
227 switch (j)
228 {
229 case 0: return rate_ * r_ * (piA_ * exp1_); // A
230 case 1: return rate_ * r_ * (piC_ * -(piR_ / piY_) * exp1_ - (piT_ / piY_) * k1_ * exp21_); // C
231 case 2: return rate_ * r_ * (piG_ * exp1_); // G
232 case 3: return rate_ * r_ * (piT_ * -(piR_ / piY_) * exp1_ + (piT_ / piY_) * k1_ * exp21_); // T, U
233 default: return 0;
234 }
235 case 2: // G
236 switch (j)
237 {
238 case 0: return rate_ * r_ * (piA_ * -(piY_ / piR_) * exp1_ + (piA_ / piR_) * k2_ * exp22_); // A
239 case 1: return rate_ * r_ * (piC_ * exp1_); // C
240 case 2: return rate_ * r_ * (piG_ * -(piY_ / piR_) * exp1_ - (piA_ / piR_) * k2_ * exp22_); // G
241 case 3: return rate_ * r_ * (piT_ * exp1_); // T, U
242 default: return 0;
243 }
244 case 3: // T, U
245 switch (j)
246 {
247 case 0: return rate_ * r_ * (piA_ * exp1_); // A
248 case 1: return rate_ * r_ * (piC_ * -(piR_ / piY_) * exp1_ + (piC_ / piY_) * k1_ * exp21_); // C
249 case 2: return rate_ * r_ * (piG_ * exp1_); // G
250 case 3: return rate_ * r_ * (piT_ * -(piR_ / piY_) * exp1_ - (piC_ / piY_) * k1_ * exp21_); // T, U
251 default: return 0;
252 }
253 default: return 0;
254 }
255}
256
257/******************************************************************************/
258
259double HKY85::d2Pij_dt2(size_t i, size_t j, double d) const
260{
261 double r_2 = rate_ * rate_ * r_ * r_;
262 l_ = rate_ * r_ * d;
263 double k1_2 = k1_ * k1_;
264 double k2_2 = k2_ * k2_;
265 exp1_ = exp(-l_);
266 exp22_ = exp(-k2_ * l_);
267 exp21_ = exp(-k1_ * l_);
268
269 switch (i)
270 {
271 case 0: // A
272 switch (j)
273 {
274 case 0: return r_2 * (piA_ * (piY_ / piR_) * exp1_ + (piG_ / piR_) * k2_2 * exp22_); // A
275 case 1: return r_2 * (piC_ * -exp1_); // C
276 case 2: return r_2 * (piG_ * (piY_ / piR_) * exp1_ - (piG_ / piR_) * k2_2 * exp22_); // G
277 case 3: return r_2 * (piT_ * -exp1_); // T, U
278 default: return 0;
279 }
280 case 1: // C
281 switch (j)
282 {
283 case 0: return r_2 * (piA_ * -exp1_); // A
284 case 1: return r_2 * (piC_ * (piR_ / piY_) * exp1_ + (piT_ / piY_) * k1_2 * exp21_); // C
285 case 2: return r_2 * (piG_ * -exp1_); // G
286 case 3: return r_2 * (piT_ * (piR_ / piY_) * exp1_ - (piT_ / piY_) * k1_2 * exp21_); // T, U
287 default: return 0;
288 }
289 case 2: // G
290 switch (j)
291 {
292 case 0: return r_2 * (piA_ * (piY_ / piR_) * exp1_ - (piA_ / piR_) * k2_2 * exp22_); // A
293 case 1: return r_2 * (piC_ * -exp1_); // C
294 case 2: return r_2 * (piG_ * (piY_ / piR_) * exp1_ + (piA_ / piR_) * k2_2 * exp22_); // G
295 case 3: return r_2 * (piT_ * -exp1_); // T, U
296 default: return 0;
297 }
298 case 3: // T, U
299 switch (j)
300 {
301 case 0: return r_2 * (piA_ * -exp1_); // A
302 case 1: return r_2 * (piC_ * (piR_ / piY_) * exp1_ - (piC_ / piY_) * k1_2 * exp21_); // C
303 case 2: return r_2 * (piG_ * -exp1_); // G
304 case 3: return r_2 * (piT_ * (piR_ / piY_) * exp1_ + (piC_ / piY_) * k1_2 * exp21_); // T, U
305 default: return 0;
306 }
307 default: return 0;
308 }
309}
310
311/******************************************************************************/
312
313const Matrix<double>& HKY85::getPij_t(double d) const
314{
315 l_ = rate_ * r_ * d;
316 exp1_ = exp(-l_);
317 exp22_ = exp(-k2_ * l_);
318 exp21_ = exp(-k1_ * l_);
319
320 // A
321 p_(0, 0) = piA_ * (1. + (piY_ / piR_) * exp1_) + (piG_ / piR_) * exp22_; // A
322 p_(0, 1) = piC_ * (1. - exp1_); // C
323 p_(0, 2) = piG_ * (1. + (piY_ / piR_) * exp1_) - (piG_ / piR_) * exp22_; // G
324 p_(0, 3) = piT_ * (1. - exp1_); // T, U
325
326 // C
327 p_(1, 0) = piA_ * (1. - exp1_); // A
328 p_(1, 1) = piC_ * (1. + (piR_ / piY_) * exp1_) + (piT_ / piY_) * exp21_; // C
329 p_(1, 2) = piG_ * (1. - exp1_); // G
330 p_(1, 3) = piT_ * (1. + (piR_ / piY_) * exp1_) - (piT_ / piY_) * exp21_; // T, U
331
332 // G
333 p_(2, 0) = piA_ * (1. + (piY_ / piR_) * exp1_) - (piA_ / piR_) * exp22_; // A
334 p_(2, 1) = piC_ * (1. - exp1_); // C
335 p_(2, 2) = piG_ * (1. + (piY_ / piR_) * exp1_) + (piA_ / piR_) * exp22_; // G
336 p_(2, 3) = piT_ * (1. - exp1_); // T, U
337
338 // T, U
339 p_(3, 0) = piA_ * (1. - exp1_); // A
340 p_(3, 1) = piC_ * (1. + (piR_ / piY_) * exp1_) - (piC_ / piY_) * exp21_; // C
341 p_(3, 2) = piG_ * (1. - exp1_); // G
342 p_(3, 3) = piT_ * (1. + (piR_ / piY_) * exp1_) + (piC_ / piY_) * exp21_; // T, U
343
344 return p_;
345}
346
347const Matrix<double>& HKY85::getdPij_dt(double d) const
348{
349 l_ = rate_ * r_ * d;
350 exp1_ = exp(-l_);
351 exp22_ = exp(-k2_ * l_);
352 exp21_ = exp(-k1_ * l_);
353
354 // A
355 p_(0, 0) = rate_ * r_ * (piA_ * -(piY_ / piR_) * exp1_ - (piG_ / piR_) * k2_ * exp22_); // A
356 p_(0, 1) = rate_ * r_ * (piC_ * exp1_); // C
357 p_(0, 2) = rate_ * r_ * (piG_ * -(piY_ / piR_) * exp1_ + (piG_ / piR_) * k2_ * exp22_); // G
358 p_(0, 3) = rate_ * r_ * (piT_ * exp1_); // T, U
359
360 // C
361 p_(1, 0) = rate_ * r_ * (piA_ * exp1_); // A
362 p_(1, 1) = rate_ * r_ * (piC_ * -(piR_ / piY_) * exp1_ - (piT_ / piY_) * k1_ * exp21_); // C
363 p_(1, 2) = rate_ * r_ * (piG_ * exp1_); // G
364 p_(1, 3) = rate_ * r_ * (piT_ * -(piR_ / piY_) * exp1_ + (piT_ / piY_) * k1_ * exp21_); // T, U
365
366 // G
367 p_(2, 0) = rate_ * r_ * (piA_ * -(piY_ / piR_) * exp1_ + (piA_ / piR_) * k2_ * exp22_); // A
368 p_(2, 1) = rate_ * r_ * (piC_ * exp1_); // C
369 p_(2, 2) = rate_ * r_ * (piG_ * -(piY_ / piR_) * exp1_ - (piA_ / piR_) * k2_ * exp22_); // G
370 p_(2, 3) = rate_ * r_ * (piT_ * exp1_); // T, U
371
372 // T, U
373 p_(3, 0) = rate_ * r_ * (piA_ * exp1_); // A
374 p_(3, 1) = rate_ * r_ * (piC_ * -(piR_ / piY_) * exp1_ + (piC_ / piY_) * k1_ * exp21_); // C
375 p_(3, 2) = rate_ * r_ * (piG_ * exp1_); // G
376 p_(3, 3) = rate_ * r_ * (piT_ * -(piR_ / piY_) * exp1_ - (piC_ / piY_) * k1_ * exp21_); // T, U
377
378 return p_;
379}
380
381const Matrix<double>& HKY85::getd2Pij_dt2(double d) const
382{
383 double r_2 = rate_ * rate_ * r_ * r_;
384 l_ = rate_ * r_ * d;
385 double k1_2 = k1_ * k1_;
386 double k2_2 = k2_ * k2_;
387 exp1_ = exp(-l_);
388 exp22_ = exp(-k2_ * l_);
389 exp21_ = exp(-k1_ * l_);
390
391 // A
392 p_(0, 0) = r_2 * (piA_ * (piY_ / piR_) * exp1_ + (piG_ / piR_) * k2_2 * exp22_); // A
393 p_(0, 1) = r_2 * (piC_ * -exp1_); // C
394 p_(0, 2) = r_2 * (piG_ * (piY_ / piR_) * exp1_ - (piG_ / piR_) * k2_2 * exp22_); // G
395 p_(0, 3) = r_2 * (piT_ * -exp1_); // T, U
396
397 // C
398 p_(1, 0) = r_2 * (piA_ * -exp1_); // A
399 p_(1, 1) = r_2 * (piC_ * (piR_ / piY_) * exp1_ + (piT_ / piY_) * k1_2 * exp21_); // C
400 p_(1, 2) = r_2 * (piG_ * -exp1_); // G
401 p_(1, 3) = r_2 * (piT_ * (piR_ / piY_) * exp1_ - (piT_ / piY_) * k1_2 * exp21_); // T, U
402
403 // G
404 p_(2, 0) = r_2 * (piA_ * (piY_ / piR_) * exp1_ - (piA_ / piR_) * k2_2 * exp22_); // A
405 p_(2, 1) = r_2 * (piC_ * -exp1_); // C
406 p_(2, 2) = r_2 * (piG_ * (piY_ / piR_) * exp1_ + (piA_ / piR_) * k2_2 * exp22_); // G
407 p_(2, 3) = r_2 * (piT_ * -exp1_); // T, U
408
409 // T, U
410 p_(3, 0) = r_2 * (piA_ * -exp1_); // A
411 p_(3, 1) = r_2 * (piC_ * (piR_ / piY_) * exp1_ - (piC_ / piY_) * k1_2 * exp21_); // C
412 p_(3, 2) = r_2 * (piG_ * -exp1_); // G
413 p_(3, 3) = r_2 * (piT_ * (piR_ / piY_) * exp1_ + (piC_ / piY_) * k1_2 * exp21_); // T, U
414
415 return p_;
416}
417
418/******************************************************************************/
419
420void HKY85::setFreq(std::map<int, double>& freqs)
421{
422 piA_ = freqs[0];
423 piC_ = freqs[1];
424 piG_ = freqs[2];
425 piT_ = freqs[3];
426 vector<string> thetas(3);
427 thetas[0] = getNamespace() + "theta";
428 thetas[1] = getNamespace() + "theta1";
429 thetas[2] = getNamespace() + "theta2";
430 ParameterList pl = getParameters().createSubList(thetas);
431 pl[0].setValue(piC_ + piG_);
432 pl[1].setValue(piA_ / (piA_ + piT_));
433 pl[2].setValue(piG_ / (piC_ + piG_));
434 setParametersValues(pl);
435}
436
437/******************************************************************************/
bpp::AbstractParameterAliasable::addParameter_
void addParameter_(Parameter *parameter)
bpp::AbstractParameterAliasable::setParametersValues
void setParametersValues(const ParameterList &parameters) override
bpp::AbstractParameterAliasable::getNamespace
std::string getNamespace() const override
bpp::AbstractParameterAliasable::getParameters
const ParameterList & getParameters() const override
bpp::AbstractParameterAliasable::getParameterValue
double getParameterValue(const std::string &name) const override
bpp::AbstractReversibleNucleotideSubstitutionModel
Specialisation abstract class for reversible nucleotide substitution model.
Definition: NucleotideSubstitutionModel.h:78
bpp::AbstractSubstitutionModel::generator_
RowMatrix< double > generator_
The generator matrix of the model.
Definition: AbstractSubstitutionModel.h:298
bpp::AbstractSubstitutionModel::eigenValues_
Vdouble eigenValues_
The vector of eigen values.
Definition: AbstractSubstitutionModel.h:320
bpp::AbstractSubstitutionModel::exchangeability_
RowMatrix< double > exchangeability_
The exchangeability matrix of the model, defined as . When the model is reversible,...
Definition: AbstractSubstitutionModel.h:310
bpp::AbstractSubstitutionModel::leftEigenVectors_
RowMatrix< double > leftEigenVectors_
The matrix made of left eigen vectors (by row) if rightEigenVectors_ is non-singular.
Definition: AbstractSubstitutionModel.h:346
bpp::AbstractSubstitutionModel::isScalable
virtual bool isScalable() const
returns if model is scalable
Definition: AbstractSubstitutionModel.h:469
bpp::AbstractSubstitutionModel::rightEigenVectors_
RowMatrix< double > rightEigenVectors_
The matrix made of right eigen vectors (by column).
Definition: AbstractSubstitutionModel.h:335
bpp::AbstractSubstitutionModel::setScale
void setScale(double scale)
Multiplies the current generator by the given scale.
Definition: AbstractSubstitutionModel.cpp:641
bpp::AbstractTransitionModel::freq_
Vdouble freq_
The vector of equilibrium frequencies.
Definition: AbstractSubstitutionModel.h:143
bpp::AbstractTransitionModel::rate_
double rate_
The rate of the model (default: 1). The generator (and all its vectorial components) is independent o...
Definition: AbstractSubstitutionModel.h:138
bpp::CanonicalStateMap
This class implements a state map where all resolved states are modeled.
Definition: StateMap.h:168
bpp::FrequencySetInterface::FREQUENCE_CONSTRAINT_SMALL
static std::shared_ptr< IntervalConstraint > FREQUENCE_CONSTRAINT_SMALL
Definition: FrequencySet.h:90
bpp::HKY85::piY_
double piY_
Definition: HKY85.h:147
bpp::HKY85::piC_
double piC_
Definition: HKY85.h:147
bpp::HKY85::piG_
double piG_
Definition: HKY85.h:147
bpp::HKY85::theta1_
double theta1_
Definition: HKY85.h:147
bpp::HKY85::dPij_dt
double dPij_dt(size_t i, size_t j, double d) const
Definition: HKY85.cpp:208
bpp::HKY85::updateMatrices_
void updateMatrices_()
Compute and diagonalize the matrix, and fill the eigenValues_, leftEigenVectors_ and rightEigenVecto...
Definition: HKY85.cpp:42
bpp::HKY85::k1_
double k1_
Definition: HKY85.h:147
bpp::HKY85::Pij_t
double Pij_t(size_t i, size_t j, double d) const
Definition: HKY85.cpp:157
bpp::HKY85::getPij_t
const Matrix< double > & getPij_t(double d) const
Definition: HKY85.cpp:313
bpp::HKY85::k2_
double k2_
Definition: HKY85.h:147
bpp::HKY85::d2Pij_dt2
double d2Pij_dt2(size_t i, size_t j, double d) const
Definition: HKY85.cpp:259
bpp::HKY85::exp21_
double exp21_
Definition: HKY85.h:148
bpp::HKY85::kappa_
double kappa_
Definition: HKY85.h:147
bpp::HKY85::exp1_
double exp1_
Definition: HKY85.h:148
bpp::HKY85::l_
double l_
Definition: HKY85.h:148
bpp::HKY85::piA_
double piA_
Definition: HKY85.h:147
bpp::HKY85::theta2_
double theta2_
Definition: HKY85.h:147
bpp::HKY85::exp22_
double exp22_
Definition: HKY85.h:148
bpp::HKY85::p_
RowMatrix< double > p_
Definition: HKY85.h:149
bpp::HKY85::theta_
double theta_
Definition: HKY85.h:147
bpp::HKY85::setFreq
void setFreq(std::map< int, double > &freqs)
This method is redefined to actualize the corresponding parameters piA, piT, piG and piC too.
Definition: HKY85.cpp:420
bpp::HKY85::r_
double r_
Definition: HKY85.h:147
bpp::HKY85::piT_
double piT_
Definition: HKY85.h:147
bpp::HKY85::getdPij_dt
const Matrix< double > & getdPij_dt(double d) const
Definition: HKY85.cpp:347
bpp::HKY85::getd2Pij_dt2
const Matrix< double > & getd2Pij_dt2(double d) const
Definition: HKY85.cpp:381
bpp::HKY85::piR_
double piR_
Definition: HKY85.h:147
bpp::ParameterList::createSubList
virtual ParameterList createSubList(const std::vector< std::string > &names) const
bpp::Parameter::R_PLUS_STAR
static const std::shared_ptr< IntervalConstraint > R_PLUS_STAR
bpp
Defines the basic types of data flow nodes.
bpp::exp
ExtendedFloat exp(const ExtendedFloat &ef)
Definition: ExtendedFloat.h:474