Version 0.17.10

.github/workflows/actions.yml

@@ -9,7 +9,7 @@ jobs:
     runs-on: ubuntu-latest
     strategy:
       matrix:
-        python-version: [3.9, '3.10', '3.11', '3.12']
+        python-version: ['3.11', '3.12', '3.13', '3.14']
 
     steps:
     - uses: actions/checkout@v4

meson.build

@@ -1,7 +1,7 @@
 project(
     'scikit-digital-health',
     'c',
-    version: '0.17.9',
+    version: '0.17.10',
     license: 'MIT',
     meson_version: '>=1.1',
 )

pyproject.toml

@@ -47,9 +47,9 @@ classifiers = [
 
 requires-python = ">=3.9"
 dependencies = [
-    "numpy>=1.25.0",
+    "numpy>=2.0.0",
     "scipy>=1.12.0",
-    "pandas>=1.0.0",
+    "pandas>=2.0.0",
     "lightgbm>=2.3.0",
     "pywavelets",
     "scikit-learn",

src/skdh/activity/core.py

@@ -27,7 +27,7 @@
     isnan,
 )
 from numpy.linalg import norm
-from pandas import Timedelta, Timestamp
+from pandas import Timedelta, Timestamp, DataFrame
 import matplotlib
 from matplotlib.backends.backend_pdf import PdfPages
 import matplotlib.lines as mlines
@@ -82,6 +82,8 @@ class ActivityLevelClassification(BaseProcess):
         Two (2) element array-like of the base and period of the window to use for
         determining days. Default is (0, 24), which will look for days starting at
         midnight and lasting 24 hours. None removes any day-based windowing.
+    save_epoch_data : bool
+        Save the epochs with calculated metrics. Default is False.
 
     Notes
     -----
@@ -137,6 +139,7 @@ def __init__(
         min_wear_time=10,
         cutpoints="migueles_wrist_adult",
         day_window=(0, 24),
+        save_epoch_data=False,
     ):
         # make sure that the short_wlen is a factor of 60, and if not send it to
         # nearest factor
@@ -165,6 +168,7 @@ def __init__(
             min_wear_time=min_wear_time,
             cutpoints=cutpoints_,
             day_window=day_window,
+            save_epoch_data=save_epoch_data,
         )
 
         self.wlen = short_wlen
@@ -181,6 +185,13 @@ def __init__(
         else:
             self.day_key = tuple(day_window)
 
+        self.save_epochs = save_epoch_data
+        self.epoch_data = {
+            "time": [],
+            "metric": [],
+            "intensity": []
+        }
+
         # enable plotting as a public method
         self.setup_plotting = self._setup_plotting
         self._update_buttons = []
@@ -228,6 +239,39 @@ def __init__(
             for lvl in self.act_levels
         ]
 
+    def save_results(self, results, file_name):
+        """
+        Save the results of the processing pipeline to a csv file. Will also
+        save per epoch metric value if `save_epoch_data` is true. The file name 
+        for per epoch results is the same as the activity endpoints file with 
+        "_per_epoch_predictions" added to the end.
+
+        Parameters
+        ----------
+        results : dict
+            Dictionary of results from the output of predict
+        file_name : str
+            File name. Can be optionally formatted (see Notes)
+
+        Notes
+        -----
+        Available format variables available:
+
+        - date: todays date expressed in yyyymmdd format.
+        - name: process name.
+        - file: file name used in the pipeline, or "" if not found.
+        """
+        file_name = super().save_results(results, file_name)
+
+        if self.save_epochs:
+            file_name = Path(file_name)
+
+            new_name = file_name.stem + "_per_epoch_predictions" + file_name.suffix
+            epoch_file = file_name.with_name(new_name)
+
+            df = DataFrame(self.epoch_data)
+            df.to_csv(epoch_file, index=False)
+
     def add(self, endpoint):
         """
         Add an endpoint to the list to be calculated.
@@ -478,7 +522,7 @@ def predict(self, *, time, accel, fs=None, wear=None, tz_name=None, **kwargs):
 
             # compute waking hours activity endpoints
             self._compute_awake_activity_endpoints(
-                res, accel, fs, iday, dwear_starts, dwear_stops, nwlen, nwlen_60, epm
+                res, time, accel, fs, iday, dwear_starts, dwear_stops, nwlen, nwlen_60, epm
             )
             # compute sleeping hours activity endpoints
             self._compute_sleep_activity_endpoints(
@@ -517,7 +561,7 @@ def _initialize_awake_values(self, results, day_n):
                 results[endpt.name][day_n] = 0.0
 
     def _compute_awake_activity_endpoints(
-        self, results, accel, fs, day_n, starts, stops, n_wlen, n_wlen_60, epm
+        self, results, time, accel, fs, day_n, starts, stops, n_wlen, n_wlen_60, epm
     ):
         # initialize values from nan to 0.0. Do this here because days with less than
         # minimum hours should have nan values
@@ -538,6 +582,14 @@ def _compute_awake_activity_endpoints(
                 accel[start:stop], n_wlen_60, fs, **self.cutpoints["kwargs"]
             )
 
+            # handle saving the epoch data
+            self._handle_epoch_data(
+                n_wlen,
+                fs,
+                time[start:stop],
+                acc_metric,
+            )
+
             for endpoint in self.wake_endpoints:
                 endpoint.predict(
                     results, day_n, acc_metric, acc_metric_60, self.wlen, epm
@@ -740,6 +792,30 @@ def _finalize_plots(self):  # pragma: no cover
 
         pp.close()
 
+    def _handle_epoch_data(self, wlen, fs, time, metric):
+        if self.save_epochs:
+            self.epoch_data["time"].extend(
+                time[:-wlen+1:wlen]
+            )
+            self.epoch_data["metric"].extend(
+                metric
+            )
+            # handle the threshold/intensity stuff
+            intensity = full(metric.size, "sedentary")
+
+            light_lt, light_ut = get_level_thresholds("light", self.cutpoints)
+            mod_lt, mod_ut = get_level_thresholds("moderate", self.cutpoints)
+            vig_lt, _ = get_level_thresholds("vig", self.cutpoints)
+
+            light_mask = (metric >= light_lt) & (metric < light_ut)
+            mod_mask = (metric >= mod_lt) & (metric < mod_ut)
+            vig_mask = metric >= vig_lt
+
+            intensity[light_mask] = "light"
+            intensity[mod_mask] = "moderate"
+            intensity[vig_mask] = "vigorous"
+
+            self.epoch_data['intensity'].extend(intensity)
 
 class StaudenmayerClassification(BaseProcess):
     """

src/skdh/features/core.py

@@ -262,7 +262,7 @@ def compute(
         # standardize the input signal
         if isinstance(signal, DataFrame):
             columns = columns if columns is not None else signal.columns
-            x = signal[columns].values.astype(float64)
+            x = signal[columns].to_numpy(copy=True).astype(float64)
         else:
             try:
                 x = asarray(signal, dtype=float64)

src/skdh/io/base.py

@@ -106,7 +106,7 @@ def handle_naive_timestamps(time, is_local, tz_name=None):
                 .total_seconds()
             )
 
-            time += offset
+            time = time + offset
         else:  # is_local, tz_name is None
             warn(
                 "Timestamps are local but naive, and no time-zone information is available. "

src/skdh/io/csv.py

@@ -20,6 +20,7 @@
     allclose,
     nonzero,
     diff,
+    int64,
 )
 from pandas import read_csv, to_datetime
 
@@ -266,7 +267,7 @@ def handle_timestamp_inconsistency_np(self, fill_dict, time, fs, data):
             t_delta = tile(arange(0, 1, 1 / n_samples), int(n_blocks))
 
             # add the time delta so that we have unique timestamps
-            time += t_delta
+            time = time + t_delta
 
         # check if we are filling gaps or not
         if self.fill_gaps:
@@ -332,29 +333,29 @@ def predict(self, *, file, tz_name=None, **kwargs):
         # update the to_datetime_kwargs based on tz_name.  tz_name==None (utc=False)
         self.to_datetime_kw.update({"utc": tz_name is not None})
 
-        # convert time column to a datetime column. Give a unique name so we shouldnt overwrite
-        raw[self.time_col_name] = to_datetime(
+        # convert time column to a datetime column.
+        time_series = to_datetime(
             raw[self.time_col_name], **self.to_datetime_kw
         )
 
         # convert timestamps if necessary
         if tz_name is not None:
             # convert, and then remove the timezone so its naive again, but now in local time
-            raw[self.time_col_name] = raw[self.time_col_name].dt.tz_convert(tz_name)
+            time_series = time_series.dt.tz_convert(tz_name)
+
 
         # now handle data gaps and second level timestamps, etc
         # raw, fs = self.handle_timestamp_inconsistency(raw, fill_values)
 
         # get the time values and convert to seconds
-        time = (
-            raw[self.time_col_name].astype(int).values / 1e9
-        )  # int gives ns, convert to s
+        time = time_series.dt.as_unit('s').astype(int64).to_numpy(copy=True)
+        # first convert to 's' representation, then to int gives correct values
 
         data = {}
         # grab the data we expect
         for dstream in self.column_names:
             try:
-                data[dstream] = raw[self.column_names[dstream]].values
+                data[dstream] = raw[self.column_names[dstream]].to_numpy(copy=True)
             except KeyError:
                 warn(
                     f"Data stream {dstream} specified in column names but all "

src/skdh/sleep/sleep.py

@@ -47,7 +47,7 @@ class Sleep(BaseProcess):
     min_rest_block : int, optional
         Number of minutes required to consider a rest period valid. Default is 30 minutes.
     max_activity_break : int, optional
-        Number of minutes of activity allowed to interrupt the major rest period. Default is 30
+        Number of minutes of activity allowed to interrupt the major rest period. Default is 60
         minutes.
     tso_min_thresh : float, optional
         Minimum allowed z-angle threshold for determining major rest period. Default is 0.1.

tests/context/conftest.py

@@ -7,7 +7,7 @@
 def ambulation_positive_data(path_tests):
     accel = read_csv(
         path_tests / "context" / "data" / "ambulation_test1.csv", header=None
-    ).values
+    ).to_numpy(copy=True)
     time = arange(0, len(accel) / 20, 1 / 20)
 
     return time, accel
@@ -17,7 +17,7 @@ def ambulation_positive_data(path_tests):
 def ambulation_negative_data(path_tests):
     accel = read_csv(
         path_tests / "context" / "data" / "ambulation_test2.csv", header=None
-    ).values
+    ).to_numpy(copy=True)
     time = arange(0, len(accel) / 20, 1 / 20)
 
     return time, accel
@@ -27,7 +27,7 @@ def ambulation_negative_data(path_tests):
 def ambulation_negative_data_50hz(path_tests):
     accel = read_csv(
         path_tests / "context" / "data" / "ambulation_test2.csv", header=None
-    ).values
+    ).to_numpy(copy=True)
     time = arange(0, len(accel) / 50, 1 / 50)
 
     return time, accel
@@ -38,7 +38,7 @@ def motion_positive_data_100hz(path_tests):
     fs = 100
     accel = read_csv(
         path_tests / "context" / "data" / "motion_test_positive.csv"
-    ).values
+    ).to_numpy(copy=True)
     accel = vstack([accel, accel, accel])
     time = arange(0, len(accel) / fs, 1 / fs)
 
