nlparci

Nonlinear regression parameter confidence intervals

Syntax

ci = nlparci(beta,resid,'covar',sigma)
ci = nlparci(beta,resid,'jacobian',J)
ci = nlparci(...,'alpha',alpha)

Description

ci = nlparci(beta,resid,'covar',sigma) returns the 95% confidence intervals ci for the nonlinear least squares parameter estimates beta. Before calling nlparci, use nlinfit to fit a nonlinear regression model and get the coefficient estimates beta, residuals resid, and estimated coefficient covariance matrix sigma.

ci = nlparci(beta,resid,'jacobian',J) is an alternative syntax that also computes 95% confidence intervals. J is the Jacobian computed by nlinfit. If the 'robust' option is used with nlinfit, use the 'covar' input rather than the 'jacobian' input so that the required sigma parameter takes the robust fitting into account.

ci = nlparci(...,'alpha',alpha) returns 100(1-alpha)% confidence intervals.

nlparci treats NaNs in resid or J as missing values, and ignores the corresponding observations.

The confidence interval calculation is valid for systems where the length of resid exceeds the length of beta and J has full column rank. When J is ill-conditioned, confidence intervals may be inaccurate.

Examples

Fit to exponential decay

Suppose you have data, and want to fit a model of the form

yi = a1 + a2exp(–a3xi) + εi.

Here the ai are the parameters you want to estimate, xi are the data points, the yi are the responses, and the εi are noise terms.

  1. Write a function handle that represents the model:

    mdl = @(a,x)(a(1) + a(2)*exp(-a(3)*x));
  2. Generate synthetic data with parameters a = [1;3;2], with the x data points distributed exponentially with parameter 2, and normally distributed noise with standard deviation 0.1:

    rng(9845,'twister') % for reproducibility
    a = [1;3;2];
    x = exprnd(2,100,1);
    epsn = normrnd(0,0.1,100,1);
    y = mdl(a,x) + epsn;
  3. Fit the model to data starting from the arbitrary guess a0 = [2;2;2]:

    a0 = [2;2;2];
    [ahat,r,J,cov,mse] = nlinfit(x,y,mdl,a0);
    ahat
    
    ahat =
        1.0153
        3.0229
        2.1070
  4. Check whether [1;3;2] is in a 95% confidence interval using the Jacobian argument in nlparci:

    ci = nlparci(ahat,r,'Jacobian',J)
    
    ci =
        0.9869    1.0438
        2.9401    3.1058
        1.9963    2.2177
  5. You can obtain the same result using the covariance argument:

    ci = nlparci(ahat,r,'covar',cov)
    
    ci =
        0.9869    1.0438
        2.9401    3.1058
        1.9963    2.2177

See Also

|

Was this topic helpful?