How to fit the data with 2 term power law function?

Hi there,

I have the followig data and I want to fit them in the form of z = r^a * v^b.

I tried to obtain log values of the function: log(z) =a*log(r)+b*log(v), but I'm confused as there are so many vectors.

r = [4.242640687119285e-06,5.656854249492381e-06,7.071067811865476e-06,8.485281374238570e-06,9.899494936611665e-06,1.131370849898476e-05,1.272792206135786e-05];
v = [0.3656,0.6093,0.9748,1.2185];
z = [2.20075780727420	2.34800174967587	2.40633263085150	2.46928643488236	2.55551444168423	2.66309682650254	2.78591291730793
1.20828256541234	1.57243262679317	1.78594063849408	1.94976208027465	2.10251979263655	2.25584172701771	2.41163339237766
0.472116287092821	0.855613921769617	1.14739856522263	1.38120431719641	1.58720801177653	1.78025160940097	1.96628391586011
0.243938638287455	0.563375613060710	0.851155979561504	1.09774989977420	1.31858268912292	1.52483525677172	1.72202346894929
];

the plot is also attached, where different curves are from different v values.

[peipei feng]

Answer 1

[Ameer Hamza]

taking log() of the equation make it very easy to find the parameter (because the equation is linear in term of a and b). Try the following code

r = [4.242640687119285e-06,5.656854249492381e-06,7.071067811865476e-06,8.485281374238570e-06,9.899494936611665e-06,1.131370849898476e-05,1.272792206135786e-05];
v = [0.3656,0.6093,0.9748,1.2185];
z = [2.20075780727420	2.34800174967587	2.40633263085150	2.46928643488236	2.55551444168423	2.66309682650254	2.78591291730793
1.20828256541234	1.57243262679317	1.78594063849408	1.94976208027465	2.10251979263655	2.25584172701771	2.41163339237766
0.472116287092821	0.855613921769617	1.14739856522263	1.38120431719641	1.58720801177653	1.78025160940097	1.96628391586011
0.243938638287455	0.563375613060710	0.851155979561504	1.09774989977420	1.31858268912292	1.52483525677172	1.72202346894929
];

a = zeros(1, 4);
b = zeros(1, 4);

for i=1:numel(a)
    V = repmat(v(i), size(r));
    Z = log(z(i, :).');
    X = [log(r.') log(V.')];
    sol = X\Z;
    a(i) = sol(1);
    b(i) = sol(2);
end

figure;
axes();
hold on
for i=1:numel(v)
    z_est = r.^a(i).*v(i).^b(i);
    plot(r, z(i,:), '+--', 'DisplayName', ['v = ' num2str(v(i)) ' real']);
    plot(r, z_est, 'DisplayName', ['v = ' num2str(v(i)) ' est']);
end
legend('Location', 'best');

[peipei feng]

Thank you for your answer. However, I found the results badly fitted. It seems that the impact of r have been negelected, resulting in flaten fitted curves.

[Ameer Hamza]

I misunderstood the question initially. See the updated answer.

[peipei feng]

Thank you very much for the updates. Actually you didn't misunderstood the question. My aim is to obtain one equation with the form z = r^a * v^b (with one a and one b) that gives me four different curves at four different speeds (v). In your updated codes, I found the values for a and b vary alot for different speeds(v), does this mean that my initial aim is unachieveable?

[Ameer Hamza]

In that case, the original solution was correct. As you can see, the value of a and b vary very much to fit these lines, so it might not be possible to get a good fit with single 'a' and 'b'.

Here is the original code for reference

r = [4.242640687119285e-06,5.656854249492381e-06,7.071067811865476e-06,8.485281374238570e-06,9.899494936611665e-06,1.131370849898476e-05,1.272792206135786e-05];
v = [0.3656,0.6093,0.9748,1.2185];
z = [2.20075780727420	2.34800174967587	2.40633263085150	2.46928643488236	2.55551444168423	2.66309682650254	2.78591291730793
1.20828256541234	1.57243262679317	1.78594063849408	1.94976208027465	2.10251979263655	2.25584172701771	2.41163339237766
0.472116287092821	0.855613921769617	1.14739856522263	1.38120431719641	1.58720801177653	1.78025160940097	1.96628391586011
0.243938638287455	0.563375613060710	0.851155979561504	1.09774989977420	1.31858268912292	1.52483525677172	1.72202346894929
];

[V, R] = ndgrid(v, r);

Y = log(z(:));
X = [log(R(:)) log(V(:))];
sol = X\Y; % least square solution

a = sol(1);
b = sol(2);

[peipei feng]

Thank you!