@@ -50,7 +50,7 @@ def motion_negative_data_100hz(path_tests):
     fs = 100
     accel = read_csv(
         path_tests / "context" / "data" / "motion_test_negative.csv"
-    ).values
+    ).to_numpy(copy=True)
     time = arange(0, len(accel) / fs, 1 / fs)
 
     return time, accel, fs
@@ -61,7 +61,7 @@ def motion_positive_data_20hz(path_tests):
     fs = 20
     accel = read_csv(
         path_tests / "context" / "data" / "motion_test_positive.csv"
-    ).values
+    ).to_numpy(copy=True)
     accel = vstack([accel, accel, accel])
     time = arange(0, len(accel) / fs, 1 / fs)
 
@@ -73,7 +73,7 @@ def motion_negative_data_20hz(path_tests):
     fs = 20
     accel = read_csv(
         path_tests / "context" / "data" / "motion_test_negative.csv"
-    ).values
+    ).to_numpy(copy=True)
     time = arange(0, len(accel) / fs, 1 / fs)
 
     return time, accel, fs

tests/gait/data/update_from_sheet_apcwt.py

@@ -56,8 +56,8 @@
     # get keys to update results for
     keys_turn = keys + ["Turn"]
 
-    tmp = {k: df.loc[:, k].values for k in keys}
-    tmp_turn = {k: df_turn.loc[:, k].values for k in keys_turn}
+    tmp = {k: df.loc[:, k].to_numpy(copy=True) for k in keys}
+    tmp_turn = {k: df_turn.loc[:, k].to_numpy(copy=True) for k in keys_turn}
 
     np.savez("gait_results_apcwt.npz", **tmp)
     np.savez("gait_results2_apcwt.npz", **tmp_turn)

