How to use Nlinfit for a function with two independent variables?

Question

Faizan Lali 2023-2-23

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1917445-how-to-use-nlinfit-for-a-function-with-two-independent-variables

评论： Walter Roberson 2023-2-27

Hi here is my data and code and I am trying to predict parameters for a function with two independent variables but Nlinfit is giving me error.

clc

clear %all; % Clear the workspace.

close all; % Close all figures.

format compact

%% Read in data

data =readmatrix('Flexible_BUC.xlsx');

x1=[7.50000000000000

7.50000000000000

7.70000000000000

5

7.50000000000000

6

3.75000000000000

8

7.50000000000000

8

8];

x2=[0.00153000000000000

0.00522000000000000

0.000189000000000000

3.73000000000000e-06

1.17000000000000e-05

8.66000000000000e-05

1.17000000000000e-05

1.51000000000000e-05

2.99000000000000e-05

7.00000000000000e-05

7.92000000000000e-05

6.06000000000000e-05

0.000163000000000000

0.000656000000000000

0.000818000000000000

0.00129000000000000

9.78000000000000e-06

0.000161000000000000

0.000183000000000000

0.000204000000000000

0.000297000000000000

0.000343000000000000

0.00705000000000000

0.00850000000000000

0.0233000000000000

0.0250000000000000

0.0267000000000000

0.00125000000000000

0.00245000000000000

0.00391000000000000

0.00878000000000000

0.00111000000000000

0.00122000000000000

3.23000000000000e-05

4.26000000000000e-05

4.53000000000000e-05

4.13000000000000e-05

8.24000000000000e-05

8.33000000000000e-05];

yobs=[5.95833333300000

0.300000000000000

0.0625000000000000

0.111111111000000

0.213809289000000

0.140625000000000

0.651315789000000

0.351694915000000

12.0555555600000

0.626846311000000

0.555555556000000

0.136363636000000

6.38793103400000

0.233051458000000

0.540000000000000

0.0277777780000000

1.05555555600000

0.113636364000000

0.0933323590000000

0.352272727000000

3.20833333300000

0.897435897000000

1.17046404700000

1.41666666700000

1.79545454500000

1.15384615400000

1.85576923100000

10.8333333300000

0.848684211000000

6.18835443000000

0.767441860000000

0.527777778000000

1.54872306000000

0.337691494000000

1.08333333300000

1.87477002000000

1.81654734900000

1.97222222200000

0.550000000000000

0.340020401000000];

xm=[x1 x2];

%% Initial parameter guesses

C1=0.25;

C2=0.73;

beta0(1)=C1; %initial guess beta 1

beta0(2)=C2; %initial guess beta 2

p=length(beta0); %p = # parameters

%% define function to be used for inverse problem

fINV=@Project_funcINV;

%fnameINV=@forderexpINV;

[beta,resids,J,COVB,mse] = nlinfit(xm,yobs,fINV,beta0);

beta

%% Functions

function y=Project_funcINV(beta0,x1,x2)

c2s=@(x)-2.40874-39.748*(1+x).^-2.856;

y=100./(1+exp(-beta0(1).*c2s(x1)+(beta0(2).*c2s(x1).*log10(100.*x2))));

end

Someone can please help, I would appreciate it.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Star Strider 2023-2-23

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1917445-how-to-use-nlinfit-for-a-function-with-two-independent-variables#answer_1177990

在 MATLAB Online 中打开

You are using the correct approach with:

xm=[x1 x2];

In the function, refer to ‘x1’ as ‘xm(:,1)’ and ‘x2’ as ‘xm(:,2)’ , passing ‘xm’ as the independent variable to ‘Project_funcINV’. I made those changes, and added a fitnlm call to display the statistics, and provided a plot of the data and the fit to it (as a line plot).

Try this —

x1=[7.50000000000000

7.50000000000000

7.70000000000000

5

7.50000000000000

6

3.75000000000000

8

7.50000000000000

8

8];

x2=[0.00153000000000000

0.00522000000000000

0.000189000000000000

3.73000000000000e-06

1.17000000000000e-05

8.66000000000000e-05

1.17000000000000e-05

1.51000000000000e-05

2.99000000000000e-05

7.00000000000000e-05

7.92000000000000e-05

6.06000000000000e-05

0.000163000000000000

0.000656000000000000

0.000818000000000000

0.00129000000000000

9.78000000000000e-06

0.000161000000000000

0.000183000000000000

0.000204000000000000

0.000297000000000000

0.000343000000000000

0.00705000000000000

0.00850000000000000

0.0233000000000000

0.0250000000000000

0.0267000000000000

0.00125000000000000

0.00245000000000000

0.00391000000000000

0.00878000000000000

0.00111000000000000

0.00122000000000000

3.23000000000000e-05

4.26000000000000e-05

4.53000000000000e-05

4.13000000000000e-05

8.24000000000000e-05

8.33000000000000e-05];

