gstools.covmodel.TPLStable

class gstools.covmodel.TPLStable(dim=3, var=1.0, len_scale=1.0, nugget=0.0, anis=1.0, angles=0.0, integral_scale=None, var_raw=None, hankel_kw=None, **opt_arg)[source]

Bases: gstools.covmodel.base.CovModel

Truncated-Power-Law with Stable modes.

Notes

The truncated power law is given by a superposition of scale-dependent variograms:

\gamma_{\ell_{\mathrm{low}},\ell_{\mathrm{up}}}(r) =
\intop_{\ell_{\mathrm{low}}}^{\ell_{\mathrm{up}}}
\gamma(r,\lambda) \frac{\rm d \lambda}{\lambda}

with Stable modes on each scale:

\gamma(r,\lambda) &=
\sigma^2(\lambda)\cdot\left(1-
\exp\left[- \left(\frac{r}{\lambda}\right)^{\alpha}\right]
\right)\\
\sigma^2(\lambda) &= C\cdot\lambda^{2H}

This results in:

\gamma_{\ell_{\mathrm{low}},\ell_{\mathrm{up}}}(r) &=
\sigma^2_{\ell_{\mathrm{low}},\ell_{\mathrm{up}}}\cdot\left(1-
\frac{2H}{\alpha} \cdot
\frac{\ell_{\mathrm{up}}^{2H} \cdot
E_{1+\frac{2H}{\alpha}}
\left[\left(\frac{r}{\ell_{\mathrm{up}}}\right)^{\alpha}\right]
- \ell_{\mathrm{low}}^{2H} \cdot
E_{1+\frac{2H}{\alpha}}
\left[\left(\frac{r}{\ell_{\mathrm{low}}}\right)^{\alpha}\right]}
{\ell_{\mathrm{up}}^{2H}-\ell_{\mathrm{low}}^{2H}}
\right) \\
\sigma^2_{\ell_{\mathrm{low}},\ell_{\mathrm{up}}} &=
\frac{C\cdot\left(\ell_{\mathrm{up}}^{2H}
-\ell_{\mathrm{low}}^{2H}\right)}{2H}

The “length scale” of this model is equivalent by the integration range:

\ell = \ell_{\mathrm{up}} -\ell_{\mathrm{low}}

If you want to define an upper scale truncation, you should set len_low and len_scale accordingly.

The following Parameters occure:

  • 0<\alpha\leq 2 : The shape parameter of the Stable model.

    • \alpha=1 : Exponential modes
    • \alpha=2 : Gaussian modes
  • C>0 : scaling factor from the Power-Law This parameter will be calculated internally by the given variance. You can access C directly by model.var_raw

  • 0<H<\frac{\alpha}{2} : hurst coefficient (model.hurst)

  • \ell_{\mathrm{low}}\geq 0 : lower length scale truncation of the model (model.len_low)

  • \ell_{\mathrm{up}}\geq 0 : upper length scale truncation of the model (model.len_up)

    This will be calculated internally by:

    • len_up = len_low + len_scale

    That means, that the len_scale in this model actually represents the integration range for the truncated power law.

  • E_s(x) is the exponential integral.

Other Parameters:
 
  • **opt_arg – The following parameters are covered by these keyword arguments
  • hurst (float, optional) – Hurst coefficient of the power law. Standard range: (0, 1). Default: 0.5
  • alpha (float, optional) – Shape parameter of the stable model. Standard range: (0, 2]. Default: 1.5
  • len_low (float, optional) – The lower length scale truncation of the model. Standard range: [0, 1000]. Default: 0.0
  • .
