---
title: "Getting Started with AHTauDesign"
author: "Xinyue Zhao, Haitao Pan"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Getting Started with AHTauDesign}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
library(AHTauDesign)

filename <- "getting-started-results.rds"
paths <- c(file.path("..", "inst", "extdata", filename),
           file.path("inst", "extdata", filename),
           system.file("extdata", filename, package = "AHTauDesign"))
paths <- paths[nzchar(paths) & file.exists(paths)]
if (!length(paths)) {
  stop("Missing precomputed vignette data. Run the one-time preparation chunk.")
}
example <- readRDS(paths[1])
if (!identical(example$example_version, "two-branch-202610")) {
  stop("The vignette RDS is from an older example. Regenerate it once.")
}
reference <- example$reference
alternative <- example$alternative
stopifnot(inherits(reference, "ahtau_planning_report"),
          inherits(alternative, "ahtau_planning_report"))


```

## Introduction

AHTauDesign connects a clinically preferred average-hazard (AH) horizon to the
trial design expected to support it. Its prospective workflow is
**calibrate–freeze–audit–decide–reverse-design**.

Here is a **synthetic example** illustrating the main package functions. We
consider two designs under the same working event-time model. The first shows
that a clinical horizon may be supported but have insufficient power. The
second shows how Stage 2 is activated after a support failure. This vignette
is based on **precomputed simulation results** (base seed **202610**).

The example is illustrative, not a clinical recommendation. All support-rule
thresholds require calibration or validation for the intended trial setting.

## The synthetic trial

Suppose a two-arm trial uses equal allocation, 48 months of uniform accrual,
independent exponential event times, and administrative censoring. The working
control-arm 36-month survival probability is 0.60, and the treatment hazard is
0.55 times the control hazard. The preferred AH horizon is
$\tau_{\mathrm{clin}}=36$ months; a 24-month horizon is an acceptable planning
alternative only if Stage 1 fails.

Our two illustrative information states are:

| Design | Sample size | Additional follow-up | Clinical horizon |
|:--|--:|--:|--:|
| Reference design | 160 | 36 months | 36 months |
| Short-follow-up design | 160 | 6 months | 36 months |

A user may specify a different trial generator. This one is used only for the
synthetic example:

```{r trial-generator, eval=FALSE}
trial_generator <- function(N, A, F, rate_C, rate_T) {
  arm <- rep(c(0, 1), each = N / 2)
  entry <- runif(N, 0, A)
  event_time <- rexp(N, ifelse(arm == 0, rate_C, rate_T))
  censor_time <- A + F - entry
  data.frame(arm = arm, X = pmin(event_time, censor_time),
             status = as.numeric(event_time <= censor_time))
}

null_generator <- function(N, A, F, rate_C) {
  trial_generator(N, A, F, rate_C, rate_C)
}

rate_C <- -log(0.60) / 36
rate_T <- 0.55 * rate_C
```

## 1. Calibration and freeze

Calibration evaluates null-model Type I error, 95% Wald coverage, and support
availability over a prespecified design-and-horizon region. Each simulated
trial is evaluated at both horizons. Reliability classification is based on
repeated-sampling behavior, **not** on the risk-set cutoffs being calibrated.

The example compares a relative-support-only candidate with one that also
requires adequate absolute mean risk-set size. It includes a small, immature
calibration setting deliberately, so the candidates can be distinguished.

```{r calibration-code, eval=FALSE}
calibration_designs <- data.frame(
  scenario_id = c("small_early", "medium_early", "large_mature", "large_later"),
  N = c(60, 120, 240, 360),
  A = 48, F = c(0, 6, 24, 36), rate_C = -log(0.60) / 36
)

rules <- data.frame(
  rule_id = c("relative_only", "relative_plus_count"),
  rule_family = c("A", "B"),
  p_risk_threshold = c(0.10, 0.10),
  n_risk_threshold = c(NA, 30),
  p_support_threshold = c(0.95, 0.95),
  use_absolute_n = c(FALSE, TRUE)
)

cal <- simulate_calibration_null(
  design_space = calibration_designs,
  trial_generator = null_generator,
  tau_grid = c(24, 36), n_sim = 10000, seed = 202610,
  rules = rules,
  reliability_criteria = list(type1_tolerance = 0.06,
                              coverage_tolerance = 0.94,
                              support_tolerance = 0.95)
)

