Feat/issue 102/handle idle sleep mode (#143)

* Initial support for reading idle_sleep_mode Get idle_sleep_mode flag from metadata Only available for .gt3x files. Add it to watch_Data class to access later * Create idle_sleep_mode_imputation.py Initial version of the idle_sleep_mode data imputation, resamples for unevenly spaced samples * Create test_idle_sleep_mode.py Initial unit tests to check effective sampling rate is correct * Update test_idle_sleep_mode.py Adding more unit tests * Update test_idle_sleep_mode.py Check filled value is assigned Spacing for AAA testing * Update idle_sleep_mode_imputation.py Fix mypy error with this ugly monstrosity * Update test_idle_sleep_mode.py Fix to ensure numerical result from time.diff().mean() * Update orchestrator.py Adding check for idle_sleep_mode into run_file() Might be hard to hit this line for testing * Update orchestrator.py Fix ruff format error * New sample data Added new sample data for idle_sleep_mode = true smoke test * Addressing PR comments Moving idle_sleep_mode imputation to after the calibration step Imputing acceleration to exactly (0,0,-1) * Update test_orchestrator.py Path error for the full_dir processing * Update README.md Added info about handling idle_sleep_mode * Update readers.py Added info to docstring about idle_sleep_mode_flag settings * Update idle_sleep_mode_imputation.py Fix docstring for new fill value * Update pyproject.toml Update version number * Update README.md Fix zenodo badge for DOI that resolves to latest version
childmindresearch · Dec 9, 2024 · 981adb2 · 981adb2
1 parent 7dc1057
commit 981adb2
Show file tree

Hide file tree

Showing 13 changed files with 227 additions and 12 deletions.
diff --git a/README.md b/README.md
@@ -1,4 +1,4 @@
-[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.13883191.svg)](https://doi.org/10.5281/zenodo.13883191)
+[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.13883190.svg)](https://doi.org/10.5281/zenodo.13883190)
 
 # `wristpy` <img src="https://media.githubusercontent.com/media/childmindresearch/wristpy/main/docs/wristpy_logo.png" align="right" width="25%"/>
 
@@ -25,18 +25,19 @@ The package currently supports the following formats:
 | BIN | GENEActiv | GENEActiv | ✅ |
 
 **Special Note**
-    The `idle_sleep_mode` for Actigraph watches will lead to uneven sampling rates during periods of no motion (read about this [here](https://actigraphcorp.my.site.com/support/s/article/Idle-Sleep-Mode-Explained)). Consequently, this causes issues when implementing wristpy's non-wear and sleep detection. As of this moment, the authors of this package do not take any steps to impute data during these time gaps and would caution to not use data collected with this mode enabled. Of course, users can make use of the readers within wristpy for their own analysis with this type of data.
+    The `idle_sleep_mode` for Actigraph watches will lead to uneven sampling rates during periods of no motion (read about this [here](https://actigraphcorp.my.site.com/support/s/article/Idle-Sleep-Mode-Explained)). Consequently, this causes issues when implementing wristpy's non-wear and sleep detection. As of this moment, we fill in the missing acceleration data with the assumption that the watch is perfeclty idle in the face-up position (Acceleration vector = [0, 0, -1]). The data is filled in at the same sampling rate as the raw acceleration data. In the special circumstance when acceleration samples are not evenly spaced, the data is resampled to the highest effective sampling rate to ensure linearly sampled data.
 
 ## Processing pipeline implementation
 
 The main processing pipeline of the wristpy module can be described as follows:
 
-- Data loading: sensor data is loaded using [`actfast`](https://github.com/childmindresearch/actfast), and a `WatchData` object is created to store all sensor data
-- Data calibration: A post-manufacturer calibration step can be applied, to ensure that the acceleration sensor is measuring 1*g* force during periods of no motion. There are three possible options: `None`, `gradient`, `ggir`.
-- Metrics Calculation: Calculates various metrics on the calibrated data, namely ENMO (euclidean norm , minus one) and angle-Z (angle of acceleration relative to the *x-y* axis).
-- Non-wear detection: We find periods of non-wear based on the acceleration data. Specifically, the standard deviation of the acceleration values in a given time window, along each axis, is used as a threshold to decide `wear` or `not wear`.
-- Sleep Detection: Using the HDCZ<sup>1</sup> and HSPT<sup>2</sup> algorithms to analyze changes in arm angle we are able to find periods of sleep. We find the sleep onset-wakeup times for all sleep windows detected.
-- Physical activity levels: Using the enmo data (aggreagated into epoch 1 time bins, 5 second default) we compute activity levels into the following categories: inactivity, light activity, moderate activity, vigorous activity. The default threshold values have been chosen based on the values presented in the Hildenbrand 2014 study<sup>3</sup>.
+- **Data loading**: sensor data is loaded using [`actfast`](https://github.com/childmindresearch/actfast), and a `WatchData` object is created to store all sensor data
+- **Data calibration**: A post-manufacturer calibration step can be applied, to ensure that the acceleration sensor is measuring 1*g* force during periods of no motion. There are three possible options: `None`, `gradient`, `ggir`.
+- ***Data imputation*** In the special case when dealing with the Actigraph `idle_sleep_mode == enabled`, the gaps in acceleration are filled in after calibration, to avoid biasing the calibration phase.
+- **Metrics Calculation**: Calculates various metrics on the calibrated data, namely ENMO (euclidean norm , minus one) and angle-Z (angle of acceleration relative to the *x-y* axis).
+- **Non-wear detection**: We find periods of non-wear based on the acceleration data. Specifically, the standard deviation of the acceleration values in a given time window, along each axis, is used as a threshold to decide `wear` or `not wear`.
+- **Sleep Detection**: Using the HDCZ<sup>1</sup> and HSPT<sup>2</sup> algorithms to analyze changes in arm angle we are able to find periods of sleep. We find the sleep onset-wakeup times for all sleep windows detected.
+- **Physical activity levels**: Using the enmo data (aggreagated into epoch 1 time bins, 5 second default) we compute activity levels into the following categories: inactivity, light activity, moderate activity, vigorous activity. The default threshold values have been chosen based on the values presented in the Hildenbrand 2014 study<sup>3</sup>.
 
 
 ## Installation

diff --git a/pyproject.toml b/pyproject.toml
@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "wristpy"
-version = "0.1.0"
+version = "0.1.1"
 description = "wristpy is a Python package designed for processing and analyzing wrist-worn accelerometer data."
 authors = [
   "Adam Santorelli <[email protected]>",

diff --git a/src/wristpy/core/models.py b/src/wristpy/core/models.py
@@ -113,6 +113,7 @@ class WatchData(BaseModel):
     battery: Optional[Measurement] = None
     capsense: Optional[Measurement] = None
     temperature: Optional[Measurement] = None
+    idle_sleep_mode_flag: Optional[bool] = None
 
     @field_validator("acceleration")
     def validate_acceleration(cls, v: Measurement) -> Measurement:

diff --git a/src/wristpy/core/orchestrator.py b/src/wristpy/core/orchestrator.py
@@ -11,7 +11,12 @@
 
 from wristpy.core import computations, config, exceptions, models
 from wristpy.io.readers import readers
-from wristpy.processing import analytics, calibration, metrics
+from wristpy.processing import (
+    analytics,
+    calibration,
+    idle_sleep_mode_imputation,
+    metrics,
+)
 
 logger = config.get_logger()
 
@@ -350,7 +355,13 @@ def _run_file(
                 "Calibration FAILED: %s. Proceeding without calibration.", exc_info
             )
             calibrated_acceleration = watch_data.acceleration
-
+    if watch_data.idle_sleep_mode_flag:
+        logger.debug("Imputing idle sleep mode gaps.")
+        calibrated_acceleration = (
+            idle_sleep_mode_imputation.impute_idle_sleep_mode_gaps(
+                calibrated_acceleration
+            )
+        )
     enmo = metrics.euclidean_norm_minus_one(calibrated_acceleration)
     anglez = metrics.angle_relative_to_horizontal(calibrated_acceleration)
     sleep_detector = analytics.GgirSleepDetection(anglez)

diff --git a/src/wristpy/io/readers/readers.py b/src/wristpy/io/readers/readers.py
@@ -1,5 +1,6 @@
 """Function to read accelerometer data from a file."""
 
+import os
 import pathlib
 from typing import Literal, Union
 
@@ -14,6 +15,8 @@ def read_watch_data(file_name: Union[pathlib.Path, str]) -> models.WatchData:
     """Read watch data from a file.
 
     Currently supported watch types are Actigraph .gt3x and GeneActiv .bin.
+    Assigns the idle_sleep_mode_flag to false unless the watchtype is .gt3x and
+    sleep_mode is enabled (based on watch metadata).
 
     Args:
         file_name: The filename to read the watch data from.
@@ -36,13 +39,19 @@ def read_watch_data(file_name: Union[pathlib.Path, str]) -> models.WatchData:
             measurements[sensor_name] = models.Measurement(
                 measurements=sensor_values, time=time
             )
+    idle_sleep_mode_flag = False
+    if os.path.splitext(file_name)[1] == ".gt3x":
+        idle_sleep_mode_flag = (
+            data["metadata"]["device_feature_enabled"]["sleep_mode"].lower() == "true"
+        )
 
     return models.WatchData(
         acceleration=measurements["acceleration"],
         lux=measurements.get("light"),
         battery=measurements.get("battery_voltage"),
         capsense=measurements.get("capsense"),
         temperature=measurements.get("temperature"),
+        idle_sleep_mode_flag=idle_sleep_mode_flag,
     )
 
 

diff --git a/src/wristpy/processing/idle_sleep_mode_imputation.py b/src/wristpy/processing/idle_sleep_mode_imputation.py
@@ -0,0 +1,81 @@
+"""Handle idle sleep mode special case."""
+
+import numpy as np
+import polars as pl
+
+from wristpy.core import models
+
+
+def impute_idle_sleep_mode_gaps(acceleration: models.Measurement) -> models.Measurement:
+    """This function imputes the gaps in the idle sleep mode data.
+
+    Gaps in the acceleration data are filled by assuming the watch is idle in a face up
+    position. The acceleration data is filled in at a linear sampling rate, estimated
+    based on the first 100 samples timestamps, with (0, 0, -1).
+
+    In cases when the sampling rate leads to unevenly spaced samples within one second,
+    eg. 30Hz sampling rate has samples spaced at 33333333ns and 33333343ns within one
+    second, the entire data set will be resampled at the highest effective sampling rate
+    that allows for for linearly spaced samples within one second,
+    to nanosecond precision.
+
+    Args:
+        acceleration: The raw acceleration data.
+
+    Returns:
+        A Measurement object with the modified acceleration data.
+    """
+
+    def _find_effective_sampling_rate(sampling_rate: int) -> int:
+        """Helper function to find the effective sampling rate.
+
+        This function finds the new sampling rate that allows for linearly spaced
+        samples within one second, to nanosecond precision.
+
+        Args:
+            sampling_rate: The original sampling rate.
+
+        Returns:
+            The new effective sampling rate.
+        """
+        for effective_sr in range(sampling_rate, 1, -1):
+            if 1e9 % (1e9 / effective_sr) == 0:
+                return effective_sr
+        return 1
+
+    acceleration_polars_df = pl.DataFrame(
+        {
+            "X": acceleration.measurements[:, 0],
+            "Y": acceleration.measurements[:, 1],
+            "Z": acceleration.measurements[:, 2],
+            "time": acceleration.time,
+        }
+    )
+
+    sampling_space_nanosec = np.mean(
+        acceleration.time[:100]
+        .diff()
+        .drop_nulls()
+        .dt.total_nanoseconds()
+        .to_numpy()
+        .astype(dtype=float)
+    )
+
+    sampling_rate = int(1e9 / sampling_space_nanosec)
+
+    effective_sampling_rate = _find_effective_sampling_rate(sampling_rate)
+    effective_sampling_interval = int(1e9 / effective_sampling_rate)
+
+    filled_acceleration = (
+        acceleration_polars_df.set_sorted("time")
+        .group_by_dynamic("time", every=f"{effective_sampling_interval}ns")
+        .agg(pl.exclude("time").mean())
+        .upsample("time", every=f"{effective_sampling_interval}ns", maintain_order=True)
+        .with_columns(
+            pl.col("X").fill_null(value=0),
+            pl.col("Y").fill_null(value=0),
+            pl.col("Z").fill_null(value=-1),
+        )
+    )
+
+    return models.Measurement.from_data_frame(filled_acceleration)
diff --git a/tests/conftest.py b/tests/conftest.py
@@ -11,6 +11,16 @@ def sample_data_gt3x() -> pathlib.Path:
     return pathlib.Path(__file__).parent / "sample_data" / "example_actigraph.gt3x"
 
 
+@pytest.fixture
+def sample_data_gt3x_idle_sleep_mode() -> pathlib.Path:
+    """Test data for .gt3x data file."""
+    return (
+        pathlib.Path(__file__).parent
+        / "sample_data"
+        / "example_actigraph_idle_sleep_mode.gt3x"
+    )
+
+
 @pytest.fixture
 def sample_data_bin() -> pathlib.Path:
     """Test data for .bin data file."""

diff --git a/tests/sample_data/example_actigraph_idle_sleep_mode.gt3x b/tests/sample_data/example_actigraph_idle_sleep_mode.gt3x
diff --git a/tests/smoke/test_orchestrator_smoke.py b/tests/smoke/test_orchestrator_smoke.py
@@ -47,3 +47,21 @@ def test_orchestrator_different_epoch(
     assert isinstance(results.nonwear_epoch, models.Measurement)
     assert isinstance(results.sleep_windows_epoch, models.Measurement)
     assert isinstance(results.physical_activity_levels, models.Measurement)
+
+
+def test_orchestrator_idle_sleep_mode_run(
+    tmp_path: pathlib.Path,
+    sample_data_gt3x_idle_sleep_mode: pathlib.Path,
+) -> None:
+    """Idle sleep mode path for orchestrator."""
+    results = orchestrator.run(
+        input=sample_data_gt3x_idle_sleep_mode, output=tmp_path / "good_file.csv"
+    )
+
+    assert (tmp_path / "good_file.csv").exists()
+    assert isinstance(results, models.OrchestratorResults)
+    assert isinstance(results.enmo, models.Measurement)
+    assert isinstance(results.anglez, models.Measurement)
+    assert isinstance(results.nonwear_epoch, models.Measurement)
+    assert isinstance(results.sleep_windows_epoch, models.Measurement)
+    assert isinstance(results.physical_activity_levels, models.Measurement)
diff --git a/tests/unit/test_idle_sleep_mode.py b/tests/unit/test_idle_sleep_mode.py
@@ -0,0 +1,78 @@
+"""Testing the idle_sleep_mode functions."""
+
+from datetime import datetime, timedelta
+
+import numpy as np
+import polars as pl
+import pytest
+
+from wristpy.core import models
+from wristpy.processing import idle_sleep_mode_imputation
+
+
+@pytest.mark.parametrize(
+    "sampling_rate, effective_sampling_rate", [(30, 25), (20, 20), (1, 1)]
+)
+def test_idle_sleep_mode_resampling(
+    sampling_rate: int, effective_sampling_rate: int
+) -> None:
+    """Test the idle_sleep_mode function."""
+    num_samples = 10000
+    dummy_date = datetime(2024, 5, 2)
+    dummy_datetime_list = [
+        dummy_date + timedelta(seconds=i / sampling_rate) for i in range(num_samples)
+    ]
+    test_time = pl.Series("time", dummy_datetime_list, dtype=pl.Datetime("ns"))
+    acceleration = models.Measurement(
+        measurements=np.ones((num_samples, 3)), time=test_time
+    )
+
+    filled_acceleration = idle_sleep_mode_imputation.impute_idle_sleep_mode_gaps(
+        acceleration
+    )
+
+    assert (
+        np.mean(
+            filled_acceleration.time.diff()
+            .drop_nulls()
+            .dt.total_nanoseconds()
+            .to_numpy()
+            .astype(dtype=float)
+        )
+        == 1e9 / effective_sampling_rate
+    )
+
+
+def test_idle_sleep_mode_gap_fill() -> None:
+    """Test the idle_sleep_mode gap fill functionality."""
+    num_samples = 10000
+    dummy_date = datetime(2024, 5, 2)
+    dummy_datetime_list = [
+        dummy_date + timedelta(seconds=i) for i in range(num_samples // 2)
+    ]
+    time_gap = dummy_date + timedelta(seconds=(1000))
+    dummy_datetime_list += [
+        time_gap + timedelta(seconds=i) for i in range(num_samples // 2, num_samples)
+    ]
+    test_time = pl.Series("time", dummy_datetime_list, dtype=pl.Datetime("ns"))
+    acceleration = models.Measurement(
+        measurements=np.ones((num_samples, 3)), time=test_time
+    )
+    expected_acceleration = (0, 0, -1)
+
+    filled_acceleration = idle_sleep_mode_imputation.impute_idle_sleep_mode_gaps(
+        acceleration
+    )
+
+    assert len(filled_acceleration.measurements) > len(acceleration.measurements)
+    assert (
+        np.mean(
+            filled_acceleration.time.diff()
+            .drop_nulls()
+            .dt.total_nanoseconds()
+            .to_numpy()
+            .astype(dtype=float)
+        )
+        == 1e9
+    )
+    assert np.all(filled_acceleration.measurements[5010] == expected_acceleration)
diff --git a/tests/unit/test_models.py b/tests/unit/test_models.py
@@ -42,7 +42,9 @@ def test_watchdata_model() -> None:
     lux = models.Measurement(measurements=sensor_data, time=time)
     temp = models.Measurement(measurements=sensor_data, time=time)
 
-    watch_data = models.WatchData(acceleration=acceleration, lux=lux, temperature=temp)
+    watch_data = models.WatchData(
+        acceleration=acceleration, lux=lux, temperature=temp, idle_sleep_mode_flag=True
+    )
 
     if watch_data.lux is not None:
         assert np.array_equal(watch_data.lux.measurements, sensor_data)
@@ -53,6 +55,7 @@ def test_watchdata_model() -> None:
             np.array([1, 2, 3]) * 1000000,
         )
     assert isinstance(watch_data.battery, type(None))
+    assert watch_data.idle_sleep_mode_flag
 
 
 def test_measurement_model_time_type() -> None:

diff --git a/tests/unit/test_orchestrator.py b/tests/unit/test_orchestrator.py
@@ -148,6 +148,7 @@ def test_run_dir(tmp_path: pathlib.Path, sample_data_gt3x: pathlib.Path) -> None
     expected_files = {
         tmp_path / "example_actigraph.csv",
         tmp_path / "example_geneactiv.csv",
+        tmp_path / "example_actigraph_idle_sleep_mode.csv",
     }
 
     results = orchestrator.run(input=input_dir, output=tmp_path, output_filetype=".csv")

diff --git a/tests/unit/test_readers.py b/tests/unit/test_readers.py
@@ -23,6 +23,7 @@ def test_gt3x_loader(sample_data_gt3x: pathlib.Path) -> None:
     assert isinstance(watch_data.battery, models.Measurement)
     assert isinstance(watch_data.capsense, models.Measurement)
     assert watch_data.temperature is None
+    assert watch_data.idle_sleep_mode_flag is False
 
 
 def test_geneactiv_bin_loader(sample_data_bin: pathlib.Path) -> None:
@@ -34,6 +35,7 @@ def test_geneactiv_bin_loader(sample_data_bin: pathlib.Path) -> None:
     assert isinstance(watch_data.battery, models.Measurement)
     assert isinstance(watch_data.temperature, models.Measurement)
     assert watch_data.capsense is None
+    assert watch_data.idle_sleep_mode_flag is False
 
 
 def test_nonexistent_file() -> None: