SIS / ISIS screening for high dimensional models

SIS stands for Sure Independence Screening and ISIS for Iterative Sure Independence Screening. These methods rapidly prescreen the most promising predictors, after which L0glm performs joint sparse selection and estimation. This combination gives excellent stability in high dimensional examples and scales well to very large predictor sets. SIS was introduced by Fan & Lv (2008), and Fan, Samworth & Wu (2009) developed the general likelihood-based ISIS framework. L0glm uses an offset-based screen followed by joint sparse refitting; the controls are documented in sis_glm().

ISIS is enabled by default for ultrahigh-dimensional models through screen="auto". Let n be the number of observations and p_{\mathrm{eligible}} the number of eligible predictors, excluding mandatory terms and predictors excluded by zero starts in every response channel. Automatic ISIS is enabled when

p_{\mathrm{eligible}}>10n.

The same rule applies in L0glm.tune(), independently for each full-data or CV training fit, using its training row count. There is no separate minimum sample size or predictor count.

Screening keeps the full fitting interface while concentrating joint fitting on a small, informative candidate set. Use screen="ISIS" to request it explicitly at any design size.

SIS or ISIS in one call

screen="SIS" performs one marginal screen. screen="ISIS" repeats screening against the current model, incorporating its predictions as offsets, and keeps the best joint-fitting round. screen="auto" enables ISIS for ultrahigh-dimensional designs. screen="none" fits the supplied design directly.

set.seed(31)
X <- matrix(rnorm(250*5000), 250, 5000)
beta <- numeric(5000); beta[c(5, 50, 500)] <- c(2, -1.5, 1)
y <- as.numeric(X %*% beta) + rnorm(250)
fit <- L0glm.fit(X, y, family=gaussian(), lambda="gic",
  screen="ISIS", screen.nthreads=2L, relaxed=TRUE)
L0glm.metrics(fit, beta_true=beta)[c("gic", "TP", "FP", "FN", "TN", "F1_score", "MCC")]
#>   gic TP FP FN   TN F1_score MCC
#> 1 797  3  0  0 4997        1   1

nsis controls a SIS retained set or an ISIS batch. The default ISIS batch is GIC-based. A small batch lets joint fitting remove false positives as the search progresses. Screening is implemented in Rcpp with OpenMP and supports ordinary GLMs, negative binomial, multinomial and Cox models.

fit$screening records selected predictors, screening rounds, elapsed time and thread use. By default screening uses all physical cores, with one thread as the fallback when detection is unavailable. Explicit screen.nthreads values remain available. Coefficients always have the original dimensions, with zeros for predictors excluded by screening.

A million-variable example

This simulation has 500 observations, one million standard-normal predictors and 60 nonzero effects of either sign. The dense design occupies about 4 GB.

n <- 500L
p <- 1000000L
set.seed(123)
X <- matrix(rnorm(n*p), n, p)
active <- sample.int(p, 60L)
beta <- numeric(p)
beta[active] <- runif(60L, -.8, .8)
set.seed(1)
y <- as.numeric(X[, active] %*% beta[active]) + rnorm(n)

fit <- L0glm.fit(X, y, family=gaussian(), lambda="gic",
  screen="ISIS", screen.nthreads=8L, nthreads=1L,
  intercept=FALSE, relaxed=TRUE)
L0glm.metrics(fit, beta_true=beta)[c("gic", "TP", "FP", "FN", "TN", "F1_score", "MCC")]

The default L0glm fit with ISIS screening recovers 40 true effects with only 3 false positives. It finished in just 3.49 s using all eight physical cores (the default). By comparison, five-fold-CV LASSO took 123.13 s, recovering 17 true effects with 50 false positives. The benchmark article reports all methods, including relaxed LASSO, together with their tuning settings.

GIC tuning and cross-validation

gic <- L0glm.tune(X, y, family=gaussian(), screen="ISIS",
  screen.nthreads=8L, tune="gic", nlambda=10, relaxed=TRUE)
cv <- L0glm.tune(X, y, family=gaussian(), screen="ISIS",
  screen.nthreads=8L, tune="cv", nfolds=5, nlambda=10, relaxed=TRUE)

Each fold screens its training observations independently. Validation can also compare a batch grid, for example nsis=c(2L, 5L, 10L) with validation=list(X=V, y=v), while retaining the default GIC batch when nsis=NULL.