# Review the candidate comparison before explicitly freezing a rule.
cal$rule_evaluation$performance
frozen <- freeze_support_rule(
  diagnostic = cal$diagnostic, rules = rules,
  rule_id = "relative_plus_count",
  rationale = "Selected from qualifying rules after UPR/APR review.",
  calibration_scope = list(designs = calibration_designs,
                           tau_grid = c(24, 36),
                           null_model = "equal exponential event hazards"),
  provenance = list(source = "synthetic vignette simulation",
                    inference_conditioning = "support_valid"),
  simulation_settings = cal$settings
)
```

```{r calibration-summary, echo=FALSE}
counts <- reference$frozen_rule$reliability_counts
knitr::kable(
  counts[, c("n_total", "n_reliable", "n_unreliable", "n_indeterminate")],
  row.names = FALSE, caption = "Null-calibration reliability classifications"
)
perf <- reference$frozen_rule$performance
knitr::kable(
  perf[, c("rule_id", "n_pass", "pass_reliable", "pass_unreliable",
           "pass_indeterminate", "UPR_conservative", "APR_confirmed")],
  digits = 3, row.names = FALSE,
  caption = "Candidate support-rule comparison (illustrative calibration)"
)
th <- reference$frozen_rule$thresholds
knitr::kable(data.frame(
  selected_rule = reference$frozen_rule$rule_id,
  p_risk_min = th$p_risk_min,
  n_risk_min = th$n_risk_min,
  p_support_min = th$p_support_min
), digits = 3, row.names = FALSE,
caption = "Frozen support thresholds for this synthetic example")
```

The candidate rule uses `p_support_threshold`, while calibration reliability
uses `support_tolerance` for support availability. The displayed reliability
labels and rule-performance measures come from the precomputed null
simulations; no label was assigned by hand.

## 2. Stage 1: supported horizon, insufficient power

At the reference design, Stage 1 checks whether the fixed 36-month horizon
satisfies the frozen rule. Overall empirical power is calculated separately,
with unavailable analyses counted as nonrejections. Passing Stage 1 does not
imply that a formal power target is met.

```{r stage1-code, eval=FALSE}
raw_reference <- simulate_design(
  trial_generator = trial_generator,
  n_sim = 1000,
  tau_grid = c(24, 36),
  trial_args = list(N = 160, A = 48, F = 36,
                    rate_C = rate_C, rate_T = rate_T),
  scenario_id = "reference",
  seed = 202611
)
audit_reference <- stage1_audit(
  raw_reference, thresholds = frozen$thresholds
)
```

```{r stage1-summary, echo=FALSE}
clin <- reference$stage1_clinical_summary
knitr::kable(
  clin[, c("tau", "n_sim", "min_mean_p_risk", "min_mean_n_risk",
           "p_support", "Stage1_pass", "power_overall")],
  digits = 3, row.names = FALSE,
  caption = "Reference design: support and overall power at 36 months"
)
```

In this example, the **clinical horizon is retained** because it passes the
support audit, even though its empirical power is below the illustrative
0.80 target. This is why **supportability and power must be distinguished**.
No Stage 2 selection is made for the reference design.

## 3. Stage 2: a support failure activates conditional selection

Now consider the same sample size and event-time model, but with only six
months of additional follow-up. The clinical horizon is audited first. If it
fails, the candidate set is limited to **Stage 1-supported horizons inside the
prespecified clinical window**. The candidate with the largest finite,
positive utility is selected; independent simulations then evaluate the two
fixed horizons.

```{r stage2-code, eval=FALSE}
raw_short <- simulate_design(
  trial_generator = trial_generator,
  n_sim = 1000,
  tau_grid = c(24, 36),
  trial_args = list(N = 160, A = 48, F = 6,
                    rate_C = rate_C, rate_T = rate_T),
  scenario_id = "short_followup",
  seed = 202612
)
audit_short <- stage1_audit(raw_short, thresholds = frozen$thresholds)
selection <- stage2_select(audit_short$scenario_summary,
                            tau_clin = 36, clinical_window = c(24, 36))

