SFI.inference.overdamped module

class SFI.inference.overdamped.OverdampedLangevinInference(data, *, max_memory_gb=1.0, **kwargs)[source]

Bases: BaseLangevinInference

Stochastic Force Inference concrete class for overdamped systems

This class provides tools for inferring force (drift) and diffusion tensors from stochastic trajectory data based on overdamped Langevin dynamics. It supports both linear and nonlinear basis function methods.

Notes

The dynamics are described by the 1st order autonomous stochastic differential equation (SDE):

dx/dt = F(x) + sqrt(2D(x)) dxi(t)

where:

  • F(x) is the Ito drift (force) term.

  • D(x) is the diffusion tensor, evaluated in the Ito convention.

  • dxi(t) is Gaussian white noise.

[Dynamical equations] Overdamped Langevin SDE

\[\frac{\mathrm{d}x}{\mathrm{d}t} = F(x) + \sqrt{2\,D(x)}\;\mathrm{d}\xi(t)\]

\(F(x)\) is the Itô drift, \(D(x)\) the diffusion tensor (Itô convention), and \(\mathrm{d}\xi\) is Gaussian white noise.

Here x is a 2D array of shape Nparticles x dimension. All particles are assumed to have identical properties.

This class provides tools to approximate F(x) and D(x) from a time series x(t) formatted as TrajectoryCollection.

Note that the Ito and Strato variants of the force inference routines do NOT refer to the convention in which the SDE is expressed (which is always Ito), but to the way stochastic integrals are performed to compute parameters.

Notes

  • Force Inference: - Linear combination of basis functions (infer_force_linear). - Parametric families (infer_force with a Basis or PSF).

  • Diffusion Inference: - Constant diffusion estimation (compute_diffusion_constant). - State-dependent diffusion with basis functions (infer_diffusion_linear and infer_diffusion).

  • Sparsification: - Force sparsification for linear inference sparsify_force, implementing PASTIS and other information criteria.

  • Error Estimation: - Normalized mean-squared error (MSE) prediction for both force and diffusion.

  • Comparison Tools for method benchmarking: - Evaluate inferred fields against known exact models (compare_to_exact).

  • Simulation: - Generate trajectories using inferred fields (simulate_bootstrapped_trajectory).

Notes

  1. Initialize with TrajectoryCollection containing the trajectory.

  2. Use the infer_* methods to infer force and diffusion fields.

  3. Optionally compute error estimates and/or compare with exact data for validation.

Notes

The code uses jnp.einsum for array manipulation, with a consistent index naming scheme for clarity: - t : Time index, = 0..Ntimesteps-1 - a, b, c… : Basis function indices, = 0..Nfunctions - 1. - m, n, o… : Spatial indices, = 0..dim-1. - i, j… : Particle indices. We also use these indices as shortcuts for array shapes. For instance basis_linear : im -> iam reads basis_linear has input a jnp.array of shape (Nparticles,dim) and outputs a jnp.array of shape (Nparticles,Nfunctions,dim).

Notes

Progress messages use Python logging. Enable with SFI.enable_logging() or logging.getLogger('SFI').setLevel(logging.INFO).

Example

Fully documented examples in the “examples” folder: Lorenz model, ActiveBrownianParticles, Ornstein-Uhlenbeck…

coeff_block(block, *, field='force')

Return the coefficient (and covariance) slice for a basis sub-block.

Compound bases (e.g. a multi-kernel or time-Fourier library) pack several conceptual blocks into one flat coefficient vector. This returns the slice for one block without hand-computed offsets.

Parameters:

  • block – (start, stop) indices, a slice, an int, or a Basis (located by matching its labels as a contiguous run of the fitted basis labels).

  • field ({"force", "diffusion"})

Returns:

(coeffs, cov) – coeffs is the 1-D slice; cov is the matching covariance block (or None if no covariance is available).

Return type:

tuple

compare_params_to_exact(theta_true, *, psf=None)

