## ----setup, include=FALSE-----------------------------------------------------
knitr::opts_chunk$set(echo = TRUE)
library(TKApprox)

## -----------------------------------------------------------------------------
# Define exponential distribution
pdf_exp <- function(x, param) dexp(x, rate = param)
cdf_exp <- function(x, param) pexp(x, rate = param)

# Prior specification
prior_spec <- list(rate = list(family = "gamma", hyperparameters = list(shape = 2, rate = 1)))

# Generate complete data
set.seed(123)
data <- rexp(20, rate = 1.5)

# Fit with complete data
fit_complete <- tk_fit(
  data = data,
  censoring_scheme = "complete",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel"
)

summary(fit_complete)

## -----------------------------------------------------------------------------
# Create right-censored data
# status = 1: observed, status = 0: right-censored
data <- c(1.2, 2.3, 1.8, 3.1, 0.9, 2.5, 1.5, 3.8, 2.0, 1.7)
status <- c(1, 1, 0, 1, 0, 1, 1, 0, 1, 1)

fit_right <- tk_fit(
  data = data,
  censoring_scheme = "right-censored",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  status = status
)

summary(fit_right)

## -----------------------------------------------------------------------------
# Create left-censored data
# status = 1: observed, status = 0: left-censored
data <- c(1.2, 2.3, 1.8, 3.1, 0.9, 2.5, 1.5, 3.8, 2.0, 1.7)
status <- c(1, 0, 1, 1, 0, 1, 1, 0, 1, 1)

fit_left <- tk_fit(
  data = data,
  censoring_scheme = "left-censored",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  status = status
)

summary(fit_left)

## -----------------------------------------------------------------------------
# Create interval-censored data
# Each row: [lower, upper]
# If lower == upper, it's an exact observation
data <- cbind(
  lower = c(1.0, 2.0, 1.5, 2.5, 1.2, 2.0, 1.8, 3.0, 1.5, 2.2),
  upper = c(1.5, 2.5, 1.5, 3.0, 1.8, 2.5, 2.2, 3.5, 2.0, 2.5)
)

fit_interval <- tk_fit(
  data = data,
  censoring_scheme = "interval-censored",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel"
)

summary(fit_interval)

## -----------------------------------------------------------------------------
# Generate Type-I censored data
set.seed(123)
true_rate <- 1.5
censoring_time <- 2.0

# Simulate failure times
failure_times <- rexp(20, rate = true_rate)

# Apply Type-I censoring
data <- pmin(failure_times, censoring_time)

fit_type1 <- tk_fit(
  data = data,
  censoring_scheme = "type-i",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  censoring_time = censoring_time
)

summary(fit_type1)

## -----------------------------------------------------------------------------
# Generate Type-II censored data
set.seed(123)
n <- 20  # total items
r <- 10  # number of failures to observe

# Simulate failure times
failure_times <- sort(rexp(n, rate = 1.5))

# Observe only first r failures
data <- failure_times[1:r]

fit_type2 <- tk_fit(
  data = data,
  censoring_scheme = "type-ii",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  n = n,
  r = r
)

summary(fit_type2)

## -----------------------------------------------------------------------------
# Generate progressive Type-II censored data
set.seed(123)
n <- 20
m <- 10  # number of observed failures

# Simulate failure times
failure_times <- sort(rexp(n, rate = 1.5))

# Specify removal scheme (remove 1 item at each failure)
removals <- rep(1, m)

# Adjust for remaining items
data <- failure_times[1:m]

fit_progressive <- tk_fit(
  data = data,
  censoring_scheme = "progressive-type2",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  removals = removals,
  n = n
)

summary(fit_progressive)

## -----------------------------------------------------------------------------
# Generate hybrid censored data
set.seed(123)
n <- 20
r <- 10
censoring_time <- 2.0

# Simulate failure times
failure_times <- sort(rexp(n, rate = 1.5))

# Apply hybrid censoring
if (failure_times[r] < censoring_time) {
  # Type-II censoring (r failures occur before T)
  data <- failure_times[1:r]
} else {
  # Type-I censoring (experiment ends at T)
  data <- failure_times[failure_times < censoring_time]
}

fit_hybrid <- tk_fit(
  data = data,
  censoring_scheme = "hybrid",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  censoring_time = censoring_time,
  r = r,
  n = n
)

summary(fit_hybrid)

## -----------------------------------------------------------------------------
# Create doubly censored data
# status = -1: left-censored, status = 0: observed, status = 1: right-censored
data <- c(1.2, 2.3, 1.8, 3.1, 0.9, 2.5, 1.5, 3.8, 2.0, 1.7)
status <- c(-1, 1, 0, 1, -1, 1, 0, 1, 0, 1)

fit_doubly <- tk_fit(
  data = data,
  censoring_scheme = "doubly-censored",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  status = status
)

summary(fit_doubly)

## ----warning=FALSE------------------------------------------------------------
set.seed(123)
true_rate <- 1.5
n <- 30

# Generate complete data
complete_data <- rexp(n, rate = true_rate)

# Fit complete data
fit_complete <- tk_fit(
  data = complete_data,
  censoring_scheme = "complete",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel"
)

# Create right-censored data
status_right <- c(rep(1, 20), rep(0, 10))
fit_right <- tk_fit(
  data = complete_data,
  censoring_scheme = "right-censored",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  status = status_right
)

# Create Type-I censored data
censoring_time <- median(complete_data)
data_type1 <- pmin(complete_data, censoring_time)
fit_type1 <- tk_fit(
  data = data_type1,
  censoring_scheme = "type-i",
  pdf = pdf_exp,
  cdf = cdf_exp,
  prior_spec = prior_spec,
  initial_values = c(rate = 1),
  loss_function = "sel",
  censoring_time = censoring_time
)

# Compare estimates
comparison <- data.frame(
  Scheme = c("Complete", "Right-Censored", "Type-I"),
  Estimate = c(coef(fit_complete), coef(fit_right), coef(fit_type1)),
  SE = c(fit_complete$standard_errors, fit_right$standard_errors, fit_type1$standard_errors),
  True = true_rate
)

print(comparison)