Focused ultrasound alters neural time irreversibility independently of response amplitude

Reviewer #1 (Public review):

Summary:

The authors introduce ordinal EPR (Dsym) as a novel metric for characterizing tFUS-evoked calcium responses. The concept is interesting and could offer an innovative way of looking at neural responses to neuromodulation more broadly. Their main result, that EPR carries information beyond mean GCaMP amplitude, is compelling, but the manuscript would be strengthened by testing whether it holds against other conventional GCaMP metrics beyond mean amplitude (e.g. decay time). More importantly, the EPR metric requires further validation before its central claim, that it reflects genuine dynamical reorganization of the underlying circuit rather than artifacts of the measurement pipeline (e.g. GCaMP indicator kinetics), can be accepted. A rigorous surrogate/synthetic signal validation would be a very valuable, if not essential, addition to this paper.

Strengths:

This manuscript frames tFUS-evoked activity in a way that is uncommon in the field. Moving beyond amplitude-based readouts to ask how sonication reshapes the temporal organization of neural activity is a novel contribution to the ultrasound neuromodulation literature. The finding that EPR shows a distinct dose-response profile from calcium amplitude and retains dose-related structure after controlling for amplitude, is a promising demonstration that this framework can extract information not visible to conventional measures.

Methods are very detailed and well explained for reproducibility purposes.

Weaknesses:

Major weaknesses:

EPR may carry information beyond mean GCaMP amplitude- but is it carrying information beyond other elements of the GCaMP signal? For the on-target recovery window, the calcium signal has not yet returned to baseline (original data traces, Fig. 2C, Fig. 4A), so the trace still contains a residual decay transient at the time Dsym is computed. This decay rate reflects neuronal activity returning to baseline and will be dependent on things such as calcium indicator kinetics. Would a change in Dsym always accompany this kind of decay, regardless of whether this has anything to do with endogenous dynamics, irreversibility etc.? Would any gradual return to baseline/decay of this kind produce elevated Dsym on its own, independent of underlying stochastic fluctuations/endogenous dynamics? Is it not expected that a system that is not at steady state shows entropy production? Perhaps the decay rate is simply different by dose, which explains the findings? Could this be tested i.e. with synthetic data or by removing the decay? The fact that a change in EPR is measured when there is no significant GCaMP change indicates that the measure is not purely dominated by decay kinetics but it is not clear whether decay kinetics etc. could be confounding this result.

GCaMP's rise and decay kinetics are asymmetric (fast rise, slow decay). Since Dsym is a measure of asymmetry between forward and time-reversed ordinal statistics, could this kinetic asymmetry alone have any impact on "irreversibility", independent of underlying neural dynamics? Could this be tested with synthetic data? The authors note that the GCaMP signal is filtered through an indicator with slow kinetics but don't mention the asymmetry.

When different stimulation periods are used, what impact if any, does this have on computation of the EPR metric?

What impact does having a noisier/lower SNR calcium signal have on the EPR metric, if any?

Minor weaknesses:

