analysis package

Submodules

PyHydroGeophysX.analysis.sensitivity module

Sensitivity and resolution analysis utilities.

PyHydroGeophysX.analysis.sensitivity.compute_cumulative_sensitivity(J: Any) → Any

Compute cumulative sensitivity as absolute column sums of the Jacobian.

PyHydroGeophysX.analysis.sensitivity.compute_depth_of_investigation(inv_class: Any, data: Any, mesh: Any, scale_low: float = 0.8, scale_high: float = 1.2, reference_resistivity: float | None = None) → Tuple[ndarray, Dict[str, ndarray]]

Estimate DOI using two inversions with different reference/initial models.

This follows the Oldenburg-Li style idea of quantifying model sensitivity to reference-model choice. The two starting models are reference_resistivity * scale_low and reference_resistivity * scale_high; when reference_resistivity is None it is taken as the median apparent resistivity of data so the perturbation brackets the survey’s actual background instead of an arbitrary absolute value.

The index is only as meaningful as the inversion’s use of its reference. ERTInversion regularizes with first-order smoothness, which ignores a homogeneous reference, so give it a reference_weight for this (a smallness term, as Oldenburg and Li (1999) use). Let both runs converge as well: a run stopped at its target misfit keeps its starting model wherever the data are weak, and the index then reflects where each run stopped. A cell that follows its reference completely scores (scale_high - scale_low) / (scale_high + scale_low); dividing by it gives Oldenburg and Li’s index, from 0 to 1.

PyHydroGeophysX.analysis.sensitivity.compute_resolution_matrix(J: Any, Wd: Any, Wm: Any, lam: Any) → Any

Compute model resolution matrix:

R = (J^T Wd^T Wd J + lam Wm^T Wm)^(-1) J^T Wd^T Wd J

Wd whitens the data, so the data term is (Wd J)^T (Wd J). A whitening matrix need not be symmetric (the inverse Cholesky factor of a data covariance is triangular), and J^T Wd Wd J, which this used to form, is then a different matrix; for a symmetric or diagonal Wd the two are the same.

lhs is the Gauss-Newton normal matrix, which is symmetric positive definite, so the inverse is applied as a Cholesky solve rather than formed. The two agree here to within 5e-9 across condition numbers from 1e4 to 1e8, so this is a cost change rather than an accuracy one: the pseudo-inverse it replaces is an SVD, measured at 16.5 times the Cholesky’s runtime for a matrix of order 800.

PyHydroGeophysX.analysis.sensitivity.plot_sensitivity_map(sensitivity: Any, mesh: Any, ax: Any = None) → Any

Plot a sensitivity map on a mesh (or as a simple line plot fallback).

Module contents

Sensitivity and resolution analysis utilities.

PyHydroGeophysX.analysis.compute_cumulative_sensitivity(J: Any) → Any

Compute cumulative sensitivity as absolute column sums of the Jacobian.

PyHydroGeophysX.analysis.compute_depth_of_investigation(inv_class: Any, data: Any, mesh: Any, scale_low: float = 0.8, scale_high: float = 1.2, reference_resistivity: float | None = None) → Tuple[ndarray, Dict[str, ndarray]]

Estimate DOI using two inversions with different reference/initial models.

This follows the Oldenburg-Li style idea of quantifying model sensitivity to reference-model choice. The two starting models are reference_resistivity * scale_low and reference_resistivity * scale_high; when reference_resistivity is None it is taken as the median apparent resistivity of data so the perturbation brackets the survey’s actual background instead of an arbitrary absolute value.

The index is only as meaningful as the inversion’s use of its reference. ERTInversion regularizes with first-order smoothness, which ignores a homogeneous reference, so give it a reference_weight for this (a smallness term, as Oldenburg and Li (1999) use). Let both runs converge as well: a run stopped at its target misfit keeps its starting model wherever the data are weak, and the index then reflects where each run stopped. A cell that follows its reference completely scores (scale_high - scale_low) / (scale_high + scale_low); dividing by it gives Oldenburg and Li’s index, from 0 to 1.

PyHydroGeophysX.analysis.compute_resolution_matrix(J: Any, Wd: Any, Wm: Any, lam: Any) → Any

Compute model resolution matrix:

R = (J^T Wd^T Wd J + lam Wm^T Wm)^(-1) J^T Wd^T Wd J

Wd whitens the data, so the data term is (Wd J)^T (Wd J). A whitening matrix need not be symmetric (the inverse Cholesky factor of a data covariance is triangular), and J^T Wd Wd J, which this used to form, is then a different matrix; for a symmetric or diagonal Wd the two are the same.

lhs is the Gauss-Newton normal matrix, which is symmetric positive definite, so the inverse is applied as a Cholesky solve rather than formed. The two agree here to within 5e-9 across condition numbers from 1e4 to 1e8, so this is a cost change rather than an accuracy one: the pseudo-inverse it replaces is an SVD, measured at 16.5 times the Cholesky’s runtime for a matrix of order 800.

PyHydroGeophysX.analysis.plot_sensitivity_map(sensitivity: Any, mesh: Any, ax: Any = None) → Any

Plot a sensitivity map on a mesh (or as a simple line plot fallback).