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

« back to all changes in this revision

Viewing changes to src/Bpp/Phyl/Likelihood/TreeLikelihood.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: TreeLikelihood.h
 
3
// Created by: Julien Dutheil
 
4
// Created on: Fri Oct 17 17:36:44 2003
 
5
//
 
6
 
 
7
/*
 
8
Copyright or Ā© or Copr. CNRS, (November 16, 2004)
 
9
 
 
10
This software is a computer program whose purpose is to provide classes
 
11
for phylogenetic data analysis.
 
12
 
 
13
This software is governed by the CeCILL  license under French law and
 
14
abiding by the rules of distribution of free software.  You can  use, 
 
15
modify and/ or redistribute the software under the terms of the CeCILL
 
16
license as circulated by CEA, CNRS and INRIA at the following URL
 
17
"http://www.cecill.info". 
 
18
 
 
19
As a counterpart to the access to the source code and  rights to copy,
 
20
modify and redistribute granted by the license, users are provided only
 
21
with a limited warranty  and the software's author,  the holder of the
 
22
economic rights,  and the successive licensors  have only  limited
 
23
liability. 
 
24
 
 
25
In this respect, the user's attention is drawn to the risks associated
 
26
with loading,  using,  modifying and/or developing or reproducing the
 
27
software by the user in light of its specific status of free software,
 
28
that may mean  that it is complicated to manipulate,  and  that  also
 
29
therefore means  that it is reserved for developers  and  experienced
 
30
professionals having in-depth computer knowledge. Users are therefore
 
31
encouraged to load and test the software's suitability as regards their
 
32
requirements in conditions enabling the security of their systems and/or 
 
33
data to be ensured and,  more generally, to use and operate it in the 
 
34
same conditions as regards security. 
 
35
 
 
36
The fact that you are presently reading this means that you have had
 
37
knowledge of the CeCILL license and that you accept its terms.
 
38
*/
 
39
 
 
40
#ifndef _TREELIKELIHOOD_H_
 
41
#define _TREELIKELIHOOD_H_
 
42
 
 
43
#include "../Node.h"
 
44
#include "../Tree.h"
 
45
#include "../Model/SubstitutionModel.h"
 
46
#include "TreeLikelihoodData.h"
 
47
 
 
48
#include <Bpp/Numeric/ParameterList.h>
 
49
#include <Bpp/Numeric/Parametrizable.h>
 
50
#include <Bpp/Numeric/Function/Functions.h>
 
51
#include <Bpp/Numeric/VectorTools.h>
 
52
 
 
53
// From SeqLib:
 
54
#include <Bpp/Seq/Alphabet/Alphabet.h>
 
55
#include <Bpp/Seq/Container/SiteContainer.h>
 
56
 
 
57
namespace bpp
 
