tests/testthat/test_fit.R

In IGPFrailty: Inverse Gaussian Process Degradation Models with Frailty

test_that("igp_fit works on laser dataset across all frailty specifications", {
  data(laser)
  
  # Classical IGP
  fit_none <- igp_fit(laser, time_col = "t", deg_col = "increase", unit_col = "unit", frailty = "none")
  expect_s3_class(fit_none, "igp_fit")
  expect_equal(fit_none$frailty, "none")
  expect_equal(unname(fit_none$coefficients["theta"]), 2.038, tolerance = 1e-2)
  expect_equal(unname(fit_none$coefficients["eta"]), 13.15, tolerance = 1e-1)
  expect_equal(unname(round(fit_none$aic, 1)), -146.2)
  
  # IGP-Gamma
  fit_gamma <- igp_fit(laser, time_col = "t", deg_col = "increase", unit_col = "unit", frailty = "gamma")
  expect_s3_class(fit_gamma, "igp_fit")
  expect_equal(fit_gamma$frailty, "gamma")
  expect_equal(unname(fit_gamma$coefficients["theta"]), 2.052, tolerance = 1e-2)
  expect_equal(unname(fit_gamma$coefficients["eta"]), 15.17, tolerance = 1e-1)
  expect_equal(unname(fit_gamma$coefficients["xi"]), 0.210, tolerance = 1e-2)
  expect_equal(unname(round(fit_gamma$aic, 1)), -174.7)
  
  # IGP-IG
  fit_ig <- igp_fit(laser, time_col = "t", deg_col = "increase", unit_col = "unit", frailty = "ig")
  expect_s3_class(fit_ig, "igp_fit")
  expect_equal(fit_ig$frailty, "ig")
  expect_equal(unname(fit_ig$coefficients["theta"]), 2.057, tolerance = 1e-2)
  expect_equal(unname(fit_ig$coefficients["eta"]), 15.13, tolerance = 1e-1)
  expect_equal(unname(fit_ig$coefficients["xi"]), 0.248, tolerance = 1e-2)
  expect_equal(unname(round(fit_ig$aic, 1)), -176.0)
})

test_that("igp_fit works on crack dataset and matches paper results", {
  data(crack)
  
  # IGP-Gamma on crack
  fit_crack_gam <- igp_fit(crack, time_col = "t", deg_col = "deg", unit_col = "specimen", frailty = "gamma")
  expect_s3_class(fit_crack_gam, "igp_fit")
  expect_equal(unname(fit_crack_gam$coefficients["theta"]), 0.00495, tolerance = 1e-3)
  expect_equal(unname(fit_crack_gam$coefficients["eta"]), 145.6, tolerance = 0.5)
  expect_equal(unname(fit_crack_gam$coefficients["xi"]), 0.416, tolerance = 1e-2)
  expect_equal(unname(round(fit_crack_gam$aic, 1)), -1316.7)
})

test_that("igp_fit handles invalid inputs gracefully", {
  data(laser)
  expect_error(igp_fit("not_a_df"), "data frame")
  expect_error(igp_fit(laser, time_col = "missing"), "not found")
})