Compare inferred parametric coefficients to known ground truth.

For a model fitted with a parametric family, returns a per-parameter dict of absolute and relative error. theta_true may be a flat array (compared elementwise to force_coefficients_full) or a {name: value} dict (unflattened from the fitted coefficients via psf.unflatten_params, falling back to self.force_psf).

Returns:

{name: {"true", "inferred", "abs_error", "rel_error"}}; also stored as self.parameter_comparison.

Return type:

dict

compare_to_exact(*, model_exact=None, data_exact=None, force_exact=None, diffusion_exact=None, maxpoints=1000)

Compare inferred vs exact using dt-weighted time means via the integrate() API.

This function evaluates the inferred force/diffusion against “exact” references on a (possibly exact/synthetic) dataset. It updates:

self.MSE_force / self.NMSE_force self.MSE_diffusion / self.NMSE_diffusion

Inputs: exact references

You can provide exact references in two ways:

  1. Preferred: model_exact A model object (from SFI.langevin submodule) exposing:

    • model_exact.force_sf : exact force/drift (SF/StateExpr-like)

    • model_exact.diffusion_sfexact diffusion (SF/StateExpr-like)

      OR a constant (float or (d,d) matrix) via model_exact.D

  2. Explicit: force_exact, diffusion_exact - force_exact: SF/StateExpr-like callable returning (N, d) - diffusion_exact:

    • callable returning (N, d, d), OR

    • float meaning σ·I, OR

    • (d,d) matrix constant diffusion.

    These are used if model_exact is not provided. If model_exact is provided, its members take precedence unless they are missing, in which case the explicit arguments can be used as fallback.

Velocity provisioning (underdamped)

If an evaluated expression advertises needs_v=True, this routine supplies:

v := dX/dt (secant velocity from the data stream)

i.e. it uses velocity(“dX”, “dt”) as the v=… keyword argument. This works for both exact and inferred expressions and keeps underdamped comparisons possible even when the dataset only stores positions.

Metrics

Force:

e = Fe - Fh MSE_force = < e^T A_inv e > NMSE_force = MSE_force / < Fh^T A_inv Fh >

Diffusion:

E = De - Dh MSE_diffusion = < tr(A_inv E A_inv E) >

[Error analysis] Normalized MSE metrics (force & diffusion)

\[\text{NMSE}_F = \frac{\langle (F_{\text{exact}} - \hat F)^\top A^{-1} (F_{\text{exact}} - \hat F) \rangle} {\langle \hat F^\top A^{-1} \hat F \rangle}\]
\[\text{NMSE}_D = \frac{\langle \operatorname{tr}(A^{-1} E\, A^{-1} E) \rangle} {\langle \operatorname{tr}(A^{-1} \hat D\, A^{-1} \hat D) \rangle}\]
where \(E = D_{\text{exact}} - \hat D\).

NMSE_diffusion = MSE_diffusion / < tr(A_inv Dh A_inv Dh) >

where A_inv is self.A_inv (typically (2 D̄)^{-1} from the inferred constant diffusion normalization).

Subsampling

Uses a simple subsampling heuristic so that the total number of evaluated points is ~<= maxpoints, accounting for the maximum number of particles.

Requirements

  • self.A_inv must exist (run compute_diffusion_constant() or otherwise set A_inv).

  • The dataset must provide streams X, dt, and if any evaluated expr needs v: dX as well.

Parameters:

maxpoints (int)

Return type:

None

comparison_scatter(*, model_exact=None, field='force', data=None, ax=None, maxpoints=2000, **plot_kw)

Scatter inferred-vs-exact force (or diffusion) along the trajectory.

Evaluates both fields on the data with force_comparison_arrays() / diffusion_comparison_arrays() and renders them with SFI.utils.plotting.comparison_scatter() (identity line + Pearson r + MSE). Replaces hand-rolled exact-vs-inferred scatters in demos.

