51 AbstractParameterAliasable(
"RE08."),
53 simpleModel_(simpleModel),
55 simpleExchangeabilities_(),
56 exp_(), p_(), lambda_(lambda), mu_(mu),
57 nestedPrefix_(
"model_" + simpleModel->getNamespace())
59 addParameter_(
new Parameter(
"RE08.lambda", lambda, &Parameter::R_PLUS));
60 addParameter_(
new Parameter(
"RE08.mu", mu, &Parameter::R_PLUS));
62 addParameters_(simpleModel->getParameters());
83 for(
size_t i = 0; i <
size_ - 1; i++)
92 for (
size_t i = 0; i <
size_ - 1; i++)
94 for (
size_t j = 0; j <
size_ - 1; j++)
138 return 1. - f * (1. - exp(-(
lambda_ +
mu_) * d));
217 for (
size_t i = 0; i <
size_ - 1; i++)
219 for (
size_t j = 0; j <
size_ - 1; j++)
225 for(
size_t j = 0; j <
size_ - 1; j++)
230 for(
size_t i = 0; i <
size_ - 1; i++)
242 RowMatrix<double> simpleDP =
simpleModel_->getdPij_dt(d);
244 for (
size_t i = 0; i <
size_ - 1; i++)
246 for (
size_t j = 0; j <
size_ - 1; j++)
248 p_(i, j) = simpleDP(i, j) * exp(-
mu_ * d)
253 for (
size_t j = 0; j <
size_ - 1; j++)
258 for (
size_t i = 0; i <
size_ - 1; i++)
270 RowMatrix<double> simpleDP =
simpleModel_->getdPij_dt(d);
271 RowMatrix<double> simpleD2P =
simpleModel_->getd2Pij_dt2(d);
273 for (
size_t i = 0; i <
size_ - 1; i++)
275 for (
size_t j = 0; j <
size_ - 1; j++)
277 p_(i, j) = simpleD2P(i, j) * exp(-
mu_ * d)
278 - 2 *
mu_ * simpleDP(i, j) * exp(-
mu_ * d)
283 for (
size_t j = 0; j <
size_ - 1; j++)
288 for(
size_t i = 0; i <
size_ - 1; i++)
297 double RE08::getInitValue(
size_t i,
int state)
const throw (IndexOutOfBoundsException, BadIntException)
299 if (i >= size_)
throw IndexOutOfBoundsException(
"RE08::getInitValue", i, 0, size_ - 1);
300 if (state < -1 || !getAlphabet()->isIntInAlphabet(state))
301 throw BadIntException(state,
"RE08::getInitValue. Character " + getAlphabet()->intToChar(state) +
" is not allowed in model.");
302 if (i == size_ - 1 && state == -1)
return 1.;
303 vector<int> states = getAlphabet()->getAlias(state);
304 for (
size_t j = 0; j < states.size(); j++)
305 if ((
int)i == states[j])
return 1.;
313 AbstractSubstitutionModel::setNamespace(prefix);
const Matrix< double > & getd2Pij_dt2(double d) const
virtual void updateMatrices()
Diagonalize the matrix, and fill the eigenValues_, iEigenValues_, leftEigenVectors_ and rightEigenVe...
void setNamespace(const std::string &prefix)
RowMatrix< double > exchangeability_
The exchangeability matrix of the model, defined as . When the model is reversible, this matrix is symetric.
const Matrix< double > & getPij_t(double d) const
This class implements a state map where all resolved states are modeled.
RowMatrix< double > simpleGenerator_
RowMatrix< double > rightEigenVectors_
The matrix made of right eigen vectors (by column).
Vdouble eigenValues_
The vector of eigen values.
Vdouble freq_
The vector of equilibrium frequencies.
const Matrix< double > & getdPij_dt(double d) const
void updateMatrices()
Compute and diagonalize the matrix, and fill the eigenValues_, leftEigenVectors_ and rightEigenVecto...
Partial implementation of the ReversibleSubstitutionModel interface.
double d2Pij_dt2(size_t i, size_t j, double d) const
Interface for reversible substitution models.
RowMatrix< double > leftEigenVectors_
The matrix made of left eigen vectors (by row) if rightEigenVectors_ is non-singular.
RowMatrix< double > generator_
The generator matrix of the model.
RowMatrix< double > simpleExchangeabilities_
size_t size_
The size of the generator, i.e. the number of states.
double dPij_dt(size_t i, size_t j, double d) const
RE08(ReversibleSubstitutionModel *simpleModel, double lambda=0, double mu=0)
Build a new Rivas & Eddy model from a standard substitution model.
std::string nestedPrefix_
virtual size_t getNumberOfStates() const =0
Get the number of states.
double getInitValue(size_t i, int state) const
std::auto_ptr< ReversibleSubstitutionModel > simpleModel_
double Pij_t(size_t i, size_t j, double d) const