Tremor analysis

This tutorial shows how to run the tremor pipeline on prepared IMU or gyroscope data to obtain aggregated tremor measures.

Load example data

Example IMU data (8 minutes) from a participant of the Personalized Parkinson Project is loaded. The data was prepared as explained in the data preparation tutorial. The prepared data contains both accelerometer and gyroscope data, but only gyroscope data is necessary for running the tremor pipeline.

from pathlib import Path
from paradigma.util import load_tsdf_dataframe

# Set the path to where the prepared data is saved and load the data.
# Note: the test data is stored in TSDF, but you can load your data in your own way
path_to_data =  Path('../../tests/data')
path_to_prepared_data = path_to_data / '1.prepared_data' / 'imu'

df_data, metadata_time, metadata_values = load_tsdf_dataframe(path_to_prepared_data, prefix='IMU')

df_data
time accelerometer_x accelerometer_y accelerometer_z gyroscope_x gyroscope_y gyroscope_z
0 0.00000 0.550718 0.574163 -0.273684 -115.670732 32.012195 -26.097561
1 0.01004 0.535885 0.623445 -0.254545 -110.609757 34.634146 -24.695122
2 0.02008 0.504306 0.651675 -0.251675 -103.231708 36.768293 -22.926829
3 0.03012 0.488517 0.686603 -0.265550 -96.280488 38.719512 -21.158537
4 0.04016 0.494258 0.725359 -0.278469 -92.560976 41.280488 -20.304878
... ... ... ... ... ... ... ...
72942 730.74468 0.234928 -0.516268 -0.802871 0.975610 -2.256098 2.256098
72943 730.75472 0.245455 -0.514354 -0.806699 0.304878 -1.707317 1.768293
72944 730.76476 0.243541 -0.511005 -0.807177 0.304878 -1.585366 1.890244
72945 730.77480 0.240191 -0.514354 -0.808134 0.000000 -1.280488 1.585366
72946 730.78484 0.243541 -0.511005 -0.808134 -0.060976 -1.036585 1.219512

72947 rows × 7 columns

Step 1: Preprocess data

IMU sensors collect data at a fixed sampling frequency, but the sampling rate is not uniform, causing variation in time differences between timestamps. The preprocess_imu_data function therefore resamples the timestamps to be uniformly distributed, and then interpolates IMU values at these new timestamps using the original timestamps and corresponding IMU values. By setting sensor to ‘gyroscope’, only gyroscope data is preprocessed and the accelerometer data is removed from the dataframe.

from paradigma.config import IMUConfig
from paradigma.preprocessing import preprocess_imu_data

config = IMUConfig()
print(f'The data is resampled to {config.sampling_frequency} Hz.')

df_preprocessed_data = preprocess_imu_data(df_data, config, sensor='gyroscope', watch_side='left')

df_preprocessed_data

The data is resampled to 100 Hz.
time gyroscope_x gyroscope_y gyroscope_z
0 0.00 -115.670732 32.012195 -26.097561
1 0.01 -110.636301 34.624710 -24.701537
2 0.02 -103.292766 36.753000 -22.942002
3 0.03 -96.349062 38.692931 -21.175227
4 0.04 -92.585735 41.237328 -20.311531
... ... ... ... ...
73074 730.74 1.150220 -2.561552 2.440945
73075 730.75 0.588721 -1.917765 1.948620
73076 730.76 0.270257 -1.626831 1.813725
73077 730.77 0.185022 -1.451942 1.793145
73078 730.78 -0.133645 -1.112332 1.347580

73079 rows × 4 columns

Step 2: Extract tremor features

The function extract_tremor_features extracts windows from the preprocessed gyroscope data using non-overlapping windows of length config.window_length_s. Next, from these windows the tremor features are extracted: 12 mel-frequency cepstral coefficients (MFCCs), frequency of the peak in the power spectral density, power below tremor (0.5 - 3 Hz), and power around the tremor peak. The latter is not used for tremor detection, but to compute aggregate measures of tremor power in Step 4.

from paradigma.config import TremorConfig
from paradigma.pipelines.tremor_pipeline import extract_tremor_features

config = TremorConfig(step='features')
print(f'The window length is {config.window_length_s} seconds')

df_features = extract_tremor_features(df_preprocessed_data, config)

df_features

