This is an open-access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work, first published in JMIR mHealth and uHealth, is properly cited. The complete bibliographic information, a link to the original publication on http://mhealth.jmir.org/, as well as this copyright and license information must be included.
Transcutaneous cervical vagus nerve stimulation (tcVNS) is a promising alternative to implantable stimulation of the vagus nerve. With demonstrated potential in myriad applications, ranging from systemic inflammation reduction to traumatic stress attenuation, closed-loop tcVNS during periods of risk could improve treatment efficacy and reduce ineffective delivery. However, achieving this requires a deeper understanding of biomarker changes over time.
The aim of the present study was to reveal the dynamics of relevant cardiovascular biomarkers, extracted from wearable sensing modalities, in response to tcVNS.
Twenty-four human subjects were recruited for a randomized double-blind clinical trial, for whom electrocardiography and photoplethysmography were used to measure heart rate and photoplethysmogram amplitude responses to tcVNS, respectively. Modeling these responses in state-space, we (1) compared the biomarkers in terms of their predictability and active vs sham differentiation, (2) studied the latency between stimulation onset and measurable effects, and (3) visualized the true and model-simulated biomarker responses to tcVNS.
The models accurately predicted future heart rate and photoplethysmogram amplitude values with root mean square errors of approximately one-fifth the standard deviations of the data. Moreover, (1) the photoplethysmogram amplitude showed superior predictability (
This work furthers the state of the art by modeling pertinent biomarker responses to tcVNS. Through subsequent analysis, we discovered three key findings with implications related to (1) wearable sensing devices for bioelectronic medicine, (2) the dominant mechanism of action for tcVNS-induced effects on cardiovascular physiology, and (3) the existence of dynamic biomarker signatures that can be leveraged when titrating therapy in closed loop.
ClinicalTrials.gov NCT02992899; https://clinicaltrials.gov/ct2/show/NCT02992899
RR2-10.1016/j.brs.2019.08.002
Transcutaneous cervical vagus nerve stimulation (tcVNS) devices have emerged as inexpensive and convenient alternatives to implantable devices for stimulation of the cervical vagus nerve [
In enabling such closed-loop approaches to tcVNS, initial challenges include (1) identifying noninvasively measured biomarkers of desirable effects, (2) understanding how these biomarkers change dynamically with the delivery of tcVNS for a deepened understanding of the characteristic responses, and, accordingly, (3) predicting each subject’s responses to tcVNS for improved treatment outcomes. In current psychiatric practice, the responses to interventions are primarily assessed through patient reports and qualitative judgments of symptoms [
Unfortunately, the existing literature on such objective biomarkers is limited to invasive or obtrusive measures (eg, brain imaging or blood biomarkers) [
High-level illustration of the modeling task investigated in this study. Transcutaneous cervical vagus nerve stimulation (tcVNS) is treated as the input to the underlying physiological system, while the observed cardiovascular responses found pertinent to tcVNS effects in previous work are treated as the output signals. This work focuses on the dynamics described by linear time-invariant difference equations, formulated as discrete-time state-space equations. The latent state is depicted here as variable x, while f and g are nonlinear functions. HR: heart rate; PPG: photoplethysmogram; ECG: electrocardiogram.
As part of a collaborative study approved by the Institutional Review Boards of the Georgia Institute of Technology (#H17126), Emory University School of Medicine (#IRB00091171), SPAWAR Systems Center Pacific, and the Department of Navy Human Research Protection Program, 24 healthy human subjects, including 12 women and 12 men, with a mean age of 31 (SD 9) years, height of 173.4 (SD 8.9) cm, and weight of 77.9 (SD 13.7) kg, were recruited. All subjects had a history of prior psychological trauma, but without current posttraumatic stress disorder or other major psychiatric disorder, and written informed consent was obtained. Conforming to a randomized double-blind protocol spanning 3 days (ClinicalTrials.gov NCT02992899), each subject received four administrations—two on the first day and one on each of the following two days—of either “active” tcVNS or “sham” stimulation. These four administrations were accompanied by no other form of stimulus to focus solely on the effects of tcVNS on human physiology. Overall, 96 doses of either active or sham tcVNS were administered in total to the 24 subjects (with 12 allocated to the active group).
The active and sham devices (gammaCore, electroCore, Basking Ridge, NJ, USA) were identical in both appearance and operation, differing only with respect to the stimulation parameters. The active devices administered voltage signals consisting of five 5-kHz sinusoidal pulses repeating at a rate of 25 Hz, while the sham devices delivered biphasic square pulses at a rate of 0.2 Hz, resulting in a perceptible tingling sensation. For both active and sham devices, electrical stimulation was delivered transcutaneously to the left side of the neck, targeting the cervical vagus nerve projection. At specified times, the researcher initiated the device and adjusted the stimulation amplitude (ranging from 0 to 5 arbitrary units [AU], corresponding to 0-30 V and 0-14 V for the active and sham device, respectively) to as high as the subject could comfortably endure (active: 3.0 [SD 0.8] AU; sham: 4.4 [SD 1.2] AU). The amplitude was then kept fixed for the remainder of a 120-second period, following which the device automatically stopped, reducing its stimulation amplitude to zero. This 2-minute timeframe replicates the programmed administrations onboard the tcVNS devices that are currently in use [
Electrocardiogram (ECG) and transmissive PPG signals, taken from the finger, were continuously measured at the locations displayed in
To extract the instantaneous heart rate from the ECG signal, finite impulse response bandpass filtering (passband of 0.6-40 Hz) was applied to cancel out-of-band noise while maintaining the waveform shape [
The focus on heart rate and PPG amplitude as the biomarkers of interest for this work was based on a rationale that both cardiac and vascular downstream responses to stimulus were of interest scientifically and may indicate different autonomically mediated mechanisms following tcVNS. Heart rate is a hallmark measure of the cardiac response to changes in autonomic tone and is controlled by both the sympathetic and parasympathetic branches of the autonomic nervous system. Parasympathetic decreases in heart rate are mediated by the release of acetylcholine that binds to muscarinic receptors in the heart, whereas sympathetic increases in heart rate are mediated by the release of epinephrine and norepinephrine that bind to beta-1 receptors in the heart. PPG is a measure of peripheral blood volume pulse, and represents a surrogate measure of vasodilation or vasoconstriction resulting primarily from changes in sympathetic tone. Peripheral vasoconstriction is mediated by alpha-1 receptors in the smooth muscle of the vasculature [
Physiological sensing, signal processing, and modeling preparation. Twenty-four subjects (12 active) underwent a clinical protocol, wherein the electrocardiogram (ECG) was measured with electrodes placed in a three-lead configuration and the photoplethysmogram (PPG) was measured from the fingertip in a transmissive LED-photodiode setup. Transcutaneous cervical vagus nerve stimulation (tcVNS) or sham stimulation was administered at predefined times, where the exact stimulation location was identified by locating the left carotid pulse. Following signal processing and biomarker extraction, the biomarkers were prepared for modeling via 5-point causal averaging, resampling to 1 Hz, normalizing to rest, and finally parsing into 4 separate vectors associated with the 4 tcVNS/sham administrations studied. By referencing stimulation initiation, the corresponding input amplitude waveforms were then constructed to replicate device administration. Amp.: amplitude; D1: day 1; D2: day 2; D3: day 3.
Following feature extraction, the heart rate and PPG amplitude values existed as beat-to-beat time series of nonuniform sampling rates (due to variability in heart rate). Thus, prior to any modeling, a few time-series processing steps were employed (
To then separately investigate the effects of tcVNS/sham administration on each of the two biomarkers, the resultant time series were parsed according to the 4 administrations per subject. Based on the data available, parsing was achieved by leaving 60 seconds of data prior to each 120-second administration and retaining 120 seconds of data postadministration, for a total of 300 seconds per administration. Note that missing data at the end of the 300-second interval relevant to our analysis affected 3 administrations among the 96 collected, and therefore the corresponding data vectors were shortened accordingly.
The accompanying input data were then created for each of the 192 subject-administration-biomarker combinations (2 biomarkers × 4 administrations × 24 subjects). This was accomplished by modeling the relative tcVNS/sham amplitude delivered to each subject with pulse-like trapezoidal signals that replicated the ramping and stabilization of true stimulation. These inputs were formed by passing a boxcar input of 120-second width and unit amplitude through a 5-point moving average filter. Since stimulation amplitude remains the only variable modulated during tcVNS/sham administration, stimulation amplitude was specified as the input variable, as done in related work [
Considering the established diversity in VNS outcomes, which itself remains an active area of research [
The forthcoming equations and corresponding explanations will thus depict single-input single-output systems. The discrete-time state-space model structure in innovations form is governed by the following two difference equations:
where
In this study, we initialized model estimates using subspace identification [
Modeling, optimization, and cross-validation; quality assurance; and analysis. (Left) Modeling, optimization, and testing were performed in a leave-one-out cross-validation (CV) scheme, where, out of the four administrations, each administration was considered once as the unseen test set. The model order (M) and input delay (τ) were optimized by minimizing the small sample size-corrected Akaike information criterion (AICc). For each model order-input delay combination, a state-space model was trained by first initializing the parameter estimates using subspace estimation (N4SID); this was then followed by prediction error minimization (PEM) to refine the parameter estimates. Ridge regression was then performed for this specific model configuration by iterating lambda (λ) logarithmically over a specified interval and minimizing the mean square error. The 1-step-ahead prediction performance was evaluated using a fit percentage formula based on the root mean square error (RMSE) normalized by the standard deviation of the data (σ). (Middle) For quality assurance beyond subjective satisfaction in visual results, the models were evaluated against two objective baseline tests from the literature: the mean test and the naïve test. (Right) To extract the model information pertinent to a deepened dynamic understanding of biomarker responses, (1) optimal input delays were compiled to assess the expected response latency following stimulation onset, (2) the biomarkers were compared against each other to identify superiority in monitoring dynamic changes following transcutaneous cervical vagus nerve stimulation (tcVNS), and (3) the responses were quantifiably visualized via controlled simulation to produce population trajectories of expected dynamic changes following transcutaneous cervical vagus nerve stimulation (tcVNS) administration, in comparison to sham.
The small sample size-corrected Akaike information criterion (AICc) was employed to select the optimal model configuration (M*, τ*); the AICc was used instead of the standard AIC due to the ratio between training data points,
To optimize τ, based on a previous effort to subjectively annotate the delays between tcVNS and biomarker changes [
Once the optimal state-space model was selected for each training set of 3 administrations, the model was then regularized separately via ridge regression [
The models were evaluated against two established baselines for dynamic modeling tasks: the mean test and the naïve test. The mean test involves comparing out-of-sample 1-step prediction performance vs mean predictors, and the naïve test involves comparing 1-step prediction performance vs the naïve predictor [
The metric used herein for evaluation is the fit %, defined as:
where
demonstrating that the fit % simply replaces the mean square error and variance with the RMSE and standard deviation, respectively. It thereby exhibits improved spread for RMSE < σ/2.
To further investigate the optimal model configurations, the automatically optimized input delays, τ*, for each of the two biomarkers were compiled. These two sets of latencies were then compared against the manually annotated delays of previous work. As detailed in Gurel et al [
With 4 fit % values, model orders, and input delays obtained per subject-biomarker combination (corresponding to the 4 models produced via cross-validation), all quantities were first averaged across all 4 models prior to comparison/compilation. To examine biomarker differences in fit % and model order, pairwise
As a final analysis step, we investigated the tcVNS-induced response dynamics captured through modeling by simulating both the active and sham models forward from a controlled state. To leverage the previously resampled biomarker time series, we assembled a second set of plots corresponding to the true experimental responses. To facilitate qualitative comparisons, both sets were constructed by compiling/simulating the true/artificial biomarker time series such that 10 seconds existed prior to the true/simulated stimulus administration and 120 seconds remained afterward. The 120 seconds corresponds to the 2-minute poststimulus period used during modeling, and the 10 seconds were included to better understand the true biomarker values prior to stimulation.
For each of the 4 administrations for the 24 subjects, the heart rate and PPG amplitude time series were extracted by simply considering the values produced as a result of the modeling preparation steps. With resampled and normalized time series in hand, each biomarker’s overall response was constructed by first averaging the biomarker responses across all 4 administrations, followed by additional average and SEM calculations across all 12 subjects in each device group.
To visualize the biomarker dynamics captured by the models, each model was simulated forward by (i) setting the initial conditions to zero (ie, initializing the system at its equilibrium) to guarantee equivalence between active and sham time series (setting
The naïve test results across all subjects are shown on the right side of
Example model predictions vs true data, along with naïve test box plot results for the entire sample. (Left) Both model vs data plots shown originate from one administration of a single active transcutaneous cervical vagus nerve stimulation (tcVNS) subject. Each biomarker’s model was trained on the three other tcVNS administration datasets available for this subject; the test results on the remaining unseen dataset in a 1-step-ahead prediction task are shown. The gray dashed lines demarcate the time frame in which tcVNS was administered, and the fit percentages are calculated as previously defined. The y-axis represents relative changes from rest in percent form. (Right) The regularized, optimal state space models strongly satisfied both the naïve test and the mean test. The statistically significant (
Mean (SD) fit % and model orders for each biomarker’s state-space model.
Metric | Heart rate | PPGa amplitude |
Fit % | 76.67 (6.96) | 81.64 (7.07) |
Model order | 9.52 (0.38) | 9.39 (0.34) |
aPPG: photoplethysmogram.
Delayed effects of transcutaneous cervical vagus nerve stimulation (tcVNS) in downstream cardiovascular biomarkers. By including input delay as a free parameter necessitating optimization, the objective, quantitative state-space model optimization process described in this paper reproduced the tcVNS delay findings of the manual, subjective annotations of [
Dynamic responses to transcutaneous cervical vagus nerve stimulation (tcVNS) for heart rate (HR) and photoplethysmogram (PPG) amplitude. The curves themselves, along with their accompanying shaded regions, represent average (SEM). (Left) True responses, resampled for population averaging and subsequent dynamic modeling. The y-axis values represent relative changes from rest in percent form, and the orange dashed lines demarcate the period of active/sham stimulus administration (t ∈ (10, 130) seconds). (Right) Simulated responses to tcVNS, produced by solving the state-space model difference equations forward in time. The y-axis values represent relative changes from rest in percent form, and the orange dashed lines demarcate the period of simulated active/sham stimulus administration (t ∈ (10, 130) seconds).
Based on the results in
Biomarker characteristics that could explain this distinction involve measurement location and physiological origin. Heart rate is a cardiac measure originating centrally from the heart, whereas the PPG signals measured in this study were sensed at the periphery—in particular, transmissively at the fingertip. As described in the Methods section, a clear physiological distinction exists; these differences may in fact explain the relative contrasts observed in this study, as well as the previous difficulties encountered in identifying a trustworthy biomarker derived from the heart (eg, heart rate variability measures) [
Cardiovascular biomarkers extracted from peripheral processes, mediated solely by the sympathetic nervous system (eg, vasoconstrictor sympathetic nerve activity [
In providing quantitative evidence for the likely dominant mechanism of action for tcVNS effects on downstream physiology, we here highlight the two prevalent hypotheses for tcVNS at either the cervical or auricular branch of the vagus nerve [
The results presented herein seem to align with the latter hypothesis: namely, that the delayed effects modeled—and subjectively observed in previous work—are a byproduct of afferent vagal activity, processing in the brain, and resultant efferent-mediated autonomic effects. Synthesizing previous sensing and measurement studies of tcVNS, afferent signatures (P1-N1 vagal somatosensory evoked potentials) tend to occur within 1 second of amplitude stabilization [
Heart rate and PPG amplitude responses were visualized both for resampled data averaged across administrations and subjects, and for simulated responses produced by solving the estimated difference equations forward in time. As shown in
This issue of Rubin causality (see
A few limitations are to be noted for this study. Although nonlinear approaches to difference equation modeling are generally discouraged at the outset unless expert knowledge or sufficient evidence seems to suggest otherwise [
Although all tcVNS clinical protocols reported to date have used the gammaCore device to which the present results readily apply, future control approaches will need to venture beyond this standardized waveform and vary parameters during stimulation to achieve desired regulation goals, while simultaneously maintaining user safety. This would additionally help in satisfying the stringent input requirements for system identification [
Finally, we note that biomarkers other than PPG amplitude and heart rate were found to be statically relevant in quantifying tcVNS effects in previous work [
In this work, we studied heart rate and PPG amplitude responses to tcVNS and derived three key findings by approaching this question from a dynamic modeling perspective. First, PPG amplitude demonstrates preeminence in both modeling amenability and active vs sham response separation, suggesting its superiority as a digital biomarker for real-time response prediction and tcVNS effect quantification. Second, a consistent delay of greater than 5 seconds exists between tcVNS onset and downstream cardiovascular biomarker responses. Latency must therefore be considered and accounted for appropriately during clinical monitoring and closed-loop system design. Moreover, this delay may provide measurable evidence for the dominance of the hypothesized vagal afferent-to-brain-to-autonomic efferent pathway in downstream cardiovascular modulation. Lastly, state-space models can successfully predict heart rate and PPG amplitude dynamics in response to tcVNS, and can help to identify the characteristic dynamic signatures that separate these biomarker responses to tcVNS from sham stimulation. This dynamic modeling and analysis thereby deepens our understanding of tcVNS effects and lays the groundwork for future closed-loop approaches in time-sensitive applications.
Supplementary material.
small sample size-corrected Akaike information criterion
arbitrary units
electrocardiogram
Food and Drug Administration
photoplethysmogram
root mean square error
transcutaneous cervical vagus nerve stimulation
This work was sponsored in part by the Defense Advanced Research Projects Agency (DARPA) Biological Technologies Office (BTO) Targeted Neuroplasticity Training (TNT) program through the Naval Information Warfare Center (NIWC) Cooperative Agreement No. N66001-16-4054, and by DARPA BTO ElectRx through NIWC Cooperative Agreement No. N66001-19-2-4002. AS was funded by the National Institutes of Health (K23HL127251). JB receives funding from an electroCore Investigator Initiated Research Grant and from a Distinguished Investigator Award from the Brain and Behavior Foundation. The gammaCore stimulation devices used in this study were provided by electroCore.
JB has received funding from an electroCore Investigator Initiated Research Grant, and the gammaCore stimulation devices used in this study were provided by electroCore.