TransportMaps.Maps.AffineTransportMapBase
¶
Module Contents¶
Classes¶
Linear map \(T({\bf x})={\bf c} + {\bf L}{\bf x}\) |
|
Linear map \(T({\bf x})={\bf c} + {\bf L}{\bf x}\) |
- class TransportMaps.Maps.AffineTransportMapBase.AffineTransportMap(**kwargs)[source]¶
Bases:
TransportMaps.Maps.AffineMapBase.AffineMap
,TransportMaps.Maps.ParametricTransportMapBase.ParametricTransportMap
Linear map \(T({\bf x})={\bf c} + {\bf L}{\bf x}\)
Note
This class supports only regular matrices. Ad-hoc implemnetations for sparse matrices should be implemented by the user.
- static build_from_Normal(pi, typeMap='sym')[source]¶
Build a linear transport map from a standard normal to a Gaussian distribution pi
- Parameters:
pi (
GaussianDistribution
) – constant term of the linear maptypeMap (str) – the linear term \(L\) is obtained as the square root of the covarinace \(\Sigma\) or precision \(\Sigma^{-1}\) matrix. For
typeMap=='sym'
, \(L=U\Lambda^{\frac{1}{2}}U^T\) where \(\Sigma = U\Lambda U^T\) is the eigenvalue decomposition of \(\Sigma\). FortypeMap=='tri'
, :maht:`L=C` where \(\Sigma=CC^T\) is the Cholesky decomposition of \(\Sigma\). FortypeMap=='kl'
, \(L=U\Lambda^{\frac{1}{2}}\) where \(\Sigma = U\Lambda U^T\) is the eigenvalue decomposition of \(\Sigma\) (this corresponds to the Karuenen-Loeve expansion). The eigenvalues and eigenvectors are ordered with \(\lambda_i\geq\lambda_{i+1}\).
- Raises:
ValueError – if the shape of linear and constant term are inconsistent.
- log_det_grad_x(x, *args, **kwargs)[source]¶
Compute: \(\log \det \nabla_{\bf x} \hat{T}({\bf x}, {\bf a})\).
- Parameters:
x (
ndarray
[\(m,d\)]) – evaluation points- Returns:
(
float
) – \(\log \det \nabla_{\bf x} \hat{T}({\bf x}, {\bf a})\) (constant at every evaluation point)- Raises:
ValueError – if \(d\) does not match the dimension of the transport map.
- grad_x_log_det_grad_x(x, *args, **kwargs)[source]¶
Compute: \(\nabla_{\bf x} \log \det \nabla_{\bf x} T({\bf x}, {\bf a})\)
- Parameters:
- Returns:
(
ndarray
[\(m,d\)]) – \(\nabla_{\bf x} \log \det \nabla_{\bf x} T({\bf x}, {\bf a})\) at every evaluation point- Raises:
ValueError – if \(d\) does not match the dimension of the transport map.
See also
- hess_x_log_det_grad_x(x, *args, **kwargs)[source]¶
Compute: \(\nabla^2_{\bf x} \log \det \nabla_{\bf x} T({\bf x}, {\bf a})\)
- Parameters:
- Returns:
(
ndarray
[\(m,d\)]) – \(\nabla^2_{\bf x} \log \det \nabla_{\bf x} T({\bf x}, {\bf a})\) at every evaluation point- Raises:
ValueError – if \(d\) does not match the dimension of the transport map.
See also
- action_hess_x_log_det_grad_x(x, dx, *args, **kwargs)[source]¶
Compute: \(\langle\nabla^2_{\bf x} \log \det \nabla_{\bf x} T({\bf x}, {\bf a}),\delta{\bf x}\rangle\)
- Parameters:
- Returns:
(
ndarray
[\(m,d\)]) – \(\nabla^2_{\bf x} \log \det \nabla_{\bf x} T({\bf x}, {\bf a})\) at every evaluation point- Raises:
ValueError – if \(d\) does not match the dimension of the transport map.
See also
- grad_x_inverse(x, *args, **kwargs)[source]¶
Compute \(\nabla_{\bf x} \hat{T}^{-1}({\bf x},{\bf a})\).
- Parameters:
x (
ndarray
[\(m,d\)]) – evaluation points- Returns:
(
ndarray
[\(d,d\)]) – gradient matrix (constant at every evaluation point).- Raises:
ValueError – if \(d\) does not match the dimension of the transport map.
- log_det_grad_x_inverse(x, *args, **kwargs)[source]¶
Compute: \(\log \det \nabla_{\bf x} T^{-1}({\bf x}, {\bf a})\).
- Parameters:
x (
ndarray
[\(m,d\)]) – evaluation points- Returns:
(float) – \(\log \det \nabla_{\bf x} T^{-1}({\bf x}, {\bf a})\) (constant at every evaluation point)
- grad_x_log_det_grad_x_inverse(x, *args, **kwargs)[source]¶
[Abstract] Compute: \(\nabla_{\bf x} \log \det \nabla_{\bf x} T^{-1}({\bf x}, {\bf a})\)
- Parameters:
- Returns:
(
ndarray
[\(m,d\)]) – \(\nabla_{\bf x} \log \det \nabla_{\bf x} T^{-1}({\bf x}, {\bf a})\) at every evaluation point
See also
- hess_x_log_det_grad_x_inverse(x, *args, **kwargs)[source]¶
[Abstract] Compute: \(\nabla^2_{\bf x} \log \det \nabla_{\bf x} T^{-1}({\bf x}, {\bf a})\)
- Parameters:
- Returns:
(
ndarray
[\(m,d,d\)]) – \(\nabla^2_{\bf x} \log \det \nabla_{\bf x} T^{-1}({\bf x}, {\bf a})\) at every evaluation point
See also
- action_hess_x_log_det_grad_x_inverse(x, dx, *args, **kwargs)[source]¶
[Abstract] Compute: \(\langle\nabla^2_{\bf x} \log \det \nabla_{\bf x} T^-1({\bf x}), \delta{\bf x}\rangle\)
- Parameters:
x (
ndarray
[\(m,d\)]) – evaluation pointsdx (
ndarray
[\(m,d\)]) – directions on which to evaluate the Hessianprecomp (
dict
) – dictionary of precomputed valuesidxs_slice (slice) – if precomputed values are present, this parameter indicates at which of the points to evaluate. The number of indices represented by
idxs_slice
must matchx.shape[0]
.
- Returns:
(
ndarray
[\(m,d\)]) – \(\langle\nabla^2_{\bf x} \log \det \nabla_{\bf x} T^-1({\bf x}), \delta{\bf x}\rangle\) at every evaluation point
See also
- class TransportMaps.Maps.AffineTransportMapBase.LinearTransportMap(c, L, Linv=None)[source]¶
Bases:
AffineTransportMap
Linear map \(T({\bf x})={\bf c} + {\bf L}{\bf x}\)
Note
This class supports only regular matrices. Ad-hoc implemnetations for sparse matrices should be implemented by the user.