Package {normalblockr}

Type: Package
Title: Gaussian Graphical Models with Latent Clustering Structure
Version: 0.3.0
Description: Implements the Normal-Block model, a Gaussian graphical model with a latent clustering structure for the multivariate analysis of continuous data. The model clusters variables and, building on the graphical lasso, infers a network of statistical dependencies between clusters rather than between individual variables, for known or unknown clusterings, with an optional zero-inflation extension for data with an excess of exact zeros. A complementary family clusters variables by their regression response to covariates rather than by their covariance, sharing one profile per cluster. See Tous & Chiquet (2026) <doi:10.1016/j.csda.2026.108347> for the model itself and its variational expectation-maximization estimation procedure.
License: GPL (≥ 3)
URL: https://github.com/jchiquet/normalblockr, https://jchiquet.github.io/normalblockr/
BugReports: https://github.com/jchiquet/normalblockr/issues
Encoding: UTF-8
LazyData: true
Depends: R (≥ 4.1.0)
Imports: Matrix, dplyr, tidyr, tibble, purrr, R6, ggplot2, igraph, sbm, corrplot, scales, MASS, stats, Rcpp
LinkingTo: Rcpp, RcppArmadillo
Suggests: testthat (≥ 3.0.0), aricode, covr, Metrics, rmarkdown, pheatmap, paletteer, ggthemes, knitr
Config/testthat/edition: 3
VignetteBuilder: knitr
Config/roxygen2/version: 8.1.0
NeedsCompilation: yes
Packaged: 2026-09-11 17:43:32 UTC; jchiquet
Author: Jeanne Tous [aut], Nestor Ngalala Manguitini [ctb], Julien Chiquet ORCID iD [aut, cre]
Maintainer: Julien Chiquet <julien.chiquet@inrae.fr>
Repository: CRAN
Date/Publication: 2026-09-11 18:10:02 UTC

normalblockr: Gaussian Graphical Models with Latent Clustering Structure

Description

Implements the Normal-Block model, a Gaussian graphical model with a latent clustering structure for the multivariate analysis of continuous data. The model clusters variables and, building on the graphical lasso, infers a network of statistical dependencies between clusters rather than between individual variables, for known or unknown clusterings, with an optional zero-inflation extension for data with an excess of exact zeros. A complementary family clusters variables by their regression response to covariates rather than by their covariance, sharing one profile per cluster. See Tous & Chiquet (2026) doi:10.1016/j.csda.2026.108347 for the model itself and its variational expectation-maximization estimation procedure.

Author(s)

Maintainer: Julien Chiquet julien.chiquet@inrae.fr (ORCID)

Authors:

Other contributors:

See Also

Useful links:

Bayesian Information Criterion for a Normal-Block Model

Description

Extracts the (variational) BIC of a fitted normal-block model, computed as 'deviance + log(n) * nb_param' (lower is better).

Usage

## S3 method for class 'NormalBlockBase'
BIC(object, ...)

Arguments

object

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

A scalar: the (variational) BIC.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
model <- normal_block(data, blocks = 3, control = NB_control(verbose = FALSE))
BIC(model)

Bayesian Information Criterion for a Collection of Normal-Block Models

Description

Returns the (variational) BIC of every model in a collection of normal-block models.

Usage

## S3 method for class 'NormalBlockCollection'
BIC(object, ...)

Arguments

object

An object inheriting from NormalBlockCollection.

...

not used, only here for S3 compatibility

Value

A numeric vector of BIC values, one per model in the collection.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
models <- normal_block(data, blocks = 2:5, control = NB_control(verbose = FALSE))
BIC(models)

NB_control

Description

Control the model settings and various optimization parameters

Usage

NB_control(
  niter = 500,
  threshold = 1e-04,
  fixed_point_niter = 5,
  sparsity_weights = NULL,
  sparsity_penalties = NULL,
  n_sparsity_penalties = 30,
  min_ratio = 0.01,
  fixed_tau = FALSE,
  clustering_init = NULL,
  verbose = TRUE,
  heuristic = FALSE,
  noise_covariance = NULL,
  refine = FALSE
)

Arguments

niter

number of iterations in model optimization

threshold

loglikelihood / elbo threshold under which optimization stops

fixed_point_niter

number of sweeps of the tau update for Normal-Block-Mean with unknown clusters. Each sweep visits the rows of tau sequentially and maximizes the ELBO exactly in each, so it can never decrease the ELBO.

sparsity_weights

weights with which the penalty should be applied in case sparsity is required, non-0 values on the diagonal mean diagonal shall be penalized too (default is non-penalized diagonal and 1s off-diagonal)

sparsity_penalties

list of penalties the user wants to test, other parameters are only used if penalties is not specified

n_sparsity_penalties

number of penalties to test.

min_ratio

ratio for sparsity between max penalty (0 edge penalty) and min penalty to test

fixed_tau

whether tau should be fixed at clustering_init during optimization useful for calls to fixed_q models in stability_selection

clustering_init

