x = [1 4 9 15 25 35 45 55 65 75 85 95 150 250 350 450 900];
y = [0 335 28 37 9 4 3.5 2 1 0.7 0.6 0.5 0.3 0.1 0.01 0.005 0.001];
Look CAREFULLY at your data!
What is y(1)? Do you see that it is tremendously out of character, inconsistent with the rest of your data?
y(1)
Of course, when you plot it on a log scale, MATLAB does not plot that point, because log(0) is undefined. (-inf if you insist.) HOWEVER, when you do a fit, fit sees it!
And what that does is completely drag your curve around. y(1) = 1000 or larger might be more consistent. But not zero.
ft = fittype('a*(x^n)','dependent','y','independent','x');
Next, starting ANY optimization at [0,0] is a HIGHLY, HUGELY dangerous thing to do.
[fitresult, ~] = fit(x(2:end)',y(2:end)',ft,'Startpoint', [1 -1])
p1 = plot(x,y,'o-','LineWidth',1.5,'color','b','DisplayName','inputs');
hold on
p2 = plot(x,(fitresult.a .* (x.^fitresult.n)),'-','LineWidth',1.5,'color','k','DisplayName','best fit');
set(gca, 'XScale', 'log', 'YScale', 'log');
legend([p1 p2])
hold off
Is the fit very good, when plotted on a log scale? Well, not really. But it never is a good idea to do this sort of thing. Why not? Because MATLAB will not worry about TINY deviations in those lower points, even though in proportional terms, it will be very obvious on a log scale.
What you should see is the fit fits the first two points quite well, and then misses the rest. But again, think about what a log scale does!!!!!
fit assumes an additive error model, with the errors presumed to be homogeneous, and additive. That is a bad thing when your dats spans 4 or 5 powers of 10. It would be far smarter to log your data and your model, THEN do the fit. Now instead of proportional error, you will have additive error. See how nice it gets?
When you log your model, the problem becomes trivial.
Your model: y = a*x^n * (proportional noise)
Logged model: log(y) = log(a) + n*log(x) + (additive homogeneous noise)
We still need to drop the bogus first data point of course. But now your model is just a linear 'poly1' model.
mdl = fit(log(x(2:end)'),log(y(2:end)'),'poly1')
Don't forget to exponentiate the coefficient out front!
a = exp(mdl.p2)
n = mdl.p1
Or, we could use fittype
ft = fittype('log(a) + n*log(x)','dependent','y','independent','x');
Note that I logged x inside fittype, but not y.
[fitresult, ~] = fit(x(2:end)',log(y(2:end)'),ft,'Startpoint', [1 -1])
And you see this generates the same result for the parameters.
