>  Docs Center  >  Libraries  >  Markwardt  >  MPPROPERR






  Craig B. Markwardt, NASA/GSFC Code 662, Greenbelt, MD 20770
  UPDATED VERSIONs can be found on my WEB PAGE:


  Propagate fitted model uncertainties to measured data points

Major Topics

  Curve and Surface Fitting

Calling Sequence



  MPPROPERR propagates the parameter uncertainties of a fitted
  model to provide estimates of the model uncertainty at each
  measurement point.
  When fitting a model to data with uncertainties, the parameters
  will have estimated uncertainties. In fact, the parameter
  variance-covariance matrix indicates the estimated uncertainties
  and correlations between parameters. These uncertainties and
  correlations can, in turn, be used to estimate the "error in the
  model" for each measurement point. In a sense, this quantity also
  reflects the sensitivity of the model to each data point.
  The algorithm used by MPPROPERR uses standard propagation of error
  techniques, assuming that errors are small. The input values of
  MPPROPERR should be found from the output keywords of MPFIT or
  MPFITFUN, as documented below.
  The user has a choice whether to compute the *full*
  variance-covariance matrix or not, depending on the setting of the
  DIAGONAL keyword. The full matrix is large, and indicates the
  correlation the sampled model function between each measurement
  point and every other point. The variance terms lie on the
  diagonal, and the covariance terms are on the off-diagonal.
  Usually however, the user will want to set /DIAGONAL, which only
  returns the "diagonal" or variance terms, which represent the
  model "uncertainty" at each measurement point. The /DIAGONAL
  setting only controls the amount of data returned to the user.
  the full *parameter* covariance matrix is always used to compute
  the output regardless of the setting for /DIAGONAL.
  When using MPPROPERR, keep in mind the following dimensions of
  the problem:
    NPOINTS - number of measurement points
    NPAR - total number of fit parameters
    NFREE - number of *free* fit parameters
  The inputs to this function are:
    BEST_FJAC - the partial derivative matrix, or Jacobian matrix,
                as estimated by MPFIT or MPFITFUN (see below),
                which has dimensions of ARRAY(NPOINTS,NFREE).
    PCOVAR - the parameter covariance matrix, as estimated by MPFIT
              or MPFITFUN (see below), which has dimensions of
    PFREE_INDEX - an index array which describes which of the
                  parameter set were variable, as returned by MPFIT
                  or MPFITFUN. Of the total parameter set PARMS,
                  only PARMS[PFREE_INDEX] were varied by MPFIT.
  There are special considerations about the values returned by
  MPPROPERR. First, if a parameter is touching a boundary
  limit when the fit is complete, then it will be marked as having
  variance and covariance of zero. To avoid this situation, one can
  re-run MPFIT or MPFITFUN with MAXITER=0 and boundary limits
  disabled. This will permit MPFIT to estimate variance and
  covariance for all parameters, without allowing them to actually
  vary during the fit.
  Also, it is important to have a quality parameter covariance
  matrix PCOVAR. If the matrix is singular or nearly singular, then
  the measurement variances and covariances will not be meaningful.
  It helps to parameterize the problem to minimize parameter
  covariances. Also, consider fitting with double precision
  quantities instead of single precision to minimize the chances of
  round-off error creating a singular covariance matrix.
  IMPORTANT NOTE: the quantities returned by this function are the
  *VARIANCE* and covariance. If the user wishes to compute
  estimated standard deviation, then one should compute
  SQRT(VARIANCE). (see example below)


  BEST_FJAC - the Jacobian matrix, as estimated by MPFIT/MPFITFUN
              (returned in keyword BEST_FJAC). This should be an
              array ARRAY(NPOINTS,NFREE) where NFREE is the number
              of free parameters.
  PCOVAR - the full parameter covariance matrix, as returned in the
            COVAR keyword of MPFIT/MPFITFUN. This should be an array
            ARRAY(NPAR,NPAR) where NPAR is the *total* number of


  The estimated uncertainty at each measurement point, due to
  propagation of errors. The dimensions depend on the value of the
  DIAGONAL keyword.
    DIAGONAL=1: returned value is ARRAY(NPOINTS)
                corresponding to the *VARIANCE* of the model
                function sampled at each measurment point
                **NOTE**: the propagated standard deviation would
                then be SQRT(RESULT).
    DIAGONAL=0: returned value is ARRAY(NPOINTS,NPOINTS)
                corresponding to the variance-covariance matrix of
                the model function, sampled at the measurement

Keyword Parameters

  DIAGONAL - if set, then compute only the "diagonal" (variance)
              terms. If not set, then propagate the full covariance
              matrix for each measurement point.
  NAN - if set, then ignore NAN values in BEST_FJAC or PCOVAR
        matrices (they would be set to zero).
  PFREE_INDEX - index list of free parameters, as returned in the
                PFREE_INDEX keyword of MPFIT/MPFITFUN. This should
                be an integer array ARRAY(NFREE), such that
                parameters PARMS[PFREE_INDEX] were freely varied during
                the fit, and the remaining parameters were not.
                Thus it should also be the case that PFREE_INDEX
                indicates the rows and columns of the parameter
                covariance matrix which were allowed to vary freely.
                Default: All parameters will be considered free.


  ; First, generate some synthetic data
  npts = 200
  x = dindgen(npts) * 0.1 - 10. ; Independent variable
  yi = gauss1(x, [2.2D, 1.4, 3000.]) ; "Ideal" Y variable
  y = yi + randomn(seed, npts) * sqrt(1000. + yi); Measured, w/ noise
  sy = sqrt(1000.D + y) ; Poisson errors
  ; Now fit a Gaussian to see how well we can recover
  p0 = [1.D, 1., 1000.] ; Initial guess (cent, width, area)
  p = mpfitfun('GAUSS1', x, y, sy, p0, $ ; Fit a function
                best_fjac=best_fjac, pfree_index=pfree_index, /calc_fjac, $
  ; Above statement calculates best Jacobian and parameter
  ; covariance matrix
  ; Propagate errors from parameter covariance matrix to estimated
  ; measurement uncertainty. The /DIAG call returns only the
  ; "diagonal" (variance) term for each measurement.
  ycovar = mpproperr(best_fjac, pcovar, pfree_index=pfree_index, /diagonal)
  sy_prop = sqrt(ycovar) ; Estimated sigma


  MINPACK-1, Jorge More', available from netlib (www.netlib.org).
  "Optimization Software Guide," Jorge More' and Stephen Wright,
    SIAM, *Frontiers in Applied Mathematics*, Number 14.

Modification History

  Written, 2010-10-27, CM
  Updated documentation, 2011-06-26, CM

© 2020 Harris Geospatial Solutions, Inc. |  Legal
My Account    |    Store    |    Contact Us