## Abstract

A model of two coupled phase oscillators is studied, where both oscillators are subject to random forces but only one oscillator is repetitively stimulated with a pulsatile stimulus. A pulse causes a reset, which is transmitted to the other oscillator via the coupling. The transmission time of the cross-trial (CT) averaged responses, i.e. the difference in time between the maxima of the CT averaged responses of both oscillators differs from the time difference between the maxima of the oscillators' resets. In fact, the transmission time of the CT averaged responses directly corresponds to the phase difference in the stable synchronized state with integer multiples of the oscillators' mean period added to it. With CT averaged responses it is impossible to reliably estimate the time elapsing, owing to the stimulus' action being transmitted between the two oscillators.

## 1. Introduction

Numerous approaches are used to study brain connectivity (see, for example, Breakspear & Stam 2005; Eichler 2005; Harrison *et al*. 2005; Kaminski 2005; Koenig *et al*. 2005; Riera *et al*. 2005; Robinson *et al*. 2005). Cross-trial (CT) averaging, averaging over an ensemble of post-stimulus responses relative to stimulus onset, is widely used for noise reduction of biological signals, such as electro-encephalogram (Dawson 1954; Chiappa 1983) and magnetoencephalography (MEG) signals (Hämäläinen *et al*. 1993; Hari & Salmelin 1997) as well as local field potentials (Steriade 1990). The CT averaged signal of the signal *x*_{j} of the *j*th oscillator reads(1.1)where the stimulus is repetitively administered at *l* different onset times *τ*_{1},*τ*_{2},…,*τ*_{l}. The assumption behind the triggered averaging is that the response *x*_{j} can be decomposed into a stereotypical evoked response *e*_{j}, which follows the stimulus with a constant delay, plus additive Gaussian noise *ξ*_{j}, so that(1.2)is fulfilled (Dawson 1954; Hämäläinen *et al*. 1993). In this case, averaging improves the signal-to-noise ratio by , where the number of responses *l* typically equals 20–300, and (Dawson 1954; Hämäläinen *et al*. 1993).

To study neural information processing and, in particular, the flow of information in interacting neural populations, it is important to analyse the transmission of stimulus-locked responses within networks of interacting neural populations. Transmission processes are typically assessed by determining the difference in time between the occurrence of marker events (like maxima or minima) of the individual averaged responses belonging to different neuronal populations (Dawson 1954; Steriade 1990; Hämäläinen *et al*. 1993).

The averaging assumption from equation (1.2) is violated by stimulated brain activity:

On-going oscillations are abundant in the brain (Hari & Salmelin 1997). Such oscillators are not ‘silent’ during the pre-stimulus period.

Evoked responses, detected with CT averaging, result from reorganizing part of such ongoing oscillations, especially by resetting their phase dynamics (Makeig

*et al*. 2002).Noise is inevitably inherent in neuronal action and, hence, cannot be modelled appropriately by simply adding it to the deterministic signal (Gielen & Moss 2001).

Coupled phase oscillators, subject to pulsatile stimulation and noise, may display an anti-phase CT response clustering, occurring after a stereotypical reset (Tass 2002, 2003

*a*,*b*). CT response clustering means that the oscillators switch between qualitatively different responses across trials. In a recent MEG study, it has been shown that simple visual stimuli make neuronal populations generate a late anti-phase CT response clustering after an early stereotypical reset (Tass*et al*. submitted).

Phase oscillators (see below) may be used as a minimal model for stimulated oscillatory brain activity since they exhibit all of these basic dynamical features (Tass 2003*a*,*b*). This paper is dedicated to studying the transmission of stimulus-locked responses in two coupled phase oscillators. The motivation behind this approach is to investigate the dynamical processes that are essential for neuronal information processing, in a reasonably simple oscillator model. Insights into the transmission dynamics can be used for the analysis of experimental data.

## 2. Stochastic model

### (a) Two phase oscillators