Parameters:

  • model_exact – Object exposing force_sf / diffusion_sf / D.

  • field ({"force", "diffusion"}) – Which field to compare.

  • data – Collection to evaluate on (default: training data).

  • ax – Target axes (default: current axes).

  • maxpoints (int) – Approximate number of points to evaluate.

  • **plot_kw – Forwarded to SFI.utils.plotting.comparison_scatter().

Return type:

matplotlib.axes.Axes

compute_diffusion_constant(method='auto')[source]

Estimate a constant (spatially uniform) diffusion matrix.

Parameters:

method ({"auto", "noisy", "WeakNoise", "MSD"}) – Estimator to use. "noisy" is the noise-robust Vestergaard–Blainey–Flyvbjerg estimator. "auto" selects "noisy" when the measurement-noise trace Tr(Λ) > 0 (localization noise detected), and "WeakNoise" otherwise.

Return type:

None

Notes

Sets diffusion_average, diffusion_inferred, A, A_inv, sqrtA, sqrtA_inv, Lambda, Lambda_trace, and metadata["diffusion_constant_method"].

compute_diffusion_error()

Estimate sampling error for diffusion inference.

Mirrors compute_force_error() for the diffusion field. Uses the diffusion Gram matrix (normal matrix) and its inverse.

For linear diffusion inference the moments covariance is proportional to the Gram matrix, giving Cov(θ_D) = cov_factor * G_D⁻¹.

Note

This error estimate is approximate. The diffusion inference is more complex than the force inference: diffusion coefficients are inferred from force residuals, a positive-definiteness constraint applies, and the simple covariance formula Cov(θ_D) = cov_factor * G_D⁻¹ may not capture all sources of uncertainty. Treat the result as a rough guide rather than a rigorous confidence interval.

Notes

self.diffusion_coefficients_covariancejnp.ndarray

Covariance matrix of the diffusion coefficients.

self.diffusion_coefficients_stderrjnp.ndarray

Standard error for each diffusion coefficient.

self.diffusion_informationfloat

Estimated information content of the inferred diffusion field.

self.diffusion_predicted_MSEfloat

Predicted normalized mean squared error.

compute_force_error()

Estimate sampling error for force inference.

[Error analysis] Force coefficient covariance & predicted error

\[\operatorname{Cov}(C) = G^{-1}, \qquad \mathbb{E}\!\left[\langle \delta F^\top A^{-1} \delta F \rangle\right] = \operatorname{Tr}\!\left(G\,\operatorname{Cov}(C)\right), \qquad I_F = \tfrac{1}{2}\,C^\top M, \qquad \text{NMSE}_{F,\text{pred}} = \frac{\operatorname{Tr}(G\cdot\operatorname{Cov}(C))}{C^\top M} = \frac{\operatorname{Tr}(G\cdot\operatorname{Cov}(C))}{2 I_F}\]

Assumes the sampling error dominates; measurement noise and discretization biases are not addressed.

This method evaluates the covariance of the inferred force coefficients, the standard error, and computes the predicted normalized mean squared error (MSE) of the inferred force field. This analysis assumes that the sampling error dominates, and measurement noise or discretization biases are not explicitly addressed. It is common to OLI and ULI (by construction of the normal matrix G).

Notes

self.force_coefficients_covariance (jnp.ndarray): Covariance matrix of the force coefficients. self.force_coefficients_stderr (jnp.ndarray): Standard error for each force coefficient. self.force_information (float): Estimated information content of the inferred force field. self.force_predicted_MSE (float): Predicted normalized mean squared error of the inferred force field.

diagnose(*, level='standard', **kwargs)

Run the consistency-check suite from SFI.diagnostics.

Convenience wrapper for SFI.diagnostics.assess(). See its docstring for the available level presets.

Parameters:

level (str)

diffusion_comparison_arrays(*, model_exact=None, diffusion_exact=None, data=None, maxpoints=2000)

