Optimizing curving fitting for speed

Share icon

I am attempting to speed up a user defined fit function. The function that I am trying to fit to looks something like:

f(x) = Integral( g(x,v)p(v) dv)

The fitting function was originally written using the point-by-point method and a for loop to numerically integrate over v. I sped up the fitting by rewriting the code to use the all-at-once routine, but I still have the slow and imprecise for loop. I recently stumbled onto the Integrate1D function and thought that it might provide further speed/precision enhancements, but I'd like to know if anybody has any experience using it with curve fitting. I did see from a previous post that it looks like I'd have to use global variables to get the fit parameters into the integrate function.

Any other ideas would also be appreciated. Thanks.
I think you've understood the problem very well. Integrate1D might help in the speed department; you may need to fiddle a bit to get the right combination of options.

You are correct that globals are required to pass fit coefficient to Integrate1D.

The fastest possible solution would be to write the function in C using the XOP Toolkit.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

July 7, 2008 at 12:55 pm - Permalink

Thanks for the reply. I am now trying to implement this and I'm actually finding that Integrate1D takes longer than a simple for loop! This doesn't seem right. My for loop for integration is this:

Wave output=0 //This is where the best fit function lives
For( v=1000; v<4000; v+=50 )
pv = (1/(2 * pi)^.5/sv) * exp( - (v-v0)^2 / (2 * sv^2) ) //a probability from a gaussian
x1 =   .41675 *1000/v //diffraction order spacing for a given atom velocity

output += pv* exp(-(pos-x0)^2 / (2* sb^2))   //0th order gaussian
output += pv* p1 * ( exp(-(pos-x0-1*x1)^2 / (2* sb^2)) +  exp(-(pos-x0+1*x1)^2 / (2* sb^2))  )  //+ and - 1st orders
output += // repeat up to 6th orders
EndFor


There were (4000-1000)/50=60 steps in the for loop. The help for integrate1D is slightly unclear as to what, exactly, the counts parameters is but it sounds like I should set it to 60 to achieve the same precision as the for loop. So, the integrate1D version is this:

output = Integrate1D(IntegrateDiffractionFunc,1000,4000,0,60)

Function IntegrateDiffractionFunc(v)
Variable v
Variable amp //the best fit diffraction amplitude at a single point
NVAR p1, p2, p3, p4, p5, p6, p7, x0, sb, sv, v0 //weighting factors and other fit parameters
Wave pos = testWave //independent variable wave
   
Variable pv = (1/(2 * pi)^.5/sv) * exp( - (v-v0)^2 / (2 * sv^2) )
Variable x1 =   .41675 *1000/v

amp += pv* exp(-(pos-x0)^2 / (2* sb^2))   //0th order
amp += pv* p1 * ( exp(-(pos-x0-1*x1)^2 / (2* sb^2)) +  exp(-(pos-x0+1*x1)^2 / (2* sb^2))  ) //1st orders
amp += // same as the above up to 6 orders


Am I doing something wrong with the Integrate1D function? My guess is that in the for loop we can compute the integral for all points in the output wave basically simultaneously, but with the Integrate1D method the integral must be computed separately for each point. Thanks.

July 14, 2008 at 05:16 pm - Permalink

Maybe you should send me an example to support@wavemetrics.com. Send an experiment file that includes both your original procedure and the new one that uses Integrate1D, plus data and instructions or a command line to make it run. I'll see what I can figure out.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

July 16, 2008 at 11:11 am - Permalink