The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.
antedep supports three model families:
Packaged datasets include bolus_inad,
cattle_growth, cochlear_implant,
labor_force_cat, and race_100km.
| 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 |
data("bolus_inad")
plot_profile(
bolus_inad$y,
blocks = bolus_inad$blocks,
block_labels = c("Group 1", "Group 2"),
title = "Bolus Longitudinal Profiles"
)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: 1y_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.8335These 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.01854549y_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: 1y_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.7852656y_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: 0y_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)# 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)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").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.ci_gau, ci_inad, ci_cat require
complete-data fits).test_timeinvariance_cat,
test_stationarity_cat,
run_stationarity_tests_cat) are not implemented yet.Use the local one-command workflow before pushing:
These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.