how to obtain the initial clustering of the q unknown blocks: a heuristic name ("ward2", "kmeans", "sbm" or "spectral"), an actual clustering (a vector of labels or a p x q indicator matrix, or a list of either per q for a collection), or "best_of_inits" to try several heuristics per model and keep the best-ELBO fit (see [NormalBlockVarBase]'s 'best_of_inits()'; not supported with 'sparsity = TRUE'). Default 'NULL', resolved per model family at fit time: "ward2" for variance-block models, "kmeans" for mean-block models ("ward2" was benchmarked substantially worse there – see [NormalBlockMeanBase]). See 'inst/methods_initialization_and_refine.md' for the heuristics' rationale, why no single one dominates, and how this interacts with 'refine' (below).

verbose

telling if information should be printed during optimization

heuristic

whether to use the heuristic approach (moment-based, no (V)EM recursion) instead of the full (V)EM. Default is FALSE. In heuristic mode, no likelihood/ELBO is computed, so 'entropy', 'loglik', 'BIC', 'ICL' and 'EBIC' are all 'NA' on the resulting model.

noise_covariance

shape of the residual covariance. Variance-block models accept "diagonal" (variable-specific) or "spherical" (common); mean-block models, whose Sigma is the full p x p residual covariance, also accept "full". Default 'NULL', resolved per model family at fit time to "diagonal" – except for a mean-block model with 'sparsity > 0', which implies "full" since a penalty on a diagonal precision matrix would have nothing to act on (explicitly asking for both is an error). The mean-block "diagonal"/"spherical" variants need no matrix inversion, hence no 'n > p' requirement, and they select the number of clusters markedly better than a full Sigma once p approaches n: the p(p+1)/2 covariance parameters otherwise drown the mean structure BIC/ICL are weighing. Use "full" when the residual associations are themselves of interest.

refine

for [NormalBlockVarCollectionClusters] only: whether 'optimize()' should automatically call 'refine()' afterwards. Default 'FALSE' since it adds real cost; call 'collection$refine()' directly at any point afterwards for the same effect without setting this.

Value

A named list of parameters to pass to [normal_block()]'s 'control' argument.

Root Base Class for Normal-Block Models

Description

R6 abstract class shared by the variance-block ([NormalBlockVarBase]) and mean-block ([NormalBlockMeanBase]) model families.

Public fields

data

object of NormalBlockData class, with responses and design matrix

Active bindings

inference_method

inference procedure used (heuristic or integrated with EM)

n

number of samples

p

number of responses per sample

d

number of variables (dimensions in X)

q

number of blocks

n_edges

number of edges of the network (non null coefficient of the sparse precision matrix Omega)

objective

evolution of the objective function during (V)EM algorithm

loglik

(or its variational lower bound)

deviance

(or its variational lower bound)

BIC

(or its variational lower bound)

entropy

Entropy of the conditional distribution when applicable

ICL

variational lower bound of the ICL

EBIC

variational lower bound of the EBIC

criteria

a vector with loglik, BIC and number of parameters

sparsity

(overall sparsity parameter)

sparsity_term

(sparsity_term term in log-likelihood due to sparsity)

memberships

cluster memberships

clustering

given as the list of elements contained in each cluster

cluster_sizes

given as a vector of cluster sizes

elements_per_cluster

given as the list of elements contained in each cluster

Methods

Public methods

NormalBlockBase$new()

Create a new ['NormalBlockBase'] object.

Usage
NormalBlockBase$new(data, q, sparsity = 0, control, zero_inflation = FALSE)
Arguments
data

object of NormalBlockData class, with responses and design matrix

q

number of block/cluster

sparsity

sparsity penalty on the network density

control

structured list of more specific parameters, to generate with NB_control

zero_inflation

whether the concrete subclass models zero-inflation; set by the ZI subclasses themselves, not meant to be set by the end user. When 'FALSE', the (costly) zero-inflation probability fit ('kappa'/'B0') is skipped entirely, since it would otherwise never be used downstream.

Returns

A new ['NormalBlockBase'] object

NormalBlockBase$update()

Update a ['NormalBlockBase'] object

All possible parameters of the child classes

Usage
NormalBlockBase$update(
  B = NA,
  dm1 = NA,
  C = NA,
  Omega = NA,
  gamma = NA,
  mu = NA,
  kappa = NA,
  alpha = NA,
  M = NA,
  S = NA,
  Psi = NA,
  Phi = NA,
  Lambda = NA,
  ll_list = NA,
  warm_started = NA,
  clustering_init = NA
)
Arguments
B

regression matrix [all]

dm1

diagonal vector of inverse variance matrix (variables level) [NBVar]

C

the matrix of groups memberships (posterior probabilities) [all]

Omega

inverse variance matrix (cluster-level for Normal Block models, variable-level for Normal Mean Block models) [all]

gamma

variance of posterior distribution of W [NBVar - known]

mu

mean for posterior distribution of W [NBVar - known]

kappa

vector of zero-inflation probabilities [ZINBVar, ZINBMean]

alpha

vector of groups probabilities [NBVar]

M

variational mean for posterior distribution of W [NBVar - unknown]

S

variational diagonal of variances for posterior distribution of W [NBVar - unknown]

Psi

variational expectation of C'Omega C, intermediary term in calculations [NBMean - unknown]

Phi

variational correction term used in the Psi/ELBO computations [NBMean - unknown]

Lambda

variational correction term used in the Sigma-hat update [NBMean - unknown]

ll_list

list of log-lik (elbo) values

warm_started

whether 'optim_initialize()' should treat the model as already initialized (reuse B/Omega/dm1/C/alpha/M/S as they stand) rather than recomputing a fresh heuristic initialization. Set by [warm_start_from()] and by [split()]/[merge()].

clustering_init

initial clustering

Returns

Update the current ['normal'] object

NormalBlockBase$optimize()

calls optimization (EM or heuristic) and updates relevant fields

Usage
NormalBlockBase$optimize(
  control = list(niter = 500, threshold = 1e-04, fixed_point_niter = 5),
  warn = TRUE
)
Arguments
control

a list for controlling the optimization process

warn

whether to warn when the (V)EM stops at the 'niter' cap without reaching 'threshold' (see 'private$warn_if_not_converged()'). Set to 'FALSE' for deliberately-truncated trial fits (cheap candidate scoring in 'candidates_split()'/'candidates_merge()', the sparsity-path warm-start probe in [NormalBlockCollectionSparsity]) where stopping at the cap is expected and not a sign of trouble.

Returns

optimizes the model and updates its parameters

NormalBlockBase$best_of_inits()

Try several clustering-initialization heuristics and keep the best-ELBO converged fit (see 'NB_control(clustering_init = )' and 'inst/methods_initialization_and_refine.md' for the rationale). Every candidate is first screened with a short 'trial_niter' run (same idea as 'candidates_split()'/'candidates_merge()'), and only the 'max_training' best-screened ones are fully retrained with 'control'.

Usage
NormalBlockBase$best_of_inits(
  inits = private$default_inits,
  trial_niter = 10,
  max_training = 2,
  control = list(niter = 500, threshold = 1e-04, fixed_point_niter = 5)
)
Arguments
inits

vector of clustering-heuristic names to try; defaults to the model family's own preferred order ('private$default_inits')

trial_niter

number of (V)EM iterations used to cheaply screen every candidate in 'inits' before fully retraining the best few

max_training

how many of the screened candidates (best 'loglik' after 'trial_niter' iterations) get fully retrained with 'control'

control

'optimize()' control list ('niter'/'threshold') used for the final full retraining of the 'max_training' best candidates

Returns

a new, already-optimized model. Does not mutate the current object; reassign the result ('model <- model$best_of_inits()').

NormalBlockBase$candidates_split()

generate and select a set of candidate models by splitting the clusters of the current model

Usage
NormalBlockBase$candidates_split(trial_niter = 5)
Arguments
trial_niter

number of (V)EM iterations used to cheaply score each candidate before the caller fully re-optimizes the best few – kept short on purpose.

NormalBlockBase$candidates_merge()

generate and select a set of candidate models by merging the clusters of the current model

Usage
NormalBlockBase$candidates_merge(max_candidates = 30, trial_niter = 2)
Arguments
max_candidates

merge candidates are, unlike split's, quadratic in q ('choose(q, q-2)' pairs): beyond 'max_candidates' pairs, only the most promising ones are actually built and trial-optimized, ranked by the family's own 'private$merge_score()'. Set to 'Inf' to always try every pair.

trial_niter

see [candidates_split()]

NormalBlockBase$predict()

Predicts observations Y for new covariates X, in Y's original units (like '$fitted', so that predicting on the training X reproduces it).

Usage
NormalBlockBase$predict(new_X)
Arguments
new_X

new set of covariates.

Returns

A n*p prediction matrix for new observations

NormalBlockBase$latent_network()

Extract interaction network in the latent space, as a matrix rather than a plot – see '$plot_network()' to plot it instead.

Usage
NormalBlockBase$latent_network(type = c("partial_cor", "support", "precision"))
Arguments
type

edge value in the network. Can be "support" (binary edges), "precision" (coefficient of the precision matrix) or "partial_cor" (partial correlation between species)

Returns

a square matrix of size 'self$q'

NormalBlockBase$plot_loglik()

plots the evolution of the objective (log-likelihood or ELBO) across the (V)EM iterations of the last call to 'optimize()'.

Usage
NormalBlockBase$plot_loglik(show_increment = TRUE)
Arguments
show_increment

whether to add a second panel with the (log10) absolute increment between iterations and the convergence 'threshold' (distinguishes true convergence from a flat-looking objective trace).

Returns

a ['ggplot2::ggplot'] graph

NormalBlockBase$plot_network()

plot the latent network. To extract the network as a matrix instead of plotting it, use '$latent_network()'.

Usage
NormalBlockBase$plot_network(
  type = c("partial_cor", "support"),
  output = c("igraph", "corrplot"),
  edge.color = c("#F8766D", "#00BFC4"),
  remove.isolated = FALSE,
  node.labels = NULL,
  layout = igraph::layout_in_circle,
  plot = TRUE
)
Arguments
type

edge value in the network. Either "precision" (coefficient of the precision matrix) or "partial_cor" (partial correlation between species).

output

Output type. Either 'igraph' (for the network) or 'corrplot' (for the adjacency matrix)

edge.color

Length 2 color vector. Color for positive/negative edges. Default is 'c("#F8766D", "#00BFC4")'. Only relevant for igraph output.

remove.isolated

if 'TRUE', isolated node are remove before plotting. Only relevant for igraph output.

node.labels

vector of character. The labels of the nodes. The default will use the column names ot the response matrix.

layout

an optional igraph layout. Only relevant for igraph output.

plot

logical. Should the final network be displayed or only sent back to the user. Default is 'TRUE'.

NormalBlockBase$plot()

plots the evolution of the objective during model optimization (see 'plot_loglik()')

Usage
NormalBlockBase$plot()

NormalBlockBase$print()

User friendly print method

Usage
NormalBlockBase$print(model = paste("A", self$who_am_I, ".\n"))
Arguments
model

First line of the print output

NormalBlockBase$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockBase$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly -- see
# normal_block() for how concrete models are created and fitted.

Base Class for a Collection of Normal-Block Models

Description

Shared scaffolding for the collections explored by [get_model()]/ [normal_block()]: a sweep over sparsity penalties (['NormalBlockVarCollectionSparsity']), over the number of clusters (['NormalBlockVarCollectionClusters']), or over both (['NormalBlockVarCollectionClustersSparsity']). Concrete subclasses set 'private$progress_field'/'private$progress_label' in their 'initialize()' and provide their own 'get_best_model()', delegating the (row of 'self$criteria' minimizing a criterion) lookup to 'private$best_id()'.

Public fields

models

list of models (or sub-collections) explored by the collection

control

store the list of user-defined model settings and optimization parameters

Active bindings

criteria

a data frame with the values of some criteria for the collection of models

loglik

not defined for a collection: accessing it raises an informative error. Use ‘logLik()' for every model’s log-likelihood, or '$get_best_model()$loglik' for a single one.

Methods

Public methods

NormalBlockCollection$optimize()

optimizes every model (or sub-collection) in the collection

Usage
NormalBlockCollection$optimize(
  control = list(niter = 500, threshold = 1e-04, verbose = TRUE)
)
Arguments
control

optimization parameters (niter and threshold). When 'control$clustering_init' is '"best_of_inits"', each leaf model is fit via its own 'best_of_inits()' instead of a plain 'optimize()'.

NormalBlockCollection$print()

User-friendly print method: model type and the range of q/sparsity explored. See 'summary()' for the full criteria table.

Usage
NormalBlockCollection$print()

NormalBlockCollection$summary()

Summarize the collection: model type, full criteria table, and the range of q/sparsity explored.

Usage
NormalBlockCollection$summary()
Returns

An object of class 'summary.NormalBlockCollection', printed with a dedicated [print.summary.NormalBlockCollection()] method.

NormalBlockCollection$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockCollection$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly -- see
# normal_block() for how collections (NormalBlockVarCollectionClusters,
# NormalBlockVarCollectionSparsity, NormalBlockVarCollectionClustersSparsity)
# are actually created and fitted.

Base Class for a Collection of Models over a Range of Cluster Counts

Description

Shared scaffolding for [NormalBlockVarCollectionClusters] and [NormalBlockMeanCollectionClusters]: everything that does not depend on the model family (model lookup, model selection, the criteria plot and the split/merge 'refine()' search). Concrete subclasses only build 'self$models' in their 'initialize()' and name themselves through 'who_am_I'.

Super class

NormalBlockCollection -> NormalBlockCollectionClusters

Active bindings

q_list

number of blocks

Methods

Public methods

Inherited methods

NormalBlockCollectionClusters$get_model()

returns the unknown-clusters model corresponding to given q

Usage
NormalBlockCollectionClusters$get_model(q)
Arguments
q

number of blocks asked by user

Returns

A unknown-clusters object with given value q

NormalBlockCollectionClusters$get_best_model()

Extract best model in the collection

Usage
NormalBlockCollectionClusters$get_best_model(
  crit = c("ICL", "BIC", "EBIC", "deviance")
)
Arguments
crit

a character for the criterion used to performed the selection. Either "ICL" or "BIC". "ICL" is the default criterion

Returns

a ['unknown-clusters'] object

NormalBlockCollectionClusters$plot()

Display various outputs (goodness-of-fit criteria, robustness, diagnostic) associated with a collection of network fits (a ['Networkfamily'])

Usage
NormalBlockCollectionClusters$plot(
  criteria = c("deviance", "ICL", "BIC", "EBIC")
)
Arguments
criteria

vector of characters. The criteria to plot in 'c("deviance", "BIC", "ICL")'. Defaults to all of them.

Returns

a ['ggplot2::ggplot'] graph

NormalBlockCollectionClusters$optimize()

optimizes every model in the collection, then – if 'control$refine' is 'TRUE' (see [NB_control()], default 'FALSE') – calls [refine()] automatically.

Usage
NormalBlockCollectionClusters$optimize(
  control = list(niter = 500, threshold = 1e-04, verbose = TRUE)
)
Arguments
control

optimization parameters (niter, threshold, verbose)

NormalBlockCollectionClusters$refine()

Tries to improve every model in the collection with a short split-and-reoptimize trial seeded from its smaller-q neighbor ('"split"'), a short merge-and-reoptimize trial seeded from its larger-q neighbor ('"merge"'), or both (the default); a candidate replaces the original only if it strictly lowers the deviance, so this can only improve (or leave unchanged) each model it touches. Only contiguous q pairs ('q' and 'q -/+ 1', both present in the collection) are refined. See 'inst/methods_initialization_and_refine.md' for the rationale and empirical evidence.

Usage
NormalBlockCollectionClusters$refine(
  trial_niter = 2,
  max_candidates = 30,
  directions = c("split", "merge"),
  verbose = self$control$verbose
)
Arguments
trial_niter

number of EM iterations used for the cheap trial candidates (passed to 'candidates_split()'/'candidates_merge()') before fully re-optimizing only the best one.

max_candidates

passed to 'candidates_merge()' (ignored for '"split"', which is never combinatorial in q) – see its documentation.

directions

which neighbor(s) to seed refinement candidates from: '"split"' (smaller-q neighbor), '"merge"' (larger-q neighbor), or both (the default).

verbose

whether to print, for each q attempted, whether the candidate from that neighbor improved on it. Defaults to 'control$verbose' (the value set at construction, see [NB_control()]).

Returns

invisibly returns 'self'; improved models replace the originals in '$models' in place.

NormalBlockCollectionClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockCollectionClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly -- see
# normal_block() for how collections are created and fitted.

Base Class for a Collection over Cluster Counts and Sparsity Levels

Description

Shared scaffolding for [NormalBlockVarCollectionClustersSparsity] and [NormalBlockMeanCollectionClustersSparsity]: a collection of sparsity sub-collections, one per q. Everything family-agnostic (two-key model lookup, model selection over both axes, the criteria heatmap) lives here; subclasses only build 'self$models' and name themselves.

Super class

NormalBlockCollection -> NormalBlockCollectionClustersSparsity

Active bindings

q_list

number of blocks

sparsity

list of penalties used for each q

Methods

Public methods

Inherited methods

NormalBlockCollectionClustersSparsity$get_model()

returns a collection of models corresponding to given q or one single model if penalty is also given

Usage
NormalBlockCollectionClustersSparsity$get_model(q, sparsity = NA)
Arguments
q

number of blocks asked by user.

sparsity

sparsity penalty asked by user

Returns

either a sparsity sub-collection or a single model object

NormalBlockCollectionClustersSparsity$get_best_model()

Extract best model in the collection

Usage
NormalBlockCollectionClustersSparsity$get_best_model(
  crit = c("ICL", "BIC", "EBIC")
)
Arguments
crit

a character for the criterion used to performed the selection. Either "BIC", "EBIC" or "ICL". "ICL" is the default criterion

Returns

a single fitted model object

NormalBlockCollectionClustersSparsity$plot()

Display various outputs (goodness-of-fit criteria, robustness, diagnostic) associated with a collection of network fits

Usage
NormalBlockCollectionClustersSparsity$plot(
  criterion = c("deviance", "ICL", "BIC", "EBIC"),
  n_intervals = NULL
)
Arguments
criterion

The criteria to plot in 'c("deviance", BIC", "EBIC", "ICL")'. Defaults deviance.

n_intervals

number of intervals into which the penalties range should be split

Returns

a ['ggplot2::ggplot'] heatmap

NormalBlockCollectionClustersSparsity$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockCollectionClustersSparsity$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly -- see
# normal_block() for how collections are created and fitted.

Base Class for a Collection of Models over a Sparsity Path

Description

Shared scaffolding for [NormalBlockVarCollectionSparsity] and [NormalBlockMeanCollectionSparsity]: the warm-started path traversal, the penalty lookup and the criteria plot. Concrete subclasses derive the penalty grid in their 'initialize()' and provide their own 'get_best_model()' (only the variance-block family offers StARS).

Super class

NormalBlockCollection -> NormalBlockCollectionSparsity

Public fields

data

object of NormalBlockData class, with responses and design matrix

Active bindings

q

number of blocks

blocks

group matrix or number of blocks

sparsity

list of sparsity penalties

Methods

Public methods

Inherited methods

NormalBlockCollectionSparsity$optimize()

optimizes every model in the sparsity path, warm-starting each one (after the first) from the previous, adjacent penalty's converged parameters (see the family's base class's 'warm_start_from()') instead of re-deriving everything from the heuristic clustering, the way the generic [NormalBlockCollection] 'optimize()' would. 'blocks' (hence q) is fixed across the whole path, only the sparsity penalty changes, so the warm start is always between models of matching shape.

Usage
NormalBlockCollectionSparsity$optimize(
  control = list(niter = 500, threshold = 1e-04, verbose = TRUE)
)
Arguments
control

optimization parameters (niter and threshold)

NormalBlockCollectionSparsity$get_model()

returns the NormalBlockVarKnownClusters model corresponding to given penalty

Usage
NormalBlockCollectionSparsity$get_model(sparsity)
Arguments
sparsity

sparsity penalty asked by user

Returns

A NormalBlockVarKnownClusters (sparse) object with given value penalty

NormalBlockCollectionSparsity$plot()

Display various outputs (goodness-of-fit criteria, robustness, diagnostic) associated with a collection of network fits (a collection)

Usage
NormalBlockCollectionSparsity$plot(
  criteria = c("deviance", "BIC", "EBIC", "ICL"),
  log.x = TRUE
)
Arguments
criteria

vector of characters. The criteria to plot in 'c("deviance", BIC", "EBIC", "ICL")'. Defaults to all of them.

log.x

logical: should the x-axis be represented in log-scale? Default is 'TRUE'.

Returns

a ['ggplot2::ggplot'] graph

NormalBlockCollectionSparsity$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockCollectionSparsity$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly -- see
# normal_block() for how collections are created and fitted.

Data Container for Normal-Block Models

Description

R6 class holding the responses and design matrix used to fit a normal-block model.

Public fields

Y

the matrix of responses (rescaled column-wise if 'scale = TRUE')

Y_scale

the per-column standard deviation Y was divided by (all 1's if 'scale = FALSE')

X

the matrix of covariates

X0

the matrix of zero-inflation covariates, if applicable

formula

describes the relationship between Y and X, and X0 if applicable, useful if not all of X's or X0's covariates should be used, should be formatted ~ X1 + X2... | Z1 + Z2... with the Normal formula before the | and the ZI formula after the |

n

sample size

d

number of covariates

d0

number of zero-inflation covariates, if applicable

p

number of variables

XtX

useful for inference in some cases

XtXm1

inverse of XtX, useful for inference

XtY

useful for inference

npY

total number of non zeros in Y

nY

total number of non zeros for each column/variable in Y

zeros

where are the zero in Y

zeros_bar

where are the non-zeros in Y

Methods

Public methods

NormalBlockData$new()

Create a new ['NormalBlockData'] object.

Usage
NormalBlockData$new(
  Y,
  X,
  X0 = NULL,
  formula = NULL,
  scale = TRUE,
  zeros = NULL
)
Arguments
Y

the matrix of responses (called Y in the model).

X

design matrix (called X in the model).

X0

zero-inflation design matrix, if applicable.

formula

describes the relationship between Y and X, useful if not all of X's covariates should be used.

scale

whether to rescale each column of Y by its own standard deviation (no centering). Default TRUE, see the class-level documentation for the rationale and its limits.

zeros

an optional 0/1 matrix of structural zeros, overriding the default 'Y == 0'.

NormalBlockData$ols_fit()

Ordinary-least-squares fit of Y on X, with its residuals and their covariance. Computed once and memoized.

Usage
NormalBlockData$ols_fit()
Returns

a list with 'B' (d x p), 'R' (n x p residuals) and 'Sigma' (p x p residual covariance)

NormalBlockData$zi_ols_fit()

Masked counterpart of 'ols_fit()': a per-variable weighted least-squares fit of B under the zero-inflation mask, with its inverse residual variances and residuals (see 'zi_weighted_fit()'). Computed once and memoized.

Usage
NormalBlockData$zi_ols_fit()
Returns

a list with 'B' (d x p), 'dm1' (p) and 'R' (n x p masked residuals)

NormalBlockData$zi_fit()

Zero-inflation component: 'p' independent logistic regressions of each variable's zero pattern on 'X0', and the fixed contribution they make to the log-likelihood. The (V)EM never revisits these, so they are a property of the data rather than of a model.

Usage
NormalBlockData$zi_fit()
Returns

a list with 'B0' (d0 x p), 'kappa' (n x p zero-inflation probabilities) and 'ZI_cond_mean' (a scalar)

NormalBlockData$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockData$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
c(n = data$n, p = data$p, d = data$d)

Base Class for Mean-Block Models

Description

R6 abstract class for the Normal-Block models where the clustering structures the mean (mu_i = C B' X_i).

Super class

NormalBlockBase -> NormalBlockMeanBase

Active bindings

model_par

a list with the matrices of the model parameters: B (covariates), dm1 (species variance), Omega (groups precision matrix)). On the internal fitting scale (‘self$data$Y', possibly column-rescaled by 'NormalBlockData(scale = TRUE)'): use '$B_original'/'$dm1_original' for the same quantities converted back to Y’s original units.

nb_param

number of parameters in the model

sparsity_weights

(weights associated to each pair of groups)

Methods

Public methods

Inherited methods

NormalBlockMeanBase$new()

Create a new ['NormalBlockMeanBase'] object.

Usage
NormalBlockMeanBase$new(
  data,
  q,
  sparsity = 0,
  control = NB_control(),
  zero_inflation = FALSE
)
Arguments
data

object of NormalBlockData class, with responses and design matrix

q

number of block/cluster

sparsity

sparsity penalty on the network density

control

structured list of more specific parameters, to generate with NB_Mean_control

zero_inflation

whether the concrete subclass models zero-inflation; set by the ZI subclasses themselves, not meant to be set by the end user.

Returns

A new ['NormalBlockMeanBase'] object

NormalBlockMeanBase$predict()

Predicts observations Y for new covariates X, in Y's original units. The mean-block mean is mu_i = C B' X_i, so the cluster-level predictor has to be mapped back to the variables through C.

Usage
NormalBlockMeanBase$predict(new_X)
Arguments
new_X

new set of covariates.

Returns

A n*p prediction matrix for new observations

NormalBlockMeanBase$warm_start_from()

Seed this model's starting parameters from another, already-optimized model with the same q, instead of the heuristic clustering-derived values set at construction time. Used by [split()]/[merge()].

Usage
NormalBlockMeanBase$warm_start_from(other)
Arguments
other

a [NormalBlockMeanBase] object, already optimized

Returns

Update the current object in place with ‘other'’s parameters

NormalBlockMeanBase$split()

Create a clone of the current ['NormalBlockMeanBase'] object after splitting cluster 'index'. Unlike the variance-block family, Omega and the sparsity weights are p x p here and do not depend on q, so they carry over unchanged; only C (tau) and B (one column per cluster) are affected. Variables are split by their current noise variance (1 / diag(Omega)) around its within-cluster median, the same criterion [NormalBlockVarBase]'s 'split()' uses via 'dm1', since 'diag(Omega)' plays the same per-variable-precision role here.

Usage
NormalBlockMeanBase$split(index, in_place = FALSE)
Arguments
index

index (integer) of the cluster to split

in_place

should the split be applied to the object itself, or should a copy be sent? Default FALSE (send a copy)

Returns

A new ['NormalBlockMeanBase'] object

NormalBlockMeanBase$merge()

Create a clone of the current ['NormalBlockMeanBase'] object after merging clusters 'indices'

Usage
NormalBlockMeanBase$merge(indices, in_place = FALSE)
Arguments
indices

indices (couple of integer) of the clusters to merge

in_place

should the merge be applied to the object itself, or should a copy be sent? Default FALSE (send a copy)

Returns

A new ['NormalBlockMeanBase'] object

NormalBlockMeanBase$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockMeanBase$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly -- use
# NormalBlockMeanKnownClusters / NormalBlockMeanUnknownClusters.

Collection of Mean-Block Models over a Range of Cluster Counts

Description

R6 class for a collection of mean-block models ([NormalBlockMeanBase]) with different numbers of clusters (q). Inherits its scaffolding ('print()'/'summary()'/'plot()'/'optimize()', the 'criteria' table) from [NormalBlockCollection]. Unlike [NormalBlockVarCollectionClusters], there is no SBM-path shortcut here: the shared clustering-heuristic registry's cov()/correlation-based methods are ill-suited to the mean-block family (see [NormalBlockMeanBase]'s own default), so each q is fit independently.

Super classes

NormalBlockCollection -> NormalBlockCollectionClusters -> NormalBlockMeanCollectionClusters

Active bindings

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockMeanCollectionClusters$new()

Create a new ['NormalBlockMeanCollectionClusters'] object.

Usage
NormalBlockMeanCollectionClusters$new(
  mydata,
  q_list,
  zero_inflation = FALSE,
  sparsity = 0,
  control = NB_control()
)
Arguments
mydata

object of NormalBlockData class, with responses and design matrix

q_list

list of q values (number of groups) in the collection

zero_inflation

whether Y carries structural zeros; every model in the collection is then zero-inflated (Sigma diagonal or spherical only, see [NormalBlockMeanBase])

sparsity

sparsity penalty on the network density

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['NormalBlockMeanCollectionClusters'] object

NormalBlockMeanCollectionClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockMeanCollectionClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
models <- normal_block(data, blocks = 2:5, model = "mean")
models$plot(c("BIC", "ICL"))
models$get_best_model()

Collection of Mean-Block Models over Cluster Counts and Sparsity Levels

Description

R6 class for a collection of mean-block models ([NormalBlockMeanBase]) over both a range of cluster counts (q) and a sparsity path, i.e. one [NormalBlockMeanCollectionSparsity] per q. Mirrors [NormalBlockVarCollectionClustersSparsity], minus its SBM-path shortcut for the initial clustering (ill-suited to this family, see [NormalBlockMeanCollectionClusters]).

Super classes

NormalBlockCollection -> NormalBlockCollectionClustersSparsity -> NormalBlockMeanCollectionClustersSparsity

Active bindings

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockMeanCollectionClustersSparsity$new()

Create a new ['NormalBlockMeanCollectionClustersSparsity'] object.

Usage
NormalBlockMeanCollectionClustersSparsity$new(
  mydata,
  q_list,
  control = NB_control()
)
Arguments
mydata

object of NormalBlockData class, with responses and design matrix

q_list

list of q values (number of groups) in the collection

control

structured list of parameters to handle sparsity control

Returns

A new ['NormalBlockMeanCollectionClustersSparsity'] object

NormalBlockMeanCollectionClustersSparsity$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockMeanCollectionClustersSparsity$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 60, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
models <- normal_block(data, blocks = 2:4, sparsity = TRUE, model = "mean",
                       control = NB_control(n_sparsity_penalties = 4))
models$plot("BIC")
models$get_best_model("BIC")

Collection of Mean-Block Models over a Sparsity Path

Description

R6 class for a collection of mean-block models ([NormalBlockMeanBase]) with a fixed clustering (or a fixed number of blocks) and different sparsity levels applied to the p x p precision matrix of the variables. Mirrors [NormalBlockVarCollectionSparsity], minus the StARS/stability selection path, which relies on 'fixed_tau', not supported by the mean-block VEM.

Super classes

NormalBlockCollection -> NormalBlockCollectionSparsity -> NormalBlockMeanCollectionSparsity

Active bindings

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockMeanCollectionSparsity$new()

Create a new ['NormalBlockMeanCollectionSparsity'] object.

Usage
NormalBlockMeanCollectionSparsity$new(mydata, blocks, control = NB_control())
Arguments
mydata

object of NormalBlockData class, with responses and design matrix

blocks

either a clustering matrix (known, fixed clustering) or a single integer (number of blocks to infer)

control

structured list of parameters to handle sparsity control

Returns

A new ['NormalBlockMeanCollectionSparsity'] object

NormalBlockMeanCollectionSparsity$get_best_model()

Extract best model in the collection

Usage
NormalBlockMeanCollectionSparsity$get_best_model(
  crit = c("BIC", "EBIC", "ICL")
)
Arguments
crit

a character for the criterion used to perform the selection, either "BIC", "EBIC" or "ICL". Default is BIC

Returns

a ['NormalBlockMeanBase'] object

NormalBlockMeanCollectionSparsity$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockMeanCollectionSparsity$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 60, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
models <- normal_block(data, blocks = 3, sparsity = TRUE, model = "mean",
                       control = NB_control(n_sparsity_penalties = 5))
models$plot(c("BIC", "EBIC"))

Mean-Block Model with Known Clustering

Description

R6 class for a Normal-Block-Mean model with a known clustering.

Super classes

NormalBlockBase -> NormalBlockMeanBase -> NormalBlockMeanKnownClusters

Active bindings

fitted

Y values predicted by the model, in Y's original units

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockMeanKnownClusters$new()

Create a new ['NormalBlockMeanKnownClusters'] object.

Usage
NormalBlockMeanKnownClusters$new(data, C, sparsity = 0, control = NB_control())
Arguments
data

object of NormalBlockData class, with responses and design matrix

C

clustering matrix C_jk = 1 if species j belongs to cluster k

sparsity

to apply on variance matrix when calling GLASSO

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['NormalBlockMeanKnownClusters'] object

NormalBlockMeanKnownClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockMeanKnownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
model <- normal_block(data, blocks = ex$parameters$C, model = "mean")
model$plot()

Mean-Block Model with Unknown Clustering

Description

R6 class for a Normal-Block-Mean model with a fixed number of clusters (but unknown clustering), inferred by variational EM.

Super classes

NormalBlockBase -> NormalBlockMeanBase -> NormalBlockMeanUnknownClusters

Active bindings

fitted

Y values predicted by the model, in Y's original units

var_par

a list with the variational parameter: tau (posterior group probabilities)

nb_param

number of parameters in the model

entropy

Entropy of the conditional distribution The only latent variable is the clustering, so the entropy of the variational distribution reduces to -sum(tau * log(tau)).

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockMeanUnknownClusters$new()

Create a new ['NormalBlockMeanUnknownClusters'] object.

Usage
NormalBlockMeanUnknownClusters$new(
  data,
  q,
  sparsity = 0,
  control = NB_control()
)
Arguments
data

object of NormalBlockData class, with responses and design matrix

q

number of clusters

sparsity

to apply on variance matrix when calling GLASSO

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['NormalBlockMeanUnknownClusters'] object

NormalBlockMeanUnknownClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockMeanUnknownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
model <- normal_block(data, blocks = 3, model = "mean")
model$clustering

Base Class for Variance-Block Models

Description

R6 abstract class for the sparse Normal-Block models, where the clustering structures the latent covariance.

Super class

NormalBlockBase -> NormalBlockVarBase

Active bindings

B_original

regression coefficients (d x p), converted back to Y's original units (undoing ‘NormalBlockData(scale = TRUE)'’s column-wise rescaling, if any). Use 'model_par$B' instead for the coefficients on the internal fitting scale.

d0

number of zi variables (dimensions in X0)

model_par

a list with the matrices of the model parameters: B (covariates), dm1 (species variance), Omega (groups precision matrix)). On the internal fitting scale (‘self$data$Y', possibly column-rescaled by 'NormalBlockData(scale = TRUE)') – use '$B_original'/'$dm1_original' for the same quantities converted back to Y’s original units.

nb_param

number of parameters in the model

sparsity_weights

(weights associated to each pair of groups)

dm1_original

inverse residual variance per variable (1 / Var(Y_j)), converted back to Y's original units. Use 'model_par$dm1' instead for the internal fitting scale. With 'noise_covariance = "spherical"', 'model_par$dm1' is a single value repeated p times; once converted back per-variable, the p values returned here generally differ from one another whenever Y's columns were rescaled by different factors.

Methods

Public methods

Inherited methods

NormalBlockVarBase$new()

Create a new ['NormalBlockVarBase'] object.

Usage
NormalBlockVarBase$new(
  data,
  q,
  sparsity = 0,
  control = NB_control(),
  zero_inflation = FALSE
)
Arguments
data

object of NormalBlockData class, with responses and design matrix

q

number of block/cluster

sparsity

sparsity penalty on the network density

control

structured list of more specific parameters, to generate with NB_control

zero_inflation

whether the concrete subclass models zero-inflation; set by the ZI subclasses themselves, not meant to be set by the end user. When 'FALSE', the (costly) zero-inflation probability fit ('kappa'/'B0') is skipped entirely, since it would otherwise never be used downstream.

Returns

A new ['NormalBlockVarBase'] object

NormalBlockVarBase$warm_start_from()

Seed this model's starting parameters from another, already-optimized model with the same q, instead of a fresh heuristic clustering. Used by [NormalBlockVarCollectionSparsity] to warm-start each penalty in a sparsity path from the previous one's solution.

Usage
NormalBlockVarBase$warm_start_from(other)
Arguments
other

a [NormalBlockVarBase] object, already optimized

Returns

Update the current object in place with ‘other'’s parameters

NormalBlockVarBase$split()

Create a clone of the current ['NormalBlockVarBase'] object after splitting cluster 'cl' We split the cluster according to the species variances

Usage
NormalBlockVarBase$split(index, in_place = FALSE)
Arguments
index

index (integer) of the cluster to split

in_place

should the split applied to the object itself, or should a copy be sent? default FALSE (send a copy)

Returns

A new ['NormalBlockVarBase'] object

NormalBlockVarBase$merge()

Create a clone of the current ['NormalBlockVarBase'] object after merging clusters 'cl1' and 'cl2'

Usage
NormalBlockVarBase$merge(indices, in_place = FALSE)
Arguments
indices

indices (couple of integer) of the clusters to merge

in_place

should the split applied to the object itself, or should a copy be sent? default FALSE (send a copy)

Returns

A new ['NormalBlockVarBase'] object

NormalBlockVarBase$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockVarBase$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

# An internal abstract base class, never instantiated directly. See
# normal_block() for how concrete models (NormalBlockVarKnownClusters,
# NormalBlockVarUnknownClusters, and their zero-inflated variants) are
# actually created and fitted.

Collection of Normal-Block Models over a Range of Cluster Counts

Description

R6 class for a collection of normal-block models with different number of clusters (q) and a fixed sparsity level.

Super classes

NormalBlockCollection -> NormalBlockCollectionClusters -> NormalBlockVarCollectionClusters

Active bindings

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockVarCollectionClusters$new()

Create a new ['NormalBlockVarCollectionClusters'] object.

Usage
NormalBlockVarCollectionClusters$new(
  mydata,
  q_list,
  zero_inflation = FALSE,
  sparsity = 0,
  control = NB_control()
)
Arguments
mydata

object of NormalBlockData class, with responses and design matrix

q_list

list of q values (number of groups) in the collection

zero_inflation

whether the models in the collection should be zero-inflated or not

sparsity

sparsity penalty on the network density

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['NormalBlockVarCollectionClusters'] object

NormalBlockVarCollectionClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockVarCollectionClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
models <- normal_block(data, blocks = 2:4, control = NB_control(verbose = FALSE))
models$get_best_model("ICL")$q

Collection of Normal-Block Models over Cluster Counts and Sparsity Levels

Description

R6 class for a collection of normal-block models with different number of clusters (q) and different sparsity levels.

Super classes

NormalBlockCollection -> NormalBlockCollectionClustersSparsity -> NormalBlockVarCollectionClustersSparsity

Active bindings

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockVarCollectionClustersSparsity$new()

Create a new ['NormalBlockVarCollectionClustersSparsity'] object.

Usage
NormalBlockVarCollectionClustersSparsity$new(
  mydata,
  q_list,
  zero_inflation = FALSE,
  control = NB_control()
)
Arguments
mydata

object of NormalBlockData class, with responses and design matrix

q_list

list of q values (number of groups) in the collection

zero_inflation

boolean to specify whether data is zero-inflated

control

structured list of parameters to handle sparsity control

Returns

A new ['NormalBlockVarCollectionClustersSparsity'] object

NormalBlockVarCollectionClustersSparsity$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockVarCollectionClustersSparsity$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
models <- normal_block(data, blocks = 2:3, sparsity = TRUE,
                       control = NB_control(verbose = FALSE, n_sparsity_penalties = 3))
models$get_best_model("BIC")$q

Collection of Normal-Block Models over a Sparsity Path

Description

R6 class for a collection of normal-block models with a fixed clustering (blocks) and different sparsity levels.

Super classes

NormalBlockCollection -> NormalBlockCollectionSparsity -> NormalBlockVarCollectionSparsity

Active bindings

sparsity_details

list of information about model's penalties

criteria

a data frame with the values of some criteria ((approximated) log-likelihood, BIC) for the collection of models

stability_path

measure of edges stability based on StARS method

stability

mean edge stability along the sparsity penalties path

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockVarCollectionSparsity$new()

Create a new ['NormalBlockVarCollectionSparsity'] object.

Usage
NormalBlockVarCollectionSparsity$new(
  mydata,
  blocks,
  zero_inflation = FALSE,
  control = NB_control()
)
Arguments
mydata

object of NormalBlockData class, with responses and design matrix

blocks

either a clustering matrix (known, fixed clustering) or a single integer (number of blocks to infer)

zero_inflation

boolean to specify whether data is zero-inflated

control

structured list of parameters to handle sparsity control

Returns

A new ['NormalBlockVarCollectionSparsity'] object

NormalBlockVarCollectionSparsity$get_best_model()

Extract best model in the collection

Usage
NormalBlockVarCollectionSparsity$get_best_model(
  crit = c("BIC", "EBIC", "ICL", "StARS"),
  stability = 0.9
)
Arguments
crit

a character for the criterion used to performed the selection.

stability

if criterion = "StARS" gives level of stability required. Either "BIC", "EBIC", "ICL" or "StARS". Default is BIC

Returns

a ['NormalBlockVarUnknownClusters'] object

NormalBlockVarCollectionSparsity$stability_selection()

Compute the stability path by stability selection

Usage
NormalBlockVarCollectionSparsity$stability_selection(
  subsamples = NULL,
  n_subsamples = 10
)
Arguments
subsamples

a list of vectors describing the subsamples. The number of vectors (or list length) determines the number of subsamples used in the stability selection. Automatically set to 20 subsamples with size '10*sqrt(n)' if 'n >= 144' and '0.8*n' otherwise following Liu et al. (2010) recommendations.

n_subsamples

number of subsamples to create if the subsamples are not given

NormalBlockVarCollectionSparsity$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockVarCollectionSparsity$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
models <- normal_block(data, blocks = ex$parameters$C, sparsity = TRUE,
                       control = NB_control(verbose = FALSE, n_sparsity_penalties = 5))
models$get_best_model("BIC")$sparsity

Normal-Block Model with Known Clustering

Description

R6 class for a normal-block model with known clustering.

Super classes

NormalBlockBase -> NormalBlockVarBase -> NormalBlockVarKnownClusters

Active bindings

posterior_par

a list with the parameters of posterior distribution W | Y

entropy

Entropy of the conditional distribution

fitted

Y values predicted by the model, in Y's original units

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockVarKnownClusters$new()

Create a new ['NormalBlockVarKnownClusters'] object.

Usage
NormalBlockVarKnownClusters$new(data, C, sparsity = 0, control = NB_control())
Arguments
data

object of NormalBlockVarData class, with responses and design matrix

C

clustering matrix C_jk = 1 if species j belongs to cluster k

sparsity

to apply on variance matrix when calling GLASSO

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['NormalBlockVarKnownClusters'] object

NormalBlockVarKnownClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockVarKnownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
model <- normal_block(data, blocks = ex$parameters$C, control = NB_control(verbose = FALSE))
model$clustering

Normal-Block Model with Unknown Clustering

Description

R6 class for a normal-block model with a fixed number of clusters (but unknown clustering).

Super classes

NormalBlockBase -> NormalBlockVarBase -> NormalBlockVarUnknownClusters

Public fields

fixed_tau

whether tau should be fixed at clustering_init during optimization, useful for stability selection

Active bindings

model_par

a list with the matrices of the model parameters: B (covariates), dm1 (species variance), Omega (groups precision matrix))

