MurrellGroup / MolecularEvolution.jl

A Julia framework for developing phylogenetic models
MIT License
9 stars 4 forks source link

Explore optimization algorithms #14

Open murrellb opened 5 months ago

murrellb commented 5 months ago

... for rates, and for discrete distributions (ie. with a sum to 1 constraint). Start here: https://github.com/libprima/PRIMA.jl

nossleinad commented 1 month ago

Switching from NLopt to PRIMA could look like this for a typical BOBYQA use case

NLopt

opt = Opt(:LN_BOBYQA, num_params)

min_objective!(opt, (x,y) -> (objective ∘ unflatten)(x))
lower_bounds!(opt, [-5.0 for i in 1:num_params])
upper_bounds!(opt, [5.0 for i in 1:num_params])
xtol_rel!(opt, 1e-12)
score,mini,did_it_work = NLopt.optimize(opt, flat_initial_params)

final_params = unflatten(mini)
optimized_model = build_model_vec(final_params)
LL = log_likelihood!(tree,optimized_model)

println(did_it_work)
println("Opt LL:",LL)

PRIMA

lower_bounds = [-5.0 for i in 1:num_params]
upper_bounds = [5.0 for i in 1:num_params]
tr_radius = 1e-12 #trust-region radius

mini, info = bobyqa(objective ∘ unflatten, flat_initial_params, xl=lower_bounds, xu=upper_bounds, rhoend=tr_radius)

final_params = unflatten(mini)
optimized_model = build_model_vec(final_params)
LL = log_likelihood!(tree,optimized_model)

println("SUCCESS: $(issuccess(info))")
println("Opt LL:",LL)