bpp-phyl3  3.0.0
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
AlignedPhyloLikelihoodAutoCorrelation.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: The Bio++ Development Group
2 //
3 // SPDX-License-Identifier: CECILL-2.1
4 
5 #ifndef BPP_PHYL_LIKELIHOOD_PHYLOLIKELIHOODS_AUTOCORRELATIONOFALIGNEDPHYLOLIKELIHOOD_H
6 #define BPP_PHYL_LIKELIHOOD_PHYLOLIKELIHOODS_AUTOCORRELATIONOFALIGNEDPHYLOLIKELIHOOD_H
7 
8 
11 
12 // From Numeric
13 
14 #include "HmmLikelihood_DF.h"
16 
17 
18 namespace bpp
19 {
29 {
30 private:
31  std::shared_ptr<HmmPhyloAlphabet> hma_;
32 
33  std::shared_ptr<AutoCorrelationTransitionMatrix> htm_;
34 
35  std::shared_ptr<HmmPhyloEmissionProbabilities> hpep_;
36 
37  mutable std::shared_ptr<HmmLikelihood_DF> hmm_;
38 
39 public:
42  std::shared_ptr<PhyloLikelihoodContainer> pC,
43  const std::vector<size_t>& nPhylo,
44  bool inCollection = true);
45 
46 protected:
53  hma_(mlc.hma_),
54  htm_(mlc.htm_),
55  hpep_(mlc.hpep_),
56  hmm_(mlc.hmm_)
57  {}
58 
60  {
62  hma_ = mlc.hma_;
63  htm_ = mlc.htm_;
64  hpep_ = mlc.hpep_;
65  hmm_ = mlc.hmm_;
66  return *this;
67  }
68 
70  {
71  return new AlignedPhyloLikelihoodAutoCorrelation(*this);
72  }
73 
74 public:
76 
77 public:
78  void setNamespace(const std::string& nameSpace);
79 
80  void fireParameterChanged(const ParameterList& parameters);
81 
88  {
89  return *hmm_;
90  }
91 
92  std::shared_ptr<LikelihoodCalculation> getLikelihoodCalculation () const
93  {
94  return hmm_;
95  }
96 
98  {
99  return *hmm_;
100  }
101 
102  std::shared_ptr<AlignedLikelihoodCalculation> getAlignedLikelihoodCalculation () const
103  {
104  return hmm_;
105  }
106 
107  /*
108  *@brief return the posterior probabilities of subprocess on each site.
109  *
110  *@return MatrixXd sites x states
111  */
113  {
114  VVdouble pp;
115  auto mat = hmm_->getHiddenStatesPosteriorProbabilities().transpose();
116  copyEigenToBpp(mat, pp);
117  return pp;
118  }
119 
121  {
122  Vdouble pp;
123  auto mat = hmm_->getHiddenStatesPosteriorProbabilitiesForASite(site);
124  copyEigenToBpp(mat, pp);
125  return pp;
126  }
127 
128  const Eigen::MatrixXd& getHmmTransitionMatrix() const
129  {
130  return hmm_->getHmmTransitionMatrix();
131  }
132 };
133 } // end of namespace bpp.
134 #endif // BPP_PHYL_LIKELIHOOD_PHYLOLIKELIHOODS_AUTOCORRELATIONOFALIGNEDPHYLOLIKELIHOOD_H
The AlignedPhyloLikelihoodSet abstract class.
AbstractAlignedPhyloLikelihoodSet & operator=(const AbstractAlignedPhyloLikelihoodSet &soap)
The PhyloLikelihoodSet class, to manage a subset of PhyloLikelihoods from a given PhyloLikelihoodCont...
const Context & context() const override
Likelihood framework based on a hmm of simple likelihoods.
AlignedPhyloLikelihoodAutoCorrelation & operator=(const AlignedPhyloLikelihoodAutoCorrelation &mlc)
AlignedLikelihoodCalculation & alignedLikelihoodCalculation() const
std::shared_ptr< AutoCorrelationTransitionMatrix > htm_
AlignedPhyloLikelihoodAutoCorrelation(const AlignedPhyloLikelihoodAutoCorrelation &mlc)
std::shared_ptr< LikelihoodCalculation > getLikelihoodCalculation() const
std::shared_ptr< AlignedLikelihoodCalculation > getAlignedLikelihoodCalculation() const
AlignedPhyloLikelihoodAutoCorrelation(Context &context, std::shared_ptr< PhyloLikelihoodContainer > pC, const std::vector< size_t > &nPhylo, bool inCollection=true)
std::shared_ptr< HmmPhyloEmissionProbabilities > hpep_
AlignedPhyloLikelihoodAutoCorrelation * clone() const
Context for dataflow node construction.
Definition: DataFlow.h:527
Defines the basic types of data flow nodes.
std::vector< double > Vdouble
template void copyEigenToBpp(const ExtendedFloatMatrixXd &eigenMatrix, std::vector< std::vector< double >> &bppMatrix)
std::vector< Vdouble > VVdouble