nb_param

number of parameters in the model

var_par

a list with the matrices of the variational parameters: M (means), S (variances), tau (posterior group probabilities)

entropy

Entropy of the conditional distribution

fitted

Y values predicted by the model, in Y's original units

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

NormalBlockVarUnknownClusters$new()

Create a new ['NormalBlockVarUnknownClusters'] object.

Usage
NormalBlockVarUnknownClusters$new(
  data,
  q,
  sparsity = 0,
  control = NB_control()
)
Arguments
data

contains the matrix of responses (Y) and the design matrix (X).

q

required number of groups

sparsity

sparsity penalty to add on blocks precision matrix for sparsity

control

structured list for specific parameters

Returns

A new ['NormalBlockVarUnknownClusters'] object

NormalBlockVarUnknownClusters$clone()

The objects of this class are cloneable with this method.

Usage
NormalBlockVarUnknownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
model <- normal_block(data, blocks = 3, control = NB_control(verbose = FALSE))
model$clustering

Zero-Inflated Mean-Block Model with Known Clustering

Description

R6 class for a zero-inflated Normal-Block-Mean model with a known clustering. Sigma is diagonal or spherical here. See [NormalBlockMeanBase] for why a full one is out of reach under a mask.

Super classes

NormalBlockBase -> NormalBlockMeanBase -> ZINormalBlockMeanKnownClusters

