API Reference

API Reference#

Quick reference#

Function	Description
`compute_scaledRI`	Compute Scaled RI (RI₀₋₃ by default) from six isoGDGT abundances
`predict_proxy_from_T`	Forward: temperature → proxy percentiles (pure Python)
`predict_T_from_proxyObs`	Inverse: proxy → temperature with full uncertainty (runs Stan)
`compute_quality_flags`	Flag reconstructed samples the calibration cannot support
`download_posteriors`	Download forward posteriors from Zenodo
`download_training_data`	Download training CSVs + CMEMS NO₃ field
`list_posteriors`	Print and return `.nc` stems in the local cache
`MahalanobisOutlierDetector`	Flag samples outside the published calibration domain (no fitting required)
`InvTConfig`	Control the number of forward draws M a reconstruction marginalises over
`build_fwd_data`	Build validated Stan data dict for forward calibration
`get_posterior`	Run forward calibration Stan sampling
`save_posterior`	Persist forward posterior as compressed NetCDF
`load_posterior`	Load a forward or invT posterior from the cache
`summarize_sampler_diagnostics`	Divergences, R-hat, ESS, E-BFMI

Prediction#

Compute Scaled Ring Index#

compute_scaledRI(gdgt0, gdgt1, gdgt2, gdgt3, cren, cren_prime, *, cren_weight=3, cren_rings=None)[source]#

Compute Scaled Ring Index from six isoGDGT abundances.

Accepts raw LC/MS peak areas or fractional abundances — both give identical results because the formula divides by the total sum of all six GDGTs, so any common scale factor drops out. Default cren_weight=3 produces scaledRI_cren3 (RI₀₋₃), the canonical proxy used in TEXAS calibration posteriors.

Parameters:

gdgt0 (float or array-like) – isoGDGT abundances — GDGT-0, GDGT-1, GDGT-2, GDGT-3, crenarchaeol, crenarchaeol regioisomer (cren’). Raw LC/MS peak areas and fractional abundances give the same result (see above).
gdgt1 (float or array-like) – isoGDGT abundances — GDGT-0, GDGT-1, GDGT-2, GDGT-3, crenarchaeol, crenarchaeol regioisomer (cren’). Raw LC/MS peak areas and fractional abundances give the same result (see above).
gdgt2 (float or array-like) – isoGDGT abundances — GDGT-0, GDGT-1, GDGT-2, GDGT-3, crenarchaeol, crenarchaeol regioisomer (cren’). Raw LC/MS peak areas and fractional abundances give the same result (see above).
gdgt3 (float or array-like) – isoGDGT abundances — GDGT-0, GDGT-1, GDGT-2, GDGT-3, crenarchaeol, crenarchaeol regioisomer (cren’). Raw LC/MS peak areas and fractional abundances give the same result (see above).
cren (float or array-like) – isoGDGT abundances — GDGT-0, GDGT-1, GDGT-2, GDGT-3, crenarchaeol, crenarchaeol regioisomer (cren’). Raw LC/MS peak areas and fractional abundances give the same result (see above).
cren_prime (float or array-like) – isoGDGT abundances — GDGT-0, GDGT-1, GDGT-2, GDGT-3, crenarchaeol, crenarchaeol regioisomer (cren’). Raw LC/MS peak areas and fractional abundances give the same result (see above).
cren_weight (float) – Weight carried by both crenarchaeol and its regioisomer, which also sets the normalisation (see Notes). 3 → scaledRI_cren3 / RI₀₋₃ (default, recommended). 4 → scaledRI / RI₀₋₄ (Zhang et al. 2016 convention).
cren_rings (int, optional) – Deprecated alias for cren_weight. Emits a DeprecationWarning.

Returns:

Scaled Ring Index, dimensionless, nominally in [0, 1].

Return type:

numpy.ndarray or float

Notes

The formula is:

