Separate Mixscape min_cells from min_de_genes to match Seurat

src/pertpy/tools/_perturbation_efficacy/_mixscape.py

@@ -52,6 +52,7 @@ def mixscape(
         new_class_name: str = "mixscape_class",
         layer: str | None = None,
         min_de_genes: int = 5,
+        min_cells: int = 5,
         logfc_threshold: float = 0.25,
         de_layer: str | None = None,
         test_method: RankGenesMethod = "wilcoxon",
@@ -75,6 +76,7 @@ def mixscape(
             new_class_name: Name of mixscape classification to be stored in `.obs`.
             layer: Key from adata.layers whose value will be used to perform tests on. Default is using `.layers["X_pert"]`.
             min_de_genes: Required number of genes that are differentially expressed for method to separate perturbed and non-perturbed cells.
+            min_cells: Minimum number of cells in a target gene class. If fewer cells remain assigned during classification, all are assigned NP.
             logfc_threshold: Limit testing to genes which show, on average, at least X-fold difference (log-scale) between the two groups of cells.
             de_layer: Layer to use for identifying differentially expressed genes. If `None`, adata.X is used.
             test_method: Method to use for differential expression testing.
@@ -252,7 +254,7 @@ def mixscape(
                         classes[orig_guide_rows[~ko_mask]] = f"{gene} NP"
 
                         current_classes = classes[orig_guide_rows]
-                        if gene_rows_elsewhere.size + (current_classes == gene).sum() < min_de_genes:
+                        if gene_rows_elsewhere.size + (current_classes == gene).sum() < min_cells:
                             classes[gene_rows_elsewhere] = "NP"
                             classes[guide_rows] = "NP"
                             converged = True

tests/tools/_perturbation_efficacy/test_mixscape.py

@@ -26,6 +26,46 @@ def test_mixscape(adata):
     assert sum(ko_result_correct) > ACCURACY_THRESHOLD * NUM_CELLS_PER_GROUP
 
 
+def test_mixscape_min_cells_is_independent_of_min_de_genes():
+    rng = np.random.default_rng(7)
+    n_cells = 6
+    n_genes = 20
+    control = rng.normal(1.0, 0.1, size=(n_cells, n_genes))
+    perturbed = rng.normal(5.0, 0.1, size=(n_cells, n_genes))
+    adata = ad.AnnData(
+        X=np.vstack([control, perturbed]),
+        obs=pd.DataFrame(
+            {"gene_target": ["NT"] * n_cells + ["target_gene"] * n_cells},
+            index=[f"cell_{i}" for i in range(2 * n_cells)],
+        ),
+        var=pd.DataFrame(index=[f"gene_{i}" for i in range(n_genes)]),
+    )
+    adata.layers["X_pert"] = adata.X.copy()
+
+    result = pt.tl.Mixscape().mixscape(
+        adata,
+        pert_key="gene_target",
+        control="NT",
+        test_method="t-test",
+        min_de_genes=10,
+        min_cells=3,
+        copy=True,
+    )
+    target = result.obs["gene_target"] == "target_gene"
+    assert (result.obs.loc[target, "mixscape_class_global"] == "KO").all()
+
+    result = pt.tl.Mixscape().mixscape(
+        adata,
+        pert_key="gene_target",
+        control="NT",
+        test_method="t-test",
+        min_de_genes=1,
+        min_cells=7,
+        copy=True,
+    )
+    assert (result.obs.loc[target, "mixscape_class_global"] == "NP").all()
+
+
 def test_mixscape_sparse_without_scaling_matches_dense(adata):
     adata.layers["X_pert"] = adata.X.toarray()
     dense = pt.tl.Mixscape().mixscape(