From 6d5bc7d0781148ccf4713968007bf63d5b8eda11 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Tue, 30 Jun 2026 16:08:48 +0000 Subject: [PATCH 01/14] Initial plan From 721352e4a557acb69f81e9bfa9465aeece6fdd65 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Tue, 30 Jun 2026 16:14:52 +0000 Subject: [PATCH 02/14] Add run_gsea function using gseapy.prerank for GSEA enrichment analysis --- docs/api_examples/enrichment_analysis.py | 38 ++++++ src/acore/enrichment_analysis/__init__.py | 131 +++++++++++++++++++ tests/enrichment_analysis/test_enrichment.py | 75 +++++++++++ 3 files changed, 244 insertions(+) diff --git a/docs/api_examples/enrichment_analysis.py b/docs/api_examples/enrichment_analysis.py index cce99e5d..eca4044f 100644 --- a/docs/api_examples/enrichment_analysis.py +++ b/docs/api_examples/enrichment_analysis.py @@ -357,6 +357,44 @@ figure = go.Figure(data=figure["data"], layout=figure["layout"]) figure +# %% [markdown] +# # GSEA (Gene Set Enrichment Analysis) +# GSEA ranks all measured features (e.g. proteins) by a continuous metric such as the +# log2 fold-change from differential expression analysis and tests whether the members of +# a gene set are enriched at the top or bottom of this ranked list. +# +# Unlike overrepresentation analysis (ORA), GSEA does **not** require a significance +# cutoff to define a foreground set: it uses the full ranked list and is therefore less +# sensitive to the choice of an arbitrary threshold. +# +# Here we use `gseapy.prerank` under the hood, which implements the classic weighted +# Kolmogorov–Smirnov enrichment statistic described in Subramanian et al. (2005). +# +# See [`gseapy.prerank`](https://gseapy.readthedocs.io/en/latest/run.html#gseapy.prerank) +# for more details. + +# %% +gsea_results = acore.enrichment_analysis.run_gsea( + regulation_data=diff_reg, + annotation=annotations, + identifier="identifier", + ranking_col="log2FC", + annotation_col="annotation", + identifier_col="identifier", + min_size=1, + permutation_num=100, + seed=42, +) +gsea_results + +# %% [markdown] +# The normalised enrichment score (NES) summarises the enrichment: positive values +# indicate that the gene set is enriched among the **up**-regulated features and negative +# values indicate enrichment among the **down**-regulated features. + +# %% +ax = gsea_results["NES"].plot.hist(title="Distribution of NES", xlabel="NES") + # %% [markdown] # # Compare two distributions - KS test # diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index b3ae590c..5e18b41b 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -34,6 +34,7 @@ __all__ = [ "run_site_regulation_enrichment", "run_up_down_regulation_enrichment", + "run_gsea", "run_fisher", "run_kolmogorov_smirnov", ] @@ -594,3 +595,133 @@ def run_ssgsea( # potentially return wide format in separate format # result = {"es": enrichment_es, "nes": enrichment_nes} return result + + +def run_gsea( + regulation_data: pd.DataFrame, + annotation: pd.DataFrame, + identifier: str = "identifier", + ranking_col: str = "log2FC", + annotation_col: str = "annotation", + identifier_col: str = "identifier", + min_size: int = 15, + max_size: int = 500, + permutation_num: int = 1000, + seed: int = 42, + weight: float = 1.0, + ascending: bool = False, +) -> pd.DataFrame: + """ + Run Gene Set Enrichment Analysis (GSEA) on a pre-ranked gene list using + `gseapy.prerank `_. + + Genes are ranked by `ranking_col` (e.g. log2 fold-change from differential + expression analysis) and enrichment is computed for each gene set defined in + `annotation`. Unlike overrepresentation analysis (ORA), GSEA does not require + a significance cutoff to split genes into foreground/background; it uses the + full ranked list instead. + + :param pd.DataFrame regulation_data: DataFrame with differential regulation + results. Must contain the column specified by `identifier` and the + column specified by `ranking_col`. Each identifier must appear at most + once (duplicates are not allowed). + :param pd.DataFrame annotation: DataFrame with functional annotations. + Must contain `annotation_col` and `identifier_col`. + :param str identifier: Column in `regulation_data` containing feature + identifiers. These must match the identifiers in `annotation_col` of + the `annotation` DataFrame. + :param str ranking_col: Column in `regulation_data` used to rank the + features (e.g. ``'log2FC'``). Features are ranked in descending order + by default so that strongly up-regulated features appear at the top. + :param str annotation_col: Column in `annotation` containing annotation + term names (e.g. pathway names). + :param str identifier_col: Column in `annotation` containing feature + identifiers that map to `identifier` in `regulation_data`. + :param int min_size: Minimum number of features in a gene set that are also + present in the ranking. Gene sets smaller than this threshold are + excluded. + :param int max_size: Maximum number of features in a gene set that are also + present in the ranking. Gene sets larger than this threshold are + excluded. + :param int permutation_num: Number of permutations used to estimate nominal + p-values. Higher values give more precise p-values but increase runtime. + Default: 1000. + :param int seed: Random seed for reproducibility. Default: 42. + :param float weight: Exponent for the enrichment score weighting. + ``weight=1`` (default) corresponds to weighted GSEA; ``weight=0`` + gives classic (unweighted) GSEA. + :param bool ascending: Sort order for the ranking. ``False`` (default) + means descending (highest ranking metric first, suited for log2FC). + :return: DataFrame with one row per gene set. Columns are: ``'Term'``, + ``'ES'`` (enrichment score), ``'NES'`` (normalised enrichment score), + ``'NOM p-val'`` (nominal p-value), ``'FDR q-val'`` (FDR-corrected + q-value), ``'FWER p-val'``, ``'Tag %'``, ``'Gene %'``, + ``'Lead_genes'`` (semi-colon-separated leading-edge genes). + :rtype: pd.DataFrame + + Example:: + + result = run_gsea( + regulation_data=diff_reg, + annotation=annotations, + identifier='identifier', + ranking_col='log2FC', + annotation_col='annotation', + identifier_col='identifier', + min_size=5, + permutation_num=100, + ) + """ + if regulation_data is None or regulation_data.empty: + raise ValueError("regulation_data is empty.") + + if identifier not in regulation_data.columns: + raise ValueError( + f"Column '{identifier}' not found in regulation_data. " + "Specify the correct column name via `identifier`." + ) + if ranking_col not in regulation_data.columns: + raise ValueError( + f"Column '{ranking_col}' not found in regulation_data. " + "Specify the correct column name via `ranking_col`." + ) + if not regulation_data[identifier].is_unique: + raise ValueError( + f"Column '{identifier}' in regulation_data contains duplicate values. " + "Each identifier must appear at most once." + ) + if annotation_col not in annotation.columns: + raise ValueError( + f"Missing annotation column '{annotation_col}' in annotation DataFrame. " + "Specify the correct column name via `annotation_col`." + ) + if identifier_col not in annotation.columns: + raise ValueError( + f"Missing identifier column '{identifier_col}' in annotation DataFrame. " + "Specify the correct column name via `identifier_col`." + ) + + # Build the pre-ranked Series: index = feature identifier, value = ranking metric + rnk = regulation_data.set_index(identifier)[ranking_col].dropna() + + # Build gene sets dict from the annotation DataFrame + gene_sets = ( + annotation.groupby(annotation_col)[identifier_col] + .apply(list) + .to_dict() + ) + + enrichment = gp.prerank( + rnk=rnk, + gene_sets=gene_sets, + min_size=min_size, + max_size=max_size, + permutation_num=permutation_num, + weight=weight, + ascending=ascending, + seed=seed, + no_plot=True, + outdir=None, + ) + result = pd.DataFrame(enrichment.res2d).drop(columns=["Name"]) + return result diff --git a/tests/enrichment_analysis/test_enrichment.py b/tests/enrichment_analysis/test_enrichment.py index 90afd7c7..dd7d0b5e 100644 --- a/tests/enrichment_analysis/test_enrichment.py +++ b/tests/enrichment_analysis/test_enrichment.py @@ -588,3 +588,78 @@ def test_run_regulation_enrichment_with_duplicates(): annotation=annotation, min_detected_in_set=1, ) + + +def test_run_gsea_basic(): + """Integration test for run_gsea: checks that the function returns + a DataFrame with the expected columns and one row per gene set.""" + regulation_data = pd.DataFrame( + { + "identifier": ["gene1", "gene2", "gene3", "gene4", "gene5", "gene6"], + "log2FC": [3.0, 2.5, 1.0, -1.0, -2.5, -3.0], + } + ) + annotation = pd.DataFrame( + { + "annotation": ["path1", "path1", "path1", "path2", "path2", "path2"], + "identifier": ["gene1", "gene2", "gene3", "gene4", "gene5", "gene6"], + } + ) + + result = ea.run_gsea( + regulation_data=regulation_data, + annotation=annotation, + identifier="identifier", + ranking_col="log2FC", + annotation_col="annotation", + identifier_col="identifier", + min_size=2, + permutation_num=10, + seed=42, + ) + + assert isinstance(result, pd.DataFrame) + assert not result.empty + expected_cols = {"Term", "ES", "NES", "NOM p-val", "FDR q-val", "Lead_genes"} + assert expected_cols.issubset(result.columns) + assert set(result["Term"]) == {"path1", "path2"} + # path1 genes are up-regulated -> positive NES + nes_path1 = result.loc[result["Term"] == "path1", "NES"].iloc[0] + nes_path2 = result.loc[result["Term"] == "path2", "NES"].iloc[0] + assert nes_path1 > 0, "path1 (up-regulated genes) should have a positive NES" + assert nes_path2 < 0, "path2 (down-regulated genes) should have a negative NES" + + +def test_run_gsea_raises_on_empty_data(): + """run_gsea should raise ValueError when regulation_data is empty.""" + with pytest.raises(ValueError, match="regulation_data is empty"): + ea.run_gsea( + regulation_data=pd.DataFrame(), + annotation=pd.DataFrame({"annotation": [], "identifier": []}), + ) + + +def test_run_gsea_raises_on_missing_identifier_col(): + """run_gsea should raise ValueError when identifier column is missing.""" + regulation_data = pd.DataFrame({"gene": ["g1"], "log2FC": [1.0]}) + annotation = pd.DataFrame({"annotation": ["path1"], "identifier": ["g1"]}) + with pytest.raises(ValueError, match="identifier"): + ea.run_gsea( + regulation_data=regulation_data, + annotation=annotation, + identifier="identifier", # not a column in regulation_data + ) + + +def test_run_gsea_raises_on_duplicate_identifiers(): + """run_gsea should raise ValueError when identifier column has duplicates.""" + regulation_data = pd.DataFrame( + {"identifier": ["gene1", "gene1"], "log2FC": [1.0, 2.0]} + ) + annotation = pd.DataFrame({"annotation": ["path1"], "identifier": ["gene1"]}) + with pytest.raises(ValueError, match="duplicate"): + ea.run_gsea( + regulation_data=regulation_data, + annotation=annotation, + ) + From b2f2ce3c25f252822bc7c5513acbc154df46d6b1 Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Fri, 2 Oct 2026 10:54:00 +0200 Subject: [PATCH 03/14] :art: update warnings and information for enrichment analysis --- src/acore/enrichment_analysis/__init__.py | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index 5e18b41b..3014798f 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -229,6 +229,12 @@ def run_up_down_regulation_enrichment( correction=correction, correction_alpha=correction_alpha, ) + if _enrichment.empty: + logger.warning( + "No significant enrichment found for comparison %s (%s).", + comparison_tag, + direction, + ) _enrichment["direction"] = direction _enrichment["comparison"] = comparison_tag ret.append(_enrichment) @@ -705,11 +711,7 @@ def run_gsea( rnk = regulation_data.set_index(identifier)[ranking_col].dropna() # Build gene sets dict from the annotation DataFrame - gene_sets = ( - annotation.groupby(annotation_col)[identifier_col] - .apply(list) - .to_dict() - ) + gene_sets = annotation.groupby(annotation_col)[identifier_col].apply(list).to_dict() enrichment = gp.prerank( rnk=rnk, From 8121ff2e1574683ca66c18a6e7bc6b5be09d4dc7 Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Fri, 2 Oct 2026 13:45:26 +0200 Subject: [PATCH 04/14] :memo: update checklist hint --- CONTRIBUTING.md | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md index 1c1f6ff1..78325347 100644 --- a/CONTRIBUTING.md +++ b/CONTRIBUTING.md @@ -118,8 +118,7 @@ Before you submit a pull request, check that it meets these guidelines: into a function with a docstring, and add the feature to the list in README.rst. 3. The pull request should pass the GitHub workflows. -See the PR template example: -[Add module PR template](https://github.com/Multiomics-Analytics-Group/acore/blob/main/.github/PULL_REQUEST_TEMPLATE/module.md) +See the checklist if you open a PR! ## Deploying From 2ffcad1ae18dc24aa5259fd35fb756df243ffd89 Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" <41898282+github-actions[bot]@users.noreply.github.com> Date: Fri, 2 Oct 2026 11:50:28 +0000 Subject: [PATCH 05/14] docs: update markdown reference --- .../markdown_ref/acore.enrichment_analysis.md | 147 +++++++++++++----- docs/markdown_ref/acore.md | 5 +- docs/markdown_ref/modules.md | 5 +- 3 files changed, 113 insertions(+), 44 deletions(-) diff --git a/docs/markdown_ref/acore.enrichment_analysis.md b/docs/markdown_ref/acore.enrichment_analysis.md index 245cfc99..c937f287 100644 --- a/docs/markdown_ref/acore.enrichment_analysis.md +++ b/docs/markdown_ref/acore.enrichment_analysis.md @@ -6,46 +6,6 @@ analysis. Most things in this module are covered in [https://www.youtube.com/watch?v=2NC1QOXmc5o](https://www.youtube.com/watch?v=2NC1QOXmc5o) by Lars Juhl Jensen. -### run_fisher(group1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], group2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] - -Run fisher’s exact test on two groups using [scipy.stats.fisher_exact](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.fisher_exact.html). - -Example: - -```default -# annotated not-annotated -# group1 a b -# group2 c d - - -odds, pvalue = stats.fisher_exact(group1=[a, b], - group2 =[c, d] - ) -``` - -### run_kolmogorov_smirnov(dist1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[float](https://docs.python.org/3/builtins/functions.html#float)], dist2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[float](https://docs.python.org/3/builtins/functions.html#float)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] - -Compute the Kolmogorov-Smirnov statistic on 2 samples. -See [scipy.stats.ks_2samp](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.ks_2samp.html) - -* **Parameters:** - * **dist1** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list)) – sequence of 1-D ndarray (first distribution to compare) - drawn from a continuous distribution - * **dist2** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list)) – sequence of 1-D ndarray (second distribution to compare) - drawn from a continuous distribution - * **alternative** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – defines the alternative hypothesis (default is ‘two-sided’): - \* **‘two-sided’** - \* **‘less’** - \* **‘greater’** -* **Returns:** - statistic float and KS statistic pvalue float Two-tailed p-value. - -Example: - -```default -result = run_kolmogorov_smirnov(dist1, dist2, alternative='two-sided') -``` - ### run_site_regulation_enrichment(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', groups: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[str](https://docs.python.org/3/builtins/stdtypes.html#str)] = ('group1', 'group2'), annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', rejected_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'rejected', group_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'group', method: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fisher', regex: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = '(\\\\w+~.+)_\\\\w\\\\d+\\\\-\\\\w+', correction: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fdr_bh', remove_duplicates: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False) → [DataFrame](https://pandera.readthedocs.io/en/stable/reference/generated/pandera.typing.DataFrame.html#pandera.typing.DataFrame)[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] This function runs a simple enrichment analysis for significantly @@ -139,6 +99,113 @@ result = run_up_down_regulation_enrichment( ) ``` +### run_gsea(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', ranking_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'log2FC', annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', identifier_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', min_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 15, max_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 500, permutation_num: [int](https://docs.python.org/3/builtins/functions.html#int) = 1000, seed: [int](https://docs.python.org/3/builtins/functions.html#int) = 42, weight: [float](https://docs.python.org/3/builtins/functions.html#float) = 1.0, ascending: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False) → DataFrame + +Run Gene Set Enrichment Analysis (GSEA) on a pre-ranked gene list using +[gseapy.prerank](https://gseapy.readthedocs.io/en/latest/run.html#gseapy.prerank). + +Genes are ranked by ranking_col (e.g. log2 fold-change from differential +expression analysis) and enrichment is computed for each gene set defined in +annotation. Unlike overrepresentation analysis (ORA), GSEA does not require +a significance cutoff to split genes into foreground/background; it uses the +full ranked list instead. + +* **Parameters:** + * **regulation_data** (*pd.DataFrame*) – DataFrame with differential regulation + results. Must contain the column specified by identifier and the + column specified by ranking_col. Each identifier must appear at most + once (duplicates are not allowed). + * **annotation** (*pd.DataFrame*) – DataFrame with functional annotations. + Must contain annotation_col and identifier_col. + * **identifier** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – Column in regulation_data containing feature + identifiers. These must match the identifiers in annotation_col of + the annotation DataFrame. + * **ranking_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – Column in regulation_data used to rank the + features (e.g. `'log2FC'`). Features are ranked in descending order + by default so that strongly up-regulated features appear at the top. + * **annotation_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – Column in annotation containing annotation + term names (e.g. pathway names). + * **identifier_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – Column in annotation containing feature + identifiers that map to identifier in regulation_data. + * **min_size** ([*int*](https://docs.python.org/3/builtins/functions.html#int)) – Minimum number of features in a gene set that are also + present in the ranking. Gene sets smaller than this threshold are + excluded. + * **max_size** ([*int*](https://docs.python.org/3/builtins/functions.html#int)) – Maximum number of features in a gene set that are also + present in the ranking. Gene sets larger than this threshold are + excluded. + * **permutation_num** ([*int*](https://docs.python.org/3/builtins/functions.html#int)) – Number of permutations used to estimate nominal + p-values. Higher values give more precise p-values but increase runtime. + Default: 1000. + * **seed** ([*int*](https://docs.python.org/3/builtins/functions.html#int)) – Random seed for reproducibility. Default: 42. + * **weight** ([*float*](https://docs.python.org/3/builtins/functions.html#float)) – Exponent for the enrichment score weighting. + `weight=1` (default) corresponds to weighted GSEA; `weight=0` + gives classic (unweighted) GSEA. + * **ascending** ([*bool*](https://docs.python.org/3/builtins/functions.html#bool)) – Sort order for the ranking. `False` (default) + means descending (highest ranking metric first, suited for log2FC). +* **Returns:** + DataFrame with one row per gene set. Columns are: `'Term'`, + `'ES'` (enrichment score), `'NES'` (normalised enrichment score), + `'NOM p-val'` (nominal p-value), `'FDR q-val'` (FDR-corrected + q-value), `'FWER p-val'`, `'Tag %'`, `'Gene %'`, + `'Lead_genes'` (semi-colon-separated leading-edge genes). +* **Return type:** + pd.DataFrame + +Example: + +```default +result = run_gsea( + regulation_data=diff_reg, + annotation=annotations, + identifier='identifier', + ranking_col='log2FC', + annotation_col='annotation', + identifier_col='identifier', + min_size=5, + permutation_num=100, +) +``` + +### run_fisher(group1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], group2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] + +Run fisher’s exact test on two groups using [scipy.stats.fisher_exact](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.fisher_exact.html). + +Example: + +```default +# annotated not-annotated +# group1 a b +# group2 c d + + +odds, pvalue = stats.fisher_exact(group1=[a, b], + group2 =[c, d] + ) +``` + +### run_kolmogorov_smirnov(dist1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[float](https://docs.python.org/3/builtins/functions.html#float)], dist2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[float](https://docs.python.org/3/builtins/functions.html#float)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] + +Compute the Kolmogorov-Smirnov statistic on 2 samples. +See [scipy.stats.ks_2samp](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.ks_2samp.html) + +* **Parameters:** + * **dist1** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list)) – sequence of 1-D ndarray (first distribution to compare) + drawn from a continuous distribution + * **dist2** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list)) – sequence of 1-D ndarray (second distribution to compare) + drawn from a continuous distribution + * **alternative** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – defines the alternative hypothesis (default is ‘two-sided’): + \* **‘two-sided’** + \* **‘less’** + \* **‘greater’** +* **Returns:** + statistic float and KS statistic pvalue float Two-tailed p-value. + +Example: + +```default +result = run_kolmogorov_smirnov(dist1, dist2, alternative='two-sided') +``` + ## Subpackages * [acore.enrichment_analysis.statistical_tests namespace](acore.enrichment_analysis.statistical_tests.md) diff --git a/docs/markdown_ref/acore.md b/docs/markdown_ref/acore.md index 30a7098e..4a286b52 100644 --- a/docs/markdown_ref/acore.md +++ b/docs/markdown_ref/acore.md @@ -67,10 +67,11 @@ * [`qc_rlsc_loess()`](acore.drift_correction.md#acore.drift_correction.loess_drift_correction.qc_rlsc_loess) * [`run_loess_drift_correction()`](acore.drift_correction.md#acore.drift_correction.loess_drift_correction.run_loess_drift_correction) * [acore.enrichment_analysis package](acore.enrichment_analysis.md) - * [`run_fisher()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_fisher) - * [`run_kolmogorov_smirnov()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_kolmogorov_smirnov) * [`run_site_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_site_regulation_enrichment) * [`run_up_down_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_up_down_regulation_enrichment) + * [`run_gsea()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_gsea) + * [`run_fisher()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_fisher) + * [`run_kolmogorov_smirnov()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_kolmogorov_smirnov) * [Subpackages](acore.enrichment_analysis.md#subpackages) * [acore.enrichment_analysis.statistical_tests namespace](acore.enrichment_analysis.statistical_tests.md) * [Submodules](acore.enrichment_analysis.statistical_tests.md#submodules) diff --git a/docs/markdown_ref/modules.md b/docs/markdown_ref/modules.md index e128c154..77317eb1 100644 --- a/docs/markdown_ref/modules.md +++ b/docs/markdown_ref/modules.md @@ -40,10 +40,11 @@ * [acore.drift_correction.cpca_drift_correction module](acore.drift_correction.md#module-acore.drift_correction.cpca_drift_correction) * [acore.drift_correction.loess_drift_correction module](acore.drift_correction.md#module-acore.drift_correction.loess_drift_correction) * [acore.enrichment_analysis package](acore.enrichment_analysis.md) - * [`run_fisher()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_fisher) - * [`run_kolmogorov_smirnov()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_kolmogorov_smirnov) * [`run_site_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_site_regulation_enrichment) * [`run_up_down_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_up_down_regulation_enrichment) + * [`run_gsea()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_gsea) + * [`run_fisher()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_fisher) + * [`run_kolmogorov_smirnov()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_kolmogorov_smirnov) * [Subpackages](acore.enrichment_analysis.md#subpackages) * [Submodules](acore.enrichment_analysis.md#submodules) * [acore.enrichment_analysis.annotate module](acore.enrichment_analysis.md#module-acore.enrichment_analysis.annotate) From 79a9b24e464f7a01b85ac0a85383aabe55477915 Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Fri, 2 Oct 2026 16:08:39 +0200 Subject: [PATCH 06/14] :art: format imports to be updated with a plot from vuecore --- .../enrichment_analysis.ipynb | 68 ++++++++++++++++++- .../enrichment_analysis.py | 8 +-- 2 files changed, 69 insertions(+), 7 deletions(-) diff --git a/docs/examples_single_analysis/enrichment_analysis.ipynb b/docs/examples_single_analysis/enrichment_analysis.ipynb index cdbb0326..2a47609b 100644 --- a/docs/examples_single_analysis/enrichment_analysis.ipynb +++ b/docs/examples_single_analysis/enrichment_analysis.ipynb @@ -47,6 +47,8 @@ "\n", "import dsp_pandas\n", "import pandas as pd\n", + "import plotly.graph_objects as go\n", + "from vuecore import viz\n", "\n", "import acore\n", "import acore.differential_regulation\n", @@ -646,9 +648,6 @@ "metadata": {}, "outputs": [], "source": [ - "import plotly.graph_objects as go\n", - "from vuecore import viz\n", - "\n", "args = {\"factor\": 2, \"loadings\": 1} # increase number of loadings or scaling factor\n", "#! pca_results has three items, but docstring requests only two -> double check\n", "figure = viz.get_pca_plot(data=pca_result, identifier=\"PCA enrichment\", args=args)\n", @@ -656,6 +655,69 @@ "figure" ] }, + { + "cell_type": "markdown", + "id": "eb6bd6cd", + "metadata": {}, + "source": [ + "# GSEA (Gene Set Enrichment Analysis)\n", + "GSEA ranks all measured features (e.g. proteins) by a continuous metric such as the\n", + "log2 fold-change from differential expression analysis and tests whether the members of\n", + "a gene set are enriched at the top or bottom of this ranked list.\n", + "\n", + "Unlike overrepresentation analysis (ORA), GSEA does **not** require a significance\n", + "cutoff to define a foreground set: it uses the full ranked list and is therefore less\n", + "sensitive to the choice of an arbitrary threshold.\n", + "\n", + "Here we use `gseapy.prerank` under the hood, which implements the classic weighted\n", + "Kolmogorov–Smirnov enrichment statistic described in\n", + "[Subramanian et al. (2005)](https://doi.org/10.1073/pnas.0506580102).\n", + "\n", + "See [`gseapy.prerank`](https://gseapy.readthedocs.io/en/latest/run.html#gseapy.prerank)\n", + "for more details." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "ef98e5b1", + "metadata": {}, + "outputs": [], + "source": [ + "gsea_results = acore.enrichment_analysis.run_gsea(\n", + " regulation_data=diff_reg,\n", + " annotation=annotations,\n", + " identifier=\"identifier\",\n", + " ranking_col=\"log2FC\",\n", + " annotation_col=\"annotation\",\n", + " identifier_col=\"identifier\",\n", + " min_size=1,\n", + " permutation_num=100,\n", + " seed=42,\n", + ")\n", + "gsea_results" + ] + }, + { + "cell_type": "markdown", + "id": "f5f96579", + "metadata": {}, + "source": [ + "The normalised enrichment score (NES) summarises the enrichment: positive values\n", + "indicate that the gene set is enriched among the **up**-regulated features and negative\n", + "values indicate enrichment among the **down**-regulated features." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "fd02b00c", + "metadata": {}, + "outputs": [], + "source": [ + "ax = gsea_results[\"NES\"].plot.hist(title=\"Distribution of NES\", xlabel=\"NES\")" + ] + }, { "cell_type": "markdown", "id": "4be5a8c8", diff --git a/docs/examples_single_analysis/enrichment_analysis.py b/docs/examples_single_analysis/enrichment_analysis.py index 2d574686..70c3162a 100644 --- a/docs/examples_single_analysis/enrichment_analysis.py +++ b/docs/examples_single_analysis/enrichment_analysis.py @@ -22,6 +22,8 @@ import dsp_pandas import pandas as pd +import plotly.graph_objects as go +from vuecore import viz import acore import acore.differential_regulation @@ -350,9 +352,6 @@ # for this, which is also developed by the Multiomics Analytics Group. # %% -import plotly.graph_objects as go -from vuecore import viz - args = {"factor": 2, "loadings": 1} # increase number of loadings or scaling factor # #! pca_results has three items, but docstring requests only two -> double check figure = viz.get_pca_plot(data=pca_result, identifier="PCA enrichment", args=args) @@ -370,7 +369,8 @@ # sensitive to the choice of an arbitrary threshold. # # Here we use `gseapy.prerank` under the hood, which implements the classic weighted -# Kolmogorov–Smirnov enrichment statistic described in Subramanian et al. (2005). +# Kolmogorov–Smirnov enrichment statistic described in +# [Subramanian et al. (2005)](https://doi.org/10.1073/pnas.0506580102). # # See [`gseapy.prerank`](https://gseapy.readthedocs.io/en/latest/run.html#gseapy.prerank) # for more details. From 47aefc61876f72efe79192d9cff2bd9f52979d48 Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Mon, 19 Oct 2026 16:30:32 +0200 Subject: [PATCH 07/14] :art: format and remove missing imports from all --- src/acore/enrichment_analysis/__init__.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index a8e81137..ede8df78 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -33,12 +33,10 @@ __all__ = [ "run_fisher", + "run_gsea", "run_kolmogorov_smirnov", "run_site_regulation_enrichment", "run_up_down_regulation_enrichment", - "run_gsea", - "run_fisher", - "run_kolmogorov_smirnov", ] From 211962fd2b31140f4d0de47d2073f2c52c049efd Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" <41898282+github-actions[bot]@users.noreply.github.com> Date: Mon, 19 Oct 2026 14:32:31 +0000 Subject: [PATCH 08/14] docs: update markdown reference --- .../markdown_ref/acore.enrichment_analysis.md | 202 +++++++++--------- docs/markdown_ref/acore.md | 6 +- docs/markdown_ref/acore.types.md | 12 +- docs/markdown_ref/modules.md | 6 +- 4 files changed, 113 insertions(+), 113 deletions(-) diff --git a/docs/markdown_ref/acore.enrichment_analysis.md b/docs/markdown_ref/acore.enrichment_analysis.md index c937f287..dc7b07ad 100644 --- a/docs/markdown_ref/acore.enrichment_analysis.md +++ b/docs/markdown_ref/acore.enrichment_analysis.md @@ -6,97 +6,21 @@ analysis. Most things in this module are covered in [https://www.youtube.com/watch?v=2NC1QOXmc5o](https://www.youtube.com/watch?v=2NC1QOXmc5o) by Lars Juhl Jensen. -### run_site_regulation_enrichment(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', groups: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[str](https://docs.python.org/3/builtins/stdtypes.html#str)] = ('group1', 'group2'), annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', rejected_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'rejected', group_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'group', method: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fisher', regex: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = '(\\\\w+~.+)_\\\\w\\\\d+\\\\-\\\\w+', correction: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fdr_bh', remove_duplicates: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False) → [DataFrame](https://pandera.readthedocs.io/en/stable/reference/generated/pandera.typing.DataFrame.html#pandera.typing.DataFrame)[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] - -This function runs a simple enrichment analysis for significantly -regulated protein sites in a dataset. +### run_fisher(group1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], group2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] -* **Parameters:** - * **regulation_data** – pandas.DataFrame resulting from differential - regulation analysis. - * **annotation** – pandas.DataFrame with annotations for features - (columns: ‘annotation’, ‘identifier’ (feature identifiers), and ‘source’). - * **identifier** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing - feature identifiers. - * **groups** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list)) – column names from regulation_data containing - group identifiers. - * **annotation_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing - annotation terms. - * **rejected_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from regulation_data containing - boolean for rejected null hypothesis. - * **group_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – column name for new column in annotation dataframe - determining if feature belongs to foreground or background. - * **method** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – method used to compute enrichment - (only ‘fisher’ is supported currently). - * **regex** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – how to extract the annotated identifier from the site identifier -* **Returns:** - pandas.DataFrame with columns: ‘terms’, ‘identifiers’, ‘foreground’, - ‘background’, foreground_pop, background_pop, ‘pvalue’, ‘padj’ and ‘rejected’. -* **Raises:** - [**ValueError**](https://docs.python.org/3/builtins/exceptions.html#ValueError) – if regulation_data is None or empty. +Run fisher’s exact test on two groups using [scipy.stats.fisher_exact](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.fisher_exact.html). Example: ```default -result = run_site_regulation_enrichment(regulation_data, - annotation, - identifier='identifier', - groups=['group1', 'group2'], - annotation_col='annotation', - rejected_col='rejected', - group_col='group', - method='fisher', - match="(\\w+~.+)_\\w\\d+\\-\\w+" -) -``` - -### run_up_down_regulation_enrichment(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', groups: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[str](https://docs.python.org/3/builtins/stdtypes.html#str)] = ('group1', 'group2'), annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', pval_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'pval', group_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'group', log2fc_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'log2FC', method: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fisher', min_detected_in_set: [int](https://docs.python.org/3/builtins/functions.html#int) = 2, correction: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fdr_bh', correction_alpha: [float](https://docs.python.org/3/builtins/functions.html#float) = 0.05, lfc_cutoff: [float](https://docs.python.org/3/builtins/functions.html#float) = 1) → [DataFrame](https://pandera.readthedocs.io/en/stable/reference/generated/pandera.typing.DataFrame.html#pandera.typing.DataFrame)[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] - -This function runs a simple enrichment analysis for significantly regulated proteins -distinguishing between upregulation per group defined by the groups. - -* **Parameters:** - * **regulation_data** ([*pandas.DataFrame*](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.html#pandas.DataFrame)) – pandas.DataFrame resulting from differential regulation - analysis (CKG’s regulation table). - * **annotation** ([*pandas.DataFrame*](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.html#pandas.DataFrame)) – pandas.DataFrame with annotations for features - (columns: ‘annotation’, ‘identifier’ (feature identifiers), and ‘source’). - * **identifier** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing feature identifiers. - * **groups** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list) *[*[*str*](https://docs.python.org/3/builtins/stdtypes.html#str) *]*) – - - column names from regulation_data containing group identifiers. - See [pandas.DataFrame.groupby](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.groupby.html) for more information. - * **annotation_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing annotation terms. - * **rejected_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from regulation_data containing boolean for - rejected null hypothesis. - * **group_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – column name for new column in annotation dataframe determining - if feature belongs to foreground or background. - * **method** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – method used to compute enrichment - (only ‘fisher’ is supported currently). - * **correction** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – method to be used for multiple-testing correction - * **alpha** ([*float*](https://docs.python.org/3/builtins/functions.html#float)) – adjusted p-value cutoff to define significance - * **lfc_cutoff** ([*float*](https://docs.python.org/3/builtins/functions.html#float)) – log fold-change cutoff to define practical significance -* **Returns:** - DataFrame adhering to EnrichmentAnalysisSchema -* **Return type:** - DataFrame[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] +# annotated not-annotated +# group1 a b +# group2 c d -Example: -```default -result = run_up_down_regulation_enrichment( - regulation_data, - annotation, - identifier='identifier', - groups=['group1', - 'group2'], - annotation_col='annotation', - rejected_col='rejected', - group_col='group', - method='fisher', - correction='fdr_bh', - alpha=0.05, - lfc_cutoff=1, -) +odds, pvalue = stats.fisher_exact(group1=[a, b], + group2 =[c, d] + ) ``` ### run_gsea(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', ranking_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'log2FC', annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', identifier_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', min_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 15, max_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 500, permutation_num: [int](https://docs.python.org/3/builtins/functions.html#int) = 1000, seed: [int](https://docs.python.org/3/builtins/functions.html#int) = 42, weight: [float](https://docs.python.org/3/builtins/functions.html#float) = 1.0, ascending: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False) → DataFrame @@ -166,23 +90,6 @@ result = run_gsea( ) ``` -### run_fisher(group1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], group2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[int](https://docs.python.org/3/builtins/functions.html#int)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] - -Run fisher’s exact test on two groups using [scipy.stats.fisher_exact](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.fisher_exact.html). - -Example: - -```default -# annotated not-annotated -# group1 a b -# group2 c d - - -odds, pvalue = stats.fisher_exact(group1=[a, b], - group2 =[c, d] - ) -``` - ### run_kolmogorov_smirnov(dist1: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[float](https://docs.python.org/3/builtins/functions.html#float)], dist2: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[float](https://docs.python.org/3/builtins/functions.html#float)], alternative: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'two-sided') → [tuple](https://docs.python.org/3/builtins/stdtypes.html#tuple)[[float](https://docs.python.org/3/builtins/functions.html#float), [float](https://docs.python.org/3/builtins/functions.html#float)] Compute the Kolmogorov-Smirnov statistic on 2 samples. @@ -206,6 +113,99 @@ Example: result = run_kolmogorov_smirnov(dist1, dist2, alternative='two-sided') ``` +### run_site_regulation_enrichment(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', groups: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[str](https://docs.python.org/3/builtins/stdtypes.html#str)] = ('group1', 'group2'), annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', rejected_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'rejected', group_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'group', method: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fisher', regex: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = '(\\\\w+~.+)_\\\\w\\\\d+\\\\-\\\\w+', correction: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fdr_bh', remove_duplicates: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False) → [DataFrame](https://pandera.readthedocs.io/en/stable/reference/generated/pandera.typing.DataFrame.html#pandera.typing.DataFrame)[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] + +This function runs a simple enrichment analysis for significantly +regulated protein sites in a dataset. + +* **Parameters:** + * **regulation_data** – pandas.DataFrame resulting from differential + regulation analysis. + * **annotation** – pandas.DataFrame with annotations for features + (columns: ‘annotation’, ‘identifier’ (feature identifiers), and ‘source’). + * **identifier** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing + feature identifiers. + * **groups** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list)) – column names from regulation_data containing + group identifiers. + * **annotation_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing + annotation terms. + * **rejected_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from regulation_data containing + boolean for rejected null hypothesis. + * **group_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – column name for new column in annotation dataframe + determining if feature belongs to foreground or background. + * **method** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – method used to compute enrichment + (only ‘fisher’ is supported currently). + * **regex** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – how to extract the annotated identifier from the site identifier +* **Returns:** + pandas.DataFrame with columns: ‘terms’, ‘identifiers’, ‘foreground’, + ‘background’, foreground_pop, background_pop, ‘pvalue’, ‘padj’ and ‘rejected’. +* **Raises:** + [**ValueError**](https://docs.python.org/3/builtins/exceptions.html#ValueError) – if regulation_data is None or empty. + +Example: + +```default +result = run_site_regulation_enrichment(regulation_data, + annotation, + identifier='identifier', + groups=['group1', 'group2'], + annotation_col='annotation', + rejected_col='rejected', + group_col='group', + method='fisher', + match="(\\w+~.+)_\\w\\d+\\-\\w+" +) +``` + +### run_up_down_regulation_enrichment(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', groups: [list](https://docs.python.org/3/builtins/stdtypes.html#list)[[str](https://docs.python.org/3/builtins/stdtypes.html#str)] = ('group1', 'group2'), annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', pval_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'pval', group_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'group', log2fc_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'log2FC', method: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fisher', min_detected_in_set: [int](https://docs.python.org/3/builtins/functions.html#int) = 2, correction: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'fdr_bh', correction_alpha: [float](https://docs.python.org/3/builtins/functions.html#float) = 0.05, lfc_cutoff: [float](https://docs.python.org/3/builtins/functions.html#float) = 1) → [DataFrame](https://pandera.readthedocs.io/en/stable/reference/generated/pandera.typing.DataFrame.html#pandera.typing.DataFrame)[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] + +This function runs a simple enrichment analysis for significantly regulated proteins +distinguishing between upregulation per group defined by the groups. + +* **Parameters:** + * **regulation_data** ([*pandas.DataFrame*](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.html#pandas.DataFrame)) – pandas.DataFrame resulting from differential regulation + analysis (CKG’s regulation table). + * **annotation** ([*pandas.DataFrame*](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.html#pandas.DataFrame)) – pandas.DataFrame with annotations for features + (columns: ‘annotation’, ‘identifier’ (feature identifiers), and ‘source’). + * **identifier** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing feature identifiers. + * **groups** ([*list*](https://docs.python.org/3/builtins/stdtypes.html#list) *[*[*str*](https://docs.python.org/3/builtins/stdtypes.html#str) *]*) – + + column names from regulation_data containing group identifiers. + See [pandas.DataFrame.groupby](https://pandas.pydata.org/pandas-docs/stable/reference/api/pandas.DataFrame.groupby.html) for more information. + * **annotation_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from annotation containing annotation terms. + * **rejected_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – name of the column from regulation_data containing boolean for + rejected null hypothesis. + * **group_col** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – column name for new column in annotation dataframe determining + if feature belongs to foreground or background. + * **method** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – method used to compute enrichment + (only ‘fisher’ is supported currently). + * **correction** ([*str*](https://docs.python.org/3/builtins/stdtypes.html#str)) – method to be used for multiple-testing correction + * **alpha** ([*float*](https://docs.python.org/3/builtins/functions.html#float)) – adjusted p-value cutoff to define significance + * **lfc_cutoff** ([*float*](https://docs.python.org/3/builtins/functions.html#float)) – log fold-change cutoff to define practical significance +* **Returns:** + DataFrame adhering to EnrichmentAnalysisSchema +* **Return type:** + DataFrame[[EnrichmentAnalysisSchema](acore.types.md#acore.types.enrichment_analysis.EnrichmentAnalysisSchema)] + +Example: + +```default +result = run_up_down_regulation_enrichment( + regulation_data, + annotation, + identifier='identifier', + groups=['group1', + 'group2'], + annotation_col='annotation', + rejected_col='rejected', + group_col='group', + method='fisher', + correction='fdr_bh', + alpha=0.05, + lfc_cutoff=1, +) +``` + ## Subpackages * [acore.enrichment_analysis.statistical_tests namespace](acore.enrichment_analysis.statistical_tests.md) diff --git a/docs/markdown_ref/acore.md b/docs/markdown_ref/acore.md index 4a286b52..04ca55e2 100644 --- a/docs/markdown_ref/acore.md +++ b/docs/markdown_ref/acore.md @@ -67,11 +67,11 @@ * [`qc_rlsc_loess()`](acore.drift_correction.md#acore.drift_correction.loess_drift_correction.qc_rlsc_loess) * [`run_loess_drift_correction()`](acore.drift_correction.md#acore.drift_correction.loess_drift_correction.run_loess_drift_correction) * [acore.enrichment_analysis package](acore.enrichment_analysis.md) - * [`run_site_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_site_regulation_enrichment) - * [`run_up_down_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_up_down_regulation_enrichment) - * [`run_gsea()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_gsea) * [`run_fisher()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_fisher) + * [`run_gsea()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_gsea) * [`run_kolmogorov_smirnov()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_kolmogorov_smirnov) + * [`run_site_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_site_regulation_enrichment) + * [`run_up_down_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_up_down_regulation_enrichment) * [Subpackages](acore.enrichment_analysis.md#subpackages) * [acore.enrichment_analysis.statistical_tests namespace](acore.enrichment_analysis.statistical_tests.md) * [Submodules](acore.enrichment_analysis.statistical_tests.md#submodules) diff --git a/docs/markdown_ref/acore.types.md b/docs/markdown_ref/acore.types.md index 066c636d..27905ec2 100644 --- a/docs/markdown_ref/acore.types.md +++ b/docs/markdown_ref/acore.types.md @@ -201,7 +201,7 @@ Load a schema from YAML. * **Returns:** the backend-specific schema object. -#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict) | [None](https://docs.python.org/3/builtins/constants.html#None) +#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict)[[Any](https://docs.python.org/3/library/typing.html#typing.Any), [Any](https://docs.python.org/3/library/typing.html#typing.Any)] Provide metadata for columns and schema level @@ -437,7 +437,7 @@ Load a schema from YAML. * **Returns:** the backend-specific schema object. -#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict) | [None](https://docs.python.org/3/builtins/constants.html#None) +#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict)[[Any](https://docs.python.org/3/library/typing.html#typing.Any), [Any](https://docs.python.org/3/library/typing.html#typing.Any)] Provide metadata for columns and schema level @@ -683,7 +683,7 @@ Load a schema from YAML. * **Returns:** the backend-specific schema object. -#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict) | [None](https://docs.python.org/3/builtins/constants.html#None) +#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict)[[Any](https://docs.python.org/3/library/typing.html#typing.Any), [Any](https://docs.python.org/3/library/typing.html#typing.Any)] Provide metadata for columns and schema level @@ -915,7 +915,7 @@ Load a schema from YAML. * **Returns:** the backend-specific schema object. -#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict) | [None](https://docs.python.org/3/builtins/constants.html#None) +#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict)[[Any](https://docs.python.org/3/library/typing.html#typing.Any), [Any](https://docs.python.org/3/library/typing.html#typing.Any)] Provide metadata for columns and schema level @@ -1131,7 +1131,7 @@ Load a schema from YAML. * **Returns:** the backend-specific schema object. -#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict) | [None](https://docs.python.org/3/builtins/constants.html#None) +#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict)[[Any](https://docs.python.org/3/library/typing.html#typing.Any), [Any](https://docs.python.org/3/library/typing.html#typing.Any)] Provide metadata for columns and schema level @@ -1341,7 +1341,7 @@ Load a schema from YAML. * **Returns:** the backend-specific schema object. -#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict) | [None](https://docs.python.org/3/builtins/constants.html#None) +#### *classmethod* get_metadata() → [dict](https://docs.python.org/3/builtins/stdtypes.html#dict)[[Any](https://docs.python.org/3/library/typing.html#typing.Any), [Any](https://docs.python.org/3/library/typing.html#typing.Any)] Provide metadata for columns and schema level diff --git a/docs/markdown_ref/modules.md b/docs/markdown_ref/modules.md index 77317eb1..d8ef4991 100644 --- a/docs/markdown_ref/modules.md +++ b/docs/markdown_ref/modules.md @@ -40,11 +40,11 @@ * [acore.drift_correction.cpca_drift_correction module](acore.drift_correction.md#module-acore.drift_correction.cpca_drift_correction) * [acore.drift_correction.loess_drift_correction module](acore.drift_correction.md#module-acore.drift_correction.loess_drift_correction) * [acore.enrichment_analysis package](acore.enrichment_analysis.md) - * [`run_site_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_site_regulation_enrichment) - * [`run_up_down_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_up_down_regulation_enrichment) - * [`run_gsea()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_gsea) * [`run_fisher()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_fisher) + * [`run_gsea()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_gsea) * [`run_kolmogorov_smirnov()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_kolmogorov_smirnov) + * [`run_site_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_site_regulation_enrichment) + * [`run_up_down_regulation_enrichment()`](acore.enrichment_analysis.md#acore.enrichment_analysis.run_up_down_regulation_enrichment) * [Subpackages](acore.enrichment_analysis.md#subpackages) * [Submodules](acore.enrichment_analysis.md#submodules) * [acore.enrichment_analysis.annotate module](acore.enrichment_analysis.md#module-acore.enrichment_analysis.annotate) From 51e2bda3c5a6fef8c96376a8d4b480c6fa2bfd3a Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Tue, 20 Oct 2026 11:12:21 +0200 Subject: [PATCH 09/14] :construction: add rejected column to GSEA output table --- .../enrichment_analysis.ipynb | 21 +++++++++++++++++++ .../enrichment_analysis.py | 9 ++++++++ src/acore/enrichment_analysis/__init__.py | 12 ++++++++--- tests/enrichment_analysis/test_enrichment.py | 3 +++ 4 files changed, 42 insertions(+), 3 deletions(-) diff --git a/docs/examples_single_analysis/enrichment_analysis.ipynb b/docs/examples_single_analysis/enrichment_analysis.ipynb index 2a47609b..7aa98f7b 100644 --- a/docs/examples_single_analysis/enrichment_analysis.ipynb +++ b/docs/examples_single_analysis/enrichment_analysis.ipynb @@ -694,6 +694,7 @@ " min_size=1,\n", " permutation_num=100,\n", " seed=42,\n", + " correction_alpha=0.05,\n", ")\n", "gsea_results" ] @@ -708,6 +709,26 @@ "values indicate enrichment among the **down**-regulated features." ] }, + { + "cell_type": "markdown", + "id": "fc0abd44", + "metadata": {}, + "source": [ + "gseapy reports a permutation-based `FDR q-val` for each gene set (derived from the NES\n", + "null distribution), so no further correction is needed. The `rejected` column flags\n", + "gene sets with `FDR q-val <= correction_alpha`." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "fdfaa460", + "metadata": {}, + "outputs": [], + "source": [ + "gsea_results[[\"Term\", \"NES\", \"NOM p-val\", \"FDR q-val\", \"rejected\"]]" + ] + }, { "cell_type": "code", "execution_count": null, diff --git a/docs/examples_single_analysis/enrichment_analysis.py b/docs/examples_single_analysis/enrichment_analysis.py index 70c3162a..9019559b 100644 --- a/docs/examples_single_analysis/enrichment_analysis.py +++ b/docs/examples_single_analysis/enrichment_analysis.py @@ -386,6 +386,7 @@ min_size=1, permutation_num=100, seed=42, + correction_alpha=0.05, ) gsea_results @@ -394,6 +395,14 @@ # indicate that the gene set is enriched among the **up**-regulated features and negative # values indicate enrichment among the **down**-regulated features. +# %% [markdown] +# gseapy reports a permutation-based `FDR q-val` for each gene set (derived from the NES +# null distribution), so no further correction is needed. The `rejected` column flags +# gene sets with `FDR q-val <= correction_alpha`. + +# %% +gsea_results[["Term", "NES", "NOM p-val", "FDR q-val", "rejected"]] + # %% ax = gsea_results["NES"].plot.hist(title="Distribution of NES", xlabel="NES") diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index ede8df78..784245d9 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -616,6 +616,7 @@ def run_gsea( seed: int = 42, weight: float = 1.0, ascending: bool = False, + correction_alpha: float = 0.05, ) -> pd.DataFrame: """ Run Gene Set Enrichment Analysis (GSEA) on a pre-ranked gene list using @@ -658,11 +659,14 @@ def run_gsea( gives classic (unweighted) GSEA. :param bool ascending: Sort order for the ranking. ``False`` (default) means descending (highest ranking metric first, suited for log2FC). + :param float correction_alpha: FDR q-value cutoff used to define the + ``'rejected'`` column. Default: 0.05. :return: DataFrame with one row per gene set. Columns are: ``'Term'``, ``'ES'`` (enrichment score), ``'NES'`` (normalised enrichment score), - ``'NOM p-val'`` (nominal p-value), ``'FDR q-val'`` (FDR-corrected - q-value), ``'FWER p-val'``, ``'Tag %'``, ``'Gene %'``, - ``'Lead_genes'`` (semi-colon-separated leading-edge genes). + ``'NOM p-val'`` (nominal p-value), ``'FDR q-val'`` (gseapy's + permutation-based FDR q-value), ``'FWER p-val'``, ``'Tag %'``, + ``'Gene %'``, ``'Lead_genes'`` (semi-colon-separated leading-edge genes) + and ``'rejected'`` (``'FDR q-val' <= correction_alpha``). :rtype: pd.DataFrame Example:: @@ -726,4 +730,6 @@ def run_gsea( outdir=None, ) result = pd.DataFrame(enrichment.res2d).drop(columns=["Name"]) + fdr = pd.to_numeric(result["FDR q-val"], errors="coerce") + result["rejected"] = (fdr <= correction_alpha).where(fdr.notna()).astype("boolean") return result diff --git a/tests/enrichment_analysis/test_enrichment.py b/tests/enrichment_analysis/test_enrichment.py index 8d3d73ff..6f3a206c 100644 --- a/tests/enrichment_analysis/test_enrichment.py +++ b/tests/enrichment_analysis/test_enrichment.py @@ -634,6 +634,9 @@ def test_run_gsea_basic(): assert not result.empty expected_cols = {"Term", "ES", "NES", "NOM p-val", "FDR q-val", "Lead_genes"} assert expected_cols.issubset(result.columns) + assert "rejected" in result.columns + fdr = result["FDR q-val"].astype(float) + assert (result["rejected"] == (fdr <= 0.05)).all() assert set(result["Term"]) == {"path1", "path2"} # path1 genes are up-regulated -> positive NES nes_path1 = result.loc[result["Term"] == "path1", "NES"].iloc[0] From 86c73dc19442a34fd6cd5791b3d2acd17a4369f5 Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" <41898282+github-actions[bot]@users.noreply.github.com> Date: Tue, 20 Oct 2026 09:14:38 +0000 Subject: [PATCH 10/14] docs: update markdown reference --- docs/markdown_ref/acore.enrichment_analysis.md | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/docs/markdown_ref/acore.enrichment_analysis.md b/docs/markdown_ref/acore.enrichment_analysis.md index dc7b07ad..a3c21580 100644 --- a/docs/markdown_ref/acore.enrichment_analysis.md +++ b/docs/markdown_ref/acore.enrichment_analysis.md @@ -23,7 +23,7 @@ odds, pvalue = stats.fisher_exact(group1=[a, b], ) ``` -### run_gsea(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', ranking_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'log2FC', annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', identifier_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', min_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 15, max_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 500, permutation_num: [int](https://docs.python.org/3/builtins/functions.html#int) = 1000, seed: [int](https://docs.python.org/3/builtins/functions.html#int) = 42, weight: [float](https://docs.python.org/3/builtins/functions.html#float) = 1.0, ascending: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False) → DataFrame +### run_gsea(regulation_data: DataFrame, annotation: DataFrame, identifier: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', ranking_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'log2FC', annotation_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'annotation', identifier_col: [str](https://docs.python.org/3/builtins/stdtypes.html#str) = 'identifier', min_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 15, max_size: [int](https://docs.python.org/3/builtins/functions.html#int) = 500, permutation_num: [int](https://docs.python.org/3/builtins/functions.html#int) = 1000, seed: [int](https://docs.python.org/3/builtins/functions.html#int) = 42, weight: [float](https://docs.python.org/3/builtins/functions.html#float) = 1.0, ascending: [bool](https://docs.python.org/3/builtins/functions.html#bool) = False, correction_alpha: [float](https://docs.python.org/3/builtins/functions.html#float) = 0.05) → DataFrame Run Gene Set Enrichment Analysis (GSEA) on a pre-ranked gene list using [gseapy.prerank](https://gseapy.readthedocs.io/en/latest/run.html#gseapy.prerank). @@ -66,12 +66,15 @@ full ranked list instead. gives classic (unweighted) GSEA. * **ascending** ([*bool*](https://docs.python.org/3/builtins/functions.html#bool)) – Sort order for the ranking. `False` (default) means descending (highest ranking metric first, suited for log2FC). + * **correction_alpha** ([*float*](https://docs.python.org/3/builtins/functions.html#float)) – FDR q-value cutoff used to define the + `'rejected'` column. Default: 0.05. * **Returns:** DataFrame with one row per gene set. Columns are: `'Term'`, `'ES'` (enrichment score), `'NES'` (normalised enrichment score), - `'NOM p-val'` (nominal p-value), `'FDR q-val'` (FDR-corrected - q-value), `'FWER p-val'`, `'Tag %'`, `'Gene %'`, - `'Lead_genes'` (semi-colon-separated leading-edge genes). + `'NOM p-val'` (nominal p-value), `'FDR q-val'` (gseapy’s + permutation-based FDR q-value), `'FWER p-val'`, `'Tag %'`, + `'Gene %'`, `'Lead_genes'` (semi-colon-separated leading-edge genes) + and `'rejected'` (`'FDR q-val' <= correction_alpha`). * **Return type:** pd.DataFrame From 18c0f248eda89308906490e64ea83940dfe6a4ed Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Tue, 20 Oct 2026 13:04:13 +0200 Subject: [PATCH 11/14] :art: use dictionary of pathways instead of local gmt file --- src/acore/enrichment_analysis/__init__.py | 34 +++++++++++++---------- 1 file changed, 19 insertions(+), 15 deletions(-) diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index 784245d9..376d5d5e 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -10,7 +10,6 @@ import logging import os import re -import uuid import gseapy as gp import pandas as pd @@ -571,22 +570,27 @@ def run_ssgsea( " as specified by `identifier_col`" ) - grouped_annotations = ( - annotation.groupby(annotation_col)[identifier_col].apply(list).reset_index() - ) - fid = uuid.uuid4() - file_path = os.path.join(outdir, str(fid) + ".gmt") - with open(file_path, "w", encoding="utf8") as out: - for _, row in grouped_annotations.iterrows(): - out.write( - row[annotation_col] - + "\t" - + "\t".join(list(filter(None, row[identifier_col]))) - + "\n" - ) + # grouped_annotations = ( + # annotation.groupby(annotation_col)[identifier_col].apply(list).reset_index() + # ) + # fid = uuid.uuid4() + # file_path = os.path.join(outdir, str(fid) + ".gmt") + # with open(file_path, "w", encoding="utf8") as out: + # for _, row in grouped_annotations.iterrows(): + # out.write( + # row[annotation_col] + # + "\t" + # + "\t".join(list(filter(None, row[identifier_col]))) + # + "\n" + # ) + # In case annotations were not provided by acore itself: + annotation = annotation.dropna(subset=[identifier_col]) + annotation = annotation[annotation[identifier_col] != ""] + gene_sets = annotation.groupby(annotation_col)[identifier_col].apply(list).to_dict() + enrichment = gp.ssgsea( data=df, - gene_sets=str(file_path), + gene_sets=gene_sets, outdir=outdir, min_size=min_size, max_size=max_size, From 055d8994dab05fb5629d8c6cdfaaf41ccc5751e4 Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Tue, 20 Oct 2026 13:57:02 +0200 Subject: [PATCH 12/14] :bug: decide to drop gmt file dropping as it contained a bug, verified I saw the difference in results after changing the implementation to skip the tmp-file gmt file dump and then confirmed it with AI: Description column: The GMT format is termdescriptiongene1gene2..., so the parser skips entries[1]. What the old code wrote: The code you replaced wrote termgene1gene2..., with no description column. The parser therefore treated gene1 as the description and dropped it. Consequence: The file-based run lost the first gene of every set. Sets were smaller, and some may have been removed by the min_size filter. Scores were computed on different sets. --- src/acore/enrichment_analysis/__init__.py | 10 +++++----- tests/enrichment_analysis/test_enrichment.py | 1 - 2 files changed, 5 insertions(+), 6 deletions(-) diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index 376d5d5e..24dc84d8 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -569,20 +569,20 @@ def run_ssgsea( f"Missing Identifier Column: {identifier_col}" " as specified by `identifier_col`" ) - + # ! Legacy implementation that was missing the description column of gmt files: + # ?gmt-columns: termdescriptiongene1gene2 # grouped_annotations = ( # annotation.groupby(annotation_col)[identifier_col].apply(list).reset_index() # ) - # fid = uuid.uuid4() - # file_path = os.path.join(outdir, str(fid) + ".gmt") - # with open(file_path, "w", encoding="utf8") as out: + # gene_sets = os.path.join(outdir, "ss_gsea" + ".gmt") + # with open(gene_sets, "w", encoding="utf8") as out: # for _, row in grouped_annotations.iterrows(): # out.write( # row[annotation_col] # + "\t" # + "\t".join(list(filter(None, row[identifier_col]))) # + "\n" - # ) + # ) # In case annotations were not provided by acore itself: annotation = annotation.dropna(subset=[identifier_col]) annotation = annotation[annotation[identifier_col] != ""] diff --git a/tests/enrichment_analysis/test_enrichment.py b/tests/enrichment_analysis/test_enrichment.py index 6f3a206c..9018e1db 100644 --- a/tests/enrichment_analysis/test_enrichment.py +++ b/tests/enrichment_analysis/test_enrichment.py @@ -677,4 +677,3 @@ def test_run_gsea_raises_on_duplicate_identifiers(): regulation_data=regulation_data, annotation=annotation, ) - From 8eb268845ba5de5ea6de5e2dbf8abd26acbf02cb Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Tue, 20 Oct 2026 14:00:40 +0200 Subject: [PATCH 13/14] :art: restore old 'wrong' version by setting minimum to 2 Same list of terms, but different enrichments score (as one more gene was always included in the set) - now the score distrubtion looks also more sensible (to me): normal instead of bimodal - show a second factor in the PCA plot --- docs/examples_single_analysis/enrichment_analysis.ipynb | 4 ++-- docs/examples_single_analysis/enrichment_analysis.py | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/docs/examples_single_analysis/enrichment_analysis.ipynb b/docs/examples_single_analysis/enrichment_analysis.ipynb index 7aa98f7b..ede67f7f 100644 --- a/docs/examples_single_analysis/enrichment_analysis.ipynb +++ b/docs/examples_single_analysis/enrichment_analysis.ipynb @@ -546,7 +546,7 @@ "enrichtments = acore.enrichment_analysis.run_ssgsea(\n", " data=df_omics,\n", " annotation=annotations,\n", - " min_size=1,\n", + " min_size=2,\n", ")\n", "enrichtments" ] @@ -648,7 +648,7 @@ "metadata": {}, "outputs": [], "source": [ - "args = {\"factor\": 2, \"loadings\": 1} # increase number of loadings or scaling factor\n", + "args = {\"factor\": 2, \"loadings\": 2} # increase number of loadings or scaling factor\n", "#! pca_results has three items, but docstring requests only two -> double check\n", "figure = viz.get_pca_plot(data=pca_result, identifier=\"PCA enrichment\", args=args)\n", "figure = go.Figure(data=figure[\"data\"], layout=figure[\"layout\"])\n", diff --git a/docs/examples_single_analysis/enrichment_analysis.py b/docs/examples_single_analysis/enrichment_analysis.py index 9019559b..8913ebe9 100644 --- a/docs/examples_single_analysis/enrichment_analysis.py +++ b/docs/examples_single_analysis/enrichment_analysis.py @@ -307,7 +307,7 @@ enrichtments = acore.enrichment_analysis.run_ssgsea( data=df_omics, annotation=annotations, - min_size=1, + min_size=2, ) enrichtments @@ -352,7 +352,7 @@ # for this, which is also developed by the Multiomics Analytics Group. # %% -args = {"factor": 2, "loadings": 1} # increase number of loadings or scaling factor +args = {"factor": 2, "loadings": 2} # increase number of loadings or scaling factor # #! pca_results has three items, but docstring requests only two -> double check figure = viz.get_pca_plot(data=pca_result, identifier="PCA enrichment", args=args) figure = go.Figure(data=figure["data"], layout=figure["layout"]) From e2a8d81bf0b17976ebaaf34b06d59716ad18808e Mon Sep 17 00:00:00 2001 From: Henry Webel Date: Tue, 20 Oct 2026 15:37:51 +0200 Subject: [PATCH 14/14] :art: check that numeric columns are numeric --- src/acore/enrichment_analysis/__init__.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/src/acore/enrichment_analysis/__init__.py b/src/acore/enrichment_analysis/__init__.py index 24dc84d8..6314c8a5 100644 --- a/src/acore/enrichment_analysis/__init__.py +++ b/src/acore/enrichment_analysis/__init__.py @@ -734,6 +734,13 @@ def run_gsea( outdir=None, ) result = pd.DataFrame(enrichment.res2d).drop(columns=["Name"]) - fdr = pd.to_numeric(result["FDR q-val"], errors="coerce") + # 'Tag %' and 'Gene %' are strings like "5/20" and "10.00%", so they stay as is. + numeric_cols = ["ES", "NES", "NOM p-val", "FDR q-val", "FWER p-val"] + non_numeric_cols = result[numeric_cols].select_dtypes(exclude="number").columns + if len(non_numeric_cols): + result[non_numeric_cols] = result[non_numeric_cols].apply( + pd.to_numeric, errors="coerce" + ) + fdr = result["FDR q-val"] result["rejected"] = (fdr <= correction_alpha).where(fdr.notna()).astype("boolean") return result