RI      = (1·GDGT1 + 2·GDGT2 + 3·GDGT3 + cren_weight·cren + cren_weight·cren')
          / (GDGT0 + GDGT1 + GDGT2 + GDGT3 + cren + cren')
scaledRI = RI / cren_weight

The parameter does two jobs at once: it is the coefficient on cren and cren’ in the numerator and the constant the whole index is divided by, so it fixes the scale on which every sample is expressed, not just the two crenarchaeol terms. That second job is what makes a value chosen here part of the proxy’s definition rather than a detail — a posterior calibrated at 3 cannot read an index built at 4.

It was called cren_rings until 2026-08-21, which was misleading twice over: it named only the first job, and it implied a count of cyclic moieties. Neither 3 nor 4 is that count. They are calibration conventions, which is why both are offered and why neither is “correct”.

Examples

>>> compute_scaledRI(0.45, 0.10, 0.08, 0.05, 0.30, 0.02)
array(0.45666...)

>>> import pandas as pd
>>> df = pd.read_csv("my_gdgt_data.csv")
>>> df["scaledRI_cren3"] = compute_scaledRI(
...     df["GDGT-0"], df["GDGT-1"], df["GDGT-2"], df["GDGT-3"],
...     df["cren"],   df["cren_prime"],
... )

Predict proxy from T#

predict_proxy_from_T(temperatures, posterior, *, n_draws=500, percentiles=[5, 50, 95], return_full=False, seed=42, gdgt23ratio=None, no3=None, no3_cutoff=None, suffix=None, fwd_cache_dir=None)[source]#

Forward prediction: temperature → proxy percentiles (Scaled RI, TEX86, or any fitted proxy).

Samples n_draws self-consistent parameter sets from the forward calibration posterior (all parameters drawn from the same posterior index, preserving correlations) and evaluates the calibration curve at each requested temperature. Corresponds to the forward model described in Eq. 1 / Eq. 6–7 of the manuscript.

Parameters:

temperatures (array-like) – Temperatures (°C) at which to evaluate the calibration curve.
posterior (xr.Dataset or str) – Forward calibration posterior — either a loaded xr.Dataset or a saved-file name string (looked up in the posterior cache).
n_draws (int) – Number of posterior draws to sample. Default 500.
percentiles (list of float) – Percentiles to return, e.g. [5, 50, 95].
return_full (bool) – If True, also return the full (n_draws × len(temperatures)) ensemble array and run metadata under keys "ensemble" and "metadata".
seed (int) – Random seed for reproducible draw sampling.
gdgt23ratio (array-like, optional) – GDGT-2/GDGT-3 ratio values (one per temperature point). Required only when the posterior was fitted with the multivariate model (β_{G₂/₃} correction).
no3 (array-like, optional) – Nitrate concentration values (one per temperature point). Required only when the posterior was fitted with NO₃ correction.
no3_cutoff (float, optional) – Nitrate threshold (μmol/L) below which the NO₃ correction applies. Defaults to the value stored in the posterior attributes.
suffix (str, optional) – Force a specific parameter suffix (e.g. "crtp"). Auto-detected by priority order when omitted.
fwd_cache_dir (Path or str, optional) – Directory to resolve posterior in when it is given as a name string. Defaults to the standard forward posterior cache. Ignored when a loaded Dataset is passed.

Returns:

"x_vals" — temperature array (°C) "pN" — one key per requested percentile, e.g. "p5", "p50", "p95" "ensemble" — full array, shape (n_draws, len(temperatures)), if return_full=True "metadata" — run metadata dict, if return_full=True

Return type:

dict with keys

Predict T from proxy observations#

predict_T_from_proxyObs(proxyObs, prior_mu_t, prior_sigma_t, fwd_posterior=None, *, proxy_name=None, temptype=None, site_name=None, predictors=None, no3=None, gdgt23ratio=None, site_lat=None, site_lon=None, no3_dataset=None, no3_dataset_var='no3_sf2tc_avg', flags=True, tex86=None, config=None, chains=4, iter_warmup=500, iter_sampling=1000, seed=42, constraint_type='unconstrained', min_temp=None, threads_per_chain=None, save_results=False, save_draws=False, filename_tag=None, cache_dir=None, fwd_cache_dir=None)[source]#

Inverse reconstruction: scaled RI → temperature percentiles.

Runs the TEXAS-Bay inverse Stan model to infer paleotemperature from observed scaled Ring Index values. Marginalises over M draws from the forward calibration posterior to propagate calibration uncertainty into the temperature reconstruction. Corresponds to Section 8 (Applications to Paleothermometry) of the manuscript.

Parameters:

proxyObs (array-like, shape (N,)) – Observed proxy values from downcore or coretop samples (e.g. scaledRI, TEX86).
prior_mu_t (float or array-like, shape (N,)) – Prior mean temperature (°C). Scalar applies the same prior to all N observations; array sets a site-specific prior per sample.
prior_sigma_t (float) – Prior temperature uncertainty (°C). Use a diffuse value (e.g. 10) when little prior information is available.
fwd_posterior (str or xr.Dataset, optional) –
The forward calibration posterior. When omitted, the full multivariate T0-shift calibration for temptype is used — tx.GHEB.sst.sri03.G23-N1p0 for SST, tx.GHEB.thm.sri03.G23-N1p0 for thermoT — which ships with the package, so no download is required. Accepts either:
- str — name of the saved posterior (without .nc extension) in the posterior cache directory. The file is loaded automatically.
- xr.Dataset — a pre-loaded posterior Dataset. No file I/O or Zenodo download is attempted; pass this when the cache is unavailable (e.g. Google Colab with a Drive-mounted .nc):
```
ds = xr.open_dataset("my_drive/posterior.nc")
result = predict_T_from_proxyObs(..., fwd_posterior=ds)
```
temptype (str, optional) – Temperature type: "SST" or "thermoT". Optional. It does not change the reconstruction, which follows the calibration supplied: when fwd_posterior is given, the target is read from that posterior’s own attributes and temptype only labels the metadata and output filenames (a value conflicting with the calibration raises a warning). It matters in exactly one case — when fwd_posterior is omitted, it chooses which default calibration is used, SST or thermoT.
site_name (str, optional) – Label attached to result metadata and output filenames.
predictors (dict, optional) – Non-thermal predictor arrays for the N observations, e.g. {"gdgt23ratio": array, "no3": array}. Must be provided when the forward posterior was fitted with the multivariate model. Overridden by no3 / gdgt23ratio shorthands when both are given.
no3 (float or array-like, optional) –
Nitrate concentration (µmol/L) for the N observations.
- Array (length N): per-observation values — use modern WOA23 values extracted at each sample’s location (ocean_prop_ds column "no3_sf2tc_avg").
- Scalar: broadcast to all N observations. Pass a value above no3_cutoff (e.g. no3=10.0 when no3_cutoff=1.0) to effectively disable the NO₃ correction — all observations fall outside the correction window.
Overrides any "no3" key in predictors. Ignored when site_lat / site_lon / no3_dataset are also provided (the lookup result takes priority).
gdgt23ratio (float or array-like, optional) – GDGT-2/GDGT-3 ratio for the N observations. Scalar or array, same broadcast rules as no3. Overrides any "gdgt23ratio" key in predictors.
site_lat (float or array-like, optional) – Decimal latitude(s) of the study site(s). Scalar for a single drill core; array of length N to assign a distinct location to each observation. Requires site_lon and no3_dataset.
site_lon (float or array-like, optional) – Decimal longitude(s) of the study site(s). Same shape rules as site_lat.
no3_dataset (xr.Dataset, optional) – WOA23-derived dataset with a (lat, lon) grid, typically the ocean_prop_ds generated in the preprocessing notebook (SI_code0). Must contain no3_dataset_var. When provided together with site_lat / site_lon, the NO₃ value at those coordinates is looked up via bilinear interpolation and used as the predictor. The result is a scalar (one drill site) or array (per-obs sites), and is broadcast to all N observations when scalar. Optional when site_lat / site_lon are given: if omitted, the ~20 MB ocean_prop_ds field is downloaded from Zenodo and cached automatically (see TEXAS.data.ocean_lookup.get_ocean_prop_ds() / TEXAS.download_ocean_properties()) — pass it explicitly to avoid the download (e.g. a pre-loaded copy, or a Colab session with no persistent cache).
no3_dataset_var (str) – Variable name to extract from no3_dataset. Default "no3_sf2tc_avg".
flags (bool, default True) – Attach result["flags"], a DataFrame with one row per observation marking rows the reconstruction cannot support — most importantly proxy values outside the calibration curve’s attainable range, which return a converged, plausible-looking temperature that is really a readout of the prior. See TEXAS.quality.compute_quality_flags(). Computing them costs no extra sampling.
tex86 (float or array-like, optional) – TEX86 for the same samples, used only by the flags. The published calibration-domain ellipse is two-dimensional (TEX86 × Scaled RI), so the outside_domain check needs both; without this it is reported as pd.NA rather than as passing.
config (InvTConfig, optional) – Controls number of forward-posterior draws (M), seed, etc. Defaults to InvTConfig(), which auto-selects M = min(500, max(100, available_draws // 4)) — M=500 for a typical 4-chain/1000-sample forward posterior (4000 draws), the value the published paleo reconstructions use.
chains (int) – Number of MCMC chains. Default 4.
iter_warmup (int) – Warmup iterations per chain. Default 500.
iter_sampling (int) – Sampling iterations per chain. Default 1000.
seed (int) – Random seed. Default 42.
constraint_type (str) –
Temperature constraint applied in the Stan model:
- "unconstrained" (default): no lower bound; P5 can be unrealistically cold near the calibration curve’s lower asymptote.
- "truncated_prior" (recommended when min_temp is set): proper truncated Normal prior via inverse-CDF reparameterization — P50 is data-driven and P5 is bounded at min_temp without warm bias.
"hard_constraint" (hard lower bound via <lower=min_temp>) was withdrawn in 2026-09: its Jacobian biases P50 warm for polar sites, which is what "truncated_prior" was written to fix. The model is kept in archive/submission-2026-04/stan_models/ and can be run by passing its absolute path as stan_model_path. "reparameterized" and "soft" were listed in earlier type hints but never existed as Stan models; all three now raise ValueError instead of failing later with a missing-file error.
min_temp (float, optional) – Lower temperature bound (°C). Required for "truncated_prior". Typically −1.8 (seawater freezing point). When provided without an explicit constraint_type, automatically selects "truncated_prior".
threads_per_chain (int, optional) – Enable within-chain parallelism via Stan’s reduce_sum.
save_results (bool) – If True, save the quantile posterior .nc and results .npz to the invT cache directory.
save_draws (bool) – If True, also save the raw posterior draws (pre-quantile) as a separate {base}_draws.nc file in the invT cache directory. The file contains t_est with dims (chain, draw, obs_idx) and is suitable for kernel-density plots or custom quantile calculation. Default False.
filename_tag (str or list of str, optional) – Extra tag(s) appended to the output filename.
cache_dir (Path or str, optional) – Directory where .nc and .npz files are written when save_results or save_draws is True. Defaults to the standard invT cache (~/.texas/cache/TEXAS_invT_posterior_cache/ for pip installs, or data/cache/TEXAS_invT_posterior_cache/ in the repo).
fwd_cache_dir (Path or str, optional) – Directory to resolve fwd_posterior in when it is given as a name string. Defaults to the standard forward posterior cache. This is a separate directory from cache_dir, which controls only where results are written.
proxy_name (str | None)

Returns:

"proxyObs" — input proxy array "proxy_name" — proxy type label (e.g. "scaledRI", "TEX86") "p5" — 5th percentile temperature (°C), shape (N,) "p50" — median temperature (°C), shape (N,) "p95" — 95th percentile temperature (°C), shape (N,) "metadata" — run metadata dict (model name, attrs, etc.) "flags" — DataFrame of per-observation quality flags, N rows,

when flags is True

Return type:

dict with keys

Examples

Keep only the observations the calibration can actually support:

>>> result = predict_T_from_proxyObs(ri, prior_mu_t=25, prior_sigma_t=10)
>>> keep = ~result["flags"]["any_flag"].to_numpy(dtype=bool)
>>> sst = result["p50"][keep]

Per-observation quality flags#

Returned automatically as result["flags"], and callable on its own against a saved reconstruction.

compute_quality_flags(proxyObs, fwd_posterior, *, predictors=None, prior_sigma_t=None, result=None, tex86=None, suffix=None, domain_confidence=0.9, prior_dominated_ratio=0.85)[source]#

Flag observations a reconstruction cannot support, one row per input.

Parameters:

proxyObs (array-like, shape (N,)) – The observations handed to the inverse model.
fwd_posterior (xr.Dataset) – The forward calibration the reconstruction used. Both the parameter draws and the attrs recording the calibration’s own data range are read from it.
predictors (dict, optional) – {"gdgt23ratio": ..., "no3": ...} as passed to the reconstruction. Scalars are broadcast. Absent predictors that the calibration uses are reported by predictor_missing.
prior_sigma_t (float, optional) – The temperature prior SD. Required for prior_dominated; without it that column is pd.NA.
result (dict, optional) – The reconstruction result. prior_dominated needs p16/p84 from it; without it that column is pd.NA.
tex86 (float or array-like, optional) – TEX86 for the same samples. The published calibration-domain ellipse is two-dimensional (TEX86 x Scaled RI), so outside_domain can only be evaluated when both are available; it is pd.NA otherwise.
suffix (str, optional) – Parameter suffix to read ("crtp", "culmeso", …). Chosen by the usual priority when omitted.
domain_confidence (float, default 0.90) – Confidence level of the Mahalanobis ellipse.
prior_dominated_ratio (float, default 0.85) – Posterior/prior SD ratio at or above which a row is prior-dominated.

Returns:

N rows. Boolean columns from FLAG_COLUMNS, graded columns from DIAGNOSTIC_COLUMNS, and any_flag.

A check that could not be evaluated is pd.NA, never False — “not assessed” and “passed” are different answers, and collapsing them would quietly launder the first into the second. any_flag therefore ORs only what was actually evaluated, and only the defect columns: ADVISORY_COLUMNS are reported but do not vote.

Return type:

pd.DataFrame

Download and cache#

Download posteriors#

download_posteriors(names=None, cache_dir=None, force=False)[source]#

Download forward calibration posteriors from Zenodo.

Parameters:

names (list of str, optional) – Subset of POSTERIOR_REGISTRY keys to download. When omitted, downloads every posterior of the current record version (~475 MB — dominated by the six full multivariate EIV posteriors at ~78–81 MB each; pass names= to fetch only what you need). Superseded initial-submission posteriors are pinned to the v0.2.0 record and must be requested by name.
cache_dir (Path or str, optional) – Destination directory. Defaults to the standard posterior cache.
force (bool) – Re-download files that already exist locally.

Returns:

Local paths of the downloaded .nc files.

Return type:

list of Path

Examples

Download only the univariate SST posterior (~0.3 MB):

>>> download_posteriors(["gen_logi_fixed_hier_crtp_univ_priorApprox_SST_scaledRI_cren3"])

Download all#

download_all(cache_dir=None, data_dir=None, force=False)[source]#

Download everything from Zenodo: forward posteriors + training data.

Files are downloaded individually; already-cached files are skipped unless force=True. Total download is ~490 MB (dominated by the six full EIV multivariate posteriors at ~78–81 MB each).

Parameters:

cache_dir (Path or str, optional) – Destination for .nc posteriors. Defaults to the standard posterior cache directory.
data_dir (Path or str, optional) – Destination for training data files. Defaults to data/spreadsheets/.
force (bool) – Re-download files that already exist locally.

Return type:

None

Download training data#

download_training_data(dest_dir=None, force=False)[source]#

Download GDGT training data files from Zenodo.

Downloads the coretop/culture/mesocosm training CSVs, the CMEMS NO₃ uncertainty field used in the EIV calibration, and the WOA23-derived ocean_prop_ds gridded ocean properties. The CSVs and the CMEMS field are needed only to re-run the SI preprocessing and calibration notebooks from scratch and are NOT required for inverse temperature reconstructions — use download_posteriors() for that. ocean_prop_ds is the exception: it is also used at inference time by predict_T_from_proxyObs() for the optional site_lat/site_lon NO₃ lookup — see download_ocean_properties() to fetch just that one file. It also isn’t hosted on the TEXAS record at all; it comes from the companion GRL paper’s Zenodo record (see TRAINING_DATA_REGISTRY["ocean_prop_ds"]["record"]).

Parameters:

dest_dir (Path or str, optional) – Destination directory. Defaults to data/spreadsheets/ in the repo (or ~/.texas/data/spreadsheets/ when pip-installed).
force (bool) – Re-download files that already exist locally.

Returns:

Local paths of the downloaded files.

Return type:

list of Path

List posteriors#

list_posteriors(model_type='both', cache_dir=None)[source]#

List available posterior files in the cache directory.

Prints a summary and returns a dict of stem names that can be passed directly to predict_T_from_proxyObs(fwd_posterior=...).

Parameters:

model_type ("forward", "invT", or "both") – Which cache to inspect. Default "both".
cache_dir (Path or str, optional) – Override the default cache root. When given, both forward and invT subdirectories are looked for under this path.

Returns:

{"forward": [...], "invT": [...]} — lists of stem names (no .nc).

Return type:

dict

Set cache directory#

set_cache_dir(path)[source]#

Override TEXAS cache directories at runtime.

Call this before any posterior I/O. For a persistent override, set the TEXAS_CACHE_DIR environment variable instead.

Parameters:: path (str | Path) – Root directory for all TEXAS caches. Three subdirectories will be used inside it: TEXAS_posterior_cache/, TEXAS_invT_posterior_cache/ and TEXAS_kriged_grids_cache/.
Return type:: None

Screening#

Screen a record#

Flags samples lying outside the TEX₈₆–Scaled RI domain the calibration was trained on, using the same 90 % χ² Mahalanobis ellipse used to build the training set. A flagged sample is an extrapolation, not necessarily a bad measurement — the recommendation is to flag and interpret with caution rather than exclude.

The ellipse is a fixed property of the published calibration (the core-top reference cluster of Section 5.1), so nothing needs fitting — an unfitted detector on the standard features adopts it automatically:

from TEXAS.data import MahalanobisOutlierDetector

detector = MahalanobisOutlierDetector(["TEX86", "scaledRI_cren3"], confidence=0.90)
df["flagged"] = detector.detect_outliers_manual(df)

detect_outliers_manual is the manuscript’s criterion: the ellipse plus the warm-end exception, so samples with TEX₈₆ > 0.75 and Scaled RI₀₋₃ > 0.75 are retained rather than flagged — the calibration retains them too (Section 5.1, Fig. 4). detect_outliers applies the bare ellipse and would flag warm samples the calibration actually covers.

Do not call fit() or fit_predict() on your own record: that re-centres the ellipse onto the record, so the more unusual your data the less it flags.

class MahalanobisOutlierDetector(features, confidence=0.9, method='chi2', use_pinv_if_singular=True)[source]#

Fit Mahalanobis distance parameters on training data, apply to any dataset.

Examples

>>> # Fit on coretop data
>>> detector = MahalanobisOutlierDetector(['TEX86', 'ringIndex'])  # default confidence=0.90
>>> detector.fit(coretop_df)
>>> coretop_df['mahal_dist'] = detector.transform(coretop_df)
>>> coretop_df['outliers'] = detector.detect_outliers(coretop_df)
>>>
>>> # Apply to downcore data
>>> downcore_df['mahal_dist'] = detector.transform(downcore_df)
>>> downcore_df['outliers'] = detector.detect_outliers(downcore_df)

Parameters:

features (list[str])
confidence (float)
method (Literal['chi2', 'chi2_rounddown'])
use_pinv_if_singular (bool)

classmethod from_calibration(confidence=0.9)[source]#

Return a detector fixed to TEXAS’s published calibration domain.

Explicit form of the default behaviour: a detector on the standard features screens against this domain without being fitted at all.

The ellipse is a property of the calibration – Section 5.1 fits it to core-top samples with GDGT-2/GDGT-3 <= 5 – not of the data being screened. Fitting on your own record instead re-centres the ellipse onto that record, so the more unusual the record the less it flags: an all-warm Paleogene section would report nothing.

Examples

>>> detector = MahalanobisOutlierDetector.from_calibration()
>>> df["flagged"] = detector.detect_outliers(df)

Parameters:: confidence (float)
Return type:: MahalanobisOutlierDetector

fit(df, *, columns=None, on_unscorable='warn')[source]#

Fit mean, covariance, and threshold on training data.

Parameters:

df (pd.DataFrame) – Training data
columns (dict, optional) – Mapping {logical_name: physical_column} for callers whose DataFrame uses different column names than self.features.
on_unscorable ({'warn', 'raise', 'ignore'}, default 'warn') – Policy for rows that map to a present column but are NaN/Inf.

Returns:

self – Fitted detector instance

Return type:

MahalanobisOutlierDetector

transform(df, col_name=None, *, columns=None, on_unscorable='warn')[source]#

Compute Mahalanobis distances using fitted parameters.

Parameters:

df (pd.DataFrame) – Data to transform
col_name (str, optional) – If provided, add distances to df as this column
columns (dict, optional) – Mapping {logical_name: physical_column}.
on_unscorable ({'warn', 'raise', 'ignore'}, default 'warn') – Policy for rows that map to a present column but are NaN/Inf.

Returns:

distances – Mahalanobis distances

Return type:

pd.Series

detect_outliers(df, col_name=None, *, columns=None, on_unscorable='warn')[source]#

Detect outliers using fitted threshold.

Parameters:

df (pd.DataFrame) – Data to screen
col_name (str, optional) – If provided, add outlier flags to df as this column
columns (dict, optional) – Mapping {logical_name: physical_column}.
on_unscorable ({'warn', 'raise', 'ignore'}, default 'warn') – Policy for rows that map to a present column but are NaN/Inf.

Returns:

outliers – Boolean series (True=outlier, False=inlier, NaN=invalid)

Return type:

pd.Series

detect_outliers_manual(df, col_name=None, exclude_condition=None, *, columns=None, on_unscorable='warn')[source]#

Detect outliers using the manuscript’s screening criterion.

This is the recommended call for screening a record. It applies the calibration-domain ellipse and then restores the warm-end exception: samples with TEX86 > 0.75 and Scaled RI(0-3) > 0.75 are retained, not flagged, because the calibration itself retains them (Section 5.1, Fig. 4 – “data with TEX86 and Scaled RI above 0.75 are retained in the calibration dataset”). Screening with the bare ellipse instead (detect_outliers()) would flag warm samples the calibration actually covers.

Needs no fit(): on the standard features the published calibration domain is adopted automatically.

Parameters:

df (pd.DataFrame) – Data to screen
col_name (str, optional) – If provided, add manual outlier flags to df as this column
exclude_condition (pd.Series, optional) – Boolean series indicating samples to exclude from outlier detection. If None, applies default: (ringIndex > 3) & (TEX86 > 0.7)
columns (dict, optional) – Mapping {logical_name: physical_column} for callers whose DataFrame uses different column names than self.features.
on_unscorable ({'warn', 'raise', 'ignore'}, default 'warn') – Policy for rows that map to a present column but are NaN/Inf.

Returns:

manual_outliers – Boolean series with manual exceptions applied

Return type:

pd.Series

Examples

>>> # Use default exception
>>> outliers = detector.detect_outliers_manual(df)
>>>
>>> # Custom exception
>>> custom_exclude = (df['SST'] > 30) & (df['TEX86'] > 0.8)
>>> outliers = detector.detect_outliers_manual(df, exclude_condition=custom_exclude)

fit_predict(df, col_name=None, *, columns=None, on_unscorable='warn')[source]#

Fit on df and flag its rows in one call (sklearn-style).

Equivalent to fit(df) followed by detect_outliers(df).

Warning

This fits the ellipse on df itself, which is right when df is a reference/training set and wrong for screening a record against the calibration – the domain would move with the data, so an unusual record flags less rather than more. To screen your own data, just call detect_outliers(): an unfitted detector on the standard features uses TEXAS’s published calibration domain automatically.

Parameters:

df (pd.DataFrame) – Data to fit on and screen
col_name (str, optional) – If provided, add outlier flags to df as this column
columns (dict, optional) – Mapping {logical_name: physical_column}.
on_unscorable ({'warn', 'raise', 'ignore'}, default 'warn') – Policy for rows that map to a present column but are NaN/Inf.

Returns:

outliers – Boolean series (True=outlier, False=inlier, NaN=invalid)

Return type:

pd.Series

fit_transform(df, dist_col=None, outlier_col=None, manual_outlier_col=None, *, columns=None, on_unscorable='warn')[source]#

Fit and transform in one step.

Parameters:

df (pd.DataFrame) – Training data
dist_col (str, optional) – Column name for distances
outlier_col (str, optional) – Column name for outlier flags
manual_outlier_col (str, optional) – Column name for manual outlier flags
columns (dict, optional) – Mapping {logical_name: physical_column}.
on_unscorable ({'warn', 'raise', 'ignore'}, default 'warn') – Policy for rows that map to a present column but are NaN/Inf.

Returns:

results – Dictionary containing: - ‘distances’: Mahalanobis distances - ‘outliers’: Outlier flags - ‘manual_outliers’: Manual outlier flags - ‘threshold’: Computed threshold

Return type:

dict

get_params()[source]#

Get fitted parameters.

Returns:: params – Dictionary of fitted parameters
Return type:: dict

Data builders#

Build forward data#

build_fwd_data(*, t_cul=None, proxy_cul=None, t_meso=None, proxy_meso=None, t_crtp=None, proxy_crtp=None, gdgt23ratio_crtp=None, sd_gdgt23ratio_crtp=None, no3_crtp=None, sd_no3_crtp=None, no3_cutoff=None, proxy_residuals_crtp=None, sd_proxyObs=None, R2_thermal=None, culmeso_posterior=None, prior_mean_t0=None, prior_sd_t0=None, prior_mean_k=None, prior_sd_k=None, prior_mean_b=None, prior_sd_b=None, prior_mean_v=None, prior_sd_v=None)[source]#

Build the Stan data dictionary for forward calibration models.

Handles all forward model variants:

culmeso / Q1_culmeso / v1_culmeso : pass t_cul, proxy_cul, t_meso, proxy_meso
culmesocore : add t_crtp, proxy_crtp
hier_crtp_multiv : add gdgt23ratio_crtp, no3_crtp
hier_crtp_multiv_priorApprox : add culmeso_posterior (extracts hyperpriors)
hier_crtp_univ_priorApprox : add culmeso_posterior (no predictors needed)

Parameters:

t_cul – Culture temperature and proxy arrays.
proxy_cul – Culture temperature and proxy arrays.
t_meso – Mesocosm temperature and proxy arrays.
proxy_meso – Mesocosm temperature and proxy arrays.
t_crtp – Coretop temperature and proxy arrays.
proxy_crtp – Coretop temperature and proxy arrays.
gdgt23ratio_crtp – GDGT-2/GDGT-3 ratio for coretop samples. Sets use_gdgt23ratio=1 if non-zero/non-NaN.
sd_gdgt23ratio_crtp – Per-site measurement SE of gdgt23ratio (same units, linear). Required for the _eiv model; always included in the data dict (defaults to zeros when not provided, which disables the G₂/₃ EIV measurement model).
no3_crtp – Nitrate concentration for coretop samples. Sets use_no3=1 if non-zero/non-NaN.
sd_no3_crtp – Per-site measurement SE of NO₃ (μmol/L, linear space). Required for the _eiv model. Always included (defaults to zeros; sites with sd=0 receive only the lognormal prior and skip the normal measurement model).
no3_cutoff (float | None) – NO3 threshold for the nonthermal correction. Priority: (1) this arg, (2) culmeso_posterior attrs, (3) auto-calculated via Spearman method.
proxy_residuals_crtp (ndarray | None) – Pre-computed proxy residuals for NO3 threshold calculation. If omitted, residuals are computed internally by fitting a generalized logistic curve. Warning: Stan models use generalized logistic — residuals from other functional forms may shift the threshold.
culmeso_posterior (Dataset | None) – xr.Dataset from a completed culmeso forward run. Auto-extracts prior_mean_*/prior_sd_* hyperpriors and no3_cutoff (if saved in attrs).
prior_mean_*/prior_sd_* – Manual hyperprior values. Override auto-extracted values from culmeso_posterior for individual params.
R2_thermal (float | None)

Returns:

Stan-ready data dict with proxyObs_* keys, N_* counts, use_* flags,: and hyperpriors — ready for get_posterior().

Return type:

dict

Build invT input data#

build_invT_inputData(proxyObs=None, prior_mu_t=None, prior_sigma_t=None, *, scaledRI=None, fwd_posterior_name=None, predictors=None, config=None, fwd_posterior=None, fwd_cache_dir=None)[source]#

Build the data dictionary for Stan’s inverse model and sampler configuration.

WORKFLOW: ───────── 1. Load forward calibration posterior from .nc file (or accept a pre-loaded Dataset) 2. Randomly sample M parameter sets from that posterior 3. Extract calibration curve parameters (t0, k, b, v, sigma) 4. Package optional environmental predictors (GDGT-2/3, NO3) if used 5. Return data dict (for Stan) + sampler_kwargs (for CmdStanPy)

Parameters:

proxyObs (ndarray | List[float]) – Observed proxy values to predict temperature from (length N). Any proxy is accepted: scaledRI, TEX86, ringIndex, etc.
prior_mu_t (ndarray | float) – Prior mean temperature (scalar or array of length N)
prior_sigma_t (float) – Prior temperature uncertainty (e.g., 10°C)
fwd_posterior_name (str | None) – Name of saved forward calibration (without .nc extension). Not required when fwd_posterior is supplied directly.
predictors (Dict[str, ndarray] | None) – Optional environmental covariates {‘gdgt23ratio’: array, ‘no3’: array}
config (InvTConfig | None) – Configuration object controlling M, seed, etc.
fwd_posterior (Dataset | None) – Pre-loaded forward posterior xr.Dataset. When provided, fwd_posterior_name is ignored and no file I/O is performed. Useful when running from Google Colab or any pip-install context where the posterior cache is not available.
fwd_cache_dir (str | Path | None) – Directory to resolve fwd_posterior_name in. Defaults to the standard forward posterior cache. Use this to read posteriors from a project-local bundle without copying them into the cache.
scaledRI (ndarray | List[float])

Returns:

Dictionary for Stan’s data block sampler_kwargs: Dictionary for CmdStanPy sampling configuration

Return type:

data

InvT configuration#

class InvTConfig(mode='ensemble', n_draws=None, seed=42, no3_cutoff=0.0, suffix=None)[source]#

Configuration for the inverse-T Stan data builder.

Parameters:

mode (str)
n_draws (int | None)
seed (int)
no3_cutoff (float | None)
suffix (str | None)

WOA23 NO₃ lookup#

lookup_no3_from_woa(lat, lon, woa_dataset, variable='no3_sf2tc_avg', method='linear')[source]#

Look up modern NO₃ at one or more lat/lon coordinates from a WOA23-derived xarray Dataset.

The dataset is typically the preprocessed ocean_prop_ds generated in SI_code1, which contains thermocline-depth-integrated WOA23 climatology on a regular (lat, lon) grid. The returned value(s) are time-invariant (climatological mean) and intended as a modern-ocean proxy for the NO₃ correction in paleo reconstructions.

Parameters:

lat (float or array-like) – Latitude(s) in decimal degrees (−90 to 90). Pass a scalar for a single drill site; pass an array of length N to match N observations.
lon (float or array-like) – Longitude(s) in decimal degrees. Both −180–180 and 0–360 conventions are accepted — the function normalises to match the dataset’s convention automatically.
woa_dataset (xr.Dataset) – WOA23-derived dataset with a (lat, lon) grid containing variable. Dimensions must be named "lat" and "lon".
variable (str) – Name of the NO₃ variable to extract. Default "no3_sf2tc_avg" (thermocline depth-integrated annual average from SI_code1).
method ({"linear", "nearest"}) – Interpolation method. "linear" (default) performs bilinear interpolation and is preferred for smooth fields. "nearest" snaps to the closest grid cell and is useful when the dataset is sparse or has NaN-masked shelves.

Returns:

NO₃ value(s) in µmol/L. Shape matches the scalar/array input: a 0-d array for scalar inputs, 1-d array of length N for array inputs. NaN is returned for locations outside the dataset’s valid range (e.g. continental shelves masked in WOA23).

Return type:

np.ndarray

Raises:

KeyError – If variable is not found in woa_dataset.
ValueError – If woa_dataset does not have "lat" and "lon" dimensions.

Examples

Single drill site:

>>> no3_val = lookup_no3_from_woa(15.3, -23.7, ocean_prop_ds)
>>> # returns scalar-equivalent float; broadcasts to all N obs automatically
>>> result = predict_T_from_proxyObs(..., no3=no3_val)

Multi-site stack (per-obs lookup):

>>> no3_arr = lookup_no3_from_woa(core_df["lat"].values,
...                                core_df["lon"].values,
...                                ocean_prop_ds)
>>> result = predict_T_from_proxyObs(..., no3=no3_arr)

Forward calibration#

Get posterior#

get_posterior(data, stan_file, temptype, proxy_name, *, iter_warmup=None, iter_sampling=None, threads_per_chain=None, chains=None, parallel_chains=None, adapt_delta=None, max_treedepth=None, **kwargs)[source]#

Run forward calibration Stan sampling and return the posterior.

Wraps StanSampler with automatic predictor detection, CPU configuration, and metadata attachment. The returned dataset can be passed directly to predict_proxy_from_T or saved with save_posterior.

Parameters:

data (dict) – Stan data dict built by build_fwd_data(). Predictor flags (use_gdgt23ratio, use_no3) are auto-detected from the arrays present; you do not need to set them manually.
stan_file (str) – Stan model name (without .stan), e.g. "gen_logi_fixed_hier_crtp_multiv_priorApprox_eiv".
temptype (str) – Temperature variable type, e.g. "SST" or "thermoT". Stored in the posterior metadata.
proxy_name (str) – Proxy type, e.g. "scaledRI_cren3". Required — stored in the .nc attrs and validated downstream when the posterior is used for inverse reconstruction.
iter_warmup (int, optional) – HMC warmup iterations per chain (default: CmdStan default, 1000).
iter_sampling (int, optional) – Post-warmup sampling iterations per chain (default: 1000).
chains (int, optional) – Number of independent chains (default: 4).
parallel_chains (int, optional) – Chains to run simultaneously (auto-detected from CPU count).
threads_per_chain (int, optional) – Threads per chain for reduce_sum models (auto-enabled for models whose filename contains reduce_sum).
adapt_delta (float, optional) – Target acceptance rate (default: 0.8). Increase toward 0.99 to reduce divergences at the cost of more leapfrog steps.
max_treedepth (int, optional) – Maximum tree depth for HMC (default: 10).
**kwargs – Additional keyword arguments forwarded to CmdStanModel.sample.

Returns:

posterior (xr.Dataset) – Forward calibration posterior with parameter draws and metadata attrs (model name, temptype, proxy_name, priors, diagnostics).
diagnostics (str) – Human-readable sampler diagnostic summary (divergences, R-hat, ESS, E-BFMI).

Raises:

ValueError – If active predictors are present but a univariate model is requested, or if an EIV model is requested without R2_thermal.

Return type:

Tuple[Dataset, str]

Examples

>>> data = build_fwd_data(t_crtp=..., proxy_crtp=..., ...)
>>> posterior, diag = get_posterior(
...     data,
...     stan_file="gen_logi_fixed_hier_crtp_univ_priorApprox",
...     temptype="SST",
...     proxy_name="scaledRI_cren3",
... )
>>> save_posterior(posterior)

Save posterior#

save_posterior(posterior, cache_dir=None, overwrite=True, filename_suffix='', layout='auto', run='auto')[source]#

Save a forward-model posterior to disk as compressed NetCDF.

The filename is auto-generated from the posterior’s metadata attrs: {model}_{temptype}[_gdgt23ratio][_no3_{cutoff}][_{proxy_name}]{suffix}.nc

Parameters:

posterior (xr.Dataset) – Forward calibration posterior returned by get_posterior(). Must have stan_model_name, temptype, and proxy_name attrs set (proxy_name is required — a warning is raised if missing).
cache_dir (str or Path, optional) – Directory to write the file. Defaults to the standard forward posterior cache (data/cache/TEXAS_posterior_cache/ for source installs, ~/.texas/cache/TEXAS_posterior_cache/ for pip installs).
overwrite (bool) – If False, raise FileExistsError when the output path already exists. Default True.
filename_suffix (str, optional) – Extra tag appended before .nc under the legacy layout, and the run/member token under the case layout. Leading/trailing underscores are stripped. Do not pass a date here. Filenames no longer carry date stamps; the run date is recorded in the run_timestamp attr, where it survives a rename and does not have to be parsed back out of a path. Prefer run= for an explicit member.
run (str or int, optional) – Optional run/member token under the case layout. "auto" (default) writes no member, so a configuration has one canonical path and a re-run replaces it; pass a value only to keep a run deliberately apart. An explicit run wins over filename_suffix for the case path. With overwrite=False an existing file raises FileExistsError.
layout ({"auto", "case", "legacy"}) – Where to write. "case" uses the CESM-style case directory (tx.v026.GHEB.sst.ri3.G23-N1p0/fwd.nc); "legacy" uses the historical long flat filename; "auto" (default) prefers the case layout and falls back to legacy with a warning if no case id can be derived. See TEXAS.utils.naming.

Returns:

Absolute path of the saved .nc file.

Return type:

Path

Load posterior#

load_posterior(model_name, model_type='forward', cache_dir=None)[source]#

Load a posterior from disk: {model_name}.nc in the appropriate cache directory.

Parameters:

model_name (str) – Name of the model file (without .nc extension)
model_type (Literal['forward', 'invT']) – Type of posterior (“forward” or “invT”)
cache_dir (str | Path | None) – Custom cache directory (overrides default locations)

Returns:

xarray.Dataset containing the posterior

Raises:

FileNotFoundError – If the posterior file doesn’t exist

Return type:

Dataset

Ensemble#

Generate ensemble (auto)#

generate_ensemble_auto(post_ds, x_vals, model_type='auto', gdgt23ratio=None, no3=None, no3_cutoff=None, return_full_ensemble=False, suffix=None, **kwargs)[source]#

Sample draws from a forward posterior and compute calibration-curve percentiles.

Inspects the posterior’s stan_model_name attr and data_vars to determine the model function, parameter names, and optional-predictor flags automatically, then delegates to generate_ensemble.

Parameters:

post_ds (xr.Dataset) – Forward calibration posterior returned by get_posterior() or loaded with load_posterior().
x_vals (np.ndarray) – Temperature values (°C) at which to evaluate the calibration curve.
model_type ({"auto", "forward", "inverse"}) – Force forward or inverse dispatch; "auto" (default) infers from the posterior. InvT posteriors are not supported — use predict_T_from_proxyObs() instead.
gdgt23ratio (np.ndarray, optional) – GDGT-2/3 ratio values; required when the posterior was fitted with a multivariate (GDGT-2/3) model.
no3 (np.ndarray, optional) – NO₃ concentrations (µmol/L); required when the posterior uses the NO₃ correction.
no3_cutoff (float, optional) – Override the NO₃ cutoff from the posterior attrs.
return_full_ensemble (bool) – If True, return the full M × N draw matrix in addition to percentiles. Default False.
suffix (str, optional) – Force a specific parameter suffix (e.g. "crtp"); overrides auto-detection.
**kwargs – Forwarded to generate_ensemble.

Returns:

Keys "p1" … "p99" (and optionally "ensemble") — each a numpy array of length len(x_vals).

Return type:

dict

Raises:

NotImplementedError – If called with an invT posterior.

Detect model and params#

detect_model_and_params(posterior_ds, suffix=None)[source]#

Auto-detect which logistic model and parameters to use. Uses a shared suffix priority via choose_suffix().

Parameters:

posterior_ds (Dataset)
suffix (str)

Diagnostics#

Sampler diagnostics#

summarize_sampler_diagnostics(fit)[source]#

Extract divergent__, treedepth__, E-BFMI, R_hat, and ESS_bulk from a CmdStanPy fit and return them as stan_diag_* attrs.

Return type:: dict

Summary table#

create_summary_table(datasets)[source]#

Build a DataFrame summarizing the stan_diag_* attrs from each xarray.Dataset.

Parameters:: datasets (list)
Return type:: DataFrame

Plotting#

Plot prior distributions#

plot_prior_distributions(priors_list=None, posterior_datasets=None, posterior_labels_list=None, show_suptitle=True, kde_bw=0.3, focus_on_posterior=True, include_groups=('t0', 'k', 'b', 'v', 'a', 'beta_G23', 'beta_NO3', 'gamma_G23', 'gamma_NO3', 'betaLogit_G23', 'betaLogit_NO3', 'sigma_proxyObs'), suffix_include=None, zoomin_suffix=None, zoomin_dataset_idx=None, use_linestyle_by_param=False, show_histogram=True, show_annotation=False, set_linewidth=1.5, set_fig_width_factor=3, set_fig_height_factor=3.5, set_leg_max_ncol=3, color_list=None, param_source_map=None, annotation_style='ci95', show_subplot_legend=True, show_figure_legend=True, show_prior_expression=True, cache_dir=None)[source]#

Plot priors + any number of posterior distributions in a grid, split by parameter group (t0, k, b, etc.).

Parameters:

param_source_map (Dict[str, int] | None) –
Optional dict mapping a param group name to the index of the dataset in posterior_datasets that should be used as the sole source for that group. All other datasets are skipped for that group.

Use this when different parameters come from different posteriors — e.g. logistic params (t0, k, b…) from a culmeso run and beta coefficients from a multivariate crtp run:
```
plot_prior_distributions(
    posterior_datasets=[culmeso_ds, crtp_multiv_ds],
    param_source_map={"beta_G23": 1, "beta_NO3": 1},
)
```
When a group is not in param_source_map, all datasets are searched as usual.
include_groups (Sequence[str]) – Parameter-group prefixes to draw, in panel order. The default covers the thermal parameters plus every coefficient naming scheme in PREDICTOR_GROUPS (beta_* for the parent model, gamma_* for t0shift, betaLogit_* for boundedCeil), so the same call works for any of them and for figures mixing several. A group with no matching variable in any dataset is skipped, as is a coefficient group whose predictor no dataset switched on.
cache_dir (str | Path | None) – Directory to resolve any name strings in posterior_datasets. Defaults to the standard forward posterior cache. Loaded Datasets are passed through untouched.
priors_list (List[str] | Dict[str, str] | None)
posterior_datasets (List[xr.Dataset] | None)
posterior_labels_list (List[str] | None)
show_suptitle (bool)
kde_bw (float)
focus_on_posterior (bool)
suffix_include (List[str] | None)
zoomin_suffix (str | List[str] | None)
zoomin_dataset_idx (int | None)
use_linestyle_by_param (bool)
show_histogram (bool)
show_annotation (bool)
set_linewidth (float)
set_fig_width_factor (float)
set_fig_height_factor (float)
set_leg_max_ncol (int)
color_list (Sequence[str] | None)
annotation_style (Literal['ci95', 'ci68', 'sigma'])
show_subplot_legend (bool)
show_figure_legend (bool)
show_prior_expression (bool)

API Reference

Contents

API Reference#

Quick reference#

Prediction#

Compute Scaled Ring Index#

Predict proxy from T#

Predict T from proxy observations#

Per-observation quality flags#

Download and cache#

Download posteriors#

Download all#

Download training data#

List posteriors#

Set cache directory#

Screening#

Screen a record#

Data builders#

Build forward data#

Build invT input data#

InvT configuration#

WOA23 NO₃ lookup#

Forward calibration#

Get posterior#

Save posterior#

Load posterior#

Ensemble#

Generate ensemble (auto)#

Detect model and params#

Diagnostics#

Sampler diagnostics#

Summary table#

Plotting#

Plot prior distributions#