ARTS built-in documentation server

Workspace Method OEM


Inversion by the so called optimal estimation method (OEM).

Work in progress ...

The cost function to minimise, including a normalisation with lengthof y, is:
   cost = cost_y + cost_x
   cost_y = 1/m * [y-yf]' * covmat_se_inv * [y-yf]
   cost_x = 1/m * [x-xa]' * covmat_sx_inv * [x-xa]

 The current implementation provides 3 methods for the minimization of
the cost functional: Linear, Gauss-Newton and Levenberg-Marquardt.
The Gauss-Newton minimizer attempts to find a minimum solution by 
fitting a quadratic function to the cost functional. The linear minimizer
is a special case of the Gauss-Newton method, since for a linear forward
model the exact solution of the minimization problem is obtained after
the first step. The Levenberg-Marquardt method adaptively constrains the
search region for the next iteration step by means of the so-called gamma-factor.
This makes the method more suitable for strongly non-linear problems.
If the gamma-factor is 0, Levenberg-Marquardt and Gauss-Newton method
are identical. Each minimization method (li,gn,lm) has an indirect
variant (li_cg,gn_cg,lm_cg), which uses the conjugate gradient solver
for the linear system that has to be solved in each minimzation step.
This of advantage for very large problems, that would otherwise require
the computation of expensive matrix products.

Description of the special input arguments:

  "li": A linear problem is assumed and a single iteration is performed.
  "li_cg": A linear problem is assumed and solved using the CG solver.
  "gn": Non-linear, with Gauss-Newton iteration scheme.
  "gn_cg": Non-linear, with Gauss-Newton and conjugate gradient solver.
  "lm": Non-linear, with Levenberg-Marquardt (LM) iteration scheme.
  "lm_cg": Non-linear, with Levenberg-Marquardt (LM) iteration scheme and conjugate gradient solver.
  No inversion is done if the cost matching the a priori state is above
  this value. If set to a negative value, all values are accepted.
  This argument also controls if the start cost is calculated. If
  set to <= 0, the start cost in oem_diagnostics is set to NaN
  when using "li" and "gn".
  A normalisation vector for x. A normalisation of x can be needed
  due to limited numerical precision. If this vector is set to be empty
  no normalisation is done (defualt case). Otherwise, this must be a
  vector with same length as x, just having values above zero.
  Elementwise division between x and *x_norm* (x./x_norm) shall give
  a vector where all values are in the order of unity. Maybe the best
  way to set *x_norm* is x_norm = sqrt( diag( Sx ) ).
  Maximum number of iterations to perform. No effect for "li".
  Iteration stop criterion. The criterion used is the same as given
  in Rodgers' "Inverse Methods for Atmospheric Sounding"
  Settings controlling the gamma factor, part of the "LM" method.
  This is a vector of length 6, having the elements (0-based index):
    0: Start value.
    1: Fractional decrease after succesfull iteration.
    2: Fractional increase after unsuccessful iteration.
    3: Maximum allowed value. If the value is passed, the inversion
       is halted.
    4: Lower treshold. If the threshold is passed, gamma is set to zero.
       If gamma must be increased from zero, gamma is set to this value.
    5: Gamma limit. This is an additional stop criterion. Convergence
       is not considered until there has been one succesful iteration
       having a gamma <= this value.
  The default setting triggers an error if "lm" is selected.
*clear matrices*
   With this flag set to 1, jacobian and dxdy are returned as empty
   Controls if there is any screen output. The overall report level
   is ignored by this WSM.

Authors: Patrick Eriksson


OEM( x, yf, jacobian, dxdy, oem_diagnostics, lm_ga_history, oem_errors, xa, covmat_sx, y, covmat_se, jacobian_quantities, inversion_iterate_agenda, method, max_start_cost, x_norm, max_iter, stop_dx, lm_ga_settings, clear_matrices, display_progress )


OUT+INx(Vector)The state vector.
OUT+INyf(Vector)A fitted measurement vector.
OUT+INjacobian(Matrix)The Jacobian matrix.
OUTdxdy(Matrix)Contribution function (or gain) matrix.
OUToem_diagnostics(Vector)Basic diagnostics of an OEM type inversion.
OUTlm_ga_history(Vector)The series of gamma values for a Marquardt-levenberg inversion.
OUToem_errors(ArrayOfString)Errors encountered during OEM execution.
INxa(Vector)The a priori state vector.
INcovmat_sx(CovarianceMatrix)Covariance matrix of a priori distribution This covariance matrix describes the Gaussian a priori distribution for an OEM retrieval.
INy(Vector)The measurement vector.
INcovmat_se(CovarianceMatrix)Covariance matrix for observation uncertainties.
INjacobian_quantities(ArrayOfRetrievalQuantity)The retrieval quantities in the Jacobian matrix.
INinversion_iterate_agenda(Agenda)Agenda recalculating spectra and Jacobian for iterative inversion methods.
GINmethod(String)Iteration method. For this and all options below, see further above.
GINmax_start_cost(Numeric, Default: Inf)Maximum allowed value of cost function at start.
GINx_norm(Vector, Default: [])Normalisation of Sx.
GINmax_iter(Index, Default: 10)Maximum number of iterations.
GINstop_dx(Numeric, Default: 0.01)Stop criterion for iterative inversions.
GINlm_ga_settings(Vector, Default: [])Settings associated with the ga factor of the LM method.
GINclear_matrices(Index, Default: 0)An option to save memory.
GINdisplay_progress(Index, Default: 0)Flag to control if inversion diagnostics shall be printed on the screen.