yobs=[5.95833333300000

0.300000000000000

0.0625000000000000

0.111111111000000

0.213809289000000

0.140625000000000

0.651315789000000

0.351694915000000

12.0555555600000

0.626846311000000

0.555555556000000

0.136363636000000

6.38793103400000

0.233051458000000

0.540000000000000

0.0277777780000000

1.05555555600000

0.113636364000000

0.0933323590000000

0.352272727000000

3.20833333300000

0.897435897000000

1.17046404700000

1.41666666700000

1.79545454500000

1.15384615400000

1.85576923100000

10.8333333300000

0.848684211000000

6.18835443000000

0.767441860000000

0.527777778000000

1.54872306000000

0.337691494000000

1.08333333300000

1.87477002000000

1.81654734900000

1.97222222200000

0.550000000000000

0.340020401000000];

xm=[x1 x2];

%% Initial parameter guesses

C1=0.25;

C2=0.73;

beta0(1)=C1; %initial guess beta 1

beta0(2)=C2; %initial guess beta 2

p=length(beta0); %p = # parameters

%% define function to be used for inverse problem

fINV=@Project_funcINV;

%fnameINV=@forderexpINV;

[beta,resids,J,COVB,mse] = nlinfit(xm,yobs,fINV,beta0);

beta

beta = 1×2

1.3585 0.1749

mdl = fitnlm(xm,yobs,fINV,beta0) % ADDED

mdl =

Nonlinear regression model: y ~ Project_funcINV(b,X) Estimated Coefficients: Estimate SE tStat pValue ________ _______ ______ __________ b1 1.3585 0.13278 10.231 1.8021e-12 b2 0.17495 0.11064 1.5812 0.12212 Number of observations: 40, Error degrees of freedom: 38 Root Mean Squared Error: 2.81 R-Squared: -0.0164, Adjusted R-Squared -0.0431 F-statistic vs. zero model: 7.46, p-value = 0.00185

figure % ADDED

stem3(x1, x2, yobs, 'filled')

hold on

plot3(x1, x2, Project_funcINV(beta,xm), '-r')

hold off

%% Functions

function y=Project_funcINV(beta0,xm)

c2s=@(x)-2.40874-39.748*(1+x).^-2.856;

y=100./(1+exp(-beta0(1).*c2s(xm(:,1))+(beta0(2).*c2s(xm(:,1)).*log10(100.*xm(:,2)))));

end

The fit is reasonably good, although ‘beta(2)’ may not be significnatly different from zero.

.

13 个评论
显示 11更早的评论隐藏 11更早的评论

Faizan Lali 2023-2-24

Hi,

I have a quick question extension to previous one. Once I calculate the 95% confidence and prediction bands for depndent variable it gives 40x1 column vector for each. So I am curious how to plot confidence interval and y_pred ? and agaisnt which "x" I should plot them ?

Is there any way out to plot them in 3 D ?

here is my code I am using to calculate confidence and prediction bands.

You help will be beneficial. Thank you!

ci=nlparci(beta, resids,J)

%R is the correlation matrix for the parameters, sigma is the standard error vector

[R,sigma]=corrcov(COVB);

R

R_2=R.*R;

%Rsq1 = 1 - sum((yobs - ypred).^2)/sum((yobs - mean(yobs)).^2);

sigma

relerr=sigma./beta' %relative error (coefficient of variance) for each parameter

%% Confidence and prediction intervals for the dependent variable

[ypred, delta] = nlpredci(fnameINV,xm,beta,resids,J,0.05,'on','curve'); %confidence band for regression line

[ypred, deltaob] =nlpredci(fnameINV,xm,beta,resids,J,0.05,'on','observation');%prediction band for individual points

yspan=range(ypred)% total span of ypred

relrmse=rmse/yspan % ratio of rmse vs. yspan

%simultaneous confidence bands for regression line

CBu=ypred+delta;

CBl=ypred-delta;

Star Strider 2023-2-24

在 MATLAB Online 中打开

