From 24f4644c3d76cfa3f453540defbd052b8c3ec56d Mon Sep 17 00:00:00 2001 From: Marc Balestreri Date: Sat, 28 Mar 2026 22:40:49 -0700 Subject: [PATCH] feat: add clinical layer to GeoData and load NHANES/FHS through the library - Add a clinical layer on GeoData (samples-as-rows, biomarkers-as-columns, same orientation as metadata) - Add load_nhanes_as_geodata and load_fhs_as_geodata so both data sources flow through the library and return a GeoData - Add a biomarker registry with unit conversions, validated end-to-end against real FHS Period 1 data (glucose mg/dL -> mmol/L) - Drop the unverified UK Biobank source preset; we have not validated it against real UK Biobank data - Add an example showing the NHANES through-the-library pattern The required_features() interface and the consumer-facing missing-data error path will be a separate PR. Addresses #194 --- biolearn/clinical/__init__.py | 4 + biolearn/clinical/convert.py | 99 ++++++++ biolearn/clinical/registry.py | 216 ++++++++++++++++++ biolearn/data/library.yaml | 37 +++ biolearn/data_library.py | 118 ++++++++++ biolearn/load.py | 30 ++- biolearn/test/test_clinical_layer.py | 167 ++++++++++++++ biolearn/test/test_load.py | 49 ++++ biolearn/test/test_registry.py | 130 +++++++++++ .../plot_load_nhanes_through_library.py | 48 ++++ 10 files changed, 892 insertions(+), 6 deletions(-) create mode 100644 biolearn/clinical/__init__.py create mode 100644 biolearn/clinical/convert.py create mode 100644 biolearn/clinical/registry.py create mode 100644 biolearn/test/test_clinical_layer.py create mode 100644 biolearn/test/test_registry.py create mode 100644 examples/01_composite_biomarkers/plot_load_nhanes_through_library.py diff --git a/biolearn/clinical/__init__.py b/biolearn/clinical/__init__.py new file mode 100644 index 0000000..380b060 --- /dev/null +++ b/biolearn/clinical/__init__.py @@ -0,0 +1,4 @@ +from biolearn.clinical.registry import BIOMARKER_REGISTRY +from biolearn.clinical.convert import convert_units + +__all__ = ["BIOMARKER_REGISTRY", "convert_units"] diff --git a/biolearn/clinical/convert.py b/biolearn/clinical/convert.py new file mode 100644 index 0000000..57729dd --- /dev/null +++ b/biolearn/clinical/convert.py @@ -0,0 +1,99 @@ +"""Unit conversion utilities for clinical biomarker data.""" + +import warnings +import pandas as pd +from biolearn.clinical.registry import BIOMARKER_REGISTRY + + +def convert_units(df, source_units=None, units=None): + """Convert biomarker columns to canonical units. + + Parameters + ---------- + df : DataFrame + DataFrame with biomarker columns (samples as rows). + source_units : str, optional + Named source preset (e.g. 'ukbiobank'). Applies preset unit + mappings for all biomarkers from that source. + units : dict, optional + Per-biomarker unit overrides. Keys are biomarker names, values + are source unit strings (e.g. ``{"creatinine": "umol/L"}``). + Overrides any preset from ``source_units``. + + Returns + ------- + DataFrame + Copy of df with converted columns. + + Raises + ------ + ValueError + If a specified unit has no known conversion. + """ + unit_map = {} + if source_units is not None: + unit_map.update(BIOMARKER_REGISTRY.get_source_preset(source_units)) + if units is not None: + unit_map.update(units) + + if not unit_map: + return df.copy() + + result = df.copy() + for biomarker, source_unit in unit_map.items(): + if biomarker not in result.columns: + continue + if biomarker not in BIOMARKER_REGISTRY: + warnings.warn( + f"Biomarker '{biomarker}' not in registry, skipping conversion." + ) + continue + + entry = BIOMARKER_REGISTRY.get(biomarker) + if source_unit == entry["unit"]: + continue # already in canonical units + + if source_unit not in entry["conversions"]: + raise ValueError( + f"No conversion from '{source_unit}' to '{entry['unit']}' " + f"for biomarker '{biomarker}'. " + f"Known source units: {list(entry['conversions'].keys())}" + ) + + converter = entry["conversions"][source_unit] + result[biomarker] = result[biomarker].apply(converter) + + return result + + +def validate_ranges(df, warn=True): + """Check biomarker values against expected ranges. + + Parameters + ---------- + df : DataFrame + DataFrame with biomarker columns (samples as rows). + warn : bool + If True, emit warnings for out-of-range values. + + Returns + ------- + dict + Mapping of biomarker name to count of out-of-range values. + """ + out_of_range = {} + for col in df.columns: + if col not in BIOMARKER_REGISTRY: + continue + lo, hi = BIOMARKER_REGISTRY.valid_range(col) + mask = (df[col] < lo) | (df[col] > hi) + count = mask.sum() + if count > 0: + out_of_range[col] = int(count) + if warn: + warnings.warn( + f"Biomarker '{col}': {count} values outside expected " + f"range [{lo}, {hi}] (unit: {BIOMARKER_REGISTRY.canonical_unit(col)}). " + f"Check units." + ) + return out_of_range diff --git a/biolearn/clinical/registry.py b/biolearn/clinical/registry.py new file mode 100644 index 0000000..9d117e0 --- /dev/null +++ b/biolearn/clinical/registry.py @@ -0,0 +1,216 @@ +"""Biomarker registry defining canonical names, units, valid ranges, and conversions. + +The canonical units match NHANES conventions established in biolearn.load. +All clinical clocks expect data in these units. +""" + +_REGISTRY = { + "albumin": { + # NHANES exposes albumin via LBDSALSI in g/L, so we use that as the + # canonical unit. Typical human serum albumin is ~35-50 g/L. + "unit": "g/L", + "range": (10.0, 60.0), + "description": "Serum albumin", + "conversions": { + "g/dL": lambda x: x * 10.0, + }, + }, + "creatinine": { + # NHANES exposes creatinine via LBDSCRSI in umol/L. Typical adult range + # is ~50-110 umol/L. + "unit": "umol/L", + "range": (10.0, 1500.0), + "description": "Serum creatinine", + "conversions": { + "mg/dL": lambda x: x * 88.42, + }, + }, + "glucose": { + "unit": "mmol/L", + "range": (1.0, 40.0), + "description": "Fasting glucose", + "conversions": { + "mg/dL": lambda x: x * 0.05551, + }, + }, + "c_reactive_protein": { + "unit": "mg/dL", + "range": (0.01, 30.0), + "description": "C-reactive protein", + "conversions": { + "mg/L": lambda x: x / 10.0, + "nmol/L": lambda x: x / 95.24, + }, + }, + "white_blood_cell_count": { + "unit": "1000 cells/uL", + "range": (1.0, 50.0), + "description": "White blood cell count", + "conversions": {}, + }, + "lymphocyte_percent": { + "unit": "%", + "range": (1.0, 80.0), + "description": "Lymphocyte percentage", + "conversions": {}, + }, + "red_blood_cell_distribution_width": { + "unit": "%", + "range": (8.0, 30.0), + "description": "Red blood cell distribution width", + "conversions": {}, + }, + "mean_cell_volume": { + "unit": "fL", + "range": (50.0, 130.0), + "description": "Mean corpuscular volume", + "conversions": {}, + }, + "alkaline_phosphate": { + "unit": "U/L", + "range": (10.0, 500.0), + "description": "Alkaline phosphatase", + "conversions": {}, + }, + "hdl_cholesterol": { + "unit": "mmol/L", + "range": (0.2, 5.0), + "description": "HDL cholesterol", + "conversions": { + "mg/dL": lambda x: x / 38.67, + }, + }, + "hemoglobin": { + "unit": "g/dL", + "range": (4.0, 22.0), + "description": "Hemoglobin", + "conversions": { + "g/L": lambda x: x / 10.0, + }, + }, + "platelet_count": { + "unit": "1000 cells/uL", + "range": (10.0, 1000.0), + "description": "Platelet count", + "conversions": {}, + }, + "mean_cell_hemoglobin": { + "unit": "pg", + "range": (15.0, 45.0), + "description": "Mean corpuscular hemoglobin", + "conversions": {}, + }, + "basophil_percent": { + "unit": "%", + "range": (0.0, 10.0), + "description": "Basophil percentage", + "conversions": {}, + }, + "lymphocyte_number": { + "unit": "1000 cells/uL", + "range": (0.1, 20.0), + "description": "Lymphocyte count", + "conversions": {}, + }, + "red_blood_cell_count": { + "unit": "million cells/uL", + "range": (1.0, 10.0), + "description": "Red blood cell count", + "conversions": {}, + }, +} + +# Source unit presets for common data sources. +# +# Only sources we have actually loaded and tested against are listed here. Do +# not add a preset unless you can validate it against real data from that source +# (see ``test_load_fhs_as_geodata`` for the FHS validation pattern). +_SOURCE_PRESETS = { + # NHANES values are already in canonical units after ``load_nhanes`` runs, + # so no conversions are needed. + "nhanes": {}, + # FHS (frmgham2.csv) ships glucose in mg/dL; canonical is mmol/L. Other + # FHS biomarkers (HDL/LDL/total cholesterol) aren't populated in Period 1, + # so they're not part of the preset until we expand load_fhs to other + # periods. + "fhs": { + "glucose": "mg/dL", + }, +} + + +class BiomarkerRegistry: + """Registry of canonical biomarker definitions for clinical clocks. + + Provides lookup of units, valid ranges, and conversion functions + for all biomarkers used by clinical aging clocks. + """ + + def __init__(self, registry=None): + self._registry = registry or _REGISTRY + + def get(self, name): + """Get biomarker definition by canonical name. + + Parameters + ---------- + name : str + Canonical biomarker name (e.g. 'albumin', 'creatinine'). + + Returns + ------- + dict + Biomarker definition with keys: unit, range, description, conversions. + + Raises + ------ + KeyError + If the biomarker name is not in the registry. + """ + if name not in self._registry: + raise KeyError( + f"Unknown biomarker: '{name}'. " + f"Known biomarkers: {', '.join(sorted(self._registry.keys()))}" + ) + return self._registry[name] + + def canonical_unit(self, name): + """Return the canonical unit for a biomarker.""" + return self.get(name)["unit"] + + def valid_range(self, name): + """Return the (min, max) valid range for a biomarker.""" + return self.get(name)["range"] + + def known_biomarkers(self): + """Return sorted list of all known biomarker names.""" + return sorted(self._registry.keys()) + + def get_source_preset(self, source_name): + """Return unit mapping for a named data source. + + Parameters + ---------- + source_name : str + Source preset name (e.g. 'nhanes', 'ukbiobank'). + + Returns + ------- + dict + Mapping of biomarker name to source unit string. + """ + if source_name not in _SOURCE_PRESETS: + raise ValueError( + f"Unknown source preset: '{source_name}'. " + f"Known presets: {', '.join(sorted(_SOURCE_PRESETS.keys()))}" + ) + return _SOURCE_PRESETS[source_name] + + def __contains__(self, name): + return name in self._registry + + def __len__(self): + return len(self._registry) + + +BIOMARKER_REGISTRY = BiomarkerRegistry() diff --git a/biolearn/data/library.yaml b/biolearn/data/library.yaml index 6f4bc26..387549f 100644 --- a/biolearn/data/library.yaml +++ b/biolearn/data/library.yaml @@ -930,6 +930,43 @@ items: parse: numeric datalinks: - https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE246337 +- id: NHANES_2010 + title: NHANES 2010 Clinical Biomarkers + summary: Blood biomarker panel and linked mortality follow-up for the 2009-2010 + NHANES cycle. Loaded directly from CDC public data files. + format: clinical + organism: human + path: https://wwwn.cdc.gov/nchs/nhanes/continuousnhanes/default.aspx?BeginYear=2009 + parser: + type: nhanes + year: 2010 + datalinks: + - https://wwwn.cdc.gov/nchs/nhanes/continuousnhanes/default.aspx?BeginYear=2009 +- id: NHANES_2012 + title: NHANES 2012 Clinical Biomarkers + summary: Blood biomarker panel and linked mortality follow-up for the 2011-2012 + NHANES cycle. Loaded directly from CDC public data files. Note that C-reactive + protein is not available in this cycle. + format: clinical + organism: human + path: https://wwwn.cdc.gov/nchs/nhanes/continuousnhanes/default.aspx?BeginYear=2011 + parser: + type: nhanes + year: 2012 + datalinks: + - https://wwwn.cdc.gov/nchs/nhanes/continuousnhanes/default.aspx?BeginYear=2011 +- id: FHS + title: Framingham Heart Study (Period 1) + summary: Period 1 (1948-1953) clinical biomarkers for the Framingham Heart Study + Original Cohort. Glucose is converted from the source mg/dL into biolearn + canonical mmol/L via the "fhs" unit-conversion preset. + format: clinical + organism: human + path: https://raw.githubusercontent.com/singator/bdah/master/data/frmgham2.csv + parser: + type: fhs + datalinks: + - https://www.framinghamheartstudy.org/ - id: GSE167998 title: 'FlowSorted.BloodExtended.EPIC: An Enhanced DNA Methylation Library for Deconvoluting Peripheral Blood' diff --git a/biolearn/data_library.py b/biolearn/data_library.py index dc3eeae..b0c62f9 100644 --- a/biolearn/data_library.py +++ b/biolearn/data_library.py @@ -183,6 +183,7 @@ def __init__( rna=None, protein_alamar=None, protein_olink=None, + clinical=None, ): """ Initializes the GeoData instance. @@ -190,12 +191,15 @@ def __init__( Args: metadata (DataFrame): Metadata associated with genomic samples. dnam (DataFrame): Methylation data associated with genomic samples. + clinical (DataFrame): Clinical biomarker data with samples as rows + and biomarkers as columns. Same orientation as ``metadata``. """ self.metadata = metadata self.dnam = dnam self.rna = rna self.protein_alamar = protein_alamar self.protein_olink = protein_olink + self.clinical = clinical def _validate_metadata_omics_consistency(self): """Validate that metadata exists for all omics samples and vice versa.""" @@ -216,6 +220,10 @@ def _validate_metadata_omics_consistency(self): if self.protein_olink is not None: omics_samples.update(self.protein_olink.columns) omics_types.append("protein_olink") + if self.clinical is not None: + # clinical uses samples-as-rows, so samples live in the index + omics_samples.update(self.clinical.index) + omics_types.append("clinical") if not omics_samples: return @@ -266,6 +274,11 @@ def copy(self): if self.protein_olink is not None else None ), + clinical=( + self.clinical.copy(deep=True) + if self.clinical is not None + else None + ), ) def quality_report(self, sites=None): @@ -358,6 +371,58 @@ def from_methylation_matrix(cls, matrix): return cls(metadata, dnam) + @classmethod + def from_clinical_matrix(cls, df, source_units=None, units=None): + """Creates a GeoData instance from a clinical biomarker DataFrame. + + Splits ``df`` into a metadata frame (age, sex, mortality fields) and a + clinical frame (everything else). Both keep the standard tabular + orientation: one row per sample, one column per field. + + Parameters + ---------- + df : DataFrame + DataFrame with samples as rows and biomarkers/metadata as columns. + Index should be sample identifiers. + source_units : str, optional + Named source preset for unit conversion (e.g. ``"fhs"``). + units : dict, optional + Per-biomarker unit overrides (e.g. ``{"glucose": "mg/dL"}``). + + Returns + ------- + GeoData + Instance with clinical and metadata layers populated. Both layers + use samples as rows and fields as columns. + """ + from biolearn.clinical.convert import convert_units, validate_ranges + + df = df.copy() + + # Separate metadata columns from biomarker columns + metadata_cols = ["age", "sex", "is_dead", "months_until_death"] + existing_meta = [c for c in metadata_cols if c in df.columns] + biomarker_cols = [c for c in df.columns if c not in metadata_cols] + + metadata = ( + df[existing_meta] + if existing_meta + else pd.DataFrame(index=df.index) + ) + + clinical = df[biomarker_cols] + + # Convert units if requested + if source_units is not None or units is not None: + clinical = convert_units( + clinical, source_units=source_units, units=units + ) + + # Warn about out-of-range values + validate_ranges(clinical, warn=True) + + return cls(metadata=metadata, clinical=clinical) + def save_csv(self, folder_path, name): """ Saves the GeoData instance to CSV files according to the DNA Methylation Array Data Standard V-2410. @@ -410,6 +475,9 @@ def save_csv(self, folder_path, name): folder_path, f"{name}_protein_olink.csv" ) self.protein_olink.to_csv(protein_file) + if self.clinical is not None: + clinical_file = os.path.join(folder_path, f"{name}_clinical.csv") + self.clinical.to_csv(clinical_file) @classmethod def load_csv(cls, folder_path, name, series_part="all", validate=True): @@ -509,12 +577,20 @@ def load_csv(cls, folder_path, name, series_part="all", validate=True): else None ) + clinical_file = os.path.join(folder_path, f"{name}_clinical.csv") + clinical_df = ( + pd.read_csv(clinical_file, index_col=0, skipinitialspace=True) + if os.path.exists(clinical_file) + else None + ) + geodata = cls( metadata_df, dnam=dnam_df, rna=rna_df, protein_alamar=protein_alamar_df, protein_olink=protein_olink_df, + clinical=clinical_df, ) if validate and metadata_df is not None: @@ -1093,6 +1169,44 @@ def __init__(self, message): super().__init__(message) +class NhanesParser: + """Parser for NHANES clinical biomarker data. + + Calls ``biolearn.load.load_nhanes(year)`` and returns a GeoData with the + clinical and metadata layers populated. The ``year`` is read from the + parser config in library.yaml. + """ + + def __init__(self, data): + self.year = data.get("year") + if self.year is None: + raise ValueError("NhanesParser requires 'year' in parser config") + + def parse(self, _): + from biolearn.load import load_nhanes + + df = load_nhanes(self.year) + return GeoData.from_clinical_matrix(df) + + +class FhsParser: + """Parser for Framingham Heart Study Period 1 clinical biomarker data. + + Loads the raw FHS CSV and runs it through ``GeoData.from_clinical_matrix`` + with ``source_units="fhs"`` so glucose lands in biolearn canonical mmol/L. + """ + + def __init__(self, data): + # No config needed for Period 1; future periods can be added here. + self.data = data + + def parse(self, _): + from biolearn.load import _load_fhs_raw + + df = _load_fhs_raw() + return GeoData.from_clinical_matrix(df, source_units="fhs") + + class DataSource: """ Represents a single data source in the DataLibrary. @@ -1205,6 +1319,10 @@ def _create_parser(self, parser_data): return GisbyOlinkParser(parser_data) if parser_type == "filbin-olink": return FilbinOlinkParser(parser_data) + if parser_type == "nhanes": + return NhanesParser(parser_data) + if parser_type == "fhs": + return FhsParser(parser_data) raise ValueError(f"Unknown parser type: {parser_type}") def _show_work_needed_warning(self): diff --git a/biolearn/load.py b/biolearn/load.py index a044c82..b830202 100644 --- a/biolearn/load.py +++ b/biolearn/load.py @@ -6,14 +6,18 @@ MG_PER_DL_TO_MMOL_PER_L = 0.05551 -def load_fhs(): - """ - Loads data from the Framingham Heart Study +def _load_fhs_raw(): + """Load Framingham Heart Study Period 1 data in its native units. + + Loads the public frmgham2.csv, filters to Period 1 rows, and returns the + columns biolearn uses. Glucose is left in mg/dL (the FHS source unit) so + callers can apply the canonical conversion themselves. Returns ------- - df: Pandas.Dataframe - A pandas dataframe where each row represents an individual and each column represents a measurement about that individual + pandas.DataFrame + Indexed by ``id`` (RANDID) with columns: ``age``, ``sex``, + ``glucose`` (mg/dL), ``is_dead``, ``months_until_death``. """ public_link = "https://raw.githubusercontent.com/singator/bdah/master/data/frmgham2.csv" df = pd.read_csv( @@ -43,10 +47,24 @@ def load_fhs(): }, axis=1, ) + return df + +def load_fhs(): + """Loads data from the Framingham Heart Study (Period 1). + + Glucose is converted from the FHS source unit (mg/dL) to biolearn's + canonical mmol/L. + + Returns + ------- + pandas.DataFrame + Each row represents an individual; columns are ``age``, ``sex``, + ``glucose`` (mmol/L), ``is_dead``, ``months_until_death``. + """ + df = _load_fhs_raw() # standardize glucose units df["glucose"] = df["glucose"].apply(lambda g: g * MG_PER_DL_TO_MMOL_PER_L) - return df diff --git a/biolearn/test/test_clinical_layer.py b/biolearn/test/test_clinical_layer.py new file mode 100644 index 0000000..3a9dacc --- /dev/null +++ b/biolearn/test/test_clinical_layer.py @@ -0,0 +1,167 @@ +import os +import numpy as np +import pandas as pd +import pytest +from biolearn.data_library import GeoData + + +def _make_clinical_df(): + """Build a small DataFrame in the canonical samples-as-rows shape. + + Values are in biolearn canonical units (g/L for albumin, umol/L for + creatinine, mmol/L for glucose, etc.). + """ + return pd.DataFrame( + { + "age": [45, 62, 38], + "sex": [1, 0, 1], + "albumin": [42.0, 38.0, 45.0], + "creatinine": [80.0, 97.0, 71.0], + "glucose": [5.1, 6.2, 4.8], + "white_blood_cell_count": [6.5, 8.0, 5.5], + "lymphocyte_percent": [30.0, 25.0, 35.0], + "mean_cell_volume": [88.0, 92.0, 86.0], + "red_blood_cell_distribution_width": [12.5, 14.0, 12.0], + "alkaline_phosphate": [65.0, 80.0, 55.0], + }, + index=pd.Index(["P1", "P2", "P3"], name="id"), + ) + + +def test_geodata_accepts_clinical_layer(): + """GeoData stores clinical data with samples-as-rows.""" + clinical = pd.DataFrame( + {"albumin": [4.2, 3.8], "creatinine": [0.9, 1.1]}, + index=["P1", "P2"], + ) + metadata = pd.DataFrame({"age": [45, 62]}, index=["P1", "P2"]) + geo = GeoData(metadata=metadata, clinical=clinical) + + assert geo.clinical is not None + assert list(geo.clinical.index) == ["P1", "P2"] + assert list(geo.clinical.columns) == ["albumin", "creatinine"] + assert geo.dnam is None + + +def test_clinical_defaults_to_none(): + """The clinical layer is opt-in and defaults to None.""" + metadata = pd.DataFrame({"age": [45]}, index=["P1"]) + geo = GeoData(metadata=metadata) + assert geo.clinical is None + + +def test_from_clinical_matrix_splits_metadata_and_biomarkers(): + """from_clinical_matrix moves age/sex/mortality into metadata.""" + df = _make_clinical_df() + geo = GeoData.from_clinical_matrix(df) + + # Metadata grabs the demographic and mortality columns when present + assert "age" in geo.metadata.columns + assert "sex" in geo.metadata.columns + assert len(geo.metadata) == 3 + + # Clinical keeps samples-as-rows + assert geo.clinical is not None + assert set(geo.clinical.index) == {"P1", "P2", "P3"} + assert "albumin" in geo.clinical.columns + assert "creatinine" in geo.clinical.columns + + # Metadata fields must not leak into clinical + assert "age" not in geo.clinical.columns + assert "sex" not in geo.clinical.columns + + +def test_from_clinical_matrix_preserves_values(): + """Values stay attached to the correct sample after the split.""" + df = _make_clinical_df() + geo = GeoData.from_clinical_matrix(df) + + assert geo.clinical.loc["P1", "albumin"] == 42.0 + assert geo.clinical.loc["P2", "creatinine"] == 97.0 + assert geo.metadata.loc["P1", "age"] == 45 + + +def test_from_clinical_matrix_without_metadata_columns(): + """Works when the input is biomarkers only.""" + df = pd.DataFrame( + # Albumin in g/L, creatinine in umol/L (canonical units) + {"albumin": [42.0, 38.0], "creatinine": [80.0, 97.0]}, + index=["P1", "P2"], + ) + geo = GeoData.from_clinical_matrix(df) + + assert len(geo.metadata.columns) == 0 + assert geo.clinical is not None + assert "albumin" in geo.clinical.columns + + +def test_from_clinical_matrix_per_biomarker_unit_override(): + """The ``units`` argument applies per-biomarker conversions.""" + # Source mg/dL glucose becomes mmol/L (canonical) using factor 0.05551 + df = pd.DataFrame( + {"glucose": [90.0, 180.0]}, + index=["P1", "P2"], + ) + geo = GeoData.from_clinical_matrix(df, units={"glucose": "mg/dL"}) + + assert abs(geo.clinical.loc["P1", "glucose"] - 90.0 * 0.05551) < 1e-6 + assert abs(geo.clinical.loc["P2", "glucose"] - 180.0 * 0.05551) < 1e-6 + + +def test_from_clinical_matrix_fhs_source_preset(): + """The ``fhs`` preset converts glucose from mg/dL to mmol/L.""" + df = pd.DataFrame( + {"glucose": [90.0]}, # mg/dL in raw FHS Period 1 + index=["P1"], + ) + geo = GeoData.from_clinical_matrix(df, source_units="fhs") + + # glucose: 90 mg/dL * 0.05551 ≈ 4.9959 mmol/L + assert abs(geo.clinical.loc["P1", "glucose"] - 90.0 * 0.05551) < 1e-6 + + +def test_copy_preserves_clinical(): + """GeoData.copy() deep-copies the clinical layer.""" + df = _make_clinical_df() + geo = GeoData.from_clinical_matrix(df) + geo_copy = geo.copy() + + # Mutate the original; the copy should not change + geo.clinical.iloc[0, 0] = -999 + + assert geo_copy.clinical.iloc[0, 0] != -999 + + +def test_save_load_roundtrip_with_clinical(tmp_path): + """Clinical data survives save_csv / load_csv unchanged.""" + df = _make_clinical_df() + geo = GeoData.from_clinical_matrix(df) + + folder = str(tmp_path) + geo.save_csv(folder, "test") + + assert os.path.exists(os.path.join(folder, "test_clinical.csv")) + + loaded = GeoData.load_csv(folder, "test", validate=False) + assert loaded.clinical is not None + assert set(loaded.clinical.index) == set(geo.clinical.index) + assert set(loaded.clinical.columns) == set(geo.clinical.columns) + + pd.testing.assert_frame_equal( + loaded.clinical.sort_index(axis=0).sort_index(axis=1), + geo.clinical.sort_index(axis=0).sort_index(axis=1), + atol=1e-10, + ) + + +def test_validate_metadata_omics_includes_clinical(): + """Consistency check picks up samples that live only in the clinical layer.""" + clinical = pd.DataFrame( + {"albumin": [4.2, 3.8, 4.5]}, + index=["P1", "P2", "P3"], + ) + metadata = pd.DataFrame({"age": [45, 62]}, index=["P1", "P2"]) + geo = GeoData(metadata=metadata, clinical=clinical) + + with pytest.warns(UserWarning, match="without metadata"): + geo._validate_metadata_omics_consistency() diff --git a/biolearn/test/test_load.py b/biolearn/test/test_load.py index 0008e57..837344b 100644 --- a/biolearn/test/test_load.py +++ b/biolearn/test/test_load.py @@ -1,4 +1,5 @@ from biolearn import load +from biolearn.data_library import DataLibrary, GeoData import pytest load_columns = ["sex", "age", "glucose", "is_dead", "months_until_death"] @@ -23,6 +24,54 @@ def test_expected_error_when_loading_unsupported_year_nhanes(): df = load.load_nhanes(1913) +def test_nhanes_2010_via_datalibrary_returns_clinical_layer(): + geo = DataLibrary().get("NHANES_2010").load() + assert isinstance(geo, GeoData) + assert geo.clinical is not None + # Samples-as-rows orientation: each row is a participant + assert "age" in geo.metadata.columns + # Biomarkers live as columns on the clinical frame + assert "glucose" in geo.clinical.columns + + +def test_fhs_via_datalibrary_applies_fhs_unit_conversion(): + """End-to-end check that the FHS source preset produces canonical units.""" + geo = DataLibrary().get("FHS").load() + assert isinstance(geo, GeoData) + assert geo.clinical is not None + assert "glucose" in geo.clinical.columns + + # Raw FHS glucose is mg/dL (~70-200). After the fhs preset converts to + # canonical mmol/L the median should be in the typical 4-12 range. + glucose = geo.clinical["glucose"].dropna() + assert ( + glucose.median() < 20 + ), "Glucose median above 20 suggests conversion didn't run (still mg/dL?)" + assert ( + glucose.median() > 2 + ), "Glucose median below 2 suggests an over-conversion bug" + + +def test_load_fhs_matches_fhs_datalibrary_glucose(): + """``load_fhs`` and ``DataLibrary().get('FHS').load()`` agree on glucose.""" + df = load.load_fhs() + geo = DataLibrary().get("FHS").load() + + common_ids = df.index.intersection(geo.clinical.index) + assert len(common_ids) > 0 + + # Both paths convert mg/dL to mmol/L. They should agree to floating-point + # precision since they use the same conversion factor. + sample_id = common_ids[0] + assert ( + abs( + df.loc[sample_id, "glucose"] + - geo.clinical.loc[sample_id, "glucose"] + ) + < 1e-6 + ) + + def verify_expected_columns(df): actual_columns = set(df.columns.to_list()) missing_columns = set(load_columns) - actual_columns diff --git a/biolearn/test/test_registry.py b/biolearn/test/test_registry.py new file mode 100644 index 0000000..1aca504 --- /dev/null +++ b/biolearn/test/test_registry.py @@ -0,0 +1,130 @@ +import pandas as pd +import pytest +from biolearn.clinical.registry import BIOMARKER_REGISTRY +from biolearn.clinical.convert import convert_units, validate_ranges + + +class TestBiomarkerRegistry: + def test_known_biomarkers_not_empty(self): + assert len(BIOMARKER_REGISTRY) > 0 + + def test_get_albumin(self): + entry = BIOMARKER_REGISTRY.get("albumin") + # Canonical unit matches what NHANES exposes via LBDSALSI + assert entry["unit"] == "g/L" + assert "range" in entry + assert "conversions" in entry + + def test_get_unknown_raises(self): + with pytest.raises(KeyError, match="Unknown biomarker"): + BIOMARKER_REGISTRY.get("nonexistent_biomarker") + + def test_canonical_unit(self): + # Both glucose and creatinine canonical units match NHANES SI columns + assert BIOMARKER_REGISTRY.canonical_unit("glucose") == "mmol/L" + assert BIOMARKER_REGISTRY.canonical_unit("creatinine") == "umol/L" + + def test_valid_range(self): + lo, hi = BIOMARKER_REGISTRY.valid_range("albumin") + assert lo < hi + + def test_contains(self): + assert "albumin" in BIOMARKER_REGISTRY + assert "fake_marker" not in BIOMARKER_REGISTRY + + def test_known_biomarkers_list(self): + names = BIOMARKER_REGISTRY.known_biomarkers() + assert isinstance(names, list) + assert "albumin" in names + assert names == sorted(names) + + def test_source_preset_nhanes(self): + preset = BIOMARKER_REGISTRY.get_source_preset("nhanes") + assert isinstance(preset, dict) + # NHANES values are pre-canonicalized in load_nhanes, so no conversions + assert preset == {} + + def test_source_preset_fhs(self): + preset = BIOMARKER_REGISTRY.get_source_preset("fhs") + # FHS Period 1 ships glucose in mg/dL; canonical is mmol/L + assert preset.get("glucose") == "mg/dL" + + def test_source_preset_unknown_raises(self): + with pytest.raises(ValueError, match="Unknown source preset"): + BIOMARKER_REGISTRY.get_source_preset("fake_source") + + def test_ukbiobank_preset_is_not_registered(self): + """We removed the UK Biobank preset because we have not validated it.""" + with pytest.raises(ValueError, match="Unknown source preset"): + BIOMARKER_REGISTRY.get_source_preset("ukbiobank") + + +class TestConvertUnits: + def test_glucose_mg_per_dl_to_mmol_per_l(self): + df = pd.DataFrame({"glucose": [90.0]}, index=["P1"]) + result = convert_units(df, units={"glucose": "mg/dL"}) + assert abs(result.loc["P1", "glucose"] - 90.0 * 0.05551) < 1e-6 + + def test_hdl_cholesterol_mg_per_dl_to_mmol_per_l(self): + df = pd.DataFrame({"hdl_cholesterol": [50.0]}, index=["P1"]) + result = convert_units(df, units={"hdl_cholesterol": "mg/dL"}) + assert abs(result.loc["P1", "hdl_cholesterol"] - 50.0 / 38.67) < 1e-3 + + def test_creatinine_mg_per_dl_to_umol_per_l(self): + # Canonical is umol/L; the conversion from mg/dL multiplies by 88.42 + df = pd.DataFrame({"creatinine": [1.0]}, index=["P1"]) + result = convert_units(df, units={"creatinine": "mg/dL"}) + assert abs(result.loc["P1", "creatinine"] - 88.42) < 1e-3 + + def test_albumin_g_per_dl_to_g_per_l(self): + # Canonical is g/L; the conversion from g/dL multiplies by 10 + df = pd.DataFrame({"albumin": [4.2]}, index=["P1"]) + result = convert_units(df, units={"albumin": "g/dL"}) + assert abs(result.loc["P1", "albumin"] - 42.0) < 1e-6 + + def test_no_conversion_returns_copy(self): + df = pd.DataFrame({"albumin": [4.2]}, index=["P1"]) + result = convert_units(df) + pd.testing.assert_frame_equal(result, df) + assert result is not df # should be a different object + + def test_already_canonical_no_change(self): + # Canonical creatinine is umol/L; passing umol/L should be a no-op + df = pd.DataFrame({"creatinine": [60.0]}, index=["P1"]) + result = convert_units(df, units={"creatinine": "umol/L"}) + assert result.loc["P1", "creatinine"] == 60.0 + + def test_unknown_unit_raises(self): + df = pd.DataFrame({"creatinine": [1.0]}, index=["P1"]) + with pytest.raises(ValueError, match="No conversion"): + convert_units(df, units={"creatinine": "fake_unit"}) + + def test_missing_column_skipped(self): + df = pd.DataFrame({"albumin": [42.0]}, index=["P1"]) + result = convert_units(df, units={"creatinine": "mg/dL"}) + assert "albumin" in result.columns + + def test_fhs_source_preset(self): + df = pd.DataFrame({"glucose": [90.0]}, index=["P1"]) + result = convert_units(df, source_units="fhs") + assert abs(result.loc["P1", "glucose"] - 90.0 * 0.05551) < 1e-6 + + +class TestValidateRanges: + def test_in_range_no_warnings(self): + # Albumin canonical unit is g/L; 40 is mid-range + df = pd.DataFrame({"albumin": [40.0]}, index=["P1"]) + result = validate_ranges(df, warn=False) + assert len(result) == 0 + + def test_out_of_range_detected(self): + # 0.1 g/L is well below the (10, 60) g/L range + df = pd.DataFrame({"albumin": [0.1]}, index=["P1"]) + result = validate_ranges(df, warn=False) + assert "albumin" in result + assert result["albumin"] == 1 + + def test_unknown_columns_ignored(self): + df = pd.DataFrame({"unknown_col": [999]}, index=["P1"]) + result = validate_ranges(df, warn=False) + assert len(result) == 0 diff --git a/examples/01_composite_biomarkers/plot_load_nhanes_through_library.py b/examples/01_composite_biomarkers/plot_load_nhanes_through_library.py new file mode 100644 index 0000000..b434c36 --- /dev/null +++ b/examples/01_composite_biomarkers/plot_load_nhanes_through_library.py @@ -0,0 +1,48 @@ +""" +Loading NHANES Data Through the biolearn Library +================================================= + +This example shows the recommended way to load NHANES blood exam data into the +library. NHANES is registered in :class:`~biolearn.data_library.DataLibrary`, +so it loads the same way as every other dataset and returns a +:class:`~biolearn.data_library.GeoData` you can pass straight to a model. +""" + +############################################################################# +# Load NHANES 2010 through DataLibrary +# --------------------------------------- +from biolearn.data_library import DataLibrary + +data_source = DataLibrary().get("NHANES_2010") +data = data_source.load() + +############################################################################# +# Inspect the clinical biomarkers +# --------------------------------------- +# ``data.clinical`` has one row per participant and one column per biomarker. +print("Clinical biomarkers available:") +print(list(data.clinical.columns)) + +print("\nFirst few rows:") +print(data.clinical.head()) + +############################################################################# +# Inspect the metadata +# --------------------------------------- +# Age, sex, and mortality fields live on ``data.metadata`` so they stay +# separate from the biomarker measurements. +print("\nMetadata columns:") +print(list(data.metadata.columns)) + +print("\nFirst few rows of metadata:") +print(data.metadata.head()) + +############################################################################# +# Plot the age distribution +# --------------------------------------- +import matplotlib.pyplot as plt + +data.metadata["age"].hist(bins=30) +plt.xlabel("Age") +plt.ylabel("Number of participants") +plt.title("NHANES 2010: age distribution")