tests/gait/data/update_from_sheet_vcwt.py

@@ -88,12 +88,12 @@
     # old = np.load("gait_results_vcwt.npz")
     # old_turn = np.load("gait_results2_vcwt.npz")
 
-    tmp = {k: df.loc[:, k].values for k in keys}
+    tmp = {k: df.loc[:, k].to_numpy(copy=True) for k in keys}
     # tmp.update(
     #     {k: np.full(tmp["Bout N"].size, np.nan) for k in old.files if k not in tmp}
     # )
 
-    tmp_turn = {k: df_turn.loc[:, k].values for k in keys_turn}
+    tmp_turn = {k: df_turn.loc[:, k].to_numpy(copy=True) for k in keys_turn}
     # tmp_turn.update(
     #     {
     #         k: np.full(tmp_turn["Bout N"].size, np.nan)

tests/gait_old/data/update_from_sheet.py

@@ -62,12 +62,12 @@
     old = np.load("gait_results.npz")
     old_turn = np.load("gait_results2.npz")
 
-    tmp = {k: df.loc[:, k].values for k in keys}
+    tmp = {k: df.loc[:, k].to_numpy(copy=True) for k in keys}
     # tmp.update(
     #     {k: np.full(tmp["Bout N"].size, np.nan) for k in old.files if k not in tmp}
     # )
 