You can use any ‘x’ values you want, providing that they span the same bounds as the original ‘x’ values. This means that you can use the linspace function to create them (although I do not recommend that with this data set). If you use fitnlm, you can use the appropriate predict function.

x1=[7.50000000000000

7.50000000000000

7.70000000000000

5

7.50000000000000

6

3.75000000000000

8

7.50000000000000

8

8];

x2=[0.00153000000000000

0.00522000000000000

0.000189000000000000

3.73000000000000e-06

1.17000000000000e-05

8.66000000000000e-05

1.17000000000000e-05

1.51000000000000e-05

2.99000000000000e-05

7.00000000000000e-05

7.92000000000000e-05

6.06000000000000e-05

0.000163000000000000

0.000656000000000000

0.000818000000000000

0.00129000000000000

9.78000000000000e-06

0.000161000000000000

0.000183000000000000

0.000204000000000000

0.000297000000000000

0.000343000000000000

0.00705000000000000

0.00850000000000000

0.0233000000000000

0.0250000000000000

0.0267000000000000

0.00125000000000000

0.00245000000000000

0.00391000000000000

0.00878000000000000

0.00111000000000000

0.00122000000000000

3.23000000000000e-05

4.26000000000000e-05

4.53000000000000e-05

4.13000000000000e-05

8.24000000000000e-05

8.33000000000000e-05];

yobs=[5.95833333300000

0.300000000000000

0.0625000000000000

0.111111111000000

0.213809289000000

0.140625000000000

0.651315789000000

0.351694915000000

12.0555555600000

0.626846311000000

0.555555556000000

0.136363636000000

6.38793103400000

0.233051458000000

0.540000000000000

0.0277777780000000

1.05555555600000

0.113636364000000

0.0933323590000000

0.352272727000000

3.20833333300000

0.897435897000000

1.17046404700000

1.41666666700000

1.79545454500000

1.15384615400000

1.85576923100000

10.8333333300000

0.848684211000000

6.18835443000000

0.767441860000000

0.527777778000000

1.54872306000000

0.337691494000000

1.08333333300000

1.87477002000000

1.81654734900000

1.97222222200000

0.550000000000000

0.340020401000000];

xm=[x1 x2];

%% Initial parameter guesses

C1=0.25;

C2=0.73;

beta0(1)=C1; %initial guess beta 1

beta0(2)=C2; %initial guess beta 2

p=length(beta0); %p = # parameters

%% define function to be used for inverse problem

fINV=@Project_funcINV;

%fnameINV=@forderexpINV;

[beta,resids,J,COVB,mse] = nlinfit(xm,yobs,fINV,beta0);

beta

beta = 1×2

1.3585 0.1749

mdl = fitnlm(xm,yobs,fINV,beta0) % ADDED

mdl =

Nonlinear regression model: y ~ Project_funcINV(b,X) Estimated Coefficients: Estimate SE tStat pValue ________ _______ ______ __________ b1 1.3585 0.13278 10.231 1.8021e-12 b2 0.17495 0.11064 1.5812 0.12212 Number of observations: 40, Error degrees of freedom: 38 Root Mean Squared Error: 2.81 R-Squared: -0.0164, Adjusted R-Squared -0.0431 F-statistic vs. zero model: 7.46, p-value = 0.00185

[ypredc, deltac] = nlpredci(@Project_funcINV,xm,beta,resids,'Jacobian',J,'Alpha',0.05,'SimOpt','on','PredOpt','curve'); %confidence band for regression line

% [ypredob, deltaob] =nlpredci(fnameINV,xm,beta,resids,J,0.05,'on','observation');%prediction band for individual points

CBu=ypredc+deltac;

CBl=ypredc-deltac;

figure % ADDED

stem3(x1, x2, yobs, 'filled')

hold on

plot3(x1, x2, Project_funcINV(beta,xm), '-r')

plot3(x1, x2, [CBu CBl], '--r')

hold off

%% Functions

function y=Project_funcINV(beta0,xm)

c2s=@(x)-2.40874-39.748*(1+x).^-2.856;

y=100./(1+exp(-beta0(1).*c2s(xm(:,1))+(beta0(2).*c2s(xm(:,1)).*log10(100.*xm(:,2)))));

end

.

Faizan Lali 2023-2-27

How you parameter estimates are so high?

Which code did you use?

Walter Roberson 2023-2-27

Alex uses a commercial program named 1stOpt that does some very nice optimization. Sometimes I am able to improve a little on his results, but not usually, and when I do manage then it is only after a couple of days of continuous computations.

请先登录，再进行评论。