dubfi.inversion.inversion_lognormal

Inversion for log-normal a priori probability distribution.

Added in version 0.1.3.

Changed in version 0.1.4: add hybrid log-normal + normal solver

Classes

InvertorOptimizerLogNormal	Numerical solver for flux inversion problem with log-normal prior PDF.
InvertorOptimizerHybrid	Numerical solver for flux inversion problem with normal and log-normal components.

Functions

lognormal2normal(→ tuple[numpy.ndarray, numpy.ndarray])	Given mu and B of a log-normal distribution, find Gaussian approximation at its maximum.
normal2lognormal(→ tuple[numpy.ndarray, numpy.ndarray])	Compute inverse of lognormal2normal().
lognormal2meanvar(→ tuple[numpy.ndarray, numpy.ndarray])	Given mu and B of a log-normal distribution, compute mean and variance (for Gaussian approximation).
meanvar2lognormal(→ tuple[numpy.ndarray, numpy.ndarray])	Given mean and variance of a random vector, construct parameters mu and B of log-normal distribution.

Module Contents

class dubfi.inversion.inversion_lognormal.InvertorOptimizerLogNormal(y: dubfi.linalg.types.AbstractVector, prior: numpy.ndarray, h: dubfi.linalg.generic.ParametrizedVector, b: numpy.ndarray, r: dubfi.linalg.generic.ParametrizedOperator, norm_prefactor: float = 0.0, min_s: float = 1e-06)

Bases: dubfi.inversion.inversion.InvertorOptimizer

Numerical solver for flux inversion problem with log-normal prior PDF.

Solve inversion problem by numerically minimizing the cost function, assuming that uncertainties may depend on the optimization parameters. In this case, the control parameter ls is the logarithm of the scaling factors.

Solve inverse problem using full ensemble-estimated model uncertainty.

Parameters:

• y (AbstractVector) – vector of observations (usually minus model far-field or background)

• prior (np.ndarray) – a priori estimate for scaling factors

• h (ParametrizedVector) – model equivalents depending on scaling factors

• b (np.ndarray) – error covariance matrix of prior scaling factors

• r (ParametrizedOperator) – error covariance matrix of model-observation comparison, depends on scaling factors

• norm_prefactor (float, default=0.0) – prefactor of normalization term in cost function.

• min_s (float, default=1e-6) – minimum allowed scaling factor

property b: numpy.ndarray

Prior error covariance matrix of scaling factors.

property prior: numpy.ndarray

Prior vector of scaling factors.

_check_before_solving() → None

Assert that prior and initial scaling factors are positive.

cost_prior(s: numpy.ndarray) → numpy.float64

Prior contribution to cost function.

\[J_\text{prior}(s) = \tfrac{1}{2} \sum_{lk} [\log(s_l)-\mu_l] (\tilde{B}^{-1})_{lk} [\log(s_k) - \mu_k] + \sum_k [\log(s_k) - \mu_k]\]

This defines a log-normal probability distribution \(P(s) \propto e^{-J_\text{prior}(s)}\) such that the element-wise logarithm of the scaling factors, \(z_k\equiv\log(s_k)\), has a probability density \(z \sim \mathcal{N}(\mu, \tilde{B})\).

The parameters \(\mu\) and \(\tilde{B}\) are defined such that the a priori log-normal probability distribution \(P(s) \propto e^{-J_\text{prior}(s)}\) has the given most likely scaling factors (\(s_0\)) and the local Gaussian approximation (\(B\)):

\[ \begin{align}\begin{aligned}(s_0)_i &= e^{\mu_i+\tfrac{1}{2}\tilde{B}_{ii}},\\B_{ij} &= e^{\mu_i + \mu_j + \tfrac{1}{2}(\tilde{B}_{ii}+\tilde{B}_{jj})} \left(e^{\tilde{B}_{ij}} - 1\right)\end{aligned}\end{align} \]

See also

lognormal2normal()

cost_prior_grad(s: numpy.ndarray) → numpy.ndarray

Prior contribution to gradient of cost function.

\[\partial_k J_\text{prior}(s) = \sum_l [\log(s_l) - \mu_l] (\tilde{B}^{-1})_{lk} s_k^{-1} + s_k^{-1}\]

cost_prior_hess(s: numpy.ndarray) → numpy.ndarray

Prior contribution to Hessian of cost function.

\[\partial_k \partial_l J_\text{prior}(s) = \frac{(\tilde{B}^{-1})_{lk}}{s_l s_k} - \delta_{lk} \left( 1 + \sum_m [\log(s_m) - \mu_m] (\tilde{B}^{-1})_{ml} \right) s_l^{-2}\]

_approximate_cost_hess(s: numpy.ndarray) → bool

Avoid computing expensive term in Hessian, use if prior PDF is dominant.

class dubfi.inversion.inversion_lognormal.InvertorOptimizerHybrid(y: dubfi.linalg.types.AbstractVector, prior: numpy.ndarray, h: dubfi.linalg.generic.ParametrizedVector, b: numpy.ndarray, r: dubfi.linalg.generic.ParametrizedOperator, lognorm_mask: numpy.ndarray, norm_prefactor: float = 0.0, min_s: float = 1e-06)

Bases: dubfi.inversion.inversion.InvertorOptimizer