The window length is 4 seconds
time mfcc_1 mfcc_2 mfcc_3 mfcc_4 mfcc_5 mfcc_6 mfcc_7 mfcc_8 mfcc_9 mfcc_10 mfcc_11 mfcc_12 freq_peak below_tremor_power tremor_power
0 0.0 11.407993 1.854440 0.335866 0.319960 0.026209 0.071381 0.118471 0.116478 -0.029843 0.052588 -0.053722 0.012292 1.00 13277.614531 1020.094627
1 4.0 10.043115 2.011087 0.462899 0.468260 -0.110635 0.151994 -0.078799 0.112665 0.081837 0.037202 -0.034632 -0.042053 1.25 1310.703842 87.538668
2 8.0 8.877666 2.008945 0.485860 0.274370 0.064770 -0.149113 -0.130340 0.050427 0.023477 0.096070 0.074539 0.020284 3.50 215.567337 56.805816
3 12.0 10.554189 2.206457 0.528599 0.410778 0.012854 0.101787 0.085519 0.102281 0.073039 0.094114 -0.044519 0.073097 1.00 7808.868582 566.219495
4 16.0 11.051087 1.628075 0.159574 0.268279 0.110523 0.146418 0.176130 0.087930 -0.038643 0.035997 -0.017965 0.020223 1.00 3945.127377 174.899788
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
177 708.0 6.180462 1.558660 0.520806 0.407299 -0.126574 0.131888 0.040739 -0.007852 -0.034759 -0.012127 -0.048380 0.056300 2.50 4.952011 0.397532
178 712.0 8.949686 1.759538 0.146314 0.258097 0.012388 0.210405 0.096675 0.086521 0.094405 0.086109 0.081392 0.072370 1.00 570.735779 32.044727
179 716.0 6.618968 1.476363 0.332030 -0.039835 -0.085920 0.072537 0.069817 0.130373 -0.043866 0.003158 -0.052160 0.038509 1.00 3.990732 1.261435
180 720.0 6.274967 1.592159 0.146655 0.111758 -0.087665 -0.046542 0.111900 0.018244 -0.041603 -0.097095 -0.061765 0.045770 2.50 3.638188 1.217930
181 724.0 6.384021 1.438339 0.099644 0.061810 -0.184136 -0.016639 -0.052482 0.006177 0.099365 -0.059942 -0.016524 -0.010869 3.00 1.664200 1.054432

182 rows × 16 columns

Step 3: Detect tremor

The function detect_tremor uses a pretrained logistic regression classifier to predict the tremor probability (pred_tremor_proba) for each window, based on the MFCCs. Using the prespecified threshold, a tremor label of 0 (no tremor) or 1 (tremor) is assigned (pred_tremor_logreg). Furthermore, the detected tremor windows are checked for rest tremor in two ways. First, the frequency of the peak should be between 3-7 Hz. Second, we want to exclude windows with significant arm movements. We consider a window to have significant arm movement if below_tremor_power exceeds config.movement_threshold. The final tremor label is saved in pred_tremor_checked. A label for predicted arm at rest (pred_arm_at_rest, which is 1 when at rest and 0 when not at rest) was also saved, to control for the amount of arm movement during the observed time period when aggregating the amount of tremor time in Step 4 (if a person is moving their arm, they cannot have rest tremor).

from importlib.resources import files
from paradigma.pipelines.tremor_pipeline import detect_tremor

print(f'A threshold of {config.movement_threshold} deg\u00b2/s\u00b2 is used to determine whether the arm is at rest or in stable posture.')

# Load the pre-trained logistic regression classifier
tremor_detection_classifier_package_filename = 'tremor_detection_clf_package.pkl'
full_path_to_classifier_package = files('paradigma') / 'assets' / tremor_detection_classifier_package_filename

# Use the logistic regression classifier to detect tremor and check for rest tremor
df_predictions = detect_tremor(df_features, config, full_path_to_classifier_package)

df_predictions[['time','pred_tremor_proba','pred_tremor_logreg','pred_arm_at_rest','pred_tremor_checked']]

A threshold of 50 deg²/s² is used to determine whether the arm is at rest or in stable posture.
time pred_tremor_proba pred_tremor_logreg pred_arm_at_rest pred_tremor_checked
0 0.0 0.001765 0 0 0
1 4.0 0.001402 0 0 0
2 8.0 0.001146 0 0 0
3 12.0 0.001356 0 0 0
4 16.0 0.001618 0 0 0
... ... ... ... ... ...
177 708.0 0.000136 0 1 0
178 712.0 0.003822 0 0 0
179 716.0 0.001298 0 1 0
180 720.0 0.001217 0 1 0
181 724.0 0.005037 0 1 0

182 rows × 5 columns

Step 4: Compute aggregated tremor measures

The final step is to compute the amount of tremor time and tremor power with the function aggregate_tremor, which aggregates over all windows in the input dataframe. Tremor time is calculated as the number of detected tremor windows, as percentage of the number of windows while the arm is at rest or in stable posture (when below_tremor_power does not exceed config.movement_threshold). This way the tremor time is controlled for the amount of time the arm is at rest or in stable posture, when rest tremor and re-emergent tremor could occur. For tremor power the following aggregates are derived: the mode, median and 90th percentile of tremor power (specified in config.aggregates_tremor_power). The median and modal tremor power reflect the typical tremor severity, whereas the 90th percentile reflects the maximal tremor severity within the observed timeframe. The aggregated tremor measures and metadata are stored in a json file.

import pprint
from paradigma.pipelines.tremor_pipeline import aggregate_tremor

print(f'Besides the modal tremor power the following tremor power aggregates are derived: {config.aggregates_tremor_power}.')

d_tremor_aggregates = aggregate_tremor(df_predictions, config)

pprint.pprint(d_tremor_aggregates)

Besides the modal tremor power the following tremor power aggregates are derived: ['mode', 'median', '90p'].
{'aggregated_tremor_measures': {'90p_tremor_power': 0.8001855655603451,
                                'median_tremor_power': 0.4242462829851676,
                                'modal_tremor_power': 0.42,
                                'perc_windows_tremor': 16.853932584269664},
 'metadata': {'nr_windows_rest': 89, 'nr_windows_total': 182}}