The major advantage of TUS as a non-invasive brain stimulation technique is its capacity for spatially focused deep brain stimulation. However, it carries a significant auditory confound, which leads to activation of widespread (not spatially restricted) neuronal networks (Kop et al., 2024; Sato et al., 2018). The data used in this paper does not adequately control for this confound (i.e. deafened animals (Guo et al., 2023; Sato et al., 2018). Whilst this is not the responsibility of the authors of this manuscript, it should be mentioned in the discussion that any of the tFUS-calcium evoked responses could be due to indirect auditory stimulation rather than a direct pressure-mediated effect. The authors may just be measuring EPR in response to the auditory confound, which does not undermine the overall impact of this paper (as it is more focused on an approach to looking at this kind of data), but should be mentioned. I do think that presence of the auditory confound could impact interpretability of the results in the manuscript. Could the authors comment on this? Could auditory mediated arousal or state shift (perhaps activating brain regions not captured by the fiber) explain any of the presented EPR results or affect interpretation?

Guo, H., Salahshoor, H., Wu, D., Yoo, S., Sato, T., Tsao, D. Y., & Shapiro, M. G. (2023). Effects of focused ultrasound in a "clean" mouse model of ultrasonic neuromodulation. iScience, 26(12). https://doi.org/10.1016/j.isci.2023.108372
Kop, B. R., Shamli Oghli, Y., Grippe, T. C., Nandi, T., Lefkes, J., Meijer, S. W., Farboud, S., Engels, M., Hamani, M., Null, M., Radetz, A., Hassan, U., Darmani, G., Chetverikov, A., den Ouden, H. E., Bergmann, T. O., Chen, R., & Verhagen, L. (2024). Auditory confounds can drive online effects of transcranial ultrasonic stimulation in humans. eLife, 12, RP88762. https://doi.org/10.7554/eLife.88762
Sato, T., Shapiro, M. G., & Tsao, D. Y. (2018). Ultrasonic Neuromodulation Causes Widespread Cortical Activation via an Indirect Auditory Mechanism. Neuron, 98(5), 1031-1041.e5. https://doi.org/10.1016/j.neuron.2018.05.009

Reviewer #2 (Public review):

Summary:

The present paper studies the entropy production rate (EPR) before, during, and after thalamic sonication in freely moving mice from a previously published dataset. The motivation is that EPR contains information about the dynamics of neural activity that is not present in the average calcium amplitude. It was observed that the change in EPR reaches a maximum at an intermediate stimulation dosage, in contrast to the change in overall calcium amplitude, which increases monotonically with dose in the on-target condition. The paper further observes a statistically significant relation between the baseline EPR and the change in EPR during stimulation, but not for recovery, as well as a similar relation for the calcium response. The manuscript draws an analogy between the non-monotonic dependence and stochastic resonance.

Overall, these results present an interesting analysis of how acoustic stimulation affects the dynamics of neural activity. In my view, these tools would be quite useful for quantifying changes in neural activity under different perturbations.

Strength:

The study is well motivated by the argument that quantifying neural response requires information about dynamics that is not contained in the average calcium amplitude. EPR or irreversibility has been shown to be a powerful tool in describing out-of-equilibrium biological processes. The paper successfully quantifies irreversibility using an ordinal surrogate of the entropy production rate. The non-monotonic dependence is an interesting result, which demonstrates (with caveats stated in Weakness) that EPR contains more information about neural response than the average calcium amplitude. Furthermore, the paper also clearly states limitations of the approach and performs a robustness check by varying the parameters used in the ordinal EPR estimation.

Weakness:

(1) My main concern, as the authors have touched upon in the Discussion, is that although irreversibility clearly provides a readout of the neuronal dynamics, it remains unclear what its biological significance is. While it is related to susceptibility, the mechanistic link remains weak. I suspect that this interpretability issue will limit the impact of this approach.

(2) I wonder how the measured irreversibility is affected by the asymmetric response of GCaMP, which typically rises quickly and decays slowly. It seems possible that this temporal asymmetry, combined with the average calcium response, already leads to a non-monotonic irreversibility curve. The paper should investigate this and other potential contributors to the temporal irreversibility in the readout.

(3) As shown in Fig. 2AB, the baseline EPR already varies considerably with acoustic intensity, by an amount comparable to the change in EPR after stimulation (Fig. 2E). Since the baseline window precedes stimulation, this variation cannot be caused by the sonication. It suggests either that the uncertainty in EPR quantification is larger than the error bars indicate, or that trials at different nominal intensities differ systematically in some other respect. It should be examined whether the non-monotonic trend is statistically significant in light of this baseline fluctuation.

Author response:

We thank the reviewers and the editor for providing a thorough and constructive review of the manuscript. Both reviewers indicated that the ordinal entropy production rate should be validated with a synthetic (surrogate) model of the calcium signal. We agree with this, and it will feature prominently in the revised version of the manuscript.

As described (see Author response image 1), we have already implemented a synthetic GCaMP model. In the simulations, the strength of the calcium intensity varies monotonically with the intensity of tFUS. Applying our D_sym measure to the synthetic data produces a monotonic dose-response for irreversibility. It does not reproduce the non-monotonic, moderate-dose peak that we found in our study. We report this preliminary finding here because it speaks directly to the reviewers’ central critique.

Below we provide details of the planned revision.

(1) A surrogate model of GCaMP’s response to tFUS does not reproduce the observed dose structure

We built a generative model of the measured calcium signal according to the underlying physics (Author response image 1A):

(1) A latent Poisson spike train whose rate during and after sonication is monotone in tFUS dose by construction, such that any non-monotonicity appearing downstream must follow from the measurement of irreversibility;

(2) Convolution with an indicator kernel having fast rise and slow decay, matching the kinetic asymmetry mentioned by the reviewers;

(3) A static saturating nonlinearity representing GCaMP6s saturation;

(4) Additive measurement noise applied after the nonlinearity, with standard deviation calibrated to the variance of the real data’s baseline periods;

(5) baseline z-scoring, windowing, and estimation of D_sym using the same parameters as the manuscript (m = 3, τ = 8, f_s = 32 Hz; 5 s baseline, 5 s stimulation, 10 s recovery).

We ran four separate versions of the simulation (with and without saturation, noise before and after the nonlinearity). Simulated GCaMP traces for the variant with saturation and post-nonlinearity noise are drawn in Author response image 1B.

Result. All four simulation variants produce ∆D_sym that increases monotonically with dose and is largest at the highest dose (Author response image 1C). During the recovery window, the difference between moderate (2.3 to 3.7 W/cm²) and maximum (7.4 W/cm²) dose is −0.297±0.066 for the full version with the nonlinearity, −0.222±0.070 without the saturating nonlinearity, and −0.006±0.003 for the two versions in which noise precedes the nonlinearity (mean ± SD across seeds). In contrast, the D_sym of the real data (reproduced for convenience in Author response image 1D) peaks at 3.7 W/cm².

Author response image 1.

A synthetic model of the calcium measurement chain does not reproduce the observed dose structure. (A) Generative model. A Poisson spike train whose rate is proportional to acoustic dose is convolved with a fast-rise, slow-decay indicator kernel, passed through a static saturating nonlinearity, corrupted by additive noise, and then analysed with the same windowing and ordinal estimator used in the manuscript. (B) Simulated traces at low, moderate, and high dose for the variant with saturation and post-nonlinearity noise. (C) Simulated change in irreversibility during recovery, for all four variants of the simulation. Every variant increases monotonically with dose. Markers denote the mean across 20 seeds (error bars show standard deviation). (D) Real data as shown in the manuscript: observed change in irreversibility during recovery for on-target sonication. Markers denote individual animals; boxes show the mean ± 1 SEM. Irreversibility peaks at moderate dose.

(1.1) Invariance of ordinal statistics

Ordinal pattern methods are invariant under any strictly monotonic pointwise transformation of the signal, because such a transformation preserves rank order. Two consequences follow, and we will state them in the Methods of the revision:

- The static component of the indicator nonlinearity, including saturation and the relationship between calcium concentration and fluorescence, cannot alter D_sym.

- The scaling of amplitude with dose does not affect D_sym.

Any potential confound of GCaMP kinetics on D_sym enters through the dynamics, which is what the simulations above address.

(1.2) Planned analyses

In the revision, we will add:

- A sweep over rise and decay constants, Hill coefficient, noise level, and spike statistics, reporting explicitly the region of parameter space, if any, in which a moderate-dose peak can be produced;

- Burstiness and adaptation in the latent source;

- Calibrating the amplitude of the simulated responses to the real data;

- Finite-sample bias and variability of D_sym as a function of window length, pattern count, and signal-to-noise ratio, addressing Reviewer 1’s questions about window duration and noise.

(2) Aspects of the data that argue against a confound of kinetics

Independent of simulation, certain results in the manuscript are not explainable from an effect of calcium indicator kinetics.

(2.1) Non-monotonic irreversibility is reported for off-target sonication – without an accompanying amplitude response. Off-target sonication produces no appreciable fluorescence response, meaning that there is no evoked transient whose kinetics could generate ordinal asymmetry, yet irreversibility still varies non-monotonically with dose. We will elaborate on the off-target EPR effects in the revision.

(2.2) The dose structure remains after adjustment for amplitude. Our nested model comparison shows that acoustic dose explains variance in D_sym beyond baseline irreversibility and the windowmatched change in fluorescence. This is inconsistent with the view of D_sym as a deterministic readout of the amplitude trajectory. In the revision, we will also include decay time constant as a covariate, as Reviewer 1 suggests.

(2.3) Prediction from the prestimulation window. Baseline D_sym is computed on a window that contains no stimulus and no evoked transient, and it predicts the subsequent calcium response after adjustment for baseline fluorescence and dose. The kinetics of a response that has not yet occurred cannot modulate the value of a quantity measured before it.

(3) Reviewer 2, point 3: variation in baseline EPR across nominal intensity

Reviewer 2 correctly observes that baseline irreversibility varies across nominal intensity by an amount comparable to the stimulation-induced change. This implies two candidate explanations, which we address separately.

Systematic differences between intensity groups. If trials at different nominal intensities were not evenly distributed through a recording session, then drifts in bleaching, arousal, or habituation would produce spurious baseline differences between intensity groups. We tested this directly in the 408 analysed CMT 2.5 Hz trials. Session position, defined as the index of a trial within its recording, does not differ across the six intensity levels: Kruskal-Wallis H(5) = 3.51, p = 0.62 on-target and H(5) = 1.70, p = 0.89 off-target, with mean session position varying only between 18.3 and 23.3 across intensities. Within individual recordings, the rank correlation between session position and intensity has a median of +0.03 and does not exceed 0.34 in magnitude in any of the 10 recordings. Acoustic intensity is therefore not confounded with position in the session, and we will report this in the revision.

Uncertainty larger than the plotted error bars. This is a fair concern and we will address it by reporting trial-resampling and pattern-count bootstrap confidence intervals for D_sym. We also note that our primary inference already includes baseline D_sym as a covariate in an ANCOVA specification.

(4) Reviewer 1, remaining points

Other conventional GCaMP metrics. The revision will add decay time constant and additional waveform descriptors to the amplitude-adjustment models described above in 2.2.

Stimulation periods of differing duration. We will report the dependence of D_sym on window length from the finite-sample bias analysis using synthetic data (see 1.2 above).

Lower signal-to-noise recordings. We will measure the dependence of D_sym on SNR with synthetic data.

The auditory confound. We will add a Discussion paragraph citing Sato et al. (2018), Guo et al. (2023) and Kop et al. (2024), stating that the responses analysed here may reflect an indirect auditory mechanism.