Speed up HC calculations 2 to 5x using Newton's method

[joethorley]

Email from David Fox

I’ve been thinking about the ssd_hc() function and the use of uniroot(). While uniroot() is no doubt optimized for performance, it nevertheless uses numerical methods to approximate the gradient function which can be time consuming. It occurred to me that I could potentially speed up the HC calculations because an explicit formula for the gradient function is readily available.

I’ve written up my work in the attached Word document and coded the method both as R-code and C++ code for the iterative bit.

I used a simple example by fitting a model-averaged SSD comprised of the lognorm and gamma distributions to the ccme_boron data. As you can see from the timings provided in the attached my R-code gave an HC5 in 0.02 seconds while ssd_hc() took twice that long. Performance improvements were achieved with the C++ coding which resulted in my code computing an HC5 in 0.01 seconds and ssd_hc() taking 0.05 seconds. I would expect the differences to be even greater when fitting the full bcanz distribution set and running simulations involving 1,000s of HC calculations.

HC_Newton_approx.pdf