Overview

The AlphaPowerHazard package implements the flexible Alpha-Power Hazard probability distribution and regression models proposed by Pal et al. (2026) (Journal of the Indian Society for Probability and Statistics, DOI: 10.1007/s41096-026-00297-5).

The baseline hazard rate function is: \[h(x; \alpha, \beta) = \alpha^x + x^{\beta-1}, \quad \alpha > 1, \beta > 0, x > 0\]

The baseline cumulative hazard function is: \[H(x; \alpha, \beta) = \frac{\alpha^x - 1}{\log \alpha} + \frac{x^\beta}{\beta}\]

The baseline survival function is: \[S(x; \alpha, \beta) = \exp\left(-H(x; \alpha, \beta)\right)\]

The probability density function is: \[f(x; \alpha, \beta) = (\alpha^x + x^{\beta-1}) \exp\left(-\left[\frac{\alpha^x - 1}{\log \alpha} + \frac{x^\beta}{\beta}\right]\right)\]

Statistical Properties

Distribution Functions

x_vals <- c(0.5, 1.0, 1.5, 2.0)
alpha <- 1.5
beta <- 0.8

# PDF, CDF, Survival, Hazard, Cumulative Hazard
dalpha_power(x_vals, alpha, beta)
#> [1] 0.66505616 0.20869926 0.06212618 0.01622428
palpha_power(x_vals, alpha, beta)
#> [1] 0.7197927 0.9165203 0.9774842 0.9948008
salpha_power(x_vals, alpha, beta)
#> [1] 0.280207319 0.083479705 0.022515808 0.005199173
halpha_power(x_vals, alpha, beta)
#> [1] 2.373443 2.500000 2.759225 3.120551
Halpha_power(x_vals, alpha, beta)
#> [1] 1.272226 2.483152 3.793538 5.259256

# Quantiles
qalpha_power(c(0.25, 0.50, 0.75), alpha, beta)
#> [1] 0.0952583 0.2572982 0.5480327

Distributional Moments & Quantile Statistics

mom <- moments_alpha_power(k = 4, alpha = 1.5, beta = 1.1)
cat("Mean:", mom$mean, "\n")
#> Mean: 0.4804367
cat("Variance:", mom$variance, "\n")
#> Variance: 0.1775296
cat("Skewness:", mom$skewness, "\n")
#> Skewness: 1.44417
cat("Kurtosis:", mom$kurtosis, "\n")
#> Kurtosis: 5.528196

qs <- quantile_stats_alpha_power(alpha = 1.5, beta = 1.1)
cat("Bowley Skewness:", qs$bowley_skewness, "\n")
#> Bowley Skewness: 0.2191885
cat("Moors Kurtosis:", qs$moors_kurtosis, "\n")
#> Moors Kurtosis: 1.243827

Hazard Minimum & Order Statistics

# Unique minimum location for beta < 1 (bathtub shape)
hrf_m <- hrf_min_alpha_power(alpha = 1.5, beta = 0.8)
cat("HRF Minimum x:", hrf_m$xmin, "h(x):", hrf_m$hmin, "\n")
#> HRF Minimum x: 0.4729624 h(x): 2.372933

# Order statistics
order_stats_alpha_power(c(0.5, 1.0), r = 2, n = 5, alpha = 1.5, beta = 1.1, type = "pdf")
#> [1] 0.537817662 0.008380725

Classical Baseline Estimation Methods

The package implements five classical parameter estimation methods for the baseline distribution: 1. Maximum Likelihood Estimation (MLE) 2. Least Squares Estimation (LSE) 3. Weighted Least Squares Estimation (WLSE) 4. Maximum Product of Spacings Estimation (MPSE) 5. Cramér-von Mises Estimation (CME)

set.seed(123)
sim_x <- ralpha_power(50, alpha = 1.5, beta = 0.8)

