42 #include <Bpp/Numeric/Matrix/MatrixTools.h> 43 #include <Bpp/Numeric/VectorTools.h> 44 #include <Bpp/Numeric/Matrix/EigenValue.h> 47 #include <Bpp/Seq/Alphabet/WordAlphabet.h> 48 #include <Bpp/Seq/Container/SequenceContainerTools.h> 60 const std::vector<SubstitutionModel*>& modelVector,
61 const std::string& prefix) :
62 AbstractParameterAliasable((prefix ==
"") ?
"Word." : prefix),
64 (prefix ==
"") ?
"Word." : prefix)
69 for (i = 0; i < nbmod - 1; i++)
71 addParameter_(
new Parameter(
"Word.relrate" + TextTools::toString(i + 1), 1.0 / static_cast<int>(nbmod - i), &Parameter::PROP_CONSTRAINT_EX));
80 const std::string& prefix) :
81 AbstractParameterAliasable((prefix ==
"") ?
"Word." : prefix),
88 const std::string& prefix) :
89 AbstractParameterAliasable((prefix ==
"") ?
"Word." : prefix),
92 (prefix ==
"") ?
"Word." : prefix)
97 for (i = 0; i < num - 1; i++)
99 addParameter_(
new Parameter(
"Word.relrate" + TextTools::toString(i + 1), 1.0 / static_cast<int>(num - i ), &Parameter::PROP_CONSTRAINT_EX));
111 for (i = 0; i < nbmod - 1; i++)
113 y = getParameterValue(
"relrate" + TextTools::toString(i + 1));
130 for (i = 0; i < salph; i++)
134 for (p = nbmod; p > 0; p--)
136 m =
VSubMod_[p-1]->getNumberOfStates();
145 vector<const Matrix<double>*> vM;
149 for (i = 0; i < nbmod; i++)
160 for (i = 0; i < nbStates; i++)
162 for (j = 0; j < nbStates; j++)
167 for (p = nbmod; p > 0; p--)
169 t =
VSubMod_[p - 1]->getNumberOfStates();
170 x *= (*vM[p - 1])(i2 % t, j2 % t);
182 vector<const Matrix<double>*> vM, vdM;
186 for (i = 0; i < nbmod; i++)
198 for (i = 0; i < nbStates; i++)
200 for (j = 0; j < nbStates; j++)
203 for (q = 0; q < nbmod; q++)
208 for (p = nbmod; p > 0; p--)
210 t =
VSubMod_[p - 1]->getNumberOfStates();
212 x *= (*vM[p - 1])(i2 % t, j2 % t);
214 x *=
rate_ *
Vrate_[p - 1] * (*vdM[p - 1])(i2 % t, j2 % t);
229 vector<const Matrix<double>*> vM, vdM, vd2M;
233 for (i = 0; i < nbmod; i++)
248 for (i = 0; i < nbStates; i++)
250 for (j = 0; j < nbStates; j++)
253 for (q = 1; q < nbmod; q++)
255 for (b = 0; b < q; b++)
260 for (p = nbmod; p > 0; p--)
262 t =
VSubMod_[p - 1]->getNumberOfStates();
263 if ((p - 1 == q) || (p - 1 == b))
264 x *=
rate_ *
Vrate_[p - 1] * (*vdM[p - 1])(i2 % t, j2 % t);
266 x *= (*vM[p - 1])(i2 % t, j2 % t);
277 for (q = 0; q < nbmod; q++)
282 for (p = nbmod; p > 0; p--)
284 t =
VSubMod_[p - 1]->getNumberOfStates();
286 x *= (*vM[p - 1])(i2 % t, j2 % t);
304 string s =
"WordSubstitutionModel model: " +
VSubMod_[0]->getName();
305 for (
size_t i = 1; i < nbmod - 1; i++)
Interface for all substitution models.
RowMatrix< double > d2pijt_
RowMatrix< double > pijt_
These ones are for bookkeeping:
virtual const RowMatrix< double > & getPij_t(double d) const
std::vector< double > Vrate_
Vdouble freq_
The vector of equilibrium frequencies.
virtual size_t getNumberOfStates() const
Get the number of states.
virtual void completeMatrices()
Called by updateMatrices to handle specific modifications for inheriting classes. ...
virtual const RowMatrix< double > & getd2Pij_dt2(double d) const
Abstract Basal class for words of substitution models.
virtual const RowMatrix< double > & getdPij_dt(double d) const
virtual std::string getName() const
Get the name of the model.
double rate_
The rate of the model (default: 1). The generator (and all its vectorial components) is independent o...
const Alphabet * getAlphabet() const
RowMatrix< double > dpijt_
virtual void updateMatrices()
Diagonalize the matrix, and fill the eigenValues_, iEigenValues_, leftEigenVectors_ and rightEigenVe...
WordSubstitutionModel(const std::vector< SubstitutionModel *> &modelVector, const std::string &prefix="")
Build a new WordSubstitutionModel object from a Vector of pointers to SubstitutionModels.
void updateMatrices()
Diagonalize the matrix, and fill the eigenValues_, iEigenValues_, leftEigenVectors_ and rightEigenVe...
Map the states of a given alphabet which have a model state.
std::vector< SubstitutionModel * > VSubMod_