Provides a user-friendly interface for fitting Gaussian mixture models with optional covariates on class membership probabilities and various covariance structures.
Usage
rtmb_mixture(
formula,
k = 2,
data = NULL,
covariance = c("diagonal", "diagonal_equal", "full", "full_equal", "full_equal_corr"),
prior = prior_flat(),
y_range = NULL,
fixed = NULL,
WAIC = FALSE,
...
)Arguments
- formula
A formula specifying the response variable(s). For multivariate, use
cbind(y1, y2) ~ 1.- k
Number of mixture components.
- data
A data frame containing the variables in the model.
- covariance
Covariance structure: "diagonal" (default), "diagonal_equal", "full", "full_equal", or "full_equal_corr".
- prior
Prior configuration: `prior_flat()`, `prior_normal()`, `prior_weak()`, `prior_rhs()`, or `prior_ssp()`. Default is `prior_flat()`. If `y_range` is supplied with the default flat prior, the wrapper automatically switches to `prior_weak()`.
- y_range
Optional theoretical range of the response variables. Use a numeric vector `c(min, max)` for a common range, or a P x 2 matrix/list for response-specific ranges. Specifying this automatically enables `prior_weak()` when `prior` is `prior_flat()` and calibrates the component mean and residual standard-deviation priors.
- fixed
Optional named list of fixed values for specific parameters.
- WAIC
Logical; if TRUE, add pointwise `log_lik` to the generate block for WAIC.
- ...
Reserved; unused arguments are rejected.
Examples
# Simulate 1D mixture data (2 components)
set.seed(123)
N <- 100
group <- rbinom(N, 1, 0.4)
y <- ifelse(group == 1, rnorm(N, 5, 1), rnorm(N, 0, 1))
dat <- data.frame(y)
# Fit a 1D Gaussian mixture model with 2 components
fit_mix <- rtmb_mixture(y ~ 1, k = 2, data = dat)
#> Pre-checking model code...
#> Checking RTMB setup...
# MAP estimation
map_mix <- fit_mix$optimize()
#> Starting RTMB optimization...
map_mix$summary()
#>
#> Call:
#> MAP Estimation via RTMB
#>
#> Negative Log-Posterior: 200.57
#> Approx. Log Marginal Likelihood (Laplace): -207.18
#>
#> Point Estimates and 95% Wald CI:
#> variable Estimate Std. Error Lower 95% Upper 95%
#> prob_mean[1] 0.41451 0.05103 0.31450 0.51452
#> prob_mean[2] 0.58549 0.05103 0.48548 0.68550
#> mu[C1] 4.81787 0.17776 4.46946 5.16628
#> mu[C2] 0.04764 0.12571 -0.19874 0.29402
#> sigma[C1] 1.00445 0.14877 0.75137 1.34278
#> sigma[C2] 0.88706 0.09688 0.71612 1.09880
#> theta[1] 0.41451 0.05103 0.31920 0.51668
#> theta[2] 0.58549 0.05103 0.48332 0.68080
#>