Skip to contents

Find Predictive Estimate of Shared Calendar Age Density from Bayesian Non-Parametric DPMM Output

Source: R/FindPredictiveCalendarAgeDensity.R
FindPredictiveCalendarAgeDensity.Rd

Given output from one of the Bayesian non-parametric summarisation functions (either PolyaUrnBivarDirichlet or WalkerBivarDirichlet) calculate the predictive (summarised/shared) calendar age density and probability intervals on a given calendar age grid (provided in cal yr BP).

Note: If you want to calculate and plot the result, use PlotPredictiveCalendarAgeDensity instead.

Usage

FindPredictiveCalendarAgeDensity(
  output_data,
  calendar_age_sequence,
  n_posterior_samples = 5000,
  interval_width = "2sigma",
  bespoke_probability = NA,
  n_burn = NA,
  n_end = NA
)

Arguments

output_data

The return value from one of the Bayesian non-parametric DPMM functions, e.g. PolyaUrnBivarDirichlet or WalkerBivarDirichlet, or a list, each item containing one of these return values. Optionally, the output data can have an extra list item named label which is used to set the label on the plot legend.

calendar_age_sequence

A vector containing the calendar age grid (in cal yr BP) on which to calculate the predictive (summarised/shared) density.

n_posterior_samples

Number of samples it will draw, after having removed n_burn, from the (thinned) realisations stored in the DPMM outputs to estimate the predictive calendar age density. These samples may be repeats if the number of, post burn-in, realisations is less than n_posterior_samples. If not given, 5000 is used.

interval_width

The confidence intervals to show for both the calibration curve and the predictive density. Choose from one of "1sigma" (68.3%), "2sigma" (95.4%) and "bespoke". Default is "2sigma".

bespoke_probability

The probability to use for the confidence interval if "bespoke" is chosen above. E.g., if 0.95 is chosen, then the 95% confidence interval is calculated. Ignored if "bespoke" is not chosen.

n_burn

The number of MCMC iterations that should be discarded as burn-in (i.e., considered to be occurring before the MCMC has converged). This relates to the number of iterations (n_iter) when running the original update functions (not the thinned output_data). Any MCMC iterations before this are not used in the calculations. If not given, the first half of the MCMC chain is discarded. Note: The maximum value that the function will allow is n_iter - 100 * n_thin (where n_iter and n_thin are the arguments given to PolyaUrnBivarDirichlet or WalkerBivarDirichlet) which would leave only 100 of the (thinned) values in output_data.

n_end

The last iteration in the original MCMC chain to use in the calculations. Assumed to be the total number of iterations performed, i.e. n_iter, if not given.

Value

A data frame of the calendar_age_BP, the density_mean and the confidence intervals for the density density_ci_lower and density_ci_upper.

See also

PlotPredictiveCalendarAgeDensity

Examples

# NOTE: All these examples are shown with a small n_iter and n_posterior_samples
# to speed up execution.
# Try n_iter and n_posterior_samples as the function defaults.

# First generate output data
polya_urn_output <- PolyaUrnBivarDirichlet(
    two_normals$c14_age,
    two_normals$c14_sig,
    intcal20,
    n_iter = 100,
    show_progress = FALSE)

# Find results for example output, 2-sigma confidence interval (default)
FindPredictiveCalendarAgeDensity(
    polya_urn_output, seq(3600, 4700, length=12), n_posterior_samples = 500)
#>    calendar_age_BP density_mean density_ci_lower density_ci_upper
#> 1             3600 7.555309e-04     5.927833e-04     8.815453e-04
#> 2             3700 5.966265e-04     3.593694e-04     7.484916e-04
#> 3             3800 3.831980e-04     1.741983e-04     5.140313e-04
#> 4             3900 2.014323e-04     6.768471e-05     2.777326e-04
#> 5             4000 8.757076e-05     2.129438e-05     1.189600e-04
#> 6             4100 3.204235e-05     5.677310e-06     4.973936e-05
#> 7             4200 1.028096e-05     1.558970e-06     1.781786e-05
#> 8             4300 3.289998e-06     6.904306e-07     5.737949e-06
#> 9             4400 1.520525e-06     5.295331e-07     2.026711e-06
#> 10            4500 2.767692e-06     5.785619e-07     7.625474e-06
#> 11            4600 2.172633e-05     4.368764e-06     5.706649e-05
#> 12            4700 1.601356e-04     6.924276e-05     2.736710e-04