How can I loop a curve fitting for varying parameters?

I would like to vary some input parameters (ccBHe and cc0He) in the curve fitting procedure but up to now it doesn't work because of a matrix miss match.

po=0.2;
T=2.4;
DHe=5.7*10^-8*po/T; % diffusion coefficient
%ccBHe=3.224*10^-5;  % initial Mantle concentration

xdata=[21200,18300,2000,700,84,33,23600,23600,15900,16800,48,7,25000,26000,21200,29800,21200,30800,27800];
ydata=[2.55E-05,2.77E-05,4.70E-05,1.28E-04,6.28E-05,2.41E-04,3.54E-05,4.81E-05,2.71E-05,2.68E-05,8.51E-05,4.13E-04,2.50E-05,2.6e-5,3e-5,2.95e-5,1.89e-5,8.83e-5,2.66e-5];

options = optimoptions('lsqcurvefit','Algorithm','trust-region-reflective','DiffMinChange',1e-8,'TolFun',1e-66,'MaxFunEvals',6000,'TolX',0.001);

%monte carlo simulation for error progagation calculation

%mantle
ccBHe=normrnd(4.25*10^-5,4.27*10^-6);          %normal distribution

%groundwater 
cc0He=normrnd(4.3*10^-4,3.6*10^-5);

fun=@(x,xdata)cc0He-(cc0He-ccBHe)*erf(xdata/(2*(DHe*x)^(1/2)));

x0=2.5e11;

x=lsqcurvefit(fun,x0,xdata,ydata,[],[],options);

figure(1)
distance=linspace(-1000,32000);
plot(xdata,ydata,'ko',distance,fun(x,distance),'b','linewidth',2)
xlim([-1000 35000])
ylim([0.1*10^-4 5*10^-4])

Answer 1

[Star Strider]

If you want to vary ‘ccBHe’ and ‘cc0He’ randomly in a loop, try this:

%monte carlo simulation for error progagation calculation

distance=linspace(-1000,32000);         % Be Careful Here — The Mapping Toolbox Has A ‘distance’ Function

for k1 = 1:10;
    %mantle
    ccBHe=normrnd(4.25*10^-5,4.27*10^-6);          %normal distribution

      %groundwater
      cc0He=normrnd(4.3*10^-4,3.6*10^-5);

      fun=@(x,xdata) cc0He-(cc0He-ccBHe).*erf(xdata./(2*(DHe*x).^(1/2)));

      x0=2.5e11;

      x(k1)=lsqcurvefit(fun,x0,xdata,ydata,[],[],options);                % Store Parameter Estimates
      ccHe(k1,:) = [ccBHe,  cc0He];                                       % Store ‘ccHe’ Values

      sim_fun(k1,:) = fun(x(k1),distance);                                % Simulate Function, Save Results
  end

figure(1)
plot(xdata,ydata,'ko',distance,sim_fun,'b','linewidth',2)
xlim([-1000 35000])
ylim([0.1*10^-4 5*10^-4])

It does the simulation in each iteration and stores the results, and plots all of them at the end. It also stores the ‘ccHe’ values in their own matrix, with the rows corresponding to the colums of ‘x’.

[Carmen]

Perfect. thank you very much. juhe

[Star Strider]

My pleasure.

The sincerest form of thanks here on MATLAB Answers is to Accept the Answer that solves your problem.