We consider a model of two phase oscillators with phases *ψ*_{1} and *ψ*_{2} obeying(2.1)(2.2)where *j*, *k*=1, 2 and *j*≠*k*. The eigenfrequencies *ω*_{1} and *ω*_{2} fulfil *ω*_{1}−*ω*_{2}=*γ* with detuning *γ*. *K* is the positive coupling constant. In biological systems the effect of a stimulus is typically phase dependent (Winfree 1980; Tass 1999). Accordingly, the stimulus of oscillator 1 is modelled by a 2*π*-periodic, time-independent function , which, for the sake of simplicity, is chosen to be:(2.3)

The results presented here are also obtained in the case of more complex type of stimuli and asymmetric coupling (see §5). Switching on and off the stimulus of oscillator 1 is modelled by:(2.4)

The random forces *F*_{1} and *F*_{2} are Gaussian white noise fulfilling and with constant noise amplitude *D*. *δ*_{jk} is Kronecker's delta, i.e. *δ*_{jk}=1 if *j*=*k* and 0 else. Dirac's delta function is defined by: for and for arbitrary *ϵ*>0. Equations (2.1) and (2.2) serve as a minimal model for two rhythmically active cortical populations of neurons, one of them being directly stimulated via the sensory input pathway (like V1), the other one being stimulated via transmission from the first area (like V3; for more details see Tass 2004*a*,*b*,*c*). We set the amplitude of both oscillators equal to 1, so that the signal of the *j*th phase oscillator reads:(2.5)

## 3. Cross-trial analysis based on stochastic phase resetting

We introduce normalized phases(3.1)and the *normalized cyclic n* : *m phase difference*(3.2)