Numerical solver for flux inversion problem with normal and log-normal components.

The state space is split into two subspaces, defined by the parameter lognorm_mask. One part has a log-normal prior PDF, the other one has a normal (Gaussian) prior PDF. Both parts must be uncorrelated in the prior PDF. The equations for the prior PDF of these two parts are take from (and documented in) InvertorOptimizerLogNormal and InvertorOptimizer.

Added in version 0.1.4.

Solve inverse problem using full ensemble-estimated model uncertainty.

Parameters:

• y (AbstractVector) – vector of observations (usually minus model far-field or background)

• prior (np.ndarray) – a priori estimate for scaling factors

• h (ParametrizedVector) – model equivalents depending on scaling factors

• b (np.ndarray) – error covariance matrix of prior scaling factors

• r (ParametrizedOperator) – error covariance matrix of model-observation comparison, depends on scaling factors

• lognorm_mask (np.ndarray) – boolean, 1d numpy array. Log-normal prior PDF is assumed for all entries marked as True in this array, normal PDF is assumed otherwise.

• norm_prefactor (float, default=0.0) – prefactor of normalization term in cost function.

• min_s (float, default=1e-6) – minimum allowed scaling factor

property b: numpy.ndarray

Prior error covariance matrix of scaling factors.

property prior: numpy.ndarray

Prior vector of scaling factors.

_check_before_solving() → None

Assert that prior and initial scaling factors are positive.

cost_prior(s: numpy.ndarray) → numpy.float64

Prior contribution to cost function.

cost_prior_grad(s: numpy.ndarray) → numpy.ndarray

Prior contribution to gradient of cost function.

cost_prior_hess(s: numpy.ndarray) → numpy.ndarray

Prior contribution to Hessian of cost function.

_approximate_cost_hess(s: numpy.ndarray) → bool

Avoid computing expensive term in Hessian, use if prior PDF is dominant.

dubfi.inversion.inversion_lognormal.lognormal2normal(mu: numpy.ndarray, b: numpy.ndarray) → tuple[numpy.ndarray, numpy.ndarray]

Given mu and B of a log-normal distribution, find Gaussian approximation at its maximum.

Consider a random vector \(x\) and denote by \(z\) its element-wise logarithm, \(z_i=\log(x_i)\). Assume that \(z\sim\mathcal{N}(\mu, B)\). Then \(P(x) \propto \exp\left(-\sum_k z_k - \tfrac{1}{2} [z-\mu]^\top B^{-1} [z-\mu] \right)\). Approximate \(P(x)\) locally at its maximum:

\[ \begin{align}\begin{aligned}x^\text{max}_i &= e^{\mu_i - \sum_j B_{ij}},\\(B_\text{Gauss}^{-1})_{ij} &= -\partial_i \partial_j \log(P(x))|_{x^\text{max}} = \frac{B_{ij}}{x^\text{max}_i x^\text{max}_j}\end{aligned}\end{align} \]

Return \(x^\text{max}\) and \(B_\text{Gauss}^{-1}\).

Note

This only makes sense if \(B\) is positive definite and real-symmetric.

dubfi.inversion.inversion_lognormal.normal2lognormal(xmax: numpy.ndarray, b_gauss: numpy.ndarray) → tuple[numpy.ndarray, numpy.ndarray]

Compute inverse of lognormal2normal().

Given a local Gaussian approximation of a log-normal distribution at its maximum, return parameters mu and B of the log-normal distribution.

Note

All entries of xmax must be positive. b_gauss must be a positive definite, real-symmetric matrix.

dubfi.inversion.inversion_lognormal.lognormal2meanvar(mu: numpy.ndarray, b: numpy.ndarray) → tuple[numpy.ndarray, numpy.ndarray]

Given mu and B of a log-normal distribution, compute mean and variance (for Gaussian approximation).

Consider a random vector \(x\) and denote by \(z\) its element-wise logarithm, \(z_i=\log(x_i)\). Assume that \(z\sim\mathcal{N}(\mu, B)\). Then \(P(x) \propto \exp\left(-\sum_k z_k - \tfrac{1}{2} [z-\mu]^\top B^{-1} [z-\mu] \right)\) and \(\langle e^{t^\top z} \rangle = e^{t^\top\mu + \tfrac{1}{2} t^\top B t}\). By choosing \(t=e_i\) or \(t=e_i+e_j\), one can show that the mean and variance of \(x\) are:

\[ \begin{align}\begin{aligned}\langle x_i \rangle &= e^{\mu_i+\tfrac{1}{2}B_{ii}},\\\langle\langle x_i x_j \rangle\rangle &= e^{\mu_i + \mu_j + \tfrac{1}{2}(B_{ii}+B_{jj})} \left(e^{B_{ij}} - 1\right)\end{aligned}\end{align} \]

This function returns the mean and variance of \(x\).

Note

This only makes sense if \(B\) is positive definite and real-symmetric.

dubfi.inversion.inversion_lognormal.meanvar2lognormal(mean: numpy.ndarray, var: numpy.ndarray) → tuple[numpy.ndarray, numpy.ndarray]

Given mean and variance of a random vector, construct parameters mu and B of log-normal distribution.

This is the inverse of lognormal2normal().

Note

All entries of mean must be positive. var must be a positive definite, real-symmetric matrix.