Code for "Device-embedded accelerometry complements neural signals for tracking parkinsonian motor states"

This folder contains the analysis code accompanying the manuscript *Liu, Yao, Abdi-Sargezeh, Sharma, Lasbareilles, Ho, Cheung, Denison, Tan, Neumann, Zhu, Liu, Starr, Little & Oswal — "Device-embedded accelerometry complements neural signals for tracking parkinsonian motor states".

The pipeline ingests RC+S neural and IPG-accelerometer recordings together with PKG bradykinesia/dyskinesia scores from 11 PD patients implanted with the Medtronic Summit RC+S system, extracts spectral and movement features in non-overlapping 2-minute windows, and evaluates the per-hemisphere relationship between cortico-STN neural activity, device accelerometry, and PKG-derived motor severity, both before and during continuous STN stimulation (cDBS).

Software prerequisites

Package	Tested version	Required for
MATLAB	R2020a or later	Everything
Signal Processing TB	(R2020a)	Filtering, multitaper PSD
Statistics & ML TB	(R2020a)	Random Forest, SVM, ElasticNet, cross-validation
Deep Learning TB	(R2020a)	Fully Connected Neural Network (FNN) regressor
FieldTrip	20240110	ft_freqanalysis, ft_connectivityanalysis, FOOOF
FOOOF (FieldTrip plug)	1.0+	Periodic / aperiodic spectral decomposition
SPM	SPM12	spm_eeg_ft2spm, used in spectral pre-processing
Lead-DBS	v3+	ea_space, only anal_electrodeVisualization.m
Python	3.8+	preprocessing/BehaviorAnalysis.py, preprocessing/update_segmented_PKG_files.py

See DEPENDENCIES.md for the full list of external function calls and which toolbox provides each one.

Data availability

Raw RC+S streams, accelerometer data, and PKG exports are not redistributed in this repository. Data are available from the corresponding authors (tao.liu@ndcn.ox.ac.uk, ashwini.oswal@ndcn.ox.ac.uk) upon reasonable request, subject to a data-use agreement consistent with the institutional review board approvals at UCSF and the University of Oxford.

