---
title: "scmix: Bayesian sparse conditional mixture clustering"
author: "Aqi Dong, Yang-Li Liao, and Danhyang Lee"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{scmix: Bayesian sparse conditional mixture clustering}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>", fig.width = 7,
                      fig.height = 3)
set.seed(1)
```

`scmix` implements the Bayesian sparse conditional mixture model of
Dong, Liao, and Lee: each mixture component factorizes into a chain of univariate
polynomial regressions, per-component/per-equation variable selection is
*sampled* (spike-and-slab under a centered Zellner $g$-prior), and an
overfitted sparse mixture selects the number of clusters -- all in one
run of a blocked Gibbs sampler. This vignette reproduces, at reduced
scale, the two headline analyses of the accompanying paper.

## 1. Curved clusters that defeat Gaussian mixtures

Two opposite parabolas: a Gaussian mixture needs many elliptical
components to trace curvature, so BIC over-selects $K$; the conditional
mixture models the curvature directly.

```{r curved}
library(scmix)
n1 <- 400; n2 <- 300
x1 <- rnorm(n1 + n2)
x2 <- c( 0.5 + 1.0 * x1[1:n1]^2,
        -2.0 - 0.8 * x1[(n1 + 1):(n1 + n2)]^2) + rnorm(n1 + n2, 0, 0.6)
Y <- cbind(x1, x2)
truth <- rep(1:2, c(n1, n2))

fit <- scmix(Y, K = 7, m = 2, n_iter = 800, burn = 300, seed = 42)
fit
table(consensus = fit$cluster, truth)
```

A single run with an overfitted cap of $K = 7$ recovers the two curved
clusters (the chain lengths here are shortened for vignette build time;
the paper uses 1500 sweeps).

```{r plot, fig.alt = "Scatter plot of the two curved clusters colored by the fitted consensus partition"}
plot(fit, Y)
```

What did it learn? Per-cluster regression structure with uncertainty --
the interpretable output that distinguishes sampled selection from
shrinkage:

```{r summary}
summary(fit)
```

## 2. A replication-setting analysis

The paper's simulation study replicates Melnykov and Wang (2023) exactly.
Here is one draw from their $p = 3$ supplementary Table S-1 setting
(three components, quadratic conditionals), analyzed the same way:

```{r mw}
sim <- scmix_sim_mw3(n = 600, seed = 7)   # bundled S-1 generator
fit2 <- scmix(sim$Y, K = 7, m = 2, n_iter = 800, burn = 300, seed = 11)
c(K_consensus = fit2$K, ARI_vs_truth = round(scmix_ari(fit2$cluster, sim$z), 3))
```

The consensus partition and $K^*$ come from the same single run; no
sweep over $K$, no search over conditioning orders. The paper's Section
4 reports the full 100-replicate version of this analysis (consensus ARI
0.967, $K^* = 3$ in 100% of replicates), its order-sensitivity study,
and the comparisons with `cmbClust`, `mclust`, and the Bayesian
cluster-weighted models of Papastamoulis and Perrakis (2026).

## Options worth knowing

- `g = "hyper"`: hyper-$g$ prior (Liang et al., 2008); removes the
  $\sigma^2$ inflation of the default unit-information prior at high
  signal-to-noise, with clustering conclusions unchanged.
- `gprior_intercept = FALSE` (default in `scmix()`): centered
  convention; location-invariant selection in the linear case.
- `alive_frac` in the fit object: per-label occupancy diagnostics for
  component-specific summaries.
```