Active bindings

fitted

Y values predicted by the model, in Y's original units

model_par

a list with model parameters: B, Omega and kappa (zero-inflation probabilities)

nb_param

number of parameters in the model

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

ZINormalBlockMeanKnownClusters$new()

Create a new ['ZINormalBlockMeanKnownClusters'] object.

Usage
ZINormalBlockMeanKnownClusters$new(
  data,
  C,
  sparsity = 0,
  control = NB_control()
)
Arguments
data

object of NormalBlockData class, with responses and design matrix

C

clustering matrix C_jk = 1 if species j belongs to cluster k

sparsity

unused here, kept for signature symmetry (must be 0)

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['ZINormalBlockMeanKnownClusters'] object

ZINormalBlockMeanKnownClusters$clone()

The objects of this class are cloneable with this method.

Usage
ZINormalBlockMeanKnownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 50, p = 20, d = 1, q = 3)
Y <- ex$Y; Y[runif(length(Y)) < 0.2] <- 0
data <- NormalBlockData$new(Y, ex$X)
model <- normal_block(data, blocks = ex$parameters$C, model = "mean",
                      zero_inflation = TRUE)
model$clustering

Zero-Inflated Mean-Block Model with Unknown Clustering

Description

R6 class for a zero-inflated Normal-Block-Mean model with a fixed number of clusters (but unknown clustering), inferred by variational EM. Sigma is diagonal or spherical here. See [NormalBlockMeanBase] for why a full one is out of reach under a mask.

