Introduction to UniLindleyApprox: Bayesian Estimation via Lindley’s Approximation

Shikhar Tyagi, Arvind Pandey, Bhupendra Singh, Vrijesh Tripathi

2026-07-29

Overview

The UniLindleyApprox package provides a generalized computational framework for performing Bayesian parameter point estimation using Lindley’s Approximation (1980) for arbitrary univariate probability distributions under complete, censored, and truncated data.

Lindley’s approximation evaluates posterior expectations of arbitrary smooth functions \(g(\boldsymbol{\theta})\): \[\mathrm{E}[g(\boldsymbol{\theta}) \mid \mathbf{x}] = \frac{\int g(\boldsymbol{\theta}) e^{\ell(\boldsymbol{\theta}) + \rho(\boldsymbol{\theta})} d\boldsymbol{\theta}}{\int e^{\ell(\boldsymbol{\theta}) + \rho(\boldsymbol{\theta})} d\boldsymbol{\theta}}\] using Taylor series expansions around the posterior mode (MAP) or MLE.

Features

Example: Exponential Distribution under Complete Data

library(UniLindleyApprox)

# Define probability functions for Exponential(rate)
dexp_custom <- function(x, theta) dexp(x, rate = theta[1])
pexp_custom <- function(x, theta) pexp(x, rate = theta[1])
sexp_custom <- function(x, theta) 1 - pexp(x, rate = theta[1])

# Gamma prior for rate parameter
logprior <- function(theta) dgamma(theta[1], shape = 2, rate = 1, log = TRUE)

# Simulated data
set.seed(42)
x <- rexp(50, rate = 2)

# Fit model using Lindley's approximation
fit <- lindley_fit(
  data = x,
  pdf = dexp_custom,
  cdf = pexp_custom,
  survival = sexp_custom,
  log_prior = logprior,
  theta0 = c(1.5),
  scheme = "complete",
  loss = "SELF"
)
#> Constructing log-likelihood for scheme: complete
#> Computing posterior...
#> Finding posterior mode...
#> Computing derivatives...
#> Computing Bayes estimate under loss: SELF
#> Estimation complete. Elapsed time: 0.02 seconds

# Print results
print(fit)
#> Bayesian Estimation Using Lindley's Approximation
#> ================================================
#> 
#> Censoring Scheme: complete 
#> Loss Function: SELF 
#> 
#> Posterior Mode (MAP):
#> [1] 1.735768
#> 
#> Bayes Estimate (SELF):
#> [1] 1.769802
#> 
#> Posterior Mean:
#> [1] 1.769802
#> 
#> Log-Likelihood: -21.6917 
#> Log-Posterior: -22.8761 
#> Convergence: Successful 
#> Iterations: 24 
#> Elapsed Time: 0.02 seconds
summary(fit)
#> Summary of Bayesian Estimation Using Lindley's Approximation
#> ============================================================
#> 
#> Model Information:
#>   Censoring Scheme: complete 
#>   Loss Function: SELF 
#>   Number of Parameters: 1 
#> 
#> Parameter Estimates:
#>  Parameter      MAP Bayes_Estimate Posterior_Mean
#>     theta1 1.735768       1.769802       1.769802
#> 
#> Model Fit:
#>   Log-Likelihood: -21.6917 
#>   Log-Posterior: -22.8761 
#>   AIC: 45.3835 
#>   BIC: 47.2955 
#> 
#> Optimization:
#>   Convergence: Successful 
#>   Iterations: 24 
#>   Elapsed Time: 0.02 seconds
#> 
#> Information Matrix:
#>          [,1]
#> [1,] 16.92727

Bayes Estimation Under Multiple Loss Functions

# SELF Bayes estimate
b_self <- bayes_estimate(fit, loss = "SELF")

# LINEX Bayes estimate
b_linex <- bayes_estimate(fit, loss = "LINEX", a = 0.5)

# GELF Bayes estimate
b_gelf <- bayes_estimate(fit, loss = "GELF", q = 2)

cat("SELF Bayes Estimate: ", b_self, "\n")
#> SELF Bayes Estimate:  1.769802
cat("LINEX Bayes Estimate: ", b_linex, "\n")
#> LINEX Bayes Estimate:  1.755127
cat("GELF Bayes Estimate: ", b_gelf, "\n")
#> GELF Bayes Estimate:  1.718997

Diagnostics & Model Information Criteria

cat("AIC: ", AIC(fit), "\n")
#> AIC:  45.38349
cat("BIC: ", BIC(fit), "\n")
#> BIC:  47.29552
cat("AICc: ", AICc(fit), "\n")
#> AICc:  45.46683
cat("HQIC: ", HQIC(fit), "\n")
#> HQIC:  46.1116