Package {FitVerse}

Type:	Package
Title:	Parametric Distribution Fitting and Analysis
Version:	1.0-3
Date:	2026-11-13
Description:	Provides a unified, user-friendly interface for fitting parametric probability distributions to continuous univariate data. 'FitVerse' supports 52 distribution families spanning symmetric, right-skewed, heavy-tailed, bounded, and extreme-value shapes, and three estimation methods: Maximum Likelihood Estimation (MLE), Method of Moments (MOM), and L-Moments (L-MOM). Automatic best-fit selection is performed using AIC, BIC, and goodness-of-fit tests (Kolmogorov-Smirnov, Anderson-Darling, Cramer-von Mises (CvM)). Every fitted model produces a publication-quality diagnostic plot: a histogram overlaid with the fitted density curve and the estimated PDF formula annotated directly on the figure. An optional interactive version is produced via 'plotly'. Additional tools include bootstrap confidence intervals for parameter estimates and return levels, batch fitting across multiple columns for automated workflows and web-upload use cases, JSON serialisation for integration with 'Shiny' web applications, and automated HTML/PDF report generation. 'FitVerse' is designed to support data characterisation in survey sampling, hydrology, and actuarial workflows, where identifying the underlying distribution of a variable is a prerequisite for downstream modelling and inference. L-moment estimation follows Hosking (1990) <doi:10.1111/j.2517-6161.1990.tb01775.x> and Hosking and Wallis (1997, ISBN:9780521430456). Model selection via AIC follows Akaike (1974) <doi:10.1109/TAC.1974.1100705> and via BIC follows Schwarz (1978) <doi:10.1214/aos/1176344136>. Bootstrap confidence intervals follow Efron and Hastie (2016, ISBN:9781107149892).
License:	GPL (≥ 3)
Encoding:	UTF-8
Depends:	R (≥ 4.1.0)
Imports:	stats, graphics, grDevices, utils, tools, actuar (≥ 3.1.0), mc2d (≥ 0.1.18), evd (≥ 2.3.3), sn (≥ 2.1.0), fitdistrplus (≥ 1.1.0), lmomco (≥ 2.3.7), ggplot2 (≥ 3.4.0), goftest (≥ 1.2.3), gridExtra (≥ 2.3)
Suggests:	shiny (≥ 1.7.0), bslib (≥ 0.5.0), DT, plotly (≥ 4.10.0), testthat (≥ 3.0.0), knitr (≥ 1.40), rmarkdown (≥ 2.20), kableExtra, jsonlite (≥ 1.8.0), base64enc, MASS, GeneralizedHyperbolic, truncnorm
VignetteBuilder:	knitr
NeedsCompilation:	no
Config/roxygen2/version:	8.0.0
Packaged:	2026-11-13 06:47:51 UTC; karunareddy
Author:	Karuna G. Reddy [aut, cre], M. G. M. Khan [aut]
Maintainer:	Karuna G. Reddy <karuna.reddy@auckland.ac.nz>
Repository:	CRAN
Date/Publication:	2026-11-13 07:30:08 UTC

FitVerse: Parametric Distribution Fitting and Analysis

Description

FitVerse provides a unified, user-friendly interface for fitting parametric probability distributions to continuous univariate data. The package supports 52 distribution families spanning symmetric, right-skewed, heavy-tailed, bounded, and extreme-value shapes, and three estimation methods: Maximum Likelihood Estimation (MLE), Method of Moments (MOM), and L-Moments (LMOM).

Automatic best-fit selection is performed using AIC and BIC, supplemented by four goodness-of-fit tests: Kolmogorov-Smirnov (KS), Anderson-Darling (AD), Cramer-von Mises (CvM), and Chi-Squared.

Every fitted model produces a publication-quality diagnostic plot: a histogram overlaid with the fitted density curve and the estimated PDF formula annotated directly on the figure. Additional diagnostic panels (Q-Q, P-P, CDF, Survival, Hazard, and Cumulative Hazard) are available via plot.fitverse().

FitVerse is designed as the data-characterisation layer of the OptiStrata ecosystem, with seamless handoff to stratifyR for optimal stratification of survey populations.

Main fitting functions

Probability tools

Inference and reporting

Sample datasets

S3 methods for "fitverse" objects

Supported distributions (52 families)

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz (University of Auckland, New Zealand)
M. G. M. Khan mgm.khan@usp.ac.fj (University of the South Pacific, Fiji)

Maintainer: Karuna G. Reddy karuna.reddy@auckland.ac.nz

References

Hosking, J. R. M. (1990). L-moments: Analysis and estimation of distributions using linear combinations of order statistics. Journal of the Royal Statistical Society, Series B, 52(1), 105–124. doi:10.1111/j.2517-6161.1990.tb01775.x

Delignette-Muller, M. L. and Dutang, C. (2015). fitdistrplus: An R package for fitting distributions. Journal of Statistical Software, 64(4), 1–34. doi:10.18637/jss.v064.i04

Akaike, H. (1974). A new look at the statistical model identification. IEEE Transactions on Automatic Control, 19(6), 716–723. doi:10.1109/TAC.1974.1100705

See Also

fitverse, fitverse_distributions, fitverse_grouped, FitVerse_datasets, fitdistrplus, actuar, lmomco

Examples

## Quick start: fit three distributions to simulated Gamma data
set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)

fit <- fitverse(x,
                dists   = c("gamma", "lognormal", "weibull"),
                method  = "MLE",
                xlab    = "Waiting time (s)",
                plot    = FALSE,
                verbose = FALSE)

print(fit)

## See all 52 supported distributions
fitverse_distributions()

Build a standardised L-moment fit result list

Description

Computes the log-likelihood of the estimated parameters on the data, then packages everything into the standard FitVerse result format.

Usage

.build_lmom_result(x, dist_key, params, n)

Build a CDF overlay plot for the report (top_n distributions)

Description

Build a CDF overlay plot for the report (top_n distributions)

Usage

.build_report_cdf_plot(fv_obj, top_n = 3L, xlab = "x")

Build a P-P plot for the report (Hazen plotting positions)

Description

Build a P-P plot for the report (Hazen plotting positions)

Usage

.build_report_pp_plot(fv_obj, xlab = "x")

Arrange two ggplot objects side by side with an optional formula strip below.

Description

Arrange two ggplot objects side by side with an optional formula strip below.

Usage

.combine_panels(p1, p2, fv_obj = NULL)

Arguments

Expose d/p/q functions from actuar or mc2d into a caller's environment

Description

fitdistrplus finds density functions via get() on the search path. When actuar/mc2d are in Imports (not Depends), their functions are invisible to fitdistrplus. Calling this inside a fitting function assigns the d/p/q functions into that function's own frame, making them visible to fitdistrplus via parent.frame() traversal.

Usage

.expose_dist_fns(pkg, r_name, envir = parent.frame())

Arguments

Return an all-NA GoF data frame (four rows, one per test)

Description

Return an all-NA GoF data frame (four rows, one per test)

Usage

.gof_na_frame()

Build a single GoF result row

Description

Build a single GoF result row

Usage

.gof_row(test_name, statistic, p_value)

Fit a single distribution by grouped MLE

Description

Maximises sum_i n_i * log[ F(upper_i) - F(lower_i) ] using optim() with Nelder-Mead. Starting values are taken from the standard MLE fit on the pseudo-observations (midpoint x count).

Usage

.grouped_mle(lower_bd, upper_bd, counts, dist_key, x_pseudo, verbose = FALSE)

Arguments

L-Moments: Exponential parexp returns: para = c(xi, alpha) where xi is location, alpha is scale R's dexp uses rate = 1/alpha (valid when xi ~= 0)

Description

L-Moments: Exponential parexp returns: para = c(xi, alpha) where xi is location, alpha is scale R's dexp uses rate = 1/alpha (valid when xi ~= 0)

Usage

.lmom_exponential(x, lmom, n)

L-Moments: Gamma pargam returns: para = c(alpha, beta) where beta is the SCALE (not rate) R's dgamma uses shape and rate; we convert: rate = 1/beta

Description

L-Moments: Gamma pargam returns: para = c(alpha, beta) where beta is the SCALE (not rate) R's dgamma uses shape and rate; we convert: rate = 1/beta

Usage

.lmom_gamma(x, lmom, n)

L-Moments: Generalised Extreme Value (GEV) pargev returns: para = c(xi, alpha, kappa) IMPORTANT: lmomco kappa has OPPOSITE sign to evd shape. evd shape = -lmomco_kappa

Description

L-Moments: Generalised Extreme Value (GEV) pargev returns: para = c(xi, alpha, kappa) IMPORTANT: lmomco kappa has OPPOSITE sign to evd shape. evd shape = -lmomco_kappa

Usage

.lmom_gev(x, lmom, n)

L-Moments: Generalised Logistic (GLO) - full shape parameter parglo returns: para = c(xi, alpha, kappa) No sign flip for GLO: lmomco kappa = our shape. (The existing .lmom_logistic drops kappa to get standard logistic; .lmom_glo keeps it for the full GLO.)

Description

L-Moments: Generalised Logistic (GLO) - full shape parameter parglo returns: para = c(xi, alpha, kappa) No sign flip for GLO: lmomco kappa = our shape. (The existing .lmom_logistic drops kappa to get standard logistic; .lmom_glo keeps it for the full GLO.)

Usage

.lmom_glo(x, lmom, n)

L-Moments: Generalised Pareto Distribution (GPD) pargpa returns: para = c(xi, alpha, kappa) IMPORTANT: lmomco kappa has OPPOSITE sign to our shape (same convention as GEV). our shape = -lmomco_kappa

Description

L-Moments: Generalised Pareto Distribution (GPD) pargpa returns: para = c(xi, alpha, kappa) IMPORTANT: lmomco kappa has OPPOSITE sign to our shape (same convention as GEV). our shape = -lmomco_kappa

Usage

.lmom_gpd(x, lmom, n)

L-Moments: Gumbel pargum returns: para = c(xi, alpha) -> loc = xi, scale = alpha

Description

L-Moments: Gumbel pargum returns: para = c(xi, alpha) -> loc = xi, scale = alpha

Usage

.lmom_gumbel(x, lmom, n)

L-Moments: Inverse Gaussian (Wald)

Description

lmomco::parwig fits the Wald distribution. Returns: para = c(mu, lambda)

Usage

.lmom_invgaussian(x, lmom, n)

L-Moments: Kappa-4 (Hosking 1994) parkap returns: para = c(xi, alpha, kappa, h) These map directly to our registry params: xi, alpha, kappa, h. Requires nmom >= 4 (already computed above with nmom=5).

Description

L-Moments: Kappa-4 (Hosking 1994) parkap returns: para = c(xi, alpha, kappa, h) These map directly to our registry params: xi, alpha, kappa, h. Requires nmom >= 4 (already computed above with nmom=5).

Usage

.lmom_kappa4(x, lmom, n)

L-Moments: Logistic (via Generalised Logistic with kappa fixed to 0) parglo returns: para = c(xi, alpha, kappa) Standard logistic: kappa = 0, location = xi, scale = alpha

Description

L-Moments: Logistic (via Generalised Logistic with kappa fixed to 0) parglo returns: para = c(xi, alpha, kappa) Standard logistic: kappa = 0, location = xi, scale = alpha

Usage

.lmom_logistic(x, lmom, n)

L-Moments: Log-Normal (3-parameter, with location fixed near 0) parln3 returns: para = c(zeta, mu, sigma) We use the 2-parameter form by ignoring zeta (lmomco sets it close to 0).

Description

L-Moments: Log-Normal (3-parameter, with location fixed near 0) parln3 returns: para = c(zeta, mu, sigma) We use the 2-parameter form by ignoring zeta (lmomco sets it close to 0).

Usage

.lmom_lognormal(x, lmom, n)

L-Moments: Log-Pearson Type III Apply Pearson III L-moment fitting to log(x), then convert.

Description

L-Moments: Log-Pearson Type III Apply Pearson III L-moment fitting to log(x), then convert.

Usage

.lmom_lpearson3(x, n)

L-Moments: Normal parnor returns: para = c(mu, sigma) -> R names: mean, sd

Description

L-Moments: Normal parnor returns: para = c(mu, sigma) -> R names: mean, sd

Usage

.lmom_normal(x, lmom, n)

L-Moments: Pearson Type III parpe3 returns: para = c(mu, sigma, gamma) - mean, std dev, skewness Convert to our (shape, scale, loc) parameterisation.

Description

L-Moments: Pearson Type III parpe3 returns: para = c(mu, sigma, gamma) - mean, std dev, skewness Convert to our (shape, scale, loc) parameterisation.

Usage

.lmom_pearson3(x, lmom, n)

L-Moments: Wakeby (five-parameter quantile distribution) parwak returns: para = c(xi, alpha, beta, gamma, delta) These map directly to our registry params: xi, alpha, beta, gamma, delta. Requires nmom >= 5 (already computed above with nmom=5).

Description

L-Moments: Wakeby (five-parameter quantile distribution) parwak returns: para = c(xi, alpha, beta, gamma, delta) These map directly to our registry params: xi, alpha, beta, gamma, delta. Requires nmom >= 5 (already computed above with nmom=5).

Usage

.lmom_wakeby(x, lmom, n)

L-Moments: Weibull parwei returns: para = c(zeta, beta, delta) where zeta = location (~=0), beta = scale, delta = shape R's dweibull: shape = delta, scale = beta

Description

L-Moments: Weibull parwei returns: para = c(zeta, beta, delta) where zeta = location (~=0), beta = scale, delta = shape R's dweibull: shape = delta, scale = beta

Usage

.lmom_weibull(x, lmom, n)

L-Moments: 3-Parameter Weibull parwei returns: para = c(zeta, beta, delta) where zeta = threshold, beta = scale, delta = shape We pass all three through, unlike the 2P handler which drops zeta.

Description

L-Moments: 3-Parameter Weibull parwei returns: para = c(zeta, beta, delta) where zeta = threshold, beta = scale, delta = shape We pass all three through, unlike the 2P handler which drops zeta.

Usage

.lmom_weibull3p(x, lmom, n)

