Type:	Package
Title:	Distribution Fitting with Calibrating Priors for Commonly Used Distributions
Version:	0.1.1
Maintainer:	Stephen Jewson <stephen.jewson@gmail.com>
Imports:	stats, mev, extraDistr, gnorm, fdrtool, pracma, rust, actuar, fExtremes
Depends:	R (≥ 3.5.0)
Description:	Generates predictive distributions based on calibrating priors for various commonly used statistical models, including models with predictors. Routines for densities, probabilities, quantiles, random deviates and the parameter posterior are provided. The predictions are generated from the Bayesian prediction integral, with priors chosen to give good reliability (also known as calibration). For homogeneous models, the prior is set to the right Haar prior, giving predictions which are exactly reliable. As a result, in repeated testing, the frequencies of out-of-sample outcomes and the probabilities from the predictions agree. For other models, the prior is chosen to give good reliability. Where possible, the Bayesian prediction integral is solved exactly. Where exact solutions are not possible, the Bayesian prediction integral is solved using the Datta-Mukerjee-Ghosh-Sweeting (DMGS) asymptotic expansion. Optionally, the prediction integral can also be solved using posterior samples generated using Paul Northrop's ratio of uniforms sampling package ('rust'). Results are also generated based on maximum likelihood, for comparison purposes. Various model selection diagnostics and testing routines are included. Based on "Reducing reliability bias in assessments of extreme weather risk using calibrating priors", Jewson, S., Sweeting, T. and Jewson, L. (2024); <doi:10.5194/ascmo-11-1-2025>.
License:	MIT + file LICENSE
BugReports:	https://github.com/stephenjewson/fitdistcp/issues
URL:	https://fitdistcp.info
Encoding:	UTF-8
LazyData:	true
RoxygenNote:	7.3.1
Suggests:	knitr, rmarkdown
NeedsCompilation:	no
Packaged:	2025-04-19 11:12:16 UTC; steph
Author:	Stephen Jewson [aut, cre]
Repository:	CRAN
Date/Publication:	2025-04-23 09:40:02 UTC

Generates a comment about the method

Description

Generates a comment about the method

Usage

adhoc_dmgs_cpmethod()

Value

String

Generates a comment about the method

Description

Generates a comment about the method

Usage

analytic_cpmethod()

Value

String

Evaluate DMGS equation 3.3

Description

Evaluate DMGS equation 3.3

Usage

bayesian_dq_4terms_v1(lddi, lddd, mu1, pidopi1, pidopi2, mu2, dim)

Arguments

Value

Vector

determine revert2ml or not

Description

determine revert2ml or not

Usage

calc_revert2ml(v5h, v6h, t3)

Arguments

Value

Logical

Cauchy Distribution Predictions Based on a Calibrating Prior

Description

The fitdistcp package contains functions that generate predictive distributions for various statistical models, with and without parameter uncertainty. Parameter uncertainty is included by using Bayesian prediction with a type of objective prior known as a calibrating prior. Calibrating priors are chosen to give predictions that give good reliability (i.e., are well calibrated), for any underlying true parameter values.

There are five functions for each model, each of which uses training data x. For model **** the five functions are as follows:

• q****_cp returns predictive quantiles at the specified probabilities p, and various other diagnostics.

• r****_cp returns n random deviates from the predictive distribution.

• d****_cp returns the predictive density function at the specified values y

• p****_cp returns the predictive distribution function at the specified values y

• t****_cp returns n random deviates from the posterior distribution of the model parameters.

The q, r, d, p routines return two sets of results, one based on maximum likelihood, and the other based on a calibrating prior. The prior used depends on the model, and is given under Details below.

Where possible, the Bayesian prediction integral is solved analytically. Otherwise, DMGS asymptotic expansions are used. Optionally, a third set of results is returned that integrates the prediction integral by sampling the parameter posterior distribution using the RUST rejection sampling algorithm.

Usage

qcauchy_cp(
  x,
  p = seq(0.1, 0.9, 0.1),
  d1 = 0.01,
  fd2 = 0.01,
  means = FALSE,
  waicscores = FALSE,
  logscores = FALSE,
  dmgs = TRUE,
  rust = FALSE,
  nrust = 1e+05,
  debug = FALSE,
  aderivs = TRUE
)

rcauchy_cp(
  n,
  x,
  d1 = 0.01,
  fd2 = 0.01,
  rust = FALSE,
  mlcp = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

dcauchy_cp(
  x,
  y = x,
  d1 = 0.01,
  fd2 = 0.01,
  rust = FALSE,
  nrust = 1000,
  debug = FALSE,
  aderivs = TRUE
)

pcauchy_cp(
  x,
  y = x,
  d1 = 0.01,
  fd2 = 0.01,
  rust = FALSE,
  nrust = 1000,
  debug = FALSE,
  aderivs = TRUE
)

tcauchy_cp(n, x, d1 = 0.01, fd2 = 0.01, debug = FALSE)

Arguments

Value

q**** returns a list containing at least the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_value: the value of the log-likelihood at the maximum.

• standard_errors: estimates of the standard errors on the parameters, from the inverse observed information matrix.

• ml_quantiles: quantiles calculated using maximum likelihood.

• cp_quantiles: predictive quantiles calculated using a calibrating prior.

• maic: the AIC score for the maximum likelihood model, times -1/2.

• cp_method: a comment about the method used to generate the cp prediction.

For models with predictors, q**** additionally returns:

• predictedparameter: the estimated value for parameter, as a function of the predictor.

• adjustedx: the detrended values of x

r**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_deviates: random deviates calculated using maximum likelihood.

• cp_deviates: predictive random deviates calculated using a calibrating prior.

• cp_method: a comment about the method used to generate the cp prediction.

d**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_pdf: density function from maximum likelihood.

• cp_pdf: predictive density function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

p*** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_cdf: distribution function from maximum likelihood.

• cp_cdf: predictive distribution function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

t*** returns a list containing the following:

• theta_samples: random samples from the parameter posterior.

Details of the Model

The Cauchy distribution has probability density function

f(x;\mu,\sigma)=\frac{1}{\pi\sigma}\left(1+\left(\frac{x-\mu}{\sigma}\right)^2\right)^{-1}

where x is the random variable and \mu,\sigma>0 are the parameters.

The calibrating prior is given by the right Haar prior, which is

\pi(\sigma) \propto \frac{1}{\sigma}

as given in Jewson et al. (2025).

Optional Return Values

q**** optionally returns the following:

If rust=TRUE:

• ru_quantiles: predictive quantiles calculated using a calibrating prior, using posterior sampling with the RUST algorithm, based on inverting an empirical CDF based on nrust samples.

If waicscores=TRUE:

• waic1: the WAIC1 score for the calibrating prior model.

• waic2: the WAIC2 score for the calibrating prior model.

If logscores=TRUE:

• ml_oos_logscore: the leave-one-out logscore for the maximum likelihood prediction (not available in EVT routines, for mathematical reasons)

• cp_oos_logscore: the leave-one-out logscore for the parameter uncertainty model available in EVT routines, for mathematical reasons)

If means=TRUE:

• ml_mean: analytic estimate of the mean of the MLE predictive distribution, where possible

• cp_mean: analytic estimate of the mean of the calibrating prior predictive distribution, where mathematically possible. Can be compared with the mean estimated from random deviates.

r**** optionally returns the following:

If rust=TRUE:

• ru_deviates: nrust predictive random deviatives calculated using a calibrating prior, using posterior sampling with RUST.

d**** optionally returns the following:

If rust=TRUE:

• ru_pdf: predictive density calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust density functions.

p**** optionally returns the following:

If rust=TRUE:

• ru_cdf: predictive probability calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust distribution functions.

Selecting these additional outputs increases runtime. They are optional so that runtime for the basic outputs is minimised. This facilitates repeated experiments that evaluate reliability over many thousands of repeats.

Details (non-homogeneous models)

This model is not homogeneous, i.e. it does not have a transitive transformation group, and so there is no right Haar prior and no method for generating exactly reliable predictions. The cp outputs are generated using a prior that has been shown in tests to give reasonable reliability. See Jewson et al. (2024) for discussion of the prior and test results. For non-homogeneous models, reliability is generally poor for small sample sizes (<20), but is still much better than maximum likelihood. For small sample sizes, it is advisable to check the level of reliability using the routine reltest.