We want to detect whether in an ensemble of responses to the stimulus there are epochs during which the phases *ϕ*_{1}, *ϕ*_{2} and/or the phase difference *φ*_{n,m} display a stereotypical, tightly stimulus-locked time-course. To this end, we deliver a series of *l* identical stimuli to oscillator 1 at random times *τ*_{1},*τ*_{2},…,*τ*_{l}. The length of the inter-stimulus intervals (ISI) is randomized according to(3.3)where the minimal ISI, *t*_{win}, is constant and large compared with the stimulation duration and the time-scale of the transient dynamics. *ζ*_{k} is uniformly distributed in [0,*k*2*π*/*ω*], where *k* is a small integer. We attach an identical time window [*t*_{a},*t*_{b}] to each stimulus (*t*_{a}<0, *t*_{b}>0), where the onset of the stimulus in each window lies in *t*=0. The window length *t*_{b}−*t*_{a} is smaller than the length of the ISI (*t*_{b}−*t*_{a}<*t*_{win}), but large compared with the time-scale of the transient dynamics.

To evaluate the dynamics of the ensemble of stimulus-locked responses statistically, we collect the values for *ϕ*_{j} and *φ*_{n,m} across all trials for each time *t* relative to stimulus onset. Accordingly, for each time *t*∈[*t*_{a},*t*_{b}] we introduce the time-dependent *CT distributions* of the normalized phases from equation (3.1) and the cyclic *n* : *m* phase difference from equation (3.2) by:(3.4)

The time-course of *ϕ*_{j} and *φ*_{n,m} is perfectly stimulus-locked at time *t* if the corresponding CT distributions from equation (3.4) are Dirac-like distributions, i.e. and for all *i*, *k*=1,…,*l*. If *ϕ*_{j} and *φ*_{n,m} are not at all stimulus-locked at time *t*, these distributions are uniform. The extent of stimulus locking of *ϕ*_{j} and *φ*_{n,m} is quantified for each time *t* by means of the time-dependent *resetting index ρ*_{j}(*t*) of *ϕ*_{j} given by(3.5)(Tallon-Baudry *et al.* 1996; Tass 2002, 2003*a*,*b*) and the *n* : *m* *synchronization index σ*_{n,m}(*t*) of *φ*_{n,m} given by(3.6)(Tass 1999, 2002, 2003*a*,*b*) where |*y*| denotes the modulus of *y*, and *ν* is an integer. *ρ*_{j}(*t*) and *σ*_{n,m}(*t*) detect whether the CT distribution of *ϕ*_{j} and *φ*_{n,m} from equation (3.4) at time *t* has one prominent peak, respectively. 0≤*ρ*_{j}(*t*)≤1, 0≤*σ*_{n,m}(*t*)≤1 are fulfilled for *t*∈[*t*_{a},*t*_{b}].

The first and the 99th percentile of the pre-stimulus distributions of the indices and serve as confidence levels in order to determine whether a stimulus causes a significant increase or decrease of the corresponding index (Tass 2002, 2003*a*,*b*).

## 4. Transmission of CT averaged responses

### (a) Transmission of resets

To investigate how a reset of oscillator 1 is transmitted to oscillator 2, we assume that the coupling is strong enough compared with the noise amplitude *D*. Without stimulation the two oscillators spontaneously synchronize in-phase (figure 1*g*,*h*) because of their in-phase coupling (equations (2.1) and (2.2)). The stimulus of oscillator 1 is modelled by , where the stimulation intensity *I* is large compared with both coupling strength *K* and noise amplitude *D* (*K*≪*I*, *D*≪*I*). Therefore, oscillator 1 is quickly reset by the strong stimulus: *ϕ*_{1} is shifted close to (figure 1*a*). This reset is reflected by an increase of the resetting index *ρ*_{1} from equation (3.5) (figure 1*b*). The CT averaged signal from equation (1.1) vanishes during the pre-stimulus period because of the randomized stimulus administration (figure 1*c*). The reset leads to an oscillatory CT averaged signal , which relaxes to zero owing to noise.

The stimulus of oscillator 1 perturbs the strong synchronization between both oscillators. Owing to the coupling between the two oscillators, the phase of oscillator 2 becomes adapted to the phase of oscillator 1 within roughly a period. Before the reset of oscillator 2 reaches its maximum, the two oscillators undergo a transient epoch of desynchronization and resynchronization (figure 1*g*,*h*): as soon as both oscillators are fully resynchronized, the reset of oscillator 2 is maximal (figure 1*e*). In this way the reset of oscillator 1 is transmitted onto oscillator 2 via the coupling. The reset of oscillator 2 corresponds to an increase of the resetting index *ρ*_{2} (figure 1*e*) and causes an oscillatory CT averaged signal , which vanishes in the course of the post-stimulus period because of noise (figure 1*f*). Since the reset of oscillator 2 occurs by transmission via the coupling, it is delayed and less pronounced compared with the reset of oscillator 1.

### (b) Transmission time

Let us consider how the transmission time of the reset, that is, the time elapsing between the maximal reset of oscillator 1 and the maximal reset of oscillator 2, depends on the coupling strength *K*. We compare it with the transmission time of the CT average: the time between the maximum of the CT averaged signal of oscillator 1 and the maximum of the CT averaged signal of oscillator 2. We normalize both quantities by dividing by the mean period of the oscillators, , and obtain the *normalized transmission time of the reset*(4.1)where is the time at which the resetting index *ρ*_{j} from equation (3.5) is maximal. The *normalized transmission time of the CT average* reads(4.2)with being the time at which the CT averaged signal from equation (1.1) is maximal. corresponds to the standard method for the detection of transmission times, as used in medicine and biology (Dawson 1954; Chiappa 1983; Steriade 1990; Hämäläinen *et al*. 1993).

With increasing coupling strength *K*, the normalized transmission time of the reset, , decreases gradually (figure 2*a*). For values of *K* greater than the intensity parameter *I* from equation (2.3), finally converges to zero; the two oscillators then behave like one composite oscillator (*I*=40 in figure 2). This convergence is not shown in figure 2; rather we focus on phenomena occurring for weaker coupling, since the latter is more realistic for applications to biology and medicine. Unlike over the whole range of the coupling constant (for 0≤*K*≤15), takes only discrete values (except for minor fluctuations): 0, 1, and 2. With an increase of *K*, jumps from 2 to 1 and finally to 0, with an overlap of the different levels of .

To understand the discrete dependence of on *K*, let us consider the time-course of for intermediate coupling (*K*=2.02) and strong coupling (*K*=11.3) (figure 2*b*,*c*). According to their coupling mechanism (equations (2.1) and (2.2)), both oscillators spontaneously tend to synchronize in-phase. A maximum of the CT averaged response requires that the phase of oscillator 2 is tightly time-locked to the stimulus across trials. Hence, a maximum of occurs, when the reset of oscillator 2 is maximal. After being reset, the first oscillator acts like a reference oscillator during several cycles. A reset of oscillator 2 means that oscillator 2 has adapted its phase to the reference oscillator. In other words, the reset of oscillator 2 requires a complete in-phase resynchronization of the two oscillators (figure 2*d*,*f*). As a consequence of the particular form of the coupling, which induces an in-phase synchronization, the maxima of and are separated by multiples of their period . The quicker the resynchronization takes place, the stronger the coupling is. Accordingly, the maxima of and coincide for strong coupling, and (figure 2*c*). In contrast, for intermediate coupling it takes one more period to attain its maximum of (figure 2*b*).

## 5. Discussion

In this paper, a model is presented which allows the study of basic transmission properties of stimulus-locked responses in two coupled phase oscillators, where only one oscillator is stimulated. Furthermore, I have explained how to detect such dynamics reliably with data analysis techniques based on stochastic phase resetting (§3). The transmission time of the CT averaged responses, that is, the difference in time between the maxima of the CT averaged responses of both oscillators, directly corresponds with the phase difference in the stable synchronized state with integer multiples of the oscillators' mean period added to it, where the integer value of added periods depends on the coupling strength (figure 2).

Thus, the transmission time of the averaged responses primarily corresponds with features of the coupling mechanism and the detuning. Put otherwise, the transmission time of the CT averaged responses basically represents spontaneous synchronization properties of two oscillators, that is, the outcome of the interaction of two oscillators in the absence of stimulation. In particular, the transmission time of the CT averaged responses is not a quantity that reflects the time elapsing owing to the stimulus action being transmitted between the two oscillators. This contradicts the assumption used in the evoked response literature (Dawson 1954; Steriade 1990; Hämäläinen *et al*. 1993; Hari & Salmelin 1997). In contrast, with the stochastic phase resetting analysis from §3, the transmission time of the stimulus-locked responses is reliably assessed by detecting the length of time between the maximal resets of both oscillators.

The results presented are also valid in general phase oscillator models, especially in the case of more complex coupling and stimulation terms obtained, for example, by introducing a constant shift term, such as (Tass 2004*a*,*b*,*c*). Furthermore, the results are also obtained in the case of non-identical coupling strength in equations (2.1) and (2.2). Obviously, the only necessary condition is that oscillator 2 is coupled to the directly stimulated oscillator 1. Otherwise there is no way for the stimulus' effect to be transmitted.

The detuning *γ*=*ω*_{1}−*ω*_{2} influences the stable synchronized state. Given a coupling strength that is strong enough, a non-vanishing detuning will shift the phase difference between the two oscillators in the stable synchronized state. Let us recall that the transmission time of the CT averaged responses directly corresponds with the phase difference in the stable synchronized state plus integer multiples of the oscillators' mean period. Accordingly, the detuning changes the transmission time of the CT averaged responses.

Unlike misleading estimates based on CT averaging, the stochastic phase resetting analysis (§3) enables a reliable detection of transmission times of stimulus-locked responses. For the technical details necessary for an experimental application refer to Tass (1999, 2002, 2003*a*,*b*, 2004*a*,*b*,*c*).

## Acknowledgments

This study was supported by the Volkswagen Foundation (76761) and the German Israeli Foundation (I-667-81.1/2000).

## Footnotes

One contribution of 21 to a Theme Issue ‘Multimodal neuroimaging of brain connectivity’.

- © 2005 The Royal Society