Evaluate the log-likelihood for a given set of parameters

Description

Used by manual MOM estimators that produce parameters without an embedded log-likelihood (because they bypass fitdistrplus).

Usage

.loglik_from_params(x, dist_key, params)

Arguments

Value

Numeric log-likelihood, or NA on error.

Resolve a d/p/q/r function by distribution key string

Description

Checks the closed-form override table first, then match.fun, then getExportedValue. Used by both probability_tools.R and plot_fitverse.R.

Usage

.lookup_dpqr(prefix, dist_key)

Arguments

Value

A function, or NULL.

Build a two-column footer grob: PDF formula (left) and parameter estimates (right). Placed below the main panels by .combine_panels().

Description

Build a two-column footer grob: PDF formula (left) and parameter estimates (right). Placed below the main panels by .combine_panels().

Usage

.make_formula_strip(fv_obj)

MLE for the Birnbaum-Saunders (fatigue-life) distribution

Description

Correct closed-form MLEs (Birnbaum & Saunders 1969; Ng et al. 2003): nu_tilde = mean(sqrt(x)) mu_tilde = mean(1 / sqrt(x)) beta_hat = nu_tilde / mu_tilde alpha_hat = sqrt(mean(x)/beta_hat + beta_hat*mean(1/x) - 2)

Usage

.mle_birnbaum_saunders(x)

MLE for the Burr XII distribution with multi-start optimisation

Description

The generic .mle_custom_optim uses a single starting point that can fail to converge on Burr XII's three-parameter likelihood surface. This dedicated fitter tries a 3x3 grid of (shape1, shape2) starting values with scale fixed at the sample median, returning the best result.

Usage

.mle_burr(x, entry)

Generic MLE via stats::optim using the internal d/p/q/r override table

Description

Used for distributions whose d-functions are in .FITVERSE_DPQR_OVERRIDES and cannot be passed to fitdistrplus directly.

Usage

.mle_custom_optim(x, dist_key, entry)

MLE for the Generalised Extreme Value distribution using evd::fgev

Description

MLE for the Generalised Extreme Value distribution using evd::fgev

Usage

.mle_evd_gev(x, entry)

MLE for the Gumbel distribution via analytical profile likelihood. The MLE scale beta solves the 1-D fixed-point equation beta = mean(x) - sum(x * w) / sum(w), w = exp(-x/beta) (Prescott & Walden, 1980). Given beta, mu = -beta * log(mean(exp(-x/beta))). This avoids evd::fgev's optim iteration-limit false non-convergence.

Description

MLE for the Gumbel distribution via analytical profile likelihood. The MLE scale beta solves the 1-D fixed-point equation beta = mean(x) - sum(x * w) / sum(w), w = exp(-x/beta) (Prescott & Walden, 1980). Given beta, mu = -beta * log(mean(exp(-x/beta))). This avoids evd::fgev's optim iteration-limit false non-convergence.

Usage

.mle_evd_gumbel(x, entry)

MLE for Exponential: rate = 1 / mean(x) (exact, no optimizer)

Description

Bypasses fitdistrplus entirely to avoid try()-generated console noise when fitdistrplus's internal finite-difference gradient computation fails for small rate values (large mean).

Usage

.mle_exponential(x)

MLE via fitdistrplus::fitdist

Description

This covers all distributions whose d/p/q/r functions are available on the search path via stats or actuar: normal, lognormal, gamma, exponential, weibull, beta, logistic, cauchy, uniform, loglogistic, pareto, burr, invgamma, invweibull, paralogistic (15 distributions). All others have dedicated fitters and are NOT routed through this function.

Usage

.mle_fitdistrplus(x, dist_key, entry)

For fixed thres, Y = x - thres ~ Gamma(shape, scale). The Gamma MLE on Y (Newton-Raphson) gives shape and scale analytically. A 1-D search over thres finds the global MLE robustly.

Description

For fixed thres, Y = x - thres ~ Gamma(shape, scale). The Gamma MLE on Y (Newton-Raphson) gives shape and scale analytically. A 1-D search over thres finds the global MLE robustly.

Usage

.mle_gamma3(x)

Profile-likelihood MLE for the Generalised Gamma distribution

Description

Reduces the 3-parameter optimisation to a 1-D problem: for fixed p, let Y = X^p ~ Gamma(d/p, scale = b^p). The Gamma MLE for Y gives d and b analytically via Newton-Raphson. A grid search + 1-D refinement over p finds the global optimum reliably.

Usage

.mle_gengamma(x)

MLE for GHyp: uses ghypFit() if GeneralizedHyperbolic is available, otherwise falls back to reparametrized pure-R optimizer.

Description

MLE for GHyp: uses ghypFit() if GeneralizedHyperbolic is available, otherwise falls back to reparametrized pure-R optimizer.

Usage

.mle_ghyp(x, entry)

MLE for the Generalised Logistic distribution

Description

The generic '.mle_custom_optim' route starts with shape = 0 (logistic). Because the '.dglo_internal' density switches to the logistic formula when |shape| < 1e-8, the numerical gradient in the shape direction is effectively zero at the starting point, and L-BFGS-B never moves shape away from 0.

Usage

.mle_glo(x, entry)

Details

Fix: use L-moments (via lmomco::parglo) as starting values, which give an accurate initial shape estimate. Fall back to a small positive shape (0.1) if L-moments fail.

MLE for the Gompertz distribution

Description

Starting values via log-linear regression on the empirical CDF: log(-log(1-F(x))) = log(alpha/beta) + beta*x so slope = beta, intercept = log(alpha/beta) => alpha = beta*exp(intercept). This reliably recovers small beta values that trip up generic starting heuristics.

Usage

.mle_gompertz(x)

sigma_hat = sqrt(mean(x^2))

Description

sigma_hat = sqrt(mean(x^2))

Usage

.mle_halfnormal(x)

MLE for the Inverse Gaussian (Wald) distribution

Description

Closed-form MLEs: mu_hat = mean(x) lambda_hat = n / sum(1/x - 1/mean(x))

Usage

.mle_invgaussian(x)

MLE for the Laplace (double-exponential) distribution

Description

Closed-form MLEs: mu_hat = median(x) b_hat = mean(|x - mu_hat|)

Usage

.mle_laplace(x)

MLE for Lindley: solve quadratic score equation Score: 2/theta - 1/(1+theta) = mean(x) => theta^2*xbar + theta*(xbar-1) - 2 = 0

Description

MLE for Lindley: solve quadratic score equation Score: 2/theta - 1/(1+theta) = mean(x) => theta^2*xbar + theta*(xbar-1) - 2 = 0

Usage

.mle_lindley(x)

For fixed thres, Y = x - thres, then meanlog = mean(log Y), sdlog = sd(log Y).

Description

For fixed thres, Y = x - thres, then meanlog = mean(log Y), sdlog = sd(log Y).

Usage

.mle_lnorm3(x)

MLE for the Log-Pearson III distribution

Description

Strategy: fit Pearson III directly on log(x) and adjust the log-likelihood for the Jacobian: loglik_LP3 = loglik_P3(log x) - [U+03A3]log(x).

Usage

.mle_lpearson3(x, entry)

Details

Root cause of past failures: the MOM skewness estimate of log(x) can be near-zero for some random seeds, producing shape = 4/gam^2 -> Inf. The Gamma density at data points is then effectively 0 everywhere, optim gets a flat 1e300 objective from the first step, cannot compute a gradient, and errors out. L-moments (via lmomco::parpe3) are robust to this.

MLE for Nakagami: closed-form estimators Omega_hat = mean(x^2); m_hat from moment equation solved numerically. Actually lmomco-free profile MLE: use optim on (m, Omega) jointly.

Description

MLE for Nakagami: closed-form estimators Omega_hat = mean(x^2); m_hat from moment equation solved numerically. Actually lmomco-free profile MLE: use optim on (m, Omega) jointly.

Usage

.mle_nakagami(x)

MLE for NIG: uses nigFit() if GeneralizedHyperbolic is available, otherwise falls back to reparametrized pure-R optimizer.

Description

MLE for NIG: uses nigFit() if GeneralizedHyperbolic is available, otherwise falls back to reparametrized pure-R optimizer.

Usage

.mle_nig(x, entry)

MLE for the Pearson Type III distribution

Description

Pearson III is a shifted/scaled Gamma: X = loc + scale * Gamma(shape). When scale > 0 (right-skewed), the support is x > loc. The generic .mle_custom_optim path can produce starting values where loc >= min(x), giving log-likelihood = -Inf at the first step; L-BFGS-B cannot move.

Usage

.mle_pearson3(x, entry)

Details

Fix: use L-moments for starting values (robust to near-zero skewness) and clamp loc strictly below min(x) before entering the optimizer.

MLE for the PERT distribution

Description

The practical MLE profiles over mode using stats::optimize(), with min and max expanded beyond the sample extremes (shape parameter fixed at the standard PERT value of 4).

Usage

.mle_pert(x)

Details

WHY A LARGER BUFFER THAN TRIANGULAR: The PERT density approaches 0 at its boundaries as (x - min)^(alpha1 - 1) where alpha1 > 1 (e.g. alpha1 ~= 2.6 for mode near the midpoint). This power-law decay is faster than the linear decay of the Triangular density. If min/max are set too close to the sample extremes, every observation near the true boundary gets an artificially tiny log-density, collapsing the total log-likelihood. A buffer of rng / sqrt(n) scales the expansion with sample size: large n needs less padding because sample extremes are already close to the true bounds, while small n benefits from wider coverage.

sigma_hat = sqrt(mean(x^2) / 2)

Description

sigma_hat = sqrt(mean(x^2) / 2)

Usage

.mle_rayleigh(x)

MLE for the Right Triangular distribution

Description

The right triangular PDF is f(x) = 2(b-x)/(b-a)^2 on [a,b]. The MLE of a is min(x) and b is max(x) (boundary estimators). Only 1 effective free parameter (k=1) for AIC purposes, same logic as Triangular and PERT.

Usage

.mle_rtriangle(x)

MLE for the Skew-Normal distribution using sn::selm

Description

MLE for the Skew-Normal distribution using sn::selm

Usage

.mle_sn(x, entry)

MLE for the Student-t (Location-Scale) distribution with multi-start

Description

The generic .mle_custom_optim uses a fixed starting value of nu=10. For data with heavy tails (true nu = 3–5), the optimizer can converge to a poor local solution when started near the Normal. This dedicated fitter tries five candidate starting values for nu (derived from the sample excess kurtosis plus a fixed grid) and returns the result with the highest log-likelihood.

Usage

.mle_student_t(x, entry)

MLE for the Triangular distribution

Description

The MLE for min and max are the sample extremes. We profile over mode using stats::optimize() on the interval (min, max), calling mc2d::dtriang directly (bypassing fitdistrplus and search-path lookup issues).

Usage

.mle_triangular(x)

MLE for Truncated Normal with fixed bounds at data min/max

Description

MLE for Truncated Normal with fixed bounds at data min/max

Usage

.mle_truncnorm(x)

MLE for the 3-Parameter Weibull distribution

Description

Uses stats::optim with L-BFGS-B directly so that the critical constraint \mu < \min(x) can be enforced via a dynamic upper bound on the threshold parameter - something fitdistrplus cannot do because its bounds are static.

Usage

.mle_weibull3p(x, entry)

Details

The log-likelihood is evaluated through dweibull3p (defined in utils.R) so the same density is used consistently across MLE, GoF, and plotting.

MOM for Birnbaum-Saunders: match E[X]=beta(1+alpha^2/2), Var[X]. beta = mean(x)/(1 + alpha^2/2); solve alpha from variance equation numerically.

Description

MOM for Birnbaum-Saunders: match E[X]=beta(1+alpha^2/2), Var[X]. beta = mean(x)/(1 + alpha^2/2); solve alpha from variance equation numerically.

Usage

.mom_birnbaum_saunders(x)

MOM for Exponentiated Exponential: use digamma-based moment equations. E[X] = (digamma(alpha+1) + euler_gamma) / lambda Var[X] = (pi^2/6 - trigamma(alpha+1)) ... actually use trigamma. Match observed mean and variance numerically to find (alpha, lambda).

Description

MOM for Exponentiated Exponential: use digamma-based moment equations. E[X] = (digamma(alpha+1) + euler_gamma) / lambda Var[X] = (pi^2/6 - trigamma(alpha+1)) ... actually use trigamma. Match observed mean and variance numerically to find (alpha, lambda).

Usage

.mom_expexp(x)

MOM via fitdistrplus::fitdist with method = "mme"

Description

Covers: normal, lognormal, gamma, exponential, beta, logistic, uniform. fitdistrplus matches theoretical moments to sample moments analytically for these families.

Usage

.mom_fitdistrplus(x, dist_key, entry)

Match mean, variance, and skewness: skewness = 2/sqrt(shape) => shape = 4/skew^2 scale = sqrt(var/shape) thres = mean - shape*scale

Description

Match mean, variance, and skewness: skewness = 2/sqrt(shape) => shape = 4/skew^2 scale = sqrt(var/shape) thres = mean - shape*scale

Usage

.mom_gamma3(x)

MOM for Generalised Gamma via numerical optimisation. Matches mean, variance, and skewness to theoretical moments. E[X] = b * Gamma((d+1)/p) / Gamma(d/p) Var[X] = b^2 * [Gamma((d+2)/p)/Gamma(d/p) - (Gamma((d+1)/p)/Gamma(d/p))^2]

Description

MOM for Generalised Gamma via numerical optimisation. Matches mean, variance, and skewness to theoretical moments. E[X] = b * Gamma((d+1)/p) / Gamma(d/p) Var[X] = b^2 * [Gamma((d+2)/p)/Gamma(d/p) - (Gamma((d+1)/p)/Gamma(d/p))^2]

Usage

.mom_gengamma(x)

GPD MOM: closed-form, valid when shape < 0.5

Description

Fix loc = min(x) - [U+03B5]; then for the excess Y = x - loc: E[Y] = scale/(1-shape), Var[Y] = scale^2/((1-shape)^2(1-2*shape)) From these: shape = (1 - E[Y]^2/Var[Y]) / 2, scale = E[Y]*(1-shape)

