KarhunenLoeveSVDAlgorithm

(Source code, png, hires.png, pdf)

../../_images/KarhunenLoeveSVDAlgorithm.png
class KarhunenLoeveSVDAlgorithm(*args)

Computation of Karhunen-Loeve decomposition using SVD approximation.

Available constructors:

KarhunenLoeveSVDAlgorithm(sample, s, centeredFlag)

KarhunenLoeveSVDAlgorithm(sample, verticesWeights, s, centeredFlag)

KarhunenLoeveSVDAlgorithm(sample, verticesWeights, sampleWeights, s, centeredFlag)

Parameters:

sample : ProcessSample

The sample containing the observations.

verticesWeights : sequence of float

The weights associated to the vertices of the mesh defining the sample.

sampleWeights : sequence of float

The weights associated to the fields of the sample.

s : float, \geq0

The threshold used to select the most significant eigenmodes, defined in KarhunenLoeveAlgorithm.

centeredFlag : logical

Flag to tell if the sample is drawn according to a centered process or if it has to be centered using the empirical mean. The default value is False.

Notes

The Karhunen-Loeve SVD algorithm solves the Fredholm problem associated to the covariance function C: see KarhunenLoeveAlgorithm to get the notations.

The SVD approach is a particular case of the quadrature approach (see KarhunenLoeveQuadratureAlgorithm) where we consider the functional basis ((\vect{\theta}_p^j(\vect{s}))_{1 \leq j \leq d, 1 \leq p \leq P} of L^2(\cD, \Rset^d) defined on \cD by:

\vect{\theta}_p^j(\vect{s})= \left[\mat{C}(\vect{s}, \vect{s}_p) \right]_{:, j}

The SVD approach is used when the covariance function is not explicitely known but only through K fields of the associated stochastic process \vect{X}: (\vect{X}_1, \dots, \vect{X}_K).

It consists in :

  • approximating \mat{C} by its empirical estimator \dfrac{1}{\tilde{K}} \mat{X}\, \Tr{\mat{X}} where \mat{X}~=~(\vect{X}_1 | \dots | \vect{X}_K) and \tilde{K}=K if the process is centered and \tilde{K}=K-1 otherwise;
  • taking the L vertices of the mesh of \vect{X} as the L quadrature points.

We suppose now that K < dL, and we note \mat{Y} = \sqrt{\mat{W}} \,\mat{X}.

As the matrix \mat{\Theta} = \mat{C} is invertible, the Galerkin and collocation approaches are equivalent and both lead to the following singular value problem for \mat{Y}:

(1)\mat{Y}\,\Tr{\mat{Y}}\,\mat{\Psi}_k  & = \tilde{K} \lambda_k \mat{\Psi}_k

The SVD decomposition of \mat{Y}\in \mathcal{M}_{dL,\tilde{K}}(\Rset) writes:

\mat{Y} = \mat{U}\, \mat{\Sigma} \, \Tr{\mat{V}}

where we have \mat{U} \in \mathcal{M}_{dL,\tilde{K}}(\Rset), \mat{\Sigma}\in \mathcal{M}_{\tilde{K},\tilde{K}}(\Rset), \mat{V} \in \mathcal{M}_{\tilde{K},\tilde{K}}(\Rset) such that :

  • \Tr{\mat{V}}\mat{V} =\mat{V}\Tr{\mat{V}}= \mat{I}_{\tilde{K}},
  • \Tr{\mat{U}}\mat{U} = \mat{I}_{\tilde{K}} ,
  • \mat{\Sigma} = diag(\sigma_1, \dots, \sigma_{\tilde{K}}).

Then the columns of \mat{U} are the eigenvectors of \mat{Y}\Tr{\mat{Y}} associated to the eigenvalues \sigma_k^2.

We deduce the modes and eigenvalues of the Fredholm problem for k \leq \tilde{K}:

\begin{align*}
  \mat{\Phi}_k = \dfrac{1}{\lambda_k} \sqrt{\mat{W}}\, \mat{U}_k
  \lambda_k = \dfrac{\sigma_k^2}{\tilde{K}}
\end{align*}

We have:

\tilde{\vect{\varphi}}_k(\vect{t})= \sum_{\ell=1}^L  C(\vect{t}, \vect{\xi}_\ell) \vect{\phi}_{\ell,k} \quad \mbox{pour} \quad k \leq \tilde{K}

Examples

Create a Karhunen-Loeve SVD algorithm:

>>> import openturns as ot
>>> mesh = ot.IntervalMesher([10]*2).build(ot.Interval([-1.0]*2, [1.0]*2))
>>> s = 0.01
>>> model = ot.AbsoluteExponential([1.0]*2)
>>> sample = ot.GaussianProcess(model, mesh).getSample(8)
>>> algorithm = ot.KarhunenLoeveSVDAlgorithm(sample, s)

Run it!

>>> algorithm.run()
>>> result = algorithm.getResult()

Methods

getClassName() Accessor to the object’s name.
getCovarianceModel() Accessor to the covariance model.
getId() Accessor to the object’s id.
getName() Accessor to the object’s name.
getResult() Get the result structure.
getSample() Accessor to the process sample.
getSampleWeights() Accessor to the weights of the sample.
getShadowedId() Accessor to the object’s shadowed id.
getThreshold() Accessor to the threshold used to select the most significant eigenmodes.
getVerticesWeights() Accessor to the weights of the vertices.
getVisibility() Accessor to the object’s visibility state.
hasName() Test if the object is named.
hasVisibleName() Test if the object has a distinguishable name.
run() Computation of the eigenvalues and eigenfunctions values at nodes.
setCovarianceModel(covariance) Accessor to the covariance model.
setName(name) Accessor to the object’s name.
setShadowedId(id) Accessor to the object’s shadowed id.
setThreshold(threshold) Accessor to the limit ratio on eigenvalues.
setVisibility(visible) Accessor to the object’s visibility state.
__init__(*args)

x.__init__(…) initializes x; see help(type(x)) for signature

getClassName()

Accessor to the object’s name.

Returns:

class_name : str

The object class name (object.__class__.__name__).

getCovarianceModel()

Accessor to the covariance model.

Returns:

covModel : CovarianceModel

The covariance model.

getId()

Accessor to the object’s id.

Returns:

id : int

Internal unique identifier.

getName()

Accessor to the object’s name.

Returns:

name : str

The name of the object.

getResult()

Get the result structure.

Returns:

resKL : KarhunenLoeveResult

The structure containing all the results of the Fredholm problem.

Notes

The structure contains all the results of the Fredholm problem.

getSample()

Accessor to the process sample.

Returns:

sample : ProcessSample

The process sample containing the observations of the process.

getSampleWeights()

Accessor to the weights of the sample.

Returns:

weights : Point

The weights associated to the fields of the sample.

Notes

The fields might not have the same weight, for example if they come from importance sampling.

getShadowedId()

Accessor to the object’s shadowed id.

Returns:

id : int

Internal unique identifier.

getThreshold()

Accessor to the threshold used to select the most significant eigenmodes.

Returns:

s : float, positive

The threshold s.

Notes

OpenTURNS truncates the sequence (\lambda_k,  \vect{\varphi}_k)_{k \geq 1} at the index K defined in (3).

getVerticesWeights()

Accessor to the weights of the vertices.

Returns:

weights : Point

The weights associated to the vertices of the mesh defining the sample field.

Notes

The vertices might not have the same weight, for example if the mesh is not regular.

getVisibility()

Accessor to the object’s visibility state.

Returns:

visible : bool

Visibility flag.

hasName()

Test if the object is named.

Returns:

hasName : bool

True if the name is not empty.

hasVisibleName()

Test if the object has a distinguishable name.

Returns:

hasVisibleName : bool

True if the name is not empty and not the default one.

run()

Computation of the eigenvalues and eigenfunctions values at nodes.

Notes

Runs the algorithm and creates the result structure KarhunenLoeveResult.

setCovarianceModel(covariance)

Accessor to the covariance model.

Parameters:

covModel : CovarianceModel

The covariance model.

setName(name)

Accessor to the object’s name.

Parameters:

name : str

The name of the object.

setShadowedId(id)

Accessor to the object’s shadowed id.

Parameters:

id : int

Internal unique identifier.

setThreshold(threshold)

Accessor to the limit ratio on eigenvalues.

Parameters:

s : float, \geq 0

The threshold s defined in (3).

setVisibility(visible)

Accessor to the object’s visibility state.

Parameters:

visible : bool

Visibility flag.