Parameters:
  • dim (int, optional) – dimension of the model. Default: 3
  • var (float, optional) – variance of the model (the nugget is not included in “this” variance) Default: 1.0
  • len_scale (float or list, optional) – length scale of the model. If a single value is given, the same length-scale will be used for every direction. If multiple values (for main and transversal directions) are given, anis will be recalculated accordingly. If only two values are given in 3D, the latter one will be used for both transversal directions. Default: 1.0
  • nugget (float, optional) – nugget of the model. Default: 0.0
  • anis (float or list, optional) –

    anisotropy ratios in the transversal directions [e_y, e_z].

    • e_y = l_y / l_x
    • e_z = l_z / l_x

    If only one value is given in 3D, e_y will be set to 1. This value will be ignored, if multiple len_scales are given. Default: 1.0

  • angles (float or list, optional) –

    angles of rotation (given in rad):

    • in 2D: given as rotation around z-axis
    • in 3D: given by yaw, pitch, and roll (known as Tait–Bryan angles)

    Default: 0.0

  • integral_scale (float or list or None, optional) – If given, len_scale will be ignored and recalculated, so that the integral scale of the model matches the given one. Default: None
  • var_raw (float or None, optional) – raw variance of the model which will be multiplied with CovModel.var_factor to result in the actual variance. If given, var will be ignored. (This is just for models that override CovModel.var_factor) Default: None
  • hankel_kw (dict or None, optional) – Modify the init-arguments of hankel.SymmetricFourierTransform used for the spectrum calculation. Use with caution (Better: Don’t!). None is equivalent to {"a": -1, "b": 1, "N": 1000, "h": 0.001}.
Attributes:
angles

numpy.ndarray: Rotation angles (in rad) of the model.

anis

numpy.ndarray: The anisotropy factors of the model.

arg

list of str: Names of all arguments.

arg_bounds

dict: Bounds for all parameters.

dim

int: The dimension of the model.

dist_func

tuple of callable: pdf, cdf and ppf.

do_rotation

bool: State if a rotation is performed.

hankel_kw

dict: hankel.SymmetricFourierTransform kwargs.

has_cdf

bool: State if a cdf is defined by the user.

has_ppf

bool: State if a ppf is defined by the user.

integral_scale

float: The main integral scale of the model.

integral_scale_vec

numpy.ndarray: The integral scales in each direction.

is_isotropic

bool: State if a model is isotropic.

len_scale

float: The main length scale of the model.

len_scale_bounds

list: Bounds for the lenght scale.

len_scale_vec

numpy.ndarray: The length scales in each direction.

len_up

float: Upper length scale truncation of the model.

name

str: The name of the CovModel class.

nugget

float: The nugget of the model.

nugget_bounds

list: Bounds for the nugget.

opt_arg

list of str: Names of the optional arguments.

opt_arg_bounds

dict: Bounds for the optional arguments.

pykrige_angle

2D rotation angle for pykrige.

pykrige_angle_x

3D rotation angle around x for pykrige.

pykrige_angle_y

3D rotation angle around y for pykrige.

pykrige_angle_z

3D rotation angle around z for pykrige.

pykrige_anis

2D anisotropy ratio for pykrige.

pykrige_anis_y

3D anisotropy ratio in y direction for pykrige.

pykrige_anis_z

3D anisotropy ratio in z direction for pykrige.

pykrige_kwargs

Keyword arguments for pykrige routines.

sill

float: The sill of the variogram.

var

float: The variance of the model.

var_bounds

list: Bounds for the variance.

var_raw

float: The raw variance of the model without factor.

Methods