Usage

.mom_gpd(x)

Gumbel MOM: closed-form (Euler-Mascheroni constant gamma_E ~= 0.5772)

Description

E[X] = mu + gamma_E . beta Var[X] = [U+03C0]^2 beta^2 / 6 => beta = sd(x).sqrt(6)/[U+03C0], mu = mean(x) - gamma_E.beta

Usage

.mom_gumbel(x)

E[X] = sigma * sqrt(2/pi) => sigma_hat = mean(x) * sqrt(pi/2)

Description

E[X] = sigma * sqrt(2/pi) => sigma_hat = mean(x) * sqrt(pi/2)

Usage

.mom_halfnormal(x)

Inverse Gamma MOM: closed-form, valid when shape alpha > 2

Description

E[X] = beta/(alpha-1) Var[X] = beta^2/((alpha-1)^2(alpha-2)) => alpha = E[X]^2/Var[X] + 2, beta = E[X](alpha-1)

Usage

.mom_invgamma(x)

MOM for Inverse Gaussian: mu=mean(x), lambda=mu^3/var(x)

Description

MOM for Inverse Gaussian: mu=mean(x), lambda=mu^3/var(x)

Usage

.mom_invgaussian(x)

Kumaraswamy MOM: numerically match first two raw moments. E[X^k; a, b] = b * B(1 + k/a, b) = b * exp(lbeta(1 + k/a, b)) Optimise (a, b) on log scale for numerical stability.

Description

Kumaraswamy MOM: numerically match first two raw moments. E[X^k; a, b] = b * B(1 + k/a, b) = b * exp(lbeta(1 + k/a, b)) Optimise (a, b) on log scale for numerical stability.

Usage

.mom_kumaraswamy(x)

MOM for Laplace: mu=mean(x), b=sd(x)/sqrt(2)

Description

MOM for Laplace: mu=mean(x), b=sd(x)/sqrt(2)

Usage

.mom_laplace(x)

MOM for Lindley: solve for theta from mean equation E[X] = (theta+2)/(theta*(theta+1)) Rearranged: xbar*theta^2 + (xbar-1)*theta - 2 = 0 (same as MLE score)

Description

MOM for Lindley: solve for theta from mean equation E[X] = (theta+2)/(theta*(theta+1)) Rearranged: xbar*theta^2 + (xbar-1)*theta - 2 = 0 (same as MLE score)

Usage

.mom_lindley(x)

Log-Logistic MOM: closed-form solution valid when shape > 2

Description

E[X] = theta [U+03C0]/beta / sin([U+03C0]/beta) Var[X] = theta^2 [2[U+03C0]/beta / sin(2[U+03C0]/beta) - ([U+03C0]/beta)^2 / sin^2([U+03C0]/beta)]

Usage

.mom_loglogistic(x)

Details

Solve for shape beta from the CV, then back-solve for scale theta.

MOM for Lomax (Pareto Type II): alpha and lambda from mean/variance. E[X] = lambda/(alpha-1), Var[X] = alpha*lambda^2/((alpha-1)^2*(alpha-2)) Requires alpha > 2.

Description

MOM for Lomax (Pareto Type II): alpha and lambda from mean/variance. E[X] = lambda/(alpha-1), Var[X] = alpha*lambda^2/((alpha-1)^2*(alpha-2)) Requires alpha > 2.

Usage

.mom_lomax(x)

Log-Pearson III MOM: apply Pearson III MOM to log(x)

Description

Log-Pearson III MOM: apply Pearson III MOM to log(x)

Usage

.mom_lpearson3(x)

Dispatcher for manual MOM estimators

Description

Dispatcher for manual MOM estimators

Usage

.mom_manual(x, dist_key, entry)

MOM for Nakagami: E[X^2] = Omega; m from mean equation. E[X] = Gamma(m+0.5)/Gamma(m) * sqrt(Omega/m) Solve numerically for m given ratio r = mean(x)/sqrt(mean(x^2)).

Description

MOM for Nakagami: E[X^2] = Omega; m from mean equation. E[X] = Gamma(m+0.5)/Gamma(m) * sqrt(Omega/m) Solve numerically for m given ratio r = mean(x)/sqrt(mean(x^2)).

Usage

.mom_nakagami(x)

Pareto MOM: closed-form, valid when shape alpha > 2

Description

E[X] = alphatheta/(alpha-1) Var[X] = alphatheta^2/((alpha-1)^2(alpha-2)) => CV^2 = 1/(alpha-2) => alpha = 1/CV^2 + 2

Usage

.mom_pareto(x)

Pearson III MOM: closed-form via mean, sd, and skewness

Description

Let mu, sigma, gamma be the sample mean, sd, and skewness coefficient. Then: shape = 4/gamma^2, scale = sigmagamma/2, loc = mu - 2sigma/gamma

Usage

.mom_pearson3(x)

E[X] = sigma * sqrt(pi/2) => sigma_hat = mean(x) / sqrt(pi/2)

Description

E[X] = sigma * sqrt(pi/2) => sigma_hat = mean(x) / sqrt(pi/2)

Usage

.mom_rayleigh(x)

MOM for Right Triangular: min=min(x), max=max(x) (closed-form boundaries)

Description

MOM for Right Triangular: min=min(x), max=max(x) (closed-form boundaries)

Usage

.mom_rtriangle(x)

MOM for Student-t location-scale. mu = mean(x), sigma = sd(x)*sqrt((nu-2)/nu) - requires nu > 2. We first estimate nu from the excess kurtosis: kurtosis = 6/(nu-4), nu>4, falling back to nu=5 when kurtosis is too small or not estimable.

Description

MOM for Student-t location-scale. mu = mean(x), sigma = sd(x)*sqrt((nu-2)/nu) - requires nu > 2. We first estimate nu from the excess kurtosis: kurtosis = 6/(nu-4), nu>4, falling back to nu=5 when kurtosis is too small or not estimable.

Usage

.mom_student_t(x)

Triangular MOM

Description

E[X] = (a + b + c) / 3 Var[X] = (a^2 + b^2 + c^2 - ab - ac - bc) / 18

Usage

.mom_triangular(x)

Details

We fix a = min(x), b = max(x) and solve for mode c from the mean equation: c = 3.mean(x) - a - b

Variance constraint is used as a sanity check.

Uniform MOM: a = mean - sqrt(3)*sd, b = mean + sqrt(3)*sd Expanded slightly to ensure all observed data lies within [a, b].

Description

Uniform MOM: a = mean - sqrt(3)*sd, b = mean + sqrt(3)*sd Expanded slightly to ensure all observed data lies within [a, b].

Usage

.mom_uniform(x)

Weibull MOM: solve numerically because moments involve the Gamma function

Description

We solve for shape k such that: CV^2 = Gamma(1+2/k) / Gamma(1+1/k)^2 - 1

Usage

.mom_weibull(x)

Details

Then scale = mean(x) / Gamma(1+1/k).

3-Parameter Weibull MOM

Description

Matches the first three central moments. Because the standardised skewness of a Weibull depends only on the shape parameter k, we solve numerically for k from the sample skewness, then recover scale lambda from the variance and threshold mu from the mean.

Usage

.mom_weibull3p(x)

Details

Moment formulas (letting Gamma_r = Gamma(1 + r/k)): E[X] = mu + lambda Gamma_1 Var[X] = lambda^2 (Gamma_2 - Gamma_1^2) Skew = (Gamma_3 - 3 Gamma_2 Gamma_1 + 2 Gamma_1^3) / (Gamma_2 - Gamma_1^2)^(3/2)

Build the histogram + density overlay panel (publication quality)

Description

Build the histogram + density overlay panel (publication quality)

Usage

.plot_density_panel(fv_obj, top_n = 3L, xlab = "x")

Convert a ggplot2 object to a base64-encoded PNG string

Description

Convert a ggplot2 object to a base64-encoded PNG string

Usage

.plot_to_base64(p, w = 1000, h = 600, res = 120)

Look up a d/p/q/r function for a fit-result list

Description

For GEV and Gumbel, returns FitVerse's own closed-form implementations to avoid evd package-version signature issues. For all other distributions, tries match.fun() first (works for stats and auto-attached packages), then falls back to getExportedValue() for packages in Imports/Suggests (e.g. actuar, sn, mc2d).

Usage

.resolve_dpqr(fit, prefix)

Arguments

Value

A function, or NULL if not found.

Anderson-Darling test (goftest package)

Description

Anderson-Darling test (goftest package)

Usage

.run_ad(x, pcdf)

Chi-squared goodness-of-fit test (Pearson, equiprobable bins)

Description

Uses k equiprobable bins defined by the theoretical quantile function (where available) or equal-width bins as a fallback. Degrees of freedom are adjusted for the number of estimated parameters: df = k - 1 - n_params (minimum 1). Bins with expected count < 1 are merged with their neighbours before the test statistic is computed.

Usage

.run_chisq(x, pcdf, n_params = 0L, q_fun = NULL, params = NULL)

Arguments

Cramer-von Mises test (goftest package)

Description

Cramer-von Mises test (goftest package)

Usage

.run_cvm(x, pcdf)

Kolmogorov-Smirnov test (base R)

Description

Kolmogorov-Smirnov test (base R)

Usage

.run_ks(x, pcdf)

Sample skewness (type 1: biased)

Description

Sample skewness (type 1: biased)

Usage

.skewness(x)

Fast Wakeby CDF via Q-grid interpolation

Description

Pre-computes the quantile function on a fine grid and uses linear interpolation to invert it. This replaces n calls to uniroot (inside lmomco::cdfwak) with a single vectorised quawak call + approx(), giving ~100x speed-up for n = 1000 and making the GoF 10-second window safe.

Usage

.wakeby_cdf_fast(x, para, n_grid = 2001L)

Warn about distributions excluded because they explicitly block a fitting method.

Description

Uses isFALSE() so it only fires for registry entries where the support flag is explicitly FALSE – distributions that simply lack an LMOM fitter are silently skipped, as usual. Up to three excluded distributions get individual messages; more than three get a single summary.

Usage

.warn_unsupported_method(keys, support_field, method_label)

Arguments

Launch the FitVerse Interactive Shiny Application

Description

Opens the FitVerse web interface in the default browser. The app provides a point-and-click interface for uploading data, fitting distributions, inspecting diagnostic plots, and downloading results — without writing any R code.

Usage

FitVerseApp(...)

Arguments

Details

Requires the shiny package (listed in Suggests). An informative error is raised if it is not installed.

The app files are stored in inst/app/ and located at runtime via system.file("app", package = "FitVerse").

Value

Called for its side effect (launches a Shiny app); returns NULL invisibly.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse for the programmatic fitting interface.

Examples

if (interactive()) {
  FitVerseApp()
}

Sample Datasets Bundled with FitVerse

Description

FitVerse includes several comma-separated value (CSV) files in inst/extdata/ that can be used to explore the package's fitting functions. All files are accessed with system.file().

Individual-observation datasets (for fitverse())

annual_max_flows.csv

Annual maximum daily peak flows (m^3/s) for the Rewa River at Navolau, Fiji, covering a 75-year record (1950–2024). The data are right-skewed; Gumbel, GEV, Log-Pearson III, and Log-Normal are natural candidates.

wind_speeds.csv

Daily maximum wind gust speeds (m/s) recorded at Nadi Airport, Fiji, over 200 days (2020–2024). Weibull and Rayleigh distributions are commonly used for wind data.

failure_times.csv

Failure times (hours) from a life test of 150 aluminium electrolytic capacitors under 85 {}^\circC stress. Weibull, Gamma, Log-Normal, and Birnbaum-Saunders are typical reliability choices.

insurance_claims.csv

Residential property insurance claim amounts (in thousands of USD) from 180 claims settled by Pacific Insurance Ltd between 2018 and 2024. Log-Normal, Gamma, and Pareto distributions are common in loss modelling.

concrete.csv

Compressive strength (MPa) and mix-design variables for 1030 concrete specimens (UCI Machine Learning Repository). The Strength column is suitable for fitting; Normal, Gamma, and Weibull are plausible choices.

Grouped / interval-censored datasets (for fitverse_grouped())

rainfall_grouped.csv

Annual maximum daily rainfall (mm) at Suva, Fiji, summarised into a frequency table with columns midpoint and count (80 years of record).

failure_times_grouped.csv

Failure times for 120 electronic components given as an interval frequency table with columns lower, upper, and count.

insurance_losses_grouped.csv

Property insurance claim amounts for 250 claims given as an interval frequency table with columns lower, upper, and count.

Format

Each file is a plain-text CSV with a column-name header row. Lines beginning with # are comments describing the dataset. All numeric values use the standard decimal point.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse for fitting to individual observations.
fitverse_grouped for fitting to grouped / interval data.

Examples

## Annual maximum flows
path <- system.file("extdata", "annual_max_flows.csv", package = "FitVerse")
x <- read.csv(path, comment.char = "#")$peak_flow_m3s
summary(x)


fit <- fitverse(x, xlab = "Peak flow (m3/s)", plot = FALSE)
print(fit)


## Grouped rainfall data
gpath <- system.file("extdata", "rainfall_grouped.csv", package = "FitVerse")
dat   <- read_grouped_csv(gpath)

gfit  <- fitverse_grouped(midpoints = dat$midpoints, counts = dat$counts,
                          xlab = "Rainfall (mm)")
print(gfit)

Extract All Individual Fit Results from a Fitted Object

Description

Retrieves the list of per-distribution fit results stored inside a fitverse object. Each element corresponds to one model that was attempted during the fitting process.

Usage

all_fits(object, ...)
## S3 method for class 'fitverse'
all_fits(object, ...)

Arguments

Value

A named list of fit results, one per distribution attempted.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse, ranking, summary.fitverse

Examples

  set.seed(1)
  x    <- rgamma(200, shape = 2, rate = 0.5)
  fit  <- fitverse(x)
  fits <- all_fits(fit)
  length(fits)

Run GoF tests for a single FitVerse fit-result list and attach results

Description

Run GoF tests for a single FitVerse fit-result list and attach results

Usage

attach_gof(fit_result, x)

Arguments

Value

The same list with an added gof element.

Minimal base64 encoder (avoids base64enc dependency when possible)

Description

Minimal base64 encoder (avoids base64enc dependency when possible)

Usage

base64enc_simple(raw_vec)

Bootstrap Confidence Intervals for Distribution Parameters

Description

Computes non-parametric percentile bootstrap confidence intervals for the parameters of the best-fit distribution in a fitverse object. Optionally also bootstraps confidence intervals for T-year return levels.

Usage

bootstrap_ci(fit,
             B              = 999L,
             conf           = 0.95,
             return_periods = NULL,
             progress       = NULL,
             seed           = NULL,
             verbose        = TRUE)

## S3 method for class 'fitverse_bootstrap'
print(x, ...)

Arguments

Details

Algorithm

The non-parametric percentile bootstrap proceeds as follows:

1. Draw B bootstrap samples of size n with replacement from the observed data.

2. Re-fit the same distribution and estimation method to each sample, using the same internal engine as fitverse().

3. Collect the B \times k matrix of parameter estimates (where k is the number of parameters).

4. For each parameter, report the \alpha/2 and 1 - \alpha/2 empirical percentiles of the bootstrap distribution as the lower and upper confidence bounds.

Replicates that fail to converge are excluded; a warning is issued if more than 10% fail.

Return-level CIs

When return_periods is supplied, the T-year return level is computed for each bootstrap replicate using return_level(), and the percentile CI is reported for each T.

Value

An S3 object of class "fitverse_bootstrap" with the following elements:

params_ci

Data frame with one row per parameter and columns: parameter, estimate, lower, upper, se_boot (bootstrap standard error).

return_level_ci

Data frame (or NULL if return_periods = NULL) with columns: T, return_level, lower, upper.

boot_params

Numeric matrix of dimension B_{eff} \times k: all successful bootstrap parameter estimates. Useful for plotting bootstrap distributions.

B

Integer: number of replicates requested.

B_eff

Integer: number of successful replicates.

conf

Numeric: confidence level.

dist_key

Character: distribution key.

method

Character: estimation method.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

References

Efron, B. and Tibshirani, R. J. (1993). An Introduction to the Bootstrap. Chapman and Hall, New York.

Davison, A. C. and Hinkley, D. V. (1997). Bootstrap Methods and Their Application. Cambridge University Press.

See Also

return_level for point estimates without CIs.
fitverse for the initial fitting step.

Examples

set.seed(42)
x <- rgamma(200, shape = 2, rate = 0.5)
fit <- fitverse(x, dists = "gamma", plot = FALSE, verbose = FALSE)

boot <- bootstrap_ci(fit, B = 499, conf = 0.95, seed = 1, verbose = FALSE)
print(boot)

# Parameter CI data frame
boot$params_ci

# With return-level CIs
boot2 <- bootstrap_ci(fit, B = 499, return_periods = c(10, 50, 100),
                      seed = 1, verbose = FALSE)
boot2$return_level_ci

Summarise all fit results into a ranked comparison table

Description

Summarise all fit results into a ranked comparison table

Usage

build_ranking_table(fits, criterion = "AIC")

Arguments

Value

A data frame ranked by GoF tier then AIC/BIC (GoF composite score as final tiebreaker), with columns: Rank, Distribution, Method, N_params, LogLik, AIC, Delta_AIC, BIC, Delta_BIC, KS_p, AD_p, ChiSq_p, GoF_group, Converged, dist_key

Validate that input is a suitable numeric vector

Description

Stops with an informative message if the data are not numeric, have fewer than 4 finite observations, or violate a stated support constraint.

Usage

check_data(x, dist_name = NULL, require_pos = FALSE, require_01 = FALSE)

Arguments

Value

Invisibly TRUE; stops on failure.

Remove non-finite values from a vector

Description

Issues a message if any values are removed so the user is aware.

Usage

clean_data(x)

Arguments

Value

Numeric vector containing only finite values.

Compare Distribution Fits and Return a Ranking Table

Description

Fits one or more parametric distributions to a numeric vector and returns an ordered ranking table without producing a plot or a full fitverse object. Use this function when you want a quick tabular comparison and do not need the fitted-object infrastructure.

Usage

compare_dists(x,
              dists     = "all",
              method    = "MLE",
              criterion = "AIC",
              xlab      = "x",
              verbose   = TRUE)

Arguments

Details

compare_dists() calls the same internal fitting engine as fitverse() but skips GoF test computation, data-summary computation, and plotting. It is therefore faster and better suited to exploratory screening of many distributions or methods.

The returned data frame can be used directly with kable or datatable for reporting.

Value

A data frame with one row per distribution-method combination, ordered by the chosen information criterion (ascending), with columns:

Rank

Integer rank (1 = best).

Distribution

Full distribution name.

Method

Estimation method ("MLE", "MOM", or "LMOM").

N_params

Number of free parameters.

LogLik

Maximised log-likelihood.

AIC

Akaike Information Criterion.

Delta_AIC

AIC difference relative to the best model.

BIC

Bayesian Information Criterion.

Delta_BIC

BIC difference relative to the best model.

Converged

Logical: did the optimiser converge?

dist_key

Internal distribution key (character).

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse for a full fitted object plus diagnostic plot.

Examples

set.seed(1)
x <- rgamma(200, shape = 3, rate = 0.5)

tbl <- compare_dists(x,
                     dists   = c("gamma", "lognormal", "weibull", "exponential"),
                     method  = "MLE",
                     verbose = FALSE)
head(tbl)

Compute AIC from log-likelihood and parameter count

Description

Compute AIC from log-likelihood and parameter count

Usage

compute_aic(loglik, k)

Arguments

Value

AIC value, or NA if log-likelihood is not finite.

Compute BIC from log-likelihood, parameter count, and sample size

Description

Compute BIC from log-likelihood, parameter count, and sample size

Usage

compute_bic(loglik, k, n)

Arguments

Value

BIC value, or NA if log-likelihood is not finite.

Compute a Descriptive Statistics Summary

Description

compute_data_summary() computes a comprehensive set of descriptive statistics for a numeric vector and returns them as a named list. It is called automatically by fitverse() and summary.fitverse(), but can also be used independently for exploratory data analysis.

Usage

compute_data_summary(x, xlab = "x", n_raw = NULL)

Arguments

Details

Skewness and kurtosis are computed using the bias-corrected (unbiased) estimators:

g_1 = \frac{n}{(n-1)(n-2)} \sum_{i=1}^n \left(\frac{x_i - \bar{x}}{s}\right)^3

g_2 = \frac{n(n+1)}{(n-1)(n-2)(n-3)} \sum_{i=1}^n \left(\frac{x_i - \bar{x}}{s}\right)^4 - \frac{3(n-1)^2}{(n-2)(n-3)}

These correspond to the skewness and kurtosis functions in the moments package (type 3 in e1071) and are consistent with the convention in Applied Statistics textbooks. The Normal distribution has skewness 0 and excess kurtosis 0.

Value

compute_data_summary() returns a named list with the following elements:

n

Integer: number of finite (non-missing) observations used.

n_missing

Integer or NA: number of non-finite values removed.

mean

Arithmetic mean.

sd

Standard deviation (denominator n-1).

cv

Coefficient of variation (sd / mean). NA when mean == 0.

median

Sample median.

iqr

Interquartile range (Q_{0.75} - Q_{0.25}).

skewness