Super classes

NormalBlockBase -> NormalBlockMeanBase -> ZINormalBlockMeanUnknownClusters

Active bindings

fitted

Y values predicted by the model, in Y's original units

var_par

a list with the variational parameter: tau (posterior group probabilities)

model_par

a list with model parameters: B, Omega and kappa (zero-inflation probabilities)

nb_param

number of parameters in the model

entropy

Entropy of the conditional distribution

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

ZINormalBlockMeanUnknownClusters$new()

Create a new ['ZINormalBlockMeanUnknownClusters'] object.

Usage
ZINormalBlockMeanUnknownClusters$new(
  data,
  q,
  sparsity = 0,
  control = NB_control()
)
Arguments
data

object of NormalBlockData class, with responses and design matrix

q

number of clusters

sparsity

unused here, kept for signature symmetry (must be 0)

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['ZINormalBlockMeanUnknownClusters'] object

ZINormalBlockMeanUnknownClusters$clone()

The objects of this class are cloneable with this method.

Usage
ZINormalBlockMeanUnknownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_mean_data(n = 50, p = 20, d = 1, q = 3)
Y <- ex$Y; Y[runif(length(Y)) < 0.2] <- 0
data <- NormalBlockData$new(Y, ex$X)
model <- normal_block(data, blocks = 3, model = "mean", zero_inflation = TRUE)
model$clustering

