Introduction to poco

Source: vignettes/introduction.Rmd

poco (POsterior COmpression, “little” in Italian) provides a small, unified API for turning a posterior draws matrix into a compact density approximation, then drawing fresh samples or evaluating the density on demand.

This vignette shows the matrix-in / samples-out workflow on a synthetic two-component example. See vignette("compress-brms"), vignette("compress-funnel"), and vignette("compress-sccomp") for end-to-end workflows on real brms, cmdstanr, and sccomp fits, and vignette("methods-benchmark") for a side-by-side comparison of the available compression backends.

A toy posterior

We build a synthetic two-blob “posterior” so we can see compression preserve multimodality.

set.seed(1)
n_per <- 2000
draws <- rbind(
  mvtnorm::rmvnorm(n_per, mean = c(-2, -1),
                   sigma = matrix(c(0.5, 0.1, 0.1, 0.3), 2)),
  mvtnorm::rmvnorm(n_per, mean = c( 2,  1),
                   sigma = matrix(c(0.4, -0.2, -0.2, 0.6), 2))
)
colnames(draws) <- c("alpha", "beta")
head(draws)
#>          alpha       beta
#> [1,] -2.425341 -0.9508481
#> [2,] -2.458845 -0.2028409
#> [3,] -1.834450 -1.4180444
#> [4,] -1.598233 -0.5608073
#> [5,] -1.620042 -1.1191850
#> [6,] -0.906606 -0.6672904

Compress

Pick a method (the rest of the API does not change):

comp <- compress_posterior(
  draws,
  method = "mclust",
  n_components = 2
)
#> ℹ mclust: trying all 14 covariance models c(EII, VII, EEI, VEI, EVI, VVI, EEE, VEE, EVE, VVE, EEV, VEV, EVV, VVV) and picking the best by BIC (n = 4000, d = 2).
#> ℹ mclust: selected model 'VEV' with G = 2 (BIC = -21,469.2) out of 14 candidate models: EII, VII, EEI, VEI, EVI, VVI, EEE, VEE, EVE, VVE, EEV, VEV, EVV, VVV.
comp
#> <posterior_compressed: mclust >
#>  parameters: 2
#>  components: 2
#>  original draws: 4000
#>  mclust model: VEV
#>  BIC: -21469.2

compression_methods() lists the currently supported choices.

compression_methods()
#> [1] "mclust"     "mvdens_gmm" "mvdens_kde"

Regenerate samples 

new_draws <- sample_posterior(comp, n_draws = 5000)
dim(new_draws)
#> [1] 5000    2
head(new_draws)
#>          alpha       beta
#> [1,]  1.130526  0.6315805
#> [2,] -1.268359 -1.1603848
#> [3,] -1.974560 -1.0396438
#> [4,]  1.264395  0.7683689
#> [5,] -1.866382 -0.7604268
#> [6,]  0.698642  2.8331890

A quick sanity check on means and covariances:

round(colMeans(draws),     3); round(colMeans(new_draws),     3)
#>  alpha   beta 
#> -0.018  0.014
#>  alpha   beta 
#> -0.049  0.000
round(cov(draws),          3); round(cov(new_draws),          3)
#>       alpha  beta
#> alpha 4.398 1.971
#> beta  1.971 1.481
#>       alpha  beta
#> alpha 4.422 2.014
#> beta  2.014 1.487

Evaluate the density

density_posterior() accepts a matrix of points (one per row).

test_pts <- rbind(
  c(-2, -1),
  c( 0,  0),
  c( 2,  1)
)
density_posterior(comp, test_pts)
#> [1] 0.19803763 0.00326898 0.17712720
density_posterior(comp, test_pts, log = TRUE)
#> [1] -1.619298 -5.723277 -1.730887

Save / load

The compressed object is a plain list, so saveRDS() (with xz compression) works out of the box.

tf <- tempfile(fileext = ".rds")
saveRDS(comp, tf, compress = "xz")
file.size(tf)
#> [1] 452

comp_loaded <- readRDS(tf)
identical(class(comp_loaded), class(comp))
#> [1] TRUE

s <- sample_posterior(tf, n_draws = 100)  # path also works
dim(s)
#> [1] 100   2

unlink(tf)

Where to next

Session info 

sessionInfo()
#> R version 4.6.0 (2026-04-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.4 LTS
#> 
#> Matrix products: default
#> BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so;  LAPACK version 3.12.0
#> 
#> locale:
#>  [1] LC_CTYPE=C.UTF-8       LC_NUMERIC=C           LC_TIME=C.UTF-8       
#>  [4] LC_COLLATE=C.UTF-8     LC_MONETARY=C.UTF-8    LC_MESSAGES=C.UTF-8   
#>  [7] LC_PAPER=C.UTF-8       LC_NAME=C              LC_ADDRESS=C          
#> [10] LC_TELEPHONE=C         LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C   
#> 
#> time zone: UTC
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] mvtnorm_1.3-7 poco_0.2.0   
#> 
#> loaded via a namespace (and not attached):
#>  [1] instantiate_0.2.3 cli_3.6.6         knitr_1.51        rlang_1.2.0      
#>  [5] xfun_0.57         otel_0.2.0        processx_3.9.0    mclust_6.1.2     
#>  [9] textshaping_1.0.5 jsonlite_2.0.0    htmltools_0.5.9   ps_1.9.3         
#> [13] ragg_1.5.2        sass_0.4.10       rmarkdown_2.31    evaluate_1.0.5   
#> [17] jquerylib_0.1.4   fastmap_1.2.0     yaml_2.3.12       lifecycle_1.0.5  
#> [21] compiler_4.6.0    fs_2.1.0          htmlwidgets_1.6.4 systemfonts_1.3.2
#> [25] digest_0.6.39     R6_2.6.1          callr_3.7.6       bslib_0.10.0     
#> [29] tools_4.6.0       pkgdown_2.2.0     cachem_1.1.0      desc_1.4.3