Skip to contents

Least Squares Separate (LSS) Analysis

lss.Rd

Computes trial-wise beta estimates using the Least Squares Separate approach of Mumford et al. (2012). This method fits a separate GLM for each trial, with the trial of interest and all other trials as separate regressors.

Usage

lss(
  Y,
  X,
  Z = NULL,
  Nuisance = NULL,
  method = c("r_optimized", "cpp_optimized", "r_vectorized", "cpp", "naive"),
  block_size = 96
)

Arguments

Y

A numeric matrix of size n × V where n is the number of timepoints and V is the number of voxels/variables

X

A numeric matrix of size n × T where T is the number of trials. Each column represents the design for one trial

Z

A numeric matrix of size n × F representing experimental regressors to include in all trial-wise models. These are regressors we want to model and get beta estimates for, but not trial-wise (e.g., intercept, condition effects, block effects). If NULL, an intercept-only design is used. Defaults to NULL

Nuisance

A numeric matrix of size n × N representing nuisance regressors to be projected out before LSS analysis (e.g., motion parameters, physiological noise). If NULL, no nuisance projection is performed. Defaults to NULL

method

Character string specifying which implementation to use. Options are:

  • "r_optimized" - Optimized R implementation (recommended, default)

  • "cpp_optimized" - Optimized C++ implementation with parallel support

  • "r_vectorized" - Standard R vectorized implementation

  • "cpp" - Standard C++ implementation

  • "naive" - Simple loop-based R implementation (for testing)

block_size

An integer specifying the voxel block size for parallel processing, only applicable when `method = "cpp_optimized"`. Defaults to 96.

Value

A numeric matrix of size T × V containing the trial-wise beta estimates. Note: Currently only returns estimates for the trial regressors (X). Beta estimates for the experimental regressors (Z) are computed but not returned.

Details

The LSS approach fits a separate GLM for each trial, where each model includes:

  • The trial of interest (from column i of X)

  • All other trials combined (sum of all other columns of X)

  • Experimental regressors (Z matrix) - these are modeled to get beta estimates but not trial-wise

If Nuisance regressors are provided, they are first projected out from both Y and X using standard linear regression residualization.

References

Mumford, J. A., Turner, B. O., Ashby, F. G., & Poldrack, R. A. (2012). Deconvolving BOLD activation in event-related designs for multivoxel pattern classification analyses. NeuroImage, 59(3), 2636-2643.

Examples

# Generate example data
n_timepoints <- 100
n_trials <- 10
n_voxels <- 50

# Create trial design matrix
X <- matrix(0, n_timepoints, n_trials)
for(i in 1:n_trials) {
  start <- (i-1) * 8 + 1
  if(start + 5 <= n_timepoints) {
    X[start:(start+5), i] <- 1
  }
}

# Create data with some signal
Y <- matrix(rnorm(n_timepoints * n_voxels), n_timepoints, n_voxels)
true_betas <- matrix(rnorm(n_trials * n_voxels, 0, 0.5), n_trials, n_voxels)
for(i in 1:n_trials) {
  Y <- Y + X[, i] %*% matrix(true_betas[i, ], 1, n_voxels)
}

# Run LSS analysis
beta_estimates <- lss(Y, X)

# With experimental regressors (intercept + condition effects)
Z <- cbind(1, scale(1:n_timepoints))
beta_estimates_with_regressors <- lss(Y, X, Z = Z)

# With nuisance regression (motion parameters)
Nuisance <- matrix(rnorm(n_timepoints * 6), n_timepoints, 6)
beta_estimates_clean <- lss(Y, X, Z = Z, Nuisance = Nuisance)