Return (D_exact, D_inferred) evaluated along the trajectory.

Like force_comparison_arrays() but for the diffusion field; a constant exact/inferred diffusion is broadcast to (n_points, d, d).

Parameters:

maxpoints (int)

force_comparison_arrays(*, model_exact=None, force_exact=None, data=None, maxpoints=2000)

Return (F_exact, F_inferred) evaluated along the trajectory.

Evaluates the exact and inferred force on the (subsampled, masked) trajectory points, supplying finite-difference velocities for underdamped fields. Feeds comparison_scatter(); also handy for custom diagnostics.

Parameters:

  • model_exact – Object exposing force_sf (e.g. an OverdampedProcess).

  • force_exact – Explicit callable exact force (overrides model_exact).

  • data – Collection to evaluate on (default: the training data).

  • maxpoints (int) – Approximate number of points to evaluate.

Returns:

(F_exact, F_inferred)

Return type:

tuple of ndarray, shape (n_points, d)

get_diffusion_timeop(method)[source]

Return the requested overdamped diffusion estimator as a TimeOperand. Output is (N, d, d). Required streams are declared on the wrapped TimeOp.

Parameters:

method (str)

Return type:

TimeOperand

holdout_score(data, *, require_error=False)

Held-out NMSE of the fitted force on an independent collection.

A side feature for data-abundant scenarios: SFI estimates its own accuracy from the training data (force_predicted_MSE) and validates fits through the diagnostics suite, neither of which costs any data. Reach for an explicit train/test split (TrajectoryCollection.split_time) only when data is plentiful, or to confirm a suspected bias floor: a ratio near 1 means the fit is sampling-limited, a ratio ≫ 1 means a bias floor (often measurement noise — see the noise-and-sampling guide).

The score is the residual-based normalised mean-square error of force_inferred on data, with the diffusion noise floor subtracted (a bias detector, not a precision instrument: its resolution is set by the χ² fluctuations of the residuals).

Parameters:

  • data (TrajectoryCollection) – Independent test data (e.g. the second half of coll.split_time(0.8)).

  • require_error (bool) – If True, run compute_force_error() first when the predicted error is missing, so ratio is always defined.

Returns:

{"holdout_NMSE", "predicted_NMSE", "ratio", "excess_z", "n_obs"}. Also stored as self.force_holdout_NMSE.

Return type:

dict

Notes

Bases that read time-dependent extras are not supported on the held-out path (the residual builders pass extras unsliced).

infer_diffusion(basis=None, *, theta_D0=None, integrator='rk4', n_substeps=1, maxiter=100)[source]

Infer state-dependent diffusion D(x) from parametric residuals.

Requires a prior parametric infer_force() call. Holds the fitted force fixed and minimises the windowed conditional NLL over the diffusion parameters (the log-det term makes the diffusion level identifiable), reusing the same flow residuals and integrate engine as the force solve.

[Inference] State-dependent diffusion inference (overdamped)

With the force \(\hat F\) held fixed, the state-dependent diffusion \(D(x;\theta_D)\) is optimised by minimising the windowed conditional negative log-likelihood; \(\Lambda\) from the force inference is held fixed. A rank-2 basis gives \(D(x) = \sum_j (\theta_D)_j\, d_j(x)\); a PSF is evaluated directly.

Parameters:

  • basis (Basis (rank 2), PSF, or None) – Diffusion model. A rank-2 Basis gives the linear parameterisation; a PSF is used directly (D(x) = PSF(x; θ_D)). When None (default), symmetric_matrix_basis(d) — all constant symmetric d×d diffusion matrices.

  • theta_D0 (dict, array, or None) – Initial diffusion parameters (default zeros).

  • integrator ({"rk4", "euler"}) – Flow predictor (default "rk4", matching infer_force()).

  • n_substeps (int) – Integrator micro-steps per Δt (default 1).

  • maxiter (int) – L-BFGS maximum iterations (default 100).

Return type:

None

Notes

Sets diffusion_inferred, diffusion_coefficients (and diffusion_basis when basis is a Basis), plus metadata.

See also

Parametric windowed estimators — algorithm and parameters

Full algorithm description.

infer_diffusion_linear(basis=None, *, M_mode='auto', G_mode='rectangle')[source]

Fit the diffusion field as a linear combination of basis functions.

This method computes the coefficients of the diffusion tensor field (self.diffusion_coefficients) using the provided basis functions. The diffusion tensor is represented as:

diffusion_inferred(x, mask) = sum_a basis_linear(x, mask)[:,a] * diffusion_coefficients[a]

Parameters:

  • basis (Basis with rank = 2 or None) – the fitting functions. When None (default), symmetric_matrix_basis(d) is used, spanning all constant symmetric d×d diffusion matrices. Requires a prior compute_diffusion_constant() call to determine d.

  • M_mode (str) – The method used for local diffusion tensor estimation and moments computation. See _diffusion_estimator documentation for additional information.

  • G_mode (str) – The method used to compute the normalization matrix G. Not investigated extensively yet for diffusion inference.

Return type:

None

Notes

self.diffusion_coefficients: The inferred coefficients for the diffusion basis functions. self.diffusion_inferred: Callable representing the inferred diffusion tensor field. self.diffusion_G: The normalization matrix used in the inference process.

Note

This inferred tensor field is not guaranteed to be nonnegative.

infer_force(F, theta0=None, *, D=None, Lambda=None, integrator='rk4', n_substeps=1, inner='auto', eiv='auto', max_outer=5, inner_maxiter=80, extra_radius=1)[source]

Infer the force field with the minimal parametric estimator.

Built on SFI.inference.parametric_core: a single RK4 flow step per observation interval defines the residual, the residual covariance gives a windowed-precision NLL, and the parameters are found by direct Gauss–Newton (linear-in-θ Basis, with the skip-trick errors-in-variables instrument) or frozen-precision L-BFGS (nonlinear-in-θ PSF). (D, Λ) are profiled natively: moment-estimator init, then one windowed conditional-NLL refinement at the fitted θ.

[Inference] Parametric windowed force inference (overdamped)

The observed positions follow \(y_t = x_t + \eta_t\) where \(\eta \sim \mathcal{N}(0, \Lambda)\). The deterministic flow \(\Phi(x;\theta) = z(\Delta t) - x\) (one RK4 step by default) defines the residual \(r_t = y_{t+1} - y_t - \Phi(y_t;\theta)\). Residuals follow a banded Gaussian whose local precision weights the Gauss–Newton normal equations; under measurement noise the left factor is replaced by the η-clean skip instrument \(\psi_{\rm inst} = \partial\Phi/\partial\theta\) evaluated at the lagged clean point (eiv=True), giving a consistent estimating equation.

Parameters:

  • F (PSF or Basis) – Parametric drift model. A Basis is converted to a PSF internally (coefficients initialised to zero) and PASTIS sparsification is enabled; it runs the fast direct-GN path. A PSF (possibly nonlinear in θ) runs the L-BFGS path.

  • theta0 (dict, array, or None) – Initial drift parameters (default: zeros).

  • D (array (d, d), optional) – Fixed diffusion matrix. If both D and Lambda are given, noise profiling is skipped entirely (fast path).

  • Lambda (array (d, d), optional) – Fixed measurement-noise covariance.

  • integrator ({"rk4", "euler"}) – Flow predictor (default "rk4", a single 4th-order step).

  • n_substeps (int) – Integrator micro-steps per observation interval (default 1 — the single-step minimal estimator).

  • inner ({"auto", "gn", "lbfgs"}) – Inner solver. "auto" → direct Gauss–Newton for a linear Basis, L-BFGS for a PSF.

  • eiv ({"auto", True, False, float}) – Measurement-noise errors-in-variables instrument. "auto" (default) → True for all models (interacting models use the same N-body flow for the instrument as for the residual); True forces the η-clean skip instrument (consistent under noise); False is the plain MLE; a float in [0, 1] blends. Active on the GN path only.

  • max_outer (int) – Outer Gauss–Newton / IRLS iterations (default 5).

  • inner_maxiter (int) – Inner L-BFGS iterations per outer step on the PSF path (default 80; raise for large nonlinear families, e.g. NNs).

  • extra_radius (int) – Precision-window padding beyond the covariance bandwidth (default 1). Raise to 2–3 in the noise-dominated regime β = Tr(Λ)/Tr(2DΔt) ≫ 1, where the windowed precision decays slowly and the default window under-resolves it.

