--- title: "Introduction to BivLaplaceRL" author: "Mahesh Divakaran" date: "`r Sys.Date()`" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Introduction to BivLaplaceRL} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", fig.width = 6, fig.height = 4 ) ``` ## Overview **BivLaplaceRL** is an R package that implements statistical methods from two research papers on bivariate Laplace transforms in reliability theory, together with a set of univariate residual life tools: 1. Bivariate Laplace transform of residual lives (Jayalekshmi, Rajesh & Nair, 2022) 2. Bivariate Laplace transform order of reversed residual lives (Jayalekshmi & Rajesh) ```{r setup} library(BivLaplaceRL) ``` ## Parametric Distributions The package provides four bivariate distributions commonly used in reliability modelling. ### Gumbel Bivariate Exponential ```{r gumbel} set.seed(42) dat <- rgumbel_biv(200, k1 = 1, k2 = 1.5, theta = 0.3) head(dat) sgumbel_biv(1, 1, k1 = 1, k2 = 1.5, theta = 0.3) ``` ### Farlie-Gumbel-Morgenstern (FGM) ```{r fgm} set.seed(1) dat_fgm <- rfgm_biv(100, theta = 0.5) pfgm_biv(0.5, 0.6, theta = 0.5) ``` ### Schur-Constant ```{r schur} set.seed(2) dat_sc <- rschur_biv(100, lambda = 1) summary(rowSums(dat_sc)) ``` ## Bivariate Laplace Transform of Residual Lives ### Closed-Form (Gumbel Distribution) $$L^*_{X_{t_1|t_2}}(s_1) = \frac{1}{k_1 + s_1 + \theta t_2}, \quad L^*_{X_{t_2|t_1}}(s_2) = \frac{1}{k_2 + s_2 + \theta t_1}.$$ ```{r closed_form} blt_residual_gumbel(s1 = 1, s2 = 1, t1 = 0.5, t2 = 0.5, k1 = 1, k2 = 1, theta = 0.5) ``` ### Numerical Computation ```{r numerical} my_surv <- function(x1, x2) exp(-x1 - x2 - 0.3 * x1 * x2) blt_residual(s1 = 1, s2 = 1, t1 = 0.5, t2 = 0.5, surv_fn = my_surv, star = TRUE) ``` ### Nonparametric Estimation ```{r np_estim} set.seed(123) dat <- rgumbel_biv(500, k1 = 1, k2 = 1, theta = 0.5) np_blt_residual(dat, s1 = 1, s2 = 1, t1 = 0.3, t2 = 0.3) # True value blt_residual_gumbel(s1 = 1, s2 = 1, t1 = 0.3, t2 = 0.3, theta = 0.5) ``` ### Simulation Study ```{r sim_study} sim_blt_residual(n_obs = 200, n_sim = 50, s1 = 1, s2 = 1, t1 = 0.3, t2 = 0.3, k1 = 1, k2 = 1, theta = 0.5) ``` ## Bivariate Hazard Gradient and MRL ```{r hazard_mrl} biv_hazard_gradient(t1 = 1, t2 = 1, k1 = 1, k2 = 1, theta = 0.3) biv_mean_residual(t1 = 0.5, t2 = 0.5, k1 = 1, k2 = 1, theta = 0) ``` ## NBUHR / NWUHR Aging Classification ```{r aging} res <- nbuhr_test(t2 = 0.5, k1 = 1, k2 = 1, theta = 0.3, t1_grid = seq(0.1, 3, 0.3)) cat("Component 1:", res$class1, "\n") cat("Component 2:", res$class2, "\n") ``` ## Bivariate Laplace Transform of Reversed Residual Lives ```{r reversed} blt_reversed(s1 = 1, s2 = 1, t1 = 0.5, t2 = 0.5, theta = 0.3) biv_rhazard_gradient(x1 = 0.5, x2 = 0.5, theta = 0.3) biv_rmrl(x1 = 0.5, x2 = 0.5, theta = 0.3) ``` ## Bivariate Stochastic Orders ```{r blt_order} sX <- function(x1, x2) sgumbel_biv(x1, x2, k1 = 1, k2 = 1, theta = 0.2) sY <- function(x1, x2) sgumbel_biv(x1, x2, k1 = 2, k2 = 2, theta = 0.2) res <- blt_order_residual(sX, sY, s1 = 1, s2 = 1, t1_grid = c(0.5, 1), t2_grid = c(0.5, 1)) cat("X <=_BLt-rl Y:", res$order_holds, "\n") res2 <- biv_whr_order(sX, sY, t1_grid = c(0.5, 1), t2_grid = c(0.5, 1)) cat("X <=_whr Y:", res2$order_holds, "\n") ``` ## Univariate Residual Life Analysis The package now provides a complete set of univariate residual life tools that complement the bivariate framework. ### Laplace Transform of Residual Life For a non-negative continuous random variable $X$ with density $f$ and survival function $\bar{F}$: $$L_X(s,t) = E[e^{-sX} \mid X > t] = \frac{1}{\bar{F}(t)}\int_t^\infty e^{-sx} f(x)\,dx.$$ ```{r uni_lt} f <- function(x) dexp(x, 1) Fb <- function(x) pexp(x, 1, lower.tail = FALSE) # At t=0: L_X(1, 0) = 1/(1+1) = 0.5 lt_residual(f, Fb, s = 1, t = 0) # At t=0.5 lt_residual(f, Fb, s = 1, t = 0.5) ``` ### Hazard Rate and MRL ```{r uni_hr_mrl} # Exp(1): constant hazard rate = 1 hazard_rate(f, Fb, t = c(0.5, 1, 2)) # Exp(1): constant MRL = 1 (memoryless) mean_residual(Fb, t = c(0, 0.5, 1, 2)) # Gamma(2,1): increasing hazard, decreasing MRL fG <- function(x) dgamma(x, shape = 2, rate = 1) FbG <- function(x) pgamma(x, shape = 2, rate = 1, lower.tail = FALSE) hazard_rate(fG, FbG, t = c(0.5, 1, 2)) mean_residual(FbG, t = c(0, 0.5, 1)) ``` ### Nonparametric Estimation ```{r uni_np} set.seed(7) x <- rexp(500, rate = 1) # Nonparametric estimate at s=1, t=0: true value = 0.5 np_lt_residual(x, s = 1, t = 0) # At t=0.5 np_lt_residual(x, s = 1, t = 0.5) ``` ### Univariate Stochastic Orders ```{r uni_orders} f1 <- function(x) dexp(x, 1) Fb1 <- function(x) pexp(x, 1, lower.tail = FALSE) f2 <- function(x) dexp(x, 2) Fb2 <- function(x) pexp(x, 2, lower.tail = FALSE) # LT-rl order: Exp(1) <=_Lt-rl Exp(2)? lt_rl_order(f1, Fb1, f2, Fb2, s_grid = c(0.5, 1, 2), t_grid = c(0, 0.5, 1))$order_holds # Hazard rate order: Exp(1) <=_hr Exp(2)? hr_order(f1, Fb1, f2, Fb2, t_grid = c(0.5, 1, 2))$order_holds # MRL order: Exp(2) <=_mrl Exp(1)? mrl_order(Fb2, Fb1, t_grid = c(0, 0.5, 1, 2))$order_holds ``` ## Entropy and Information Generating Functions ```{r entropy} # Shannon entropy of Exp(1) = 1 shannon_entropy(function(x) dexp(x, 1)) # Golomb IGF of Exp(1) at alpha=2: 0.5 info_gen_function(function(x) dexp(x, 1), alpha = 2) ``` ## References Jayalekshmi S., Rajesh G., Nair N.U. (2022). Bivariate Laplace transform of residual lives and their properties. *Communications in Statistics---Theory and Methods*. Jayalekshmi S., Rajesh G., Nair N.U. (2022). Bivariate Laplace transform order and ordering of reversed residual lives. *International Journal of Reliability, Quality and Safety Engineering*. Belzunce F., Ortega E., Ruiz J.M. (1999). The Laplace order and ordering of residual lives. *Statistics & Probability Letters*, 42(2), 145--156.