# Generate fresh trials; the selected horizon is held fixed.
raw_eval <- simulate_design(
  trial_generator = trial_generator,
  n_sim = 1000,
  tau_grid = c(36, selection$final_decision_summary$tau_star),
  trial_args = list(N = 160, A = 48, F = 6,
                    rate_C = rate_C, rate_T = rate_T),
  scenario_id = "independent_evaluation",
  seed = 202613
)
audit_eval <- stage1_audit(raw_eval, thresholds = frozen$thresholds)
independent_eval <- stage2_evaluate(selection, audit_eval$scenario_summary)
independent_eval$power_difference
```

```{r stage2-summary, echo=FALSE}
cand <- alternative$stage2_result$candidate_summary
knitr::kable(
  cand[, c("tau", "min_mean_n_risk", "p_support",
           "Stage1_pass", "within_clinical_window", "utility_U", "selected")],
  digits = 3, row.names = FALSE,
  caption = "Supported candidates, selection utility, and chosen horizon"
)
decision <- alternative$stage2_result$final_decision_summary
knitr::kable(
  data.frame(clinical_horizon = decision$tau_clin,
             chosen_horizon = decision$tau_star,
             decision = decision$decision_type),
  row.names = FALSE, caption = "Conditional Stage 2 decision"
)
knitr::kable(
  alternative$stage2_evaluation$eval_summary[
    , c("role", "tau", "Stage1_pass", "p_support", "power_overall")],
  digits = 3, row.names = FALSE,
  caption = "Independent evaluation of the fixed horizons"
)
knitr::kable(
  data.frame(power_difference = alternative$stage2_evaluation$power_difference),
  digits = 3, row.names = FALSE,
  caption = "Power difference: selected alternative minus clinical horizon"
)
```

The reported `power_difference` retains its sign and compares power at two
**different** AH estimands; it is not a precision comparison for one estimand.

Changing $\tau$ changes the estimand, not just the precision. The Stage 2
selected horizon is a **planning alternative**, not an automatic clinical
replacement for the preferred horizon.

## 4. Part II: preserve the original clinical horizon

Part II addresses the resource question for the **same 36-month estimand**:
what sample size would meet both Stage 1 and an illustrative power target of
0.80? The sample-size search uses a coarse grid, a finer local grid, and
independent validation. Only validated overall power receives isotonic
stabilization; Stage 1 support diagnostics remain unsmoothed.

```{r part2-code, eval=FALSE}
part2_12 <- find_n_req(
  F_followup = 12, tau = 36, thresholds = frozen$thresholds,
  trial_generator = trial_generator,
  trial_args = list(A = 48, rate_C = rate_C, rate_T = rate_T),
  N_lower = 160, N_upper = 600, coarse_step = 20, fine_step = 10,
  n_coarse = 200, n_fine = 400, n_validate = 2000,
  power_target = 0.80,
  search_seed = 202614, validation_seed = 202615
)
part2_36 <- find_n_req(
  F_followup = 36, tau = 36, thresholds = frozen$thresholds,
  trial_generator = trial_generator,
  trial_args = list(A = 48, rate_C = rate_C, rate_T = rate_T),
  N_lower = 160, N_upper = 600, coarse_step = 20, fine_step = 10,
  n_coarse = 200, n_fine = 400, n_validate = 2000,
  power_target = 0.80,
  search_seed = 202616, validation_seed = 202617
)
```

```{r part2-summary, echo=FALSE}
req <- reference$part2_requirements
knitr::kable(
  req[, c("F_followup", "tau", "N_req", "status", "power_target",
          "raw_power", "isotonic_power")],
  digits = 3, row.names = FALSE,
  caption = "Validated sample-size requirements for the fixed 36-month horizon"
)
```

These are **simulation-based planning estimates on a 10-participant grid**,
not exact analytical minima. The illustrative simulation counts are too
modest for definitive sample-size recommendations in real trials.

## 5. Planning report

`planning_report()` combines completed results and their provenance. It does
**not** simulate trials, rerank horizons, or rerun the sample-size search.

```{r planning-report-code, eval=FALSE}
report_reference <- planning_report(
  trial_spec = list(scenario_id = "reference", N = 160, A = 48, F = 36,
                    tau_clin = 36, clinical_window = c(24, 36),
                    power_target = 0.80, survival_model = "exponential"),
  frozen_rule = frozen,
  stage1_result = audit_reference,
  part2_results = list(F12 = part2_12, F36 = part2_36),
  simulation_diagnostics = list(base_seed = 202610,
                                stage1_seed = 202611)
)
report_alternative <- planning_report(
  trial_spec = list(scenario_id = "short_followup", N = 160, A = 48,
                    F = 6, tau_clin = 36, clinical_window = c(24, 36),
                    survival_model = "exponential"),
  frozen_rule = frozen,
  stage1_result = audit_short,
  stage2_result = selection,
  stage2_evaluation = independent_eval,
  simulation_diagnostics = list(selection_seed = 202612,
                                independent_evaluation_seed = 202613)
)
```

```{r planning-report-summary, echo=FALSE}
a <- reference$decision_summary
b <- alternative$decision_summary
knitr::kable(data.frame(
  design = c("Reference", "Short follow-up"),
  Stage1_pass = c(a$Stage1_pass, b$Stage1_pass),
  power_overall = c(a$power_overall, b$power_overall),
  Stage2_decision = c(a$stage2_decision, b$stage2_decision),
  final_horizon = c(a$tau_star, b$tau_star)
), digits = 3, row.names = FALSE,
caption = "Major decisions from the two synthetic planning reports")
```

The full report objects contain frozen calibration evidence, support and
power diagnostics, conditional selection records, independent evaluation,
Part II validation traces, and simulation metadata. The tables above display
only the key planning results rather than thousands of replicate-level rows.

<!-- Developer-only preparation. Run this chunk ONCE interactively in RStudio,
     after devtools::load_all(). It is never run when the vignette is knitted.
     This deliberately uses no external helper script and adds no package R file. -->

```{r prepare-vignette-data, eval=FALSE, include=FALSE}
package_root <- if (file.exists("DESCRIPTION")) "." else ".."
stopifnot(file.exists(file.path(package_root, "DESCRIPTION")))
base_seed <- 202610

