Lucile
Lucile

Reputation: 45

Fitting a power law: "NaN computed by model function"

I try to fit a power law function with matlab (y=ax^b)

Here are my x and y matrices

I simply compute the fitting with

fit(x,y,'power1')

I get this error:

Error using fit>iFit (line 415) NaN computed by model function, fitting cannot continue. Try using or tightening upper and lower bounds on coefficients. Error in fit (line 109) [fitobj, goodness, output, convmsg] = iFit( xdatain, ydatain, fittypeobj, ...

WHYYYY!? There is no 0 in my x and y matrices, nothing that would return a NaN value I think, and I can compute without any problem the inverse relationship fit(y,x,'power1').

Thanks for any help/suggestion!!

EDIT: (just as a precision) Excel does find a power law fit for (x,y)!

EDIT2 : code, once x and y are stored as variables.:

[p_powerlaw,results_powerlaw] = fit(x,y,'power1');

EDIT3: I have changed the link. Now there, in my dropbox, you will find the .mat of x and y...try to fit them with the power1... :P doesn't work! why? I don't get it...

And try to copy paste the values of the 2 matrices x and y into other matrices (not by assigning, really by copy-pasting the values)...no problem with the fitting..!

Upvotes: 2

Views: 3031

Answers (1)

am304
am304

Reputation: 13876

There is no equivalent to fit in Octave, however, I can do some nonlinear curve fitting (the MATLAB equivalent would be either lsqcurvefit or lsqnonlin the optimization toolbox):

% Define function handle for optimisation (equivalent to power1 fit in MATLAB)
f = @(p, x) p(1) * x.^p(2);

% Initial guess for power law coefficients
init = [1 1]';

% Run optimisation - data is contained in X matrix
[p,y_power,cvg,op] = nonlin_curvefit(f,init,X(:,1),X(:,2));

% Sort out the fitted data to be monotonically increasing
Y = zeros(size(X));
[Y(:,1),idx] = sort(X(:,1),1);
Y(:,2) = y_power(idx);

% Plot the data
plot(X(:,1),X(:,2),'x',Y(:,1),Y(:,2))

% Display the power law coefficients
clc
disp('y = a * x^b');
disp(['a = ' num2str(p(1))])
disp(['b = ' num2str(p(2))])

This gives me the following coefficients (the same as @Ander Biguri) and graph:

y = a * x^b
a = 13.416
b = 0.94642

enter image description here

Upvotes: 0

Related Questions