BMCI

class typhon.retrieval.bmci.BMCI(y, x, s_o)[source]

Bayesian Monte Carlo Integration

This class implements methods for solving Bayesian inverse problems using Monte Carlo integration. The method uses a data base of atmospheric states \(x_i\) and corresponding observations \(\mathbf{y}_i\), which are used to compute integrals of the posterior distribution by means of importance sampling:

\[\int_a^b f(x') p(x' | \mathbf{y}) \: dx' \approx \sum_{a \leq x_i \leq b} \frac{w_i(\mathbf{y}) f(x)} {\sum_j w_j(\mathbf{y})}.\]

The method assumes that measurement uncertainties can be described by a multivariate Gaussian distribution with covariance matrix $mathbf{S}_o$, so that the weights \(w_i(\mathbf{y})\) are given by

\[w_i(\mathbf{y}) = \exp \{ -\frac{1}{2} (\mathbf{y} - \mathbf{y}_i)^T \mathbf{S}_o^{-1} (\mathbf{y} - \mathbf{y}_i) \}.\]

The measurements in the database y is assumed to be given as an array of shape (n, m), where n is the number of cases in the database. Currently only scalar retrievals are supported, which means that x is assumed to a n-element vector containing the retrieval quantities corresponding to the observations in \(\mathbf{y}\).

Attributes

x: numpy.array, shape = (n,)

The retrieval quantity corresponding to the atmospheric states represented in the data base.

y: numpy.array, shape = (n, m)

The measured or simulated brightness temperatures corresponding to the atmospheric states represented in the data base.

s_o: numpy.array, shape = (m, m)

The covariance matrix describing the measurement uncertainty.

n: int

The number of entries in the retrieval database.

m: int

The number of channels in a single measurement \(\mathbf{y}_i\).

pc1: numpy.array, shape = (m, )

The eigenvector corresponding to the smallest eigenvalue of the observation uncertainty covariance matrix. Along this vector the entries in the database will be ordered, which can be used to accelerate the retrieval.

pc1_e: float

The smallest eigenvalue of the observation uncertainty covariance matrix.

pc1_proj: numpy.array, shape (n, )

The projections of the measurements in y onto pc1 along which the database entries are ordered.

__init__(y, x, s_o)[source]

Create a QRNN instance from a given training data base y, x and measurement uncertainty given by the covariance matrix s_o.

Args

y(numpy.array): 2D array containing the measured or simulated

brightness temperatures corresponding to the atmospheric states represented in the data base. These are sorted along the eigenvector of the observation uncertainty covariance matrix with the smallest eigentvector.

x: 1D array The retrieval quantity corresponding to the atmospheric states represented in the data base. These will be sorted along the eigenvector of the observation uncertainty covariance matrix with the smallest eigentvector.

s_o: 2D array The covariance matrix describing the measurement uncertainty.

Methods

__init__(y, x, s_o)

Create a QRNN instance from a given training data base y, x and measurement uncertainty given by the covariance matrix s_o.

cdf(y_obs[, x2_max])

A posteriori cumulative distribution function (CDF).

crps(y_obs, x_true[, x2_max])

Compute the Continuous Ranked Probability Score.

pdf(y_obs[, x2_max, n_points])

A posteriori probability density function (PDF).

predict(y_obs[, x2_max])

This performs the BMCI integration to approximate the mean and variance of the posterior distribution:

predict_quantiles(y_obs, quantiles[, x2_max])

This estimates the quantiles given in quantiles by approximating the CDF of the posterior as

weights(y_obs[, x2_max])

Compute the importance sampling weights for a given observation y.