## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 7,
  fig.height = 4
)

## ----quick-start--------------------------------------------------------------
library(serieshaz)

# Define three failure mode components
disk    <- dfr_weibull(shape = 2, scale = 500)
memory  <- dfr_exponential(0.001)
psu     <- dfr_gompertz(a = 0.0001, b = 0.02)

# Compose into a series system
server <- dfr_dist_series(list(disk, memory, psu))
print(server)

## ----layout-------------------------------------------------------------------
param_layout(server)

## ----hazard-surv--------------------------------------------------------------
h <- hazard(server)
S <- surv(server)

# Evaluate at t = 100 hours
cat(sprintf("System hazard at t=100:  %.6f\n", h(100)))
cat(sprintf("System survival at t=100: %.4f\n", S(100)))

## ----cdf-density--------------------------------------------------------------
F_sys <- cdf(server)
f_sys <- density(server)

cat(sprintf("P(T <= 100) = %.4f\n", F_sys(100)))
cat(sprintf("f(100)      = %.6f\n", f_sys(100)))

## ----sampling-----------------------------------------------------------------
samp <- sampler(server)
set.seed(42)
times <- samp(5)
times

## ----loglik-------------------------------------------------------------------
ll <- loglik(server)
df <- data.frame(t = c(100, 200, 150, 300, 250), delta = c(1, 1, 0, 1, 0))
ll(df)

## ----fit-demo-----------------------------------------------------------------
# Simpler 2-component model for a clean fit demo
fit_model <- dfr_dist_series(list(
    dfr_weibull(shape = 2, scale = 100),
    dfr_exponential(0.01)
))
true_par <- c(2, 100, 0.01)

set.seed(42)
fit_samp <- sampler(fit_model)
train <- data.frame(t = fit_samp(500), delta = rep(1L, 500))

solver <- fit(fit_model)
result <- suppressWarnings(solver(train, par = c(1.5, 80, 0.02)))

cat("True parameters:  ", true_par, "\n")
cat("Fitted parameters:", round(coef(result), 3), "\n")

## ----introspection------------------------------------------------------------
ncomponents(server)  # 3

# Extract the Weibull (disk) component
disk_comp <- component(server, 1)
params(disk_comp)

## ----hazard-decomp------------------------------------------------------------
# Get per-component hazard closures
h1 <- component_hazard(server, 1)  # disk
h2 <- component_hazard(server, 2)  # memory
h3 <- component_hazard(server, 3)  # PSU

# At t = 200, which component contributes most?
t0 <- 200
hazards <- c(Disk = h1(t0), Memory = h2(t0), PSU = h3(t0))
cat(sprintf("Disk:   %.6f (%.1f%%)\n", hazards[1], 100 * hazards[1] / sum(hazards)))
cat(sprintf("Memory: %.6f (%.1f%%)\n", hazards[2], 100 * hazards[2] / sum(hazards)))
cat(sprintf("PSU:    %.6f (%.1f%%)\n", hazards[3], 100 * hazards[3] / sum(hazards)))
cat(sprintf("System: %.6f\n", h(t0)))
cat(sprintf("Sum:    %.6f (matches system)\n", sum(hazards)))

## ----failure-attr-------------------------------------------------------------
set.seed(42)
mat <- sample_components(server, n = 10000)

# System lifetime = min across components
t_sys <- apply(mat, 1, min)

# Which component caused each failure?
failing <- apply(mat, 1, which.min)
cat("Failure proportions:\n")
round(table(failing) / length(failing), 3)