---
title: "Semiparametric Bayesian Regression for Dependent Current Status Data"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Semiparametric Bayesian Regression for Dependent Current Status Data}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

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

## Introduction

The **`SemiParamBernsteinDepCS`** package provides a semiparametric Bayesian regression framework using Bernstein polynomial baseline models for analyzing dependent current status data. The methodology accommodates proportional hazards (PH) and proportional odds (PO) regression models paired with Archimedean copulas (Gumbel, Frank, and Clayton) to capture dependence between event time and observation/censoring time.

For full mathematical details, see Sharma and Balakrishnan (2026) <doi:10.1080/02664763.2026.2701921>.

## Simulating Dependent Current Status Data

We generate a sample dataset under a Proportional Hazards (PH) model with a Gumbel copula:

```{r sim_data}
library(SemiParamBernsteinDepCS)

set.seed(2026)
sim_data <- sim_bernstein_depcs(
  n = 100,
  model_type = "PH",
  copula = "gumbel",
  beta_Y = c(0.5, -0.3),
  beta_T = c(-0.2, 0.4)
)

head(sim_data)
```

## Model Fitting

We fit the semiparametric Bayesian regression model using `fit_semiparam_bernstein_depcs()`:

```{r fit_model}
fit_ph <- fit_semiparam_bernstein_depcs(
  formula_Y = delta ~ x1 + x2,
  formula_T = time ~ x1 + x2,
  data = sim_data,
  model_type = "PH",
  copula = "gumbel",
  order_m = 2,
  n_iter = 400,
  n_burn = 100,
  seed = 42
)

print(fit_ph)
```

## Summary & Model Diagnostics

Extract posterior means, standard deviations, and 95% HPD credible intervals:

```{r summary_model}
summary(fit_ph)
```

## Prediction & Forecasting

Compute predicted marginal survival curves and HPD bounds for new covariate values:

```{r predict_model}
pred_res <- predict(fit_ph, times = seq(0.1, 2.5, length.out = 20))
head(pred_res$predictions_Y)
```

## Plotting

Visualize predicted marginal survival curves:

```{r plot_model, fig.width=6, fig.height=4}
plot(fit_ph, type = "survival")
```

## Real Data Example (Primary Biliary Cirrhosis)

The package includes the benchmark `pbc_depcs` dataset:

```{r pbc_example}
data(pbc_depcs)
head(pbc_depcs)
```