krigekit.kriging_indicator#

Python wrapper for Multiple Indicator Kriging (MIK) and Sequential Indicator Simulation (SIS) via the Fortran t_kriging_indicator type.

Each of the K categories (or threshold classes) is treated as one indicator variable (ivar = 1..K). Estimation produces K probability values per block; simulation produces a single drawn category (encoded as a one-hot binary vector).

Indicator encoding and the K x K coregionalization are built with IndicatorVariogramSystem and transferred with its apply(); this wrapper only allocates the engine, solves, and post-processes probabilities.

Typical workflow — estimation#

>>> system = IndicatorVariogramSystem(categories=[1, 2, 3])
>>> system.set_categorical_obs(obs_coord, obs_cat)
>>> system.set_indicator_vgm(vtype="sph", a_major=1000,
...                          sill_strategy="theoretical", cross_strategy="closure")
>>> ik = IndicatorKriging(ncat=3, ndim=2)
>>> system.apply(ik)                # transfers indicators + K x K structures
>>> ik.set_grid(coord=grid_coord)
>>> for k in range(1, 4): ik.set_search(ivar=k, nmax=20)
>>> ik.solve()
>>> probs, var = ik.get_results()   # probs.shape == (ngrid, ncat)
>>> del ik

Typical workflow — SIS#

>>> system = IndicatorVariogramSystem(categories=[1, 2, 3])
>>> system.set_categorical_obs(obs_coord, obs_cat)
>>> system.set_indicator_vgm(vtype="sph", sill=0.2, a_major=500,
...                          sill_strategy="uniform", cross_strategy="uniform")
>>> ik = IndicatorKriging(ncat=3, ndim=2, nsim=100)
>>> system.apply(ik)
>>> ik.set_grid(coord=grid_coord)
>>> ik.set_sim()    # generates U(0,1) draws in Python and passes to Fortran
>>> for k in range(1, 4): ik.set_search(ivar=k, nmax=20)
>>> ik.solve()
>>> sims, _ = ik.get_results()     # sims.shape == (ngrid, ncat, nsim)
...                                 # each [:, :, i] is a one-hot encoded realisation

Classes#

IndicatorKriging

Multiple Indicator Kriging and Sequential Indicator Simulation,

Module Contents#

class krigekit.kriging_indicator.IndicatorKriging(ncat: int, nvar: int | None = None, ndim: int = 2, ndrift: int = 0, unbias: int = 1, nsim: int = 0, anisotropic_search: bool = False, weight_correction: bool = False, use_old_weight: bool = False, store_weight: bool = False, cross_validation: bool = False, write_mat: bool = False, neglect_error: bool = True, varying_vgm: bool = False, std_ck: bool = False, verbose: bool = False, pf_cache: bool = False, weight_file: str = '', bounds: tuple | None = None, seed: int | None = None)#

Bases: krigekit.kriging.Kriging

Multiple Indicator Kriging and Sequential Indicator Simulation, with optional secondary co-variate support.

Extends Kriging — all setup, solve, and results methods are inherited unchanged. The differences are:

The Fortran object is a t_kriging_indicator (created via krige_ind_create), which overrides prepare, sim_draw, and post_solve to implement indicator-specific behaviour.
ncat names the K indicator categories. nvar defaults to ncat (pure MIS) but can be set larger to add secondary continuous co-variates for co-kriging MIS (see below).
Indicator encoding and the K x K coregionalization are built with IndicatorVariogramSystem and transferred via its apply(); this engine wrapper no longer owns that construction.

Parameters:

ncat (int) – Number of categories K. Indicator variables occupy ivar = 1..ncat.
nvar (int, optional) – Total number of co-kriging variables. Defaults to ncat (pure MIS). Set nvar = ncat + M to include M secondary continuous variables (ivar = ncat+1 .. nvar). Secondary variables contribute to the kriging weights but are excluded from the CDF draw and probability normalisation.
ndim (int) – Number of spatial dimensions (2 or 3).
nsim (int) – 0 = estimation (returns probabilities); >0 = SIS (returns one-hot draws).
**kwargs – All other keyword arguments are passed through to Kriging.

Notes

For SIS (nsim > 0), call set_sim() after set_grid().

Co-kriging MIS example (K=3 categories + 1 secondary variable):

system = IndicatorVariogramSystem(categories=[1, 2, 3])
system.set_categorical_obs(coord, cats)
system.set_indicator_vgm(vtype="sph", a_major=1000,
                         sill_strategy="theoretical",
                         cross_strategy="closure")
ik = IndicatorKriging(ncat=3, nvar=4, ndim=2)
system.apply(ik)                            # indicator block (ivar 1..3)
ik.set_obs(ivar=4, coord=sec_coord, value=sec_val)   # secondary
ik.set_vgm(ivar=4, jvar=4, ...)             # secondary auto/cross models
ik.set_grid(coord=grid_coord)
for k in range(1, 5):
    ik.set_search(ivar=k, nmax=20)
ik.solve()
probs, var = ik.get_results()   # shape (ngrid, 3) — secondary excluded

Parameters:

ndim (int) – Number of spatial dimensions (2 or 3).
nvar (int) – Number of variables (1 for ordinary/simple kriging, >1 for cokriging).
ndrift (int) – Number of external drift functions (0 = no drift).
unbias (int) – 1 = ordinary kriging (sum-of-weights = 1 constraint); 0 = simple kriging (no constraint, uses sk_mean).
nsim (int) – Number of simulations. 0 = kriging only; >0 = SGSIM.
anisotropic_search (bool) – Use anisotropic search ellipse for neighbour lookup.
weight_correction (bool) – Force kriging weights to be non-negative and sum to 1.
use_old_weight (bool) – Read pre-computed weights from weight_file instead of solving.
store_weight (bool) – Write computed weights to weight_file while also estimating blocks.
cross_validation (bool) – Leave-one-out cross-validation mode.
write_mat (bool) – Write matrix for debugging.
neglect_error (bool) – Ignore solver errors and set failed block to NaN instead of aborting.
varying_vgm (bool) – Use a different variogram per estimation block (spatially varying anisotropy). When True, call set_vgm_block() for each block after set_grid(). Defaults to False (single global model).
std_ck (bool) –
Co-kriging unbiasedness formulation (only relevant when nvar > 1 and unbias=1).
- False (default) — Isaaks & Srivastava: single combined constraint (Σw₁ + Σw₂ = 1) plus a local-mean correction applied post-solve. Matches the GSLIB/legacy behaviour.
- True — standard cokriging: separate per-variable constraints (Σwᵢ = 1 for own variable, Σwⱼ = 0 for others), equivalent to the gstat/ISATIS formulation. Use this to match gstat output.
verbose (bool) – Print progress messages.
pf_cache (bool) – Enable the persistent between-solve factorization cache. When True, the Cholesky factor of K is stored after the first solve() and reused on subsequent calls when the neighbour set and variogram are unchanged (speeds up repeated solves on the same observation grid). Defaults to False; enable only when you plan to call solve() multiple times and need the speedup.
weight_file (str) – Path to the weight file (required when use_old_weight or store_weight).
bounds (tuple(float, float) or None) – (lower, upper) clipping bounds for the estimate. None means no clipping (uses Fortran defaults: [-huge, +huge]).
seed (int, optional) – Random seed.

get_results(copy: bool = False, squeeze: bool = True)#

Return indicator results, excluding secondary covariate channels.

The kriging engine stores all nvar estimates internally because secondary variables participate in the cokriging system. Public indicator results contain only the first ncat variables, matching the probability or one-hot category array documented by this class.

set_sim(randpath: numpy.ndarray | None = None, sample: numpy.ndarray | None = None)#

Set up Sequential Indicator Simulation parameters.

Overrides set_sim(). When sample is None, delegates to the parent with no sample so that the Fortran set_sim_indicator override generates U(0, 1) draws directly; no sample array is created in Python. When sample is supplied, validates that every value lies in [0, 1] then delegates to the parent.

Parameters:

randpath (ndarray of int, shape (nblocks,), optional) – Random visiting order (1-based). Generated with a random permutation if omitted.
sample (ndarray, shape (nblocks, nvar, nsim), optional) – Pre-drawn U(0, 1) samples. Every value must lie in [0, 1]. When None, Fortran generates U(0, 1) via set_sim_indicator.