diff --git a/docs/release-notes/4348.fix.md b/docs/release-notes/4348.fix.md new file mode 100644 index 0000000000..66262bfe2b --- /dev/null +++ b/docs/release-notes/4348.fix.md @@ -0,0 +1 @@ +Correct the direction of binary logistic-regression marker scores in `tl.rank_genes_groups` when the reported group is the classifier’s first class {smaller}`D O’Toole` diff --git a/src/scanpy/tools/_rank_genes_groups.py b/src/scanpy/tools/_rank_genes_groups.py index 4f45761e69..1d399c5a65 100644 --- a/src/scanpy/tools/_rank_genes_groups.py +++ b/src/scanpy/tools/_rank_genes_groups.py @@ -597,11 +597,14 @@ def logreg(self, **kwds) -> Generator[_TestResult, None, None]: # not all codes necessarily appear in data existing_codes = np.unique(self.grouping.cat.codes) for igroup, cat in enumerate(self.groups_order): + cat_code: int = np.argmax(self.grouping.cat.categories == cat) if len(self.groups_order) <= 2: # binary logistic regression + # Binary coefficients point toward classes_[1]; orient them + # toward the group under which the scores will be reported. scores = scores_all[0] + if cat_code == clf.classes_[0]: + scores = -scores else: - # cat code is index of cat value in .categories - cat_code: int = np.argmax(self.grouping.cat.categories == cat) # index of scores row is index of cat code in array of existing codes scores_idx: int = np.argmax(existing_codes == cat_code) scores = scores_all[scores_idx] diff --git a/tests/test_rank_genes_groups_logreg.py b/tests/test_rank_genes_groups_logreg.py index 5fde09be3b..ac8857a60b 100644 --- a/tests/test_rank_genes_groups_logreg.py +++ b/tests/test_rank_genes_groups_logreg.py @@ -62,3 +62,66 @@ def test_rank_genes_groups_with_unsorted_groups(): bdata.uns["rank_genes_groups"]["scores"]["Three"] ).to_numpy() np.testing.assert_equal(array_ad, array_bd) + + +@pytest.mark.parametrize("categories", [["A", "unused", "B"], ["B", "unused", "A"]]) +@pytest.mark.parametrize("target", ["A", "B"]) +@pytest.mark.parametrize("representation", ["dense", "sparse", "raw", "layer"]) +def test_binary_logreg_scores_point_toward_requested_group( + categories, target, representation +): + from anndata import AnnData + from scipy.sparse import csr_matrix # noqa: TID251 + + x = np.vstack([ + np.tile([10.0, 0.0, 1.0], (30, 1)), + np.tile([0.0, 10.0, 1.0], (30, 1)), + ]) + adata = AnnData( + x, + obs=pd.DataFrame( + {"group": pd.Categorical(["A"] * 30 + ["B"] * 30, categories=categories)}, + index=[f"cell_{i}" for i in range(60)], + ), + var=pd.DataFrame(index=["A_marker", "B_marker", "shared"]), + ) + kwargs = {"use_raw": False} + if representation == "sparse": + adata.X = csr_matrix(adata.X) + elif representation == "raw": + adata.raw = adata.copy() + adata.X = np.zeros_like(x) + kwargs = {"use_raw": True} + elif representation == "layer": + adata.layers["expression"] = adata.X.copy() + adata.X = np.zeros_like(x) + kwargs["layer"] = "expression" + reference = "B" if target == "A" else "A" + sc.tl.rank_genes_groups( + adata, "group", groups=[target], reference=reference, method="logreg", **kwargs + ) + result = adata.uns["rank_genes_groups"] + names = result["names"][target] + scores = dict(zip(names, result["scores"][target], strict=True)) + assert names[0] == f"{target}_marker" + assert scores[f"{target}_marker"] > 0 + assert scores[f"{reference}_marker"] < 0 + + +@pytest.mark.parametrize("categories", [["A", "B"], ["B", "A"]]) +def test_binary_logreg_default_group_direction(categories): + from anndata import AnnData + + adata = AnnData( + np.vstack([np.tile([10.0, 0.0], (20, 1)), np.tile([0.0, 10.0], (20, 1))]), + obs=pd.DataFrame( + {"group": pd.Categorical(["A"] * 20 + ["B"] * 20, categories=categories)}, + index=[f"cell_{i}" for i in range(40)], + ), + var=pd.DataFrame(index=["A_marker", "B_marker"]), + ) + sc.tl.rank_genes_groups(adata, "group", method="logreg", use_raw=False) + result = adata.uns["rank_genes_groups"] + target = categories[0] + assert result["names"][target][0] == f"{target}_marker" + assert result["scores"][target][0] > 0