MarginalAncestralReconstruction.h

Go to the documentation of this file.

// SPDX-FileCopyrightText: The Bio++ Development Group
//
// SPDX-License-Identifier: CECILL-2.1
 
#ifndef BPP_PHYL_LIKELIHOOD_MARGINALANCESTRALRECONSTRUCTION_H
#define BPP_PHYL_LIKELIHOOD_MARGINALANCESTRALRECONSTRUCTION_H
 
 
#include "../AncestralStateReconstruction.h"
#include "DataFlow/DataFlowCWise.h"
#include "DataFlow/LikelihoodCalculationSingleProcess.h"
 
// From bpp-seq:
#include <Bpp/Seq/Alphabet/Alphabet.h>
#include <Bpp/Seq/Container/AlignedSequenceContainer.h>
#include <Bpp/Seq/Sequence.h>
 
// From the STL:
#include <vector>
 
namespace bpp
{
class MarginalAncestralReconstruction :
  public virtual AncestralStateReconstruction
{
private:
  std::shared_ptr<LikelihoodCalculationSingleProcess> likelihood_;
  std::shared_ptr<const ParametrizablePhyloTree> tree_;
  mutable std::shared_ptr<const Alphabet> alphabet_;
  size_t nbSites_;
  size_t nbDistinctSites_;
  // size_t nbClasses_;
  size_t nbStates_;
  PatternType rootPatternLinks_;
 
public:
  MarginalAncestralReconstruction(std::shared_ptr<LikelihoodCalculationSingleProcess> drl) :
    likelihood_      (drl),
    tree_            (drl->substitutionProcess().getParametrizablePhyloTree()),
    alphabet_        (drl->stateMap().getAlphabet()),
    nbSites_         (drl->getNumberOfSites()),
    nbDistinctSites_ (drl->getNumberOfDistinctSites()),
    nbStates_        (drl->stateMap().getNumberOfModelStates()),
    rootPatternLinks_(drl->getRootArrayPositions())
  {
    if (!tree_)
      throw Exception("MarginalAncestralReconstruction::MarginalAncestralReconstruction: missing ParametrizablePhyloTree.");
  }
 
  MarginalAncestralReconstruction(const MarginalAncestralReconstruction& masr) :
    likelihood_      (masr.likelihood_),
    tree_            (masr.tree_),
    alphabet_        (masr.alphabet_),
    nbSites_         (masr.nbSites_),
    nbDistinctSites_ (masr.nbDistinctSites_),
    nbStates_        (masr.nbStates_),
    rootPatternLinks_(masr.rootPatternLinks_)
  {}
 
  MarginalAncestralReconstruction& operator=(const MarginalAncestralReconstruction& masr)
  {
    likelihood_       = masr.likelihood_;
    tree_             = masr.tree_;
    alphabet_         = masr.alphabet_;
    nbSites_          = masr.nbSites_;
    nbDistinctSites_  = masr.nbDistinctSites_;
    nbStates_         = masr.nbStates_;
    rootPatternLinks_ = masr.rootPatternLinks_;
    return *this;
  }
 
 
  MarginalAncestralReconstruction* clone() const { return new MarginalAncestralReconstruction(*this); }
 
  virtual ~MarginalAncestralReconstruction() {}
 
public:
  std::shared_ptr<const Alphabet> getAlphabet() const
  {
    return alphabet_;
  }
 
  std::vector<size_t> getAncestralStatesForNode(uint nodeId, VVdouble& probs, bool sample) const;
 
  std::vector<size_t> getAncestralStatesForNode(uint nodeId) const override
  {
    VVdouble probs(nbSites_);
    return getAncestralStatesForNode(nodeId, probs, false);
  }
 
  std::map<uint, std::vector<size_t>> getAllAncestralStates() const override;
 
  std::unique_ptr<Sequence> getAncestralSequenceForNode(uint nodeId, VVdouble* probs, bool sample) const;
 
  std::unique_ptr<Sequence> getAncestralSequenceForNode(uint nodeId) const override
  {
    return getAncestralSequenceForNode(nodeId, 0, false);
  }
 
  std::unique_ptr<AlignedSequenceContainer> getAncestralSequences() const override
  {
    return getAncestralSequences(false);
  }
 
  std::unique_ptr<AlignedSequenceContainer> getAncestralSequences(bool sample) const;
 
private:
  void recursiveMarginalAncestralStates(
      const std::shared_ptr<PhyloNode> node,
      std::map<uint, std::vector<size_t>>& ancestors,
      AlignmentDataInterface& data) const;
};
} // end of namespace bpp.
#endif // BPP_PHYL_LIKELIHOOD_MARGINALANCESTRALRECONSTRUCTION_H