-    tmp_turn = {k: df_turn.loc[:, k].values for k in keys_turn}
+    tmp_turn = {k: df_turn.loc[:, k].to_numpy(copy=True) for k in keys_turn}
     # tmp_turn.update(
     #     {
     #         k: np.full(tmp_turn["Bout N"].size, np.nan)

tests/io/test_read_base.py

@@ -28,21 +28,19 @@ def test_in_pipeline(dummy_reader_class, dummy_process):
 def test_handle_naive_timestamps():
     # get a timestamp array of naive timestamps
     ts = (
-        date_range("2023-11-05 00:00:00", "2023-11-05 05:00:00", freq="0.1h").view(
-            "int64"
-        )
-        / 1e9
+        date_range(
+            "2023-11-05 00:00:00", "2023-11-05 05:00:00", freq="0.1h"
+        ).as_unit('s').view("int64")
     )
     # get the true timestamps - note that with the DST change there is an extra hour here
     # hence the earlier end time by 1 hour
     ts_true = (
         date_range(
-            "2023-11-05 00:00:00", "2023-11-05 04:00:00", freq="0.1h", tz="US/Eastern"
-        ).view("int64")
-        / 1e9
+            "2023-11-05 00:00:00", "2023-11-05 04:00:00", freq="0.1h", tz="America/New_York"
+        ).as_unit('s').view("int64")
     )
 
-    ts_pred = handle_naive_timestamps(ts, is_local=True, tz_name="US/Eastern")
+    ts_pred = handle_naive_timestamps(ts, is_local=True, tz_name="America/New_York")
 
     assert allclose(ts_pred - ts_true[0], ts_true - ts_true[0])