The scripts assume the following on-disk layout (paths are configurable inside each prepro_*.m / figure script — most scripts currently use a hard-coded Z:\Tao\… Windows root that will need to be edited to match the reproducer's environment):

<DATA_ROOT>/
  ucsfData/                           # off-stim 2-min epoched neural + accel
  ucsfData_reformat/                  # reformatted off-stim sessions
  ucsfData_reformat_stim/cDBS/rawTD   # cDBS-on raw TD
  results/
    Corr2min/                         # per-hemisphere correlation .mat / .fig
    stimOnOff/                        # cDBS-on correlation .mat / .fig
    AccelFeat/                        # accel feature comparison .mat / .fig
    regressionResults/                # CV regression .mat outputs
    22Session_fooof_psd_BK.fig        # produced by anal_spectrumScores.m
    22Session_fooof_psd_DK.fig

Run order

The numbered stages below take raw RC+S data through to every figure and table in the final draft. Within a stage, scripts marked (BK) and (DK) must both be run.

1. Raw data ingest & reformat

• preprocessing/MATLAB_clean_sync/reformat_preprocess.m — off-stim TD + accel
• preprocessing/MATLAB_clean_sync/reformat_preprocess_stimTD.m — cDBS-on TD
• preprocessing/MATLAB_clean_sync/reformat_preprocess_accel.m — accel
• preprocessing/MATLAB_clean_sync/PKG_matcher.m — align PKG outputs

2. Pre-processing & 2-minute epoching

• preprocessing/prepro_getTD_beforeScore.m — neural TD
• preprocessing/prepro_getAccel_beforeScore.m — accelerometer
• preprocessing/prepro_syncContraAccelTD.m — pair contralateral accel
• preprocessing/prepro_syncContraAccelTD_stim.m — cDBS-on counterpart
• preprocessing/update_segmented_PKG_files.py — segment PKG into 2-min bins

3. Feature extraction

• preprocessing/prepro_fooofPowerCohExtraction.m — neural FOOOF + coherence (off-stim)
• preprocessing/prepro_fooofPowerCohExtraction_stim.m — same, cDBS-on
• preprocessing/prepro_aperiodicContortions_individual_full.m — aperiodic only
• preprocessing/prepro_aperiodicContortions_individual_stim_full.m — cDBS-on
• accel_featureExtractionBKDK.m — 29 movement features per 2-min window
• accel_featureExtraction_amp.m — accelerometer amplitude (Euclidean norm)
• accel_optimizeThreshold.m / accel_optimizeThreshold_ipsi.m — grid search for the time-below-threshold feature
• anal_fooofFittingR2.m / anal_fooofFittingR2_cDBS.m — FOOOF fit-quality QC (Methods §"Signal processing")

4. Correlation & feature comparison

• anal_stimOnOffCorr2min.m (BK) / anal_stimOnOffCorr2minDK.m (DK) — per-hemisphere Spearman correlations between neural features and PKG scores; produces inputs to Figure 2 (off-stim) and Figure 6A (cDBS-on).
• anal_StimOnOffCorr2minAccel_ipsi.m — Pearson correlation + KL divergence of accel features vs PKG; produces inputs to Figure 3 and Figure 6B.
• anal_stimOnOnlyCorr2minAccel_ipsi.m (BK) / anal_stimOnOnlyCorr2minAccelDK_ipsi.m (DK) — accel-vs-PKG under cDBS-on only.
• accel_featureComparison.m (BK) / accel_featureComparison_DK.m (DK) / accel_featureComparison_stim.m (cDBS) — boxplot/violin sources for Figure 3A, Figure 6B.
• anal_stimOnOffAccelComparison.m — paired stim-off vs stim-on accel feature comparison.
• anal_spectrumScores.m — high vs low symptom-severity PSD groupings used in Figure 2A.

5. Validation analysis ("acc-derived surrogate of motor state")

• anal_vaildateCorr2min.m (BK) / anal_vaildateCorr2minDK.m (DK) — off-stim, sources for Figure 5.
• anal_vaildateCorr2minAccel_ipsi.m (BK) / anal_vaildateCorr2minAccelDK_ipsi.m (DK) — same, additional accel-feature pairings. (Filename typo "vaild" preserved for reproducibility — please do not rename.)

6. Regression — off-stim (Figure 4, Table 2 cDBS-off rows)

For each combination of regressor × symptom, run the matching script:

• Random Forest: anal_contraNeural_ipsiAccelRegress_individual.m (BK), anal_contraNeural_ipsiAccelRegress_individual_DK.m (DK)
• SVM: ..._SVR.m, ..._DK_SVR.m
• ElasticNet: ..._ElasticNet.m, ..._DK_ElasticNet.m
• FNN: ..._FNN.m, ..._DK_FNN.m

Feature importance (OOB permutation, RF only):

• anal_contraNeuralIpsiAccelRegressShapley_individual.m (BK)
• anal_contraNeuralIpsiAccelRegressShapley_individual_DK.m (DK)

7. Regression — cDBS-on (Figure 6C, Table 2 cDBS-on rows)

Same set with the CDBS infix:

• Random Forest: anal_contraNeuralIpsiAccelCDBS_individual.m / ..._DK.m
• SVM: ..._SVR.m / ..._DK_SVR.m
• ElasticNet: ..._ElasticNet.m / ..._DK_ElasticNet.m
• FNN: ..._FNN.m / ..._DK_FNN.m
• Feature importance: anal_contraNeuralIpsiAccelCDBSShapley_individual.m / ..._DK.m

8. Figures (read from <DATA_ROOT>/results/)

• figure/Figure2CorrelationsPowCohBKDK.m → Figure 2
• figure/Figure3ComparisonCorrelationsAccBKDK.m → Figure 3
• figure/Figure4RegressionNeuralAccBKDK.m → Figure 4
• figure/Figure5ValidationAccBKDK_ipsi.m → Figure 5
• figure/Figure6StimOnOffCorrelationPowCohAccBKDK_cDBS.m → Figure 6
• figure/FigureS1OptimalThresholdForAccelFeat.m → Supp. Figure 1
• figure/FigureS2RegressionNeuralAccBKDK_STN_CTX.m → Supp. Figure 2
• figure/FigureS3RegressionNeuralAccBKDK_cDBS.m → Supp. Figure 3
• figure/FigureS4ValidationAccBKDK_ipsi_cDBS.m → Supp. Figure 4
• figure/RegressionAccBKDKStimOn.m, figure/RegressionFeaturesImportance.m — auxiliary panels embedded in Figures 4/6.
• Figure 1 (electrode locations) is drawn from anal_electrodeVisualization.m (requires Lead-DBS).

Script → Figure / Table mapping

Output	Producing analysis script(s)	Figure/Table script
Figure 1 (electrodes)	anal_electrodeVisualization.m	—
Figure 2 (PSD + neural-vs-PKG correlations)	anal_spectrumScores.m, anal_stimOnOffCorr2min.m, anal_stimOnOffCorr2minDK.m	figure/Figure2CorrelationsPowCohBKDK.m
Figure 3 (accel feature comparison + scatter)	accel_featureComparison.m, accel_featureComparison_DK.m, anal_StimOnOffCorr2minAccel_ipsi.m	figure/Figure3ComparisonCorrelationsAccBKDK.m
Figure 4 (RF predictions, Wilcoxon, OOB importance)	anal_contraNeural_ipsiAccelRegress_individual{,_DK}.m, anal_contraNeuralIpsiAccelRegressShapley_individual{,_DK}.m	figure/Figure4RegressionNeuralAccBKDK.m, figure/RegressionFeaturesImportance.m
Figure 5 (acc-derived surrogate vs neural)	anal_vaildateCorr2min.m, anal_vaildateCorr2minDK.m, anal_vaildateCorr2minAccel_ipsi.m, anal_vaildateCorr2minAccelDK_ipsi.m	figure/Figure5ValidationAccBKDK_ipsi.m
Figure 6 (cDBS on vs off)	anal_stimOnOffCorr2min{,DK}.m (cDBS branch), anal_stimOnOffAccelComparison.m, anal_stimOnOnlyCorr2minAccel{,DK}_ipsi.m, anal_contraNeuralIpsiAccelCDBS_individual{,_DK}.m	figure/Figure6StimOnOffCorrelationPowCohAccBKDK_cDBS.m, figure/RegressionAccBKDKStimOn.m
Supp. Figure 1 (optimal threshold)	accel_optimizeThreshold{,_ipsi}.m	figure/FigureS1OptimalThresholdForAccelFeat.m
Supp. Figure 2 (RF predictions split STN vs CTX)	anal_contraNeural_ipsiAccelRegress_individual{,_DK}.m	figure/FigureS2RegressionNeuralAccBKDK_STN_CTX.m
Supp. Figure 3 (cDBS-on RF predictions + importance)	anal_contraNeuralIpsiAccelCDBS_individual{,_DK}.m, anal_contraNeuralIpsiAccelCDBSShapley_individual{,_DK}.m	figure/FigureS3RegressionNeuralAccBKDK_cDBS.m
Supp. Figure 4 (validation under cDBS-on)	anal_vaildateCorr2minAccel{,DK}_ipsi.m (cDBS branch)	figure/FigureS4ValidationAccBKDK_ipsi_cDBS.m
Table 1 (demographics)	— (manual)	—
Table 2 (R² across regressors / feature sets / stim conditions)	All eight anal_contraNeural_ipsi* and anal_contraNeuralIpsiAccelCDBS* scripts	aggregated manually
Supp. Table 1 (stim parameters)	— (manual)	—
Supp. Table 2 (29 accel feature definitions)	— (definitions tabulated alongside accel_featureExtractionBKDK.m)	—
Supp. Table 3 (regression hyper-parameters)	— (declared inline in each regression script)	—

Notes for reviewers / re-implementers

• All correlation analyses are run per hemisphere; group-level p-values use the Benjamini–Hochberg FDR procedure implemented in fdr_bh.m.
• 5-fold stratified cross-validation with 2 repeats (10 folds total) is implemented inline inside each regression script — there is no separate CV helper.
• Some legacy filenames preserve the typo vaildate (rather than validate); these are retained intentionally so that intermediate .mat filenames continue to match.
• See PIPELINE.md for a higher-level walk-through of the pipeline stages.