calc_integral_scale(self) Calculate the integral scale of the isotrope model.
check_arg_bounds(self) Check arguments to be within the given bounds.
check_opt_arg(self) Check the optional arguments.
cor(self, h) Normalziled correlation function taking a normalized range.
cor_spatial(self, pos) Spatial correlation respecting anisotropy and rotation.
correlation(self, r) Truncated-Power-Law with Stable modes - correlation function.
cov_nugget(self, r) Covariance of the model respecting the nugget at r=0.
cov_spatial(self, pos) Spatial covariance respecting anisotropy and rotation.
covariance(self, r) Covariance of the model.
default_arg_bounds(self) Provide default boundaries for arguments.
default_opt_arg(self) Defaults for the optional arguments.
default_opt_arg_bounds(self) Defaults for boundaries of the optional arguments.
fit_variogram(self, x_data, y_data[, maxfev]) Fiting the isotropic variogram-model to given data.
fix_dim(self) Set a fix dimension for the model.
ln_spectral_rad_pdf(self, r) Log radial spectral density of the model.
percentile_scale(self[, per]) Calculate the percentile scale of the isotrope model.
plot(self[, func]) Plot a function of a the CovModel.
pykrige_vario(self[, args, r]) Isotropic variogram of the model for pykrige.
set_arg_bounds(self, \*\*kwargs) Set bounds for the parameters of the model.
spectral_density(self, k) Spectral density of the covariance model.
spectral_rad_pdf(self, r) Radial spectral density of the model.
spectrum(self, k) Spectrum of the covariance model.
var_factor(self) Factor for C (Power-Law factor) to result in variance.
vario_nugget(self, r) Isotropic variogram of the model respecting the nugget at r=0.
vario_spatial(self, pos) Spatial variogram respecting anisotropy and rotation.
variogram(self, r) Isotropic variogram of the model.
calc_integral_scale(self)

Calculate the integral scale of the isotrope model.

check_arg_bounds(self)

Check arguments to be within the given bounds.

check_opt_arg(self)[source]

Check the optional arguments.

Warns:alpha – If alpha is < 0.3, the model tends to a nugget model and gets numerically unstable.
cor(self, h)

Normalziled correlation function taking a normalized range.

Given by: \mathrm{cor}\left(r/\ell\right) = \rho(r)

cor_spatial(self, pos)

Spatial correlation respecting anisotropy and rotation.

correlation(self, r)[source]

Truncated-Power-Law with Stable modes - correlation function.

If len_low=0 we have a simple representation:

\rho(r) =
\frac{2H}{\alpha} \cdot
E_{1+\frac{2H}{\alpha}}
\left[
\left(\frac{r}{\ell}\right)^{\alpha}
\right]

The general case:

\rho(r) =
\frac{2H}{\alpha} \cdot
\frac{\ell_{\mathrm{up}}^{2H} \cdot
E_{1+\frac{2H}{\alpha}}
\left[\left(\frac{r}{\ell_{\mathrm{up}}}\right)^{\alpha}\right]
- \ell_{\mathrm{low}}^{2H} \cdot
E_{1+\frac{2H}{\alpha}}
\left[\left(\frac{r}{\ell_{\mathrm{low}}}\right)^{\alpha}\right]}
{\ell_{\mathrm{up}}^{2H}-\ell_{\mathrm{low}}^{2H}}

cov_nugget(self, r)

Covariance of the model respecting the nugget at r=0.

Given by: C\left(r\right)=
\sigma^2\cdot\rho\left(r\right)

Where \rho(r) is the correlation function.

cov_spatial(self, pos)

Spatial covariance respecting anisotropy and rotation.

covariance(self, r)

Covariance of the model.

Given by: C\left(r\right)=
\sigma^2\cdot\rho\left(r\right)

Where \rho(r) is the correlation function.

default_arg_bounds(self)

Provide default boundaries for arguments.

Given as a dictionary.

default_opt_arg(self)[source]

Defaults for the optional arguments.

  • {"hurst": 0.5, "alpha": 1.5, "len_low": 0.0}
Returns:Defaults for optional arguments
Return type:dict
default_opt_arg_bounds(self)[source]

Defaults for boundaries of the optional arguments.

  • {"hurst": [0, 1, "oo"], "alpha": [0, 2, "oc"], "len_low": [0, 1000, "cc"]}
Returns:Boundaries for optional arguments
Return type:dict
fit_variogram(self, x_data, y_data, maxfev=1000, **para_deselect)

Fiting the isotropic variogram-model to given data.

