Introduction to antedep

Overview

antedep supports three model families:

AD (Gaussian) for continuous longitudinal outcomes
INAD for count longitudinal outcomes
CAT (categorical AD) for discrete state longitudinal outcomes

Packaged datasets include bolus_inad, cattle_growth, cochlear_implant, labor_force_cat, and race_100km.

data("race_100km")
dim(race_100km$y)
colMeans(race_100km$y)

Production-Readiness Matrix

Model	Data type	Complete-data fit/logLik	Missing-data fit/logLik	Notes
AD	Continuous	Ready	Ready (`fit_gau`, `logL_gau`)	Missing-data fit supports EM / observed-data likelihood modes
INAD	Counts	Ready	Ready (`fit_inad`, `logL_inad`)	Missing-data fit supports `na_action = "marginalize"`
CAT	Categorical states	Ready	Ready (`fit_cat`, `logL_cat`)	Missing-data fit supports orders 0, 1, 2

Visualizing Longitudinal Data

data("bolus_inad")
plot_profile(
  bolus_inad$y,
  blocks = bolus_inad$blocks,
  block_labels = c("Group 1", "Group 2"),
  title = "Bolus Longitudinal Profiles"
)

# More detailed dependence diagnostic
plot_prism(bolus_inad$y)

AD Workflow (Continuous Outcomes)

Complete-data fit and order comparison

y_gau <- simulate_gau(n_subjects = 20, n_time = 4, order = 1, phi = 0.5)
fit_gau1 <- fit_gau(y_gau, order = 1)
fit_gau2 <- fit_gau(y_gau, order = 2)

logLik(fit_gau1)
#> 'log Lik.' -105.672 (df=7)
BIC(fit_gau1, fit_gau2)
#>          df      BIC
#> fit_gau1  7 232.3141
#> fit_gau2  9 235.8026
bic_order_gau(y_gau, max_order = 2)
#> BIC-based Order Selection for Gaussian AD
#> ==========================================
#> 
#>  order     log_l      bic selected
#>      0 -112.6549 249.2756    FALSE
#>      1 -105.6720 244.2971     TRUE
#>      2 -104.4205 247.7855    FALSE
#> 
#> Best order by BIC: 1

# Example AD order test
test_order_gau(y_gau, p = 0, q = 1)$p_value
#> [1] 0.002952214

Missing-data AD fit

y_gau_miss <- y_gau
y_gau_miss[sample(length(y_gau_miss), 8)] <- NA

fit_gau_em <- fit_gau(y_gau_miss, order = 1, na_action = "em", em_max_iter = 10)
fit_gau_cc <- fit_gau(y_gau_miss, order = 1, na_action = "complete")
#> Warning: Only 13/20 subjects have complete data (65%).
#> Consider using na_action = 'em' for less information loss.

c(em_logLik = fit_gau_em$log_l, complete_case_logLik = fit_gau_cc$log_l)
#>            em_logLik complete_case_logLik 
#>            -94.00107            -63.11882
summary(fit_gau_em)
#> Gaussian Antedependence Model -- Summary
#> ========================================
#> 
#>   Order            : 1
#>   Subjects (n)     : 20
#>   Time points (T)  : 4
#>   Parameters       : 7
#>   Missing values   : 8 (10.0%)
#> 
#> Fit statistics:
#>   Log-likelihood   : -94.0011
#>   AIC              : 202.0021
#>   BIC              : 208.9723
#> 
#> Estimated mean (mu):
#> [1]  0.0961  0.1565  0.1244 -0.1773
#> 
#> Antedependence parameters (phi):
#> [1] 0.6772 0.1696 0.4342
#> 
#> Innovation std dev (sigma):
#> [1] 0.9905 0.6902 1.0489 0.8335

AD mean and covariance structure tests (complete data)

These AD inference tests currently require complete data.

# One-sample and two-sample mean-profile tests
mu0 <- colMeans(y_gau)
blocks_gau <- rep(1:2, each = nrow(y_gau) / 2)

test_one_sample_gau(y_gau, mu0 = mu0, p = 1)$p_value
#> [1] 1
test_two_sample_gau(y_gau, blocks = blocks_gau, p = 1)$p_value
#> [1] 0.01854549

# Wald contrast test: H0 that adjacent means are equal
C <- diff(diag(ncol(y_gau)))
test_contrast_gau(y_gau, C = C, p = 1)$p_value
#> [1] 0.8183341

# Covariance homogeneity across groups
test_homogeneity_gau(y_gau, blocks = blocks_gau, p = 1)$p_value
#>        g1 
#> 0.3237498

INAD Workflow (Count Outcomes)

Complete-data INAD fit

y_inad <- simulate_inad(
  n_subjects = 20,
  n_time = 4,
  order = 1,
  thinning = "binom",
  innovation = "pois",
  alpha = 0.4,
  theta = 3
)

fit_inad1 <- fit_inad(y_inad, order = 1, thinning = "binom", innovation = "pois")
logLik(fit_inad1)
#> 'log Lik.' -158.3541 (df=7)
bic_order_inad(y_inad, max_order = 2, thinning = "binom", innovation = "pois")
#> 
#> BIC Model Order Comparison for INAD
#> ====================================
#> 
#>  order logLik n_params   BIC
#>      0 -164.7        4 341.3
#>      1 -158.4        7 337.7
#>      2 -158.4        9 343.7
#> 
#> Best order by BIC: 1

Missing-data INAD fit

y_inad_miss <- y_inad
y_inad_miss[sample(length(y_inad_miss), 10)] <- NA