Details (DMGS integration)

For this model, the Bayesian prediction equation cannot be solved analytically, and is approximated using the DMGS asymptotic expansions given by Datta et al. (2000). This approximation seems to work well for medium and large sample sizes, but may not work well for small sample sizes (e.g., <20 data points). For small sample sizes, it may be preferable to use posterior sampling by setting rust=TRUE and looking at the ru outputs. The performance for any sample size, in terms of reliability, can be tested using reltest.

Details (RUST)

The Bayesian prediction equation can also be integrated using ratio-of-uniforms-sampling-with-transformation (RUST), using the option rust=TRUE. fitdistcp then calls Paul Northrop's rust package (Northrop, 2023). The RUST calculations are slower than the DMGS calculations.

For small sample sizes (e.g., n<20), and the very extreme tail, the DMGS approximation is somewhat poor (although always better than maximum likelihood) and it may be better to use RUST. For medium sample sizes (30+), DMGS is reasonably accurate, except for the very far tail.

It is advisable to check the RUST results for convergence versus the number of RUST samples.

It may be interesting to compare the DMGS and RUST results.

Author(s)

Stephen Jewson stephen.jewson@gmail.com

References

If you use this package, we would be grateful if you would cite the following reference, which gives the various calibrating priors, and tests them for reliability:

• Jewson S., Sweeting T. and Jewson L. (2024): Reducing Reliability Bias in Assessments of Extreme Weather Risk using Calibrating Priors; ASCMO Advances in Statistical Climatology, Meteorology and Oceanography), https://ascmo.copernicus.org/articles/11/1/2025/.

See Also

An introduction to fitdistcp, with more examples, is given on this webpage.

The fitdistcp package currently includes the following models (in alphabetical order):

• Cauchy (cauchy),

• Cauchy with linear predictor on the mean (cauchy_p1),

• Exponential (exp),

• Exponential with log-linear predictor on the scale (exp_p1),

• Frechet with known location parameter (frechet_k1),

• Frechet with log-linear predictor on the scale and known location parameter (frechet_p2k1),

• Gamma (gamma),

• Generalized normal (gnorm),

• GEV (gev),

• GEV with linear predictor on the location (gev_p1),

• GEV with linear predictor on the location and log-linear prediction on the scale (gev_p12),

• GEV with linear predictor on the location, log-linear prediction on the scale, and linear predictor on the shape (gev_p123),

• GEV with linear predictor on the location and known shape (gev_p1k3),

• GEV with known shape (gev_k3),

• GPD with known location (gpd_k1),

• Gumbel (gumbel),

• Gumbel with linear predictor on the mean(gumbel_p1),

• Half-normal (halfnorm),

• Inverse gamma (invgamma),

• Inverse Gaussian (invgauss),

• t distribution with unknown location and scale and known DoF (lst_k3),

• t distribution with unknown location and scale, linear predictor on the location, and known DoF (lst_p1k3),

• Logistic (logis),

• Logistic with linear predictor on the location (logis_p1),

• Log-normal (lnorm),

• Log-normal with linear predictor on the location (lnorm_p1),

• Normal (norm),

• Normal with linear predictor on the mean (norm_p1),

• Pareto with known scale (pareto_k2),

• Pareto with log-linear predictor on the shape and known scale (pareto_p1k2),

• Uniform (unif),

• Weibull (weibull),

• Weibull with linear predictor on the scale (weibull_p2),

The level of predictive probability matching achieved by the maximum likelihood and calibrating prior quantiles, for any model, sample size and true parameter values, can be demonstrated using the routine reltest.

Model selection among models can be demonstrated using the routines ms_flat_1tail, ms_flat_2tail, ms_predictors_1tail, and ms_predictors_2tail,

Examples

#
# example 1
x=fitdistcp::d42cauchy_example_data_v1
p=c(1:9)/10
q=qlogis_cp(x,p,rust=TRUE,nrust=1000)
xmin=min(q$ml_quantiles,q$cp_quantiles);
xmax=max(q$ml_quantiles,q$cp_quantiles);
plot(q$ml_quantiles,p,xlab="quantile estimates",xlim=c(xmin,xmax),
	sub="(from qcauchy_cp)",
	main="Cauchy: quantile estimates");
points(q$cp_quantiles,p,col="red",lwd=2)
points(q$ru_quantiles,p,col="blue")

DMGS equation 3.3, f1 term

Description

DMGS equation 3.3, f1 term

Usage

cauchy_f1f(y, v1, d1, v2, fd2)

Arguments

Value

Matrix

The first derivative of the density

Description

The first derivative of the density

Usage

cauchy_f1fa(x, v1, v2)

Arguments

Value

Vector

DMGS equation 3.3, f2 term

Description

DMGS equation 3.3, f2 term

Usage

cauchy_f2f(y, v1, d1, v2, fd2)

Arguments

Value

3d array

The second derivative of the density

Description

The second derivative of the density

Usage

cauchy_f2fa(x, v1, v2)

Arguments

Value

Matrix

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_fd(x, v1, v2)

Arguments

Value

Vector

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_fdd(x, v1, v2)

Arguments

Value

Matrix

Second derivative matrix of the normalized log-likelihood

Description

Second derivative matrix of the normalized log-likelihood

Usage

cauchy_ldd(x, v1, d1, v2, fd2)

Arguments

Value

Square scalar matrix

The second derivative of the normalized log-likelihood

Description

The second derivative of the normalized log-likelihood

Usage

cauchy_ldda(x, v1, v2)

Arguments

Value

Matrix

Third derivative tensor of the normalized log-likelihood

Description

Third derivative tensor of the normalized log-likelihood

Usage

cauchy_lddd(x, v1, d1, v2, fd2)

Arguments

Value

Cubic scalar array

The third derivative of the normalized log-likelihood

Description

The third derivative of the normalized log-likelihood

Usage

cauchy_lddda(x, v1, v2)

Arguments

Value

3d array

One component of the second derivative of the normalized log-likelihood

Description

One component of the second derivative of the normalized log-likelihood

Usage

cauchy_lmn(x, v1, d1, v2, fd2, mm, nn)

Arguments

Value

Scalar value

One component of the third derivative of the normalized log-likelihood

Description

One component of the third derivative of the normalized log-likelihood

Usage

cauchy_lmnp(x, v1, d1, v2, fd2, mm, nn, rr)

Arguments

Value

Scalar value

Logf for RUST

Description

Logf for RUST

Usage

cauchy_logf(params, x)

Arguments

Value

Scalar value.

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_logfdd(x, v1, v2)

Arguments

Value

Matrix

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_logfddd(x, v1, v2)

Arguments

Value

3d array

log-likelihood function

Description

log-likelihood function

Usage

cauchy_loglik(vv, x)

Arguments

Value

Scalar value.

Log scores for MLE and RHP predictions calculated using leave-one-out

Description

Log scores for MLE and RHP predictions calculated using leave-one-out

Usage

