Abstract
The dorsal anterior cingulate cortex (dACC) is a critical brain area for pain and autonomic processing, making it a promising noninvasive therapeutic target. We leverage the high spatial resolution and deep focal lengths of low-intensity focused ultrasound (LIFU) to noninvasively modulate the dACC for effects on behavioral and cardiac autonomic responses using transient heat pain stimuli. A N = 16 healthy human volunteers (6 M/10 F) received transient contact heat pain during either LIFU to the dACC or Sham stimulation. Continuous electroencephalogram (EEG), electrocardiogram (ECG), and electrodermal response (EDR) were recorded. Outcome measures included pain ratings, heart rate variability, EDR response, blood pressure, and the amplitude of the contact heat-evoked potential (CHEP).
LIFU reduced pain ratings by 1.09 ± 0.20 points relative to Sham. LIFU increased heart rate variability indexed by the standard deviation of normal sinus beats (SDNN), low-frequency (LF) power, and the low-frequency/high-frequency (LF/HF) ratio. There were no effects on the blood pressure or EDR. LIFU resulted in a 38.1% reduction in the P2 CHEP amplitude. Results demonstrate LIFU to the dACC reduces pain and alters autonomic responses to acute heat pain stimuli. This has implications for the causal understanding of human pain and autonomic processing in the dACC and potential future therapeutic options for pain relief and modulation of homeostatic signals.
Significance Statement
New lines of inquiry demonstrate autonomic signals like heart rate variability (HRV) are aberrant in chronic pain and mental health disorders, which may contribute to their underlying etiology. The dorsal anterior cingulate cortex (dACC) is a key center in pain processing with direct influences on autonomic function, but its depth precludes direct access without invasive surgery. For the first time in humans, we demonstrate a low-intensity focused ultrasound (LIFU) can noninvasively and selectively modulate the dACC to reduce acute pain perception, autonomic responses, and pain processing signals. This work further establishes a causal role of the dACC in pain and autonomic processing with potential future clinical applications in chronic pain and neuropsychological populations.
Introduction
The spinothalamic tract conveys information on the physiologic state of all tissues in the body through the thalamus to cortical regions including the insula and cingulate cortex (Melzack, 1999; Craig, 2002, 2004; Dum et al., 2009). Following a painful stimulus, nociceptive signals transmit along this pathway and integrate with autonomic information in brainstem regions like the ventrolateral medulla, solitary nucleus, parabrachial nucleus, and periaqueductal gray (Yezierski, 1988; Craig, 1995; Dum et al., 2009). This integrated information is sent through the mediodorsal thalamus to the dorsal anterior cingulate cortex (dACC), making it a critical region for pain and autonomic processing and a promising noninvasive therapeutic target (Beissner et al., 2013; Wager et al., 2013; Vogt, 2016; Seeley, 2019). Increased activity in the dACC is associated with pain and emotional processing (Hutchison et al., 1999; Shackman et al., 2011; Hayes and Northoff, 2012) and encodes pain intensity and affect (Rainville et al., 1997; Wager et al., 2013; Xiao and Zhang, 2018). Overactivity in the dACC appears to be involved in the onset and maintenance of chronic pain conditions (Bliss et al., 2016; Vanneste and De Ridder, 2021).
dACC and autonomic function
Pain is a potent autonomic challenge, and the dACC integrates cognitive processes, pain, and autonomic control (Shackman et al., 2011; Hayes and Northoff, 2012). Direct electrical stimulation of the dACC generates autonomic and behavioral responses like tachycardia and heat sensations (Xue et al., 2023). Activity in the dACC is associated with low-frequency (LF) power (Critchley et al., 2003; Rebollo et al., 2018), a measure of heart rate variability that may index baroreflexes activity (Moak et al., 2007; Goldstein et al., 2011; Reyes Del Paso et al., 2013). Baroreflexes look to play a role in pain suppression (Bruehl et al., 1992; Dworkin et al., 1994; Bruehl and Chung, 2004; Chung et al., 2008; Suarez-Roca et al., 2019), and the dACC is implicated in baroreflex activity (Gianaros et al., 2012; Ginty et al., 2013). The dACC is thus a critical hub for autonomic function and pain processing, suggesting these roles are inextricably linked in chronic pain populations (Medford and Critchley, 2010; Seifert et al., 2013).
Noninvasive dACC modulation and low-intensity focused ultrasound
Noninvasive neuromodulation techniques like transcranial magnetic stimulation (TMS) and transcranial electrical stimulation (TES) have previously attempted to target the dACC with mixed results (Galhardoni et al., 2019; Khan et al., 2020), likely because they do not confer the spatial resolution or depth penetration to selectively reach targets deeper than approximately 3.5 cm in the brain (Deng et al., 2013; Gomez-Tames et al., 2020). Low-intensity focused ultrasound (LIFU) is a neuromodulatory approach that uses mechanical energy to nondestructively and reversibly modulate neuronal activity with a millimeter-scale spatial resolution and adjustable depth of focus suitable for targeting deeper brain targets with high spatial precision (Legon et al., 2014, 2018a; Blackmore et al., 2019; Darmani et al., 2022). This is in contrast to high-intensity focused ultrasound (HIFU), which uses high energy levels to thermally ablate neural tissue or the use of ultrasound with contrast agents to open the blood–brain barrier (Lipsman et al., 2013; Meng et al., 2021). LIFU has been used for cortical and sub-cortical neuromodulation in humans (Legon et al., 2014, 2018a,b; Cain et al., 2021; Nakajima et al., 2022) and the side effect profile has previously been shown to be similar to other forms of noninvasive neuromodulation (Legon et al., 2020; Sarica et al., 2022).
This paper examines how modulation of the dACC with LIFU affects subjective report of pain, autonomic responses, and the contact heat evoked potential (CHEP) that is a widely studied brain marker of pain processing (Chen et al., 2001, 2006; Granovsky et al., 2016; De Schoenmacker et al., 2021) sourced to the ACC (Valeriani et al., 2002; Shenoy et al., 2011; Jin et al., 2018). We hypothesize LIFU will confer an inhibitory effect (Legon et al., 2014, 2023) to reduce pain, increase HRV, and attenuate the CHEP amplitude.
Materials and Methods
Participants
A total of 16 healthy participants (6 M/10 F, aged 28 ± 5.46 years) met inclusion criteria and were enrolled in the study. Preliminary data investigating LIFU to the anterior insula on pain and heart rate variability using the CHEP (Legon et al., 2023) suggested that, at a ß = 0.8 and p = 0.05, we needed an N = 15 to detect a difference in subjective pain report of 1 point and an N = 16 to detect an effect size of 0.75 in HRV. As both were primary outcome measures, we collected an N = 16. Inclusion criteria were males and females aged 18–65. Exclusion criteria were in accordance with contraindications to noninvasive neuromodulation as outlined by Rossi et al. (2009) for transcranial magnetic stimulation in addition to contraindications to MRI and CT; pregnancy, an active medical disorder or treatment with potential CNS effects or a history of neurologic disorder (e.g., multiple sclerosis, Parkinson’s disease, epilepsy, etc.), a history of head injury resulting in loss of consciousness for >10 min, and/or a history of alcohol or drug dependence. All procedures were reviewed and approved by the Virginia Tech Institutional Review Board.
Overall study design
Data were collected in three sessions on three separate days with a minimum of 24 h between sessions. Session 1 was an imaging visit consisting of an anatomical computerized tomography (CT) and structural magnetic resonance imaging (MRI). Sessions 2 and 3 were formal testing days where participants were counterbalanced for Sham or real LIFU to the dACC.
Imaging acquisition
MRI data were acquired on a Siemens 3 T Prisma scanner (Siemens Medical Solutions,) at the Fralin Biomedical Research Institute’s Human Neuroimaging Laboratory. Anatomical scans were acquired using a T1-weighted MPRAGE sequence with a TR = 1,400 ms, TI = 600 ms, TE = 2.66 ms, flip angle = 12°, voxel size = 0.5 × 0.5 × 1.0 mm, FoV read = 245 mm, FoV phase of 87.5%, 192 slices, and ascending acquisition.
CT scans were collected with a Kernel = Hr60 in the bone window, FoV = 250 mm, kilovolts (kV) = 120, rotation time = 1 s, delay = 2 s, pitch = 0.55, caudocranial image acquisition order, 1.0 mm image increments for a total of 121 images and total scan time of 13.14 s.
Experimental design
Overview of testing procedures and timeline
Prior to formal testing for sessions 2 and 3, participants first completed a baseline report of symptoms (ROS) questionnaire that asks about the presence and severity of symptoms (sleepiness, headache, etc.) as previously used in our human LIFU studies (Legon et al., 2020) (see Appendix for Questionnaire). They were then seated in a comfortable chair with arm and neck support and connected to continuous electroencephalogram (EEG), electrocardiogram (ECG), and electrodermal response (EDR) monitoring (see below for details). Pain thresholding (see below), baseline blood pressure, and a 5 minute (min) baseline EEG, ECG, and EDR were then performed. After the 5 min baseline, 10 min of formal contact heat pain testing with LIFU time locked to each stimulus was performed. After pain testing, participants rested for 5 min. Continuous EEG, ECG, and EDR data were collected throughout all time windows. After 30 min post LIFU stimulation, participants completed another ROS questionnaire and an Auditory Masking Questionnaire (AMQ). Details of each procedure are below. The study outline is depicted in Figure 1.
Outline of acute contact heat pain testing. A, Participants completed two arms on separate days: inactive Sham and low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate cortex (dACC). For each session, data collection was divided into three time windows: A 5 min baseline, 10 min of acute contact heat testing, and a 5 min rest period. Contact heat pain testing consisted of 40 stimuli (see B) delivered to the dorsum of the right hand with a random inter-stimulus interval (ISI) of 10–20 s. Continuous electroencephalogram (EEG), electrocardiogram (ECG), and electrodermal response (EDR) were recorded throughout all time windows. Auditory masking was delivered continuously over the contact heat testing period for both conditions. B, Blow-up of orange bar in A showing details of contact heat stimulus and timing and parameters of low-intensity focused ultrasound (LIFU). The contact heat stimulus lasted 300 msec and was titrated to a 5/9 on a numerical pain rating scale. One second of LIFU was time locked to the contact heat stimulus, beginning 200 milliseconds (ms) prior to contact heat (0 ms) and ending 800 ms after. LIFU was delivered with a fundamental frequency of 500 kHz, a pulse repetition frequency (PRF) of 1 kHz, a duty cycle (DC) of 36%, and a total burst duration of 1 s. Refer to Figure 9 for individual subject results of the Auditory Masking Questionnaire (AMQ). Refer to Figure 10 for the Report of Symptoms questionnaire that illustrates reported adverse events for both LIFU and Sham conditions.
Data acquisition
Electroencephalogram (EEG)
Surface EEG was collected from channels CZ and PZ. The vertex (CZ) has previously been shown as the best site for recording the CHEP (Valeriani et al., 2002; Granovsky et al., 2005; Chen et al., 2006; Kramer et al., 2012; De Schoenmacker et al., 2022). Data were continuously acquired at a 1 kHz sampling rate using a DC amplifier (GES 400, Magstim EGI, OR, USA) and Net Station 5.4 EEG software with two 10 mm silver-over-silver chloride cup electrodes, referenced to the bilateral mastoid processes. The scalp for each site was first prepared with a mild abrasive gel (Nuprep; Weaver and Company) and rubbing alcohol. Cup electrodes were filled with a conductive paste (Ten20 Conductive; Weaver and Company) and held in place with medical tape. Electrode impedances were verified (<50 kΩ) before recording. Data were stored on a PC for offline data analysis.
Electrocardiogram (ECG)
Two latex-free ECG electrodes (MedGel, MDSM611903) were attached symmetrically to the anterior surface of the bilateral forearms immediately distal to the antecubital fossa. ECG data were continuously collected and sampled at 1 kHz using the same 64-channel EEG recording system with integrated grounding (GES 400, Magstim EGI) with a physiologic data acquisition box (Physio16, Magstim, EGI) and Net Station 5.4 EEG software and stored on a PC for offline data analysis.
Electrodermal response (EDR)
EDR data were collected using the Consensys GSR Development Kit (Shimmer). It was placed on the right wrist with a photoplethysmogram sensor on the proximal fifth digit and skin conductance electrodes on the proximal second and third digits. EDR data were continuously collected and sampled at 128 Hz using the Consensys software and stored on a PC for offline data analysis.
Blood Pressure (BP)
Systolic and diastolic BPs were collected using the QardioArm wireless portable blood pressure monitor placed on the upper left arm. BP measurements were taken twice, once at the beginning of each session prior to any intervention or task (pre-session BP), and once at the end of each session after all session data were collected (post-session BP).
Contact heat pain stimulus
The contact heat pain stimulus was delivered using the TCS (thermal cutaneous stimulator) system (QST.Lab) and a contact 3 × 3.2 × 2.4 mm peltier device (T-03). The peltier device used for CHEP stimulus delivery is capable of delivering rapid ramping heat stimuli that activates peripheral Aδ-fibers but not C-fibers (Truini et al., 2007; Rosner et al., 2018) and is a validated way of studying nociception and pain perception in both healthy (Chen et al., 2001, 2006; Granovsky et al., 2016) and chronic pain populations (Lenoir et al., 2020). This QST device is a validated tool to elicit the CHEP (De Schoenmacker et al., 2022).
Pain thresholding
Prior to formal pain thresholding, participants were familiarized with the QST probe and the heat stimulus. Formal thresholding consisted of delivering three consecutive stimuli (as above starting at 55°C) at a fixed ISI of 5 s. Participants were required to rate the pain intensity of each stimulus on a 0–9 numerical rating scale where 0 = no pain; 1–3 = mild pain; 4–6 = moderate pain; 7–9 = severe pain. The average of the three stimuli was used as their threshold. The temperature was raised in 1°C increments until a 5/9 average rating was achieved. Thresholding was performed this way because pilot data indicated that thresholding using a continuous graded stimulus did not correlate with thresholds using the discrete 300 ms stimulus. We found using the actual CHEP stimulus as the thresholding stimulus produced more consistent pain ratings during formal testing.
Contact heat testing
Six points were first identified in a 3 × 2 cm grid (1 cm spacing) on the dorsum of the right hand to guide the location of the heat probe on the hand. The grid was positioned in the radial nerve distribution for all participants. Each CHEP stimulation consisted of a stimulus duration of 0.3 s (50 ms ramp up, 200 ms plateau, 50 ms ramp down) at a ramping speed of 300°C/s. In between each stimulus there was a random inter-stimulus interval (ISI) of 10–20 s. The probe was moved to the next spot on the 3 × 2 grid following each stimulus in order to minimize repeated stimuli at a single site and prevent habituation or sensitization to the CHEP stimulus at any one hand site. A total of 40 stimuli were delivered. Participants were required to rate the pain intensity on the 0–9 numerical rating scale using a numerical keypad in their left hand. There was no time restriction placed on the response. Pain ratings were recorded using a custom MATLAB script using MATLAB R2022a (The MathWorks, Inc.).
LIFU transducer and waveform
We used a single-element focused ultrasound transducer (Sonic Concepts, model H-104) with a center frequency of 0.5 MHz, 64 mm aperture, 52 mm focal length from the exit plane, and an f# of 0.81. Transcranial ultrasonic neuromodulation waveforms were generated using a 30 MHz dual channel function generator (4078B; BK Precision Instruments). Channel 1 was a 5Vp-p square wave burst of 1 kHz (N = 1,000) with a 36% duty cycle used to gate channel 2 that was a 500 kHz sine wave. This resulted in a 1 s waveform of 1,000 pulses of 360 microseconds (µs). The output of channel 2 was sent through a 100 W linear RF amplifier (E&I 2100L; Electronics & Innovation) before being sent to the transducer through its corresponding matching network. LIFU delivery was time-locked to the CHEP stimulus, beginning at −200 ms (prior to stimulus) and ending at +800 ms after heat stimulus delivery. This time-locked duration was chosen due to sex and age variability in CHEP latency (Granovsky et al., 2016), ensuring the CHEP was generated within the LIFU delivery.
Empirical acoustic field mapping
The acoustic intensity profile of the ultrasound waveform was measured in an acoustic test tank filled with deionized, degassed, and filtered water (Precision Acoustics Ltd.). A calibrated hydrophone (HNR-0500, Onda Corp.) mounted on a motorized stage was used to measure the pressure from the ultrasound transducer in the acoustic test tank. The ultrasound transducer was positioned in the tank using a custom setup and leveled to ensure the hydrophone was perpendicular to the surface of the transducer. XY and YZ planar scans were performed at an isotropic 0.25 × 0.25 mm resolution. A voltage sweep from input voltages of 20–250 mVpp in 10 mVpp increments was also performed to determine the necessary input voltage to obtain the desired extracranial pressure.
LIFU targeting
The Montreal Neurologic Institute (MNI) XYZ coordinates (0,18,30) were chosen based upon Neurosynth reverse inference software for the term “pain” as well as prior literature on the functional magnetic resonance imaging (fMRI) signature of heat pain (Yarkoni et al., 2011; Wager et al., 2013). This target was sufficiently anterior to allow placement of the EEG CZ electrode along with the transducer so that they did not overlap. LIFU targeting was aided using a neuronavigation system (Brainsight, Rogue Research). Individual subject T1-weighted structural MRIs were registered in the Brainsight neuronavigation software. Participant dACC targets were transformed from standard MNI space to individual space using the software’s native transformation so that the coordinates applied to each subject’s individual anatomy. In one subject (male subject 4 in Table 1), the native transformation placed the original MNI coordinate on the corpus callosum, so a secondary site (MNI, 0,26,26) was chosen as this allowed direct targeting of gray matter and is also a common activation site for thermal pain in Neurosynth (Yarkoni et al., 2011; Lieberman and Eisenberger, 2015; Wager et al., 2016). The LIFU transducer was coupled to the head by first parting the hair using standard ultrasound gel and then using an acoustically transparent gel puck (Strohman et al., 2023). Puck thickness was selected for each participant based on the distance from the scalp to the dACC target (Table 1) in order to individualize axial (depth) targeting between participants. The optimal positioning of the transducer was midline at the top of the head anterior to the vertex for all participants. The transducer was secured to the head and positioning was verified online continuously throughout the experiment using the neuronavigation system’s native real-time tracking software.
Depth to dACC target and estimated intracranial pressure
Auditory masking
LIFU is known to produce audible sounds that may confound the data, but these can be effectively masked (Braun et al., 2020; Liang et al., 2023). For auditory masking, participants were given headphones connected to a tablet with a multitoned white noise generator. The sound options were a series of environmental sounds such as thunderstorms, traffic, wind chimes, etc. Participants were asked to select their sounds, including layering them on top of each other, and instructed to turn the volume up until they could not hear normal conversation. This was verified by querying them by standing behind them and verbally asking “can you hear me?” at a normal conversational level. A lack of response indicated an adequate volume level. The mask was constantly applied throughout the entire period of LIFU or Sham application. After each session, participants were queried on the effectiveness of auditory masking using the AMQ.
Auditory Masking Questionnaire
The AMQ was administered on both LIFU and Sham sessions after all testing was completed. It consisted of three questions to evaluate the success of auditory masking: “I could hear the LIFU stimulation,” “I could feel the LIFU stimulation,” and “I believe I experienced LIFU stimulation.” Each question was answered on a 7-point Likert scale where 1 = “Strongly Disagree,” 2 = “Disagree,” 3 = Slightly Disagree,”4 = “Unsure,” 5 = “Slightly Agree,” 6 = “Agree,” and 7= “Strongly Agree.”
Blinding
This was a single-blind study. It was not possible (at the time of collection) to blind the experimenters administering the intervention. However, the researchers responsible for data analysis were blinded which follows recommendations to reduce bias for investigations of noninvasive brain stimulation on autonomic function (Schmaußer et al., 2022). The AMQ was used to verify blinding of participants.
Report of symptoms questionnaire
The report of symptoms (ROS) asks about the presence of symptoms including severity (absent, mild, moderate, severe) and has been previously used in our human LIFU studies. All symptoms are scored on a 0–3 scale, with a 0 indicating absence of the symptom up to a 3 indicating a severe symptom. The ROS questionnaire collected before and after LIFU or Sham was summed across participants to provide a total count of absent, mild, moderate, and severe intensity reports for each symptom. Counts of symptoms after LIFU/Sham were adjusted by taking the pre–post LIFU/Sham difference to provide an indication of any change of symptoms.
Data preprocessing
Electroencephalogram (EEG)
EEG data were band-pass filtered (2–100 Hz) using a third-order Butterworth filter. Data were manually inspected for artifact (eye blink, muscle activity), and contaminated epochs were removed. Two time windows were then extracted from the data for analysis: A 5 min resting baseline and the CHEP pain testing (∼10 min). These time windows are represented in Figure 1 as the baseline and contact heat testing periods. To examine LIFU effects on the CHEP, the data within the contact heat testing time window were epoched from −2000 to 2000 msec around each of the 40 trials and baseline corrected (−1,500 to −500 ms).
Electrocardiogram
ECG data were filtered from 10–30 Hz using a third order Butterworth filter. Data were then manually inspected for movement artifact and removed. R-R peaks were found using the automated detection program findpeaks in MATLAB and confirmed manually. Data were then divided into the same two epochs as the EEG data (baseline, contact heat testing).
Electrodermal response
EDR data were zero meaned and filtered using a bandpass (0.1–25 Hz) third-order Butterworth filter. Data were manually inspected for movement artifact and removed. Data were divided into the same two epochs as the EEG data (baseline, contact heat testing). Within the contact heat testing window, EDR data were further epoched around the onset of the CHEP stimulus (−5 to 10 s), mean subtracted, and averaged for the 40 CHEPs trials for each condition for each individual.
Statistical analysis
For all variables of interest, Bartlett’s test was used to assess homoscedasticity. Variables that violated homogeneity of variances were evaluated using nonparametric statistical tests.
Pain ratings
Group average pain ratings
Pain ratings were averaged across the 40 trials for each participant for each condition (LIFU, Sham). A total of 11 pain ratings (spread across six participants and both LIFU and Sham conditions) were missing, representing 0.86% of the data. There was no discernible pattern to the missing data indicative of experimental error, suggesting the missing data were due to either an erroneous keypress, multiple keypresses, or a lack of response on some trials. The maximum number of missing pain ratings within a condition for any participant was three. To account for this, missing values were replaced with the average of the two nearest adjacent nonmissing values. A paired t-test was used to statistically compare conditions (LIFU, Sham) at a significance level of p < 0.05.
Autonomic metrics
Heart rate variability
Common heart rate variability (HRV) metrics were calculated according to prior literature (Shaffer and Ginsberg, 2017). Time domain metrics included the mean difference between normal sinus beats (MNN), the standard deviation of the inter-beat interval for normal sinus beats (SDNN), root mean squared of successive differences (RMSSD), percentage of successive intervals that differ by more than 50 ms (pNN50), and the coefficient of variation for successive differences between normal sinus beats (σ/μ) (CV). Time domain metrics were normalized to their respective mean of the 5 min baseline data and compared across Sham and LIFU conditions using paired t-tests at a Bonferroni-corrected significance level of p < 0.005. Data are reported as the mean ± SEM normalized units where 1 represents the mean of the 5 min baseline time window.
For the frequency domain, R-R interval data were cubic spline interpolated and transformed using the fast Fourier transform. Data during contact heat testing were normalized to the 5 min pre-task baseline data, and the normalized powers were averaged across the LF power (0.04–0.15 Hz) and HF power (0.15–0.4 Hz) bands for each subject for each condition (Shaffer and Ginsberg, 2017). The normalized LF/HF ratio was calculated by dividing the LF by the HF power. Normalized frequency domain powers were compared across Sham and LIFU conditions using a paired t-test at a Bonferroni-corrected significance level of p < 0.005. To investigate changes across specific HRV frequencies with each of the LF and HF bands, we performed nonparametric permutation testing of each frequency across conditions (Sham, LIFU) from 0.04 to 0.4 Hz. Permutation statistics were generated using 1,000 iterations at each frequency interval (0.001 Hz). Results were then cluster thresholded to 10 intervals, requiring a consecutive 0.01 Hz to be significant at p < 0.05.
Electrodermal response
Mean peak-to-peak (absolute of min–max) of the EDR data in the interval 0–10 s after the CHEP stimulus were compared using a Wilcoxon signed rank test with a Bonferroni-corrected significance level of p < 0.005.
Blood pressure
Post-session minus pre-session differences in systolic (SBP) and diastolic (DBP) blood pressures and mean arterial pressure (MAP) were calculated and compared between conditions. MAP was calculated according to the common equation used for BP acquired from a standard blood pressure cuff:
EEG
CHEP peak-to-peak analysis
The pre-processed data were first averaged across the 40 trials for each condition and each participant from the time window −2,000 ms to +2,000 ms around the CHEP stimulus. Classically defined CHEP components (De Schoenmacker et al., 2021) were assessed by manually identifying the N2 and P2 peaks in the trial-averaged data for each condition and participant. A distinct inflection of the waveform was necessary for inclusion in statistical analyses. The maximum number of missing trials within a single participant was two. To account for this, two random trials were removed from each participant, leaving 38 stimuli for each participant for the peak-to-peak analysis. The N2-P2 peak-to-peak amplitude as well as individual N2 and P2 amplitudes were then compared for differences between LIFU and Sham conditions. Wilcoxon signed rank test was used for all CHEP amplitude testing with a bonferroni-corrected significance set at p < 0.016 for each of the three metrics of the CHEP amplitude (N2-P2, N2 only, and P2 only).
Correlation between CHEP amplitudes and pain ratings
We next aimed to investigate what the relationship was between CHEP amplitudes and pain ratings and if LIFU changed this relationship. To test this, automated detection of minimum and maximum peak was performed for each trial for each participant and condition using the window 200–700 ms following the onset of the CHEP stimulus. Two participants were missing the 40th trial for the LIFU condition and one participant was missing the 40th trial for the Sham condition. To account for these three missing data points, the 39th stimulation was duplicated, creating 40 stimulations for all participants. The CHEP peak-to-peak amplitudes were then averaged across participants, creating 40 group-averaged amplitudes for each condition. The N2-P2 peak-to-peak CHEP amplitude, the N2 amplitude, and the P2 amplitude were each correlated with pain ratings across trials using Pearson’s correlations for both Sham and LIFU conditions. Rho values were then transformed into z-scores using Fisher’s r-to-z transformation (Meng et al., 1992). Z-scores of the correlation coefficients were then compared between the LIFU and Sham conditions using Fisher’s z-test at a Bonferroni-corrected significance level of p < 0.016.
Quantitative modeling of ultrasound wave propagation
Simulations were performed using the k-wave MATLAB toolbox (Treeby and Cox, 2010), which uses a pseudo-spectral time domain method to solve discretized wave equations on a spatial grid. Acoustic simulations were performed for each participant. Each participants’ dataset consisted of their structural MRI and CT of the head. CT images were used to construct the acoustic model of the skull, while MR images were used to target LIFU to the dACC coordinates. CT and MR images were first co-registered using custom scripts written in MATLAB, and the images were then up-sampled to a finer resolution for use in the acoustic simulations. The skull was extracted manually using a threshold intensity value and the intracranial space was assumed to be homogenous as ultrasound reflections between soft tissues are small (Mueller et al., 2016). Acoustic parameters were calculated from CT data assuming a linear relationship between skull porosity and the acoustic parameters (Aubry et al., 2003; Marquet et al., 2009).
The computational model of the ultrasound transducer used in simulations was constructed to recreate empirical acoustic pressure maps of focused ultrasound transmitted in an acoustic test tank, similar to a previous work (Mueller et al., 2016; Legon et al., 2018a). The modeled ultrasound transducer and waveform was designed to match the transducer and waveform used during pain testing. The shape of the waveform produced by the transducer can be visualized in Figure 2. The skull properties used in the model are previously reported by Legon et al. (2018a). Estimated intracranial pressure values were extracted from each participants’ dACC brain volume, expressed in kilopascals (kPa) ± (SEM). Intensity parameters were derived from measured values of pressure using the approximation of plane progressive acoustic radiation waves (Preston, 1986).
Transducer characteristics and acoustic modeling. A, Pseudo-color free water XY lateral beam pressure at Z maximum of the 500 kHz single-element transducer measured in kilopascals (kPa). B, Pseudo-color free water YZ axial beam pressure. The transducer face is at Z = 0 mm. C, Full-width half maximum (FWHM) overlay of the empirically measured YZ axial beam in acoustic test tank (blue) and the modeled beam (red) used for acoustic modeling. D, Acoustic modeling results using individual MRI and CT scans from a single representative subject. From left to right: Coronal view depicting left (L) and right (R) directions. Sagittal view depicting posterior (P) and anterior (A) directions. Transverse cross-section taken from the white line in the coronal and sagittal views.
LIFU targeting data analysis
The distance from the scalp marker set in the neuronavigation system to the predetermined dACC target was calculated and expressed in millimeters as an indicator of transducer placement (targeting) error. Targeting error was averaged across the 40 trials for each participant in both Sham and LIFU conditions and expressed as mean ± standard deviation. Group-averaged data were compared using a paired t-test to evaluate differences in targeting between conditions at a significance level of p < 0.05.
Analysis of auditory confounds
Data from the auditory query for each question (all scored from 1–7) were first collapsed into a 1–3 score in order to reduce the dimensions for comparison between conditions. An original response of a 1–3 (“Strongly Disagree,” “Disagree,” “Slightly Disagree”) was given a score of one indicating a response of “Disagree.” A score of four was assigned a value of two indicating a response of “Unsure.” An original score of 5–7 (“Slightly Agree,” “Agree,” and “Strongly Agree”) was assigned a score of three indicating a response of “Agree.” Data were averaged across participants for each condition and compared using a signed rank test at a significance level of p < 0.05.
Results
Empirical acoustic measurements
Measurements of the transducer’s beam profile in free water revealed a full-width half maximum (FWHM) lateral (XY) resolution of ±1.5 mm at the axial (Z) beam focus (Fig. 2A) and a FWHM axial (YZ) resolution of ±10 mm with a focal depth of 52 mm from the exit plane of the transducer (Fig. 2B). The axial FWHM comparison between the empirical measurements in free water and the modeled waveform demonstrated good agreement validating use in the acoustic models (Fig. 2C).
Acoustic modeling
The extracranial peak negative pressure was 368.5 kPa for all participants. This translated to an intensity outside of the skull of 4.34 W/cm2 spatial-peak, pulsed averaged (Isppa) or 1.56 W/cm2 spatial-peak, temporal averaged (Ispta). Acoustic modeling estimated a mean ± SD intracranial pressure at the dACC target of 115.2 ± 53.1 kPa with a range of 66.3–287.4 kPa (Table 1). This resulted in an average attenuation through the skull of 68.4 ± 13.8%. Modeling results of the intracranial beam profile in a single representative subject in the coronal, sagittal, and transverse views show targeting of the dACC in Figure 2D. Table 1 shows individual subject estimated intracranial pressure.
LIFU targeting
The mean ± SD depth from the scalp to the dACC target across all participants was 47.34 ± 3.62 mm (range, 41.5–54.5 mm). For females the mean ± SD depth was 45.73 ± 2.76 mm (range, 41.5–51.6 mm) and for males it was 50.02 ± 3.44 mm (range, 45.6–54.5 mm). Given the FWHM of the axial beam profile at ± 10 mm and focus at 52 mm, the beam overlapped well with the dACC target across participants. Table 1 shows individual subject depths of dACC targets from the scalp.
The mean ± SD targeting error across participants was 1.50 ± 0.81 mm in the Sham condition and 1.18 ± 0.40 mm in the LIFU condition (Fig. 3A). A paired t-test showed no significant difference between conditions (t(15) = −1.39, p = 0.18). A histogram of the distance from the target for each LIFU stimulation across all participants and conditions can be seen in Figure 3B.
Targeting error. A, (Top) Mean ± standard deviation (SD) targeting error distance of the transducer on the scalp in millimeters (mm) for each subject in the Sham condition. Subject 8 in the Sham condition shows high variability due to loss of data. (Bottom) Mean ± SD targeting error distance in millimeters (mm) for each subject in the low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) condition. B, Group (N = 16) histogram of the absolute targeting error (mm) for both Sham (blue) and LIFU to dACC (red) conditions. The x-axis is targeting error in 0.5 mm bins.
Pain ratings
The mean ± SEM pain ratings were 3.61 ± 0.26 for the Sham condition and 2.52 ± 0.37 for the LIFU condition, with a difference of 1.09 ± 0.20 or 30.3% (range, −2.30 to +0.23). A paired t-test revealed a significant difference in pain ratings between conditions (t(15) = −5.39, p = 7.46 × 10−5) (Fig. 4A). Group (N = 16) pain ratings across the 40 trials for both Sham and LIFU conditions are shown in Figure 4B.
Effect of LIFU to dACC on pain ratings. A, Group (N = 16) mean ± SEM pain ratings (0–9 scale). Red lines are individual subject data. Asterisk indicates significant difference between conditions at p < 0.0001. B, Group (N = 16) mean ± SEM pain ratings across 40 trials for Sham (blue) and low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) (red). Thin lines and shaded regions are the mean ± SEM for individual trials. Thick lines represent five-trial averages.
Autonomic responses
Time domain
SDNN
The mean ± SEM normalized SDNN for the Sham condition was 0.77 ± 0.044 and for the LIFU condition was 0.95 ± 0.063 with a difference of +0.17 ± 0.046 under the LIFU condition relative to Sham. A paired t-test showed a significant difference between the conditions (t(15) = 3.76, p = 0.002) (Fig. 5A). This would suggest that during painful stimuli, SDNN is reduced as compared to baseline levels (below 1) and that LIFU looks to return these levels toward baseline resting levels.
Effect of LIFU to the dACC on time domain autonomic responses. The y-axes for each normalized metric (A–E) are the respective metrics during contact heat-evoked potential (CHEP) testing normalized to values during the 5 min baseline period. Horizontal dashed black lines represent the pre-task baseline. A, Group (N = 16) mean ± SEM normalized standard deviation of normal sinus beats (SDNN). Red lines represent individual subject data. Asterisk indicates significant differences between conditions at p < 0.01. B, Group (N = 16) mean ± SEM normalized CV. Red lines represent individual subject data. C, Group (N = 16) mean ± SEM normalized mean normal sinus difference (MNN). D, Group (N = 16) mean ± SEM normalized root mean square of successive differences (RMSSD). E, Group (N = 16) mean ± SEM normalized proportion of normal sinus differences greater than 50 ms (pNN50). F, Group (N = 16) mean ± SEM electrodermal response (EDR). Y-axis is amplitude in microsiemens (μS).
Coefficient of variation (CV)
The mean ± SEM normalized CV was 0.80 ± 0.054 for the Sham condition and 0.92 ± 0.06 for the LIFU condition with a difference of +0.12 ± 0.04 under the LIFU condition relative to Sham. A paired t-test demonstrated no significant difference after correction for multiple comparisons (t(15) = 3.01, p = 0.009) (Fig. 5B).
MNN
The mean ± SEM normalized MNN was 1.03 ± 0.02 under the Sham condition and 1.03 ± 0.01 under the LIFU condition with a difference of 0.0013 ± 0.014. A paired t-test showed no significant difference between conditions (t(15) = −0.09, p = 0.93) (Fig. 5C).
RMSSD
The mean ± SEM normalized RMSSD was 1.03 ± 0.02 under the Sham condition and 1.03 ± 0.01 under the LIFU condition with a difference of 0.0002 ± 0.014. A paired t-test showed no significant difference between conditions (t(15) = −0.01, p = 0.99) (Fig. 5D).
pNN50
The mean ± SEM normalized pNN50 was 1.38 ± 0.16 under the Sham condition and 1.41 ± 0.12 under the LIFU condition for a difference of 0.031 ± 0.22. A paired t-test showed no significant difference between conditions (t(15) = 0.14, p = 0.89) (Fig. 5E).
Electrodermal response
EDR mean ± SEM amplitude was 0.014 ± 0.0036 µS in the Sham condition and 0.009 ± 0.0021 µS in the LIFU condition. A Wilcoxon sign rank revealed a nonsignificant difference of −0.0037 ± 0.0031 µS between conditions (z = −1.31, p = 0.19) (Fig. 5F).
Frequency domain
Low-frequency (LF) power
The mean ± SEM normalized LF power for the Sham condition was 0.52 ± 0.06 and 0.93 ± 0.09 for the LIFU condition with a difference of +0.41 ± 0.09 under the LIFU condition relative to Sham. A paired t-test demonstrated a significant difference between conditions (t(15) = 4.54, p = 3.89 × 10−4) (Fig. 6A).
Effect of LIFU to the dACC on frequency-domain autonomic responses. The y-axes for each normalized metric are the respective powers during contact heat-evoked potential (CHEP) testing normalized to powers during the 5 min baseline period. Horizontal dashed black lines represent the pre-task baseline. A, Group (N = 16) mean ± SEM normalized low-frequency (LF) power. Red lines represent individual subject data. Asterisk indicates significant differences between conditions at p < 0.001. B, Group (N = 16) mean ± SEM normalized high-frequency (HF) power. C, Group (N = 16) mean ± SEM normalized low-frequency/high-frequency (LF/HF) ratio. Red lines represent individual subject data. Asterisk indicates significant differences between conditions at p < 0.01. D, Group (N = 16) mean ± SEM normalized heart rate variability (HRV) power. X-axis is frequency in hertz (Hz) and y-axis is normalized power. Lines and shaded regions represent the mean ± SEM for Sham (blue) and low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) (red) conditions. Black vertical bars represent significant differences between conditions at p < 0.05. Dashed vertical black line at 0.15 on the x-axis delineates LF (0.04–0.15) and HF (0.15–0.4) power.
High-frequency (HF) power
The mean ± SEM normalized HF power for the Sham condition was 0.95 ± 0.11 and for the LIFU condition was 1.13 ± 0.01 with a difference of +0.18 ± 0.10 under the LIFU condition relative to Sham. A paired t-test showed no significant difference between conditions (t(15) = 1.79, p = 0.09) (Fig. 6B).
LF/HF ratio
The mean ± SEM normalized LF/HF ratio in the Sham condition was 0.57 ± 0.06 and 0.83 ± 0.06 in the LIFU condition with a difference of +0.26 ± 0.07 under the LIFU condition relative to Sham. A paired t-test demonstrated a significant difference between conditions (t(15) = 3.61, p = 0.003) (Fig. 6C).
To determine if effects were driven by specific frequencies, we performed permutation statistics (1,000 permutations; p < 0.05) across frequencies 0.04–0.4 Hz (0.001 resolution). Permutation statistics revealed LIFU to the dACC increased normalized LF power in the frequency range 0.04–0.13 Hz. There were no significant differences in normalized HF power in the permutation testing (Fig. 6D).
Blood pressure
Under the Sham condition, the mean ± SEM pre-session SBP was 119.27 ± 2.84 mmHg and the DBP was 70.33 ± 1.86 mmHg. The mean ± SEM post-session SBP was 119.00 ± 2.55 mmHg and the DBP was 69.6 ± 2.34 mmHg. Under the LIFU condition, the mean ± SEM pre-session SBP was 121.2 ± 3.21 mmHg and the DBP was 70.87 ± 2.15 mmHg. The mean ± SEM post-session SBP was 117.67 ± 3.75 mmHg and the DBP was 71.27 ± 2.68 mmHg. Mean ± SEM post-session minus pre-session SBP differences were −0.27 ± 2.57 mmHg for Sham and −3.53 ± 1.89 mmHg for LIFU. A paired t-test revealed no significant difference between conditions (t(14) = −1.29, p = 0.22). Mean ± SEM post-session minus pre-session DBP differences were −0.73 ± 1.42 mmHg for Sham and +0.40 ± 1.27 mmHg for LIFU. A paired t-test showed no significant difference between conditions (t(14) = 1.03, p = 0.32).
Mean arterial pressure (MAP)
Under the Sham condition, the mean ± SEM pre-session MAP was 86.64 ± 2.07 while post-session MAP was 87.64 ± 2.31. Under the LIFU condition, the mean ± SEM pre-session MAP was 86.07 ± 2.28 while the post-session MAP was 86.73 ± 2.84. The mean ± SEM post-session minus pre-session MAP differences were −0.58 ± 1.57 mmHg for Sham and −0.91 ± 1.25 mmHg for LIFU. A paired t-test revealed no significant difference between conditions (t(14) = −0.23, p = 0.82).
CHEP amplitude
We conducted peak-to-peak analysis to test whether LIFU affected the N2-P2 peak-to-peak amplitude of the CHEP as well as the N2 and P2 amplitudes separately. The mean ± SEM N2-P2 peak-to-peak amplitude was 26.73 ± 0.73 µV for the Sham condition and 20.06 ± 0.48 µV for the LIFU condition for a difference of −6.67 ± 2.81 µV or 25.1% under the LIFU condition relative to Sham. A Wilcoxon signed rank test revealed no significant difference in the mean N2-P2 amplitude between conditions after correction for multiple comparisons (z = −2.17, p = 0.030) (Fig. 7A,D).
Effect of LIFU to the dACC on CHEP amplitudes. A, Group (N = 16) mean ± SEM N2-P2 peak-to-peak contact heat-evoked potential (CHEP) amplitude. Red lines are individual subject data. B, Group (N = 16) mean ± SEM P2 CHEP amplitude. Red lines are individual subject data. Asterisk denotes a significant difference between groups at p < 0.05. C, Group (N = 16) mean ± SEM N2 CHEP amplitude. Red lines are individual subject data. D, Group (N = 16) mean ± SEM CHEP trace for Sham (blue) and low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) (red) conditions. Amplitude is in microvolts (μV). Time 0 ms represents delivery of contact heat stimulus to the dorsum of the right hand. Arrows point to N2 and P2 peaks. E, Group (N = 16) mean ± SEM P2 CHEP amplitude across the 40 trials for Sham (blue) and LIFU to dACC (red). Thin lines and shaded regions are the mean ± SEM amplitudes for each of the 40 trials. Thick lines and dots represent five-trial averages. F, (Left) Group (N = 16) mean ± SEM P2 amplitudes in orange and pain ratings in green across the 40 trials for the Sham condition. (Right) Group (N = 16) scatter plot depicting the correlation between the P2 amplitude and pain ratings across trials for the Sham condition. r = correlation coefficient, p = p-value, black line = least-squares line. G, (Left) Group (N = 16) mean ± SEM P2 amplitudes in orange and pain ratings in green across the 40 trials for the LIFU to dACC condition. (Right) Group (N = 16) scatter plot depicting the correlation between the P2 amplitude and pain ratings across trials for the LIFU to dACC condition. r = correlation coefficient, p = p-value, black line = least-squares line.
Investigation of the N2 and P2 components of the CHEP separately demonstrated the mean ± SEM amplitude of the P2 component was 17.20 ± 2.28 µV for the Sham condition and 10.65 ± 1.05 µV for the LIFU condition with a difference of −6.54 ± 2.39 µV or −38.1% under the LIFU condition relative to Sham. The Wilcoxon signed rank test showed a significant difference in the P2 component between conditions (z = −2.48, p = 0.013) (Fig. 7B,D). P2 amplitudes across trials for both Sham and LIFU conditions are shown in Figure 7E. The mean ± SEM amplitude of the N2 component was 9.53 ± 1.74 µV for the Sham condition and 9.41 ± 1.23 µV for the LIFU condition. The Wilcoxon signed rank test showed no significant difference in the N2 component between conditions (z = −0.05, p = 0.96) (Fig. 7C,D).
CHEP amplitude and pain rating correlations
To examine the relationship between CHEP amplitudes and behavior, we performed correlations between pain ratings and the N2-P2 peak-to-peak amplitude, N2 amplitude, and P2 amplitude across the 40 trials for each condition. We used Fisher’s Z test to analyze the change in correlations between conditions (Meng et al., 1992). The P2 component did not correlate with pain ratings for the Sham condition (r = 0.06, p = 0.74) but had a significant positive correlation for the LIFU condition (r = 0.64, p = 9.64 × 10−6) (Fig. 7F,G). A comparison of the Fisher z-transformed correlation coefficients revealed a significantly increased correlation for the LIFU condition (z = 3.01, p = 0.003). The N2-P2 peak-to-peak amplitude did not significantly correlate with pain ratings for the Sham condition (r = 0.31, p = 0.054) but showed a significant positive correlation for the LIFU condition (r = 0.64, p = 8.17 × 10−6) (Fig. 8A,B). Comparison of the fisher z-transformed correlation coefficients revealed a nonsignificant increase in the correlation under the LIFU condition (z = 1.91, p = 0.057). The N2 amplitude had a significant positive correlation with pain ratings for the Sham condition (r = 0.42, p = 0.007) and for the LIFU condition (r = 0.43, p = 0.005) (Fig. 8C,D). Comparison of the fisher z-transformed correlation coefficients revealed no significant difference between conditions (z = −0.06, p = 0.95).
Correlation of CHEP amplitudes and pain ratings. A, (left) Group (N = 16) mean ± SEM N2-P2 peak-to-peak contact heat-evoked potential (CHEP) amplitudes (orange) and pain ratings (green) across the 40 trials for the Sham condition. (Right) Group (N = 16) scatter plot depicting the correlation between the N2-P2 peak-to-peak CHEP amplitude and pain ratings across trials for the Sham condition. r = correlation coefficient, p = p-value, black line = least-squares line. B, (left) Group (N = 16) mean ± SEM N2-P2 peak-to-peak amplitudes (orange) and pain ratings (green) across the 40 trials for the low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) condition. (right) Group (N = 16) scatter plot depicting the correlation between the N2-P2 peak-to-peak CHEP amplitude and pain ratings across trials for the LIFU to dACC condition. r = correlation coefficient, p = p-value, black line = least-squares line. C, (left) Group (N = 16) mean ± SEM N2 amplitude (orange) and pain ratings (green) across the 40 trials for the Sham condition. (right) Group (N = 16) scatter plot depicting the correlation between the N2 amplitude and pain ratings across trials for the Sham condition. r = correlation coefficient, p = p-value, black line = least-squares line. D, (left) Group (N = 16) mean ± SEM N2 amplitude (orange) and pain ratings (green) across the 40 trials for the LIFU to dACC condition. (right) Group (N = 16) scatter plot depicting the correlation between the N2 amplitude and pain ratings across trials for the LIFU to dACC condition. r = correlation coefficient, p = p-value, black line = least-squares line.
Auditory masking
Hearing
For the question “I could hear the LIFU stimulation,” the mean ± SD response for the Sham condition was 1.56 ± 0.89 and 1.69 ± 0.95 for the LIFU condition indicating that participants were between “Disagree” and “Unsure” in their responses. A Wilcoxon sign rank test demonstrated no significant difference between conditions (z = 0.44, p = 0.66).
Believing
For the question “I believe I experienced LIFU stimulation,” the mean ± SD response for the Sham condition was 1.31 ± 0.60 and 1.63 ± 0.96 for the LIFU condition. A Wilcoxon sign rank test demonstrated no significant difference between conditions (z = 1.39, p = 0.16).
Feeling
For the question “I could feel the LIFU stimulation,” the mean ± SD response for the Sham condition was 1.11 ± 0.46 and 1.05 ± 0.22 for the LIFU condition. A Wilcoxon sign rank test demonstrated no significant difference between conditions (z = −1.00, p = 0.32).
Individual subject results for all questions can be seen in Figure 9 and a single-subject analysis of the three participants who appeared to correctly determine the condition can be found in Table 2. For these three participants, the mean ± SD difference in pain ratings, SDNN, and the P2 CHEP amplitudes were −0.80 ± 1.00 points, + 0.14 ± 0.10 n.u., and −6.00 ± 7.47 µV, respectively. Removing these participants from the data does not affect our results.
Individual subject auditory query results. A, Individual subject results between Sham (blue) and low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) (red) conditions for the question “I could hear the LIFU stimulation.” B, Individual subject results between Sham (black) and LIFU to dACC (gray) conditions for the question “I believe I experienced LIFU stimulation.” C, Individual subject results between Sham (black) and LIFU to dACC (gray) conditions for the question “I could feel the LIFU stimulation.” The participants on the x-axis are in the same order for all three questions.
Mean differences for subjects who reported correctly hearing the LIFU
Adverse events
For the Sham condition, one participant reported mild unusual feelings, one reported mild tingling, one reported mild itching, one reported mild sleepiness, two reported mild inattentiveness, one reported mild twitching, one reported a mild anxious feeling, and three reported mild balance issues (Fig. 10A). For the LIFU condition, one participant reported a mild headache, one reported mild unusual scalp sensations, two reported mild neck pain, one reported mild tooth pain, and one reported mild nausea (Fig. 10B). No moderate or serious adverse events were reported for either condition.
Report of symptoms questionnaire results. Report of symptoms was administered before and after contact heat testing for both conditions. The presence of symptoms for each participant after contact heat testing were adjusted for the presence of symptoms prior to contact heat testing by taking the pre-post difference. Each bar represents the total count for the presence of each symptom and severity across all participants (N = 16). No change or improvement in symptoms is represented as absent. A, Presence of symptoms in the Sham condition. B, Presence of symptoms in the low-intensity focused ultrasound (LIFU) to the dorsal anterior cingulate (dACC) condition.
Discussion
We evaluated the effect of single element 500 kHz LIFU to the dACC on pain, time and frequency domain metrics of heart rate variability, blood pressure, and electrodermal response, as well as CHEPs during application of brief, noxious contact heat. LIFU to the dACC reduced pain by 1.09 points on a 10-point scale and increased HRV metrics including the SDNN, LF power, and LF/HF ratio as compared to Sham stimulation. LIFU also attenuated the CHEP P2 amplitude.
Noninvasive methods to reduce acute pain perception
We observed decreases in pain levels by 1.09 points or 30.3% by selectively targeting the dACC (MNI coordinate 0,18,30), in contrast to previous reports using other noninvasive techniques that have largely shown equivocal effects. Tzabazis tested three coil configurations and showed one resulted in significant reductions in pain intensity and dACC activity using positron emission tomography (PET) (Tzabazis et al., 2013). Galhardoni targeted the dACC using a double-cone H-6 TMS coil with no reduction in pain ratings (Galhardoni et al., 2019). Modeling of the electric field in the head from various TMS coil configurations demonstrates significant depth-focality trade-offs such that even with deep H-type TMS coils the effective electric field intensity does not extend beyond 3–3.5 cm from the scalp (Deng et al., 2013; Zibman et al., 2021; Drakaki et al., 2022). The average depth of MNI coordinate used in this study (0,18,30) across male and female subjects was 4.7 cm with the most superficial being 4.15 cm.
Studies using TES to target the ACC for pain show some effects, but with limited spatial selectivity of the applied electric field in the brain and frequently lacking sham conditions (Auvichayapat et al., 2018; Magis et al., 2018; Xiong et al., 2023). While some work shows TES can reach the dACC (Mariano et al., 2019), the authors acknowledge they did not account for individual brain morphology. The applied electric field from TES is highly influenced by individual brain anatomy including CSF thickness and depth of sulci, impacting its distribution (Opitz et al., 2015). By contrast, LIFU focuses mechanical energy at adjustable depths that is unaffected by brain morphology (Mueller et al., 2016). This improvement in spatial localization and depth penetration likely contributes to the significantly reduced pain found here.
LIFU for noninvasive autonomic modulation
Brief, painful stimuli elevate heart rate and alter HRV consistent with increased sympathetic or decreased parasympathetic activity indexed by decreased SDNN and HF power (Loggia et al., 2011; Meeuse et al., 2013; Forte et al., 2022). Behavioral methods to reduce pain levels using virtual reality or deep breathing also increases HRV, especially SDNN (Chalaye et al., 2009; Colloca et al., 2020). Similarly, LIFU to the dACC increased SDNN, LF power, and LF/HF ratio toward pre-task baseline levels during painful stimuli.
The dACC is innervated by spinothalamic input from the ventrocaudal mediodorsal (MDvc) nucleus of the thalamus, which receives input from parabrachial nucleus and periaqueductal gray and projects to the fundus of the cingulate sulcus, providing an integrated homeostatic pathway to the dACC (Yezierski, 1988; Apkarian and Hodge, 1989; Ray and Price, 1993; Craig, 1995, 2004; Dum et al., 2009). The dACC is considered a homeostatic effector region due to its direct role in autonomic and behavioral control, particularly the sympathetic autonomic system (Shackman et al., 2011; Beissner et al., 2013). Direct electrical stimulation of the dACC elicits sympathetic autonomic responses including tachycardia as well as heat sensations, consistent with excitation (Parvizi et al., 2013; Oane et al., 2020; Xue et al., 2023). LIFU likely conferred inhibitory effects as, in addition to reduced autonomic responses, we observed reductions in the CHEP P2 amplitude and the parameters we employed have previously demonstrated inhibitory effects in other brain regions (Legon et al., 2014, 2018b,a).
HRV measures homeostatic regulation, indexing the ability of the autonomic nervous system to dynamically respond to task demands (Rajendra Acharya et al., 2006; ChuDuc et al., 2013). The SDNN is a validated time domain metric of overall HRV and in short-term recordings, like the present study, SDNN is primarily influenced by parasympathetic activity with a higher SDNN being consistent with reduced sympathetic or increased parasympathetic activity (Shaffer and Ginsberg, 2017). Reductions of SDNN with delivery of pain and increases toward baseline under the LIFU condition is consistent with direct inhibition of sympathetic activity or indirect facilitation of parasympathetic control and an indication we can shift homeostatic responses to painful stimuli. LF power (0.04–0.15 Hz) is influenced by sympathetic and parasympathetic systems and correlates with SDNN (Shaffer et al., 2014; Shaffer and Ginsberg, 2017). In resting conditions, like our 5 min baseline period, LF power reflects baroreflex/parasympathetic activity, not cardiac sympathetic innervation (Goldstein et al., 2011; Reyes Del Paso et al., 2013). The increased LF power under the LIFU condition thus fits our findings of increased SDNN, providing confirmation using frequency domain HRV metrics. Furthermore, our dACC coordinate for LIFU closely matches prior literature linking dACC activity and LF power (Critchley et al., 2003; Rebollo et al., 2018). Interpretations of the LF/HF ratio primarily indexing sympathovagal balance is controversial (Billman, 2013), and it is likely appropriate to interpret LF/HF ratio based on the influence of measurement conditions on LF power (Shaffer et al., 2014).
Our results point toward LIFU to dACC directly or indirectly attenuating sympathetic responses to heat pain as our autonomic metrics decreased during contact heat testing under Sham and increased toward pre-task baseline levels under LIFU. Reductions in sympathetic (or facilitation of parasympathetic) activity may have applications for chronic pain and anxiety populations, which exhibit parasympathetic dysregulation, sympathetic hyper-reactions, and reduced HRV (Chalmers et al., 2014; Tracy et al., 2016; Zhuo, 2016; Teed et al., 2022).
Contribution of dACC to nociceptive potentials
The CHEP is a well-studied biomarker of pain processing in the dACC similar to the laser-evoked potential (LEP) (Chen et al., 2001; Garcia-Larrea et al., 2003; Shenoy et al., 2011; Granovsky et al., 2016; De Schoenmacker et al., 2021). Peak-to-peak analysis showed LIFU to the dACC selectively reduced the P2 amplitude of the CHEP. Correlations between CHEP amplitudes and pain ratings showed a significant change in the relationship between the P2 amplitude and pain ratings, but not N2-P2 or N2 amplitudes, further suggesting we selectively modulated the P2 amplitude.
Intracranial EEG responses to laser stimuli show the dACC predominately encodes the positive amplitude of the LEP (Bastuji et al., 2016). Dipole sourcing corroborates the dACC in generating the LEP (Jin et al., 2018) while regions of peak cerebral blood flow during LEP stimuli matches where we delivered LIFU (Garcia-Larrea et al., 2003).
While the dACC at least partially contributes to pain evoked potentials, it is unclear precisely what N2 and P2 components of the LEP/CHEP index. Evidence points to attentional modulation of the positive (P2) component of the CHEP with a recent study in autism spectrum disorder showing the P2 strongly correlates with attention to pain, supporting prior LEP work (Ohara et al., 2004; Boyle et al., 2008; Chien et al., 2017). Research using fMRI shows reduced dACC activity during distraction from pain (Bantick et al., 2002). Distraction induced reductions in positive heat-evoked amplitudes appears to be concomitant with reduced unpleasantness ratings, or affective pain (Boyle et al., 2008). That this positive amplitude has been associated with attentional and affective components of pain matches proposed roles of the dACC in pain processing (Rainville et al., 1997; Price, 2000; Shackman et al., 2011; Xiao and Zhang, 2018) and the role of the dACC in the salience network (Seeley, 2019). The dACC is thus a likely contributor to positive heat-evoked amplitudes while playing a central role in affective and attentional aspects of pain.
Auditory masking
LIFU can produce auditory effects that may confound neuromodulatory effects (Guo et al., 2018; Sato et al., 2018; Johnstone et al., 2021). To mitigate this, we employed a multitoned mask to blind subjects to the condition, which was recently demonstrated to effectively mask the sound of LIFU at the same 1 kHz pulse repetition frequency used in the present study (Liang et al., 2023). The masking query data support effective masking as most subjects were unsure they received LIFU and those that did correctly identify the LIFU condition did not drive the observed differences between LIFU and Sham conditions.
Limitations
While participants and researchers responsible for data analysis were blinded to the condition, it was not possible at the time of study to blind the researchers delivering the LIFU. As such, it is possible the researcher(s) could have inadvertently conveyed to the subject whether they received Verum or Sham stimulation. Another potential limitation is the lack of an active Sham condition, which would directly address nonspecific auditory confounds.
Conclusions and Future Work
LIFU provides a noninvasive and selective method to causally support of the role of the dACC in pain and autonomic processing in humans. Future work should aim to evaluate the longevity and durability of LIFU’s effects, LIFU’s ability to modulate autonomic function independent of pain, as well as its efficacy as a therapeutic option in chronic pain populations.
Data Availability Statement
Data and source code is available upon reasonable request.
Footnotes
This work was funded in part by grants awarded to WL from the Seale Innovation Fund. The authors would like to thank Jessica Florig for help with data collection and analysis.
The authors declare no competing financial interests.
- Correspondence should be addressed to Andrew Strohman at andrewstrohman{at}vt.edu or Wynn Legon at wlegon{at}vt.edu.
