test-marp.R
In marp: Model-Averaged Renewal Process

test_that("marp", {
  # (optional) load the small test dataset
  # data_file <- system.file("extdata", "small.txt", package = "marp", mustWork = TRUE)
  # data <- read.table(data_file)$V1

  # set some parameters
  m <- 10 # number of iterations for MLE optimization
  t <- seq(100,200,by=10) # time intervals
  y <- 304 # cut-off year for estimating probablity
  which.model <- 2 # specify the generating model

  # fix the random seed
  set.seed(42)

  # sample data for testing
  data <- rgamma(30, 3, 0.01)

  # fit renewal model
  suppressWarnings(  # suppressing warnings from stats::nlm: NA/Inf replaced by maximum positive value
    res <- marp::marp(data, t, m, y, which.model)
  )

  # check result
  expect_equal(res$par1, c(3.4136086430979953e-03, 2.7626793657057762e+00, 2.6037039674870583e+00, 3.3080162440951688e+02, 5.4882183788378658e+00, 2.9294512422957860e+02))
  expect_equal(res$par2, c(NA, 9.4307059277139432e-03, 2.4759796859031687e+02, 1.8010183507666513e+00, 6.5084541680686814e-01, 7.1824719073918109e-01))
  expect_equal(res$logL, c(-200.39955882401648, -193.76650950646049, -195.12976531758551, -193.92377676645586, -194.33021376270909, -194.40996001614093))
  expect_equal(res$AIC, c(402.79911764803296, 391.53301901292099, 394.25953063517102, 391.84755353291172, 392.66042752541819, 392.81992003228186))
  expect_equal(res$BIC, c(404.20031502969510, 394.33541377624528, 397.06192539849530, 394.64994829623600, 395.46282228874247, 395.62231479560614))
  expect_equal(res$mu_hat, c(292.94512187913182, 292.94512912200048, 319.72017228620746, 294.16945213908519, 298.87285747700128, 292.94512422957860))
  expect_equal(res$pr_hat, c(0.60038574701819891, 0.42154974433034809, 0.53433568234281148, 0.30779792692414687, 0.56416103510057725, 0.61794524610544410), tolerance = 1e-6)

  haz_hat_expected <-
    c(
      -5.6799852941338829, -5.6799852941338829, -5.6799852941338829, -5.6799852941338829, -5.6799852941338829, -5.6799852941338829,
      -5.6799852941338829, -5.6799852941338829, -5.6799852941338829, -5.6799852941338829, -5.6799852941338829, -6.0942031084732298,
      -5.9967873794574516, -5.9117418563554684, -5.8368230853439300, -5.7703089176306639, -5.7108525626839901, -5.6573839062669986,
      -5.6090408956082456, -5.5651206740587922, -5.5250440506799734, -5.4883291920475745, -6.0990192429336094, -5.9701664705134210,
      -5.8576899644670348, -5.7593884711134971, -5.6734972529860741, -5.5985621349393231, -5.5333565788683616, -5.4768259914915305,
      -5.4280496904694857, -5.3862145095364315, -5.3505961502861927, -6.1714638710963881, -6.0951186680582552, -6.0254209583640863,
      -5.9613052806725335, -5.9019434350392981, -5.8466788789061646, -5.7949823391436279, -5.7464209045603756, -5.7006359661738628,
      -5.6573271297614109, -5.6162402596857071, -5.9235521978533958, -5.8023896004395645, -5.7047473880293342, -5.6252373537796752,
      -5.5599409055534252, -5.5059486025117375, -5.4610610586440487, -5.4235891601883868, -5.3922173604047572, -5.3659081375131672,
      -5.3438339586221275, -5.7911126719889303, -5.6765973314326752, -5.5892417143301261, -5.5216608261560411, -5.4687921205249133,
      -5.4270729562323066, -5.3939387902533049, -5.3675067327627373, -5.3463701567645607, -5.3294619641245422, -5.3159614865560094
    )
  expect_true(all.equal(res$haz_hat, haz_hat_expected))
})