## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

## -----------------------------------------------------------------------------
library(ReliaGrowR)

id   <- c(1, 1, 1, 2, 2, 2, 3, 3, 3, 3)
time <- c(100, 350, 500, 80, 300, 600, 150, 250, 400, 700)

## -----------------------------------------------------------------------------
result <- mcf(id, time)
plot(result, main = "Mean Cumulative Function",
     xlab = "Time", ylab = "MCF")

## -----------------------------------------------------------------------------
df <- data.frame(id = id, time = time)
result2 <- mcf(data = df)

## -----------------------------------------------------------------------------
id    <- c(1, 1, 1, 2, 2, 2, 3, 3, 3)
time  <- c(100, 350, 500, 80, 300, 400, 150, 250, 700)
event <- c(  1,   1,   0,  1,   1,   0,   1,   1,   1)

result <- mcf(id, time, event)

## -----------------------------------------------------------------------------
id   <- c(1, 1, 2, 2, 3, 3)
time <- c(100, 300, 150, 400, 200, 350)

# Without end_time: system observation ends at last event
mcf_basic <- mcf(id, time)

# With end_time: all systems observed until time 800
mcf_adj <- mcf(id, time, end_time = c("1" = 800, "2" = 800, "3" = 800))

## -----------------------------------------------------------------------------
par(mfrow = c(1, 2))
plot(mcf_basic, main = "MCF (inferred exposure)",
     xlab = "Time", ylab = "MCF")
plot(mcf_adj, main = "MCF (explicit exposure)",
     xlab = "Time", ylab = "MCF")

## -----------------------------------------------------------------------------
id    <- c(1, 1, 1, 2, 2, 2, 3, 3, 3)
time  <- c(100, 350, 500, 80, 300, 400, 150, 250, 700)
event <- c(  1,   1,   0,  1,   1,   0,   1,   1,   1)

exp_result <- exposure(id, time, event)
mcf_result <- mcf(id, time, event, end_time = exp_result$end_times)

## -----------------------------------------------------------------------------
id   <- c(1, 1, 1, 2, 2, 2, 3, 3, 3, 3)
time <- c(100, 350, 500, 80, 300, 600, 150, 250, 400, 700)
result <- exposure(id, time)

## ----fig.height=8-------------------------------------------------------------
plot(result)

## -----------------------------------------------------------------------------
plot(result, which = "exposure")

## -----------------------------------------------------------------------------
plot(result, which = "event_rate")

## -----------------------------------------------------------------------------
id    <- c(1, 1, 1, 2, 2, 2, 3, 3, 3)
time  <- c(100, 350, 500, 80, 300, 400, 150, 250, 700)
event <- c(  1,   1,   0,  1,   1,   0,   1,   1,   1)

result <- exposure(id, time, event)

## -----------------------------------------------------------------------------
id   <- c(1, 1, 1, 2, 2, 2, 3, 3, 3, 3)
time <- c(100, 350, 500, 80, 300, 600, 150, 250, 400, 700)
m <- mcf(id, time)

## -----------------------------------------------------------------------------
fit_mle <- nhpp(m)
plot(fit_mle, main = "Power Law NHPP (MLE)",
     xlab = "Time")

## -----------------------------------------------------------------------------
fit_ls <- nhpp(m, method = "LS")

## -----------------------------------------------------------------------------
fit_ll <- nhpp(m, model_type = "Log-Linear")
plot(fit_ll, main = "Log-Linear NHPP",
     xlab = "Time")

## -----------------------------------------------------------------------------
id2   <- c(1,1,1,1,1, 2,2,2,2,2, 3,3,3,3,3)
time2 <- c(100,200,300,400,500, 80,200,350,450,600, 150,250,400,550,700)

m2 <- mcf(id2, time2)
fit_pw <- nhpp(m2, breaks = c(350), method = "LS")
plot(fit_pw, main = "Piecewise Power Law NHPP",
     xlab = "Time")

## -----------------------------------------------------------------------------
fc <- predict_nhpp(fit_mle, time = c(800, 1000, 1500))
print(fc)

## -----------------------------------------------------------------------------
plot(fc, main = "NHPP Forecast",
     xlab = "Time")

## -----------------------------------------------------------------------------
fc_ll <- predict_nhpp(fit_ll, time = c(800, 1000))