calibration

Description of the calibration functions implemented in the Python library bnn_for_14C_calibration:

bnn_for_14C_calibration.calibration.individual_calibration(c14age, c14sig, alpha=0.05, covariables=False, mesure_likelihood=['gaussian_mixture', 'curve_gaussian_approximation'][0], compute_calage_posterior_mean_and_std=False, sample_size=1000)

Perform an individual (independent) radiocarbon calibration using precomputed BNN predictions.

This function loads precomputed midpoint predictions for the two parts of the calibration curve (recent and older parts), constructs a posterior density on a common rescaled time axis, computes the Highest Posterior Density (HPD) region for the calibrated date, and returns a dictionary with results and diagnostics. Optionally, it can also compute Monte Carlo estimates of the posterior mean and standard deviation by sampling from the approximate posterior density.

Parameters:

Name	Type	Description	Default
c14age	float	The measured radiocarbon age (conventional radiocarbon age) in years BP (Before Present).	required
c14sig	float	The 1-sigma measurement uncertainty of c14age.	required
alpha	float	Tail probability for HPD region calculation (default 0.05).	0.05
covariables	bool	If True, use models trained with covariates (Beryllium-10 and Earth's geomagnetic field paleo-intensity); otherwise use models without covariates. Default is False.	False
mesure_likelihood	str	Which measurement-likelihood model to use. Must be one of: "gaussian_mixture", "curve_gaussian_approximation". Default is "gaussian_mixture".	"gaussian_mixture"
compute_calage_posterior_mean_and_std	bool	If True, sample from the approximate posterior (via the piecewise approximation) and compute posterior mean and standard deviation. Default is False.	False
sample_size	int	Number of posterior samples to draw when compute_calage_posterior_mean_and_std is True. Default is 1000.	1000

Returns:

Name	Type	Description
results	Dict[str, Any]	A dictionary containing many entries including (but not limited to): "calage_posterior_mode" : int Posterior mode (most probable calibrated date, rounded to int). "calage_posterior_mode_density" : float Posterior density (scaled) at the mode. "connexe_HPD_intervals" : ndarray Array of connected HPD intervals in the rescaled domain. "connexe_HPD_intervals_unscaled" : ndarray Same HPD intervals converted back to original date units. "connexe_HPD_intervals_unscaled_round" : ndarray HPD intervals rounded to integers (useful as year ranges). "HPD_region_length" : int Total length of HPD regions in original date units (rounded). "middle_points" : ndarray Middle points (unscaled) used to compute the posterior density. "middle_points_density" : ndarray Posterior density evaluated at middle_points. "calage_posterior_mean" : Optional[int] Posterior mean (int) if sampling was requested, otherwise None. "calage_posterior_std" : Optional[int] Posterior std (int) if sampling was requested, otherwise None. "calage_sample" : Optional[ndarray] Posterior samples (in original date units) if sampling was requested, otherwise None. "c14age", "c14sig", "covariables", "alpha" : input parameters echoed back.

Notes

Mathematical formulation of measurement-likelihood model

The posterior density of the calibrated date \(d\) (associated to a measurement \(m\) of std \(\sigma\)) is proportional to the product of the prior and the likelihood: $$ p(d | m) \propto p(d) \times L(m | d), $$ where the likelihood term \( L(m | d) \) is given by: $$ L(m|d) = E_{BNN}[ \exp(-(m - F^{14}C(d))^2 / (2σ^2)) ]. $$

Depending on the chosen mesure_likelihood, the likelihood term is estimated by:

1. Gaussian Mixture Approximation $$ \hat{L}(m|d) = \frac{1}{N} \sum_{i=1}^N \frac{1}{\sqrt{2\pi}\sigma} \exp\left[-\frac{(m - \hat{F}^{14}_iC(d))^2}{2\sigma^2}\right] $$ where \( \hat{F}^{14}_iC(d) \) are stochastic predictions from the BNN.

2. Curve Gaussian Approximation (Central Limit Theorem) $$ \hat{L}(m|d) = \frac{1}{\sqrt{2\pi(\sigma^2 + s_d^2)}} \exp\left[-\frac{(m - \mu_d)^2}{2(\sigma^2 + s_d^2)}\right] $$ where \( \mu_d \) and \( s_d \) are the mean and std of BNN predictions.

Approximation and sampling strategy of the posterior distribution

The function relies on precomputed midpoint predictions loaded by bnn_load_predictions_ and follows the piecewise approximation / sampling strategy used elsewhere in this package (as described for example in mono_cal_date_approx_density_sample).

Examples:

>>> # radiocarbon data
>>> c14age = 10400
>>> c14sig = 18

>>> # independent calibration using BNN curve
>>> individual_calib_results = individual_calibration(
...    c14age, c14sig, 
...    sample_size=10000,
...    compute_calage_posterior_mean_and_std=True
...)
>>> print(individual_calib_results)
{'calage_posterior_mode': 12208, 'calage_posterior_mode_density': 0.004876364275474341, 
'connexe_HPD_intervals': array([[0.18939556, 0.19895774]]), 
'connexe_HPD_intervals_density': [0.9500607994817136], 
'HPD_threshold': 0.0007750276513481828, 'connexe_HPD_intervals_unscaled': array([[12051.406, 12660.058]]), 
'connexe_HPD_intervals_unscaled_round': array([[12051, 12661]]), 'HPD_region_length': 609, 'alpha': 0.05, 
'middle_points': array([-3.38430000e+00, -2.15290000e+00, -9.21500000e-01, ...,
    5.49893275e+04,  5.49935965e+04,  5.49978655e+04]), 
'middle_points_density': array([0., 0., 0., ..., 0., 0., 0.]), 
'calage_posterior_mean': 12254, 'calage_posterior_std': 140, 
'calage_sample': array([12188.33655321, 11929.49681614, 12206.41856239, ...,
   12450.48293508, 12188.77765758, 12446.75205379]), 'c14age': 10400, 'c14sig': 18, 'covariables': False}

bnn_for_14C_calibration.calibration.IntCal20_calibration(c14age, c14sig, alpha=0.05, compute_calage_posterior_mean_and_std=False, sample_size=1000)

Perform Bayesian calibration of a radiocarbon measurement using the IntCal20 calibration curve.

Parameters:

Name	Type	Description	Default
c14age	float	Radiocarbon measured age (in BP).	required
c14sig	float	Laboratory uncertainty associated with the radiocarbon age (standard deviation, in BP).	required
alpha	(float, optional(default=0.05))	Significance level used to compute Highest Posterior Density (HPD) region. The returned region covers posterior probability mass equal to (1 - alpha).	0.05
compute_calage_posterior_mean_and_std	(bool, optional(default=False))	If True, posterior samples of the calibrated age are drawn and used to compute: - posterior mean, - posterior standard deviation. If False, the sampling step is skipped.	False
sample_size	(int, optional(default=1000))	Number of posterior samples to draw when compute_calage_posterior_mean_and_std = True.	1000

Returns:

Name	Type	Description
results	Dict[str, Any]	A dictionary containing: 'connexe_HPD_intervals' : array of HPD intervals in scaled space 'connexe_HPD_intervals_unscaled' : HPD intervals transformed back into calendar years 'connexe_HPD_intervals_unscaled_round' : integer-rounded calendar HPD intervals 'calage_posterior_mode' : posterior mode of the calibrated date (calendar scale) 'HPD_region_length' : total length (in years) of the HPD region 'middle_points' : calendar ages at which posterior density was evaluated 'middle_points_density' : posterior density evaluated at these points 'calage_posterior_mean' : posterior mean (if computed) 'calage_posterior_std' : posterior standard deviation (if computed) 'calage_sample' : posterior sample values (if sampling was performed) 'alpha' : significance level used 'c14age', 'c14sig' : the original radiocarbon measurement 'covariables' : always None for IntCal20 curve

Notes

Calibration model

The posterior distribution of the calendar date d is: $$ p(d | m) \propto p(d) \times L(m | d), $$ where:
- \(m\) is the measured radiocarbon age of std \(\sigma\) (in the F14C domain),
- \(L(m | d)\) is the measurement-likelihood (Gaussian error model in the F14C domain),
- \(p(d)\) is the prior (here uniform over the range covered by IntCal20 calibration curve).
When calibration is done with IntCal20, the likelihood term is estimated by $$ \hat{L}(m|d) = \frac{1}{\sqrt{2\pi(\sigma^2 + s_{IntCal20}(d)^2)}} \exp\left[-\frac{(m - \mu_{IntCal20}(d))^2}{2(\sigma^2 + s_{IntCal20}(d)^2)}\right] $$ where \( \mu_{IntCal20}(d) \) and \( s_{IntCal20}(d)\) are the mean and std published for the IntCal20 curve (they are interpolated for points where they are not given). This likelihood estimate comes from a Gaussian approximation of the curve obtained by Central Limit Theorem.

Approximation and sampling strategy of the posterior distribution

The interval covered by IntCal20 is divided into a large number of small intervals (typically one every ≈2 years), and on each interval the posterior density is approximated by piecewise-constant values. This allows:

• fast HPD region extraction,
• deriving a sampling strategy.

Further description about this can be found in the documentation of mono_cal_date_approx_density_sample for example.

HPD computation

HPD regions are obtained by sorting discretized posterior density values and finding the smallest set of intervals that accumulates posterior mass (1 - alpha).

Optional posterior moments

If compute_calage_posterior_mean_and_std=True:

• posterior samples are drawn using piecewise-constant inverse transform sampling,
• calendar ages are recovered by inverse scaling,
• sample mean and variance are computed.

Examples:

>>> out = IntCal20_calibration(10400, 18)
>>> out['connexe_HPD_intervals_unscaled_round']
array([[12101, 12124],
[12129, 12132],
[12140, 12157],
[12164, 12401],
[12424, 12478]])

>>> # Computing posterior mean and standard deviation
>>> out = IntCal20_calibration(10400, 18, compute_calage_posterior_mean_and_std=True)
>>> out['calage_posterior_mean']
12291
>>> out['calage_posterior_std']
105

bnn_for_14C_calibration.calibration.joint_calibration(c14ages, c14sigs, alpha=0.05, covariables=False, compute_calage_posterior_mean_and_std=True, compute_calage_posterior_mode=False, chaine_length=1000)

Perform joint Bayesian calibration of multiple radiocarbon dates using the Metropolis–Hastings within Gibbs sampler defined in multi_cal_date_approx_density_MCMC_sampler_for_concatenated_curve.

This function acts as a wrapper around the MCMC sampler: it converts measurements from conventional radiocarbon ages (\(^{14}\)C) to the corresponding F\(^{14}\)C values, runs the joint sampler, and computes posterior summaries (posterior mean, standard deviation, posterior mode).

Parameters:

Name	Type	Description	Default
c14ages	ndarray	Measured radiocarbon ages (conventional radiocarbon ages) in years BP (Before Present). Shape: (n_dates,).	required
c14sigs	ndarray	Laboratory measurement standard deviations associated to each conventional radiocarbon age. Shape: (n_dates,).	required
alpha	float	Tail probability for HPD regions (default: 0.05). Note: This argument is kept for compatibility, although joint calibration does not compute connexe HPD regions directly (these can be obtained via marginal calibration functions). Marginal credible intervals can be computed using the coordinates of the MCMC chain sampled from the joint density.	0.05
covariables	bool	Whether to use geophysical covariates (10Be, paleointensity) as additional inputs to the calibration BNNs.	False
compute_calage_posterior_mean_and_std	bool	If True, compute the posterior mean and posterior standard deviation of calibrated dates (calendar scale), based on the MCMC chain.	True
compute_calage_posterior_mode	bool	If True, compute the posterior mode as the MCMC sample with highest joint log-density.	False
chaine_length	int	Length of the MCMC chain (default: 1000).	1000

Returns:

Type

Description

Dict[str, Union[ndarray, float, int, None]]

A dictionary with the following keys: 'chaine' : np.ndarray MCMC chain of calibrated dates (calendar scale). Shape: (n_dates, chaine_length). 'chaine_log_joint_density_unscaled' : np.ndarray Log of the joint target density evaluated along the MCMC chain. Shape: (chaine_length,). 'acceptance_rate' : float Proportion of iterations where at least one coordinate changes (=modified_global_acceptance_rate returned by multi_cal_date_approx_density_MCMC_sampler_for_concatenated_curve). 'marginal_acceptance_rates' : np.ndarray Marginal acceptance rate for each date. Shape: (n_dates,). 'mode_index' : int or None Index of the chain with maximum joint posterior density. 'calage_posterior_mode' : np.ndarray or None Posterior mode for each calibrated date according to the joint density (vector of length n_dates) if requested. 'calage_posterior_mean' : np.ndarray or None Posterior mean of calibrated dates (length n_dates) if requested. 'calage_posterior_std' : np.ndarray or None Posterior standard deviation for each calibrated date (length n_dates) if requested.

Notes

1. Conventional radiocarbon ages conversion
Inputs are transformed from the ¹⁴C domain into the F¹⁴C domain prior to calibration using the deterministic conversions: $$ F^{14}C = e^{-\frac{\mathrm{C14Age}}{8033}} $$ and associated propagation for uncertainties is carried out using the delta method.

2. MCMC sampler
The sampler used internally performs joint calibration by evaluating the target density constructed from the concatenated curve (two pre-trained BNNs to estimate the two parts of the calibration curve).

3. Posterior summaries
The posterior summaries are calculated on the calendar scale, using direct moment estimators from the MCMC chain.