trial_generator <- function(N, A, F, rate_C, rate_T) {
  stopifnot(N >= 2, N %% 2 == 0)
  arm <- rep(c(0, 1), each = N / 2)
  entry <- runif(N, 0, A)
  event_time <- rexp(N, ifelse(arm == 0, rate_C, rate_T))
  censor_time <- A + F - entry
  data.frame(arm = arm, X = pmin(event_time, censor_time),
             status = as.numeric(event_time <= censor_time))
}
null_generator <- function(N, A, F, rate_C) {
  trial_generator(N, A, F, rate_C, rate_C)
}
rate_C <- -log(0.60) / 36
rate_T <- 0.55 * rate_C

# Calibration includes small, immature settings and larger, mature settings.
calibration_designs <- data.frame(
  scenario_id = c("small_early", "medium_early", "large_mature", "large_later"),
  N = c(60, 120, 240, 360), A = 48,
  F = c(0, 6, 24, 36), rate_C = rate_C
)
rules <- data.frame(
  rule_id = c("relative_only", "relative_plus_count"),
  rule_family = c("A", "B"),
  p_risk_threshold = c(0.10, 0.10),
  n_risk_threshold = c(NA, 30),
  p_support_threshold = c(0.95, 0.95),
  use_absolute_n = c(FALSE, TRUE)
)
cal <- simulate_calibration_null(
  design_space = calibration_designs,
  trial_generator = null_generator, tau_grid = c(24, 36),
  n_sim = 10000, seed = base_seed,
  rules = rules,
  reliability_criteria = list(type1_tolerance = 0.06,
                              coverage_tolerance = 0.94,
                              support_tolerance = 0.95)
)
print(cal$reliability_counts)
print(cal$rule_evaluation$performance[, c("rule_id", "n_pass",
                                               "UPR_conservative", "APR_confirmed")])
if (!any(cal$classified_diagnostic$reliability_status != "Reliable")) {
  stop("Calibration did not show a non-Reliable setting; review the example.")
}
if (identical(as.vector(cal$rule_evaluation$pass_matrix[, 1]),
              as.vector(cal$rule_evaluation$pass_matrix[, 2]))) {
  stop("Candidate rules gave identical decisions; revise calibration space.")
}

# Rule must qualify under the actual Monte Carlo classifications.
frozen <- freeze_support_rule(
  diagnostic = cal$diagnostic, rules = rules,
  rule_id = "relative_plus_count",
  rationale = "Qualifying rule selected after comparing conservative UPR and confirmed APR.",
  calibration_scope = list(designs = calibration_designs,
                           tau_grid = c(24, 36),
                           null_model = "equal exponential hazards"),
  provenance = list(source = "synthetic vignette simulation",
                    inference_conditioning = "support_valid"),
  simulation_settings = cal$settings
)

# Stage 1: supported but below target power.
raw_reference <- simulate_design(
  trial_generator = trial_generator,
  n_sim = 1000,
  tau_grid = c(24, 36),
  trial_args = list(N = 160, A = 48, F = 36,
                    rate_C = rate_C, rate_T = rate_T),
  scenario_id = "reference",
  seed = base_seed + 1
)
audit_reference <- stage1_audit(raw_reference, thresholds = frozen$thresholds)
clinical_reference <- subset(audit_reference$scenario_summary, tau == 36)
if (!isTRUE(clinical_reference$Stage1_pass) ||
    !(clinical_reference$power_overall < 0.80)) {
  stop("Reference example did not exhibit Stage 1 Pass with power < 0.80.")
}