Parameters:
  • x_data (numpy.ndarray) – The radii of the meassured variogram.
  • y_data (numpy.ndarray) – The messured variogram
  • maxfev (int, optional) – The maximum number of calls to the function in scipy curvefit. Default: 1000
  • **para_deselect – You can deselect the parameters to be fitted, by setting them “False” as keywords. By default, all parameters are fitted.
Returns:

Notes

You can set the bounds for each parameter by accessing CovModel.set_arg_bounds.

The fitted parameters will be instantly set in the model.

fix_dim(self)

Set a fix dimension for the model.

ln_spectral_rad_pdf(self, r)

Log radial spectral density of the model.

percentile_scale(self, per=0.9)

Calculate the percentile scale of the isotrope model.

This is the distance, where the given percentile of the variance is reached by the variogram

plot(self, func='variogram', **kwargs)

Plot a function of a the CovModel.

Parameters:
  • func (str, optional) –

    Function to be plotted. Could be one of:

    • ”variogram”
    • ”covariance”
    • ”correlation”
    • ”vario_spatial”
    • ”cov_spatial”
    • ”cor_spatial”
    • ”spectrum”
    • ”spectral_density”
    • ”spectral_rad_pdf”
  • **kwargs – Keyword arguments forwarded to the plotting function “plot_” + func in gstools.covmodel.plot.
pykrige_vario(self, args=None, r=0)

Isotropic variogram of the model for pykrige.

Given by: \gamma\left(r\right)=
\sigma^2\cdot\left(1-\rho\left(r\right)\right)+n

Where \rho(r) is the correlation function.

set_arg_bounds(self, **kwargs)

Set bounds for the parameters of the model.

Parameters:**kwargs

Parameter name as keyword (“var”, “len_scale”, “nugget”, <opt_arg>) and a list of 2 or 3 values as value:

  • [a, b] or
  • [a, b, <type>]

<type> is one of "oo", "cc", "oc" or "co" to define if the bounds are open (“o”) or closed (“c”).

spectral_density(self, k)

Spectral density of the covariance model.

This is given by:

\tilde{S}(k) = \frac{S(k)}{\sigma^2}

Where S(k) is the spectrum of the covariance model.

Parameters:k (float) – Radius of the phase: k=\left\Vert\mathbf{k}\right\Vert
spectral_rad_pdf(self, r)

Radial spectral density of the model.

spectrum(self, k)

Spectrum of the covariance model.

This is given by:

S(k) = \left(\frac{1}{2\pi}\right)^n
\int C(r) e^{i b\mathbf{k}\cdot\mathbf{r}} d^n\mathbf{r}

Internally, this is calculated by the hankel transformation:

S(k) = \left(\frac{1}{2\pi}\right)^n \cdot
\frac{(2\pi)^{n/2}}{k^{n/2-1}}
\int_0^\infty r^{n/2} C(r) J_{n/2-1}(kr) dr

Where C(r) is the covariance function of the model.

Parameters:k (float) – Radius of the phase: k=\left\Vert\mathbf{k}\right\Vert
var_factor(self)[source]

Factor for C (Power-Law factor) to result in variance.

This is used to result in the right variance, which is depending on the hurst coefficient and the length-scale extents

\frac{\ell_{\mathrm{up}}^{2H} - \ell_{\mathrm{low}}^{2H}}{2H}

Returns:factor
Return type:float
vario_nugget(self, r)

Isotropic variogram of the model respecting the nugget at r=0.

Given by: \gamma\left(r\right)=
\sigma^2\cdot\left(1-\rho\left(r\right)\right)+n

Where \rho(r) is the correlation function.

vario_spatial(self, pos)

Spatial variogram respecting anisotropy and rotation.

variogram(self, r)

Isotropic variogram of the model.

Given by: \gamma\left(r\right)=
\sigma^2\cdot\left(1-\rho\left(r\right)\right)+n

