~ubuntu-branches/ubuntu/precise/libbpp-phyl/precise

« back to all changes in this revision

Viewing changes to src/Bpp/Phyl/Likelihood/AbstractNonHomogeneousTreeLikelihood.h

  • Committer: Bazaar Package Importer
  • Author(s): Julien Dutheil
  • Date: 2011-06-09 11:00:00 UTC
  • Revision ID: james.westby@ubuntu.com-20110609110000-yvx78svv6w7xxgph
Tags: upstream-2.0.2
ImportĀ upstreamĀ versionĀ 2.0.2

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
//
 
2
// File: AbstractNonHomogeneousTreeLikelihood.h
 
3
// Created by: Julien Dutheil
 
4
// Created on: Tue Oct 09 16:07 2007
 
5
// From file: AbstractHomogeneousTreeLikelihood.h
 
6
//
 
7
 
 
8
/*
 
9
Copyright or Ā© or Copr. CNRS, (November 16, 2004)
 
10
 
 
11
This software is a computer program whose purpose is to provide classes
 
12
for phylogenetic data analysis.
 
13
 
 
14
This software is governed by the CeCILL  license under French law and
 
15
abiding by the rules of distribution of free software.  You can  use, 
 
16
modify and/ or redistribute the software under the terms of the CeCILL
 
17
license as circulated by CEA, CNRS and INRIA at the following URL
 
18
"http://www.cecill.info". 
 
19
 
 
20
As a counterpart to the access to the source code and  rights to copy,
 
21
modify and redistribute granted by the license, users are provided only
 
22
with a limited warranty  and the software's author,  the holder of the
 
23
economic rights,  and the successive licensors  have only  limited
 
24
liability. 
 
25
 
 
26
In this respect, the user's attention is drawn to the risks associated
 
27
with loading,  using,  modifying and/or developing or reproducing the
 
28
software by the user in light of its specific status of free software,
 
29
that may mean  that it is complicated to manipulate,  and  that  also
 
30
therefore means  that it is reserved for developers  and  experienced
 
31
professionals having in-depth computer knowledge. Users are therefore
 
32
encouraged to load and test the software's suitability as regards their
 
33
requirements in conditions enabling the security of their systems and/or 
 
34
data to be ensured and,  more generally, to use and operate it in the 
 
35
same conditions as regards security. 
 
36
 
 
37
The fact that you are presently reading this means that you have had
 
38
knowledge of the CeCILL license and that you accept its terms.
 
39
*/
 
40
 
 
41
#ifndef _ABSTRACTNONHOMOGENEOUSTREELIKELIHOOD_H_
 
42
#define _ABSTRACTNONHOMOGENEOUSTREELIKELIHOOD_H_
 
43
 
 
44
#include "NonHomogeneousTreeLikelihood.h"
 
45
#include "AbstractDiscreteRatesAcrossSitesTreeLikelihood.h"
 
46
 
 
47
//From the STL:
 
48
#include <memory>
 
49
 
 
50
namespace bpp
 