Zero-Inflated Normal-Block Model with Known Clustering

Description

R6 class for a zero-inflated normal-block model with a known clustering.

Super classes

NormalBlockBase -> NormalBlockVarBase -> ZINormalBlockVarKnownClusters

Active bindings

posterior_par

a list with the parameters of posterior distribution W | Y

entropy

Entropy of the conditional distribution

nb_param

number of parameters in the model

model_par

a list with model parameters: B (covariates), dm1 (species variance), Omega (groups precision matrix), kappa (zero-inflation probabilities)

fitted

Y values predicted by the model, in Y's original units

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

ZINormalBlockVarKnownClusters$new()

Create a new ['ZINormalBlockVarKnownClusters'] object.

Usage
ZINormalBlockVarKnownClusters$new(
  data,
  C,
  sparsity = 0,
  control = NB_control()
)
Arguments
data

object of NormalBlockData class, with responses and design matrix

C

clustering matrix C_jk = 1 if species j belongs to cluster k

sparsity

to apply on variance matrix when calling GLASSO

control

structured list of more specific parameters, to generate with NB_control

Returns

A new ['ZINormalBlockVarKnownClusters'] object

ZINormalBlockVarKnownClusters$clone()

The objects of this class are cloneable with this method.

Usage
ZINormalBlockVarKnownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3, kappa = rep(0.3, 20))
data <- NormalBlockData$new(ex$Y, ex$X)
model <- normal_block(data, blocks = ex$parameters$C, zero_inflation = TRUE,
                      control = NB_control(verbose = FALSE))
model$clustering

Zero-Inflated Normal-Block Model with Unknown Clustering

Description

R6 class for a zero-inflated normal-block model with a fixed number of clusters (but unknown clustering).

Super classes

NormalBlockBase -> NormalBlockVarBase -> ZINormalBlockVarUnknownClusters

Public fields

fixed_tau

whether tau should be fixed at clustering_init during optimization, useful for stability selection

Active bindings

nb_param

number of parameters in the model

var_par

a list with variational parameters

model_par

a list with model parameters: B (covariates), dm1 (species variance), Omega (blocks precision matrix), kappa (zero-inflation probabilities)

entropy

Entropy of the conditional distribution

fitted

Y values predicted by the model, in Y's original units

who_am_I

a method to print what model is being fitted

Methods

Public methods

Inherited methods

ZINormalBlockVarUnknownClusters$new()

Create a new ['ZINormalBlockVarUnknownClusters'] object.

Usage
ZINormalBlockVarUnknownClusters$new(
  data,
  q,
  sparsity = 0,
  control = NB_control()
)
Arguments
data

object of NormalBlockVarData class, with responses and design matrix

q

required number of groups

sparsity

to apply on variance matrix when calling GLASSO

control

structured list of more specific parameters

Returns

A new ['ZINormalBlockVarUnknownClusters'] object

ZINormalBlockVarUnknownClusters$clone()

The objects of this class are cloneable with this method.

Usage
ZINormalBlockVarUnknownClusters$clone(deep = FALSE)
Arguments
deep

Whether to make a deep clone.

Examples

ex <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3, kappa = rep(0.3, 20))
data <- NormalBlockData$new(ex$Y, ex$X)
model <- normal_block(data, blocks = 3, zero_inflation = TRUE,
                      control = NB_control(verbose = FALSE))
model$clustering

Breast cancer proteomics data (TCGA, RPPA)

Description

Reverse-phase protein array (RPPA) measurements of 163 proteins across 346 breast cancer tumor samples from The Cancer Genome Atlas (TCGA), along with the PAM50 molecular subtype of each sample. Used as the running example in Tous & Chiquet (2026).

Usage

brca_rppa

Format

A list with 3 elements:

expr

a 346 x 163 numeric matrix of protein expression levels (normalized log-ratios), samples in rows (named by TCGA sample identifier), proteins in columns.

covariates

a data frame with one row per sample, in the same order as the rows of 'expr': 'sampleId' (matches 'rownames(expr)') and 11 clinical variables, factors unless noted otherwise – 'CN_CLUSTER' (copy-number cluster, 5 levels), 'CONVERTED_STAGE' (tumor stage), 'ER_STATUS' (estrogen receptor status), 'FRACTION_GENOME_ALTERED' (numeric, in \[0, 1\]), 'HER2_STATUS', 'METASTASIS_CODED', 'METHYLATION_CLUSTER' (5 levels), 'MUTATION_COUNT' (numeric), 'PAM50_SUBTYPE' (5 levels: Basal-like, HER2-enriched, Luminal A, Luminal B, Normal-like), 'RPPA_CLUSTER' (7 levels) and 'AGE' (numeric, in years). A few of these have missing values for a handful of samples ('FRACTION_GENOME_ALTERED', 'HER2_STATUS', 'METASTASIS_CODED').

gene_annotation

a data frame with one row per protein, in the same order as the columns of 'expr': 'protein' (matches 'colnames(expr)'), 'entrezGeneId', 'hugoGeneSymbol', 'type' ('"protein-coding"' or '"phosphoprotein"') and 'go_bp_term' (a Gene Ontology Biological Process classification, one term per protein, or '"unknown"' if none could be found – see Details).

Details

'go_bp_term' is derived from 'org.Hs.eg.db', queried by 'entrezGeneId' for '"protein-coding"' rows. Phospho-site antibodies ('type == "phosphoprotein"', e.g. '"EGFR_PY1068"') carry a placeholder, negative 'entrezGeneId' (there is no separate gene for a specific phosphorylation site, only the underlying gene) and are instead queried by gene symbol (the part of 'protein' before the first '"_"'). Each protein typically has dozens of GO Biological Process terms; 'go_bp_term' keeps, for each protein, the one term shared by the most *other* proteins in this dataset – giving a handful of non-trivial groups rather than one near-singleton group per protein, more useful for comparing against an inferred clustering. See 'data-raw/brca_rppa_go_annotation.R' for the full extraction code.

Source