Return type:

None

Notes

Sets standard force_* attributes: force_inferred, force_psf, force_G, force_G_pinv, force_coefficients_full, force_coefficients, force_support, force_moments.

Also sets diffusion_average, A, A_inv, Lambda from the profiled (D, Λ).

When F is a Basis, additionally sets force_basis, force_G_full, and force_scorer so that sparsify_force() can be called afterwards.

See also

Parametric windowed estimators — concepts

Mathematical foundations.

Parametric windowed estimators — algorithm and parameters

Detailed algorithm description.

infer_force_linear(basis, *, preset='auto', M_mode=None, G_mode=None)[source]

Infer the force field as a linear combination of basis functions (linear regression).

[Inference] Linear force regression (overdamped)

\[\hat F(x) = \sum_a C_a\, b_a(x) \qquad\text{where}\qquad G\,C = M\]

\(G_{ab} = \langle b_a(x_t)\, b_b(x_t) \rangle\) is the Gram matrix, \(M_a = \langle v_t \cdot A^{-1} \cdot b_a \rangle\) are the force moments, and \(A = 2\bar D\).

This method computes the force field coefficients (self.force_coefficients) using the provided Basis object. The force field is represented as:

inferred_force(x) = sum_a basis_linear(x)[:,a] * force_coefficients[a]

These coefficients are computed by solving a linear system:

G . force_coefficients = force_moments

and the different options account for the manner to compute G and force_moments. In its simplest form:

G_ab = < b_a(xt) b_b(xt) >                [G_mode = 'rectangle']
force_moments[a] = < dX[t]/dt b_a(xt) >   [mode = 'Ito' ]

but this is rarely the best choice of parameters.

Parameters:

  • basis (Basis) – Basis The fitting functions, encoded as a single callable Basis object (see SFI.statefunc submodule for doc and SFI.bases for examples).

  • preset (str) – Single-keyword convention bundle (mirrors the underdamped engine), default 'auto'. 'auto' picks 'robust' when measurement noise is detected (Tr(Λ) > 0) and 'clean' otherwise. Presets: 'robust' (Stratonovich moments + 'shift' Gram, noise-robust); 'clean' (Itô + 'trapeze', clean / fine-sampled data); 'KM' (Kramers–Moyal: Itô + 'rectangle'); 'legacy-sfi-v1.0' (Stratonovich + 'rectangle', the published SFI v1.0 convention).

  • M_mode (str, optional) – Override the preset’s moment convention: 'Ito', 'Ito-shift', or 'Strato'. None (default) uses the preset.

  • G_mode (str, optional) – Override the preset’s Gram normalization: 'rectangle' (2020 PRX), 'trapeze' (2024 Amiri et al. PRR), or 'shift' (<b_a(xt) b_b(x_t+dt)>, decorrelates measurement noise). None (default) uses the preset.

Notes

Updates the following attributes:

  • self.force_scorer: SparseScorer for model selection.

  • self.force_coefficients: The inferred coefficients for the basis functions.

  • self.force_inferred: Callable function representing the inferred force field.

  • self.force_G: The normalization matrix used in the inference process.

static load_results(path)

Reload a dict previously saved by save_results().

See SFI.inference.serialization.load_results() for details.