fit_inad_miss <- fit_inad(
  y_inad_miss,
  order = 1,
  thinning = "binom",
  innovation = "pois",
  na_action = "marginalize",
  max_iter = 10
)

summary(fit_inad_miss)
#> INAD Antedependence Model -- Summary
#> ====================================
#> 
#>   Order            : 1
#>   Thinning         : binom
#>   Innovation       : pois
#>   Subjects (n)     : 20
#>   Time points (T)  : 4
#>   Parameters       : 7
#>   Missing values   : 10 (12.5%)
#> 
#> Fit statistics:
#>   Log-likelihood   : -142.2786
#>   AIC              : 298.5573
#>   BIC              : 305.5274
#> 
#> Estimated innovation mean (theta):
#> [1] 3.3775 2.0946 2.4217 3.2136
#> 
#> Estimated thinning parameter (alpha):
#> [1] 0.0000 0.7107 0.5782 0.1964

# Missing-data LRTs (observed-data likelihood)
blocks_inad <- rep(1:2, each = nrow(y_inad_miss) / 2)
# Order comparisons are boundary problems: compare BIC (no calibrated p-value)
ord_inad <- test_order_inad(y_inad_miss, order_null = 0, order_alt = 1, blocks = blocks_inad)
c(LR = ord_inad$lrt_stat, BIC_order0 = ord_inad$bic_null, BIC_order1 = ord_inad$bic_alt)
#>         LR BIC_order0 BIC_order1 
#>   6.743751 306.205451 308.448896
test_homogeneity_inad(y_inad_miss, blocks = blocks_inad, order = 1, test = "mean")$p_value
#> [1] 0.7852656

CAT Workflow (Categorical Outcomes)

Complete-data CAT fit and inference

y_cat <- simulate_cat(n_subjects = 20, n_time = 4, order = 1, n_categories = 3)
fit_cat1 <- fit_cat(y_cat, order = 1)
logLik(fit_cat1)
#> 'log Lik.' -79.08119 (df=20)

# Wald CIs (complete data)
head(confint(fit_cat1), 3)
#>                lower     upper
#> pi_t1_cat1 0.0602273 0.4397727
#> pi_t1_cat2 0.1852967 0.6147033
#> pi_t1_cat3 0.1409627 0.5590373

# Order tests for CAT
run_order_tests_cat(y_cat, max_order = 2)
#> Sequential Order Tests for Categorical AD
#> =========================================
#> 
#>      comparison order_null order_alt
#>  AD(0) vs AD(1)          0         1
#>  AD(1) vs AD(2)          1         2
#>                                                                   method
#>  Likelihood Ratio Test for Categorical AD Order (H0: AD(0) vs H1: AD(1))
#>  Likelihood Ratio Test for Categorical AD Order (H0: AD(1) vs H1: AD(2))
#>  lrt_stat df p_value significant
#>   13.0444 12  0.3658       FALSE
#>   37.2866 24  0.0410        TRUE
#> 
#> Selected order: 0

Missing-data CAT fit

y_cat_miss <- y_cat
y_cat_miss[sample(length(y_cat_miss), 12)] <- NA

fit_cat_miss <- fit_cat(y_cat_miss, order = 1, na_action = "marginalize")
logLik(fit_cat_miss)
#> 'log Lik.' -65.91133 (df=20)
fit_cat_miss$settings$na_action_effective
#> [1] "marginalize"

# EM alternative for CAT missing data (orders 0/1)
fit_cat_miss_em <- fit_cat(y_cat_miss, order = 1, na_action = "em", em_max_iter = 10)
#> Warning in em_cat(y = y, order = p, blocks = blocks_input, homogeneous =
#> homogeneous, : EM did not converge after 10 iterations (change = 0.000809).
logLik(fit_cat_miss_em)
#> 'log Lik.' -65.91256 (df=20)

# Missing-data order/homogeneity tests are supported
blocks_cat <- rep(1:2, each = nrow(y_cat_miss) / 2)
test_order_cat(y_cat_miss, order_null = 0, order_alt = 1)$p_value
#> [1] 0.2500927
test_homogeneity_cat(y_cat_miss, blocks = blocks_cat, order = 1)$p_value
#> [1] 0.3780152

Additional Inference Helpers (Complete Data)

# INAD homogeneity/stationarity tools
blocks <- rep(1:2, each = nrow(y_inad) / 2)
run_homogeneity_tests_inad(y_inad, blocks = blocks, order = 1)
run_stationarity_tests_inad(y_inad, order = 1, blocks = blocks)

# CAT stationarity/time-invariance tools
test_timeinvariance_cat(y_cat, order = 1)
test_stationarity_cat(y_cat, order = 1)

Known Limitations

EM entry points: em_gau (Gaussian), em_inad (INAD), and em_cat (CAT, orders 0/1) are available; for CAT order 2 with missing data, use fit_cat(na_action = "marginalize").
For CAT missingness, fit_cat(na_action = "marginalize") supports orders 0/1/2, while fit_cat(na_action = "em") (or em_cat) supports orders 0/1 with optional safeguard step-halving.
Missing-data confidence intervals are not implemented yet (ci_gau, ci_inad, ci_cat require complete-data fits).
AD missing-data LRT/mean/covariance tests are not implemented yet.
CAT missing-data time-invariance/stationarity tests (test_timeinvariance_cat, test_stationarity_cat, run_stationarity_tests_cat) are not implemented yet.
Current order support is up to order 2 for AD, INAD, and CAT.
CAT missing-data marginalization currently supports orders 0, 1, and 2.