TCGA breast cancer cohort (Cancer Genome Atlas Network, 2012, doi:10.1038/nature11412), downloaded via cBioPortal ('brca_tcga_pub' study, https://www.cbioportal.org). 'go_bp_term' (see Details) was added afterwards from 'org.Hs.eg.db'.

Examples

expr <- brca_rppa$expr
X <- model.matrix(~ 0 + PAM50_SUBTYPE, data = brca_rppa$covariates)
nb_data <- NormalBlockData$new(expr, X)
table(brca_rppa$gene_annotation$go_bp_term)

Extract Model Coefficients

Description

Extract coefficients from a normal-block model.

Usage

## S3 method for class 'NormalBlockBase'
coef(object, ...)

Arguments

object

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

A matrix of coefficients extracted from the model.

Extract Fitted Values

Description

Extract fitted values from a normal-block model.

Usage

## S3 method for class 'NormalBlockBase'
fitted(object, ...)

Arguments

object

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

A matrix of fitted values extracted from the object.

Generate Normal Block Mean Data

Description

A function to draw data from the normal block model (see details). The function returns both the generated data and the corresponding model parameters, in a list.

Usage

generate_normal_block_mean_data(
  n = 100,
  p = 40,
  d = 1,
  q = 3,
  kappa = 0,
  omega_structure = "erdos-renyi",
  u_v = c(0.3, 0.1),
  SNR = 5,
  alpha = rep(1/q, q),
  range_X = c(0, 10)
)

Arguments

n

number of individuals. Default to 100.

p

number of variables. Default to 40.

d

number of covariates. Default to 1.

q

number of groups. Default to 3.

kappa

vector (or scalar) of variable-wise probability of zero inflation. Default to 0.

omega_structure

the structure of the graph on which the precision matrix between variables is built. Can be a symmetric matrix with p rows/columns or a character picked in "erdos-renyi", "preferential_attachment", "community" in which case a graph is drawn with sensible generation parameters. See generate_precision_matrix for details.

u_v

two-size vector of positive numbers v and u controlling the generation of the precision matrix Omega: v scales the off-diagonal elements of the precision matrix (magnitude of partial correlations), and u is a positive number added to the diagonal elements to ensure positive-definiteness. The default value is c(0.3, 0.1).

SNR

Signal to noise ratio: magnitude of the regression parameters B will be adjusted so that tr(var(XB)) and tr(Sigma) match the desired SNR.

alpha

the q-size vector of group proportion. Default to rep(1/q, q)

range_X

A 2-size vector defining the range of the uniform distribution used to draw values in X, the regressor matrix. Default is c(0, 10)

Value

A named list with the following element - Y a matrix of responses - X a regressor/design matrix - a list of model parameters, encompassing - B: matrix of regression coefficients - C: matrix of group membership - Omega: precision matrix of the variables - Sigma: covariance matrix of the variables - kappa: vector of ZI inflation probabilities (one per variable)

Generate Normal Block Var Data

Description

A function to draw data from the normal block model (see details). The function returns both the generated data and the corresponding model parameters, in a list.

Usage

generate_normal_block_var_data(
  n = 100,
  p = 40,
  d = 1,
  q = 3,
  kappa = 0,
  omega_structure = "erdos-renyi",
  u_v = c(0.3, 0.1),
  SNR = 0.75,
  alpha = rep(1/q, q),
  range_X = c(0, 10),
  range_D = c(0.5, 1.5)
)

Arguments

n

number of individuals. Default to 100.

p

number of variables. Default to 40.

d

number of covariates. Default to 1.

q

number of groups. Default to 3.

kappa

vector (or scalar) of variable-wise probability of zero inflation. Default to 0.

omega_structure

the structure of the graph on which the precision matrix between groups is built. Can be a symmetric matrix with q rows/columns or a character picked in "erdos-renyi", "preferential_attachment", "community" in which case a graph is drawn with sensible generation parameters. See generate_precision_matrix for details.

u_v

two-size vector of positive numbers v and u controlling the generation of the precision matrix Omega: v scales the off-diagonal elements of the precision matrix (magnitude of partial correlations), and u is a positive number added to the diagonal elements to ensure positive-definiteness. The default value is c(0.3, 0.1).

SNR

Signal to noise ratio: magnitude of the regression parameters B will be adjusted so that tr(var(XB)) and tr(Sigma) match the desired SNR.

alpha

the q-size vector of group proportion. Default to rep(1/q, q)

range_X

A 2-size vector defining the range of the uniform distribution used to draw values in X, the regressor matrix. Default is c(0, 10)

range_D

A 2-size vector defining the range of the uniform distribution used to draw values in D, the diagonal matrix of variances of variables. Default is c(0.5, 1.5)

Value

A named list with the following element - Y a matrix of responses - X a regressor/design matrix - a list of model parameters, encompassing - B: matrix of regression coefficients - C: matrix of group membership - D: diagonal matrix of variance of the variables - Omega: precision matrix of the groups - Sigma: covariance matrix of the groups - kappa: vector of ZI inflation probabilities (one per variable)

Create a Normal-Block Model Object

Description

Creates the appropriate normal-block model (or collection of models) depending on the parametrization.

Usage

get_model(
  data,
  blocks,
  sparsity = 0,
  zero_inflation = FALSE,
  control = NB_control(),
  model = c("var", "mean")
)

Arguments

data

contains the matrix of responses (Y) and the design matrix (X).

blocks

either an integer (number of blocks), a vector of integer (list of possible number of block) or a p * q matrix (for indicating block membership when its known)

sparsity

boolean to say whether the model should have a changing penalty OR float to run model with a single penalty value

zero_inflation

boolean to indicate if Y is zero-inflated and adjust fitted model as a consequence

control

a list-like structure for detailed control on parameters should be generated with NB_control() for collections of sparse models

model

which model family to fit, "var" (the default) or "mean" – see [normal_block()]

Check if an Object is a Normal-Block Model

Description

Checks if a model is of class [NormalBlockBase()], i.e. any normal-block model (variance-block or mean-block family).

Usage

isNB(object)

Arguments

object

An R object.

Value

A boolean telling whether object inherits from the NormalBlockBase class.

Extract Log-Likelihood of a Normal-Block Model

Description

Returns the (variational) log-likelihood of a fitted normal-block model as a '"logLik"' object, compatible with [stats::AIC()] and [stats::BIC()].

Usage

## S3 method for class 'NormalBlockBase'
logLik(object, ...)

Arguments

object

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

An object of class '"logLik"'. The numeric value is the log-likelihood or its variational lower bound (ELBO). Attributes 'df' and 'nobs' hold the number of parameters and observations.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
model <- normal_block(data, blocks = 3, control = NB_control(verbose = FALSE))
logLik(model)

Extract Log-Likelihood of a Collection of Normal-Block Models

Description

Returns the log-likelihood of every model in a collection of normal-block models (see [NormalBlockVarCollectionClusters], [NormalBlockVarCollectionSparsity], [NormalBlockVarCollectionClustersSparsity]).

Usage

## S3 method for class 'NormalBlockCollection'
logLik(object, ...)

Arguments

object

An object inheriting from NormalBlockCollection.

...

not used, only here for S3 compatibility

Value

A numeric vector of log-likelihood values, one per model in the collection.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
models <- normal_block(data, blocks = 2:5, control = NB_control(verbose = FALSE))
logLik(models)

Normal-block model

Description

Fit a normal-block model with a variational or heuristic algorithm

Usage

normal_block(
  data,
  blocks,
  sparsity = 0,
  zero_inflation = FALSE,
  control = NB_control(),
  model = c("var", "mean")
)

Arguments

data

NormalBlockData object, contains the matrix of responses (Y, n x p) and the design matrix (X, n x d), must be created with NormalBlockData$new.

blocks

either a integer (number of blocks), a vector of integer (list of possible number of block) or a p * q matrix (for indicating block membership when its known)

sparsity

either TRUE to run the optimization for different sparsity penalty values OR float to run model with a single sparsity penalty value

zero_inflation

boolean to indicate if Y is zero-inflated and adjust fitted model as a consequence

control

a list-like structure for detailed control on parameters should be generated with NB_control().

model

which model family to fit: "var" (the default) structures the clustering in the latent covariance (see [NormalBlockVarBase]), "mean" structures it in the mean (mu_i = C B' X_i, see [NormalBlockMeanBase]). The mean-block family covers the same collections as the variance-block one ([NormalBlockMeanCollectionClusters], [NormalBlockMeanCollectionSparsity] and their crossing, [NormalBlockMeanCollectionClustersSparsity]) – each q simply fit independently, without the variance-block family's SBM-path shortcut. Zero-inflation is supported ([ZINormalBlockMeanKnownClusters], [ZINormalBlockMeanUnknownClusters], and collections of those over a range of q), with a diagonal or spherical Sigma only, hence not along a sparsity path – see [NormalBlockMeanBase].

Value

an R6 object with one of the model classes (or a collection of model objects).

Examples

## Normal Data
ex_data <- generate_normal_block_var_data(n=100, p=30, d=1, q=3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
my_normal_block <- normal_block(data, blocks = 1:6)
my_normal_block$plot(c("deviance", "BIC", "ICL"))
Y_hat <- my_normal_block$get_best_model()$fitted
plot(data$Y, Y_hat, log = "xy"); abline(0,1)
## Normal Data with Zero Inflation
ex_data_zi <- generate_normal_block_var_data(n=50, p=50, d=1, q=3, kappa = rep(0.5,50))
zidata <- NormalBlockData$new(ex_data_zi$Y, ex_data_zi$X)
my_normal_block <- normal_block(zidata, blocks = 1:6, zero_inflation = TRUE)
## Mean-Block model (clustering in the mean rather than the covariance)
ex_mean <- generate_normal_block_mean_data(n=50, p=20, d=1, q=3)
mean_data <- NormalBlockData$new(ex_mean$Y, ex_mean$X)
my_mean_block <- normal_block(mean_data, blocks = 3, model = "mean")

Cluster variables in the mean, then in the residual covariance

Description

Fits the two model families one after the other: a mean-block model ([NormalBlockMeanBase]) groups the variables by how they respond to the covariates, then a variance-block model ([NormalBlockVarBase]) groups the residuals of that fit by how they co-vary. The two answer different questions and generally return unrelated partitions, so running both is often more informative than choosing one.

Usage

normal_block_sequential(
  data,
  blocks_mean,
  blocks_var,
  crit = c("ICL", "BIC"),
  zero_inflation = FALSE,
  control_mean = NB_control(verbose = FALSE),
  control_var = NB_control(verbose = FALSE)
)

Arguments

data

a [NormalBlockData] object

blocks_mean

number of clusters for the mean-block stage: an integer, a vector of integers to explore, or a p x q indicator matrix

blocks_var

idem for the variance-block stage, run on the residuals

crit

criterion used to pick a model when a range is explored, "ICL" (the default) or "BIC"

zero_inflation

whether Y carries structural zeros. Both stages are then zero-inflated, sharing the same mask: the second one would otherwise re-derive it from residuals, which are never exactly zero.

control_mean

control list for the mean-block stage, see [NB_control()]

control_var

control list for the variance-block stage

Details

The second stage uses an intercept-only design on purpose: the covariate effects have already been removed by the first stage.

This is a heuristic two-stage estimator, not a joint model. On simulated data carrying two genuinely distinct structures it recovers both exactly (see 'inst/mean_block_analyses/sequential_mean_then_variance.R').

Value

an object of class 'normal_block_sequential', a list with the fitted 'mean' and 'var' models and the residual matrix 'residuals' handed from one stage to the other.

Examples

ex   <- generate_normal_block_mean_data(n = 80, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex$Y, ex$X)
fit  <- normal_block_sequential(data, blocks_mean = 3, blocks_var = 2)
fit

French stream fish community data (ONEMA / OFB electrofishing surveys)

Description

Fish biomass per species, aggregated by sampling station, together with environmental covariates for each station. Derived from the fish community monitoring of stream sections across metropolitan France (1995-2018) by the French Office for Biodiversity (formerly ONEMA, "Office National de l'Eau et des Milieux Aquatiques"), using electrofishing.

Usage

onema

Format

A list with 2 elements:

biomass

a 399 x 46 numeric matrix of total biomass (grams) per species (3-letter species code, columns) and station (row names).

covariates

a data frame with one row per station, in the same order as the rows of 'biomass': 'station' (matches 'rownames(biomass)') and 14 environmental variables – 'slope', 'alt' (altitude), 'd_source' (distance to source), 'strahler' (Strahler stream order), 'width_river_mean'/'width_river_cv', 'avg_depth_station_mean'/ 'avg_depth_station_cv', 'DBO_med'/'DBO_cv' (biological oxygen demand), 'flow_med'/'flow_cv' and 'temperature_med'/'temperature_cv' ('_med' = median, '_cv' = coefficient of variation, over the monitoring period).

Details

The source data is organized around fishing *operations* ('opcod'): several electrofishing operations can be carried out at the same *station* over time. 'onema' aggregates this to one row per station: every operation is first mapped to its station ('fishing_protocol.csv'), then biomass (in grams) is summed over every operation recorded at a given station, separately for each species (i.e. 'biomass[s, ]' is the total biomass of each species ever caught at station 's', not a single operation's catch). Stations without environmental data are dropped. Species observed at fewer than 10 stations are also dropped: at that level of rarity, a zero-inflated fit's iterative initialization can fit a species' handful of non-zero observations exactly, driving its estimated noise precision to infinity (these species also carry essentially no information for clustering anyway).

Source

Danet, A., Mouchet, M., Bonnaffe, W., Thebault, E., Fontaine, C. (2021) "Species richness and food-web structure jointly drive community biomass and its temporal stability in fish communities", data set, Zenodo, doi:10.5281/zenodo.5095656.

Examples

Y <- log(1 + onema$biomass)
X <- model.matrix(~ 1, data = onema$covariates)
nb_data <- NormalBlockData$new(Y, X)

out <- normal_block(nb_data, 2:15, control = NB_control(clustering_init = "ward2"))

Plot a Normal-Block Model

Description

Plots the evolution of the objective (log-likelihood or ELBO) across the (V)EM iterations of the last call to 'optimize()', see '$plot_loglik()'.

Usage

## S3 method for class 'NormalBlockBase'
plot(x, ...)

Arguments

x

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

Invisibly returns the [ggplot2::ggplot] object; called for its side effect of plotting.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
model <- normal_block(data, blocks = 3, control = NB_control(verbose = FALSE))
plot(model)

Predict Method for Variance-Block Models

Description

Predicts observations Y for new covariates X. Specific to the variance-block family: the mean-block models have their own formula (mu = C B' X).

Usage

## S3 method for class 'NormalBlockBase'
predict(object, new_X, ...)

Arguments

object

An object of class NormalBlockVarBase.

new_X

New set of covariates.

...

not used, only here for S3 compatibility

Value

A n*p prediction matrix for new observations

Print a Normal-Block Model

Description

Print a short summary of a fitted normal-block model: model type, goodness-of-fit criteria, and the useful fields/methods to explore it further.

Usage

## S3 method for class 'NormalBlockBase'
print(x, ...)

Arguments

x

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

Invisibly returns 'x'.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
model <- normal_block(data, blocks = 3, control = NB_control(verbose = FALSE))
print(model)

Print a Collection of Normal-Block Models

Description

Print a short summary of a collection of normal-block models: model type and the range of q/sparsity explored. See [summary.NormalBlockCollection()] for the full criteria table.

Usage

## S3 method for class 'NormalBlockCollection'
print(x, ...)

Arguments

x

An object inheriting from NormalBlockCollection.

...

not used, only here for S3 compatibility

Value

Invisibly returns 'x'.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
models <- normal_block(data, blocks = 2:5, control = NB_control(verbose = FALSE))
print(models)

Print a Sequential Mean-then-Variance Fit

Description

Reports both stages and, when aricode is available, how related the two partitions are.

Usage

## S3 method for class 'normal_block_sequential'
print(x, ...)

Arguments

x

an object of class 'normal_block_sequential'

...

not used, only here for S3 compatibility

Value

Invisibly returns 'x'.

Print a Normal-Block Model Summary

Description

Print method for objects returned by [summary.NormalBlockBase()].

Usage

## S3 method for class 'summary.NormalBlockBase'
print(x, ...)

Arguments

x

An object of class 'summary.NormalBlockBase'.

...

not used, only here for S3 compatibility

Value

Invisibly returns 'x'.

Print a Collection Summary

Description

Print method for objects returned by [summary.NormalBlockCollection()].

Usage

## S3 method for class 'summary.NormalBlockCollection'
print(x, ...)

Arguments

x

An object of class 'summary.NormalBlockCollection'.

...

not used, only here for S3 compatibility

Value

Invisibly returns 'x'.

Extract the Covariance Matrix

Description

Extract the covariance matrix 'Omega^-1': between latent blocks (q x q) for the variance-block models, between variables (p x p) for the mean-block models.

Usage

## S3 method for class 'NormalBlockBase'
sigma(object, ...)

Arguments

object

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

The covariance matrix, of size q x q or p x p depending on the model.

Summarize a Normal-Block Model

Description

Summarizes a fitted normal-block model: model type, goodness-of-fit criteria, cluster sizes, and the density of the inferred network between blocks.

Usage

## S3 method for class 'NormalBlockBase'
summary(object, ...)

Arguments

object

An object of class NormalBlockBase.

...

not used, only here for S3 compatibility

Value

An object of class 'summary.NormalBlockBase' (a list with the model's 'who_am_I', 'criteria', 'cluster_sizes' and network 'density'), printed with a dedicated [print.summary.NormalBlockBase()] method.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
model <- normal_block(data, blocks = 3, control = NB_control(verbose = FALSE))
summary(model)

Summarize a Collection of Normal-Block Models

Description

Summarizes a collection of normal-block models: model type, the full criteria table, and the range of q/sparsity explored.

Usage

## S3 method for class 'NormalBlockCollection'
summary(object, ...)

Arguments

object

An object inheriting from NormalBlockCollection.

...

not used, only here for S3 compatibility

Value

An object of class 'summary.NormalBlockCollection' (a list with the collection's 'who_am_I', 'criteria', 'q_range' and 'sparsity_range'), printed with a dedicated [print.summary.NormalBlockCollection()] method.

Examples

ex_data <- generate_normal_block_var_data(n = 50, p = 20, d = 1, q = 3)
data <- NormalBlockData$new(ex_data$Y, ex_data$X)
models <- normal_block(data, blocks = 2:5, control = NB_control(verbose = FALSE))
summary(models)

University webpages text data (CMU "4 Universities" / WebKB)

Description

Term-frequency data derived from the student personal webpages of the CMU "World Wide Knowledge Base" (WebKB) project's "4 Universities" dataset (Cornell, Texas, Washington, Wisconsin computer science departments, collected January 1997). Used as one of the running examples in Tous & Chiquet (2026).

Usage

university

Format

A list with 3 elements:

frequencies

a 504 x 1867 numeric matrix of term frequencies (row sums to 1): for each document (row, named by its original file path) and term (column), the fraction of that document's (post-preprocessing) word count made up of that term.

entropies

a named numeric vector of length 1867 (one value per column of 'frequencies', same names/order), the normalized Shannon entropy of each term's distribution across documents.

terms

a character vector of the 100 column names of 'frequencies' with the highest 'entropies', ordered by decreasing entropy.

Details

Each of the 504 pages (rows) is a document; each of the 1867 columns is a term retained after lower-casing, stripping URLs/HTML tags/control characters/punctuation/numbers, removing English stopwords and dropping terms occurring in fewer than 2 documents. 'entropies' ranks every term by how evenly it is spread across documents (Shannon entropy of each term's normalized document distribution, in \[0, 1\], 1 = perfectly uniform); 'terms' is the 100 highest-entropy terms, i.e. the terms that are the most informative for distinguishing documents from one another rather than just reflecting a few documents' idiosyncratic vocabulary – the transformation used by Tan et al. (2015).

Source

CMU Text Learning Group, "World Wide Knowledge Base (Web->KB) project", http://www.cs.cmu.edu/~webkb/; "4 Universities" subset, https://www.cs.cmu.edu/afs/cs/project/theo-20/www/data/. Tan, P.-N., Steinbach, M., Kumar, V. (2015) "Introduction to Data Mining" (transform used to derive 'entropies'/'terms').

Examples

Y <- log(1 + university$frequencies[, university$terms])
nb_data <- NormalBlockData$new(Y, X = matrix(1, nrow(Y), 1))

out <- normal_block(nb_data, 2:15)