---
title: "Censoring Schemes"
author: "Your Name"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Censoring Schemes}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
library(TKApprox)
```

## Introduction

TKApprox supports a wide range of censoring schemes commonly encountered in reliability, survival analysis, and lifetime data analysis. This vignette demonstrates how to use each censoring scheme with the package.

## Complete (Uncensored) Data

The simplest case is complete data with no censoring.

```{r}
# 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)
```

## Right-Censored Data

Right censoring occurs when we only know that an event occurred after a certain time.

```{r}
# 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)
```

## Left-Censored Data

Left censoring occurs when we only know that an event occurred before a certain time.

```{r}
# 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)
```

## Interval-Censored Data

Interval censoring occurs when we only know that an event occurred within a time interval.

```{r}
# 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)
```

## Type-I Censoring

Type-I censoring (time censoring) occurs when the experiment is terminated at a fixed time T.

```{r}
# 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)
```

## Type-II Censoring

Type-II censoring occurs when the experiment ends after a specified number of failures r.

```{r}
# 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)
```

## Progressive Type-II Censoring

Progressive Type-II censoring involves removing surviving items at each failure time.

```{r}
# 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)
```

## Hybrid Censoring

Hybrid censoring combines Type-I and Type-II: the experiment ends at min(T, r-th failure).

```{r}
# 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)
```

## Doubly Censored Data

Doubly censored data involves both left and right censoring.

```{r}
# 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)
```

## Comparing Censoring Schemes

Let's compare estimates from different censoring schemes using simulated data:

```{r, 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)
```

## Tips for Working with Censored Data

1. **CDF is required**: For all censored schemes, you must provide a CDF function
2. **Status coding**: Be careful with status coding conventions
   - Right-censored: 1 = observed, 0 = censored
   - Left-censored: 1 = observed, 0 = censored
   - Doubly-censored: -1 = left-censored, 0 = observed, 1 = right-censored
3. **Interval data**: Use a matrix/data.frame with columns `lower` and `upper`
4. **Initial values**: Good initial values are especially important for censored data
5. **Sample size**: Censoring reduces effective sample size; consider this when interpreting results

## Next Steps

- See "Loss Functions" for information on different Bayesian estimation methods
- See "Prior Specification" for advanced prior modeling