2
// File: OptimizationTools.h
3
// Created by: Julien Dutheil
4
// Created on: Sun Dec 14 09:43:32 2003
8
Copyright or Ā© or Copr. CNRS, (November 16, 2004)
10
This software is a computer program whose purpose is to provide classes
11
for phylogenetic data analysis.
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".
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
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.
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.
40
#ifndef _OPTIMIZATIONTOOLS_H_
41
#define _OPTIMIZATIONTOOLS_H_
43
#include "Likelihood/ClockTreeLikelihood.h"
44
#include "Likelihood/NNIHomogeneousTreeLikelihood.h"
45
#include "Likelihood/ClockTreeLikelihood.h"
46
#include "NNITopologySearch.h"
47
#include "Parsimony/DRTreeParsimonyScore.h"
48
#include "TreeTemplate.h"
49
#include "Distance/DistanceEstimation.h"
50
#include "Distance/AgglomerativeDistanceMethod.h"
52
#include <Bpp/Io/OutputStream.h>
53
#include <Bpp/App/ApplicationTools.h>
54
#include <Bpp/Numeric/Function/SimpleNewtonMultiDimensions.h>
61
* @brief Debugging wrapper: output data to a log file if a 0 likelihood is obtained during optimization.
63
class NaNWatcher: public DerivableSecondOrderWrapper
67
NaNWatcher(TreeLikelihood* tl) :
68
DerivableSecondOrderWrapper(tl) {}
70
NaNWatcher* clone() const { return new NaNWatcher(*this); }
73
* @brief This function is redefined to output an error message if there is a 0 likelihood.
75
* @return The value of the likelihood function;
77
double getValue() const throw (Exception);
84
* @brief Listener used internally by the optimizeTreeNNI method.
86
class NNITopologyListener:
87
public virtual TopologyListener
90
NNITopologySearch* topoSearch_;
91
ParameterList parameters_;
93
OutputStream* messenger_;
94
OutputStream* profiler_;
95
unsigned int verbose_;
96
unsigned int optimizeCounter_;
97
unsigned int optimizeNumerical_;
98
std::string optMethod_;
100
bool reparametrization_;
104
* @brief Build a new NNITopologyListener object.
106
* This listener listens to a NNITopologySearch object, and optimizes numerical parameters every *n* topological movements.
107
* Optimization is performed using the optimizeNumericalParameters method (see there documentation for more details).
109
* @param ts The NNITopologySearch object attached to this listener.
110
* @param parameters The list of parameters to optimize. Use tl->getIndependentParameters() in order to estimate all parameters.
111
* @param tolerance Tolerance to use during optimizaton.
112
* @param messenger Where to output messages.
113
* @param profiler Where to output optimization steps.
114
* @param verbose Verbose level during optimization.
115
* @param optMethod Optimization method to use.
116
* @param nStep The number of optimization steps to perform.
117
* @param reparametrization Tell if parameters should be transformed in order to remove constraints.
118
* This can improve optimization, but is a bit slower.
121
NNITopologySearch* ts,
122
const ParameterList& parameters,
124
OutputStream* messenger,
125
OutputStream* profiler,
126
unsigned int verbose,
127
const std::string& optMethod,
129
bool reparametrization) :
130
topoSearch_(ts), parameters_(parameters), tolerance_(tolerance),
131
messenger_(messenger), profiler_(profiler),
133
optimizeCounter_(0), optimizeNumerical_(1),
134
optMethod_(optMethod), nStep_(nStep), reparametrization_(reparametrization) {}
136
NNITopologyListener(const NNITopologyListener& tl) :
137
topoSearch_(tl.topoSearch_),
138
parameters_(tl.parameters_),
139
tolerance_(tl.tolerance_),
140
messenger_(tl.messenger_),
141
profiler_(tl.profiler_),
142
verbose_(tl.verbose_),
143
optimizeCounter_(tl.optimizeCounter_),
144
optimizeNumerical_(tl.optimizeNumerical_),
145
optMethod_(tl.optMethod_),
147
reparametrization_(tl.reparametrization_)
150
NNITopologyListener& operator=(const NNITopologyListener& tl)
152
topoSearch_ = tl.topoSearch_;
153
parameters_ = tl.parameters_;
154
tolerance_ = tl.tolerance_;
155
messenger_ = tl.messenger_;
156
profiler_ = tl.profiler_;
157
verbose_ = tl.verbose_;
158
optimizeCounter_ = tl.optimizeCounter_;
159
optimizeNumerical_ = tl.optimizeNumerical_;
160
optMethod_ = tl.optMethod_;
162
reparametrization_ = tl.reparametrization_;
166
NNITopologyListener* clone() const { return new NNITopologyListener(*this); }
168
virtual ~NNITopologyListener() {}
171
void topologyChangeTested(const TopologyChangeEvent& event) {}
172
void topologyChangeSuccessful(const TopologyChangeEvent& event);
173
void setNumericalOptimizationCounter(unsigned int c) { optimizeNumerical_ = c; }
178
* @brief Listener used internally by the optimizeTreeNNI2 method.
180
class NNITopologyListener2:
181
public TopologyListener
184
NNITopologySearch* topoSearch_;
185
ParameterList parameters_;
187
OutputStream* messenger_;
188
OutputStream* profiler_;
189
unsigned int verbose_;
190
unsigned int optimizeCounter_;
191
unsigned int optimizeNumerical_;
192
std::string optMethod_;
193
bool reparametrization_;
197
* @brief Build a new NNITopologyListener2 object.
199
* This listener listens to a NNITopologySearch object, and optimizes numerical parameters every *n* topological movements.
200
* Optimization is performed using the optimizeNumericalParameters2 method (see there documentation for more details).
202
* @param ts The NNITopologySearch object attached to this listener.
203
* @param parameters The list of parameters to optimize. Use ts->getIndependentParameters() in order to estimate all parameters.
204
* @param tolerance Tolerance to use during optimizaton.
205
* @param messenger Where to output messages.
206
* @param profiler Where to output optimization steps.
207
* @param verbose Verbose level during optimization.
208
* @param optMethod Optimization method to use.
209
* @param reparametrization Tell if parameters should be transformed in order to remove constraints.
210
* This can improve optimization, but is a bit slower.
212
NNITopologyListener2(
213
NNITopologySearch* ts,
214
const ParameterList& parameters,
216
OutputStream* messenger,
217
OutputStream* profiler,
218
unsigned int verbose,
219
const std::string& optMethod,
220
bool reparametrization) :
221
topoSearch_(ts), parameters_(parameters), tolerance_(tolerance),
222
messenger_(messenger), profiler_(profiler),
224
optimizeCounter_(0), optimizeNumerical_(1),
225
optMethod_(optMethod), reparametrization_(reparametrization) {}
227
NNITopologyListener2(const NNITopologyListener2& tl) :
228
topoSearch_(tl.topoSearch_),
229
parameters_(tl.parameters_),
230
tolerance_(tl.tolerance_),
231
messenger_(tl.messenger_),
232
profiler_(tl.profiler_),
233
verbose_(tl.verbose_),
234
optimizeCounter_(tl.optimizeCounter_),
235
optimizeNumerical_(tl.optimizeNumerical_),
236
optMethod_(tl.optMethod_),
237
reparametrization_(tl.reparametrization_)
240
NNITopologyListener2& operator=(const NNITopologyListener2& tl)
242
topoSearch_ = tl.topoSearch_;
243
parameters_ = tl.parameters_;
244
tolerance_ = tl.tolerance_;
245
messenger_ = tl.messenger_;
246
profiler_ = tl.profiler_;
247
verbose_ = tl.verbose_;
248
optimizeCounter_ = tl.optimizeCounter_;
249
optimizeNumerical_ = tl.optimizeNumerical_;
250
optMethod_ = tl.optMethod_;
251
reparametrization_ = tl.reparametrization_;
255
NNITopologyListener2* clone() const { return new NNITopologyListener2(*this); }
257
virtual ~NNITopologyListener2() {}
260
void topologyChangeTested(const TopologyChangeEvent& event) {}
261
void topologyChangeSuccessful(const TopologyChangeEvent& event);
262
void setNumericalOptimizationCounter(unsigned int c) { optimizeNumerical_ = c; }
273
* @brief Optimization methods for phylogenetic inference.
275
* This class provides optimization methods for phylogenetics.
276
* Parameters of the optimization methods are set to work with TreeLikelihood
277
* object. Some non trivial parameters are left to the user choice (tolerance, maximum
278
* number of function evaluation, output streams).
280
class OptimizationTools
284
virtual ~OptimizationTools();
288
static std::string OPTIMIZATION_GRADIENT;
289
static std::string OPTIMIZATION_NEWTON;
290
static std::string OPTIMIZATION_BRENT;
291
static std::string OPTIMIZATION_BFGS;
294
* @brief Optimize numerical parameters (branch length, substitution model & rate distribution) of a TreeLikelihood function.
296
* Uses Newton's method for branch length and Brent or BFGS one dimensional method for other parameters.
298
* A condition over function values is used as a stop condition for the algorithm.
300
* @see BrentOneDimension, BFGSMultiDimensions
302
* @param tl A pointer toward the TreeLikelihood object to optimize.
303
* @param parameters The list of parameters to optimize. Use tl->getIndependentParameters() in order to estimate all parameters.
304
* @param listener A pointer toward an optimization listener, if needed.
305
* @param nstep The number of progressive steps to perform (see NewtonBrentMetaOptimizer). 1 means full precision from start.
306
* @param tolerance The tolerance to use in the algorithm.
307
* @param tlEvalMax The maximum number of function evaluations.
308
* @param messageHandler The massage handler.
309
* @param profiler The profiler.
310
* @param reparametrization Tell if parameters should be transformed in order to remove constraints.
311
* This can improve optimization, but is a bit slower.
312
* @param verbose The verbose level.
313
* @param optMethodDeriv Optimization type for derivable parameters (first or second order derivatives).
314
* @see OPTIMIZATION_NEWTON, OPTIMIZATION_GRADIENT
315
* @param optMethodModel Optimization type for model parameters (Brent or BFGS).
316
* @see OPTIMIZATION_BRENT, OPTIMIZATION_BFGS
317
* @throw Exception any exception thrown by the Optimizer.
319
static unsigned int optimizeNumericalParameters(
320
DiscreteRatesAcrossSitesTreeLikelihood* tl,
321
const ParameterList& parameters,
322
OptimizationListener* listener = 0,
323
unsigned int nstep = 1,
324
double tolerance = 0.000001,
325
unsigned int tlEvalMax = 1000000,
326
OutputStream* messageHandler = ApplicationTools::message,
327
OutputStream* profiler = ApplicationTools::message,
328
bool reparametrization = false,
329
unsigned int verbose = 1,
330
const std::string& optMethodDeriv = OPTIMIZATION_NEWTON,
331
const std::string& optMethodModel = OPTIMIZATION_BFGS)
335
* @brief Optimize numerical parameters (branch length, substitution model & rate distribution) of a TreeLikelihood function.
337
* Uses Newton's method for all parameters, branch length derivatives are computed analytically, derivatives for other parameters numerically.
339
* @see PseudoNewtonOptimizer
341
* @param tl A pointer toward the TreeLikelihood object to optimize.
342
* @param parameters The list of parameters to optimize. Use tl->getIndependentParameters() in order to estimate all parameters.
343
* @param listener A pointer toward an optimization listener, if needed.
344
* @param tolerance The tolerance to use in the algorithm.
345
* @param tlEvalMax The maximum number of function evaluations.
346
* @param messageHandler The massage handler.
347
* @param profiler The profiler.
348
* @param reparametrization Tell if parameters should be transformed in order to remove constraints.
349
* This can improve optimization, but is a bit slower.
350
* @param verbose The verbose level.
351
* @param optMethodDeriv Optimization type for derivable parameters (first or second order derivatives).
352
* @see OPTIMIZATION_NEWTON, OPTIMIZATION_GRADIENT
353
* @throw Exception any exception thrown by the Optimizer.
355
static unsigned int optimizeNumericalParameters2(
356
DiscreteRatesAcrossSitesTreeLikelihood* tl,
357
const ParameterList& parameters,
358
OptimizationListener* listener = 0,
359
double tolerance = 0.000001,
360
unsigned int tlEvalMax = 1000000,
361
OutputStream* messageHandler = ApplicationTools::message,
362
OutputStream* profiler = ApplicationTools::message,
363
bool reparametrization = false,
364
unsigned int verbose = 1,
365
const std::string& optMethodDeriv = OPTIMIZATION_NEWTON)
369
* @brief Optimize branch lengths parameters of a TreeLikelihood function.
371
* Uses Newton's method.
373
* A condition over function values is used as a stop condition for the algorithm.
375
* @see NewtonBrentMetaOptimizer
377
* @param tl A pointer toward the TreeLikelihood object to optimize.
378
* @param parameters The list of parameters to optimize. The intersection of branch length parameters and the input set will be used. Use tl->getBranchLengthsParameters() in order to estimate all branch length parameters.
379
* @param listener A pointer toward an optimization listener, if needed.
380
* @param tolerance The tolerance to use in the algorithm.
381
* @param tlEvalMax The maximum number of function evaluations.
382
* @param messageHandler The massage handler.
383
* @param profiler The profiler.
384
* @param verbose The verbose level.
385
* @param optMethodDeriv Optimization type for derivable parameters (first or second order derivatives).
386
* @see OPTIMIZATION_NEWTON, OPTIMIZATION_GRADIENT
387
* @throw Exception any exception thrown by the Optimizer.
389
static unsigned int optimizeBranchLengthsParameters(
390
DiscreteRatesAcrossSitesTreeLikelihood * tl,
391
const ParameterList& parameters,
392
OptimizationListener* listener = 0,
393
double tolerance = 0.000001,
394
unsigned int tlEvalMax = 1000000,
395
OutputStream* messageHandler = ApplicationTools::message,
396
OutputStream* profiler = ApplicationTools::message,
397
unsigned int verbose = 1,
398
const std::string& optMethodDeriv = OPTIMIZATION_NEWTON)
402
* @brief Optimize numerical parameters assuming a global clock (branch heights, substitution model & rate distribution) of a ClockTreeLikelihood function.
404
* Uses Newton or conjugate gradient method for branch length and Brent's one dimensional method for other parameters
405
* (NewtonBrentMetaOptimizer).
406
* Derivatives are computed analytically.
408
* A condition over function values is used as a stop condition for the algorithm.
410
* @see NewtonBrentMetaOptimizer
412
* @param cl A pointer toward the ClockTreeLikelihood object to optimize.
413
* @param parameters The list of parameters to optimize. Use cl->getIndependentParameters() in order to estimate all parameters.
414
* @param listener A pointer toward an optimization listener, if needed.
415
* @param nstep The number of progressive steps to perform (see NewtonBrentMetaOptimizer). 1 means full precision from start.
416
* @param tolerance The tolerance to use in the algorithm.
417
* @param tlEvalMax The maximum number of function evaluations.
418
* @param messageHandler The massage handler.
419
* @param profiler The profiler.
420
* @param verbose The verbose level.
421
* @param optMethodDeriv Optimization type for derivable parameters (first or second order derivatives).
422
* @see OPTIMIZATION_NEWTON, OPTIMIZATION_GRADIENT
423
* @throw Exception any exception thrown by the Optimizer.
425
static unsigned int optimizeNumericalParametersWithGlobalClock(
426
DiscreteRatesAcrossSitesClockTreeLikelihood* cl,
427
const ParameterList& parameters,
428
OptimizationListener* listener = 0,
429
unsigned int nstep = 1,
430
double tolerance = 0.000001,
431
unsigned int tlEvalMax = 1000000,
432
OutputStream* messageHandler = ApplicationTools::message,
433
OutputStream* profiler = ApplicationTools::message,
434
unsigned int verbose = 1,
435
const std::string& optMethodDeriv = OPTIMIZATION_GRADIENT)
439
* @brief Optimize numerical parameters assuming a global clock (branch heights, substitution model & rate distribution) of a ClockTreeLikelihood function.
441
* Uses Newton or conjugate gradient method for all parameters, branch length derivatives are computed analytically, derivatives for other parameters numerically.
443
* @see PseudoNewtonOptimizer
445
* @param cl A pointer toward the ClockTreeLikelihood object to optimize.
446
* @param parameters The list of parameters to optimize. Use cl->getIndependentParameters() in order to estimate all parameters.
447
* @param listener A pointer toward an optimization listener, if needed.
448
* @param tolerance The tolerance to use in the algorithm.
449
* @param tlEvalMax The maximum number of function evaluations.
450
* @param messageHandler The massage handler.
451
* @param profiler The profiler.
452
* @param verbose The verbose level.
453
* @param optMethodDeriv Optimization type for derivable parameters (first or second order derivatives).
454
* @see OPTIMIZATION_NEWTON, OPTIMIZATION_GRADIENT
455
* @throw Exception any exception thrown by the Optimizer.
457
static unsigned int optimizeNumericalParametersWithGlobalClock2(
458
DiscreteRatesAcrossSitesClockTreeLikelihood* cl,
459
const ParameterList& parameters,
460
OptimizationListener* listener = 0,
461
double tolerance = 0.000001,
462
unsigned int tlEvalMax = 1000000,
463
OutputStream* messageHandler = ApplicationTools::message,
464
OutputStream* profiler = ApplicationTools::message,
465
unsigned int verbose = 1,
466
const std::string& optMethodDeriv = OPTIMIZATION_GRADIENT)
473
public virtual Function,
474
public ParametrizableAdapter
478
mutable ParameterList brLen_, lambda_;
481
ScaleFunction(TreeLikelihood* tl);
483
ScaleFunction(const ScaleFunction& sf) :
484
tl_(sf.tl_), brLen_(sf.brLen_), lambda_(sf.lambda_)
487
ScaleFunction& operator=(const ScaleFunction& sf)
491
lambda_ = sf.lambda_;
495
virtual ~ScaleFunction();
497
#ifndef NO_VIRTUAL_COV
502
clone() const { return new ScaleFunction(*this); }
505
void setParameters(const ParameterList& lambda) throw (ParameterNotFoundException, ConstraintException);
506
double getValue() const throw (ParameterException);
507
const ParameterList& getParameters() const throw (Exception) { return lambda_; }
508
const Parameter& getParameter(const std::string& name) const throw (ParameterNotFoundException)
510
if(name == "lambda") return lambda_[0];
511
else throw ParameterNotFoundException("ScaleFunction::getParameter.", name);
513
double getParameterValue(const std::string& name) const throw (ParameterNotFoundException)
515
return lambda_.getParameter(name).getValue();
517
unsigned int getNumberOfParameters() const { return 1; }
518
unsigned int getNumberOfIndependentParameters() const { return 1; }
524
* @brief Optimize the scale of a TreeLikelihood.
526
* This method only works on branch lengths parameters.
527
* It multiply all branch length by a factor 'x' which is optimized
528
* using Brent's algorithm in one dimension.
529
* This method may be usefull for scaling a tree whose branch lengths
530
* come from the Neighbor-Joining algorithm for instance.
532
* Practically, and contrarily to what one may expect, this does not
533
* speed up the optimization!
535
* A condition over parameters is used as a stop condition for the algorithm.
537
* @param tl A pointer toward the TreeLikelihood object to optimize.
538
* @param tolerance The tolerance to use in the algorithm.
539
* @param tlEvalMax The maximum number of function evaluations.
540
* @param messageHandler The massage handler.
541
* @param profiler The profiler.
542
* @param verbose The verbose level.
543
* @throw Exception any exception thrown by the optimizer.
545
static unsigned int optimizeTreeScale(
547
double tolerance = 0.000001,
548
int tlEvalMax = 1000000,
549
OutputStream* messageHandler = ApplicationTools::message,
550
OutputStream* profiler = ApplicationTools::message,
551
unsigned int verbose = 1)
555
* @brief Optimize all parameters from a TreeLikelihood object, including tree topology using Nearest Neighbor Interchanges.
557
* This function takes as input a TreeLikelihood object implementing the NNISearchable interface.
560
* A NNITopologySearch object is instanciated and is associated an additional TopologyListener.
561
* This listener is used to re-estimate numerical parameters after one or several topology change.
562
* By default, the PHYML option is used for the NNITopologySearch object, and numerical parameters are re-estimated
563
* every 4 NNI runs (as in the phyml software).
565
* The optimizeNumericalParameters method is used for estimating numerical parameters.
566
* The tolerance passed to this function is specified as input parameters.
567
* They are generally very high to avoid local optima.
569
* @param tl A pointer toward the TreeLikelihood object to optimize.
570
* @param parameters The list of parameters to optimize. Use tl->getIndependentParameters() in order to estimate all parameters.
571
* @param optimizeNumFirst Tell if we must optimize numerical parameters before searching topology.
572
* @param tolBefore The tolerance to use when estimating numerical parameters before topology search (if optimizeNumFirst is set to 'true').
573
* @param tolDuring The tolerance to use when estimating numerical parameters during the topology search.
574
* @param tlEvalMax The maximum number of function evaluations.
575
* @param numStep Number of NNI rounds before re-estimating numerical parameters.
576
* @param messageHandler The massage handler.
577
* @param profiler The profiler.
578
* @param reparametrization Tell if parameters should be transformed in order to remove constraints.
579
* This can improve optimization, but is a bit slower.
580
* @param verbose The verbose level.
581
* @param optMethod Option passed to optimizeNumericalParameters.
582
* @param nStep Option passed to optimizeNumericalParameters.
583
* @param nniMethod NNI algorithm to use.
584
* @return A pointer toward the final likelihood object.
585
* This pointer may be the same as passed in argument (tl), but in some cases the algorithm
586
* clone this object. We may change this bahavior in the future...
587
* You hence should write something like
589
* tl = OptimizationTools::optimizeTreeNNI(tl, ...);
591
* @throw Exception any exception thrown by the optimizer.
593
static NNIHomogeneousTreeLikelihood* optimizeTreeNNI(
594
NNIHomogeneousTreeLikelihood* tl,
595
const ParameterList& parameters,
596
bool optimizeNumFirst = true,
597
double tolBefore = 100,
598
double tolDuring = 100,
599
int tlEvalMax = 1000000,
600
unsigned int numStep = 1,
601
OutputStream* messageHandler = ApplicationTools::message,
602
OutputStream* profiler = ApplicationTools::message,
603
bool reparametrization = false,
604
unsigned int verbose = 1,
605
const std::string& optMethod = OptimizationTools::OPTIMIZATION_NEWTON,
606
unsigned int nStep = 1,
607
const std::string& nniMethod = NNITopologySearch::PHYML)
611
* @brief Optimize all parameters from a TreeLikelihood object, including tree topology using Nearest Neighbor Interchanges.
613
* This function takes as input a TreeLikelihood object implementing the NNISearchable interface.
616
* A NNITopologySearch object is instanciated and is associated an additional TopologyListener.
617
* This listener is used to re-estimate numerical parameters after one or several topology change.
618
* By default, the PHYML option is used for the NNITopologySearch object, and numerical parameters are re-estimated
619
* every 4 NNI runs (as in the phyml software).
621
* The optimizeNumericalParameters2 method is used for estimating numerical parameters.
622
* The tolerance passed to this function is specified as input parameters.
623
* They are generally very high to avoid local optima.
625
* @param tl A pointer toward the TreeLikelihood object to optimize.
626
* @param parameters The list of parameters to optimize. Use tl->getIndependentParameters() in order to estimate all parameters.
627
* @param optimizeNumFirst Tell if we must optimize numerical parameters before searching topology.
628
* @param tolBefore The tolerance to use when estimating numerical parameters before topology search (if optimizeNumFirst is set to 'true').
629
* @param tolDuring The tolerance to use when estimating numerical parameters during the topology search.
630
* @param tlEvalMax The maximum number of function evaluations.
631
* @param numStep Number of NNI rounds before re-estimating numerical parameters.
632
* @param messageHandler The massage handler.
633
* @param profiler The profiler.
634
* @param reparametrization Tell if parameters should be transformed in order to remove constraints.
635
* This can improve optimization, but is a bit slower.
636
* @param verbose The verbose level.
637
* @param optMethod Option passed to optimizeNumericalParameters2.
638
* @param nniMethod NNI algorithm to use.
639
* @return A pointer toward the final likelihood object.
640
* This pointer may be the same as passed in argument (tl), but in some cases the algorithm
641
* clone this object. We may change this bahavior in the future...
642
* You hence should write something like
644
* tl = OptimizationTools::optimizeTreeNNI2(tl, ...);
646
* @throw Exception any exception thrown by the optimizer.
648
static NNIHomogeneousTreeLikelihood* optimizeTreeNNI2(
649
NNIHomogeneousTreeLikelihood* tl,
650
const ParameterList& parameters,
651
bool optimizeNumFirst = true,
652
double tolBefore = 100,
653
double tolDuring = 100,
654
int tlEvalMax = 1000000,
655
unsigned int numStep = 1,
656
OutputStream* messageHandler = ApplicationTools::message,
657
OutputStream* profiler = ApplicationTools::message,
658
bool reparametrization = false,
659
unsigned int verbose = 1,
660
const std::string& optMethod = OptimizationTools::OPTIMIZATION_NEWTON,
661
const std::string& nniMethod = NNITopologySearch::PHYML)
665
* @brief Optimize tree topology from a DRTreeParsimonyScore using Nearest Neighbor Interchanges.
667
* @param tp A pointer toward the DRTreeParsimonyScore object to optimize.
668
* @param verbose The verbose level.
669
* @return A pointer toward the final parsimony score object.
670
* This pointer may be the same as passed in argument (tl), but in some cases the algorithm
671
* clone this object. We may change this bahavior in the future...
672
* You hence should write something like
674
* tp = OptimizationTools::optimizeTreeNNI(tp, ...);
678
static DRTreeParsimonyScore* optimizeTreeNNI(
679
DRTreeParsimonyScore* tp,
680
unsigned int verbose = 1);
683
* @brief Build a tree using a distance method.
685
* This method estimate a distance matrix using a DistanceEstimation object, and then builds the phylogenetic tree using a AgglomerativeDistanceMethod object.
686
* The main issue here is to estimate non-branch lengths parameters, as substitution model and rate distribution parameters.
687
* Three options are provideed here:
688
* - DISTANCEMETHOD_INIT (default) keep parameters to there initial value,
689
* - DISTANCEMETHOD_PAIRWISE estimated parameters in a pairwise manner, which is standard but not good at all...
690
* - DISTANCEMETHOD_ITERATIONS uses Ninio et al's iterative algorithm, which uses Maximum Likelihood to estimate these parameters, and then update the distance matrix.
691
* Ninio M, Privman E, Pupko T, Friedman N.
692
* Phylogeny reconstruction: increasing the accuracy of pairwise distance estimation using Bayesian inference of evolutionary rates.
693
* Bioinformatics. 2007 Jan 15;23(2):e136-41.
695
* @param estimationMethod The distance estimation object to use.
696
* @param reconstructionMethod The tree reconstruction object to use.
697
* @param parametersToIgnore A list of parameters to ignore while optimizing parameters.
698
* @param optimizeBrLen Tell if branch lengths should be optimized together with other parameters. This may lead to more accurate parameter estimation, but is slower.
699
* @param rooted Tell if rooted trees must be produced.
700
* @param param String describing the type of optimization to use.
701
* @param tolerance Threshold on likelihood for stopping the iterative procedure. Used only with param=DISTANCEMETHOD_ITERATIONS.
702
* @param tlEvalMax Maximum number of likelihood computations to perform when optimizing parameters. Used only with param=DISTANCEMETHOD_ITERATIONS.
703
* @param profiler Output stream used by optimizer. Used only with param=DISTANCEMETHOD_ITERATIONS.
704
* @param messenger Output stream used by optimizer. Used only with param=DISTANCEMETHOD_ITERATIONS.
705
* @param verbose Verbose level.
707
static TreeTemplate<Node>* buildDistanceTree(
708
DistanceEstimation& estimationMethod,
709
AgglomerativeDistanceMethod& reconstructionMethod,
710
const ParameterList& parametersToIgnore,
711
bool optimizeBrLen = false,
713
const std::string& param = DISTANCEMETHOD_INIT,
714
double tolerance = 0.000001,
715
unsigned int tlEvalMax = 1000000,
716
OutputStream* profiler = 0,
717
OutputStream* messenger = 0,
718
unsigned int verbose = 0) throw (Exception);
721
static std::string DISTANCEMETHOD_INIT;
722
static std::string DISTANCEMETHOD_PAIRWISE;
723
static std::string DISTANCEMETHOD_ITERATIONS;
727
} //end of namespace bpp.
729
#endif //_OPTIMIZATIONTOOLS_H_