# Independent Stage 2 selection scenario: short follow-up causes support failure.
raw_short <- simulate_design(
  trial_generator = trial_generator,
  n_sim = 1000,
  tau_grid = c(24, 36),
  trial_args = list(N = 160, A = 48, F = 6,
                    rate_C = rate_C, rate_T = rate_T),
  scenario_id = "short_followup",
  seed = base_seed + 2
)
audit_short <- stage1_audit(raw_short, thresholds = frozen$thresholds)
selection <- stage2_select(audit_short$scenario_summary,
                            tau_clin = 36, clinical_window = c(24, 36))
if (!identical(selection$final_decision_summary$decision_type,
               "Utility-driven alternative")) {
  stop("Stage 2 example did not select an eligible planning alternative.")
}

# Fresh simulations only for the two FIXED horizons.
chosen_tau <- selection$final_decision_summary$tau_star
raw_eval <- simulate_design(
  trial_generator = trial_generator,
  n_sim = 1000,
  tau_grid = c(36, chosen_tau),
  trial_args = list(N = 160, A = 48, F = 6,
                    rate_C = rate_C, rate_T = rate_T),
  scenario_id = "independent_evaluation",
  seed = base_seed + 3
)
audit_eval <- stage1_audit(raw_eval, thresholds = frozen$thresholds)
independent_eval <- stage2_evaluate(selection, audit_eval$scenario_summary)

# Part II keeps the 36-month AH estimand fixed. Simulation counts are illustrative.
part2_12 <- find_n_req(
  F_followup = 12, tau = 36, thresholds = frozen$thresholds,
  trial_generator = trial_generator,
  trial_args = list(A = 48, rate_C = rate_C, rate_T = rate_T),
  N_lower = 160, N_upper = 600, coarse_step = 20, fine_step = 10,
  n_coarse = 200, n_fine = 400, n_validate = 2000,
  power_target = 0.80, search_seed = base_seed + 4,
  validation_seed = base_seed + 5
)
part2_36 <- find_n_req(
  F_followup = 36, tau = 36, thresholds = frozen$thresholds,
  trial_generator = trial_generator,
  trial_args = list(A = 48, rate_C = rate_C, rate_T = rate_T),
  N_lower = 160, N_upper = 600, coarse_step = 20, fine_step = 10,
  n_coarse = 200, n_fine = 400, n_validate = 2000,
  power_target = 0.80, search_seed = base_seed + 6,
  validation_seed = base_seed + 7
)

report_reference <- planning_report(
  trial_spec = list(scenario_id = "reference", N = 160, A = 48, F = 36,
                    allocation = "1:1", tau_clin = 36,
                    clinical_window = c(24, 36), power_target = 0.80,
                    survival_model = "exponential"),
  frozen_rule = frozen, stage1_result = audit_reference,
  part2_results = list(F12 = part2_12, F36 = part2_36),
  simulation_diagnostics = list(
    base_seed = base_seed, stage1_seed = base_seed + 1,
    part2_F12_search_seed = base_seed + 4,
    part2_F12_validation_seed = base_seed + 5,
    part2_F36_search_seed = base_seed + 6,
    part2_F36_validation_seed = base_seed + 7,
    stage1_n_sim = 1000, note = "Illustrative Monte Carlo planning example")
)
report_alternative <- planning_report(
  trial_spec = list(scenario_id = "short_followup", N = 160, A = 48, F = 6,
                    allocation = "1:1", tau_clin = 36,
                    clinical_window = c(24, 36), survival_model = "exponential"),
  frozen_rule = frozen, stage1_result = audit_short,
  stage2_result = selection, stage2_evaluation = independent_eval,
  simulation_diagnostics = list(
    base_seed = base_seed, selection_seed = base_seed + 2,
    evaluation_seed = base_seed + 3,
    selection_n_sim = 1000, evaluation_n_sim = 1000)
)

# Keep only the completed, auditable summary objects. Raw null-calibration
# replicate results are deliberately not bundled with the vignette.
example <- list(
  example_version = "two-branch-202610",
  reference = report_reference,
  alternative = report_alternative
)
output_dir <- file.path(package_root, "inst", "extdata")
dir.create(output_dir, recursive = TRUE, showWarnings = FALSE)
output_file <- file.path(output_dir, "getting-started-results.rds")
saveRDS(example, output_file, compress = "xz")
stopifnot(file.exists(output_file))
cat("SUCCESS: updated precomputed vignette data:",
    normalizePath(output_file), "\n")
```