Where \rho(r) is the correlation function.

angles

Rotation angles (in rad) of the model.

Type:numpy.ndarray
anis

The anisotropy factors of the model.

Type:numpy.ndarray
arg

Names of all arguments.

Type:list of str
arg_bounds

Bounds for all parameters.

Notes

Keys are the opt-arg names and values are lists of 2 or 3 values:

  • [a, b] or
  • [a, b, <type>]

<type> is one of "oo", "cc", "oc" or "co" to define if the bounds are open (“o”) or closed (“c”).

Type:dict
dim

The dimension of the model.

Type:int
dist_func

pdf, cdf and ppf.

Spectral distribution info from the model.

Type:tuple of callable
do_rotation

State if a rotation is performed.

Type:bool
hankel_kw

hankel.SymmetricFourierTransform kwargs.

Type:dict
has_cdf

State if a cdf is defined by the user.

Type:bool
has_ppf

State if a ppf is defined by the user.

Type:bool
integral_scale

The main integral scale of the model.

Raises:ValueError – If integral scale is not setable.
Type:float
integral_scale_vec

The integral scales in each direction.

Notes

This is calculated by:
  • integral_scale_vec[0] = integral_scale
  • integral_scale_vec[1] = integral_scale*anis[0]
  • integral_scale_vec[2] = integral_scale*anis[1]
Type:numpy.ndarray
is_isotropic

State if a model is isotropic.

Type:bool
len_scale

The main length scale of the model.

Type:float
len_scale_bounds

Bounds for the lenght scale.

Notes

Is a list of 2 or 3 values:

  • [a, b] or
  • [a, b, <type>]

<type> is one of "oo", "cc", "oc" or "co" to define if the bounds are open (“o”) or closed (“c”).

Type:list
len_scale_vec

The length scales in each direction.

Notes

This is calculated by:
  • len_scale_vec[0] = len_scale
  • len_scale_vec[1] = len_scale*anis[0]
  • len_scale_vec[2] = len_scale*anis[1]
Type:numpy.ndarray
len_up

Upper length scale truncation of the model.

  • len_up = len_low + len_scale
Type:float
name

The name of the CovModel class.

Type:str
nugget

The nugget of the model.

Type:float
nugget_bounds

Bounds for the nugget.

Notes

Is a list of 2 or 3 values:

  • [a, b] or
  • [a, b, <type>]

<type> is one of "oo", "cc", "oc" or "co" to define if the bounds are open (“o”) or closed (“c”).

Type:list
opt_arg

Names of the optional arguments.

Type:list of str
opt_arg_bounds

Bounds for the optional arguments.

Notes

Keys are the opt-arg names and values are lists of 2 or 3 values:

  • [a, b] or
  • [a, b, <type>]

<type> is one of "oo", "cc", "oc" or "co" to define if the bounds are open (“o”) or closed (“c”).

Type:dict
pykrige_angle

2D rotation angle for pykrige.

pykrige_angle_x

3D rotation angle around x for pykrige.

pykrige_angle_y

3D rotation angle around y for pykrige.

pykrige_angle_z

3D rotation angle around z for pykrige.

pykrige_anis

2D anisotropy ratio for pykrige.

pykrige_anis_y

3D anisotropy ratio in y direction for pykrige.

pykrige_anis_z

3D anisotropy ratio in z direction for pykrige.

pykrige_kwargs

Keyword arguments for pykrige routines.

sill

The sill of the variogram.

Notes

This is calculated by:
  • sill = variance + nugget
Type:float
var

The variance of the model.

Type:float
var_bounds

Bounds for the variance.

Notes

Is a list of 2 or 3 values:

  • [a, b] or
  • [a, b, <type>]

<type> is one of "oo", "cc", "oc" or "co" to define if the bounds are open (“o”) or closed (“c”).

Type:list
var_raw

The raw variance of the model without factor.

(See. CovModel.var_factor)

Type:float