---
title: "An Analysis Workflow with netmatchRI"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{An Analysis Workflow with netmatchRI}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE, comment = "#>",
  fig.width = 10, fig.height = 5, out.width = "75%", fig.align = "center", dpi = 150
)
```

`netmatchRI` provides network-constrained matching and randomization-based
inference for observational network data. Dual-penalty matching minimizes total
covariate dissimilarity while enforcing penalties on network proximity between
and within treatment arms in matched sets. Randomization-based inference
combines evidence from correlated matched sets based on this design.

## Data

The data frame should contain one row per unit, a binary treatment, baseline
covariates, and an outcome. The network may be supplied as an adjacency or
graph-distance matrix, with rows and columns in the same order as the data.
The bundled example contains 300 units simulated under no treatment effect
and network dependence. It was selected to illustrate how accounting for
dependence can change the inferential conclusion.

```{r}
library(netmatchRI)

data("netmatch_example", package = "netmatchRI")
sim <- netmatch_example
dat <- sim$data
net_dist <- sim$net_dist

dim(dat)
table(dat$Z)
```

## Match

Construct a dual-penalty design with `kappa = 2`, so that every pair of units
within a matched set is separated by a network distance greater than two.
HiGHS is the default open-source solver; Gurobi and GLPK can be selected
with the `solver` argument.

```{r}
m <- netmatch(
  data = dat,
  treat = "Z",
  covariates = c("X1", "X2", "X3"),
  network = net_dist,
  method = "dual",
  kappa = 2,
  estimand = "ATT"
)

m
```

## Summary and diagnostics

The summary reports matched sample sizes, effective sample sizes, matched-set
sizes, balance, network distances, and solver status. `diagnose_match()` checks
covariate balance and network distance within matched sets.

```{r}
summary(m)

diag <- diagnose_match(m)
diag$covariate_balance
diag$network_summary
```

## Compare matching designs

Pass all three matching designs as a named list. The two plots compare
covariate balance and covariate dissimilarity across these designs. Set
`include_before = FALSE` to omit the original unmatched sample.

```{r, fig.width=7, fig.height=6}
m_covariate <- netmatch(
  dat, "Z", c("X1", "X2", "X3"), net_dist,
  method = "covariate"
)
m_single <- netmatch(
  dat, "Z", c("X1", "X2", "X3"), net_dist,
  method = "single", kappa = 2
)
designs <- list(
  "Covariate-only" = m_covariate,
  "Single-penalty" = m_single,
  "Dual-penalty" = m
)
plot_covariate_balance(designs, include_before = TRUE)
plot_covariate_similarity(designs, include_before = TRUE)
```

## Randomization-based inference

We use the dual-penalty design `m` to test the sharp null that no unit's outcome
is affected by any unit's treatment. This includes both direct effects and
treatment interference. The tests combine matched-set Mann–Whitney U statistics
and calculate two-sided p-values using a normal approximation.

`RI_unadjusted()` ignores covariance across matched sets. `RI_design()` uses
the design-based covariance bounds, while `RI_adjusted()` scales these bounds
using the specified `eta` and `rho`. These parameters describe assumed bounds
on covariance, rather than estimates of covariance from the data.

```{r}
ri_unadjusted <- RI_unadjusted(m, "Y")
ri_design <- RI_design(m, "Y")
ri_adjusted <- RI_adjusted(m, "Y", eta = 0.05, rho = 0.20)

p_values <- data.frame(
  inference = c("Unadjusted", "Design-based", "Adjusted (eta = 0.05, rho = 0.20)"),
  p_value = c(ri_unadjusted$result$p_value,
              ri_design$result$p_value,
              ri_adjusted$result$p_value)
)
p_values
```

## Sensitivity

Evaluate the p-value over a grid of eta and rho values, and compute the critical curve for a
chosen test level. Both plots below use the dual-penalty design `m`;
the dashed black line in the p-value plot is its unadjusted p-value, and
the dotted red line marks the 0.05 test level.

```{r}
sens <- sensitivity_grid(
  m,
  "Y",
  eta = c(0, 0.01, 0.05, 0.10, 0.15, 0.20),
  rho = seq(0, 0.5, by = 0.05)
)
plot_sensitivity(sens)

crit <- critical_sensitivity(m, "Y", rho = seq(0, 0.5, by = 0.01))
crit
crit$summary
plot_sensitivity(crit)
```

## Interpretation

The sensitivity grid shows how the p-value changes as `eta` scales the
design-based covariance bounds and `rho` controls their decay with the minimum
network distance between matched sets. Smaller `rho` means faster decay;
`rho = 1` means no decay. Setting `eta = 0` gives the unadjusted analysis,
while `eta = rho = 1` gives the design-based analysis.

The critical curve gives the value of `eta` at each `rho` for which the
p-value equals `alpha = 0.05`. Within `[0, 1]`, values below
the curve retain rejection and values above it do not.

The summary reports the critical eta under no decay (`rho = 1`) and the
critical ratio. The latter compares the average weighted covariance bound
across relevant pairs of matched sets with the average weighted within-set
variance at the critical boundary. If the unadjusted p-value is greater than
`alpha`, there is no critical boundary.

The inferential assumptions are described in Dong and Lee (2026), cited by
`citation("netmatchRI")`. In particular, the adjusted analysis requires the
specified sensitivity bounds to cover the cross-matched-set covariances.