# Fit - Nonlinear Regression - can this run faster? Python >20x faster?

17 views (last 30 days)
arnold on 19 Sep 2015
Answered: Ilya on 20 Sep 2015
Hi,
I'm looking for a way to speed up fitting an exponential model to a stack of eleven double images. I've been at it for days but can't seem to find a way to speed things up. Even using parfor this takes 20 minutes on a downsampled image stack (100x100), which is WAY slowwer than single thread of python (6 minutes per 1024x1024 pixel stack!). I was wondering if somebody could help me out speeding things up a bit - I guess I must clearly be missing something...
I load all 11 images into an array and permute it in order to be able to access the fitting dimension right.
imgdata2 = permute(imgdata,[3,1,2]);
size(imgdata2)
ans =
11 1024 1024
The model I want to fit is:
fun = @(j, x) j(1).*exp(x./0.0257)+j(2).*exp(x./(2*0.0257))- j(3);
pguess = [100e-15, 10e-9, 0.04];
Y is given and the same for all fits. A simple fit now looks like this:
x = imgdata2(:,i,k)./1000;
out = fitnlm(x,y,fun, pguess);
Now the fastest way I could figure out to fit all 1024*1024 pixels of the image stacks was using parfor (yes, I pre-start the workers). For this I had to create a dummy cell array and then unwrap the results afterwards in another loop.
%%preallocate
fitresults_wrapper = cell(1,(size(imgdata2,2)*size(imgdata2,3)));
%%fitting
parfor i = 1:(size(imgdata2,2)*size(imgdata2,3)); % for all 1024 by 1024 'pixel-stacks'
x = imgdata2(:,i)./1000;
if sum(isnan(x))== 0
fitresults_wrapper{i} = nlinfit(x,y,fun, pguess);
end;
end;
Now. Is there no way to speed things up? If this always takes a couple of minutes it's pretty much useles. A simple python loop can do the full images in 6 minutes using just one thread (curve_fit(func, x, y, [50e-12,1e-6, 0.04])) I was unable to find a way to vectorize the whole thing, which is usually the first thing I do in matlab in order to speed things up. Is there a way to vectorize fitting here?
Any help or hints appreciated Arnold

Ilya on 20 Sep 2015
I am no expert on image processing, but based on your code samples, you don't need to use fitnlm. When you do this
x = imgdata2(:,i,k)./1000;
or this
x = imgdata2(:,i)./1000;
you get a vector. Just form a matrix with two columns by putting
X = [exp(x./0.0257) exp(x./(2*0.0257))];
and pass this matrix to a linear utility such as say fitlm to fit for the two linear coefficients, j(1) and j(2), and intercept, -j(3). This is going to be much faster.