I'm solving a problem (page 16 is in English) in numerical analysis and this is the solution:

>> z=[0:500:4000 5000:1000:12000];
>> data=[5050 4980 4930 4890 4870 4865 4860 4860 4865 4875 4885 4905 4920 4935 4950 4970 4990];
>> fun=@(p)(4800 + p(1))*ones(size(z)) +p(2)/1000*z+p(3)*exp(p(4)/1000*z)-data;
>> x0=[0 0 0 -1];
>> opt = optimset('MaxFunEvals',1000);
>> p=lsqnonlin(fun,x0,[],[],opt);

Local minimum possible.

lsqnonlin stopped because the final change in the sum of squares relative to 
its initial value is less than the default value of the function tolerance.

<stopping criteria details>

>> fitf=@(t)(4800 + p(1))*ones(size(t)) + p(2)/1000*t+ p(3)*exp(p(4)/1000*t);
>> tt=linspace(0,12000,1000);
>> plot(z,data,'r-',tt,fitf(tt),'b-');
>> p

p =

  -20.2090   17.3368  272.9057   -0.7528


On the plot it looks alright:

enter image description here

But when I asked about it the coefficients to the polynomial in the answer were some other values:

ode45 usage in this case?

In my program run, the coeff. appear to be p4=-0.7528 and the other are p1=-20.2090, p2=17.3368, p3=272.9057.

So c(z) becomes




And then I can make the substitution z’=u and z’’=u’ so the matlab code to continue with will be:

function dZ=sys(x,Z)
    c=@(z)4800 -20.2090 + (17.3368)*z/1000+ (272.9057)*exp(0.7528z/1000); % c(z)

% Z(1):=z
% Z(2):=u
    dZ=zeros(2,1);    % a column vector
    dZ(2)=(c/cosd(7.8))^2*((1287*exp(Z(1)/1000))/5000 + 729/50000)/...
        (4800 + 5.66 + (14.58)*Z(1)/1000+ (257.4)*exp(-Z(1)/1000))^3;

Do you agree?


Both your coefficients and the coefficients in the answer are correct but yours are better.

The information of fit for the coefficients in the answer:

Chi-squared: 813.032070619
R-squared: 0.983460126607
R-squared adjusted: 0.981097287551
Root Mean Squared Error (RMSE): 6.91559223196
p1 = 5.6596316127433939E+00
     std err squared: 2.17628E+01
     t-stat: 1.21319E+00
     p-stat: 2.45129E-01
     95% confidence intervals: [-4.34594E+00, 1.56652E+01]
p2 = 1.4581975826123806E+01
     std err squared: 4.22221E-01
     t-stat: 2.24412E+01
     p-stat: 2.24487E-12
     95% confidence intervals: [1.31883E+01, 1.59756E+01]
p3 = 2.5743180704436827E+02
     std err squared: 8.25256E+01
     t-stat: 2.83379E+01
     p-stat: 9.14824E-14
     95% confidence intervals: [2.37948E+02, 2.76916E+02] 

The information of fit for your coefficients:

Chi-squared: 123.211540448
R-squared: 0.997493452776
R-squared adjusted: 0.996915018801
Root Mean Squared Error (RMSE): 2.69216226729
p1 = -2.0209016655555111E+01
     std err squared: 1.78565E+01
     t-stat: -4.78242E+00
     p-stat: 3.57957E-04
     95% confidence intervals: [-2.93381E+01, -1.10800E+01]
p2 = 1.7336770761968697E+01
     std err squared: 2.24807E-01
     t-stat: 3.65648E+01
     p-stat: 1.70974E-14
     95% confidence intervals: [1.63125E+01, 1.83611E+01]
p3 = 2.7290573753852300E+02
     std err squared: 1.98905E+01
     t-stat: 6.11913E+01
     p-stat: 0.00000E+00
     95% confidence intervals: [2.63271E+02, 2.82541E+02]
p4 = 7.5277833016177342E-01
     std err squared: 6.80991E-04
     t-stat: 2.88467E+01
     p-stat: 3.59046E-13
     95% confidence intervals: [6.96402E-01, 8.09155E-01] 

These fits calculated by zunzun

Additionally i think you missed a minus before p4 in c(z) and you should change the coefficients in dZ(2) with your coeffs.

  • $\begingroup$ Thank you for the answer. Can there be even more solution? Why can there be several solutions? How do I justify which coefficients I use? I will work through the program and run it once more to be sure that the calculations are correct. $\endgroup$
    – Niklas R.
    Sep 18 '13 at 14:53
  • 2
    $\begingroup$ Actually the two coefficient sets (in the answer and yours) are not for the same problem namely not solutions of same fitting problem. Because we actually ignore $p4$ by taking $p4=1$ and then find the other coefficients. But you, you find four coefficients. Although i am not expert in statistics as i know when you fit a curve or a line to a given data set you check the value of "R-squared" to measure the goodness of the fit. This value is between $0$ and $1$. As $R^2$ approach to $1$ your fit gets better. As you see the coefficients you have are better for the fit. $\endgroup$
    – Ömer
    Sep 18 '13 at 15:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.