58
{
 
59
 
 
60
/**
 
61
 * @brief The TreeLikelihood interface.
 
62
 *
 
63
 * This interface defines the methods needed for computing the likelihood
 
64
 * of a phylogenetic tree, given a dataset.
 
65
 */ 
 
66
class TreeLikelihood:
 
67
  public virtual DerivableSecondOrder
 
68
{
 
69
  public:
 
70
    /**
 
71
     * @brief An iterator over a set of branches, specified by their node ids.
 
72
     */
 
73
    class BranchIterator
 
74
    {
 
75
      public:
 
76
        virtual ~BranchIterator() {}
 
77
 
 
78
      public:
 
79
        /**
 
80
         * @return The id of the next node in the set.
 
81
         */
 
82
        virtual int next() throw (Exception) = 0;
 
83
        /**
 
84
         * @return True if there is at least another node in the set.
 
85
         */
 
86
        virtual bool hasNext() const = 0;
 
87
    };
 
88
 
 
89
    /**
 
90
     * @brief An iterator over a set of sites, speicfied by their position.
 
91
     *
 
92
     * In most cases, the position will reflect the index of an inner array used for likelihood storage.
 
93
     */
 
94
    class SiteIterator
 
95
    {
 
96
      public:
 
97
        virtual ~SiteIterator() {}
 
98
 
 
99
      public:
 
100
        /**
 
101
         * @return The position of the next site in the set.
 
102
         */
 
103
        virtual unsigned int next() throw (Exception) = 0;
 
104
        /**
 
105
         * @return True is there is at least another site in the set.
 
106
         */
 
107
        virtual bool hasNext() const = 0;
 
108
    };
 
109
 
 
110
    /**
 
111
     * @brief A pair of SubstitutionModel / SiteIterator.
 
112
     */
 
113
    class ConstBranchModelDescription
 
114
    {
 
115
      public:
 
116
        virtual ~ConstBranchModelDescription() {}
 
117
 
 
118
      public:
 
119
        virtual const SubstitutionModel* getModel() const = 0;
 
120
        virtual SiteIterator* getNewSiteIterator() const = 0;
 
121
    };
 
122
 
 
123
    /**
 
124
     * @brief Iterates through all models used for all sites on a given branch.
 
125
     */
 
126
    class ConstBranchModelIterator
 
127
    {
 
128
      public:
 
129
        virtual ~ConstBranchModelIterator() {}
 
130
 
 
131
      public:
 
132
        virtual ConstBranchModelDescription* next() throw (Exception) = 0;
 
133
        virtual bool hasNext() const = 0;
 
134
    };
 
135
 
 
136
    /**
 
137
     * @brief A pair of SubstitutionModel / BranchIterator.
 
138
     */
 
139
    class ConstSiteModelDescription
 
140
    {
 
141
      public:
 
142
        virtual ~ConstSiteModelDescription() {}
 
143
 
 
144
      public:
 
145
        virtual const SubstitutionModel* getModel() const = 0;
 
146
        virtual BranchIterator* getNewBranchIterator() const = 0;
 
147
    };
 
148
 
 
149
    /**
 
150
     * @brief Iterates through all models used for all branches on a given site.
 
151
     */
 
152
    class ConstSiteModelIterator
 
153
    {
 
154
      public:
 
155
        virtual ~ConstSiteModelIterator() {}
 
156
 
 
157
      public:
 
158
        virtual ConstSiteModelDescription* next() throw (Exception) = 0;
 
159
        virtual bool hasNext() const = 0;
 
160
    };
 
161
 
 
162
  public:
 
163
    TreeLikelihood() {}
 
164
    virtual ~TreeLikelihood() {}
 
165
 
 
166
#ifndef NO_VIRTUAL_COV
 
167
    TreeLikelihood* clone() const = 0;
 
168
#endif
 
169
 
 
170
  public:
 
171
 
 
172
    /**
 
173
     * @brief Set the dataset for which the likelihood must be evaluated.
 
174
     *
 
175
     * @param sites The data set to use.
 
176
     */
 
177
    virtual void setData(const SiteContainer& sites) = 0;
 
178
    
 
179
    /**
 
180
     * @brief Get the dataset for which the likelihood must be evaluated.
 
181
     *
 
182
     * @return A pointer toward the site container where the sequences are stored.
 
183
     */
 
184
    virtual const SiteContainer* getData() const = 0;
 
185
 
 
186
    /**
 
187
     * @brief Init the likelihood object.
 
188
     *
 
189
     * This method is used to initialize all parameters.
 
190
     * It is typically called after the constructor and the setData method.
 
191
     * It contains virtual methods that can't be called in the constructor.
 
192
     * @throw Exception if something bad happened, for instance if no data are associated to the likelihood function.
 
193
     */
 
194
    virtual void initialize() throw (Exception) = 0;
 
195
 
 
196
    /**
 
197
     * @return 'true' is the likelihood function has been initialized.
 
198
     */
 
199
    virtual bool isInitialized() const = 0;
 
200
 
 
201
    /**
 
202
     * @return The underlying likelihood data structure.
 
203
     */
 
204
    virtual TreeLikelihoodData* getLikelihoodData() = 0;
 
205
 
 
206
    /**
 
207
     * @return The underlying likelihood data structure.
 
208
     */
 
209
    virtual const TreeLikelihoodData* getLikelihoodData() const = 0;
 
210
 
 
211
    /**
 
212
     * @brief Get the likelihood for a site.
 
213
     *
 
214
     * @param site The site index to analyse.
 
215
     * @return The likelihood for site <i>site</i>.
 
216
     */
 
217
    virtual double getLikelihoodForASite(unsigned int site) const = 0;
 
218
 
 
219
    /**
 
220
     * @brief Get the logarithm of the likelihood for a site.
 
221
     *
 
222
     * @param site The site index to analyse.
 
223
     * @return The logarithm of the likelihood for site <i>site</i>.
 
224
     */
 
225
    virtual double getLogLikelihoodForASite(unsigned int site) const = 0;
 
226
 
 
227
    /**
 
228
     * @brief Get the likelihood for a site and for a state.
 
229
     *
 
230
     * @param site The site index to analyse.
 
231
     * @param state The state to consider.
 
232
     * @return The likelihood for site <i>site</i> and state <i>state</i>.
 
233
     */
 
234
    virtual double getLikelihoodForASiteForAState(unsigned int site, int state) const = 0;
 
235
 
 
236
    /**
 
237
     * @brief Get the logarithm of the likelihood for a site and for a state.
 
238
     *
 
239
     * @param site The site index to analyse.
 
240
     * @param state The state to consider.
 
241
     * @return The logarithm of the likelihood for site <i>site</i> and state <i>state</i>.
 
242
     */
 
243
    virtual double getLogLikelihoodForASiteForAState(unsigned int site, int state) const = 0;
 
244
 
 
245
    /**
 
246
     * @brief Get the likelihood for each site.
 
247
     *
 
248
     * @return A vector with all likelihoods for each site.
 
249
     */
 
250
    virtual Vdouble getLikelihoodForEachSite() const = 0;
 
251
 
 
252
    /**
 
253
     * @brief Get the logarithm of the likelihood for each site.
 
254
     *
 
255
     * @return A vector with all log likelihoods for each site.
 
256
     */
 
257
    virtual Vdouble getLogLikelihoodForEachSite() const = 0;
 
258
 
 
259
    /**
 
260
     * @brief Get the likelihood for each site and for each state.
 
261
     *
 
262
     * @return A 2d vector with all likelihoods for each site and for each state.
 
263
     */
 
264
    virtual VVdouble getLikelihoodForEachSiteForEachState() const = 0;
 
265
 
 
266
    /**
 
267
     * @brief Get the logarithm of the likelihood for each site and for each state.
 
268
     *
 
269
     * @return A 2d vector with all log likelihoods for each site and for each state.
 
270
     */
 
271
    virtual VVdouble getLogLikelihoodForEachSiteForEachState() const = 0;
 
272
    
 
273
    /**
 
274
     * @brief Get the likelihood for the whole dataset.
 
275
     *
 
276
     * @return The likelihood of the dataset.
 
277
     */
 
278
    virtual double getLikelihood() const = 0;
 
279
 
 
280
    /**
 
281
     * @brief Get the logarithm of the likelihood for the whole dataset.
 
282
     *
 
283
     * @return The logarithm of the likelihood of the dataset.
 
284
     */
 
285
    virtual double getLogLikelihood() const = 0;
 
286
  
 
287
    /**
 
288
     * @brief Get the tree (topology and branch lengths).
 
289
     *
 
290
     * @return The tree of this TreeLikelihood object.
 
291
      */
 
292
    virtual const Tree& getTree() const = 0;
 
293
 
 
294
    /**
 
295
     * @brief Get the number of sites in the dataset.
 
296
     *
 
297
     * @return the number of sites in the dataset.
 
298
     */
 
299
    virtual unsigned int getNumberOfSites() const = 0;
 
300
 
 
301
    /**
 
302
     * @brief Get the number of states in the alphabet associated to the dataset.
 
303
     *
 
304
     * @return the number of states in the alphabet associated to the dataset.
 
305
     */    
 
306
    virtual unsigned int getNumberOfStates() const = 0;
 
307
    
 
308
    /**
 
309
     * @brief Get the alphabet associated to the dataset.
 
310
     *
 
311
     * @return the alphabet associated to the dataset.
 
312
     */    
 
313
    virtual const Alphabet* getAlphabet() const = 0;
 
314
   
 
315
    /**
 
316
     * @name Retrieve some particular parameters subsets.
 
317
     *
 
318
     * @{
 
319
     */
 
320
    
 
321
    /**
 
322
     * @brief Get the branch lengths parameters.
 
323
     *
 
324
     * @return A ParameterList with all branch lengths.
 
325
     */
 
326
    virtual ParameterList getBranchLengthsParameters() const = 0;
 
327
    
 
328
    /**
 
329
     * @brief Get the parameters associated to substitution model(s).
 
330
     *
 
331
     * @return A ParameterList.
 
332
     */
 
333
    virtual ParameterList getSubstitutionModelParameters() const = 0;
 
334
 
 
335
    /**
 
336
     * @brief Get the substitution model associated to a given node and alignment column.
 
337
     *
 
338
     * @param nodeId The id of the request node.
 
339
     * @param siteIndex The index of the alignment position.
 
340
     * @see getSiteIndex
 
341
     * @return A pointer toward the corresponding model.
 
342
     * @throw NodeNotFoundException This exception may be thrown if the node is not found (depending on the implementation).
 
343
     */
 
344
    virtual const SubstitutionModel* getSubstitutionModel(int nodeId, unsigned int siteIndex) const throw (NodeNotFoundException) = 0;
 
345
 
 
346
    /**
 
347
     * @brief Get the substitution model associated to a given node and alignment column.
 
348
     *
 
349
     * @param nodeId The id of the request node.
 
350
     * @param siteIndex The index of the alignment position.
 
351
     * @see getSiteIndex
 
352
     * @return A pointer toward the corresponding model.
 
353
     * @throw NodeNotFoundException This exception may be thrown if the node is not found (depending on the implementation).
 
354
     */
 
355
    virtual SubstitutionModel* getSubstitutionModel(int nodeId, unsigned int siteIndex) throw (NodeNotFoundException) = 0;
 
356
 
 
357
    /**
 
358
     * @brief Retrieves all Pij(t) for a particular branch, defined by the upper node and site.
 
359
     *
 
360
     * These intermediate results may be used by other methods.
 
361
     *
 
362
     * @param nodeId The node defining the branch of interest.
 
363
     * @param siteIndex The index of the alignment position.
 
364
     * @see getSiteIndex
 
365
     * @return An array of dimension 2, where a[x][y] is the probability of substituting from x to y.
 
366
     */
 
367
    virtual VVdouble getTransitionProbabilities(int nodeId, unsigned int siteIndex) const = 0;
 
368
 
 
369
    virtual ConstBranchModelIterator* getNewBranchModelIterator(int nodeId) const = 0;
 
370
 
 
371
    virtual ConstSiteModelIterator* getNewSiteModelIterator(unsigned int siteIndex) const = 0;
 
372
 
 
373
    /**
 
374
     * @brief Get the index (used for inner computations) of a given site (original alignment column).
 
375
     *
 
376
     * @param site An alignment position.
 
377
     * @return The site index corresponding to the given input alignment position.
 
378
     */
 
379
    virtual unsigned int getSiteIndex(unsigned int site) const throw (IndexOutOfBoundsException) = 0;
 
380
 
 
381
    /**
 
382
     * @brief Get the values of the frequencies for each state in the alphabet at the root node.
 
383
     *
 
384
     * For reversible models, these are the equilibrium frequencies.
 
385
     * For non-reversible models, these usually are distinct parameters.
 
386
     *
 
387
     * For models without site partitioning, the set of frequencies is the same for all positions.
 
388
     * For partition models, the frequencies may differ from one site to another.
 
389
     *
 
390
     * @param siteIndex The index of the alignment position.
 
391
     * @see getSiteIndex
 
392
     * @return A vector with ancestral frequencies for each state in the alphabet;
 
393
     */
 
394
    virtual const std::vector<double>& getRootFrequencies(unsigned int siteIndex) const = 0;
 
395
    
 
396
    /** @} */
 
397
 
 
398
    /**
 
399
     * @brief Tell if derivatives must be computed.
 
400
     *
 
401
     * This methods calls the enableFirstOrderDerivatives and enableSecondOrderDerivatives.
 
402
     *
 
403
     * @param yn Yes or no.
 
404
     */
 
405
    virtual void enableDerivatives(bool yn) = 0;
 
406
 
 
407
    /**
 
408
     * @brief All derivable parameters.
 
409
     *
 
410
     * Usually, this contains all branch lengths parameters.
 
411
     *
 
412
     * @return A ParameterList.
 
413
     */
 
414
    virtual ParameterList getDerivableParameters() const = 0;
 
415
 
 
416
    /**
 
417
     * @brief All non derivable parameters.
 
418
     *
 
419
     * Usually, this contains all substitution model parameters and rate distribution.
 
420
     *
 
421
     * @return A ParameterList.
 
422
     */
 
423
    virtual ParameterList getNonDerivableParameters() const = 0;
 
424
 
 
425
};
 
426
 
 
427
} //end of namespace bpp.
 
428
 
 
429
#endif  //_TREELIKELIHOOD_H_
 
430