confidence interval of non linear fit of multiple data set
3 views (last 30 days)
Show older comments
Hi everyone,
I am trying to estimate the confidence interval of my fitting parameters. I have written the code with the constraints and it works properly but now I want to estimate the confidence intervals for the fitting parameters.
I decided to use bootci but I am not able to include this function in my code.
Do you have any suggestion?
The code is attached to the question:
clear
%fit per tln+glu under Tg
R=8.314;
T1=220; %temperature reffered to y1
T2=300; %temperature reffered to y2
yfcn1 = @(b,x) b(1).*exp(-x(:,2).^2.*b(2)).*(1-2*(exp(-b(3)./(R.*x(:,1))+b(4)/R)./(1+exp(-b(3)./(R.*x(:,1))+b(4)/R)).^2).*(1-sin(x(:,2).*b(5))./(x(:,2)*b(5))));
yfcn2 = @(b,x) b(6).*exp(-x(:,2).^2.*b(7)).*(1-2*(exp(-b(3)./(R.*x(:,1))+b(4)/R)./(1+exp(-b(3)./(R.*x(:,1))+b(4)/R)).^2).*(1-sin(x(:,2).*b(5))./(x(:,2)*b(5))));
x=[0.5215 0.7756 1.2679 1.4701 1.6702 1.8680 2.0633 2.2693 2.4584 2.6442 2.8264 3.0046 3.0890 3.2611 3.4287 3.5917 3.7497 3.9309 4.0774 4.2183 4.3535 4.4827 4.5427 4.6628];
y1=[0.9936 0.9375 0.9081 0.8648 0.8568 0.8114 0.7711 0.8010 0.7884 0.7389 0.7901 0.7825 0.7903 0.7501 0.7070 0.7489 0.6441 0.7105 0.6735 0.6385 0.6357 0.6962 0.5946 0.6783];
y1_err= [ 0.0637 0.0526 0.0330 0.0235 0.0298 0.0223 0.0388 0.0223 0.0333 0.0326 0.0410 0.0282 0.0561 0.0235 0.0271 0.0218 0.0333 0.0252 0.0344 0.0261 0.0499 0.0396 0.0655 0.0901];
y2=[0.8748 0.8726 0.7922 0.7782 0.7396 0.6958 0.6603 0.6503 0.6556 0.6461 0.6021 0.5820 0.6220 0.5768 0.4950 0.5300 0.5234 0.5170 0.4369 0.4508 0.4409 0.4392 0.4100 0.6699];
y2_err=[ 0.0562 0.0480 0.0287 0.0211 0.0260 0.0194 0.0339 0.0188 0.0287 0.0289 0.0332 0.0225 0.0460 0.0191 0.0211 0.0169 0.0280 0.0198 0.0256 0.0204 0.0392 0.0283 0.0504 0.0856];
format long E
T1=220; %temperature reffered to y1
T2=300; %temperature reffered to y2
T1v = T1*ones(size(x));
T2v = T2*ones(size(x));
%T3v = T3*ones(size(x));
xm = x(:)*ones(1,2);
ym = [y1(:) y2(:)];%, y3(:)];
Tm = [T1v(:) T2v(:)];% T3v(:) ];
yerr=[y1_err(:) y2_err(:)];% y3_err(:)];
xv = xm(:);
yv = ym(:);
Tv = Tm(:);
yerrv=yerr(:);
weights=1./yerrv;
xTm = [Tv xv];
%B0 = randn(7,1)*0.1;
B0=[0.5 1e-3 0.5 1e-3 5e4 45 1.5]';
yfcn = @(b,x) yfcn1(b,x).*(x(:,1)==T1) + yfcn2(b,x).*(x(:,1)==T2);
fitfcnw = @(b) norm(weights.*(yv - yfcn(b,xTm)));
lb=[0,0,1e3,10,0,0,0];
ub=[1,1,4e5,1000,2,1,1];
A = @simple_constraint;
problem = createOptimProblem('fmincon', 'x0',B0,'nonlcon',A, 'objective',fitfcnw,'lb',lb,'ub',ub);%
gs = GlobalSearch('PlotFcns',@gsplotbestf);
[B,fval] = run(gs,problem)
B(:)
format short eng
figure(1)
for k = 1:2
idx = (1:numel(x))+numel(x)*(k-1);
subplot(2,1,k)
errorbar(x.^2, ym(:,k),yerr(:,k), '.')
hold on
plot(x.^2, yfcn(B,xTm(idx,:)), '-r')
hold off
grid
ylabel('Substance [Units]')
title(sprintf('y_{%d}, T = %d', k,xTm(idx(1),1)))
ylim([min(yv) max(yv)+0.2])
if k == 1
text(5, max(yv)+0.1, sprintf('$y = %.3f\\times e^{-x^2\\times %.3f}(1-2\\frac{e^{\\frac{%.3f}{%3d\\times R}-\\frac{%.3f}{R}}}{(1+e^{\\frac{%.3f}{%3d\\times R}-\\frac{%.3f}{R}})^2}(1-\\frac{sin(%.3fx}{%.3fx}))$',B(1:3),T1,B(4),B(3),T1,B(4:5),B(5)), 'Interpreter','latex', 'FontSize',12)
elseif k == 2
text(5, max(yv)+0.1, sprintf('$y = %.3f\\times e^{-x^2\\times %.3f}(1-2\\frac{e^{\\frac{%.3f}{%3d\\times R}-\\frac{%.3f}{R}}}{(1+e^{\\frac{%.3f}{%3d\\times R}-\\frac{%.3f}{R}})^2}(1-\\frac{sin(%.3fx}{%.3fx}))$',B(6:7),B(3),T2,B(4),B(3),T2,B(4:5),B(5)), 'Interpreter','latex', 'FontSize',12)
end
end
xlabel('Q^2')
0 Comments
Answers (0)
See Also
Categories
Find more on Linear and Nonlinear Regression in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!