Why does IntegrateLayers default to 2000 genes for StdAssay even when more HVGs are provided?

[wangbenwang123]

I noticed that when using IntegrateLayers with a StdAssay (e.g., RNA), the function defaults to using only 2000 genes (VariableFeatures(..., nfeatures = 2000L)), even if:

1. FindVariableFeatures was previously run with nfeatures > 2000 (e.g., 4000).
2. The PCA (RunPCA) was computed using the full set of HVGs (~4000 genes).

Code Example:

# Step 1: Select 4000 HVGs
sc <- FindVariableFeatures(sc, nfeatures = 4000, assay = "RNA")  
sc <- ScaleData(sc, assay = 'RNA')

# Step 2: RunPCA uses ~4000 genes
sc <- RunPCA(sc, assay = "RNA")  

# Step 3: IntegrateLayers silently uses only 2000 genes
sc <- IntegrateLayers(object = sc, method = RPCAIntegration, assay = "RNA")

 if (inherits(x = object[[assay]], what = "SCTAssay")) {
    layers <- "data"
    scale.layer <- "scale.data"
    features <- features %||% SelectSCTIntegrationFeatures(object = object, 
      assay = assay)
  }
  else if (inherits(x = object[[assay]], what = "StdAssay")) {
    layers <- Layers(object = object, assay = assay, search = layers %||% 
      "data")
    scale.layer <- Layers(object = object, search = scale.layer)
    features <- features %||% VariableFeatures(object = object, 
      assay = assay, nfeatures = 2000L)
  }

Questions:

1. Is this intentional? If so, why is there no warning when the integration uses fewer genes than the PCA?
2. How can I force IntegrateLayers to use all HVGs (e.g., 4000) without manually specifying features?
   • Passing features = VariableFeatures(sc) seems to work, but is this the recommended approach?
3. Could this lead to inconsistencies? For example, if PCA captures biological variance from 4000 genes but integration only uses 2000.

Context:

• Seurat v5
• Standard workflow (RNA assay, not SCT).

---