How to fit data on a specific region, using lsqcurvefit or any other?

Question

Alfredo Scigliani 2023-3-21

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1932955-how-to-fit-data-on-a-specific-region-using-lsqcurvefit-or-any-other

评论： Mathieu NOE 2023-4-6

Untitled2.m

I am trying to use lsqcurvefit to fit 4 parameters, I am usually able to do it but in this problem, the fitting function is only valid close to the plateau values. So I selected the values over which I apply the fitting. for example (14:37) in line 12.

The function seem to be a bit sensitive to initial guesses and that is not allowing me to achieve a good fit.

I wonder if the problem is that I am using lsqcurvefit and I should be using a different tool or if I need to create an iteration of initial guesses so that I can finally get the correct combination?

It is acceptable to increase or decrease a bit the range of the fitting function helps to get a better fit. Right now it is set for 14:37 for data set 1 and 12:37 for data set 2.

function:

code:

clear; clc; clf; close all;

x1 = [0.0005 0.0007 0.0008 0.001 0.0013895 0.0019307 0.0026827 0.0037276 0.0051795 0.0071969 0.01 0.013895 0.019307 0.026827 0.037276 0.051795 0.071969 0.1 0.13895 0.19307 0.26827 0.37276 0.51795 0.71969 1 1.3895 1.9307 2.6827 3.7276 5.1795 7.1969 10 13.895 19.307 26.827 37.276 51.795 71.969 100 138.95 193.07 268.27 372.76 517.95 719.69 1000];

y1 = [0.011171 0.015461 0.017605 0.021901 0.029975 0.041604 0.057316 0.078645 0.10726 0.14398 0.18833 0.24482 0.31533 0.38492 0.45433 0.51514 0.53269 0.54151 0.56182 0.56775 0.56635 0.56416 0.55786 0.55369 0.53971 0.52712 0.51333 0.50512 0.49684 0.48869 0.50622 0.51273 0.53057 0.55463 0.61141 0.78976 0.96969 1.1629 1.4075 1.6804 2.1478 2.6287 3.5048 4.4863 5.6525 7.3741];

x1_fit = x1(12:37);

y1_fit = y1(12:37);

x2 = [0.0005 0.0007 0.0008 0.001 0.0013895 0.0019307 0.0026827 0.0037276 0.0051795 0.0071969 0.01 0.013895 0.019307 0.026827 0.037276 0.051795 0.071969 0.1 0.13895 0.19307 0.26827 0.37276 0.51795 0.71969 1 1.3895 1.9307 2.6827 3.7276 5.1795 7.1969 10 13.895 19.307 26.827 37.276 51.795 71.969 100 138.95 193.07 268.27 372.76 517.95 719.69];

y2 = [0.0098762 0.013754 0.015562 0.019716 0.027332 0.038068 0.052772 0.073219 0.10152 0.13964 0.19046 0.25666 0.43342 0.53389 0.62881 0.7047 0.74682 0.73801 0.68675 0.64029 0.60505 0.59568 0.59975 0.61275 0.6274 0.64335 0.6483 0.64144 0.6314 0.64113 0.66358 0.70689 0.74901 0.79417 0.86603 0.93868 1.1244 1.3195 2.0893 2.6337 3.0766 3.6119 4.7272 6.0726 7.5609];

x2_fit = x2(14:37);

y2_fit = y2(14:37);

scale1 = 38.50;

scale2 = 48.27;

C1_guess = [1 1 0.001 500];

C2_guess = [1 1 0.001 500];

figure

subplot(1,2,1)

fun = @(C,x) C(1)*sinh(C(2)+(log10(x) - C(3)))+ C(4);

C_1 = lsqcurvefit(fun, C1_guess, x1_fit, y1_fit);

Local minimum possible. lsqcurvefit stopped because the final change in the sum of squares relative to its initial value is less than the value of the function tolerance.

a = loglog(x1*scale1, y1,'ks', 'MarkerFaceColor', 'c');

hold on

a1 = loglog(x1_fit*scale1, fun(C_1,x1_fit), 'r-', 'LineWidth', 1);

xlabel('x1')

ylabel('y1')

subplot(1,2,2)

C_2 = lsqcurvefit(fun, C2_guess, x2_fit, y2_fit);

Local minimum possible. lsqcurvefit stopped because the final change in the sum of squares relative to its initial value is less than the value of the function tolerance.

b = loglog(x2*scale2, y2,'ks', 'MarkerFaceColor', 'g');

hold on

b1 = loglog(x2_fit*scale2, fun(C_2,x2_fit), 'r-', 'LineWidth', 1);

xlabel('x2')

ylabel('y2')

4 个评论
显示 2更早的评论隐藏 2更早的评论

Alfredo Scigliani 2023-3-21

编辑：Alfredo Scigliani 2023-3-21

And if you extract the second output of lsqcurvefit (which is the rsb) it is about 0.20

Matt J 2023-3-22

The fitting function is supposed to pass through the data points almkst perfectly, lets say a deviation about 0.001.

That cannot happen with the sinh model that you've shown. The sinh function is unimodal, whereas the plateau in your data has oscillations considerably bigger than 0.001.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Mathieu NOE 2023-3-22

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1932955-how-to-fit-data-on-a-specific-region-using-lsqcurvefit-or-any-other#answer_1198300

在 MATLAB Online 中打开

hello

fyi, some results obtained with the poor's man solution (based on fminsearch)

NB : I changed your log10 to log in my code

