Preclustering response channels before sparse Poisson fitting

L0glm.fit(precluster=TRUE) and L0nnpois(precluster=TRUE) cluster and sum response channels, screen predictors on the reduced responses, then refine the sparse fit on every original channel. The workflow supports nonnegative identity-link Poisson models. Standalone glm_cluster() also supports scalar-family Bregman/deviance clustering.

Details

For conditionally independent Poisson channels with an identity link, $$Y^A_{ig}=\sum_{j\in G_g}Y_{ij}\sim \mathrm{Poisson}(x_i^T\Gamma_g+o^A_{ig}),\qquad \Gamma_g=\sum_{j\in G_g}\beta_j,\quad o^A_{ig}=\sum_{j\in G_g}o_{ij}.$$ Nonnegativity preserves row support under summation. The cluster screen uses raw counts and a permissive screen.factor * lambda (default 0.1). Original-channel refinement uses the requested penalty, calculated once from the original predictor universe. Its default is two grouped passes; an explicitly supplied group.maxit takes precedence unless control.precluster$refit.passes is supplied. relaxed retains its usual selected-support coefficient-refitting meaning.

With start=NULL, original-channel refinement broadcasts each candidate's summed aggregate coefficient divided by the original channel count. Explicit starts retain their original coefficient-coordinate meaning; scalar 1 is expanded in normalized coordinates before summation. Starting zeroes remain excluded in their original channels; compulsory terms remain eligible. Coefficients always retain the original predictor dimensions, with excluded rows zero. Shared observation weights remain unchanged; offsets and coefficient prior centres are summed across groups. This workflow requires weights shared across channels.

With cluster.refit=FALSE, the returned model fits summed cluster outcomes at the requested lambda. Its coefficients/predictions have one column per cluster and its metrics describe those aggregated outcomes. This differs from the permissive screening model retained inside an original-channel fit. Aggregate likelihood criteria and original-channel criteria describe different data and should be compared on their own response scales.

control.precluster is a named list. Supported entries are k.max=32L, k=NULL, linkage="ward", criterion="gic", rowbins=1L, size_fac=NULL (use training column totals), nthreads=NULL (physical cores), screen.factor=0.1, refit.passes=2L and clusterids=NULL. Supplied cluster labels reuse a partition; responses are always aggregated afresh. Clustering threads are serial when parallelize=FALSE. Five-row binning is useful for ordered chromatographic scans; the general default of one row avoids dependence on arbitrary observation order.

Clustering uses training response profiles. Observation offsets are carried into the aggregated and original regression stages. Validation uses the original channels for original-channel refinement; aggregate-only validation uses the same training-derived membership with summed responses/offsets. L0glm.tune() supports original-channel refinement, learning a partition once on the full data and separately in each CV training fold. Training column totals determine size factors in that fold. Partitions are reused across lambda values; selection coefficients are initialized independently. Externally supplied clusterids or size_fac are reused as provided. To tune an aggregate-only model, call glm_aggregate() first and tune its returned responses directly.

The workflow does not add a background basis or change compulsory terms. Include a low-df B-spline basis through always.include for a smooth background. Its nonnegative coefficients are refined separately from the selected peaks on the original channels; see compulsory_terms.

The returned ordinary L0glm_fit gains preclustering: memberships, channel names, candidates, screening coefficients, actual penalties, initialization, clustering settings, cut history and phase timings. One-channel original fits, zero-penalty original fits and starts excluding all selectable rows use the ordinary fitter and record the bypass reason.

See also

glm_cluster(), glm_aggregate(), L0glm.fit(), L0nnpois()

Examples

set.seed(8)
X <- diag(80)
B <- matrix(0,80,8)
B[c(12,35,64),] <- matrix(runif(24,2,12),3,8)
Y <- matrix(rpois(80*8, X %*% B + .1),80,8)
fit <- L0nnpois(X,Y,lambda="gic",offset=.1,precluster=TRUE,
  parallelize=FALSE,control.precluster=list(k.max=4L))
dim(coef(fit))
fit$preclustering$candidates
agg <- L0nnpois(X,Y,lambda="gic",offset=.1,precluster=TRUE,
  cluster.refit=FALSE,parallelize=FALSE,control.precluster=list(k=4L))
dim(coef(agg))