unsigned int OptimizationTools::optimizeNumericalParametersWithGlobalClock2( DiscreteRatesAcrossSitesClockTreeLikelihood* cl, const ParameterList& parameters, OptimizationListener* listener, double tolerance, unsigned int tlEvalMax, OutputStream* messageHandler, OutputStream* profiler, unsigned int verbose, const std::string& optMethodDeriv) throw (Exception) { AbstractNumericalDerivative* fun = 0; // Build optimizer: Optimizer* optimizer = 0; if (optMethodDeriv == OPTIMIZATION_GRADIENT) { fun = new TwoPointsNumericalDerivative(cl); fun->setInterval(0.0000001); optimizer = new ConjugateGradientMultiDimensions(fun); } else if (optMethodDeriv == OPTIMIZATION_NEWTON) { fun = new ThreePointsNumericalDerivative(cl); fun->setInterval(0.0001); optimizer = new PseudoNewtonOptimizer(fun); } else throw Exception("OptimizationTools::optimizeBranchLengthsParameters. Unknown optimization method: " + optMethodDeriv); // Numerical derivatives: ParameterList tmp = parameters.getCommonParametersWith(cl->getParameters()); fun->setParametersToDerivate(tmp.getParameterNames()); optimizer->setVerbose(verbose); optimizer->setProfiler(profiler); optimizer->setMessageHandler(messageHandler); optimizer->setMaximumNumberOfEvaluations(tlEvalMax); optimizer->getStopCondition()->setTolerance(tolerance); // Optimize TreeLikelihood function: optimizer->setConstraintPolicy(AutoParameter::CONSTRAINTS_AUTO); if (listener) optimizer->addOptimizationListener(listener); optimizer->init(parameters); optimizer->optimize(); if (verbose > 0) ApplicationTools::displayMessage("\n"); // We're done. unsigned int n = optimizer->getNumberOfEvaluations(); delete optimizer; // We're done. return n; }
unsigned int OptimizationTools::optimizeNumericalParameters2( DiscreteRatesAcrossSitesTreeLikelihood* tl, const ParameterList& parameters, OptimizationListener* listener, double tolerance, unsigned int tlEvalMax, OutputStream* messageHandler, OutputStream* profiler, bool reparametrization, bool useClock, unsigned int verbose, const std::string& optMethodDeriv) throw (Exception) { DerivableSecondOrder* f = tl; ParameterList pl = parameters; // Shall we use a molecular clock constraint on branch lengths? auto_ptr<GlobalClockTreeLikelihoodFunctionWrapper> fclock; if (useClock) { fclock.reset(new GlobalClockTreeLikelihoodFunctionWrapper(tl)); f = fclock.get(); if (verbose > 0) ApplicationTools::displayResult("Log-likelihood after adding clock", -tl->getLogLikelihood()); // Reset parameters to use new branch lengths. WARNING! 'old' branch parameters do not exist anymore and have been replaced by heights pl = fclock->getParameters().getCommonParametersWith(parameters); pl.addParameters(fclock->getHeightParameters()); } // Shall we reparametrize the function to remove constraints? auto_ptr<DerivableSecondOrder> frep; if (reparametrization) { frep.reset(new ReparametrizationDerivableSecondOrderWrapper(f, pl)); f = frep.get(); // Reset parameters to remove constraints: pl = f->getParameters().subList(pl.getParameterNames()); } auto_ptr<AbstractNumericalDerivative> fnum; // Build optimizer: auto_ptr<Optimizer> optimizer; if (optMethodDeriv == OPTIMIZATION_GRADIENT) { fnum.reset(new TwoPointsNumericalDerivative(f)); fnum->setInterval(0.0000001); optimizer.reset(new ConjugateGradientMultiDimensions(fnum.get())); } else if (optMethodDeriv == OPTIMIZATION_NEWTON) { fnum.reset(new ThreePointsNumericalDerivative(f)); fnum->setInterval(0.0001); optimizer.reset(new PseudoNewtonOptimizer(fnum.get())); } else if (optMethodDeriv == OPTIMIZATION_BFGS) { fnum.reset(new TwoPointsNumericalDerivative(f)); fnum->setInterval(0.0001); optimizer.reset(new BfgsMultiDimensions(fnum.get())); } else throw Exception("OptimizationTools::optimizeNumericalParameters2. Unknown optimization method: " + optMethodDeriv); // Numerical derivatives: ParameterList tmp = tl->getNonDerivableParameters(); if (useClock) tmp.addParameters(fclock->getHeightParameters()); fnum->setParametersToDerivate(tmp.getParameterNames()); optimizer->setVerbose(verbose); optimizer->setProfiler(profiler); optimizer->setMessageHandler(messageHandler); optimizer->setMaximumNumberOfEvaluations(tlEvalMax); optimizer->getStopCondition()->setTolerance(tolerance); // Optimize TreeLikelihood function: optimizer->setConstraintPolicy(AutoParameter::CONSTRAINTS_AUTO); NaNListener* nanListener = new NaNListener(optimizer.get(), tl); optimizer->addOptimizationListener(nanListener); if (listener) optimizer->addOptimizationListener(listener); optimizer->init(pl); optimizer->optimize(); if (verbose > 0) ApplicationTools::displayMessage("\n"); // We're done. return optimizer->getNumberOfEvaluations(); }
unsigned int OptimizationTools::optimizeNumericalParametersWithGlobalClock( DiscreteRatesAcrossSitesClockTreeLikelihood* cl, const ParameterList& parameters, OptimizationListener* listener, unsigned int nstep, double tolerance, unsigned int tlEvalMax, OutputStream* messageHandler, OutputStream* profiler, unsigned int verbose, const std::string& optMethodDeriv) throw (Exception) { AbstractNumericalDerivative* fun = 0; // Build optimizer: MetaOptimizerInfos* desc = new MetaOptimizerInfos(); if (optMethodDeriv == OPTIMIZATION_GRADIENT) { fun = new TwoPointsNumericalDerivative(cl); fun->setInterval(0.0000001); desc->addOptimizer("Branch length parameters", new ConjugateGradientMultiDimensions(fun), cl->getBranchLengthsParameters().getParameterNames(), 2, MetaOptimizerInfos::IT_TYPE_FULL); } else if (optMethodDeriv == OPTIMIZATION_NEWTON) { fun = new ThreePointsNumericalDerivative(cl); fun->setInterval(0.0001); desc->addOptimizer("Branch length parameters", new PseudoNewtonOptimizer(fun), cl->getBranchLengthsParameters().getParameterNames(), 2, MetaOptimizerInfos::IT_TYPE_FULL); } else throw Exception("OptimizationTools::optimizeNumericalParametersWithGlobalClock. Unknown optimization method: " + optMethodDeriv); // Numerical derivatives: ParameterList tmp = parameters.getCommonParametersWith(cl->getBranchLengthsParameters()); fun->setParametersToDerivate(tmp.getParameterNames()); ParameterList plsm = parameters.getCommonParametersWith(cl->getSubstitutionModelParameters()); if (plsm.size() < 10) desc->addOptimizer("Substitution model parameter", new SimpleMultiDimensions(cl), plsm.getParameterNames(), 0, MetaOptimizerInfos::IT_TYPE_STEP); else desc->addOptimizer("Substitution model parameters", new DownhillSimplexMethod(cl), plsm.getParameterNames(), 0, MetaOptimizerInfos::IT_TYPE_FULL); ParameterList plrd = parameters.getCommonParametersWith(cl->getRateDistributionParameters()); if (plrd.size() < 10) desc->addOptimizer("Rate distribution parameter", new SimpleMultiDimensions(cl), plrd.getParameterNames(), 0, MetaOptimizerInfos::IT_TYPE_STEP); else desc->addOptimizer("Rate dsitribution parameters", new DownhillSimplexMethod(cl), plrd.getParameterNames(), 0, MetaOptimizerInfos::IT_TYPE_FULL); MetaOptimizer optimizer(fun, desc, nstep); optimizer.setVerbose(verbose); optimizer.setProfiler(profiler); optimizer.setMessageHandler(messageHandler); optimizer.setMaximumNumberOfEvaluations(tlEvalMax); optimizer.getStopCondition()->setTolerance(tolerance); // Optimize TreeLikelihood function: optimizer.setConstraintPolicy(AutoParameter::CONSTRAINTS_AUTO); if (listener) optimizer.addOptimizationListener(listener); optimizer.init(parameters); optimizer.optimize(); if (verbose > 0) ApplicationTools::displayMessage("\n"); // We're done. return optimizer.getNumberOfEvaluations(); }
unsigned int OptimizationTools::optimizeNumericalParameters( DiscreteRatesAcrossSitesTreeLikelihood* tl, const ParameterList& parameters, OptimizationListener* listener, unsigned int nstep, double tolerance, unsigned int tlEvalMax, OutputStream* messageHandler, OutputStream* profiler, bool reparametrization, unsigned int verbose, const std::string& optMethodDeriv, const std::string& optMethodModel) throw (Exception) { DerivableSecondOrder* f = tl; ParameterList pl = parameters; // Shall we reparametrize the function to remove constraints? auto_ptr<DerivableSecondOrder> frep; if (reparametrization) { frep.reset(new ReparametrizationDerivableSecondOrderWrapper(f, parameters)); f = frep.get(); // Reset parameters to remove constraints: pl = f->getParameters().subList(parameters.getParameterNames()); } // /////////////// // Build optimizer: // Branch lengths MetaOptimizerInfos* desc = new MetaOptimizerInfos(); MetaOptimizer* poptimizer = 0; AbstractNumericalDerivative* fnum = new ThreePointsNumericalDerivative(f); if (optMethodDeriv == OPTIMIZATION_GRADIENT) desc->addOptimizer("Branch length parameters", new ConjugateGradientMultiDimensions(f), tl->getBranchLengthsParameters().getParameterNames(), 2, MetaOptimizerInfos::IT_TYPE_FULL); else if (optMethodDeriv == OPTIMIZATION_NEWTON) desc->addOptimizer("Branch length parameters", new PseudoNewtonOptimizer(f), tl->getBranchLengthsParameters().getParameterNames(), 2, MetaOptimizerInfos::IT_TYPE_FULL); else if (optMethodDeriv == OPTIMIZATION_BFGS) desc->addOptimizer("Branch length parameters", new BfgsMultiDimensions(f), tl->getBranchLengthsParameters().getParameterNames(), 2, MetaOptimizerInfos::IT_TYPE_FULL); else throw Exception("OptimizationTools::optimizeNumericalParameters. Unknown optimization method: " + optMethodDeriv); // Other parameters if (optMethodModel == OPTIMIZATION_BRENT) { ParameterList plsm = parameters.getCommonParametersWith(tl->getSubstitutionModelParameters()); desc->addOptimizer("Substitution model parameter", new SimpleMultiDimensions(f), plsm.getParameterNames(), 0, MetaOptimizerInfos::IT_TYPE_STEP); ParameterList plrd = parameters.getCommonParametersWith(tl->getRateDistributionParameters()); desc->addOptimizer("Rate distribution parameter", new SimpleMultiDimensions(f), plrd.getParameterNames(), 0, MetaOptimizerInfos::IT_TYPE_STEP); poptimizer = new MetaOptimizer(f, desc, nstep); } else if (optMethodModel == OPTIMIZATION_BFGS) { vector<string> vNameDer; ParameterList plsm = parameters.getCommonParametersWith(tl->getSubstitutionModelParameters()); vNameDer = plsm.getParameterNames(); ParameterList plrd = parameters.getCommonParametersWith(tl->getRateDistributionParameters()); vector<string> vNameDer2 = plrd.getParameterNames(); vNameDer.insert(vNameDer.begin(), vNameDer2.begin(), vNameDer2.end()); fnum->setParametersToDerivate(vNameDer); desc->addOptimizer("Rate & model distribution parameters", new BfgsMultiDimensions(fnum), vNameDer, 1, MetaOptimizerInfos::IT_TYPE_FULL); poptimizer = new MetaOptimizer(fnum, desc, nstep); } else throw Exception("OptimizationTools::optimizeNumericalParameters. Unknown optimization method: " + optMethodModel); poptimizer->setVerbose(verbose); poptimizer->setProfiler(profiler); poptimizer->setMessageHandler(messageHandler); poptimizer->setMaximumNumberOfEvaluations(tlEvalMax); poptimizer->getStopCondition()->setTolerance(tolerance); // Optimize TreeLikelihood function: poptimizer->setConstraintPolicy(AutoParameter::CONSTRAINTS_AUTO); NaNListener* nanListener = new NaNListener(poptimizer, tl); poptimizer->addOptimizationListener(nanListener); if (listener) poptimizer->addOptimizationListener(listener); poptimizer->init(pl); poptimizer->optimize(); if (verbose > 0) ApplicationTools::displayMessage("\n"); // We're done. unsigned int nb = poptimizer->getNumberOfEvaluations(); delete poptimizer; return nb; }