fit_mle  <- fit_alpha_power_base(sim_x, method = "mle")
fit_lse  <- fit_alpha_power_base(sim_x, method = "lse")
fit_wlse <- fit_alpha_power_base(sim_x, method = "wlse")
fit_mpse <- fit_alpha_power_base(sim_x, method = "mpse")
fit_cme  <- fit_alpha_power_base(sim_x, method = "cme")

fit_mle
#> Alpha-Power Hazard Base Parameter Estimation (MLE)
#> ====================================================
#>       Estimate  StdErr
#> alpha  1.38760 0.58822
#> beta   0.86409 0.11333
#> 
#> Objective Value: -6.59988 
#> Convergence Code: 0

Hazard Regression Models (M1-M4)

The package supports four hazard regression models (M1, M2, M3, M4) across uncensored data, right censoring, left censoring, interval censoring, and progressive censoring schemes:

Model M1: \(h(y|Z) = (\alpha^y + y^{\beta-1}) \exp(Z^\top \theta)\)
Model M2: \(h(y|Z) = \alpha^y + y^{\exp(Z^\top \theta)-1}\)
Model M3: \(h(y|Z) = (1 + \exp(Z^\top \theta))^y + y^{\beta-1}\)
Model M4: \(h(y|Z) = (1 + \exp(Z^\top \theta))^y + y^{\exp(Z^\top \theta)-1}\)

set.seed(42)
n <- 60
df <- data.frame(
  time   = ralpha_power(n, alpha = 1.5, beta = 0.8),
  status = sample(c(1, 1, 1, 0), n, replace = TRUE),
  age    = rnorm(n),
  bmi    = rnorm(n)
)

fit_m1 <- alpha_power_hazard(survival::Surv(time, status) ~ age + bmi, data = df, model = "M1")
summary(fit_m1)
#> Call:
#> alpha_power_hazard(formula = survival::Surv(time, status) ~ age + 
#>     bmi, data = df, model = "M1")
#> 
#> =========================================================
#>  Alpha-Power Hazard Regression Model (Model M1)
#> =========================================================
#>  Sample size           : 60 
#>  No. of covariates     : 2 
#>  No. of baseline pars  : 2 
#> 
#> ----- Parameter Estimates --------------------------------
#>       Estimate StdErr t_stat  p_value CI_lower CI_upper Signif
#> alpha   1.2723 0.6197  2.053 0.044742  1.00010  2.48684      *
#> beta    1.2455 0.3979  3.130 0.002773  0.46567  2.02530     **
#> age    -0.4556 0.2271 -2.006 0.049669 -0.90079 -0.01051      *
#> bmi     0.3870 0.1593  2.429 0.018398  0.07466  0.69925      *
#> 
#> ----- Model Fit ------------------------------------------
#>   Log-Likelihood : -20.9800
#>   AIC            : 49.9600
#>   BIC            : 58.3374
#> 
#> ----- Overall F-test (covariates) ------------------------
#>   F-stat = 5.1235,  p-value = 0.009049
#> 
#> ----- Multicollinearity Diagnostics ----------------------
#>       VIF Tolerance
#> age 1.002     0.998
#> bmi 1.002     0.998

Model Diagnostics & Predictions

res <- residuals_alpha_power(fit_m1)
head(res$cox_snell)
#> [1] 1.2701583 1.1401953 0.1364794 0.4670338 0.9322080 0.2069983
head(res$martingale)
#> [1] -0.27015830 -0.14019529  0.86352056  0.53296623  0.06779196  0.79300174

# Predictions
pred_s <- predict_alpha_power(fit_m1, type = "survival")
head(pred_s$survival[, 1:4])
#>            [,1]       [,2]      [,3]      [,4]
#> [1,] 0.28078717 0.23375690 0.8936540 0.4223338
#> [2,] 0.36920939 0.31975657 0.9155756 0.5085588
#> [3,] 0.21400373 0.17131072 0.8724243 0.3512434
#> [4,] 0.50247593 0.45496584 0.9408973 0.6268589
#> [5,] 0.06955134 0.04733927 0.7898032 0.1638179
#> [6,] 0.40806221 0.35854573 0.9237210 0.5442893

Introduction to AlphaPowerHazard