Return type:

dict

predict_time_profile(basis_block, t, *, field='force', x=None)

Evaluate a (time-dependent) basis block’s coefficient profile at t.

Contracts the fitted coefficients of basis_block with the basis’s own evaluation at times t (via the reserved time extra), returning the time profile — e.g. -k(t) for the x block of a time-Fourier trap. Avoids re-deriving the design matrix by hand.

Parameters:

field (str)

print_report()

Print a summary report of the inference results.

Provides insights into the inferred diffusion and force fields, along with error metrics such as sampling error, trajectory length, discretization bias, and measurement noise.

report_dict()

Return a structured summary of inference results as a dictionary.

This is the machine-readable counterpart of print_report(). All values are plain Python scalars or numpy arrays (no JAX arrays).

Returns:

Keys include "diffusion_average", "Lambda", "force_information", "force_predicted_MSE", "NMSE_force", "NMSE_diffusion", and others when available. Missing quantities are omitted.

Return type:

dict

save_results(path)

Save report_dict() to <path>.npz + <path>.json.

See SFI.inference.serialization.save_results() for details.

Return type:

Path

simulate_bootstrapped_trajectory(key, oversampling=1, simulate=True, dataset=0)[source]

Simulate an overdamped Langevin trajectory with the inferred force and diffusion fields.

This function generates a trajectory using the inferred force field and diffusion tensor inferred from the input data, matching the original time series and initial conditions.

Parameters:

  • key – JAX random key for generating noise in the simulation.

  • oversampling (int, optional) – Factor for oversampling (i.e. number of intermediate simulated points between two recorded points). Defaults to 1.

  • simulate – if True, performs the simulation with the first data point as initial position; if False, returns an uninitialized object which can be simulated with flexible initial position and parameters.

  • dataset (int, optional) – Which experiment of a pooled fit to reproduce. The inferred model is collapsed to this condition via specialize() (folding per_dataset_scalar / dataset_indicator at dataset); the resulting process is a standalone single-condition model that does not read dataset_index. Defaults to 0.

Returns:

Simulated Langevin process object.

Return type:

OverdampedProcess

sparsify_force(*, criterion='PASTIS', p=0.05, method='beam', max_k=None, **strategy_kwargs)

Sparsify the inferred force by selecting a subset of basis functions.

Builds a Pareto front of sparse models using the chosen method, then selects the model that maximises the given information criterion.

Parameters:

  • criterion ("PASTIS" | "AIC" | "BIC" | "EBIC" | "SIC", default "PASTIS") – Information criterion for model selection.

  • p (float, default 0.05) – Prior-scale parameter \(p_0\) for the PASTIS penalty.

  • method (str, default "beam") –

    Search strategy. One of:

    • "beam" — bidirectional beam search (PASTIS original). Extra kwargs: beam_width (int, default 3), aic_patience (int, default 2).

    • "greedy" — forward stepwise selection. Extra kwargs: direction ("forward" | "backward" | "both", default "forward").

    • "stlsq" — Sequential Thresholded Least Squares (SINDy-style). Extra kwargs: threshold (float or None), mode ("relative" | "absolute"), n_thresholds (int).

    • "lasso" — \(\ell_1\)-penalised coordinate descent. Extra kwargs: alpha (float or None), n_alphas (int).

    • "hillclimb" — stochastic hill-climbing (Gerardos & Ronceray, 2025). Extra kwargs: ic, patience (int), seed (int or None).

  • max_k (int or None) – Maximum model size. Defaults to the full basis size.

  • **strategy_kwargs – Passed to the strategy constructor.

Returns:

The full Pareto-front result, also stored as self.force_sparsity_result.

Return type:

SparsityResult

summary(field='force')

Return a formatted coefficient table for the inferred model.

Parameters:

field ("force" or "diffusion") – Which inferred field to summarize.

Returns:

Multi-line table ready for print().

Return type:

str