Sample skewness (Fisher's g_1).

kurtosis

Excess kurtosis (Fisher's g_2).

min

Sample minimum.

max

Sample maximum.

range

Range (max - min).

q25

25th percentile.

q75

75th percentile.

xlab

Character: the variable label.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse (calls this automatically), summary.fitverse.

Examples

set.seed(3)
x <- rlnorm(500, meanlog = 1, sdlog = 0.5)

ds <- compute_data_summary(x, xlab = "Claim amount", n_raw = 510)
ds$skewness
ds$cv

Compute four goodness-of-fit tests for a fitted distribution

Description

Compute four goodness-of-fit tests for a fitted distribution

Usage

compute_gof(x, dist_key, params)

Arguments

Value

A data frame with columns: Test, Statistic, P_value, Stars (or NA rows if a test cannot be computed).

Density, CDF, Quantile, and Random Sampling for a Fitted Distribution

Description

Convenience functions that evaluate the probability density function (PDF), cumulative distribution function (CDF), quantile function (inverse CDF), and random-number generator for the best-fit distribution stored in a fitverse object.

These follow the standard R d/p/q/r naming convention and delegate to the underlying d*, p*, q*, and r* functions for the best-fit distribution.

Usage

dfitverse(x,   fit, log = FALSE)
pfitverse(q,   fit, lower.tail = TRUE, log.p = FALSE)
qfitverse(p,   fit, lower.tail = TRUE, log.p = FALSE)
rfitverse(n,   fit)

Arguments

Details

The correct d/p/q/r functions are resolved via .lookup_dpqr(), which follows a two-step priority:

1. Internal closed-form implementations (highest priority): GEV, Gumbel, GPD, Pearson III, Log-Pearson III, and GLO all use FitVerse's own vectorised implementations registered in .FITVERSE_DPQR_OVERRIDES. This avoids any package-version signature issues with evd and provides full support for the new hydrological distributions that have no dedicated CRAN d/p/q/r package.

2. Package namespace lookup: For all other distributions, stats functions are found via match.fun(); Imports/Suggests packages (actuar, sn, mc2d) are resolved via getExportedValue().

The estimated parameters are extracted from fit$best_fit$params and passed as named arguments to the underlying function. This is equivalent to calling, for example, dgamma(x, shape = fit$best_fit$params["shape"], rate = fit$best_fit$params["rate"]).

Value

dfitverse()

Numeric vector of density values.

pfitverse()

Numeric vector of probabilities.

qfitverse()

Numeric vector of quantiles.

rfitverse()

Numeric vector of random samples.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

return_level for T-year return levels.
fitverse for the fitting step.

Examples

set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)
fit <- fitverse(x, dists = "gamma", plot = FALSE, verbose = FALSE)

# Density at x = 2
dfitverse(2, fit)

# Probability that X <= 5
pfitverse(5, fit)

# 90th percentile
qfitverse(0.90, fit)

# Generate 10 random samples
rfitverse(10, fit)

Custom d/p/q/r Functions for Additional Distributions

Description

Density, CDF, quantile, and random-generation functions for distribution families that are either unavailable on CRAN, parameterised differently from existing packages, or implemented from scratch for performance and consistency within FitVerse.

The following families are covered:

Student-t (location-scale)

Three-parameter Student-t with location \mu, scale \sigma > 0, and degrees of freedom \nu > 0: dstudent_t, pstudent_t, qstudent_t, rstudent_t.

Laplace

Symmetric double-exponential with location \mu and scale b > 0: dlaplace, plaplace, qlaplace, rlaplace.

Inverse Gaussian

Two-parameter (mean \mu > 0, shape \lambda > 0): dinvgaussian, pinvgaussian, qinvgaussian, rinvgaussian.

Lomax (Pareto Type II)

Shape \alpha > 0 and scale \lambda > 0: dlomax, plomax, qlomax, rlomax.

Generalised Gamma

Three-parameter generalisation of the Gamma (shape p > 0, scale d > 0, power b > 0): dgengamma, pgengamma, qgengamma, rgengamma.

Gompertz

Parameters \alpha > 0 and \beta > 0: dgompertz, pgompertz, qgompertz, rgompertz.

Birnbaum-Saunders (Fatigue-Life)

Shape \alpha > 0 and scale \beta > 0: dbirnbaum_saunders, pbirnbaum_saunders, qbirnbaum_saunders, rbirnbaum_saunders.

Rayleigh

Scale \sigma > 0: drayleigh, prayleigh, qrayleigh, rrayleigh.

Half-Normal

Scale \sigma > 0: dhalfnormal, phalfnormal, qhalfnormal, rhalfnormal.

Exponentiated Weibull

Shape parameters shape1 > 0 and shape2 > 0, scale scale > 0: dexpweibull, pexpweibull, qexpweibull, rexpweibull.

Three-Parameter Gamma

Gamma with location threshold \gamma (shape \alpha > 0, scale \beta > 0, threshold \gamma): dgamma3, pgamma3, qgamma3, rgamma3.

Three-Parameter Log-Normal

Log-Normal with threshold \gamma (meanlog \mu, sdlog \sigma > 0, threshold \gamma): dlnorm3, plnorm3, qlnorm3, rlnorm3.

Right Triangular

Special triangular with min and max, where the mode equals max: drtriangle, prtriangle, qrtriangle, rrtriangle.

Kumaraswamy

Bounded distribution with shape parameters shape1 > 0 and shape2 > 0 on (0, 1): dkumaraswamy, pkumaraswamy, qkumaraswamy, rkumaraswamy.

Generalised Beta of the Second Kind (GB2)

Four-parameter family with shape1 > 0, scale > 0, shape2 > 0, shape3 > 0: dgb2, pgb2, qgb2, rgb2.

Kappa Four-Parameter

Generalised extreme-value family with location xi, scale alpha > 0, and shape parameters kappa and h: dkap4, pkap4, qkap4, rkap4.

Lindley

One-parameter distribution with theta > 0: dlindley, plindley, qlindley, rlindley.

Generalised Normal

Three-parameter distribution with location mu, scale sigma > 0, and shape beta > 0: dgnorm, pgnorm, qgnorm, rgnorm.

Dagum

Three-parameter Burr Type III with shape1 > 0, scale > 0, shape2 > 0: ddagum, pdagum, qdagum, rdagum.

Johnson SU

Four-parameter unbounded system with gamma, delta > 0, xi, lambda > 0: djohnsonsu, pjohnsonsu, qjohnsonsu, rjohnsonsu.

Asymmetric Laplace

Three-parameter asymmetric Laplace with location mu, scale sigma > 0, and asymmetry kappa > 0: dalaplace, palaplace, qalaplace, ralaplace.

Exponentiated Exponential

Two-parameter family with shape alpha > 0 and rate lambda > 0: dexpexp, pexpexp, qexpexp, rexpexp.

Normal-Inverse Gaussian (NIG)

Four-parameter family with location mu, spread delta > 0, steepness alpha > 0, and asymmetry beta: dnig, pnig, qnig, rnig.

Wakeby

Five-parameter quantile-defined distribution with location xi, scale parameters alpha and beta, and shape parameters gamma and delta: dwakeby, pwakeby, qwakeby, rwakeby.

Generalised Hyperbolic (GHyp)

Five-parameter family with lambda, alpha > 0, beta, delta > 0, mu: dghyp, pghyp, qghyp, rghyp.

Truncated Normal

Normal distribution truncated to [a, b] with mean and sd > 0: dtruncnorm, ptruncnorm, qtruncnorm, rtruncnorm.

Nakagami

Two-parameter fading distribution with shape m >= 0.5 and spread Omega > 0: dnakagami, pnakagami, qnakagami, rnakagami.

Usage

dstudent_t(x,   nu, mu, sigma, log = FALSE)
pstudent_t(q,   nu, mu, sigma)
qstudent_t(p,   nu, mu, sigma)
rstudent_t(n,   nu, mu, sigma)

dlaplace(x,   mu, b, log = FALSE)
plaplace(q,   mu, b)
qlaplace(p,   mu, b)
rlaplace(n,   mu, b)

dinvgaussian(x,   mu, lambda, log = FALSE)
pinvgaussian(q,   mu, lambda)
qinvgaussian(p,   mu, lambda)
rinvgaussian(n,   mu, lambda)

dlomax(x,   alpha, lambda, log = FALSE)
plomax(q,   alpha, lambda)
qlomax(p,   alpha, lambda)
rlomax(n,   alpha, lambda)

dgengamma(x,   p, d, b, log = FALSE)
pgengamma(q,   p, d, b)
qgengamma(prob, p, d, b)
rgengamma(n,   p, d, b)

dgompertz(x,   alpha, beta, log = FALSE)
pgompertz(q,   alpha, beta)
qgompertz(p,   alpha, beta)
rgompertz(n,   alpha, beta)

dbirnbaum_saunders(x,   alpha, beta, log = FALSE)
pbirnbaum_saunders(q,   alpha, beta)
qbirnbaum_saunders(p,   alpha, beta)
rbirnbaum_saunders(n,   alpha, beta)

drayleigh(x,   sigma, log = FALSE)
prayleigh(q,   sigma, lower.tail = TRUE, log.p = FALSE)
qrayleigh(p,   sigma, lower.tail = TRUE, log.p = FALSE)
rrayleigh(n,   sigma)

dhalfnormal(x,   sigma, log = FALSE)
phalfnormal(q,   sigma, lower.tail = TRUE, log.p = FALSE)
qhalfnormal(p,   sigma, lower.tail = TRUE, log.p = FALSE)
rhalfnormal(n,   sigma)

dexpweibull(x,   shape1, scale, shape2, log = FALSE)
pexpweibull(q,   shape1, scale, shape2, lower.tail = TRUE, log.p = FALSE)
qexpweibull(p,   shape1, scale, shape2, lower.tail = TRUE, log.p = FALSE)
rexpweibull(n,   shape1, scale, shape2)

dgamma3(x,   shape, scale, thres, log = FALSE)
pgamma3(q,   shape, scale, thres, lower.tail = TRUE, log.p = FALSE)
qgamma3(p,   shape, scale, thres, lower.tail = TRUE, log.p = FALSE)
rgamma3(n,   shape, scale, thres)

dlnorm3(x,   meanlog, sdlog, thres, log = FALSE)
plnorm3(q,   meanlog, sdlog, thres, lower.tail = TRUE, log.p = FALSE)
qlnorm3(p,   meanlog, sdlog, thres, lower.tail = TRUE, log.p = FALSE)
rlnorm3(n,   meanlog, sdlog, thres)

drtriangle(x,   min, max, log = FALSE)
prtriangle(q,   min, max)
qrtriangle(p,   min, max)
rrtriangle(n,   min, max)

dkumaraswamy(x,   shape1, shape2, log = FALSE)
pkumaraswamy(q,   shape1, shape2)
qkumaraswamy(p,   shape1, shape2)
rkumaraswamy(n,   shape1, shape2)

dgb2(x, shape1, scale, shape2, shape3, log = FALSE)
pgb2(q, shape1, scale, shape2, shape3, lower.tail = TRUE, log.p = FALSE)
qgb2(p, shape1, scale, shape2, shape3, lower.tail = TRUE, log.p = FALSE)
rgb2(n, shape1, scale, shape2, shape3)

dkap4(x, xi, alpha, kappa, h, log = FALSE)
pkap4(q, xi, alpha, kappa, h, lower.tail = TRUE, log.p = FALSE)
qkap4(p, xi, alpha, kappa, h, lower.tail = TRUE, log.p = FALSE)
rkap4(n, xi, alpha, kappa, h)

dlindley(x, theta, log = FALSE)
plindley(q, theta, lower.tail = TRUE, log.p = FALSE)
qlindley(p, theta, lower.tail = TRUE, log.p = FALSE)
rlindley(n, theta)

dgnorm(x, mu, sigma, beta, log = FALSE)
pgnorm(q, mu, sigma, beta, lower.tail = TRUE, log.p = FALSE)
qgnorm(p, mu, sigma, beta, lower.tail = TRUE, log.p = FALSE)
rgnorm(n, mu, sigma, beta)

ddagum(x, shape1, scale, shape2, log = FALSE)
pdagum(q, shape1, scale, shape2, lower.tail = TRUE, log.p = FALSE)
qdagum(p, shape1, scale, shape2, lower.tail = TRUE, log.p = FALSE)
rdagum(n, shape1, scale, shape2)

djohnsonsu(x, gamma, delta, xi, lambda, log = FALSE)
pjohnsonsu(q, gamma, delta, xi, lambda, lower.tail = TRUE, log.p = FALSE)
qjohnsonsu(p, gamma, delta, xi, lambda, lower.tail = TRUE, log.p = FALSE)
rjohnsonsu(n, gamma, delta, xi, lambda)

dalaplace(x, mu, sigma, kappa, log = FALSE)
palaplace(q, mu, sigma, kappa, lower.tail = TRUE, log.p = FALSE)
qalaplace(p, mu, sigma, kappa, lower.tail = TRUE, log.p = FALSE)
ralaplace(n, mu, sigma, kappa)

dexpexp(x, alpha, lambda, log = FALSE)
pexpexp(q, alpha, lambda, lower.tail = TRUE, log.p = FALSE)
qexpexp(p, alpha, lambda, lower.tail = TRUE, log.p = FALSE)
rexpexp(n, alpha, lambda)

dnig(x, mu, delta, alpha, beta, log = FALSE)
pnig(q, mu, delta, alpha, beta, lower.tail = TRUE, log.p = FALSE)
qnig(p, mu, delta, alpha, beta, lower.tail = TRUE, log.p = FALSE)
rnig(n, mu, delta, alpha, beta)

dwakeby(x, xi, alpha, beta, gamma, delta, log = FALSE)
pwakeby(q, xi, alpha, beta, gamma, delta, lower.tail = TRUE, log.p = FALSE)
qwakeby(p, xi, alpha, beta, gamma, delta, lower.tail = TRUE, log.p = FALSE)
rwakeby(n, xi, alpha, beta, gamma, delta)

dghyp(x, lambda, alpha, beta, delta, mu, log = FALSE)
pghyp(q, lambda, alpha, beta, delta, mu, lower.tail = TRUE, log.p = FALSE)
qghyp(p, lambda, alpha, beta, delta, mu, lower.tail = TRUE, log.p = FALSE)
rghyp(n, lambda, alpha, beta, delta, mu)

dtruncnorm(x, mean, sd, a, b, log = FALSE)
ptruncnorm(q, mean, sd, a, b, lower.tail = TRUE, log.p = FALSE)
qtruncnorm(p, mean, sd, a, b, lower.tail = TRUE, log.p = FALSE)
rtruncnorm(n, mean, sd, a, b)

dnakagami(x, m, Omega, log = FALSE)
pnakagami(q, m, Omega, lower.tail = TRUE, log.p = FALSE)
qnakagami(p, m, Omega, lower.tail = TRUE, log.p = FALSE)
rnakagami(n, m, Omega)

Arguments

Value

A numeric vector of densities, probabilities, quantiles, or random samples, consistent with the standard R d/p/q/r convention.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse to fit these distributions to data; fitverse_distributions for the full distribution catalogue; dfitverse for evaluating the fitted best-fit distribution.

Examples

## Rayleigh
curve(drayleigh(x, sigma = 2), from = 0, to = 10,
      main = "Rayleigh(sigma=2)", ylab = "Density")

## Kumaraswamy (bounded to (0,1))
curve(dkumaraswamy(x, shape1 = 2, shape2 = 5), from = 0, to = 1,
      main = "Kumaraswamy(2, 5)", ylab = "Density")

## Birnbaum-Saunders round-trip
p <- c(0.1, 0.5, 0.9)
q <- qbirnbaum_saunders(p, alpha = 0.5, beta = 100)
pbirnbaum_saunders(q, alpha = 0.5, beta = 100)   # recovers p

## Inverse Gaussian random samples
set.seed(1)
x <- rinvgaussian(200, mu = 5, lambda = 10)
hist(x, breaks = 20, main = "Inverse Gaussian(5, 10)")

Three-Parameter Weibull Distribution

Description

Density, cumulative distribution function, quantile function, and random generation for the three-parameter Weibull distribution with shape \beta > 0, scale \alpha > 0, and threshold (location) \gamma. The support is x > \gamma.

Usage

dweibull3p(x,     shape, scale, thres, log = FALSE)
pweibull3p(q,     shape, scale, thres, lower.tail = TRUE, log.p = FALSE)
qweibull3p(p,     shape, scale, thres, lower.tail = TRUE, log.p = FALSE)
rweibull3p(n,     shape, scale, thres)

Arguments

Details

The three-parameter Weibull is a location-shifted version of the standard two-parameter Weibull. Let Y = X - \gamma; then Y \sim \mathrm{Weibull}(\beta, \alpha). The PDF for x > \gamma is

f(x) = \frac{\beta}{\alpha} \left(\frac{x - \gamma}{\alpha}\right)^{\beta - 1} \exp\!\left[-\left(\frac{x - \gamma}{\alpha}\right)^{\beta}\right].

Internally these functions shift x by thres and delegate to dweibull, pweibull, qweibull, and rweibull.

Value

dweibull3p

Numeric vector of density values.

pweibull3p

Numeric vector of probabilities.

qweibull3p

Numeric vector of quantiles.

rweibull3p

Numeric vector of random samples.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

dweibull for the two-parameter Weibull; fitverse with dists = "weibull3" to fit this distribution by MLE.

Examples

## PDF and CDF
curve(dweibull3p(x, shape = 2, scale = 10, thres = 5),
      from = 5, to = 35, ylab = "Density",
      main = "Weibull(2, 10) shifted by 5")

## Quantile / CDF round-trip
p  <- c(0.1, 0.5, 0.9)
q  <- qweibull3p(p, shape = 2, scale = 10, thres = 5)
pweibull3p(q, shape = 2, scale = 10, thres = 5)   # should recover p

## Random samples
set.seed(1)
x <- rweibull3p(200, shape = 1.8, scale = 50, thres = 10)
hist(x, breaks = 20, main = "Weibull(1.8, 50) with threshold 10")

Filter distribution keys to those compatible with the observed data

Description

Distributions whose support is violated (e.g. positive-only distribution applied to data with negatives) are silently excluded.

Usage

filter_by_support(x, keys)

Arguments

Value

Filtered character vector.

Create a unique key for a fit result: "dist_key_METHOD"

Description

Create a unique key for a fit result: "dist_key_METHOD"

Usage

fit_list_key(dist_key, method)

Fit a single distribution by L-Moments

Description

Fit a single distribution by L-Moments

Usage

fit_lmom(x, dist_key)

Arguments

Value

Named list with fit results, or NULL if unsupported/failed.

Fit all eligible distributions by L-Moments

Description

Fit all eligible distributions by L-Moments

Usage

fit_lmom_all(x, keys)

Arguments

Value

Named list of fit results (NULLs dropped).

Fit a single distribution by Maximum Likelihood Estimation

Description

Fit a single distribution by Maximum Likelihood Estimation

Usage

fit_mle(x, dist_key)

Arguments

Value

Named list with fit results (see file header), or NULL on failure.

Fit a list of distributions by MLE and collect results

Description

Fit a list of distributions by MLE and collect results

Usage

fit_mle_all(x, keys)

Arguments

Value

Named list of fit-result lists (NULL entries for failures).

Fit a single distribution by Method of Moments

Description

Fit a single distribution by Method of Moments

Usage

fit_mom(x, dist_key)

Arguments

Value

Named list with fit results, or NULL if MOM is not available.

Fit all eligible distributions by Method of Moments

Description

Fit all eligible distributions by Method of Moments

Usage

fit_mom_all(x, keys)

Arguments

Value

Named list of fit results (NULLs dropped).

Fit Probability Distributions to Univariate Continuous Data

Description

The primary function of the FitVerse package. Given a numeric data vector, fitverse() fits one or more parametric probability distributions, ranks them by an information criterion, selects the best-fitting model, computes goodness-of-fit statistics, and optionally produces a publication-quality diagnostic plot.

Three estimation methods are supported: Maximum Likelihood Estimation ("MLE"), Method of Moments ("MOM"), and L-Moments ("LMOM"). When multiple methods are requested, all results are combined into a single ranking table and the global best model is returned.

Tip: For a point-and-click experience, launch the built-in interactive Shiny application with FitVerseApp().

Usage

fitverse(x,
         dists       = "all",
         method      = "MLE",
         criterion   = "AIC",
         xlab        = "x",
         plot        = TRUE,
         interactive = FALSE,
         top_n       = 3L,
         progress    = NULL,
         verbose     = TRUE,
         ...)

Arguments

Details

Workflow

fitverse() executes the following steps:

1. Validate and clean – removes non-finite values; checks that at least 4 observations remain.

2. Filter distributions – distributions whose support is incompatible with the data (e.g. positive-only Gamma when data contains negative values) are silently excluded.

3. Fit – calls fit_mle_all(), fit_mom_all(), and/or fit_lmom_all() according to method.

4. Goodness-of-fit – attaches Kolmogorov-Smirnov (KS), Anderson-Darling (AD), and Cramer-von Mises (CvM) test results to each fit result.

5. Rank – builds a table ordered by criterion (ascending AIC or BIC).

6. Select – chooses the model with the lowest AIC (or BIC).

7. Summarise – computes compute_data_summary() and stores it in the returned object.

8. Plot (optional) – calls plot.fitverse().

Estimation methods

MLE

Maximises the log-likelihood \ell(\boldsymbol\theta) = \sum_{i=1}^n \log f(x_i; \boldsymbol\theta). Implemented via fitdistrplus::fitdist() with bounded L-BFGS-B optimisation where parameter constraints exist.

MOM

Equates the first k sample moments to their theoretical counterparts and solves for the k parameters. Closed-form solutions are used where available (Normal, Gamma, Log-Normal, Exponential, Beta, Logistic, Uniform, Gumbel, Pareto, Inverse Gamma); numerical solution via stats::optimise() for Weibull and Log-Logistic. fitdistrplus MME is used for standard stats distributions.

LMOM

Equates sample L-moments (linear combinations of order statistics) to their theoretical counterparts. L-moments are more robust to outliers than conventional moments and are widely used in hydrology. Implemented via lmomco; available for Normal, Gamma, Exponential, Log-Normal, Weibull, GEV, GLO, Gumbel, Logistic, and Log-Logistic. Reference: Hosking (1990).

Model selection

The best-fit model is the one with the lowest AIC (or BIC) across all fitted distributions and methods. Delta AIC values in the ranking table quantify evidence against each model relative to the best: models with Delta AIC < 2 have "substantial support", 4–7 "considerably less", and > 10 "essentially none" (Burnham and Anderson, 2002).

Goodness-of-fit tests

Three non-parametric tests are applied to the best-fit model:

Kolmogorov-Smirnov (KS)

Measures the maximum absolute difference between the empirical and fitted CDFs. Conservative; not recommended as the primary criterion.

Anderson-Darling (AD)

Gives extra weight to the distribution tails. More powerful than KS for detecting tail departures. Implemented via goftest.

Cramer-von Mises (CvM)

Integrates the squared difference between empirical and fitted CDFs. Implemented via goftest.

A large p-value (> 0.05) indicates insufficient evidence to reject the fitted distribution.

Value

An S3 object of class "fitverse". It is a named list with the following elements:

best_fit

Named list describing the best-fitting model, with elements:

dist_key

Character: distribution key (e.g. "gamma").

method

Character: estimation method ("MLE", "MOM", or "LMOM").

params

Named numeric vector of parameter estimates.

loglik

Numeric: maximised log-likelihood.

aic

Numeric: Akaike Information Criterion.

bic

Numeric: Bayesian Information Criterion.

n_params

Integer: number of estimated parameters.

n

Integer: sample size used.

converged

Logical: did the optimiser converge?

gof

Data frame of goodness-of-fit test results (Test, Statistic, P_value, Stars), or NULL if unavailable.

all_fits

Named list of all fit-result lists, one per distribution-method combination. Keys are of the form "distkey_METHOD" (e.g. "gamma_MLE").

ranking

Data frame with columns: Rank, Distribution, Method, N_params, LogLik, AIC, Delta_AIC, BIC, Delta_BIC, KS_p, AD_p, Converged, dist_key. Ordered by the chosen criterion (ascending).

data

Numeric vector: the cleaned input data.

n

Integer: number of finite observations used.

xlab

Character: the variable label.

data_summary

Named list from compute_data_summary(): n, n_missing, mean, sd, cv, median, iqr, skewness, kurtosis, min, max, range, q25, q75.

criterion

Character: "AIC" or "BIC".

method

Character vector: method(s) used.

call

The matched function call.

The object is returned invisibly. Use print(), summary(), plot(), coef(), logLik(), AIC(), and BIC() for standard S3 method access.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

References

Burnham, K. P. and Anderson, D. R. (2002). Model Selection and Multimodel Inference (2nd ed.). Springer.

Delignette-Muller, M. L. and Dutang, C. (2015). fitdistrplus: An R package for fitting distributions. Journal of Statistical Software, 64(4), 1–34. doi:10.18637/jss.v064.i04

Hosking, J. R. M. (1990). L-moments: Analysis and estimation of distributions using linear combinations of order statistics. Journal of the Royal Statistical Society, Series B, 52(1), 105–124.

Akaike, H. (1974). A new look at the statistical model identification. IEEE Transactions on Automatic Control, 19(6), 716–723.

See Also

compare_dists for a ranking table only (no plot).
fitverse_distributions to list all supported distributions.
plot.fitverse, summary.fitverse, bootstrap_ci, return_level, generate_report.

Examples

## -----------------------------------------------------------
## Example 1: Basic MLE fit -- three distributions
## -----------------------------------------------------------
set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)