51
{
 
52
 
 
53
/**
 
54
 * @brief Partial implementation for branch non-homogeneous models of the TreeLikelihood interface.
 
55
 *
 
56
 * This class provides a pointer toward a single substitution model + several utilitary variables.
 
57
 */
 
58
class AbstractNonHomogeneousTreeLikelihood:
 
59
  public virtual NonHomogeneousTreeLikelihood,
 
60
  public AbstractDiscreteRatesAcrossSitesTreeLikelihood
 
61
{
 
62
  public:
 
63
 
 
64
    class ConstNonHomogeneousSiteModelIterator :
 
65
      public ConstSiteModelIterator
 
66
    {
 
67
      private:
 
68
        std::vector<ConstNoPartitionSiteModelDescription> siteModelDescriptions_;
 
69
        unsigned int index_;
 
70
        unsigned int nbModels_;
 
71
 
 
72
      public:
 
73
        ConstNonHomogeneousSiteModelIterator(const SubstitutionModelSet* modelSet) :
 
74
          siteModelDescriptions_(), index_(0), nbModels_(modelSet->getNumberOfModels())
 
75
        {
 
76
          for (unsigned int i = 0; i < nbModels_; ++i)
 
77
            siteModelDescriptions_.push_back(ConstNoPartitionSiteModelDescription(modelSet->getModel(i), modelSet->getNodesWithModel(i)));        
 
78
        }
 
79
 
 
80
      public:
 
81
        ConstSiteModelDescription* next() throw (Exception)
 
82
        {
 
83
          if (!hasNext())
 
84
            throw Exception("AbstractNonHomogeneousTreeLikelihood::ConstHomogeneousSiteModelIterator::next(). No more site in the set.");
 
85
          return &siteModelDescriptions_[index_++];
 
86
        }
 
87
 
 
88
        bool hasNext() const { return index_ < nbModels_; }
 
89
    };
 
90
 
 
91
  protected:
 
92
    SubstitutionModelSet* modelSet_;
 
93
    ParameterList brLenParameters_;
 
94
    
 
95
    mutable std::map<int, VVVdouble> pxy_;
 
96
 
 
97
    mutable std::map<int, VVVdouble> dpxy_;
 
98
 
 
99
    mutable std::map<int, VVVdouble> d2pxy_;
 
100
        
 
101
    std::vector<double> rootFreqs_;
 
102
        
 
103
    /**
 
104
     * @brief Pointer toward all nodes in the tree.
 
105
     *
 
106
     * The position in the array is the number used in the parameter name.
 
107
     * This may be different from the node id, unless you used the resetNodeId method on the input tree.
 
108
     */
 
109
    std::vector<Node*> nodes_;
 
110
 
 
111
    /**
 
112
     * @brief An index linking nodes to their id, for faster access than the getNode() method.
 
113
     */
 
114
    mutable std::map<int, const Node*> idToNode_;
 
115
 
 
116
    //some values we'll need:
 
117
    unsigned int nbSites_,         //the number of sites in the container
 
118
                 nbDistinctSites_, //the number of distinct sites
 
119
                 nbClasses_,       //the number of rate classes
 
120
                 nbStates_,        //the number of states in the alphabet
 
121
                 nbNodes_;         //the number of nodes in the tree
 
122
 
 
123
    bool verbose_;
 
124
 
 
125
    double minimumBrLen_;
 
126
    std::auto_ptr<Constraint> brLenConstraint_;
 
127
 
 
128
    int root1_, root2_;
 
129
 
 
130
  public:
 
131
    AbstractNonHomogeneousTreeLikelihood(
 
132
      const Tree& tree,
 
133
      SubstitutionModelSet* modelSet,
 
134
      DiscreteDistribution* rDist,
 
135
      bool verbose = true)
 
136
      throw (Exception);
 
137
 
 
138
    /**
 
139
     * @brief Copy constructor
 
140
     *
 
141
     * This constructor is to be called by the derived class copy constructor.
 
142
     */
 
143
    AbstractNonHomogeneousTreeLikelihood(const AbstractNonHomogeneousTreeLikelihood& lik);
 
144
    
 
145
    /**
 
146
     * @brief Assignation operator
 
147
     *
 
148
     * This operator is to be called by the derived class operator.
 
149
     */
 
150
    AbstractNonHomogeneousTreeLikelihood& operator=(const AbstractNonHomogeneousTreeLikelihood& lik);
 
151
 
 
152
    virtual ~AbstractNonHomogeneousTreeLikelihood() {}
 
153
    
 
154
  private:
 
155
 
 
156
    /**
 
157
     * @brief Method called by constructor.
 
158
     */
 
159
    void init_(
 
160
        const Tree& tree,
 
161
        SubstitutionModelSet* modelSet,
 
162
        DiscreteDistribution* rDist,
 
163
        bool verbose) throw (Exception);
 
164
 
 
165
  public:
 
166
    
 
167
    /**
 
168
     * @name The TreeLikelihood interface.
 
169
     *
 
170
     * Other methods are implemented in the AbstractTreeLikelihood class.
 
171
     *
 
172
     * @{
 
173
     */
 
174
    void initialize() throw(Exception);
 
175
    
 
176
    ParameterList getBranchLengthsParameters() const;
 
177
    
 
178
    ParameterList getSubstitutionModelParameters() const;
 
179
    
 
180
    ParameterList getRateDistributionParameters() const
 
181
    {
 
182
      return AbstractDiscreteRatesAcrossSitesTreeLikelihood::getRateDistributionParameters();
 
183
    }
 
184
 
 
185
    const SubstitutionModel* getSubstitutionModelForNode(int nodeId) const throw (NodeNotFoundException) 
 
186
    {
 
187
      return modelSet_->getModelForNode(nodeId);
 
188
    }
 
189
 
 
190
    SubstitutionModel* getSubstitutionModelForNode(int nodeId) throw (NodeNotFoundException)
 
191
    {
 
192
      return modelSet_->getModelForNode(nodeId);
 
193
    }
 
194
 
 
195
    const std::vector<double>& getRootFrequencies(unsigned int siteIndex) const { return rootFreqs_; }
 
196
    
 
197
    VVVdouble getTransitionProbabilitiesPerRateClass(int nodeId, unsigned int siteIndex) const { return pxy_[nodeId]; }
 
198
 
 
199
    ConstBranchModelIterator* getNewBranchModelIterator(int nodeId) const
 
200
    {
 
201
      return new ConstNoPartitionBranchModelIterator(modelSet_->getModelForNode(nodeId), nbDistinctSites_);
 
202
    }
 
203
 
 
204
    ConstSiteModelIterator* getNewSiteModelIterator(unsigned int siteIndex) const
 
205
    {
 
206
      return new ConstNonHomogeneousSiteModelIterator(modelSet_);
 
207
    }
 
208
       
 
209
    /** @} */
 
210
 
 
211
    /**
 
212
     * @name The NonHomogeneousTreeLikelihood interface.
 
213
     *
 
214
     * Other methods are implemented in the AbstractTreeLikelihood class.
 
215
     *
 
216
     * @{
 
217
     */
 
218
    const SubstitutionModelSet* getSubstitutionModelSet() const { return modelSet_; }
 
219
    
 
220
    SubstitutionModelSet* getSubstitutionModelSet() { return modelSet_; }
 
221
    
 
222
    void setSubstitutionModelSet(SubstitutionModelSet* modelSet) throw (Exception);
 
223
 
 
224
    ParameterList getRootFrequenciesParameters() const
 
225
    {
 
226
      return modelSet_->getRootFrequenciesParameters();
 
227
    }
 
228
    /** @} */
 
229
    
 
230
  public: //Specific methods:
 
231
 
 
232
    /**
 
233
     * @brief This builds the <i>parameters</i> list from all parametrizable objects,
 
234
     * <i>i.e.</i> substitution model, rate distribution and tree.
 
235
     */
 
236
    virtual void initParameters();
 
237
 
 
238
    /**
 
239
     * @brief All parameters are stored in a parameter list.
 
240
     * This function apply these parameters to the substitution model,
 
241
     * to the rate distribution and to the branch lengths.
 
242
     */
 
243
    virtual void applyParameters() throw (Exception);  
 
244
 
 
245
    virtual void initBranchLengthsParameters();
 
246
 
 
247
    virtual void setMinimumBranchLength(double minimum)
 
248
    {
 
249
      minimumBrLen_ = minimum;
 
250
      if (brLenConstraint_.get()) brLenConstraint_.release();
 
251
      brLenConstraint_.reset(new IncludingPositiveReal(minimumBrLen_));
 
252
      initBranchLengthsParameters();
 
253
    }
 
254
 
 
255
    virtual double getMinimumBranchLength() const { return minimumBrLen_; }
 
256
 
 
257
  protected:
 
258
    /**
 
259
     * @brief Fill the pxy_, dpxy_ and d2pxy_ arrays for all nodes.
 
260
     */
 
261
    virtual void computeAllTransitionProbabilities();
 
262
    /**
 
263
     * @brief Fill the pxy_, dpxy_ and d2pxy_ arrays for one node.
 
264
     */
 
265
    virtual void computeTransitionProbabilitiesForNode(const Node * node);
 
266
 
 
267
};
 
268
 
 
269
} //end of namespace bpp.
 
270
 
 
271
#endif //_ABSTRACTNONHOMOGENEOUSTREELIKELIHOOD_H_
 
272