diff --git a/.gitignore b/.gitignore index d77eadec..2e80d370 100644 --- a/.gitignore +++ b/.gitignore @@ -141,6 +141,9 @@ venv.bak/ *.code-workspace .vscode +# Zed +.zed/ + # mkdocs documentation /site diff --git a/HOW_TO_RELEASE.md b/HOW_TO_RELEASE.md index 0a051930..c84e38ce 100644 --- a/HOW_TO_RELEASE.md +++ b/HOW_TO_RELEASE.md @@ -9,62 +9,10 @@ The document describes how to release `gpsea` to *PyPi*. - remove deprecated methods targeted for removal in this version. The `TODO` markers are labeled using the target version (e.g. `TODO[v0.3.0]`) - bump versions to a release: - - `src/gpsea/__init__.py` + - `pyproject.toml` - ensure the CI passes -- deploy to PyPi (described below) - merge to `main` -- create a GitHub release from the latest `main` commit, including a new tag +- create a GitHub release from the latest `main` commit, including a new tag. The `release.yml` workflow takes over and deploys the new release to PyPi. - merge `main` to `develop` - bump versions to a `dev` version to begin next development iteration -## Deploy to PyPi - -### Virtual environment for deployment - -As an optional one-time step, consider creating a dedicated virtual environment with the packages required -for testing, building, and deployment: - -```shell -# Create and activate the virtual environment -python3 -m venv build -source build/bin/activate - -# Install the build packages -python3 -m pip install build twine -``` - -### Setup PyPi credentials - -As another one-time action, you must create a profile on PyPi and generate an access token. -The token is used to upload the packages. - -First, create an account (e.g. associated with your GitHub account), configure 2FA, store recovery codes, etc. -Then, generate a token to upload the packages. You can generate a token per project or for all projects. -Store the token into `$HOME/.pypirc` file with `-rw-------` permissions. The file should look like: - -``` -[pypi] - username = __token__ - password = -``` - -Now we're ready to publish packages! - -### Deploy -Run the following to deploy `gpsea` to PyPi: - -```bash -# Ensure you're on the release branch -cd gpsea - -# Build the package -python3 -m build - -# Deploy -python3 -m twine upload dist/* - -# Clear the built and deployed files -rm -rf build dist -``` - -The commands will build source distribution and a wheel, and deploy the source distribution and wheel to PyPi. diff --git a/docs/conf.py b/docs/conf.py index df80de54..7cd76577 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -53,7 +53,7 @@ # General information about the project. project = u'GPSEA' -copyright = u'2025' +copyright = u'2026' author = u'Lauren Rekerle, Daniel Danis, Peter N Robinson' # The version info for the project you're documenting, acts as replacement for diff --git a/docs/tutorial.rst b/docs/tutorial.rst index 061961a2..256322aa 100644 --- a/docs/tutorial.rst +++ b/docs/tutorial.rst @@ -332,11 +332,12 @@ Statistical testing Now we can perform the testing and evaluate the results. ->>> result = analysis.compare_genotype_vs_phenotypes( +>>> result = analysis.compare_genotype_vs_phenotypes( # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE ... cohort=cohort, ... gt_clf=gt_clf, ... pheno_clfs=pheno_clfs, ... ) +HPO terms processed: ... >>> result.total_tests 30 diff --git a/docs/user-guide/analyses/phenotype-classes.rst b/docs/user-guide/analyses/phenotype-classes.rst index 41785f31..27c838cb 100644 --- a/docs/user-guide/analyses/phenotype-classes.rst +++ b/docs/user-guide/analyses/phenotype-classes.rst @@ -287,11 +287,12 @@ Analysis We can now test associations between the genotype classes and the HPO terms: ->>> result = analysis.compare_genotype_vs_phenotypes( +>>> result = analysis.compare_genotype_vs_phenotypes( # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE ... cohort=cohort, ... gt_clf=gt_clf, ... pheno_clfs=pheno_clfs, ... ) +HPO terms processed: ... >>> len(result.phenotypes) 369 >>> result.total_tests diff --git a/pyproject.toml b/pyproject.toml index fbf34d56..de4c8c78 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,10 +4,10 @@ build-backend = "setuptools.build_meta" [project] name = "gpsea" -version = "0.9.13" +version = "0.9.14" authors = [ { name = "Lauren Rekerle", email = "lauren.rekerle@jax.org" }, - { name = "Daniel Danis", email = "daniel.danis@bih-charite.de" }, + { name = "Daniel Danis", email = "daniel.gordon.danis@protonmail.com" }, { name = "Peter Robinson", email = "peter.robinson@bih-charite.de" }, ] description = "Discover genotype-phenotype correlations with GA4GH phenopackets" @@ -54,6 +54,7 @@ test = [ "phenopacket-store-toolkit>=0.1.2", "pytest>=7.0.0,<8.0.0", "ruff==0.12.1", + "basedpyright==v1.39.10", ] docs = [ "sphinx>=7.0.0", @@ -70,6 +71,3 @@ bugtracker = "https://github.com/P2GX/gpsea/issues" [tool.setuptools] package-dir = { "" = "src" } - -[tool.setuptools.dynamic] -version = { attr = "gpsea.__version__" } diff --git a/src/gpsea/analysis/clf/_pheno.py b/src/gpsea/analysis/clf/_pheno.py index 6da89e43..e184adf2 100644 --- a/src/gpsea/analysis/clf/_pheno.py +++ b/src/gpsea/analysis/clf/_pheno.py @@ -76,9 +76,6 @@ def test( ) -> typing.Optional[PhenotypeCategorization[hpotk.TermId]]: self._check_patient(patient) - if len(patient.phenotypes) == 0: - return None - for phenotype in patient.phenotypes: if phenotype.is_present: if self._query == phenotype.identifier or any( @@ -94,7 +91,10 @@ def test( ): return self._phenotype_excluded - return None + if self._missing_implies_phenotype_excluded: + return self._phenotype_excluded + else: + return None def __eq__(self, value: object) -> bool: return ( diff --git a/src/gpsea/analysis/clf/_util.py b/src/gpsea/analysis/clf/_util.py index 7235252a..38720985 100644 --- a/src/gpsea/analysis/clf/_util.py +++ b/src/gpsea/analysis/clf/_util.py @@ -4,7 +4,8 @@ import hpotk -from ._pheno import PhenotypeClassifier, HpoClassifier +from ._api import PhenotypeClassifier +from ._pheno import HpoClassifier from gpsea.model import Patient diff --git a/src/gpsea/analysis/pcats/_impl.py b/src/gpsea/analysis/pcats/_impl.py index 2bcbc382..64c5dd00 100644 --- a/src/gpsea/analysis/pcats/_impl.py +++ b/src/gpsea/analysis/pcats/_impl.py @@ -1,24 +1,22 @@ import abc +import collections.abc import os +import sys import typing - from collections import Counter import hpotk import numpy as np import pandas as pd - +import tqdm from statsmodels.stats import multitest from gpsea.model import Patient -from ..clf import GenotypeClassifier -from ..clf import P, PhenotypeClassifier +from .._base import MultiPhenotypeAnalysisResult, StatisticResult +from ..clf import GenotypeClassifier, P, PhenotypeClassifier from ..mtc_filter import PhenotypeMtcFilter, PhenotypeMtcResult - from .stats import CountStatistic -from .._base import MultiPhenotypeAnalysisResult, StatisticResult - DEFAULT_MTC_PROCEDURE = "fdr_bh" """ @@ -27,12 +25,12 @@ def apply_classifiers_on_individuals( - individuals: typing.Iterable[Patient], + individuals: collections.abc.Iterable[Patient], gt_clf: GenotypeClassifier, - pheno_clfs: typing.Sequence[PhenotypeClassifier[P]], -) -> typing.Tuple[ - typing.Sequence[int], - typing.Sequence[pd.DataFrame], + pheno_clfs: collections.abc.Sequence[PhenotypeClassifier[P]], +) -> tuple[ + collections.abc.Sequence[int], + collections.abc.Sequence[pd.DataFrame], ]: """ Classify individuals with the genotype and phenotype classifiers. @@ -57,7 +55,12 @@ def apply_classifiers_on_individuals( # Apply genotype and phenotype predicates count_dict = {} - for ph_predicate in pheno_clfs: + for ph_predicate in tqdm.tqdm( + pheno_clfs, + desc="HPO terms processed", + file=sys.stdout, + unit=" terms", + ): if ph_predicate.phenotype not in count_dict: # Make an empty frame for keeping track of the counts. count_dict[ph_predicate.phenotype] = pd.DataFrame( @@ -83,7 +86,7 @@ def apply_classifiers_on_individuals( # Convert dicts to numpy arrays n_usable_patients = [n_usable_patient_counter[ph_predicate.phenotype] for ph_predicate in pheno_clfs] - counts = [count_dict[ph_predicate.phenotype] for ph_predicate in pheno_clfs] + counts: list[pd.DataFrame] = [count_dict[ph_predicate.phenotype] for ph_predicate in pheno_clfs] return n_usable_patients, counts diff --git a/tests/analysis/clf/test_disease.py b/tests/analysis/clf/test_disease.py new file mode 100644 index 00000000..58725a73 --- /dev/null +++ b/tests/analysis/clf/test_disease.py @@ -0,0 +1,36 @@ +import hpotk +import pytest + +from gpsea.analysis.clf import DiseasePresenceClassifier +from gpsea.model import Cohort, Patient + + +class TestDiseasePresencePredicate: + @pytest.mark.parametrize( + "patient_id, patient_category", + [ + ("HetSingleVar", "Yes"), + ("HomoVar", "No"), + ], + ) + def test_disease_predicate( + self, + patient_id: str, + patient_category: str, + toy_cohort: Cohort, + ): + patient = find_patient(patient_id, toy_cohort) + disease_id = hpotk.TermId.from_curie("OMIM:148050") + predicate = DiseasePresenceClassifier(disease_id) + actual = predicate.test(patient) + + assert actual is not None + assert actual.phenotype == disease_id + assert actual.category.name == patient_category + + +def find_patient(pat_id: str, cohort: Cohort) -> Patient: + for pat in cohort.all_patients: + if pat.patient_id == pat_id: + return pat + raise ValueError(f"Could not find patient {pat_id}") diff --git a/tests/analysis/clf/test_phenotype.py b/tests/analysis/clf/test_phenotype.py index 35f44f67..6dfdd54d 100644 --- a/tests/analysis/clf/test_phenotype.py +++ b/tests/analysis/clf/test_phenotype.py @@ -1,14 +1,18 @@ import hpotk +import pytest -from gpsea.model import Cohort - -from gpsea.analysis.clf import PhenotypeClassifier -from gpsea.analysis.clf import prepare_hpo_terms_of_interest, prepare_classifiers_for_terms_of_interest +from gpsea.analysis.clf import ( + HpoClassifier, + PhenotypeClassifier, + prepare_classifiers_for_terms_of_interest, + prepare_hpo_terms_of_interest, +) +from gpsea.model import Cohort, Patient def test_prepare_hpo_terms_of_interest( suox_cohort: Cohort, - hpo: hpotk.MinimalOntology, + hpo: hpotk.MinimalOntology[hpotk.TermId, hpotk.MinimalTerm], ): terms = prepare_hpo_terms_of_interest( cohort=suox_cohort, @@ -20,7 +24,7 @@ def test_prepare_hpo_terms_of_interest( def test_prepare_predicates_for_terms_of_interest( suox_cohort: Cohort, - hpo: hpotk.MinimalOntology, + hpo: hpotk.MinimalOntology[hpotk.TermId, hpotk.MinimalTerm], ): predicates = prepare_classifiers_for_terms_of_interest( cohort=suox_cohort, @@ -29,3 +33,123 @@ def test_prepare_predicates_for_terms_of_interest( assert len(predicates) == 71 assert all(isinstance(p, PhenotypeClassifier) for p in predicates) + + +class TestHpoPredicate: + @pytest.mark.parametrize( + "curie, patient_id, expected", + # Patient "HetSingleVar" has Phenotypes: + # Measured and present - 'HP:0001166;HP:0002266', # Arachnodactyly;Focal clonic seizure + # Measured but excluded - 'HP:0001257', # Spasticity + [ + # Test exact match + ( + "HP:0001166", # Arachnodactyly + "HetSingleVar", + "Yes", + ), + # Test inferred annotations + ( + "HP:0001250", # Seizure + "HetSingleVar", + "Yes", + ), + # Test excluded feature + ( + "HP:0001257", # Spasticity + "HetSingleVar", + "No", + ), + ], + ) + def test_phenotype_predicate__present_or_excluded( + self, + toy_cohort: Cohort, + hpo: hpotk.MinimalOntology[hpotk.TermId, hpotk.MinimalTerm], + curie: str, + patient_id: str, + expected: str, + ): + patient = find_patient(patient_id, toy_cohort) + term_id = hpotk.TermId.from_curie(curie) + predicate = HpoClassifier(hpo=hpo, query=term_id) + actual = predicate.test(patient) + + assert actual is not None + assert actual.phenotype == term_id + assert actual.category.name == expected + + def test_phenotype_predicate__unknown( + self, + toy_cohort: Cohort, + hpo: hpotk.MinimalOntology[hpotk.TermId, hpotk.MinimalTerm], + ): + # Not Measured and not Observed - 'HP:0006280', # Chronic pancreatitis + patient = find_patient("HetSingleVar", toy_cohort) + term_id = hpotk.TermId.from_curie("HP:0006280") + predicate = HpoClassifier(hpo=hpo, query=term_id) + actual = predicate.test(patient) + + assert actual is None + + @pytest.mark.parametrize( + "curie,patient_id,expected", + [ + # Test exact match + ( + "HP:0001166", # Arachnodactyly + "HetDoubleVar1", + "Yes", + ), + # An explicitly excluded feature is categorized + # even when missing implies excluded + ( + "HP:0001257", # Spasticity + "HetSingleVar", + "No", + ), + # An unannotated feature is implied to be excluded. + # This works in an individual with an excluded feature. + # Not Measured and not Observed - 'HP:0006280' + ( + "HP:0006280", # Chronic pancreatitis + "HetSingleVar", + "No", + ), + # An unannotated feature is implied to be excluded. + # This works even in an individual with no excluded features. + # Not Measured and not Observed - 'HP:0006280' + ( + "HP:0006280", # Chronic pancreatitis + "HetDoubleVar1", + "No", + ), + ], + ) + def test_phenotype_predicate__missing_implies_excluded( + self, + toy_cohort: Cohort, + hpo: hpotk.MinimalOntology[hpotk.TermId, hpotk.MinimalTerm], + curie: str, + patient_id: str, + expected: str, + ): + patient = find_patient(patient_id, toy_cohort) + term_id = hpotk.TermId.from_curie(curie) + predicate = HpoClassifier( + hpo=hpo, + query=term_id, + missing_implies_phenotype_excluded=True, + ) + actual = predicate.test(patient) + + assert actual is not None + assert actual.phenotype == term_id + assert actual.category.name == expected + + +def find_patient(pat_id: str, cohort: Cohort) -> Patient: + for pat in cohort.all_patients: + if pat.patient_id == pat_id: + return pat + raise ValueError(f"Could not find patient {pat_id}") diff --git a/tests/test_predicates.py b/tests/test_predicates.py deleted file mode 100644 index 1ab9858f..00000000 --- a/tests/test_predicates.py +++ /dev/null @@ -1,94 +0,0 @@ -import hpotk -import pytest - -from gpsea.analysis.clf import HpoClassifier, DiseasePresenceClassifier -from gpsea.model import Cohort, Patient - - -def find_patient(pat_id: str, cohort: Cohort) -> Patient: - for pat in cohort.all_patients: - if pat.patient_id == pat_id: - return pat - raise ValueError(f"Could not find patient {pat_id}") - - -class TestHpoPredicate: - @pytest.mark.parametrize( - "curie, patient_id, expected", - # Patient "HetSingleVar" has Phenotypes: - # Measured and present - 'HP:0001166;HP:0002266', # Arachnodactyly;Focal clonic seizure - # Measured but excluded - 'HP:0001257', # Spasticity - [ - # Test exact match - ( - "HP:0001166", # Arachnodactyly - "HetSingleVar", - "Yes", - ), - # Test inferred annotations - ( - "HP:0001250", # Seizure - "HetSingleVar", - "Yes", - ), - # Test excluded feature - ( - "HP:0001257", # Spasticity - "HetSingleVar", - "No", - ), - ], - ) - def test_phenotype_predicate__present_or_excluded( - self, - toy_cohort: Cohort, - hpo: hpotk.MinimalOntology, - curie: str, - patient_id: str, - expected: str, - ): - patient = find_patient(patient_id, toy_cohort) - term_id = hpotk.TermId.from_curie(curie) - predicate = HpoClassifier(hpo=hpo, query=term_id) - actual = predicate.test(patient) - - assert actual is not None - assert actual.phenotype == term_id - assert actual.category.name == expected - - def test_phenotype_predicate__unknown( - self, - toy_cohort: Cohort, - hpo: hpotk.MinimalOntology, - ): - # Not Measured and not Observed - 'HP:0006280', # Chronic pancreatitis - patient = find_patient("HetSingleVar", toy_cohort) - term_id = hpotk.TermId.from_curie("HP:0006280") - predicate = HpoClassifier(hpo=hpo, query=term_id) - actual = predicate.test(patient) - - assert actual is None - - -class TestDiseasePresencePredicate: - @pytest.mark.parametrize( - "patient_id, patient_category", - [ - ("HetSingleVar", "Yes"), - ("HomoVar", "No"), - ], - ) - def test_disease_predicate( - self, - patient_id: str, - patient_category: str, - toy_cohort: Cohort, - ): - patient = find_patient(patient_id, toy_cohort) - disease_id = hpotk.TermId.from_curie("OMIM:148050") - predicate = DiseasePresenceClassifier(disease_id) - actual = predicate.test(patient) - - assert actual is not None - assert actual.phenotype == disease_id - assert actual.category.name == patient_category