fit <- fitverse(x,
                dists   = c("gamma", "lognormal", "weibull"),
                method  = "MLE",
                xlab    = "Waiting time (s)",
                plot    = FALSE,
                verbose = FALSE)

print(fit)
coef(fit)        # parameter estimates
AIC(fit)         # AIC of the best-fit model
logLik(fit)      # log-likelihood

## -----------------------------------------------------------
## Example 2: All three estimation methods, ranked by BIC
## -----------------------------------------------------------

fit2 <- fitverse(x,
                 method    = c("MLE", "MOM", "LMOM"),
                 criterion = "BIC",
                 plot      = FALSE,
                 verbose   = FALSE)
summary(fit2)


## -----------------------------------------------------------
## Example 3: Automatic best-fit selection across all distributions
## -----------------------------------------------------------

set.seed(7)
x3 <- rlnorm(500, meanlog = 1.5, sdlog = 0.6)

fit3 <- fitverse(x3,
                 xlab    = "Insurance claim size",
                 method  = "MLE",
                 plot    = FALSE,
                 verbose = FALSE)
fit3$ranking[1:5, c("Distribution","Method","AIC","Delta_AIC","AD_p")]


## -----------------------------------------------------------
## Example 4: Accessing the data summary
## -----------------------------------------------------------

fit4 <- fitverse(x, plot = FALSE, verbose = FALSE)
fit4$data_summary$mean
fit4$data_summary$skewness

Fit Distributions to All Columns of a Data Frame

Description

Applies fitverse() to every numeric column of a data frame (or matrix, or named list), collects the results, and returns a tidy summary table alongside the full list of fitted objects.

Use fitverse_batch() for automated data characterisation workflows, survey-population profiling, or as a pre-processing step before optimal stratification with stratifyR.

Usage

fitverse_batch(data,
               cols      = NULL,
               dists     = "all",
               method    = "MLE",
               criterion = "AIC",
               top_n     = 3L,
               plot      = FALSE,
               progress  = NULL,
               verbose   = TRUE)

## S3 method for class 'fitverse_batch'
print(x, ...)

## S3 method for class 'fitverse_batch'
summary(object, ...)

Arguments

Details

fitverse_batch() iterates over columns with a tryCatch() wrapper so that a fitting failure for one column does not abort the entire batch. Failed columns are collected in the failed element of the returned object.

Column names containing spaces or special characters are handled correctly via make.names() internally.

Value

An S3 object of class "fitverse_batch" with elements:

fits

Named list of "fitverse" objects, one per successfully fitted column.

summary

Data frame with one row per column, columns: variable, n, best_dist, method, AIC, BIC, loglik, KS_p, AD_p.

failed

Character vector of column names for which fitting failed, or character(0) if all succeeded.

criterion

Character: the criterion used.

call

The matched call.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse for single-column fitting.
generate_report for producing a report from batch results.

Examples

# Fit all numeric columns of the airquality dataset
data(airquality)
batch <- fitverse_batch(airquality,
                        method  = "MLE",
                        verbose = FALSE)
print(batch)
summary(batch)

# Access the fit for the Ozone column
batch$fits[["Ozone"]]

Fit Parametric Distributions to Censored or Truncated Data

Description

Estimates distribution parameters by maximising the censored log-likelihood for data that may be right-censored, left-censored, interval-censored, and/or left-/right-truncated, or any combination of these. Seven distribution families are supported. Fits are ranked by AIC or BIC, and a nonparametric Kaplan–Meier estimate is overlaid automatically on the survival or CDF plot.

Usage

fitverse_censored(x,
                  status      = NULL,
                  x_upper     = NULL,
                  trunc_lower = -Inf,
                  trunc_upper =  Inf,
                  dists       = c("normal", "lognormal", "weibull",
                                  "exponential", "gamma",
                                  "loglogistic", "gumbel"),
                  criterion   = c("AIC", "BIC"),
                  n_starts    = 3L,
                  xlab        = "x",
                  draw        = TRUE)

## S3 method for class 'fitverse_censored'
print(x, ...)

## S3 method for class 'fitverse_censored'
plot(x, top_n = 3L, ...)

Arguments

Value

An object of class "fitverse_censored", which is a named list with the following elements:

fits

List of fit objects ordered by criterion (best first). Each element is a list with: dist_key (character), display (character display name), par (named numeric vector of parameter estimates), log_lik (numeric), aic (numeric), bic (numeric), n_par (integer), converged (logical).

rankings

Data frame with columns Rank, Distribution, LogLik, AIC, BIC, Converged, sorted by criterion.

x, status, x_upper

Cleaned copies of the input data vectors (non-finite values removed).

trunc_lower, trunc_upper

Truncation bounds as supplied.

xlab

Variable label.

criterion

Ranking criterion used ("AIC" or "BIC").

References

Meeker, W. Q. and Escobar, L. A. (1998). Statistical Methods for Reliability Data. Wiley, New York.

Lawless, J. F. (2003). Statistical Models and Methods for Lifetime Data, 2nd ed. Wiley, Hoboken.

Kaplan, E. L. and Meier, P. (1958). Nonparametric estimation from incomplete observations. Journal of the American Statistical Association, 53(282), 457–481. doi:10.1080/01621459.1958.10501452

See Also

fitverse for uncensored data fitting; fitverse_distributions for the full distribution library.

Examples

## Right-censored survival data (Weibull true model)
set.seed(1)
times  <- rweibull(100, shape = 1.8, scale = 50)
censor <- rexp(100, rate = 1/60)
x      <- pmin(times, censor)
status <- as.integer(times <= censor)   # 1 = event, 0 = right-censored

fit_s <- fitverse_censored(x, status,
                           dists = c("weibull", "lognormal", "gamma"),
                           xlab  = "Survival time (hours)")
print(fit_s)
plot(fit_s, top_n = 3)


## Left-censored environmental data (below detection limit)
set.seed(2)
conc       <- rlnorm(80, meanlog = 1.5, sdlog = 0.7)
dl         <- 3.0
x_env      <- ifelse(conc >= dl, conc, dl)
status_env <- ifelse(conc >= dl, 1L, 2L)

fit_e <- fitverse_censored(x_env, status_env,
                           dists = c("lognormal", "gamma"),
                           xlab  = "Concentration (mg/L)")
fit_e$rankings


## Left-truncated occupational exposure data
set.seed(3)
exposure <- rlnorm(300, meanlog = 2.5, sdlog = 0.5)
exposure <- exposure[exposure > 5]   # only workers above LOQ recorded

fit_t <- fitverse_censored(exposure, trunc_lower = 5,
                           dists = c("lognormal", "gamma", "weibull"),
                           xlab  = "Exposure (mg/m3)")
fit_t$fits[[1L]]$par

List All Distributions Supported by FitVerse

Description

Returns a data frame describing all 52 probability distributions that FitVerse can fit, including their distribution keys, full names, number of parameters, support, available estimation methods, and the underlying R function source.

This page also documents the mathematical details of the four hydrological / extreme-value distributions: Generalised Pareto (GPD), Pearson Type III (P3), Log-Pearson Type III (LP3), and Generalised Logistic (GLO).

Usage

fitverse_distributions(method = NULL)

Arguments

Details

Complete distribution table (52 families)

Hydrological and extreme-value distributions (v1.0 additions)

The four distributions below – GPD, Pearson III, Log-Pearson III, and GLO – are central to flood-frequency analysis and Peaks-Over-Threshold modelling. They share a common implementation strategy: closed-form d/p/q/r functions are registered in the internal .FITVERSE_DPQR_OVERRIDES table, bypassing any external package dependency and guaranteeing consistent behaviour across platforms.

Generalised Pareto Distribution (GPD) – key: "gpd"

Parameters: \mu (location / threshold), \sigma > 0 (scale), \xi (shape).

PDF (\xi \neq 0):

f(x\,;\,\mu,\sigma,\xi)=\frac{1}{\sigma} \Bigl(1+\xi\,\frac{x-\mu}{\sigma}\Bigr)^{-1/\xi-1}

Support: x \ge \mu when \xi \ge 0; \mu \le x \le \mu - \sigma/\xi when \xi < 0. Reduces to \mathrm{Exp}(1/\sigma) when \xi = 0.

Quantile function:

Q(p)=\mu+\frac{\sigma}{\xi}\bigl[(1-p)^{-\xi}-1\bigr],\quad\xi\neq 0 \qquad Q(p)=\mu-\sigma\log(1-p),\quad\xi=0

Method of Moments: Fix \hat\mu=\min(x)-\varepsilon; let y_i=x_i-\hat\mu. Then \hat\xi=(1-\bar y^2/s_y^2)/2, \hat\sigma=\bar y(1-\hat\xi). Valid only when \hat\xi<1/2.

L-Moments: via lmomco::pargpa(). Sign convention: FitVerse \xi=-\kappa_{\text{lmomco}}.

Use case: Peaks-Over-Threshold analysis; insurance excess-of-loss modelling (Pickands, 1975).

Pearson Type III (P3) – key: "pearson3"

Parameters: \alpha>0 (shape), \beta (scale; may be negative), \zeta (location).

Relationship to Gamma:

X\sim\mathrm{P3}(\alpha,\beta,\zeta) \;\Leftrightarrow\; (X-\zeta)/\beta\sim\mathrm{Gamma}(\alpha,\,\mathrm{rate}=1)

PDF:

f(x\,;\,\alpha,\beta,\zeta)= \frac{|\beta|^{-1}}{\Gamma(\alpha)} \Bigl(\frac{x-\zeta}{\beta}\Bigr)^{\!\alpha-1} \exp\!\Bigl(-\frac{x-\zeta}{\beta}\Bigr), \quad\frac{x-\zeta}{\beta}>0

Positive \beta gives right-skewed; negative \beta gives left-skewed data.

Method of Moments: Let \bar x,s,g be the sample mean, standard deviation, and skewness coefficient. Then:

\hat\alpha=4/g^2,\qquad\hat\beta=sg/2,\qquad\hat\zeta=\bar x-2s/g

L-Moments: via lmomco::parpe3(), which returns (\mu,\sigma,\gamma). FitVerse converts using the MOM formulae above.

Use case: Annual peak-flow frequency analysis (US Bulletin 17C; England et al., 2019); meteorological extremes.

Log-Pearson Type III (LP3) – key: "lpearson3"

Parameters: \alpha>0 (shape), \beta (scale), \zeta (location) – all defined in log-space.

X\sim\mathrm{LP3}(\alpha,\beta,\zeta) \;\Leftrightarrow\; \log X\sim\mathrm{P3}(\alpha,\beta,\zeta)

PDF:

f(x\,;\,\alpha,\beta,\zeta) =\frac{f_{\mathrm{P3}}(\log x\,;\,\alpha,\beta,\zeta)}{x},\quad x>0

CDF and quantile:

F(x)=F_{\mathrm{P3}}(\log x),\qquad Q(p)=\exp\bigl[Q_{\mathrm{P3}}(p)\bigr]

Method of Moments and L-Moments: Pearson III estimators applied to \log x_i. L-Moments use lmomco::parpe3() on the log-transformed sample.

Use case: The US federal standard for flood-frequency estimation (Bulletin 17C; England et al., 2019).

Generalised Logistic Distribution (GLO) – key: "glo"

Parameters: \xi (location), \alpha>0 (scale), \kappa (shape).

CDF (\kappa\neq 0):

F(x\,;\,\xi,\alpha,\kappa) =\frac{1}{1+\bigl[1-\kappa(x-\xi)/\alpha\bigr]^{1/\kappa}}

Reduces to \mathrm{Logistic}(\xi,\alpha) when \kappa=0.

Quantile function:

Q(p)=\xi+\frac{\alpha}{\kappa}\Bigl[1-\Bigl(\frac{1-p}{p}\Bigr)^{\!\kappa}\Bigr], \quad\kappa\neq 0 \qquad Q(p)=\xi+\alpha\log\frac{p}{1-p},\quad\kappa=0

Support: (-\infty,\,\xi+\alpha/\kappa) for \kappa>0; (\xi+\alpha/\kappa,\,\infty) for \kappa<0; all real for \kappa=0.

L-Moments: via lmomco::parglo(). No sign flip: FitVerse \kappa=\kappa_{\text{lmomco}}.

MOM is not implemented for GLO (complex moment system; use LMOM or MLE).

Use case: The recommended flood-frequency distribution in the United Kingdom (WINFAP; Robson & Reed, 1999).

L-moment shape-parameter sign conventions

lmomco uses the opposite sign for the shape parameter of GEV and GPD:

Value

A data frame with one row per distribution and columns:

Key

Character: the string to pass in the dists argument of fitverse() and related functions.

Name

Full distribution name.

Family

Broad family tag (e.g. "Extreme-value", "Right-skewed").

Package

The R package supplying d/p/q/r functions, or NA for distributions with FitVerse-internal implementations.

MLE

Character "Yes" or "-": MLE available?

MOM

Character "Yes" or "-": MOM available?

LMOM

Character "Yes" or "-": L-Moments available?

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

References

England, J. F., Jr., et al. (2019). Guidelines for Determining Flood Flow Frequency – Bulletin 17C. USGS Techniques and Methods, Book 4, Chapter B5.

Hosking, J. R. M. and Wallis, J. R. (1997). Regional Frequency Analysis: An Approach Based on L-moments. Cambridge University Press.

Pickands, J. (1975). Statistical inference using extreme order statistics. The Annals of Statistics, 3(1), 119–131.

Robson, A. and Reed, D. (1999). Flood Estimation Handbook, Volume 3. Institute of Hydrology, Wallingford, UK.

See Also

fitverse for fitting any of the 52 distributions.
compare_dists for quick multi-distribution screening.
return_level for T-year return levels.
bootstrap_ci for parametric bootstrap confidence intervals.

Examples

# All 52 distributions
fitverse_distributions()

# Only those with L-moment support
fitverse_distributions(method = "LMOM")

# Extreme-value / hydrological families
ev <- fitverse_distributions()
ev[ev$Key %in% c("gev","gumbel","gpd","pearson3","lpearson3","glo"), ]

## ---- GPD (Peaks-Over-Threshold) ----------------------------------------
set.seed(42)
x_excess <- evd::rgpd(300, loc = 0, scale = 2, shape = 0.1)
fit_gpd   <- fitverse(x_excess, dists = "gpd", method = "LMOM",
                      plot = FALSE, verbose = FALSE)
coef(fit_gpd)
return_level(fit_gpd, return_periods = c(10, 50, 100))

## ---- Pearson III (Bulletin 17C style) ----------------------------------
## Pearson III with shape=5, scale=3, loc=10 is equivalent to
## a shifted Gamma: X = loc + Gamma(shape, scale)
set.seed(42)
x_p3  <- rgamma(200, shape = 5, scale = 3) + 10
fit_p3 <- fitverse(x_p3, dists = "pearson3", method = "LMOM",
                   plot = FALSE, verbose = FALSE)
coef(fit_p3)

## ---- Log-Pearson III ----------------------------------------------------
## Log-Pearson III: log(X) ~ Pearson III(shape, scale, loc)
## so X = exp(loc + Gamma(shape, scale))
set.seed(42)
x_lp3  <- exp(rgamma(200, shape = 5, scale = 0.4) + 2)
fit_lp3 <- fitverse(x_lp3, dists = "lpearson3", method = "LMOM",
                    plot = FALSE, verbose = FALSE)
coef(fit_lp3)

## ---- GLO (UK WINFAP style) ----------------------------------------------
## Generate GLO samples via the lmomco quantile function
set.seed(42)
x_glo  <- lmomco::qlmomco(runif(200),
             lmomco::vec2par(c(50, 8, 0.12), type = "glo"))
fit_glo <- fitverse(x_glo, dists = "glo", method = "LMOM",
                    plot = FALSE, verbose = FALSE)
coef(fit_glo)
return_level(fit_glo, return_periods = c(2, 10, 50, 100, 200))

Fit Distributions to Grouped or Interval-Censored Data

Description

fitverse_grouped() fits one or more parametric distributions to data that are available only in grouped (frequency-table) or interval-censored form. It accepts either class midpoints with counts, or explicit lower and upper interval boundaries with counts.

read_grouped_csv() reads a two- or three-column CSV file into the format expected by fitverse_grouped().

Usage

fitverse_grouped(midpoints = NULL,
                 counts,
                 lower   = NULL,
                 upper   = NULL,
                 dists   = "auto",
                 method  = "All",
                 xlab    = "x",
                 top_n   = 5L,
                 verbose = FALSE)

read_grouped_csv(file, sep = ",")

Arguments

Details

Input formats

fitverse_grouped() accepts two input formats.

Midpoint format: Supply midpoints and counts. Interval boundaries are reconstructed by assuming equal-width classes.

Boundary format: Supply lower, upper, and counts. Classes need not be equal-width. Open-ended classes (e.g., an upper bound of Inf) are handled by replacing the infinite bound with a large finite value for midpoint computation.

Estimation

For each candidate distribution, parameters are estimated by maximising the interval-censored log-likelihood

\ell(\theta) = \sum_{j} n_j \log\bigl[F(u_j;\theta) - F(l_j;\theta)\bigr]

where n_j is the count in class j, and l_j, u_j are its lower and upper boundaries.

CSV format for read_grouped_csv()

The CSV file should have either two columns (midpoint, count) or three columns (lower, upper, count). Column names are detected automatically from a set of common variants (mid, midpoint, lower, upper, etc.). Lines starting with # are treated as comments and ignored.

Value

fitverse_grouped() returns an object of class "fitverse_grouped" containing:

results

A data frame of fitted distributions ranked by AIC, with columns for distribution name, AIC, BIC, log-likelihood, estimated parameters, and goodness-of-fit statistics.

best

The name of the best-fitting distribution.

data

A list with elements midpoints, counts, lower, and upper.

call

The matched call.

read_grouped_csv() returns a named list with elements midpoints and counts (midpoint format) or lower, upper, and counts (boundary format), ready to be passed directly to fitverse_grouped().

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse for fitting to individual (ungrouped) data.
fitverse_distributions for the full distribution catalogue.

Examples

## Midpoint format
mid   <- c(5, 15, 25, 35, 45, 55, 65, 75)
count <- c(3, 12, 25, 30, 18, 8, 3, 1)
fit <- fitverse_grouped(midpoints = mid, counts = count,
                        method = "MLE", verbose = FALSE)
print(fit)

## Boundary format
low  <- c(0, 10, 20, 30, 40, 50, 60, 70)
high <- c(10, 20, 30, 40, 50, 60, 70, 80)
fit2 <- fitverse_grouped(lower = low, upper = high, counts = count,
                         method = "MLE", verbose = FALSE)
print(fit2)

## Reading from CSV
## A file with columns: lower, upper, count
tmp <- tempfile(fileext = ".csv")
writeLines(c(
  "# Example grouped data",
  "lower,upper,count",
  "0,10,5",
  "10,20,18",
  "20,30,30",
  "30,40,25",
  "40,50,12",
  "50,60,4"
), tmp)
dat <- read_grouped_csv(tmp)
fit3 <- fitverse_grouped(lower = dat$lower, upper = dat$upper,
                         counts = dat$counts)
print(fit3)

S3 Methods for fitverse Objects

Description

Standard S3 methods for objects of class "fitverse" as returned by fitverse(). These methods provide a consistent interface for extracting, displaying, and serialising fitted distribution results.

Usage

## S3 method for class 'fitverse'
print(x, ...)

## S3 method for class 'fitverse'
coef(object, ...)

## S3 method for class 'fitverse'
logLik(object, ...)

## S3 method for class 'fitverse'
AIC(object, ..., k = 2)

## S3 method for class 'fitverse'
BIC(object, ...)

## S3 method for class 'fitverse'
as.list(x, ...)

Arguments

Details

print.fitverse()

Prints a concise single-screen summary: best-fit distribution, estimation method, number of parameters, parameter estimates, log-likelihood, AIC, BIC, and up to three goodness-of-fit test results.

coef.fitverse()

Returns the named numeric vector of parameter estimates for the best-fit model. For example, c(shape = 2.14, rate = 0.43) for a Gamma fit.

logLik.fitverse()

Returns an object of class "logLik" with the maximised log-likelihood and the degrees of freedom (number of parameters). Compatible with stats::AIC() and stats::BIC().

AIC.fitverse()

Returns the Akaike Information Criterion \mathrm{AIC} = -2\ell + k \cdot p, where \ell is the log-likelihood, p is the number of parameters, and k is the penalty (default 2).

BIC.fitverse()

Returns the Bayesian Information Criterion \mathrm{BIC} = -2\ell + p \log n.

as.list.fitverse()

Converts the fitverse object to a plain list suitable for JSON serialisation (used by to_json()). Non-serialisable internals (raw fit objects, closures) are stripped. Returns: best_fit (dist_key, method, params, loglik, aic, bic), ranking, data_summary, n, xlab, criterion.

Value

print.fitverse()

Invisibly returns x.

coef.fitverse()

Named numeric vector of parameter estimates.

logLik.fitverse()

Object of class "logLik".

AIC.fitverse()

Numeric scalar.

BIC.fitverse()

Numeric scalar.

as.list.fitverse()

A named list safe for JSON serialisation.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse, to_json, plot.fitverse.

Examples

set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)
fit <- fitverse(x, plot = FALSE, verbose = FALSE)

print(fit)
coef(fit)
logLik(fit)
AIC(fit)
BIC(fit)

# Serialisation-ready list
lst <- as.list(fit)
names(lst)

Summary of a fitverse Object

Description

summary.fitverse returns an object of class "fitverse_summary" containing the best-fit results and a ranked table of all models fitted. Printing the object (via print.fitverse_summary) displays a formatted diagnostic report; assigning the result to a variable suppresses output.

Usage

## S3 method for class 'fitverse'
summary(object, top_n = 10L, ...)
## S3 method for class 'fitverse_summary'
print(x, ...)

Arguments

Value

An object of class "fitverse_summary", a list with components:

object

The original fitverse object.

top_n

The number of rows requested.

data_summary

Data summary list (or NULL).

ranking

Data frame with the top top_n rows of the ranking table.

print.fitverse_summary returns x invisibly.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse, print.fitverse, ranking, all_fits

Examples

  set.seed(1)
  x   <- rgamma(200, shape = 2, rate = 0.5)
  fit <- fitverse(x)

  # Assign without printing:
  s <- summary(fit)

  # Inspect individual components:
  head(s$ranking)
  s$data_summary$mean

  # Print explicitly:
  print(s)

Format a named numeric vector as "name (symbol) = val, ..."

Description

Format a named numeric vector as "name (symbol) = val, ..."

Usage

format_params(params, digits = 4, labels = NULL)

Arguments

Value

Single character string.

FitVerse colour palette for multi-distribution comparison plots

Description

FitVerse colour palette for multi-distribution comparison plots

Usage

fv_palette(n = 10)

Arguments

Value

Character vector of hex colour codes.

Generate an HTML or PDF Report for a Fitted Distribution

Description

Produces a standalone, self-contained HTML (or PDF) report that documents the fitted distribution, estimated parameters, goodness-of-fit tests, descriptive statistics, diagnostic plots, and (optionally) bootstrap confidence intervals.