cauchy_logscores(logscores, x, d1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

Two scalars

DMGS equation 3.3, mu1 term

Description

DMGS equation 3.3, mu1 term

Usage

cauchy_mu1f(alpha, v1, d1, v2, fd2)

Arguments

Value

Matrix

DMGS equation 3.3, mu2 term

Description

DMGS equation 3.3, mu2 term

Usage

cauchy_mu2f(alpha, v1, d1, v2, fd2)

Arguments

Value

3d array

Cauchy Distribution with a Predictor, Predictions Based on a Calibrating Prior

Description

The fitdistcp package contains functions that generate predictive distributions for various statistical models, with and without parameter uncertainty. Parameter uncertainty is included by using Bayesian prediction with a type of objective prior known as a calibrating prior. Calibrating priors are chosen to give predictions that give good reliability (i.e., are well calibrated), for any underlying true parameter values.

There are five functions for each model, each of which uses training data x. For model **** the five functions are as follows:

• q****_cp returns predictive quantiles at the specified probabilities p, and various other diagnostics.

• r****_cp returns n random deviates from the predictive distribution.

• d****_cp returns the predictive density function at the specified values y

• p****_cp returns the predictive distribution function at the specified values y

• t****_cp returns n random deviates from the posterior distribution of the model parameters.

The q, r, d, p routines return two sets of results, one based on maximum likelihood, and the other based on a calibrating prior. The prior used depends on the model, and is given under Details below.

Where possible, the Bayesian prediction integral is solved analytically. Otherwise, DMGS asymptotic expansions are used. Optionally, a third set of results is returned that integrates the prediction integral by sampling the parameter posterior distribution using the RUST rejection sampling algorithm.

Usage

qcauchy_p1_cp(
  x,
  t,
  t0 = NA,
  n0 = NA,
  p = seq(0.1, 0.9, 0.1),
  d1 = 0.01,
  d2 = 0.01,
  fd3 = 0.01,
  means = FALSE,
  waicscores = FALSE,
  logscores = FALSE,
  dmgs = TRUE,
  rust = FALSE,
  nrust = 1e+05,
  predictordata = TRUE,
  centering = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

rcauchy_p1_cp(
  n,
  x,
  t,
  t0 = NA,
  n0 = NA,
  d1 = 0.01,
  d2 = 0.01,
  fd3 = 0.01,
  rust = FALSE,
  mlcp = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

dcauchy_p1_cp(
  x,
  t,
  t0 = NA,
  n0 = NA,
  y = x,
  d1 = 0.01,
  d2 = 0.01,
  fd3 = 0.01,
  rust = FALSE,
  nrust = 1000,
  centering = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

pcauchy_p1_cp(
  x,
  t,
  t0 = NA,
  n0 = NA,
  y = x,
  d1 = 0.01,
  d2 = 0.01,
  fd3 = 0.01,
  rust = FALSE,
  nrust = 1000,
  centering = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

tcauchy_p1_cp(n, x, t, d1 = 0.01, d2 = 0.01, fd3 = 0.01, debug = FALSE)

Arguments

Value

q**** returns a list containing at least the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_value: the value of the log-likelihood at the maximum.

• standard_errors: estimates of the standard errors on the parameters, from the inverse observed information matrix.

• ml_quantiles: quantiles calculated using maximum likelihood.

• cp_quantiles: predictive quantiles calculated using a calibrating prior.

• maic: the AIC score for the maximum likelihood model, times -1/2.

• cp_method: a comment about the method used to generate the cp prediction.

For models with predictors, q**** additionally returns:

• predictedparameter: the estimated value for parameter, as a function of the predictor.

• adjustedx: the detrended values of x

r**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_deviates: random deviates calculated using maximum likelihood.

• cp_deviates: predictive random deviates calculated using a calibrating prior.

• cp_method: a comment about the method used to generate the cp prediction.

d**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_pdf: density function from maximum likelihood.

• cp_pdf: predictive density function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

p*** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_cdf: distribution function from maximum likelihood.

• cp_cdf: predictive distribution function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

t*** returns a list containing the following:

• theta_samples: random samples from the parameter posterior.

Details of the Model

The Cauchy distribution with a predictor has probability density function

f(x;a,b,\sigma)=\frac{1}{\pi\sigma}\left(1+\left(\frac{x-\mu(a,b)}{\sigma}\right)^2\right)^{-1}

where x is the random variable, \mu=a+bt is the location parameter as a function of parameters a,b, and \sigma>0 is the scale parameter.

The calibrating prior is given by the right Haar prior, which is

\pi(\sigma) \propto \frac{1}{\sigma}

as given in Jewson et al. (2025).

Optional Return Values

q**** optionally returns the following:

If rust=TRUE:

• ru_quantiles: predictive quantiles calculated using a calibrating prior, using posterior sampling with the RUST algorithm, based on inverting an empirical CDF based on nrust samples.

If waicscores=TRUE:

• waic1: the WAIC1 score for the calibrating prior model.

• waic2: the WAIC2 score for the calibrating prior model.

If logscores=TRUE:

• ml_oos_logscore: the leave-one-out logscore for the maximum likelihood prediction (not available in EVT routines, for mathematical reasons)

• cp_oos_logscore: the leave-one-out logscore for the parameter uncertainty model available in EVT routines, for mathematical reasons)

If means=TRUE:

• ml_mean: analytic estimate of the mean of the MLE predictive distribution, where possible

• cp_mean: analytic estimate of the mean of the calibrating prior predictive distribution, where mathematically possible. Can be compared with the mean estimated from random deviates.

r**** optionally returns the following:

If rust=TRUE:

• ru_deviates: nrust predictive random deviatives calculated using a calibrating prior, using posterior sampling with RUST.

d**** optionally returns the following:

If rust=TRUE:

• ru_pdf: predictive density calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust density functions.

p**** optionally returns the following:

If rust=TRUE:

• ru_cdf: predictive probability calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust distribution functions.

Selecting these additional outputs increases runtime. They are optional so that runtime for the basic outputs is minimised. This facilitates repeated experiments that evaluate reliability over many thousands of repeats.

Details (homogeneous models)

This model is a homogeneous model, and the cp results are based on the right Haar prior. For homogeneous models (models with sharply transitive transformation groups), a Bayesian prediction based on the right Haar prior gives exact reliability, as shown by Severini et al. (2002), even when the true parameters are unknown. This means that probabilities in the prediction will correspond to frequencies of future outcomes in repeated trials (if the model is correct).

Maximum likelihood prediction does not give reliable predictions, even when the model is correct, because it does not account for parameter uncertainty. In particular, maximum likelihood predictions typically underestimate the tail in repeated trials.

The reliability of the maximum likelihood and the calibrating prior predictive quantiles produced by the q****_cp routines in fitdistcp can be quantified using repeated simulations with the routine reltest.

Details (DMGS integration)

For this model, the Bayesian prediction equation cannot be solved analytically, and is approximated using the DMGS asymptotic expansions given by Datta et al. (2000). This approximation seems to work well for medium and large sample sizes, but may not work well for small sample sizes (e.g., <20 data points). For small sample sizes, it may be preferable to use posterior sampling by setting rust=TRUE and looking at the ru outputs. The performance for any sample size, in terms of reliability, can be tested using reltest.

Details (RUST)

The Bayesian prediction equation can also be integrated using ratio-of-uniforms-sampling-with-transformation (RUST), using the option rust=TRUE. fitdistcp then calls Paul Northrop's rust package (Northrop, 2023). The RUST calculations are slower than the DMGS calculations.

For small sample sizes (e.g., n<20), and the very extreme tail, the DMGS approximation is somewhat poor (although always better than maximum likelihood) and it may be better to use RUST. For medium sample sizes (30+), DMGS is reasonably accurate, except for the very far tail.

It is advisable to check the RUST results for convergence versus the number of RUST samples.

It may be interesting to compare the DMGS and RUST results.

Author(s)

Stephen Jewson stephen.jewson@gmail.com

References

If you use this package, we would be grateful if you would cite the following reference, which gives the various calibrating priors, and tests them for reliability:

• Jewson S., Sweeting T. and Jewson L. (2024): Reducing Reliability Bias in Assessments of Extreme Weather Risk using Calibrating Priors; ASCMO Advances in Statistical Climatology, Meteorology and Oceanography), https://ascmo.copernicus.org/articles/11/1/2025/.

See Also

An introduction to fitdistcp, with more examples, is given on this webpage.

The fitdistcp package currently includes the following models (in alphabetical order):

• Cauchy (cauchy),

• Cauchy with linear predictor on the mean (cauchy_p1),

• Exponential (exp),

• Exponential with log-linear predictor on the scale (exp_p1),

• Frechet with known location parameter (frechet_k1),

• Frechet with log-linear predictor on the scale and known location parameter (frechet_p2k1),

• Gamma (gamma),

• Generalized normal (gnorm),

• GEV (gev),

• GEV with linear predictor on the location (gev_p1),

• GEV with linear predictor on the location and log-linear prediction on the scale (gev_p12),

• GEV with linear predictor on the location, log-linear prediction on the scale, and linear predictor on the shape (gev_p123),

• GEV with linear predictor on the location and known shape (gev_p1k3),

• GEV with known shape (gev_k3),

• GPD with known location (gpd_k1),

• Gumbel (gumbel),

• Gumbel with linear predictor on the mean(gumbel_p1),

• Half-normal (halfnorm),

• Inverse gamma (invgamma),

• Inverse Gaussian (invgauss),

• t distribution with unknown location and scale and known DoF (lst_k3),

• t distribution with unknown location and scale, linear predictor on the location, and known DoF (lst_p1k3),

• Logistic (logis),

• Logistic with linear predictor on the location (logis_p1),

• Log-normal (lnorm),

• Log-normal with linear predictor on the location (lnorm_p1),

• Normal (norm),

• Normal with linear predictor on the mean (norm_p1),

• Pareto with known scale (pareto_k2),

• Pareto with log-linear predictor on the shape and known scale (pareto_p1k2),

• Uniform (unif),

• Weibull (weibull),

• Weibull with linear predictor on the scale (weibull_p2),

The level of predictive probability matching achieved by the maximum likelihood and calibrating prior quantiles, for any model, sample size and true parameter values, can be demonstrated using the routine reltest.

Model selection among models can be demonstrated using the routines ms_flat_1tail, ms_flat_2tail, ms_predictors_1tail, and ms_predictors_2tail,

Examples

#
# example 1
x=fitdistcp::d64cauchy_p1_example_data_v1_x
tt=fitdistcp::d64cauchy_p1_example_data_v1_t
p=c(1:9)/10
n0=10
q=qcauchy_p1_cp(x,tt,n0=n0,p=p,rust=TRUE,nrust=1000)
xmin=min(q$ml_quantiles,q$cp_quantiles);
xmax=max(q$ml_quantiles,q$cp_quantiles);
plot(q$ml_quantiles,p,xlab="quantile estimates",xlim=c(xmin,xmax),
	sub="(from qcauchy_p1_cp)",
	main="Cauchy w/ p1: quantile estimates");
points(q$cp_quantiles,p,col="red",lwd=2)
points(q$ru_quantiles,p,col="blue")

DMGS equation 2.1, f1 term

Description

DMGS equation 2.1, f1 term

Usage

cauchy_p1_f1f(y, t0, v1, d1, v2, d2, v3, fd3)

Arguments

Value

Matrix

The first derivative of the density

Description

The first derivative of the density

Usage

cauchy_p1_f1fa(x, t, v1, v2, v3)

Arguments

Value

Vector

DMGS equation 2.1, f2 term

Description

DMGS equation 2.1, f2 term

Usage

cauchy_p1_f2f(y, t0, v1, d1, v2, d2, v3, fd3)

Arguments

Value

3d array

The second derivative of the density

Description

The second derivative of the density

Usage

cauchy_p1_f2fa(x, t, v1, v2, v3)

Arguments

Value

Matrix

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_p1_fd(x, t, v1, v2, v3)

Arguments

Value

Vector

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_p1_fdd(x, t, v1, v2, v3)

Arguments

Value

Matrix

Matrix

Second derivative matrix of the normalized log-likelihood

Description

Second derivative matrix of the normalized log-likelihood

Usage

cauchy_p1_ldd(x, t, v1, d1, v2, d2, v3, fd3)

Arguments

Value

Square scalar matrix

The second derivative of the normalized log-likelihood

Description

The second derivative of the normalized log-likelihood

Usage

cauchy_p1_ldda(x, t, v1, v2, v3)

Arguments

Value

Matrix

Third derivative tensor of the normalized log-likelihood

Description

Third derivative tensor of the normalized log-likelihood

Usage

cauchy_p1_lddd(x, t, v1, d1, v2, d2, v3, fd3)

Arguments

Value

Cubic scalar array

The third derivative of the normalized log-likelihood

Description

The third derivative of the normalized log-likelihood

Usage

cauchy_p1_lddda(x, t, v1, v2, v3)

Arguments

Value

3d array

One component of the second derivative of the normalized log-likelihood

Description

One component of the second derivative of the normalized log-likelihood

Usage

cauchy_p1_lmn(x, t, v1, d1, v2, d2, v3, fd3, mm, nn)

Arguments

Value

Scalar value

One component of the second derivative of the normalized log-likelihood

Description

One component of the second derivative of the normalized log-likelihood

Usage

cauchy_p1_lmnp(x, t, v1, d1, v2, d2, v3, fd3, mm, nn, rr)

Arguments

Value

Scalar value

Logf for RUST

Description

Logf for RUST

Usage

cauchy_p1_logf(params, x, t)

Arguments

Value

Scalar value.

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_p1_logfdd(x, t, v1, v2, v3)

Arguments

Value

Matrix

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

cauchy_p1_logfddd(x, t, v1, v2, v3)

Arguments

Value

3d array

Cauchy-with-p1 observed log-likelihood function

Description

Cauchy-with-p1 observed log-likelihood function

Usage

cauchy_p1_loglik(vv, x, t)

Arguments

Value

Scalar value.

Log scores for MLE and RHP predictions calculated using leave-one-out

Description

Log scores for MLE and RHP predictions calculated using leave-one-out

Usage

cauchy_p1_logscores(logscores, x, t, d1, d2, fd3, aderivs = TRUE)

Arguments

Value

Two scalars

Cauchy distribution: RHP mean

Description

Cauchy distribution: RHP mean

Usage

cauchy_p1_means(t0, ml_params, lddi, lddd, lambdad_rhp, nx, dim = 2)

Arguments

Value

Two scalars

DMGS equation 3.3, mu1 term

Description

DMGS equation 3.3, mu1 term

Usage

cauchy_p1_mu1f(alpha, t0, v1, d1, v2, d2, v3, fd3)

Arguments

Value

Matrix

DMGS equation 3.3, mu2 term

Description

DMGS equation 3.3, mu2 term

Usage

cauchy_p1_mu2f(alpha, t0, v1, d1, v2, d2, v3, fd3)

Arguments

Value

3d array

DMGS equation 2.1, p1 term

Description

DMGS equation 2.1, p1 term

Usage

cauchy_p1_p1f(y, t0, v1, d1, v2, d2, v3, fd3)

Arguments

Value

Matrix

DMGS equation 2.1, p2 term

Description

DMGS equation 2.1, p2 term

Usage

cauchy_p1_p2f(y, t0, v1, d1, v2, d2, v3, fd3)

Arguments

Value

3d array

Predicted Parameter and Generalized Residuals

Description

Predicted Parameter and Generalized Residuals

Usage

cauchy_p1_predictordata(predictordata, x, t, t0, params)

Arguments

Value

Two vectors

Waic

Description

Waic

Usage

cauchy_p1_waic(
  waicscores,
  x,
  t,
  v1hat,
  d1,
  v2hat,
  d2,
  v3hat,
  fd3,
  lddi,
  lddd,
  lambdad,
  aderivs
)

Arguments

Value

Two numeric values.

DMGS equation 3.3, p1 term

Description

DMGS equation 3.3, p1 term

Usage

cauchy_p1f(y, v1, d1, v2, fd2)

Arguments

Value

Matrix

DMGS equation 3.3, p2 term

Description

DMGS equation 3.3, p2 term

Usage

cauchy_p2f(y, v1, d1, v2, fd2)

Arguments

Value

3d array

Waic

Description

Waic

Usage

cauchy_waic(waicscores, x, v1hat, d1, v2hat, fd2, lddi, lddd, lambdad, aderivs)

Arguments

Value

Two numeric values.

Generates a comment about the method

Description

Generates a comment about the method

Usage

crhpflat_dmgs_cpmethod()

Value

String

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

This is data to be included in my package

Description

This is data to be included in my package

Cauchy-with-p1 density function

Description

Cauchy-with-p1 density function

Usage

dcauchy_p1(x, t0, ymn, slope, scale, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dcauchy_p1sub(x, t, y, t0, d1, d2, fd3, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dcauchysub(x, y, d1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Extract the results from derivatives and put them into f2

Description

Extract the results from derivatives and put them into f2

Usage

deriv_copyfdd(temp1, nx, dim)

Arguments

Value

3d array

Extract the results from derivatives and put them into ldd

Description

Extract the results from derivatives and put them into ldd

Usage

deriv_copyld2(temp1, nx, dim)

Arguments

Value

3d array

Extract the results from derivatives and put them into ldd

Description

Extract the results from derivatives and put them into ldd

Usage

deriv_copyldd(temp1, nx, dim)

Arguments

Value

Matrix

Extract the results from derivatives and put them into lddd

Description

Extract the results from derivatives and put them into lddd

Usage

deriv_copylddd(temp1, nx, dim)

Arguments

Value

3d array

Exponential-with-p1 density function

Description

Exponential-with-p1 density function

Usage

dexp_p1(x, t0, ymn, slope, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dexp_p1sub(x, t, y, t0, d1, d2, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dexpsub(x, y, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Frechet_k1-with-p2 density function

Description

Frechet_k1-with-p2 density function

Usage

dfrechet_p2k1(x, t0, ymn, slope, lambda, log = FALSE, kloc)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dfrechet_p2k1sub(x, t, y, t0, d1, d2, fd3, kloc, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dfrechetsub(x, y, kloc, fd1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dgammasub(x, y, fd1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dgev_k3sub(x, y, d1 = 0.01, fd2 = 0.01, kshape, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

GEVD-with-p1: Density function

Description

GEVD-with-p1: Density function

Usage

dgev_p1(x, t0, ymn, slope, sigma, xi, log = FALSE)

Arguments

Value

Vector

GEVD-with-p1: Density function

Description

GEVD-with-p1: Density function

Usage

dgev_p12(x, t1, t2, ymn, slope, sigma1, sigma2, xi, log = FALSE)

Arguments

Value

Vector

GEVD-with-p1: Density function

Description

GEVD-with-p1: Density function

Usage

dgev_p123(x, t1, t2, t3, ymn, slope, sigma1, sigma2, xi1, xi2, log = FALSE)

Arguments

Value

Vector

Densities for 5 predictions

Description

Densities for 5 predictions

Usage

dgev_p123sub(
  x,
  t1,
  t2,
  t3,
  y,
  t01,
  t02,
  t03,
  ics,
  d1 = 0.01,
  d2 = 0.01,
  d3 = 0.01,
  d4 = 0.01,
  d5 = 0.01,
  d6 = 0.01,
  extramodels,
  debug,
  aderivs = TRUE
)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities for 5 predictions

Description

Densities for 5 predictions

Usage

dgev_p12sub(
  x,
  t1,
  t2,
  y,
  t01,
  t02,
  ics,
  d1 = 0.01,
  d2 = 0.01,
  d3 = 0.01,
  d4 = 0.01,
  d5 = 0.01,
  minxi,
  maxxi,
  debug,
  extramodels = FALSE,
  aderivs = TRUE
)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

GEV-with-known-shape-with-p1 density function

Description

GEV-with-known-shape-with-p1 density function

Usage

dgev_p1k3(x, t0, ymn, slope, sigma, log = FALSE, kshape)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dgev_p1k3sub(x, t, y, t0, d1, d2, fd3, kshape, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities for 5 predictions

Description

Densities for 5 predictions

Usage

dgev_p1sub(
  x,
  t,
  y,
  t0,
  ics,
  d1 = 0.01,
  d2 = 0.01,
  fd3 = 0.01,
  d4 = 0.01,
  minxi,
  maxxi,
  extramodels = FALSE,
  aderivs = TRUE
)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities for 5 predictions

Description

Densities for 5 predictions

Usage

dgevsub(
  x,
  y,
  ics,
  d1 = 0.01,
  fd2 = 0.01,
  d3 = 0.01,
  customprior,
  minxi,
  maxxi,
  extramodels = FALSE,
  aderivs = TRUE
)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dgnorm_k3sub(x, y, d1 = 0.01, fd2 = 0.01, kbeta, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities for 5 predictions

Description

Densities for 5 predictions

Usage

dgpdsub(
  x,
  y,
  ics,
  fd1 = 0.01,
  d2 = 0.01,
  kloc = 0,
  dlogpi = 0,
  minxi,
  maxxi,
  extramodels = FALSE,
  aderivs = TRUE
)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Gumbel-with-p1 density function

Description

Gumbel-with-p1 density function

Usage

dgumbel_p1(x, t0, ymn, slope, sigma, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dgumbel_p1sub(x, t, y, t0, d1, d2, fd3, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dgumbelsub(x, y, d1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dhalfnormsub(x, y, fd1 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and cp

Description

Densities from MLE and cp

Usage

dinvgammasub(x, y, fd1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dinvgausssub(x, y, prior, fd1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlnorm_dmgssub(x, y, d1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Normal-with-p1 density function

Description

Normal-with-p1 density function

Usage

dlnorm_p1(x, t0, ymn, slope, sigma, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlnorm_p1sub(x, t, y, t0, debug = FALSE, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlnormsub(x, y, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlogis2sub(x, y, d1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Logistic-with-p1 density function

Description

Logistic-with-p1 density function

Usage

dlogis_p1(x, t0, ymn, slope, scale, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlogis_p1sub(x, t, y, t0, d1, d2, fd3, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlst_k3sub(x, y, d1 = 0.01, fd2 = 0.01, kdf, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

LST-with-p1 density function

Description

LST-with-p1 density function

Usage

dlst_p1k3(x, t0, ymn, slope, sigma, log = FALSE, kdf)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dlst_p1k3sub(x, t, y, t0, d1, d2, fd3, kdf, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Evaluate DMGS equation 3.3

Description

Evaluate DMGS equation 3.3

Usage

dmgs(lddi, lddd, mu1, pidopi, mu2, dim)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dnorm_dmgssub(x, y, d1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Normal-with-p1 density function

Description

Normal-with-p1 density function

Usage

dnorm_p1(x, t0, ymn, slope, sigma, log = FALSE)

Arguments

Value

Vector

Linear regression formula, densities

Description

Linear regression formula, densities

Usage

dnorm_p1_formula(y, ta, ta0, nx, muhat0, v3hat)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dnorm_p1sub(x, t, y, t0, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dnormsub(x, y, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dpareto_k2_sub(x, y, kscale, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

pareto_k1-with-p2 density function

Description

pareto_k1-with-p2 density function

Usage

dpareto_p1k2(x, t0, ymn, slope, kscale, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dpareto_p1k2sub(x, t, y, t0, d1, d2, kscale, aderivs = TRUE, debug = FALSE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Predictive PDFs

Description

Predictive PDFs

Usage

dunif_formula(x, y)

Arguments

Value

Two vectors

Weibull-with-p1 density function

Description

Weibull-with-p1 density function

Usage

dweibull_p2(x, t0, shape, ymn, slope, log = FALSE)

Arguments

Value

Vector

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dweibull_p2sub(x, t, y, t0, fd1, d2, d3, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Densities from MLE and RHP

Description

Densities from MLE and RHP

Usage

dweibullsub(x, y, fd1 = 0.01, fd2 = 0.01, aderivs = TRUE)

Arguments

Value

A vector of parameter estimates, two pdf vectors, two cdf vectors

Exponential Distribution Predictions Based on a Calibrating Prior

Description

The fitdistcp package contains functions that generate predictive distributions for various statistical models, with and without parameter uncertainty. Parameter uncertainty is included by using Bayesian prediction with a type of objective prior known as a calibrating prior. Calibrating priors are chosen to give predictions that give good reliability (i.e., are well calibrated), for any underlying true parameter values.

There are five functions for each model, each of which uses training data x. For model **** the five functions are as follows:

• q****_cp returns predictive quantiles at the specified probabilities p, and various other diagnostics.

• r****_cp returns n random deviates from the predictive distribution.

• d****_cp returns the predictive density function at the specified values y

• p****_cp returns the predictive distribution function at the specified values y

• t****_cp returns n random deviates from the posterior distribution of the model parameters.

The q, r, d, p routines return two sets of results, one based on maximum likelihood, and the other based on a calibrating prior. The prior used depends on the model, and is given under Details below.

Where possible, the Bayesian prediction integral is solved analytically. Otherwise, DMGS asymptotic expansions are used. Optionally, a third set of results is returned that integrates the prediction integral by sampling the parameter posterior distribution using the RUST rejection sampling algorithm.

Usage

qexp_cp(
  x,
  p = seq(0.1, 0.9, 0.1),
  fd1 = 0.01,
  means = FALSE,
  waicscores = FALSE,
  logscores = FALSE,
  rust = FALSE,
  nrust = 1e+05,
  debug = FALSE,
  aderivs = TRUE
)

rexp_cp(n, x, rust = FALSE, mlcp = TRUE, debug = FALSE, aderivs = TRUE)

dexp_cp(x, y = x, rust = FALSE, nrust = 1000, debug = FALSE, aderivs = TRUE)

pexp_cp(x, y = x, rust = FALSE, nrust = 1000, debug = FALSE, aderivs = TRUE)

texp_cp(n, x, debug = FALSE)

Arguments

Value

q**** returns a list containing at least the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_value: the value of the log-likelihood at the maximum.

• standard_errors: estimates of the standard errors on the parameters, from the inverse observed information matrix.

• ml_quantiles: quantiles calculated using maximum likelihood.

• cp_quantiles: predictive quantiles calculated using a calibrating prior.

• maic: the AIC score for the maximum likelihood model, times -1/2.

• cp_method: a comment about the method used to generate the cp prediction.

For models with predictors, q**** additionally returns:

• predictedparameter: the estimated value for parameter, as a function of the predictor.

• adjustedx: the detrended values of x

r**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_deviates: random deviates calculated using maximum likelihood.

• cp_deviates: predictive random deviates calculated using a calibrating prior.

• cp_method: a comment about the method used to generate the cp prediction.

d**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_pdf: density function from maximum likelihood.

• cp_pdf: predictive density function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

p*** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_cdf: distribution function from maximum likelihood.

• cp_cdf: predictive distribution function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

t*** returns a list containing the following:

• theta_samples: random samples from the parameter posterior.

Details of the Model

The exponential distribution has exceedance distribution function

S(x;\lambda)=\exp(-\lambda x)

where x \ge 0 is the random variable and \lambda >0 is the rate parameter.

The calibrating prior is given by the right Haar prior, which is

\pi(\lambda) \propto \frac{1}{\lambda}

as given in Jewson et al. (2025).

Optional Return Values

q**** optionally returns the following:

If rust=TRUE:

• ru_quantiles: predictive quantiles calculated using a calibrating prior, using posterior sampling with the RUST algorithm, based on inverting an empirical CDF based on nrust samples.

If waicscores=TRUE:

• waic1: the WAIC1 score for the calibrating prior model.

• waic2: the WAIC2 score for the calibrating prior model.

If logscores=TRUE:

• ml_oos_logscore: the leave-one-out logscore for the maximum likelihood prediction (not available in EVT routines, for mathematical reasons)

• cp_oos_logscore: the leave-one-out logscore for the parameter uncertainty model available in EVT routines, for mathematical reasons)

If means=TRUE:

• ml_mean: analytic estimate of the mean of the MLE predictive distribution, where possible

• cp_mean: analytic estimate of the mean of the calibrating prior predictive distribution, where mathematically possible. Can be compared with the mean estimated from random deviates.

r**** optionally returns the following:

If rust=TRUE:

• ru_deviates: nrust predictive random deviatives calculated using a calibrating prior, using posterior sampling with RUST.

d**** optionally returns the following:

If rust=TRUE:

• ru_pdf: predictive density calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust density functions.

p**** optionally returns the following:

If rust=TRUE:

• ru_cdf: predictive probability calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust distribution functions.

Selecting these additional outputs increases runtime. They are optional so that runtime for the basic outputs is minimised. This facilitates repeated experiments that evaluate reliability over many thousands of repeats.

Details (homogeneous models)

This model is a homogeneous model, and the cp results are based on the right Haar prior. For homogeneous models (models with sharply transitive transformation groups), a Bayesian prediction based on the right Haar prior gives exact reliability, as shown by Severini et al. (2002), even when the true parameters are unknown. This means that probabilities in the prediction will correspond to frequencies of future outcomes in repeated trials (if the model is correct).

Maximum likelihood prediction does not give reliable predictions, even when the model is correct, because it does not account for parameter uncertainty. In particular, maximum likelihood predictions typically underestimate the tail in repeated trials.

The reliability of the maximum likelihood and the calibrating prior predictive quantiles produced by the q****_cp routines in fitdistcp can be quantified using repeated simulations with the routine reltest.

Details (analytic integration)

For this model, the Bayesian prediction equation is integrated analytically.

Details (RUST)

The Bayesian prediction equation can also be integrated using ratio-of-uniforms-sampling-with-transformation (RUST), using the option rust=TRUE. fitdistcp then calls Paul Northrop's rust package (Northrop, 2023). The RUST calculations are slower than the DMGS calculations.

For small sample sizes (e.g., n<20), and the very extreme tail, the DMGS approximation is somewhat poor (although always better than maximum likelihood) and it may be better to use RUST. For medium sample sizes (30+), DMGS is reasonably accurate, except for the very far tail.

It is advisable to check the RUST results for convergence versus the number of RUST samples.

It may be interesting to compare the DMGS and RUST results.

Author(s)

Stephen Jewson stephen.jewson@gmail.com

References

If you use this package, we would be grateful if you would cite the following reference, which gives the various calibrating priors, and tests them for reliability:

• Jewson S., Sweeting T. and Jewson L. (2024): Reducing Reliability Bias in Assessments of Extreme Weather Risk using Calibrating Priors; ASCMO Advances in Statistical Climatology, Meteorology and Oceanography), https://ascmo.copernicus.org/articles/11/1/2025/.

See Also

An introduction to fitdistcp, with more examples, is given on this webpage.

The fitdistcp package currently includes the following models (in alphabetical order):

• Cauchy (cauchy),

• Cauchy with linear predictor on the mean (cauchy_p1),

• Exponential (exp),

• Exponential with log-linear predictor on the scale (exp_p1),

• Frechet with known location parameter (frechet_k1),

• Frechet with log-linear predictor on the scale and known location parameter (frechet_p2k1),

• Gamma (gamma),

• Generalized normal (gnorm),

• GEV (gev),

• GEV with linear predictor on the location (gev_p1),

• GEV with linear predictor on the location and log-linear prediction on the scale (gev_p12),

• GEV with linear predictor on the location, log-linear prediction on the scale, and linear predictor on the shape (gev_p123),

• GEV with linear predictor on the location and known shape (gev_p1k3),

• GEV with known shape (gev_k3),

• GPD with known location (gpd_k1),

• Gumbel (gumbel),

• Gumbel with linear predictor on the mean(gumbel_p1),

• Half-normal (halfnorm),

• Inverse gamma (invgamma),

• Inverse Gaussian (invgauss),

• t distribution with unknown location and scale and known DoF (lst_k3),

• t distribution with unknown location and scale, linear predictor on the location, and known DoF (lst_p1k3),

• Logistic (logis),

• Logistic with linear predictor on the location (logis_p1),

• Log-normal (lnorm),

• Log-normal with linear predictor on the location (lnorm_p1),

• Normal (norm),

• Normal with linear predictor on the mean (norm_p1),

• Pareto with known scale (pareto_k2),

• Pareto with log-linear predictor on the shape and known scale (pareto_p1k2),

• Uniform (unif),

• Weibull (weibull),

• Weibull with linear predictor on the scale (weibull_p2),

The level of predictive probability matching achieved by the maximum likelihood and calibrating prior quantiles, for any model, sample size and true parameter values, can be demonstrated using the routine reltest.

Model selection among models can be demonstrated using the routines ms_flat_1tail, ms_flat_2tail, ms_predictors_1tail, and ms_predictors_2tail,

Examples

#
# example 1
x=fitdistcp::d10exp_example_data_v1
p=c(1:9)/10
q=qexp_cp(x,p,rust=TRUE,nrust=1000)
xmin=min(q$ml_quantiles,q$cp_quantiles);
xmax=max(q$ml_quantiles,q$cp_quantiles);
plot(q$ml_quantiles,p,xlab="quantile estimates",xlim=c(xmin,xmax),
	sub="(from qexp_cp)",
	main="Exponential: quantile estimates");
points(q$cp_quantiles,p,col="red",lwd=2)
points(q$ru_quantiles,p,col="blue")

DMGS equation 2.1, f1 term

Description

DMGS equation 2.1, f1 term

Usage

exp_f1f(y, v1, fd1)

Arguments

Value

Matrix

The first derivative of the density

Description

The first derivative of the density

Usage

exp_f1fa(x, v1)

Arguments

Value

Vector

DMGS equation 2.1, f2 term

Description

DMGS equation 2.1, f2 term

Usage

exp_f2f(y, v1, fd1)

Arguments

Value

3d array

The second derivative of the density

Description

The second derivative of the density

Usage

exp_f2fa(x, v1)

Arguments

Value

Matrix

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_fd(x, v1)

Arguments

Value

Vector

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_fdd(x, v1)

Arguments

Value

Matrix

Third derivative of the normalized log-likelihood

Description

Third derivative of the normalized log-likelihood

Usage

exp_l111(x, v1, fd1)

Arguments

Value

Scalar value

The second derivative of the normalized log-likelihood

Description

The second derivative of the normalized log-likelihood

Usage

exp_ldd(x, v1, fd1)

Arguments

Value

Square scalar matrix

The second derivative of the normalized log-likelihood

Description

The second derivative of the normalized log-likelihood

Usage

exp_ldda(x, v1)

Arguments

Value

Matrix

Third derivative tensor of the log-likelihood

Description

Third derivative tensor of the log-likelihood

Usage

exp_lddd(x, v1, fd1)

Arguments

Value

Cubic scalar array

The third derivative of the normalized log-likelihood

Description

The third derivative of the normalized log-likelihood

Usage

exp_lddda(x, v1)

Arguments

Value

3d array

Logf for RUST

Description

Logf for RUST

Usage

exp_logf(params, x)

Arguments

Value

Scalar value.

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_logfdd(x, v1)

Arguments

Value

Matrix

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_logfddd(x, v1)

Arguments

Value

3d array

Log scores for MLE and RHP predictions calculated using leave-one-out

Description

Log scores for MLE and RHP predictions calculated using leave-one-out

Usage

exp_logscores(logscores, x)

Arguments

Value

Two scalars

Exponential Distribution with a Predictor, Predictions Based on a Calibrating Prior

Description

The fitdistcp package contains functions that generate predictive distributions for various statistical models, with and without parameter uncertainty. Parameter uncertainty is included by using Bayesian prediction with a type of objective prior known as a calibrating prior. Calibrating priors are chosen to give predictions that give good reliability (i.e., are well calibrated), for any underlying true parameter values.

There are five functions for each model, each of which uses training data x. For model **** the five functions are as follows:

• q****_cp returns predictive quantiles at the specified probabilities p, and various other diagnostics.

• r****_cp returns n random deviates from the predictive distribution.

• d****_cp returns the predictive density function at the specified values y

• p****_cp returns the predictive distribution function at the specified values y

• t****_cp returns n random deviates from the posterior distribution of the model parameters.

The q, r, d, p routines return two sets of results, one based on maximum likelihood, and the other based on a calibrating prior. The prior used depends on the model, and is given under Details below.

Where possible, the Bayesian prediction integral is solved analytically. Otherwise, DMGS asymptotic expansions are used. Optionally, a third set of results is returned that integrates the prediction integral by sampling the parameter posterior distribution using the RUST rejection sampling algorithm.

Usage

qexp_p1_cp(
  x,
  t,
  t0 = NA,
  n0 = NA,
  p = seq(0.1, 0.9, 0.1),
  d1 = 0.01,
  d2 = 0.01,
  means = FALSE,
  waicscores = FALSE,
  logscores = FALSE,
  dmgs = TRUE,
  rust = FALSE,
  nrust = 1e+05,
  predictordata = TRUE,
  centering = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

rexp_p1_cp(
  n,
  x,
  t,
  t0 = NA,
  n0 = NA,
  d1 = 0.01,
  d2 = 0.01,
  rust = FALSE,
  mlcp = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

dexp_p1_cp(
  x,
  t,
  t0 = NA,
  n0 = NA,
  y = x,
  d1 = 0.01,
  d2 = 0.01,
  rust = FALSE,
  nrust = 1000,
  centering = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

pexp_p1_cp(
  x,
  t,
  t0 = NA,
  n0 = NA,
  y = x,
  d1 = 0.01,
  d2 = 0.01,
  rust = FALSE,
  nrust = 1000,
  centering = TRUE,
  debug = FALSE,
  aderivs = TRUE
)

texp_p1_cp(n, x, t, d1 = 0.01, d2 = 0.01, debug = FALSE)

Arguments

Value

q**** returns a list containing at least the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_value: the value of the log-likelihood at the maximum.

• standard_errors: estimates of the standard errors on the parameters, from the inverse observed information matrix.

• ml_quantiles: quantiles calculated using maximum likelihood.

• cp_quantiles: predictive quantiles calculated using a calibrating prior.

• maic: the AIC score for the maximum likelihood model, times -1/2.

• cp_method: a comment about the method used to generate the cp prediction.

For models with predictors, q**** additionally returns:

• predictedparameter: the estimated value for parameter, as a function of the predictor.

• adjustedx: the detrended values of x

r**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_deviates: random deviates calculated using maximum likelihood.

• cp_deviates: predictive random deviates calculated using a calibrating prior.

• cp_method: a comment about the method used to generate the cp prediction.

d**** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_pdf: density function from maximum likelihood.

• cp_pdf: predictive density function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

p*** returns a list containing the following:

• ml_params: maximum likelihood estimates for the parameters.

• ml_cdf: distribution function from maximum likelihood.

• cp_cdf: predictive distribution function calculated using a calibrating prior (not available in EVT routines, for mathematical reasons, unless using RUST).

• cp_method: a comment about the method used to generate the cp prediction.

t*** returns a list containing the following:

• theta_samples: random samples from the parameter posterior.

Details of the Model

The exponential distribution with a predictor has exceedance distribution function

S(x;a,b)=\exp(-x \lambda (a,b))

where x \ge 0 is the random variable and \lambda(a,b)=e^{-a-bt} is the rate parameter, modelled as a function of the parameters a,b and a predictor t.

The calibrating prior is given by the right Haar prior, which is

\pi(a,b) \propto 1

. as given in Jewson et al. (2025).

Optional Return Values

q**** optionally returns the following:

If rust=TRUE:

• ru_quantiles: predictive quantiles calculated using a calibrating prior, using posterior sampling with the RUST algorithm, based on inverting an empirical CDF based on nrust samples.

If waicscores=TRUE:

• waic1: the WAIC1 score for the calibrating prior model.

• waic2: the WAIC2 score for the calibrating prior model.

If logscores=TRUE:

• ml_oos_logscore: the leave-one-out logscore for the maximum likelihood prediction (not available in EVT routines, for mathematical reasons)

• cp_oos_logscore: the leave-one-out logscore for the parameter uncertainty model available in EVT routines, for mathematical reasons)

If means=TRUE:

• ml_mean: analytic estimate of the mean of the MLE predictive distribution, where possible

• cp_mean: analytic estimate of the mean of the calibrating prior predictive distribution, where mathematically possible. Can be compared with the mean estimated from random deviates.

r**** optionally returns the following:

If rust=TRUE:

• ru_deviates: nrust predictive random deviatives calculated using a calibrating prior, using posterior sampling with RUST.

d**** optionally returns the following:

If rust=TRUE:

• ru_pdf: predictive density calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust density functions.

p**** optionally returns the following:

If rust=TRUE:

• ru_cdf: predictive probability calculated using a calibrating prior, using posterior sampling with RUST, averaging over nrust distribution functions.

Selecting these additional outputs increases runtime. They are optional so that runtime for the basic outputs is minimised. This facilitates repeated experiments that evaluate reliability over many thousands of repeats.

Details (homogeneous models)

This model is a homogeneous model, and the cp results are based on the right Haar prior. For homogeneous models (models with sharply transitive transformation groups), a Bayesian prediction based on the right Haar prior gives exact reliability, as shown by Severini et al. (2002), even when the true parameters are unknown. This means that probabilities in the prediction will correspond to frequencies of future outcomes in repeated trials (if the model is correct).

Maximum likelihood prediction does not give reliable predictions, even when the model is correct, because it does not account for parameter uncertainty. In particular, maximum likelihood predictions typically underestimate the tail in repeated trials.

The reliability of the maximum likelihood and the calibrating prior predictive quantiles produced by the q****_cp routines in fitdistcp can be quantified using repeated simulations with the routine reltest.

Details (DMGS integration)

For this model, the Bayesian prediction equation cannot be solved analytically, and is approximated using the DMGS asymptotic expansions given by Datta et al. (2000). This approximation seems to work well for medium and large sample sizes, but may not work well for small sample sizes (e.g., <20 data points). For small sample sizes, it may be preferable to use posterior sampling by setting rust=TRUE and looking at the ru outputs. The performance for any sample size, in terms of reliability, can be tested using reltest.

Details (RUST)

The Bayesian prediction equation can also be integrated using ratio-of-uniforms-sampling-with-transformation (RUST), using the option rust=TRUE. fitdistcp then calls Paul Northrop's rust package (Northrop, 2023). The RUST calculations are slower than the DMGS calculations.

For small sample sizes (e.g., n<20), and the very extreme tail, the DMGS approximation is somewhat poor (although always better than maximum likelihood) and it may be better to use RUST. For medium sample sizes (30+), DMGS is reasonably accurate, except for the very far tail.

It is advisable to check the RUST results for convergence versus the number of RUST samples.

It may be interesting to compare the DMGS and RUST results.

Author(s)

Stephen Jewson stephen.jewson@gmail.com

References

If you use this package, we would be grateful if you would cite the following reference, which gives the various calibrating priors, and tests them for reliability:

• Jewson S., Sweeting T. and Jewson L. (2024): Reducing Reliability Bias in Assessments of Extreme Weather Risk using Calibrating Priors; ASCMO Advances in Statistical Climatology, Meteorology and Oceanography), https://ascmo.copernicus.org/articles/11/1/2025/.

See Also

An introduction to fitdistcp, with more examples, is given on this webpage.

The fitdistcp package currently includes the following models (in alphabetical order):

• Cauchy (cauchy),

• Cauchy with linear predictor on the mean (cauchy_p1),

• Exponential (exp),

• Exponential with log-linear predictor on the scale (exp_p1),

• Frechet with known location parameter (frechet_k1),

• Frechet with log-linear predictor on the scale and known location parameter (frechet_p2k1),

• Gamma (gamma),

• Generalized normal (gnorm),

• GEV (gev),

• GEV with linear predictor on the location (gev_p1),

• GEV with linear predictor on the location and log-linear prediction on the scale (gev_p12),

• GEV with linear predictor on the location, log-linear prediction on the scale, and linear predictor on the shape (gev_p123),

• GEV with linear predictor on the location and known shape (gev_p1k3),

• GEV with known shape (gev_k3),

• GPD with known location (gpd_k1),

• Gumbel (gumbel),

• Gumbel with linear predictor on the mean(gumbel_p1),

• Half-normal (halfnorm),

• Inverse gamma (invgamma),

• Inverse Gaussian (invgauss),

• t distribution with unknown location and scale and known DoF (lst_k3),

• t distribution with unknown location and scale, linear predictor on the location, and known DoF (lst_p1k3),

• Logistic (logis),

• Logistic with linear predictor on the location (logis_p1),

• Log-normal (lnorm),

• Log-normal with linear predictor on the location (lnorm_p1),

• Normal (norm),

• Normal with linear predictor on the mean (norm_p1),

• Pareto with known scale (pareto_k2),

• Pareto with log-linear predictor on the shape and known scale (pareto_p1k2),

• Uniform (unif),

• Weibull (weibull),

• Weibull with linear predictor on the scale (weibull_p2),

The level of predictive probability matching achieved by the maximum likelihood and calibrating prior quantiles, for any model, sample size and true parameter values, can be demonstrated using the routine reltest.

Model selection among models can be demonstrated using the routines ms_flat_1tail, ms_flat_2tail, ms_predictors_1tail, and ms_predictors_2tail,

Examples

#
# example 1
x=fitdistcp::d55exp_p1_example_data_v1_x
tt=fitdistcp::d55exp_p1_example_data_v1_t
p=c(1:9)/10
n0=10
q=qexp_p1_cp(x,tt,n0=n0,p=p,rust=TRUE,nrust=1000)
xmin=min(q$ml_quantiles,q$cp_quantiles);
xmax=max(q$ml_quantiles,q$cp_quantiles);
plot(q$ml_quantiles,p,xlab="quantile estimates",xlim=c(xmin,xmax),
	sub="(from qexp_p1_cp)",
	main="Exponential w/ p1: quantile estimates");
points(q$cp_quantiles,p,col="red",lwd=2)
points(q$ru_quantiles,p,col="blue")

DMGS equation 2.1, f1 term

Description

DMGS equation 2.1, f1 term

Usage

exp_p1_f1f(y, t0, v1, d1, v2, d2)

Arguments

Value

Matrix

The first derivative of the density

Description

The first derivative of the density

Usage

exp_p1_f1fa(x, t, v1, v2)

Arguments

Value

Vector

DMGS equation 2.1, f2 term

Description

DMGS equation 2.1, f2 term

Usage

exp_p1_f2f(y, t0, v1, d1, v2, d2)

Arguments

Value

3d array

The second derivative of the density

Description

The second derivative of the density

Usage

exp_p1_f2fa(x, t, v1, v2)

Arguments

Value

Matrix

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

First derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_p1_fd(x, t, v1, v2)

Arguments

Value

Vector

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_p1_fdd(x, t, v1, v2)

Arguments

Value

Matrix

Second derivative matrix of the normalized log-likelihood

Description

Second derivative matrix of the normalized log-likelihood

Usage

exp_p1_ldd(x, t, v1, d1, v2, d2)

Arguments

Value

Square scalar matrix

The second derivative of the normalized log-likelihood

Description

The second derivative of the normalized log-likelihood

Usage

exp_p1_ldda(x, t, v1, v2)

Arguments

Value

Matrix

Third derivative tensor of the normalized log-likelihood

Description

Third derivative tensor of the normalized log-likelihood

Usage

exp_p1_lddd(x, t, v1, d1, v2, d2)

Arguments

Value

Cubic scalar array

The third derivative of the normalized log-likelihood

Description

The third derivative of the normalized log-likelihood

Usage

exp_p1_lddda(x, t, v1, v2)

Arguments

Value

3d array

One component of the second derivative of the normalized log-likelihood

Description

One component of the second derivative of the normalized log-likelihood

Usage

exp_p1_lmn(x, t, v1, d1, v2, d2, mm, nn)

Arguments

Value

Scalar value

One component of the third derivative of the normalized log-likelihood

Description

One component of the third derivative of the normalized log-likelihood

Usage

exp_p1_lmnp(x, t, v1, d1, v2, d2, mm, nn, rr)

Arguments

Value

Scalar value

Logf for RUST

Description

Logf for RUST

Usage

exp_p1_logf(params, x, t)

Arguments

Value

Scalar value.

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Second derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_p1_logfdd(x, t, v1, v2)

Arguments

Value

Matrix

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Description

Third derivative of the log density Created by Stephen Jewson using Deriv() by Andrew Clausen and Serguei Sokol

Usage

exp_p1_logfddd(x, t, v1, v2)

Arguments

Value

3d array

observed log-likelihood function

Description

observed log-likelihood function

Usage

exp_p1_loglik(vv, x, t)

Arguments

Value

Scalar value.

Log scores for MLE and RHP predictions calculated using leave-one-out

Description

Log scores for MLE and RHP predictions calculated using leave-one-out

Usage

exp_p1_logscores(logscores, x, t, d1, d2, aderivs)

Arguments

Value

Two scalars

exp distribution: RHP means

Description

exp distribution: RHP means

Usage

exp_p1_means(means, t0, ml_params, lddi, lddd, lambdad_rhp, nx, dim = 2)

Arguments

Value

Two scalars

DMGS equation 3.3, mu1 term

Description

DMGS equation 3.3, mu1 term

Usage

exp_p1_mu1f(alpha, t0, v1, d1, v2, d2)

Arguments

Value

Matrix

Minus the first derivative of the cdf, at alpha

Description

Minus the first derivative of the cdf, at alpha

Usage

exp_p1_mu1fa(alpha, t, v1, v2)

Arguments

Value

Vector

DMGS equation 3.3, mu2 term

Description

DMGS equation 3.3, mu2 term

Usage

exp_p1_mu2f(alpha, t0, v1, d1, v2, d2)

Arguments

Value

3d array

Minus the second derivative of the cdf, at alpha

Description

Minus the second derivative of the cdf, at alpha

Usage

exp_p1_mu2fa(alpha, t, v1, v2)

Arguments