for the first set of data , this is what I could obtain

sol = 0.0019 2.0232 0.0002 0.5240

plot in log log scale

x1 = [0.0005 0.0007 0.0008 0.001 0.0013895 0.0019307 0.0026827 0.0037276 0.0051795 0.0071969 0.01 0.013895 0.019307 0.026827 0.037276 0.051795 0.071969 0.1 0.13895 0.19307 0.26827 0.37276 0.51795 0.71969 1 1.3895 1.9307 2.6827 3.7276 5.1795 7.1969 10 13.895 19.307 26.827 37.276 51.795 71.969 100 138.95 193.07 268.27 372.76 517.95 719.69 1000];
y1 = [0.011171 0.015461 0.017605 0.021901 0.029975 0.041604 0.057316 0.078645 0.10726 0.14398 0.18833 0.24482 0.31533 0.38492 0.45433 0.51514 0.53269 0.54151 0.56182 0.56775 0.56635 0.56416 0.55786 0.55369 0.53971 0.52712 0.51333 0.50512 0.49684 0.48869 0.50622 0.51273 0.53057 0.55463 0.61141 0.78976 0.96969 1.1629 1.4075 1.6804 2.1478 2.6287 3.5048 4.4863 5.6525 7.3741];
x = x1(12:46);
y = y1(12:46);
% curve fit using fminsearch
% We would like to fit the function
f = @(a,b,c,d,x) a*sinh(b+(log(x) - c))+ d;
obj_fun = @(params) norm(f(params(1), params(2), params(3), params(4),x)-y);
C1_guess = [0.001 1 0.001 1];
sol = fminsearch(obj_fun, C1_guess); % 
a_sol = sol(1);
b_sol = sol(2);
c_sol = abs(sol(3));
d_sol = sol(4);
y_fit = f(a_sol, b_sol,c_sol, d_sol, x);
Rsquared = my_Rsquared_coeff(y,y_fit); % correlation coefficient
loglog(x, y_fit, '-',x,y, 'r .-', 'MarkerSize', 25)
title(['Fit equation to data : R² =  ' num2str(Rsquared) ], 'FontSize', 20)
xlabel('x data', 'FontSize', 20)
ylabel('y data', 'FontSize', 20)
toc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function Rsquared = my_Rsquared_coeff(data,data_fit)
    % R2 correlation coefficient computation
    
    % The total sum of squares
    sum_of_squares = sum((data-mean(data)).^2);
    
    % The sum of squares of residuals, also called the residual sum of squares:
    sum_of_squares_of_residuals = sum((data-data_fit).^2);
    
    % definition of the coefficient of correlation is
    Rsquared = 1 - sum_of_squares_of_residuals/sum_of_squares;
    
    
end

2 个评论
显示无隐藏无

Mathieu NOE 2023-3-22

在 MATLAB Online 中打开

And for the second data set

sol = 0.0142 0.3511 0.0112 0.6948

x2 = [0.0005 0.0007 0.0008 0.001 0.0013895 0.0019307 0.0026827 0.0037276 0.0051795 0.0071969 0.01 0.013895 0.019307 0.026827 0.037276 0.051795 0.071969 0.1 0.13895 0.19307 0.26827 0.37276 0.51795 0.71969 1 1.3895 1.9307 2.6827 3.7276 5.1795 7.1969 10 13.895 19.307 26.827 37.276 51.795 71.969 100 138.95 193.07 268.27 372.76 517.95 719.69];
y2 = [0.0098762 0.013754 0.015562 0.019716 0.027332 0.038068 0.052772 0.073219 0.10152 0.13964 0.19046 0.25666 0.43342 0.53389 0.62881 0.7047 0.74682 0.73801 0.68675 0.64029 0.60505 0.59568 0.59975 0.61275 0.6274 0.64335 0.6483 0.64144 0.6314 0.64113 0.66358 0.70689 0.74901 0.79417 0.86603 0.93868 1.1244 1.3195 2.0893 2.6337 3.0766 3.6119 4.7272 6.0726 7.5609];
x = x2(12:45);
y = y2(12:45);
% curve fit using fminsearch
% We would like to fit the function
f = @(a,b,c,d,x) a*sinh(b+(log(x) - c))+ d;
obj_fun = @(params) norm(f(params(1), params(2), params(3), params(4),x)-y);
C1_guess = [0.001 1 0.001 1];
sol = fminsearch(obj_fun, C1_guess) % 
a_sol = sol(1);
b_sol = sol(2);
c_sol = abs(sol(3));
d_sol = sol(4);
y_fit = f(a_sol, b_sol,c_sol, d_sol, x);
Rsquared = my_Rsquared_coeff(y,y_fit); % correlation coefficient
loglog(x, y_fit, '-',x,y, 'r .-', 'MarkerSize', 25)
title(['Fit equation to data : R² =  ' num2str(Rsquared) ], 'FontSize', 20)
xlabel('x data', 'FontSize', 20)
ylabel('y data', 'FontSize', 20)
toc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function Rsquared = my_Rsquared_coeff(data,data_fit)
    % R2 correlation coefficient computation
    
    % The total sum of squares
    sum_of_squares = sum((data-mean(data)).^2);
    
    % The sum of squares of residuals, also called the residual sum of squares:
    sum_of_squares_of_residuals = sum((data-data_fit).^2);
    
    % definition of the coefficient of correlation is
    Rsquared = 1 - sum_of_squares_of_residuals/sum_of_squares;
    
    
end