If rmarkdown is installed, the report is rendered from an R Markdown template. Otherwise, a pure HTML fallback is used with the diagnostic plot embedded as a base64-encoded PNG image.

Usage

generate_report(fit,
                output_file = NULL,
                title       = "FitVerse: Distribution Fitting Report",
                author      = "",
                boot        = NULL,
                top_n       = 5L,
                open        = TRUE)

Arguments

Details

R Markdown path (recommended)

When rmarkdown is available, generate_report() writes a temporary .Rmd file and calls rmarkdown::render() to produce a self-contained HTML document using the html_document format with self_contained = TRUE. The fitted object is passed via a temporary .rds file to avoid workspace pollution.

HTML fallback path

When rmarkdown is unavailable, a self-contained HTML file is built in pure R:

• The diagnostic plot is saved to a temporary PNG and encoded as a base64 data URI using base64enc (if available) or a manual hex-to-base64 conversion.

• Tables are rendered as HTML <table> elements.

• The output is a single-file HTML document with embedded CSS.

Value

Invisibly returns the normalised path to the output file (character).

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

bootstrap_ci for confidence intervals to include in the report.
save_plot to save the diagnostic plot separately.

Examples

set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)
fit <- fitverse(x, xlab = "Waiting time (s)", plot = FALSE, verbose = FALSE)

# Write report to a temporary directory so no files are left in the working dir
path <- generate_report(fit,
                        title       = "Gamma Fit Example",
                        author      = "Karuna Reddy",
                        output_file = file.path(tempdir(), "fitverse_report.html"),
                        open        = FALSE)
message("Report saved to: ", path)

Extract the AD p-value from a fit result

Description

Extract the AD p-value from a fit result

Usage

get_ad_pvalue(fit_result)

Extract the Chi-Squared p-value from a fit result

Description

Extract the Chi-Squared p-value from a fit result

Usage

get_chisq_pvalue(fit_result)

Retrieve a Distribution Registry Entry

Description

Look up a single distribution entry from the FitVerse distribution registry by its key. This function is primarily intended for use by the FitVerseApp Shiny application, which needs to access display names, parameter labels, and PDF equation strings at runtime.

Usage

get_dist_entry(dist_key)

Arguments

Value

A named list containing the registry fields for the requested distribution (display name, parameter names, support type, method availability flags, etc.). Raises an error if dist_key is not found in the registry.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse_distributions for a tabular view of all registered distributions and their supported methods.

Extract the KS p-value from a fit result (for ranking purposes)

Description

Extract the KS p-value from a fit result (for ranking purposes)

Usage

get_ks_pvalue(fit_result)

Check whether an optional package is installed

Description

Emits an informative message() (not an error) if the package is missing, so that the caller can skip optional functionality gracefully.

Usage

has_pkg(pkg, feature = NULL)

Arguments

Value

Logical TRUE if available, FALSE otherwise.

Sample excess kurtosis

Description

Sample excess kurtosis

Usage

kurtosis_stat(x)

Arguments

Value

Numeric scalar.

L-Moment Ratio Diagram and Extreme-Value Threshold Selection

Description

Two diagnostic tools for exploratory data analysis in FitVerse.

lmrd() produces an L-moment ratio diagram plotting the sample (\tau_3, \tau_4) pair against the theoretical curves of major distribution families (GLO, GEV, GNO, PE3, GPA) and their known special cases (Normal, Logistic, Gumbel, Exponential, Uniform). The diagram guides distributional identification before formal model fitting.

ev_threshold() plots three threshold-stability diagnostics for Peaks-Over-Threshold (POT) extreme-value analysis: the Mean Residual Life (MRL) plot, the shape-stability plot, and the modified-scale-stability plot. These are used to identify the lowest threshold above which the Generalised Pareto Distribution (GPD) is a valid model for exceedances.

Usage

lmrd(x = NULL,
     xlab = NULL,
     distributions = c("glo", "gev", "gno", "pe3", "gpa"),
     show_ref_points = TRUE,
     draw = TRUE)

ev_threshold(x,
             thresholds      = NULL,
             n_thresholds    = 40L,
             min_exceedances = 10L,
             conf            = 0.95,
             xlab            = "x",
             draw            = TRUE)

Arguments

Details

lmrd() computes the sample L-skewness \hat\tau_3 and L-kurtosis \hat\tau_4 via lmoms and plots the sample point against the theoretical loci for the GLO, GEV, GPA, GNO, and PE3 families. Exact closed-form curves are used for GLO, GEV, and GPA (Hosking, 1990); polynomial approximations are used for GNO and PE3 (Hosking and Wallis, 1997).

ev_threshold() fits a GPD to exceedances at each candidate threshold u and extracts the estimated shape \hat\xi and modified scale \hat\sigma^* = \hat\sigma - \hat\xi u. A stable region in both plots, combined with approximate linearity of the MRL plot, indicates a suitable threshold.

Value

Both functions return a ggplot object invisibly and (by default) display the plot.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

References

Hosking, J. R. M. (1990). L-moments: Analysis and estimation of distributions using linear combinations of order statistics. Journal of the Royal Statistical Society B, 52(1), 105–124. doi:10.1111/j.2517-6161.1990.tb01775.x

Hosking, J. R. M. and Wallis, J. R. (1997). Regional Frequency Analysis: An Approach Based on L-Moments. Cambridge University Press.

Davison, A. C. and Smith, R. L. (1990). Models for exceedances over high thresholds. Journal of the Royal Statistical Society B, 52(3), 393–442.

See Also

fitverse with dists = "gpd" for formal GPD fitting after threshold selection.

Examples

set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)
lmrd(x, xlab = "Simulated Gamma data")



set.seed(1)
x_ev <- rexp(500, rate = 0.1)
ev_threshold(x_ev, n_thresholds = 30, xlab = "Simulated exponential")

Right-pad a string to a fixed width (for aligned console output)

Description

Right-pad a string to a fixed width (for aligned console output)

Usage

pad_right(s, width = 30)

Arguments

Value

String padded with trailing spaces.

Diagnostic Plot for a fitverse Object

Description

Produces a two-panel diagnostic figure for a fitted fitverse object:

1. A histogram of the data overlaid with the fitted density curves for the top-ranked distributions. The PDF formula is annotated directly on the figure.

2. A Q-Q plot of the empirical quantiles against theoretical quantiles from the best-fit distribution.

Either or both panels can be requested. The plot can optionally be returned as an interactive plotly HTML widget.

Usage

## S3 method for class 'fitverse'
plot(x,
     which       = "both",
     top_n       = 3L,
     xlab        = NULL,
     interactive = FALSE,
     draw        = TRUE,
     ...)

Arguments

Details

The density panel uses ggplot2::geom_histogram(aes(y = after_stat(density))) with Freedman-Diaconis binning, overlaid with geom_line() curves evaluated on a fine grid of 512 points spanning the observed data range.

The Q-Q panel compares theoretical quantiles from the best-fit distribution (evaluated via qfitverse()) against the empirical quantiles of the data, with a 45-degree reference line.

When which = "both", the two ggplot2 objects are arranged side-by-side using gridExtra::grid.arrange(). When interactive = TRUE and which != "both", a single panel is converted with plotly::ggplotly().

Value

Invisibly returns the plot object(s):

• A single ggplot2 object when which is "density" or "qq".

• A grob (arranged grid object) when which = "both".

• A plotly widget when interactive = TRUE.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

save_plot to write the figure to a PNG or PDF file.
fitverse which calls this function automatically when plot = TRUE.

Examples

set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)
fit <- fitverse(x, plot = FALSE, verbose = FALSE)

# Both panels (default)
plot(fit)

# Density panel only with custom x-axis label
plot(fit, which = "density", xlab = "Value", top_n = 5)

# Q-Q plot only
plot(fit, which = "qq")

Render a star-rating for a p-value (used in GoF table output)

Description

Render a star-rating for a p-value (used in GoF table output)

Usage

pval_stars(p)

Arguments

Value

Character string "***", "**", "*", "." or " ".

Extract the Ranking Table from a Fitted Object

Description

Retrieves the full ranking table from a fitverse object without using the $ operator. The table contains one row per distribution fitted, sorted by the selection criterion (AIC or BIC).

Usage

ranking(object, ...)
## S3 method for class 'fitverse'
ranking(object, ...)

Arguments

Value

A data frame with columns including Rank, Distribution, Method, AIC, Delta_AIC, BIC, KS_p, AD_p, and others depending on the methods used.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

fitverse, all_fits, summary.fitverse

Examples

  set.seed(1)
  x   <- rgamma(200, shape = 2, rate = 0.5)
  fit <- fitverse(x)
  r   <- ranking(fit)
  head(r)

Compute T-Year Return Levels from a Fitted Distribution

Description

Computes return levels (quantiles) corresponding to specified return periods from the best-fit distribution in a fitverse object. Return levels are widely used in hydrology, climatology, and actuarial science to characterise the magnitude of rare events.

Usage

return_level(fit,
             return_periods = c(2, 5, 10, 20, 50, 100, 200, 500))

Arguments

Details

A T-year return level x_T is the value that is exceeded on average once every T years (or time units), i.e.:

P(X > x_T) = \frac{1}{T}

x_T = F^{-1}\!\left(1 - \frac{1}{T}\right)

where F^{-1} is the quantile function of the fitted distribution.

Quantiles are computed using qfitverse(), which resolves the quantile function from the internal distribution registry.

Value

A data frame with one row per return period and columns:

T

Numeric: the return period.

exceedance_prob

Numeric: 1/T, the annual exceedance probability.

return_level

Numeric: the estimated T-year return level (upper 1/T quantile of the fitted distribution).

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

References

Coles, S. (2001). An Introduction to Statistical Modeling of Extreme Values. Springer, London.

Katz, R. W., Parlange, M. B. and Naveau, P. (2002). Statistics of extremes in hydrology. Advances in Water Resources, 25(8-12), 1287–1304.

See Also

qfitverse for quantiles at arbitrary probabilities.
bootstrap_ci for confidence intervals on return levels.

Examples

# Flood frequency analysis with GEV distribution
set.seed(10)
floods <- evd::rgev(80, loc = 200, scale = 40, shape = 0.15)

fit <- fitverse(floods,
                dists   = c("gev", "gumbel", "lognormal"),
                method  = "MLE",
                xlab    = "Annual maximum flood (m3/s)",
                plot    = FALSE,
                verbose = FALSE)

rl <- return_level(fit)
print(rl)

# 100-year return level
rl[rl$T == 100, "return_level"]

Attempt a fitting call and return NULL (with a message) on failure

Description

Wraps any fitting function in tryCatch so that a single failed distribution does not abort the entire fitverse() run.

Usage

safe_fit(expr, dist_name, method)

Arguments

Value

The result of evaluating expr, or NULL on error/warning-as-error.

Save a FitVerse Diagnostic Plot to File

Description

Renders the two-panel diagnostic plot for a fitverse object and saves it as a high-resolution PNG (or PDF) file, without displaying the figure on the current graphics device. Suitable for embedding in reports or automated pipelines.

Usage

save_plot(fit,
          file        = file.path(tempdir(), "fitverse_plot.png"),
          which       = "both",
          top_n       = 3L,
          xlab        = NULL,
          width       = 14,
          height      = 5.5,
          dpi         = 600,
          ...)

Arguments

Details

Uses grDevices::png() (or grDevices::pdf() / grDevices::svg() depending on file extension) and gridExtra::grid.arrange() to render both panels in a single device call, ensuring consistent alignment.

The function opens the graphics device, renders the figure, and closes the device via on.exit(), so the current interactive device is unaffected.

Value

Invisibly returns the normalised path to the saved file (character).

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

plot.fitverse for interactive display.

Examples

set.seed(42)
x <- rgamma(300, shape = 2, rate = 0.5)
fit <- fitverse(x, plot = FALSE, verbose = FALSE)

save_plot(fit, file = tempfile(fileext = ".png"))

Identify the best-fitting model from a ranked table

Description

Returns the Rank-1 model from the ranked table produced by build_ranking_table(). Rank 1 is determined by the four-key sort: convergence, then GoF tier, then AIC/BIC, then GoF composite score.

Usage

select_best_fit(ranking, fits)

Arguments

Value

A single fit-result list for the best model.

Sample skewness (Fisher's moment coefficient)

Description

Sample skewness (Fisher's moment coefficient)

Usage

skewness_stat(x)

Arguments

Value

Numeric scalar.

Serialise a fitverse Object to JSON

Description

Converts a fitverse object to a JSON string. Useful for passing fit results to web front-ends, storing results in document databases, or integrating with REST APIs.

The serialisation uses toJSON() on the output of as.list.fitverse(), which strips non-serialisable internals (raw optimiser fit objects, closures, environments) and returns a flat, portable representation.

Usage

to_json(x, pretty = TRUE, ...)

## S3 method for class 'fitverse'
to_json(x, pretty = TRUE, ...)

Arguments

Details

Requires the jsonlite package (listed in Suggests). An informative error is raised if it is not installed.

The JSON output includes:

• best_fit: distribution key, method, parameter estimates, loglik, AIC, BIC.

• ranking: full ranking data frame.

• data_summary: descriptive statistics.

• n, xlab, criterion.

The raw data vector is deliberately not included to keep the output compact and privacy-friendly for web contexts.

Value

A single character string containing the JSON representation of x.

Author(s)

Karuna G. Reddy karuna.reddy@auckland.ac.nz
M. G. M. Khan mgm.khan@usp.ac.fj

See Also

as.list.fitverse for the underlying list conversion.
fromJSON to parse the JSON back to R.

Examples

if (requireNamespace("jsonlite", quietly = TRUE)) {
  set.seed(42)
  x <- rgamma(300, shape = 2, rate = 0.5)
  fit <- fitverse(x, plot = FALSE, verbose = FALSE)

  json_str <- to_json(fit)
  cat(substr(json_str, 1, 400))

  # Round-trip back to R list
  lst <- jsonlite::fromJSON(json_str)
  lst$best_fit$params
}