daijiang / phyr

Functions for phylogenetic analyses
https://daijiang.github.io/phyr/
GNU General Public License v3.0
30 stars 10 forks source link

Optimization in C++ using nlopt #15

Closed lucasnell closed 6 years ago

lucasnell commented 6 years ago

I setup an example optimization using nlopt. It's OLS for a univariate regression, but could be easily applied to the likelihood functions defined in here. You just adjust the likelihood function, the object type you input, and the "initial guesses" array.

For now, the only algorithms I have working are BOBYQA, COBYLA, and PRAXIS. See here for more info on these. Nelder-Mead and some others blow up the R session when used, and I believe it's an nlopt bug.

This version works using the developmental version of nloptr on github, although I could only get that version installed when nlopt was installed separately on my machine.

daijiang commented 6 years ago

Thanks @lucasnell ! But I am not sure how easy to use the function you wrote within other functions such as cor_phylo. Have you tried to use it within cor_phylo? If you can make it work, I'd be happy to merge this pull request. If this will take you too much time, then don't worry about it. I think calling optim into c++ should be fine. Thanks again!

lucasnell commented 6 years ago

It should be very simple to implement, but perhaps I'm too optimistic. And just to be clear, I see this function being especially useful for the parametric bootstrapping of cor_phylo, where going from C++ to R 2,000 times should be pretty costly.

I guess I'll add cor_phylo in a separate pull request first, then add the nlopt version?

daijiang commented 6 years ago

You can keep changing files by adding cor_phylo and these changes will be in this same pull request.

For the bootstrapping, it makes more sense to write the whole loop in c++ (not just the function within the loop), so that we only going from R to c++ once, and c++ to R at the end.

If you can manage to use nlopt within cor_phylo that would be awesome!

lucasnell commented 6 years ago

Okay, I'll just do that.

And for the bootstrapping, I totally agree. I just meant that performance-wise using nlopt would probably make the biggest difference for the bootstrapping versus the original model fit because going from C++ to R would happen more times during bootstrapping.

I'll try to get this done by Friday afternoon if possible.

lucasnell commented 6 years ago

The current version has cor_phylo using nlopt optimization in C++. I'm still figuring out why Nelder-Mead isn't working, but BOBYQA works well.

daijiang commented 6 years ago

Thanks @lucasnell !! You are awesome!

Can you write some tests for cor_phylo? See the tests folder for examples. The main purpose is to make sure that things work as expected. For example, does it give the same result as the ape version? For this purpose, it may be a good idea to keep the base::optim() option.

After you are done with this function. I wonder can you help me with moving the optimization in other functions into nlopt. At this moment, I just used the nloptr::nloptr() function. But ideally, we should do what you have done with the cor_phylo function. No problem at all if you are busy. The current code is working relatively fast. The small potential gain in speed may not worth the effort to use nlopt.

lucasnell commented 6 years ago

I believe cor_phylo should be mostly ready to go now, except for the lack of parametric bootstrapping.