diff --git a/slaf/integrations/scanpy.py b/slaf/integrations/scanpy.py
index c672edf..3bb3944 100644
--- a/slaf/integrations/scanpy.py
+++ b/slaf/integrations/scanpy.py
@@ -934,6 +934,7 @@ def highly_variable_genes(
         min_disp: float = 0.5,
         max_disp: float = np.inf,
         n_top_genes: int | None = None,
+        flavor: str = "seurat",
         inplace: bool = True,
     ) -> pd.DataFrame | None:
         """
@@ -951,6 +952,10 @@ def highly_variable_genes(
             max_disp: Maximum dispersion for genes to be considered.
             n_top_genes: Number of top genes to select by dispersion.
                         If specified, overrides min_disp and max_disp criteria.
+            flavor: HVG selection method. Only "seurat" (mean/dispersion-based
+                   selection over the full expression matrix) is currently
+                   implemented. Other scanpy flavors (e.g. "seurat_v3",
+                   "cell_ranger") are not yet supported.
             inplace: Whether to modify the adata object in place. Currently not
                     fully implemented - returns None when True.
 
@@ -959,6 +964,7 @@ def highly_variable_genes(
                                 and highly_variable column. If inplace=True, returns None.
 
         Raises:
+            NotImplementedError: If `flavor` is not "seurat".
             RuntimeError: If the SLAF array is not properly initialized.
 
         Examples:
@@ -986,6 +992,12 @@ def highly_variable_genes(
             Top genes selected: 1000
         """
 
+        if flavor != "seurat":
+            raise NotImplementedError(
+                f"highly_variable_genes(flavor={flavor!r}) is not implemented; "
+                "only flavor='seurat' (the default) is currently supported."
+            )
+
         # Calculate gene statistics via SQL using simple aggregation (no JOINs)
         stats_sql = """
         SELECT
diff --git a/tests/test_scanpy.py b/tests/test_scanpy.py
index 89dfff3..bed1af1 100644
--- a/tests/test_scanpy.py
+++ b/tests/test_scanpy.py
@@ -340,6 +340,27 @@ def test_highly_variable_genes_statistics(self, tiny_slaf):
         assert all(hvg_result["variance"] >= 0)
         assert all(hvg_result["dispersion"] >= 0)
 
+    def test_highly_variable_genes_flavor_seurat_default(self, tiny_slaf):
+        """flavor='seurat' (scanpy's default) should be accepted and behave as before"""
+        lazy_adata = LazyAnnData(tiny_slaf)
+
+        hvg_result = pp.highly_variable_genes(
+            lazy_adata, flavor="seurat", inplace=False
+        )
+
+        assert hvg_result is not None
+        assert "highly_variable" in hvg_result.columns
+
+    def test_highly_variable_genes_flavor_unsupported_raises(self, tiny_slaf):
+        """An unsupported flavor should raise a clear NotImplementedError,
+        not a TypeError about an unexpected keyword argument"""
+        lazy_adata = LazyAnnData(tiny_slaf)
+
+        with pytest.raises(NotImplementedError, match="seurat_v3"):
+            pp.highly_variable_genes(
+                lazy_adata, flavor="seurat_v3", n_top_genes=10, inplace=False
+            )
+
     def test_normalize_total_placeholder(self, tiny_slaf):
         """Test normalize_total placeholder implementation"""
         lazy_adata = LazyAnnData(tiny_slaf)