DriPP: Driven Point Processes to Model Stimuli Induced Patterns in M/EEG Signals
Authors: Cédric Allain
Alexandre Gramfort
Thomas Moreau
Venue: ICLR 2022 Poster (2022)
Keywords: Electrophysiology
Neuroscience
Temporal point processes
Convolutional Dictionary Learning

Abstract: The quantitative analysis of non-invasive electrophysiology signals from electroencephalography (EEG) and magnetoencephalography (MEG) boils down to the identification of temporal patterns such as evoked responses, transient bursts of neural oscillations but also blinks or heartbeats for data cleaning. Several works have shown that these patterns can be extracted efficiently in an unsupervised way, e.g., using Convolutional Dictionary Learning. This leads to an event-based description of the data. Given these events, a natural question is to estimate how their occurrences are modulated by certain cognitive tasks and experimental manipulations. To address it, we propose a point process approach. While point processes have been used in neuroscience in the past, in particular for single cell recordings (spike trains), techniques such as Convolutional Dictionary Learning make them amenable to human studies based on EEG/MEG signals. We develop a novel statistical point process model – called driven temporal point processes (DriPP) – where the intensity function of the point process model is linked to a set of point processes corresponding to stimulation events. We derive a fast and principled expectation-maximization algorithm to estimate the parameters of this model. Simulations reveal that model parameters can be identified from long enough signals. Results on standard MEG datasets demonstrate that our methodology reveals event-related neural responses – both evoked and induced – and isolates non-task specific temporal patterns.

Published as a conference paper at ICLR 2022
DRIPP: DRIVEN POINT PROCESSES TO MODEL STIM- ULI INDUCED PATTERNS IN M/EEG SIGNALS
Cédric Allain, Alexandre Gramfort & Thomas Moreau Université Paris- Paris-Saclay, Inria, CEA, Palaiseau, 91 91120, France (cedric. cllliin, lexandrre elexandre.gramf fort, moomassmoreuull (inriai
ABSTRACT
The quantitative analysis of non ion-invasive electrophysiology signals from electroen- cephalography (EEG) and magnetoencephallgraphh (MEG) boils down to the identification of temporal patterns such as evoked responses, transient bursts of neural oscillations but also blinks or heartbeats for data cleaning. Several works have shown that these patterns can extracted efficiently in an unsupervised way, e.g., using Convolutional Dictionary Learning. This leads an event-based description of the data. Given these events, a natural question is to estimate how their occurrences nipulations. To address it, we propose a point process approach. While point processes have been used in neuroscience in the past, in particular for single cell recordings (spike trains), techniques such as Convolutional Dictionary Learning make them amenable human studies based on EEG/MEG signals. We develop a novel statistical point process model called driven temporal point processes (DriPP) - where the intensity function of the point process model is linked to a set of point processes corresponding to stimulation events. We derive a fast and principled iopectation-maxiiizaan of this model. Simulations reveal that model parameters can be identified from long enough signals. Results on standard MEG datasets demonstrate that our methodology reveals event-related neural responses both evoked and induced and isolates non-task-specine temporal patterns.
modulated by certain cognitive tasks and experimental ma-
(EM) algorithm to estimate the parameters
INTRODUCTION
Statistical analysis of human neural recordings is at the core of modern neuroscience research. Thanks to non-invasive recording technologies such as electroencephaloorrppyy (EEG) and mag- netoencephalograply (MEG), or invasive techniques such as electrocorticography (ECoG) and stereotactic EEG (sEEG), the ambition is to obtain detailed quantitative description of neural signals at the millisecond timescale when human subjects perform different cognitive tasks (Baillet, 2017 70. During neuroscience experiments, human subjects are exposed to several external stimuli, and we interested in knowing how these stimuli influence neural activity.
After pre-processee rre-processing steps, such filtering or Independent Component Analysis (ICA: Winkler et al. 2015) to remove artifacts, common techniques rely on epoch averaging to highlight evoked responses or time-frequency analysis quantify power changes in certain frequency bands (Cohen, 2014) for induced responses. While these approaches have led to numerous neuroscience findings, has also been criticized. Indeed, averaging tends to blur out the responses due small jitters in time-locked responses, and the Fourier analysis of different frequency bands tends to neglect the harmonic structure of the signal, leading to the so-called "Fourier fallacy" (lasper, 1948; Jones, 2016). In doing, one may conclude to a spurious correlation between components that have actually the same origin. Morcover, artifact removal using ICA requires a tedious step of selecting the correct components.
Driven by these drawbacks, a recent trend of work aims to go beyond these classical tools by isolating prototypical waveforms related to the stimuli in the signal (Cole & Voytek, 2017; Dupré la Tour et al., 2018; Donoghue et al., 2020). The core idea consists in decomposing neural signals as combinations of time-invariant patterns, which typically correspond transient bursts of neural activity (Sherman

---

Published as a conference paper at ICLR 2022
et al., 2016), or artifacts such as eye blinks or heartbeats. In machine learning, various unsupervised algorithms have been historically proposed to efficiently identify patterns and their locations from multivariate temporal signals images (Lewicki & Sejnowski, 1999; Jost et al., 2006; Heide et al., 2015; Bristow et al., 2013; Wohlberg, 20), with applications such as audio classification (Grosse et al., 2007) image inpainting (Wohlherg, 2016a). For neural signals in particular, several methods have been proposed to tackle this task, such as the sliding window matching (SWM: Gips et al. 2017), the learning of recurrent waveforms (Brockmeier & Príncipe, 2016), adaptive waveform learning (AWL; Hitziger et al. 2017) or convolutional dictionary learning (CDL; Jas et al. 2017; Dupré la Tour al. 2018). Equipped with such algorithms, the multivariate neural signals are then represented by set of spatio-temporal patterns, called atoms, with their respective onsets, called activations. Out of all these methods, CDL has emerged as a convenient and efficient tool to extract patterns, in particular due to its ability to easily include physical priors for the patterns to recover. For example, for M/EEG data, Dupré la Tour et al. (2018) have proposed a CDL method which extracts atoms that appertain to electrical dipoles in the brain by imposing rank-1 I structure. While these methods output characteristic patterns and an event-based representation of the temporal dynamics, it is often tedious and requires certain domain knowledge to quantify how stimuli affect the atoms' activations. Knowing such effects allows determining whether an atom is triggered by specific type of stimulus, and if so, quantify by how much, and with what latency.
As activations are random signals that consist of discrete events, a natural statistical framework is the one of temporal point processes (PP). PP have received a surge of interest in machine learn- ing (Bompaire, 2019; Shchur et al., 2020; Mei et al., 2020) with diverse applications in fields such as healthcare (Lasko, 2014; Lian al., 2015) or modelling of communities on social networks (Long et al., 2015). In neuroscience, PP have also been studied in the past, in particular to model single cell recordings and neural spike trains (Truccolo et al., 2005; Okatan et al., 2005; Kim al., 201 1; Rad & Paninski, 2011), sometimes coupled with spatial statistics (Pillow et al., 2008) network models (Galves& & Löcherbach, 2015). However, existing models do not directly address our question, namely, the characterization of the influence of deterministic PP the stimuli onsets stochastic one the neural activations derived from M/EEG recordings.
In this paper, we propose novel method called driven point process (DriPP) to model the activation probability for CDL. This method is inspired from Hawkes processes (HP; Hawkes 1971), and models the intensity function of stochastic process conditioned the realization of set of PP. called drivers, parametrized using truncated Gaussian kernels to better model latency effects in neural responses. The resulting process can capture the surge of activations associated to external events, thus providing a direct statistical characterization of how much a stimulus impacts the neural response, as well as the mean and standard deviation of the response's 's latency. We derive an efficient -atiattiiimmmmmmiiitiiinoo (EM) based inference algorithm and show on synthetic data that it reliably estimates the model parameters, even in the context of MEG/EEG experiments with tens to hundreds of events at most, Finally, the evaluation of DriPP on the output of CDL for standard MEG datasets shows that reveals neural responses linked to stimuli that can be mapped precisely both in time and in brain space. Our methodology offers unified approach to decide if some waveforms extracted with CDL are unrelated to a cognitive task, such as artifacts or spontaneous brain activity or if they provoked by stimulus no matter if they are 'evoked' or 'induced' as more commonly described in the neuroscience literature (Tallon-Baudry et al., 1996). While these different effects are commonly extracted using different analysis pipelines, DriPP simply reveals them stimuli-induces neural responses using a single unified method, that does not require any manual tuning or selection.
2 DRIVEN TEMPORAL POINT PROCESS (DRIPP)
temporal point process (PP) is stochastic process whose realization consists of discrete events {t+} occurring in continuous time, t; E IR+ (Daley & Vere-Jones, 2003). In case where the probability that an event occurs at time only depends on the past events F, == {ta 6, < t), PP are usually
characterized through the conditional intensity function X R + R+:
P PNetde N₂ 1)
(t]F)) = lim
(1)
dr
dt
2

---

Published as a conference paper at ICLR 2022
where Nr E1>1 LesS is the counting process associated the PP. This function corresponds to the expected infinitesimal rate at which events are occurring at time i given the arrival times of past events prior (Daley Vere-Jones, 2003).
The proposed model DriPP is adapted from the Hawkes process (HP; Hawkes 1971), as the occurrence of a past event the driver increases the likelihood of occurrence activation events the near future. However, here we suppose that the stochastic point process in our model of neural activations does not have the self-exxct behavior characteristic of HP. Instead, the sources of activation in the DriPP model are either the drivers some spontaneous background activity, but not its previous activations. More specifically, in DriPP, the intensity function at time / between a stochastic process k whose set of events is denoted A and non-empty set of drivers P whose events
denoted T == (t(P) byg ( J.PE E P is composed of a baseline intensity μide 2 0 and triggering kernels Kk.pp IR R:
(L. = +
Σ PEP S
DR (1- (1-1,) (,))))
where Okkpp 0 is coefficient which controls the relative importance of the driver p the occurrence of events on the stochastic process k. Note that when the driver processes are known, the intensity function is deterministic, and thus corresponds to the intensity of an inhomogeneous Poisson process (Daley & Vere-Jones, 2003). The coefficient Og,, is set to be non-negative so that we only model excitatory effects, as events on the driver only increase the likelihood of occurrence of new events on the stochastic process. Inhibition effects are assumed non-existent. Figure illustrates how events Tp on the driver influence the intensity function after a short latency period.
A critical parametrization of this model is the choice of the triggering kernels KA.p. To best model the latency, we use a parametric truncated normal distribution of mcan myp E IR and standard deviation Kk.p 0, with support (a: b] RR,
Namely,
2 minz dk.y
(x))(x)) = Ks (2;; m..p, k..p. a,
Lasesb
(3)
K..p +
m-mma ak.p
db
a- a-mmag
CV.P
where & (resp. 4) is the probability density function (resp. cumulative distribution function) of the standard normal distribution. This parametrization differs from the usual exponential kernel usually considered in HP, that captures responses with low latency. Note that the truncation values a, b R are supposed independent of both the stochastic process and the drivers, hence they are similar for all kernel p P. Indeed, in the context of this paper, those values delimit the time interval during which a neuronal response might occur following an external stimulus. In other words, the interval [a b] denotes the range of possible latency values. In the following, we denote by T == (k) the duration of the process
PARAMETERS INFERENCE WITH AN EM-BASED ALGORITHM
We propose to infer the model parameters Oş.P (μk, Qg.P, P, bold the vector version of the parameter, i.e., Ft.P = (Tk.p) via maximum maximum-likelihood using an
R..)), M.P..Ok.P)
where we denote in
EM-based algorithm (Lewis& & Mohler, 2011; Xu et al., 2016). The pseudocode of the algorithm is presented in Algorithm The expectation- matiiiiiiimmiiiiiiiiiininn (EM) algorithm (Dempster et al., 1977) an iterative algorithm that allows to find the maximum likelihood estimates (MLE) of parameters in probabilistic model when the latter depends on non-observable latent variables. First, from (2), we derive the negative log-likelihood of the model (see details in Appendix A.
CK.P (OkP)) T
Σ up PEP
Σ log tEAA
1k +
Σ PEP 1
Σ
Kkepkkppp (
(1-((ot)
(4)
Expectation step For a given estimate, the E-step aims at computing the events' assignation, Le., the probability that event comes from either the kernel the baseline intensity. At iteration
3

---

Published as a conference paper at ICLR 2022
Signal
Atom
Activation A&
Stimulus T Intensity Ap.p
14
(K) 2
(1)
(K)
2
Convolutional dictionary learning (CDL) applied univariate signal (blue) de- composes it as the convolution of a temporal pattern (orange) and sparse activation signal (black). Intensity function Ae.p defined by its baseline μl and the stimulus events Tp (green). Intensity increases following stimulation events with certain latency.
Figure Top:
Bottom:
let p(") (t) Pk (t: k)) be the probability that the activation at time t E ;0; has been triggered
by the baseline intensity of the stochastic process k, and Pa (t) pf() (t; p, k) be the probability that the activation at time has been triggered by the driver p. By the definition of our intensity model
(2), we have:
(n)
R
(n)
Eg,000
(0) 21
KR.P
(n)
( (.7))
(1)
and E (1)
P(
(5)
(1)
(n)
where g(x) denotes the value of the parameter f at step of the algorithm, and similarly. if f is func- tion parameter 0, f(m) (x; == (a; 0(x)). Note that Vt E [0;T], (( + Eper PL (1)
Maximization step Once this assignation has been computed, one needs update the parameters of the model using MLE. To obtain the update equations, we fix the probabilities p(n) and p(a) and cancel the negative log-likelihood derivatives with respect to each parameter. For given values of
probabilities (())
(t) and ((*)
(t), we here derive succinctly the update for parameters i and a::
OLLP (x)
0pk
k.
Σ tEA.
(t)
=0
P! (t)
(t)
+
(n+1) Pk
T
1
Σ
pll)
(t)
BEAr
tEA,
(6)
OCLP
(m)
(0(1)
k.P
(n)
= 00
p(n)
(a+1) 0 + k.
= Rgl
1 E (())(t) Thp tE.A,
(7)
tE.A,
k.
where (-) denotes the projection onto the set E. These two updates amount to maximizing the probabilities that the events assigned to the driver or the baseline stay assigned to the same generation process. Then, give the update equations for m and σ, which corresponds parametric estimates of each truncated Gaussian kernel parameter with events assigned to the kernel. Detailed computations provided in Appendix A.2.
(1) ElEAe Eipese (I-1(1) ME.
(n)
(1-15))
(a+1) map
(n)2 Cm
mip
(n) Kep a,
(n)
(8)
ElEAL ()(
k.
(n)
mk,p' k.p a
1/3
( (1) (t) (1)
(n) (1n) (n)) C m k.p hip a,b)
Okp ELEdi Eigd'
(n)
(a+1) Tlettoo)
( (1)
k.p
C. (n)
(n)
EIEA, p(m) (t)
k σ k. a.
4

---

where X (-) denotes the Lebesgue measure, and where DR.p = (1- LE Ar) An) (a; the set of all empirical delays possibly linked to the driver with ((2) (t) = max(t', E Tp 1' < denoting the timestamp of the last event on the driver p that occurred before time Here, the baseline intensity (0) is set the average number of process' events that occur outside any kemne support, i..., the events that guaranteed to be exogenous or spontaneous. Similarly, the kernel intensity (9) computed as the increase in the average number of activations over the kernel support, compared to (")), The initial guess for (o) and (U) are obtained with their parametric estimates, considering Published as a conference paper at ICLR 2022
where, C (m, a, = Sa exp (- (M (u-m²²) du, Cm (,, =
ac 3m
(m, J, a, b) and
o²
Co (m,o,a, b) == ac do (m, a, a, b). Here = > 0 is predetermined to ensure that or remains strictly positive, In practice, we such that we avoid the overfitting that can occur when the kernel's mass is too concentrated. Note that once the initial values of the parameters are determined, the EM algorithm is entirely deterministic.
Also, when the estimate of parameter m is too far from the kernel's support (a b], we are in pathological case where EM is diverging due to indeterminacy between setting a 0 and pushing m to infinity due to the discrete nature of our events. Thus, consider that the stochastic process is not linked to the considered driver, and fall back to the MLE estimator. The algorithm is therefore stopped and we set ar = Oz#r.
Algorithm 1: EM-based algorithm
input :Ag. Tp. a, b, T
output: : The estimated values for parameters
Hr, at, m and
Initialize (0) (0) m()) (())
(0)
for i 0. N do
a = Og then
((+1) (MLE);
break;
It is worth noting that if Vp P, PR then the intensity is reduced to its base- line. thus the negative log-likelihood is Cx.p (0k.p)
end
Define 1 Compute
10+11 (i+1), m(+ (1+1);
P&T #AA log lk. where #A denotes the cardinality of the set A. Thus, we can termi- nate the EM algorithm by directly computing the MLE for Mk+ namely: (MLE) # Mk
end
return ( (it1)), m ((++))
Initialization strategy We propose "smart start" initialization strategy, where parameters are initialized based on their role in the model. It reads:
# 4+ # (UPED
1k T (UPEP UPET, I1' +a;' +
(0)
(10)
(0)
#DK
14
(0)
Ap EP
(11)
A (UPET [t' + a;tt' +b])
(0)
Σ and
(0) kkp
Σ dEDu.p
ld - m0) VppP (12)
mk.p
#Dk.P dEDE.P
#Due
that all event the kernel support assigned to the considered driver.
4 EXPERIMENTS
We evaluated our model on several experiments, using both synthetic and empirical MEG data. We used Python (Python Software Foundation, 2019) and its scientific libraries (Virtanen et al., 2020; Hunter, 2007; Harris et al., 2020). We relied alphaese for CDL with rank-1 constraints on MEG (Dupré la Tour et al., 2018) and we used MNE (Gramfort et al., 2013) to load and manipulate the MEG datasets. Computations were on CPU Intel(R) Xeon(R) E5-2699, with 44 physical cores.
5

---

Published as a conference paper at ICLR 2022
j4
0.8, 0.8, 0.4, 0.2
0.8, 0.8, 0.4, 0.03
2.5
Ground truth Estimated
6-
1.5
0.0
0.4
0.6
0.8
1.0
0.0
0.2
0.6
0.8
1.0
Time (s)
Time (s)
Figure 2: True and estimated intensity functions following driving event time zero for two different kernels, on synthetic data. Left: "wide" kernel with о = 0.2. Right: "sharp" kernel with 0.05. Parameters used T 10000, P/S 0.6. On synthetic data, the EM algorithm successfully retrieves the true values of parameters, for both shapes of kernels.
4.1 EVALUATION OF THE EM CONVERGENCE ON SYNTHETIC DATA
For a given number of drivers and set of corresponding parameters O we first generate the drivers' processes and then simulate the stochastic process for pre-determined duration T. Each driver's timestamps simulated as follows: given interstimuli interval (ISI), set of S Lisi equidistant timestamps is generated - where [-] denotes the floor function. Then P timestamps are uniformly sampled without replacement from this set. In all our experiments, we fixed the ISI to for the "wide" kernel, and to 1.4s for the "sharp" one. Finally, one-dimensional non-homogeneous Poisson process is simulated following Lewis' thinning algorithm (Lewis & Shedler, 1979), given the predefined intensity function À and the drivers' timestamps.
Figure illustrates the intensity function recovery with two drivers considered together: the first one has "wide" kernel with standard deviation a 0.2s, and the second one has a "sharp" kernel with = s. Both kernels have support [0.03 0.8s] and mean m = 0.4 s, the coefficients a are both set to 0.8 and the baseline intensity parameter ys 0.8. We report 8 estimated intensities obtained from independent simulations of the processes using 7 = 10000 and P/S = 0.6 - that we plot over each one of the driver's kernel's support. The EM algorithm run for 50 iterations using the "smart start" initialization strategy described in Section 2. Note that here, the randomness only comes from the data generation, as the EM algorithm uses deterministic initialization. Figures demonstrate that the EM algorithm is able to successfully recover the parameters for both shapes of kernels.
To provide quantitative evaluation of the parameters' recovery, we compute, for each driver p P. the foo norm between the intensity X computed with the true parameters θ, and the estimated intensity Ap with parameters Gp:
|-- - 1p11 max |1 + af,t((t) - app((tp(() rE(0:]]
(13)
The rationale for using the fo norm is to ensure that errors during baseline and within the kernel support are given equal importance. Figure presents the parameter recovery for the same scenario with varying P/S and T. To easily compare the EM performances on the two shapes of kernels, Figure (resp. Figure A. reports the mean (resp. standard deviation) of the relative foo norm that is the foc divided by the maximum of the true intensity X computed for each of the driver over 30 repetitions with different realizations of the process. The results show that the more data are available, either due to longer process duration (increase in 7') or due to higher event density (increase in P/S), the better are the parameter estimates. The convergence appears to be almost linear in these two cases. Moreover, the average computation time for an EM algorithm in Figure took 18.16 18.16s. s showing the efficiency of our inference method. In addition, we report in appendix the scaling of the EM computation time function of T' in Figure A.2.

---

Published as a conference paper at ICLR 2022
0.8, 0.8, 0.4, 0.2
# 0.8, 0.8, m 0.4,
0.05
100
Mean || |lle /mnar
P/S
0.1
0.3
10-1
0.6
10a
10a
10a
102
10c
104
T (s)
T (S)
Figure 3: Mean of the relative infinite norm function of process duration T and the percentage of events kept P/S, for two kernel shapes on synthetic data: wide kernel (left) and sharp kernel (right). The accuracy of the EM estimates increases with longer and denser processes.
4.2 EVOKED AND INDUCED EFFECTS CHARACTERIZATION IN MEG DATA
Datasets Experiments on MEG data were run on two datasets from MNE Python package (Gramfort et al., 2014; 2013): the sample dataset and the somatosensory (somato) dataset'. These datasets were selected as they elicit two distinet types of event-related neural activations: evoked responses which are time locked to the onset of the driver process, and induced responses which exhibit random jitters. Complementary experiments were performed on the larger Cam-CAN dataset (Shafto et al., 2014)². Presentation of the dataset, data pre- pre-processing and obtained results 3 subjects are presented in Appendix A.7. The presented results are self- refermetrrieed as they exhibit, for each subject, the atoms that have the higher ratio a/ a/μ. For all studied datasets, full results are presented in supplementary materials.
The sample dataset contains M/EEG recordings of human subject presented with audio and visual stimuli. In this experiment, checkerboard patterns are presented to the subject in the left and right visual field, interspersed by tones to the left right The experiment lasts about 4.6 min and approximately 70 stimuli per type are presented to the subject. The interval between the stimuli is on average of 750 ms, all types combined, with a minimum of 593 ms. Occasionally, a smiley face is presented at the center of the visual field. The subject was asked to press button with the right index finger as soon as possible after the appearance of the face. the following, we are only interested in the four main stimuli types: auditory left, auditory right, visual left, and visual right. For the somato dataset, a human subject is scanned with MEG during 15 min, while 111 stimulations of his left median nerve were made. The minimum ISI is 7s.
Experimental setting For both datasets, only the 204 gradiometer channels are analyzed. The signals are pre-processed using high-pass filtering at 2 Hz to remove slow drifts in the data, and are resampled to 150 Hz to limit the atom size in the CDL. CDL is computed using a al alphacsc (Dupré la Tour et al., 2018) with the GreedyCDL method. For the sample dataset, 40 atoms of duration 1s each are extracted, and for the somato dataset, 20 atoms of duration 0.53 s are estimated. The extracted atoms' activations are binarized using threshold of 6 10-11 (resp. 1 x 10-10) for sample (resp. somato), and the times of the events are shifted to make them correspond to the peak amplitude time of the atom. Then, for every atom, the intensity function is estimated using the EM-based algorithm with 400 iterations and the "smart start" initialization strategy. Kernels' truncation values are parame-arameterss for the EM and thus must be pre-determined. The upper truncation value b is chosen smaller than the minimum ISI. Here, we used in addition some previous domain knowledge to set coherent values for each dataset. Hence, for the sample (resp. somato) dataset, kernel support is fixed at [0.03 s; 0.5 (resp. ;0s; 2s). See Appendix A.4 for an analysis on how these hyperparameters influence on the obtained results presented below.
Evoked responses in sample dataset Results on the sample dataset are presented in Figure 4. We plot the spatial and temporal representations of four selected atoms, as well as the estimated intensity functions related to the two types of stimuli: auditory (blue) and visual (orange). The first two atoms
Both available at https::/ /tttpp::///n...ee http://eee.tooo...ooooooooooooos wtwpooootttele////eteeoeeeeeeeee tttp://tttteeeeeeeeeeeeetett//ee indd: index. mndex.htm html Available https: //////www /www. saaa-can.oogg Index wttp:e.ppppcooooeeeeddattttoeeee

---

Published as a conference paper at ICLR 2022
Atom 0
Atom
Atom 2
Atom
Spatial
Temporal
0.0
0.5
0.5
0.5
0.5
40
auditory
0.085
0.188
Intensity
visual
20
0
0
0.25
0.5
0.25
0.5 0
0.25
0.5
0.25
0.5
Time (s)
Time (s)
Time
Time (s)
Figure 4: Spatial and temporal patterns of atoms from sample dataset, and their respective estimated intensity functions following stimulus (cue at time = 0 s), for auditory and visual stimuli. The heartbeat and eye-blink artifacts are not linked to any stimuli. An auditory stimulus will induce neural response similar atom 2, with a mean latency of 85 ms.
are specifically handpicked to exhibit the usual artifacts, and the last two are selected as they have the bigger ratios a/1 for their respective learned intensity functions. Even though the intensity is learned with the two stimuli conjointly, we plot the two corresponding "intensities at the kernel separately, i.e., Vp P, Vt (0:).5), we plot (t), 0, 2,
Spatial and temporal representations of atom 0 (resp. atom 1) indicate that it corresponds to the heartbeat (resp. the eye blink) artifact. These two atoms thus expected not be linked any stimuli. This is confirmed by the shape of the intensities estimated with DriPP that is mostly flat, which indicates that the activation of these two atoms are independent of auditory and visual stimuli. Note that these two artifacts can also be recovered by an Independent Component Analysis (ICA), as shown Figure A.5. Indeed, the cosine similarity between the spatial maps of the eye blink (resp. the heartbeat) artifact extracted with CDL and its corresponding component in ICA analysis is 99.58 % (resp. 99.78 %0, as presented in Figure A.6. contrast, by looking at the spatial and temporal patterns of atom (resp. atom 6), it can be associated with auditory (resp. visual) evoked response. Given the spatial topography atom 2, we conclude to a bilateral auditory response and the peak transient temporal pattern suggests evoked response that is confirmed by the estimated intensity function that contains a narrow peak around 85 ms post-stimulus. This is the M100 response - here the auditory one well known in the MEG literature (its equivalent in EEG the N100) (Näätänen & Picton, 1987). The M100 is indeed a peak observed in the evoked response between 80 and 120 milliseconds after the onset of stimulus in an adult population. Regarding atom 6, topography is right lateralized in the occipital region, suggesting visual evoked response. This is confirmed by the intensity function estimated that reports a relationship between this atom and the visual stimuli. Here also, the intensity peak is narrow, which is characteristic of an evoked response. This reflects right lateralized response along the right ventral visual stream this subject. This may connected to the P200, peak of the electrical potential between 150 and 275 ms after visual onset. Moreover, the intensities estimated with DriPP for the unrelated tasks are completely flat. We have a = O, which indicates that atoms' activations exogenous spontaneous relatively to unrelated stimuli. For comparison, we present in Appendix A.5 similar results obtained with dedicated MEG data analysis tools, such evoked responses and fime-ncyyyency plots.
8

---

Published as a conference paper at ICLR 2022
Atom 0
Atom 2
Alom 7
Intensity
0.60
Atom 0
Spatial
0.55
2
0.50
0.2
0.45
Temporal
0.0
0.40
-0.2
0.35
0
0.25
0.5 0
0.25
0.5 0
0.25
0.5 0
1
2
Time (s)
Time (s)
Time (s)
Time (s)
Figure 5: Spatial and temporal patterns of atoms from somato dataset, and their respective estimated intensity functions following a somatosensory stimulus (cue at time = =0s) s). The eye-blink artifact (atom 0) is not linked to the stimulus, and neither is the a-wave (atom 7). A somatosensory stimulus will induce neural response similar to atom 2, with mean latency of 1 s.
Induced response in somato dataset Results on the somato dataset are presented in Figure 5. Similar to the results on sample, spatial and temporal patterns of handpicked atoms are plotted alongside the intensity functions obtained with DriPP. Thanks to their spatial and temporal patterns, and with some domain knowledge, it is possible to categorize these 3 atoms: atom 2 corresponds to fı-wave located in the secondary somatosensory region (S2), atom 7 corresponds to an a-wave originating in the occipital visual areas, whereas atom corresponds to the eye-blink artifact. As a-waves spontaneous brain activity, they that atom is not linked to the task, as confirmed by its estimated intensity function where & = 0. yr-aae -, its respective intensity is nonflat with broad peak close to s. which characterizes an induced response. Moreover, similar to results on the sample dataset, we recover the eye-blink artifact that also has flat intensity function. This allows us be confident in the interpretation of the obtained results. Some other μı-aave atoms atoms and 4 are presented in Figure A010 in Appendix A.6. They have an estimated intensity similar to atom i.e, non-flat with broad peak close to s. The usual time/frequency analysis reported in Figure A.9 exhibits the induced response the pa-avve.
phase-locked to the stimuli. It is thus expected
For atom 2 that corresponds
5 DISCUSSION
This work proposed a point process (PP) based approach specially designed to model how external stimuli can influence the occurrences of recurring spatio-temporal patterns, called atoms, extracted from M/EEG recordings using convolutional dictionary learning (CDL). The key advantage of the developed method that by estimating few parameters (one baseline parameter and 3 parameters per considered driver), it provides direct statistical characterization of when and how each stimulus is responsible for the occurrences of neural responses. Importantly, it can achieve this with relatively limited data which is well adapted to MEG/EEG experiments that last only few minutes, hence leading to tens or hundreds of events at most. This work proposed an EM algorithm derived for PP
novel kernel function: the truncated Gaussian, which differs from the usual parametrization models that capture immediate responses, e.g., with exponential kernels. As opposed to competing methods that can involve manual selection of task-related neural sources, DriPP offers unified approach to extract waveforms and automatically select the ones that are likely to be triggered by the considered stimuli. Note however that DriPP has been developed based on a point process framework. which is weent-baedd. When working with continuous stimuli, other techniques must be considered (e.g.. STRF, ppatio-tempor response functions; Drennan & Lalor 2019). Future work will explore the modeling of inhibitions effects (a 0) frequently observed in neuroscience data. Yet, this could potentially degrade the interpretability of the model parameters, and it also requires a new definition of the intensity function one needs ensure it stays non-negative at all times.

---

Published as a conference paper at ICLR 2022
REFERENCES
S. Baillet. Magnetoenceehhalogappyy for brain clectrophysiology and imaging. Nature Neuroscience, 20:327 EP -, 02 2017.
M. Bompaire. Muchine learning based on Hawkes processes and stochustic optimization. Theses, Université Paris Saclay (COmUE), July 2019.
H. Bristow, A. Eriksson, and S. Lucey. Fast convolutional sparse coding. In Computer Vision and Pultern Recognition (CVPR), pp. 391-398, 2013.
A. J. Brockmeier and J. C. Principe. Learning recurrent waveforms within EEGs. IEEE Transactions on Biomedical Engineering, 63(1):43-54, 2016.
M. X. Cohen. Analyzing Neural Time Series Data: Theory and Practice. The MIT Press, 2014.
S. R. Cole and B. Voytek. Brain Oscillations and the Importance of Waveform Shape. Trends in Cognitive Sciences, 21(2):137-149, 2017.
D. J. Daley and D. Vere-Jones. An introduction to the theory of point processes. Volume I: Elementary theory and methods. Probability and Its Applications. Springer- Verlag New York, 2003.
A. Dempster, N. M. Laird, and D. B. Rubin. Maximum likelihood from incomplete data via the EM algorithm. Journal of the Royal Statistical Society: Series B (Methodloogical)) 39(1): 1977.
T. Donoghue, M. Haller, E. Peterson, Varma, P. Sebastian, Gao, Noto, Lara, D. Wallis, R. T. Knight, A. Shestyuk, and B. Voytek. Parameterizing neural power spectra into periodic and aperiodic components. Nature Neuroscience, 23512):1655-1655 2020.
D. P. Drennan and E. C. Lalor. Cortical tracking of complex sound envelopes: modeling the changes in response with intensity. eneuro, 6(3), 2019.
Dupré la Tour, Moreau, M. Jas, and A. Gramfort. Multivariate convolutional sparse coding for electromagnetic brain signals. Advances in Neural Information Processing Systems, 31:3292-3302 2018.
A. Galves and E. Löcherbach, Modeling networks of spiking neurons as interacting processes with memory of variable length. arXiv preprint a50. 2015.
B. Gips, A. Bahramisharif, E. Lowet, M. Roberts, de Weerd, Jensen, and J. van der Eerden. Discovering recurring patterns in electrophysiological recordings. J. Neurosci. Methods, 275: 66-79, 2017.
A. Gramfort, M. Luessi, E. Larson, D. A. Engemann, D. Strohmeier, C. Brodbeck, Goj, M. Jas, T. Brooks, L. Parkkonen, et al. MEG and EEG data analysis with MNE-Python. Frontiers in neuroscience, 7:267, 2013.
A. Gramfort, M. Luessi, E. Larson, D. A. Engemann, D. Strohmeier, C. Brodbeck, L. Parkkonen, and M. S. Hämäläinen. MNE software for processing MEG and EEG data. Neuroimage, 86:446-460, 2014.
R. Grosse, R. Raina, H. Kwong, and A. N. Shift-invariant sparse coding for audio classification. In 23rd Conference on Uncertainty in Artificial Intelligence (UAI), pp. 149-158, AUAI Press, 2007. ISBN 0-9749039-3-0.
C. Harris, K. Millman, J. van der Walt, Gommers, P. Virtanen, D. Cournapeau, E. Wieser, J. Taylor, S. Berg, N. J. Smith, R. Kern, M. Picus, S. Hoyer, M. H. van Kerkwijk, M. Brett, A. Haldane, del Río, M. Wiebe, Peterson, P. Gérard-Marchaat,, K. Sheppard, Reddy, W. Weckesser, H. Abbasi, C. Gohlke, and T. F. Oliphant. Array programming with NumPy. Nature, 585(7825):357-362 September 2020.
A. G. Hawkes. Point spectra some mutually exciting point processes. Journal of the Royal Statistical Society: Series B (Methodological), 33(3):438 443, 1971.
10

---

Published as a conference paper at ICLR 2022
F. Heide, W. Heidrich, and G. Wetzstein. Fast and flexible convolutional sparse coding. In Computer Vision and Pattern Recognition (CVPR), pp. 5135-5143. IEEE, 2015.
S. Hitziger, M. Clerc, S. Saillet, C. Bénar, and Papadopoulo T. Adaptive waveform learning: A framework for modeling variability in neurophysiological signals. IEEE Transactions Signal Processing, 65(16):432444338 2017.
D. Hunter. Matplotlib: 2d graphics environment. Computing in Science Engineering, 9(3):90-95, 2007.
M. Jas, T. Dupré La Tour, U. Şimşekli, and A. Gramfort, Learning the morphology of brain signals using alpha-stable convolutional sparse coding. In Advances in Neural Information Processing Systems (NIPS). pp. 1-15, 2017.
H. H. Jasper. Charting the sea of brain waves. Science, 108(2805):343-4444 1948.
R. Jones. When brain rhythms aren't i 'rhythmic": implication for their mechanisms and meaning.
Curr. Opin. Neurobiol., 40:72 00772--02220 2016.
P. Jost, P. Vandergheynst, S. Lesage, and R. Gribonval. MoTIF: an efficient algorithm for learning translation invariant dictionaries. Acoustics, Speech and Signal Processing (ICASSP), volume IEEE, 2006.
S. Kim, D. Putrino, S. Ghosh, and E. N. Brown. A granger causality measure for point process models of ensemble neural spiking activity. PLOS Comput Biol, 7(3)=100111 1001110, 201
T. A. Lasko. Efficient inference of gaussianpproeess-mone to medical event data. In Uncertainty in artificial intelligence: proceedings of the... conference. Conference on Uncertainty in Artificial Intelligence, volume 2014, pp. 469. NIH Public Access, 2014.
renewal processes with application
M. S. Lewicki and T. J. Sejnowski. Coding time-varying signals using sparse, shift-invariant representations. M. Kearns, S. A. Solla, and D. A. Cohn (eds.), Advances in Neural Information Processing Systems (NIPS), pp. 730-736. MIT Press, 1999.
PA W Lewis and Gerald S Shedler. Simulation of nonhomogeneous Poisson processes by thinning. Naval research logistics quarterly, 26(3):403-413, 1979.
R. Lewis and G. Mohler. nonparametric EM algorithm for multiscale Hawkes processes. Journal of Nonparametric Statistics, 1(1):1-20, 2011.
W. Lian, R. Henao, V. Rao, J. Lucas, and L. Carin. multitask point process predictive model. In International Conference on Machine Learning. pp. 2030-2038. PMLR, 2015.
T. Long, F. Mehrdad, S. Le, and Z. Hongyuan. NerCodec: Community Detection from Individual Activities, pp. 91 91-99...1 SIAM International Conference on Data Mining (SDM), 2015.
H. Mei, T. Wan, and J. Eisner. Noise-contrastive estimation for multivariate point processes. In
H. Larochelle, M. Ranzato, R. Hadsell, M. F. Balcan, and H. Lin (eds.), Advances in Neural Information Processing Systems, volume 33, pp. 5204-5214. Curran Associates, Inc., 2020.
R. Näätänen and T. Picton. The NI wave of the human electric and magnetic response to sound: review and an analysis of the component structure. Psychophysiology, 24(4):375-425. 1987.
M. Okatan, M. A. Wilson, and E. N. Brown. Analyzing functional connectivity using a network likelihood model of ensemble neural spiking activity. Neural computation, 1799:1927-1161, 2005.
W. Pillow, J. Shlens, L. Paninski, A. Sher, A. M. Litke, E. J. Chichilnisky, and P. Simoncelli. Spatio-temporal correlations and visual signalling in complete neuronal population. Nature, 454 (7207)9995 999, 2008
Python Software Foundation. Python Language Reference, version 3.8. http://python.org/ 2019.
K. Rad and Paninski. Information rates and optimal decoding in large neural populations. Advances in neural information processing systems, 24:846-854 201 1.
11

---

Published as a conference paper at ICLR 2022
M. A. Shafto, L. K. Tyler, M. Dixon, J. R. Taylor, B. Rowe, R. Cusack, J. Calder, W. D. Marslen- Wilson, J. Duncan, Dalgleish, et al. The Cambridge Centre for Ageing and Neuroscience (Cam-CAN) study protocol: cross-scctional, lifespan, multidisciplinary examination of healthy cognitive ageing. BMC neurology, 14(1):1-25, 2014.
O. Shchur, N. Gao, M. Biloš, and S. Glinnemann. Fast and flexible temporal point processes with triangular maps. In H. Larochelle, M. Ranzato, Hadsell, M. F. Balcan, and H. Lin (eds.), Advances in Neural Information Processing Systems, volume 33, pp. 73-84. Curran Associates, Inc., 2020.
M. A. Sherman, S. Lee, R. Law, S. Haegens, C. A. Thorn, M. S. Hämäläinen, I. Moore, and S. R. Jones. Neural mechanisms of transient neocortical beta rhythms: Converging evidence from humans, computational modeling, monkeys, and mice. Proceedings of the National Academy of Sciences,
2016.
C. Tallon-Baudry O. Bertrand, C. Delpuech, and J. Pernier. Stimulus specificity of phase-locked and non non-phase-locker 40 hz visual responses in human. Journal of Neuroscience, 16(13):4240-4444 1996.
W. Truccolo, U. Eden, M. R. Fellows, Donoghue, and N. Brown. point process framework for relating neural spiking activity to spiking history, neural ensemble, and extrinsic covariate effects. Journal of neurophysiology. 93(2):0074-1089 2005.
Virtanen, R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau, E. Burovski, P. Peterson, W. Weckesser, J. Bright, S. van der Walt, M. Brett, Wilson, K. Millman, N. Mayorov, A. R. J. Nelson, E. Jones, R. Kern, E. Larson, C J Carey, i. Polat, Y. Feng, E. W. Moore, J. VanderPlas, D. Laxalde, Perktold, R. Cimrman, I. Henriksen,
A. Quintero, Harris, A. M. Archibald, A. H. Ribeiro. F. Pedregosa, P. van Mulbregt. and SciPy 1.0 Contributors.
SciPy 1.0: Fundamental Algorithms for Scientific Computing Python. Nature Methods, 17:
261 271-272, 2020.
I. Winkler, S. Debener, K.-R. Müller, and M. Tangermann. On the influence of high-pass filtering on ICA-based artifact reduction in EEG-ERP. 2015 37th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), pp. 4101-4105. IEEE, 2015.
B. Wohlberg. Convolutional sparse representation of color images. In IEEE Southwest Symposium on Image Analysis and Interpretation (SSIAI), pp. 57-60, 2016a.
B. Wohlberg. Efficient algorithms for convolutional sparse representations. Image Processing, IEEE Transactions on, 25(1):301-311 2016b.
H. Xu, M. Farajtabar, and H. Zha. Learning granger causality for Hawkes processes. International conference on machine learning, pp. 1717-1726, 2016.
12

---

Published as a conference paper at ICLR 2022
A APPENDIX
A.1 DETAILS OF THE NEGATIVE LOG-LIKELIHOOD COMPUTATION
We defined the negative log-likelihood for parameters kG, (Juk: M,
m,,P; (k,P) as follows:
(k,P(OKP)) = L/¹ At.p(t)dt -
Σ tEAx
log (t)
(14)
We show here in details that 50¹ Xn,p (t) = MAT + ΣpEp ak. Without loss of generality,
can first assume that Vp E P, (p)
mn b S T. Hence,
up, ['" KE. (t-) dt
P,Vi- = 1 ,..
Thus, we have that
St (1)( 1"/
Psk Σ
Σ
Ckppkkkpp (1- ((o) de
pEP i.
P r
Σ Okk
Σ
Kk.p (1 (1- (1e)) di
pEP
ET
P + Σ
Okpupp
PEP
Finally, we have that the negative log-likelihood writes
(x,p(Op,)) = PxTT +
Σ
Up
Σ BE
log (c,
PEP
A.2 DETAILS EM- BASED ALGORITHM COMPUTATIONS
First, recall that the negative log-likelihood for parameters OL.P (JUEs Mr Th.)) is as follows:
LLP (0kp)) = H
Σ
Σ
log
P4 +
+ Σ
m axppii (t- (/ ( oni))
PEP
teAr
PEP
Recall also that we defined p(n) and Pl as follows:
(n)
(n) (t ak.p Ei St. Kkp (t- ((-((((())
(n)
(())(1)
(1)
and (P
(15)
where ((n) denotes the value of the parameter step of the algorithm, and similarly, if f is function of parameter 0, (x; 0) = O(n).
Below are details of calculation for parameters update for m and a. We first rewrite the kernel function to have it under form that simplifies further computations.
1 Φb b-mm
2-
A(2, (1,m,,,a,0))
Lassh
(m-m²)
exp
lasesb
C C(m,o,a
13

---

Published as a conference paper at ICLR 2022
where
C (m, = {(" exp
( 1 (u m)) o²
du
(16)
Hence, we can precompute its derivatives with respect to m and a:
0 dm
n(x; (0 0, = ( ("0=
Cm (m,a,a.b))
1 a(x;
(17)
o² C (,,,o,a,b)
and
(x mm)
Co am, s(x;m,,,a,b)
x(x; m, o,a.b)
(18)
do
Cm,,o,a,b)
where Cm (m, o, b) = ac 3m
8C (m, J, a, b) and C. am,,, =
do (m,
Update equation for ma.p
akpplik (1-1(()
OL.P (0kp)) = Σ 0mk.p
- (P)
-RR
Cm (mx.p Kk.p. a, C (m..p: Kk.ps a, b)
=
Ak,p(t)
K.p
LEAR i,clp) gt
1-11 ((0) Kipp (1 (1- ((01)
mk.p
+
Cm (m Kk. 0,
Σ
P₂(1)
Okkp Σ Σ
oh
C (mapp kkpp a.
Ak,p(1)
LEA=
k tEAr (,t1 St
(19)
Hence, by canceling the previous derivative,
(1-(()) (()) (1-(())
(-(
Kk.p EIEA: Eiples
Cm (n)
k. kip a.
(n+1) mkp
ErEA, P()
(t)
k.p
(n)
C
m m kp a.
(n) (1))
(20)
Update equation for dk.p
(t-t(0) {(p))
t
2
C (t {(00)
OLLP (0tp)) -- Σ
mK.
Co (mx.ps kh b.
0ok.
k.p
C (ma.p: kh,p a,
de,p(1)
LEA: Ei
((0) (t) mk.
Kk.pp (t- {(00)
Co (m a.p
Σ Pa(1)
Olk,p
C (mp.p Kkpp U,
K.P tEAr i.t!") St
\e,
LEAr
(21)
Hence, by canceling the previous derivative,
(L-(()) (m)-mm)))))) (n) (
( ( (1-)0))
1/3
C m hips b)
(x) 0k
ElEAe Eni
(n+1)
k.p
(n) (12) kpp k.pp a,
EIEA, ((
Co
(22)
Finally, to ensure that the coefficient stays strictly positive in order to avoid computational enrors, (n+1)
we add a projection step onto |e; +o), with e > 0: o (e:+000) k.p
(n+1)
14

---

Published as a conference paper at ICLR 2022
14 0.8, 0.8, m 0.4,
1 0.8, 0.8, m 0.4, 0.05
STD || 1l Dimax
P/S
0.1
10
0.3
0.6
10-2
10a
10 T (s)
104
10a
10a T (s)
10
Figure A.I: Standard deviation of the relative infinite norm as function of the process duration T and the percentage of events kept P/S, for two kernel shapes on synthetic data: wide kernel (Ieft) and sharp kernel (right). The variance of the EM estimates continually decreases with longer and denser processes.
A.3 EXPERIMENTS SYNTHETIIC
In this section, we present figures that are complementary to the ones obtained on synthetic data, presented subsection 4. In Figure A. we report the standard deviation associated with Figure i.e, the standard deviation of the relative infinite norm for different process duration 7 and different proportion of events kept P/S. S, for the two kernel shapes.
Figure A.2 reports some details on the computation time of the experiment done synthetic data which produced the Figure and Figure A. For each value of T, the mean computation time is computed over 90 experiments (3 values of P/S times 30 random seeds). Results are presented in Figure A.2.
CPU time (s)
40
20
102
103
10
T (s)
Figure A.2: Mean and 95 % CI of computation time (in seconds) for one EM algorithm, as a function of the process duration T (in seconds). Results obtained synthetic data.
A.4 IMPACT OF MODEL HYPERPARAMETER
In this section, we dwell on the analysis of how setting hyperparameter values may impact the obtained results, with the aim of determining whether is possible to set these parameters using general rule thumb, without degrading previous results. More specifically, we will look at the impact of two hyperparameters the results we obtained on the MNA sample dataset: the threshold value applied to the atoms' activation values to binarized them, currently set at 6 x 10 10 and the kernel support, currently set [0.03s: 0.5 using previous and domain knowledge. To do so, conducted two experiments on sample where each varies one hyperparameter. We plot the intensity function learned by DriPP for the same four atoms namely, the artifacts (heartheat and eye
15

---

Published as a conference paper at ICLR 2022
blink) and the auditory and visual responses -, separately for the two stimuli (auditory and visual), similarly to Figure 4.
For the first experiment, presented in Figure A.3, we varied the threshold, expressed as percentile, between 0 and 80. A threshold of 20 means that we only keep activations whose values are above the 20 % percentile computed over all strictly positive activations, i.e., the smaller the threshold is, the more activations are kept. The value of the threshold used in all other experiments is the 60 % percentile. One can observe that for the two artifacts (atoms O and I), when the threshold gets smaller, the learned intensity functions get flatter, indicating that the stimulus has no influence on the atom activation. However, for the two others atoms, the effect of a smaller threshold is the opposite, as the intensity functions have higher peak, indicating a bigger value of the parameter, and thus strengthening the link stimulus-atom. Thus, the threshold value could be set to small percentile and therefore computed without manual intervention, without degrading the current results.
For the second experiment, we now focus on the kernel truncation values a and b. Results are presented in Figure A.4. We set = 0, we did for somato and Cam-CAN datasets, and we vary b between 0.5 - the current value - and 10, a large value compared to the ISI. One can observe that for b = 0.5 or b = 1, the results are either unchanged (atoms 1, 2, and 6) or better (atom O, as the intensity function is totally flat for this artifact), indicating that setting a = 0 and close to the average ISI of 0. 0.750s s does not hinder the results. However, when b is too large, the results degrade quickly, up to the point where all learned intensities are flat lines, indicating that our model does not find any link between the stimuli and the atoms. This is due the fact that this hyperparameter is of great importance in the initialization step, as the greater b is, the more atom's 's activations are considered being on kernel support. Thus, setting the upper truncation value to a value close the average ISI seems to give reliable results.
A.5 USUAL M/EEG DATA ANALYSIS
We present in this section some results obtained using usual M/EEG data analysis, such as Independent Component Analysis (ICA), epoch averaging, or time/frequency analysis. First, on MNE sample dataset, we proceed to an ICA to manually identify usual artifacts. To do so, similarly as the CDL pre-processing, the raw signals are filtered (high-pass pals-pass filter at 2 Hz), and 40 independent components are fitted. The two components 1 and 3, that we manually identify as corresponding to the eye blink and heartbeat artifacts, are presented in Figure A.5. In Figure A.6, we associate for each of the CDL atoms presented in Figure 4 the ICA component that has the maximum cosine similarity. One can observe that the artifact atoms and components are highly similar, suggesting that CDL and ICA have equal performance on the artifact detection. Note that this high similarity is only based on the spatial pattern. Indeed, ICA does not provide temporal waveforms for the atoms as well as their temporal onsets, contrarily to CDL.
However, for the auditory and visual response, the result is different. For the auditory one (atom 1), there is not really an ICA equivalent, as it is most correlated with the eye-blink ICA component. Regarding the visual atom (atom 6), there is ICA component that presents a high similarity. While the two related components correspond to neural sources on the occipital cortex, the atom 6 obtained with CSC is more right lateralized, suggesting a source in the right ventral visual stream. Note that unlike ICA, which recover full time courses for each source, CDL also provides the onset of the patterns, which we later use for automated identification of event related components. Finally, this demonstrates that CDL is strong competitor to ICA for artifact identification, while simultaneously enabling to reveal evoked or induced neural responses in an automated way.
As mentioned, the ICA is commonly used to remove manually identified artifacts to reconstruct the original signals free of those artifacts. However, there are two drawbacks of this method, the first one being that some Nomain-related knowledge is needed in order to correctly identify the artifacts. The second drawbacks happen when the signal is reconstructed after the removal of certain components. Indeed, such a reconstruction will lead to a loss of information across all channels, as the artifacts are shared by all the sensors. Thanks to the Convolutional dictionary learning (CDL) that extracts the different artifacts directly from the raw data, method does not suffer from these drawbacks.
Still on MNE sample dataset, we compute from raw data the epoch average following an auditory stimulus (stimuli in the left ear and in the right ear combined) and plot on Figure A.7 the obtained
16

---

Published as a conference paper at ICLR 2022
auditory
Visual
60
60
0
Intensity atom
40
60
20
20
0.0
0.2
0.4
0.0
0.2
0.4
60
60
0
Intensity atom
20
40
60
20
20
0.2
0.4
0.0
0.2
0.4
60
60
0
Intensity atom
20
60
20
20
0.2
0.4
0.0
0.2
0.4
60
60
0
Intensity atom
20
40
60
20
20
0.2
0.4
0.0
0.2
0.4
T (5)
T (s)
Figure A.3: Influence of the threshold (expressed as a percentile) on the obtained results on MNE sample dataset, for the 4 main atoms, for auditory and visual stimulus. The value of the threshold presented in Figure is 60 %. The threshold value has limited impact on obtained results, and thus could be determined with general rule of thumb, as percentile over all activations values.
17

---

Published as a conference paper at ICLR 2022
auditory
visual
40
40
0.5
Intensity atom
1
20
20
5
10
0.2
0.4
0.0
0.2
0.4
40
40
Intensity atom
1
20
20
5
10
0.0
0.2
0.4
0.0
0.2
0.4
40
40
Intensity atom 2
1
20
20
5
10
0.0
0.2
0.4
0.0
0.2
0.4
40
40
Intensity atom 6
1
20
10
0.0
0.2
0.4
0.0
0.2
0.4
T (s)
T (s)
Figure A.4: Influence of the kernel truncation upper bound b (with MNE sample dataset, for the main atoms, for auditory and visual stimulus. By setting too high, all intensity functions are completely flat, indicating total loss of information.
the obtained results
18

---

Published as a conference paper at ICLR 2022
Leye daar don sanll الم
ICA00O
ICA003
EOG o61
5.0
12.5
Time
ICA components
ICA000
ICA003
Figure A.5: The first two components alongside with electrooculography (EOG) signal (top) and their corresponding topomaps (bottom), on MNE sample dataset. ICA000 can be associated to the heartbeat artifact, whereas ICA003 can be associated to the eye-blink artifact.
19

---

Published as a conference paper at ICLR 2022
Atom
Atom
Atom 2
Atom
ICAO
ICA3
ICA3
ICA15
99.78%
99.58%
32.73%
81.46%
Figure A.6: Spatial representation of the four atoms presented in Figure 4 (top) alongside their ICA components with the maximum cosine similarity (value indicated the bottom). Atoms and ICA are computed MNE sample dataset.
evoked signals. Figure A.8 is similar but for the visual stimuli (again, stimuli on the left visual field and on the right visual field are combined).
Gradiometers channels)
0.090
0.165
0.250
120
90
60
30
No.e-116
100
fT/cm
100
Time (8)
Figure A.7: Evoked signals following auditory stimulus (cue at time 0), MNE sample dataset. Baseline correction applied from beginning the data until time point zero.
On somato dataset, in order exhibit the induced response to the somatosensory stimulus on the left median nerve the subject, we perform a time/frequenny analysis, presented on Figure A.9. In order perform this analysis, complex process including the of Morlet wavelets is performed.
20

---

Published as a conference paper at ICLR 2022
Gradiometers (203 channels)
0.092 S
0.181 8
0.361 s
Nov 112
iT/om
-50
0.1
0.0
0.1
0.2
0.3
0.4
Time (8)
Figure A.8: Evoked signals following a visual stimulus (cue at time = 0) on MNE sample dataset. Baseline correction applied from beginning of the data until time point zero.
(1.00 s. 7.7 Hz)
(1.30 s, 7.7 Hz)
1e-21
1.0
0.5
0.0
-1.0
35.00
27.21
21.15
16.44
12.78
7.72
6.00
0.0
0.5
1.0
1.5
2.0
Time (s)
Figure A.9: Time-frequency plane for epoched signals following somatosensory stimulus (cue at time 0), on MNE somato dataset, with overviews at and 1.3 seconds and 7.7 Hz. Baseline correction applied from time point -0.5 until time point zero.
Finally, note that these methods that are commonly used the M/EEG data analysis comprise thoughtfully data pre-processing well manual intervention requiring domain knowledge to identify both artifacts and evoked and induced responses. As the proposed method this paper is composed of unified pipeline, it is significant gain. Indeed, DriPP able to automatically isolate artifacts without prior removal of artifacts using ICA, as well as capture the diversity of latencies corresponding to induced responses.
Regarding the statistical significance of the link between a stimulus and a neural pattern, one can consider performing statistical test, where wishes reject the null hypothesis of independence.
21

---

Published as a conference paper at ICLR 2022
atom id
atom type
p-value auditory
p-value visual
0
heartbeat
2.25 x 10-1
1.19 x 10-1
eye-blink
5.85 x 10-1
8.48 x 10-1
auditory
2.31
10-97
6.31 10-1
6
visual
7.78 x 10-1
5.12 x 10-50
Table On MNE sample dataset, statistical univariate Student I-test, Ho: independence between an atom and the stimulus (auditory visual). The atom ids correspond to the ones presented in Figure 4.
We have performed a statistical test using student t-test to check if the mean activation probability on [a, b] segments the same as the mean estimated on the baseline [o, 7] U.114 + a, + bj. We present the results of the performed test in Table 1. We can clearly see that both artifacts (heartbeat and eye-blink) are not linked to stimuli while neural responses are linked to related stimuli (p-values very small). We would like to stress that while this test is interesting, it is much more dependent on the selection of the support interval [a, bj than the proposed method. Moreover, this test does not allow assess the latency of the responses, and whether the atom is an induced response an evoked one.
A.6 EXPERIMENTS ON MNE SOMATO DATASET
We present in this section extra results obtained on the real dataset somato, that are complementary of Figure 5. Figure A.10 shows 3 atoms that correspond all to a μ-wave located in the secondary somatosensory region (S2), with three different shapes of kernels in their estimated intensity functions.
A.7 EXPERIMENTS ON CAM-CAN DATASET
The Cam-CAN dataset contains data of M/EEG recordings of 643 human subjects submitted to audio and visual stimuli. In this experiment, 120 bimodal audio/visual trials and eight unimodal trials included to discourage strategic responding to one modality (four visual only and four auditory only) - are presented to each subject. For each bimodal trial, participants see two checkerboards presented to the left and right of central fixation (34 milliseconds duration) and simultaneously hear a 300 milliseconds binaural tone at one of three frequencies (300, 600, or 1200 Hz, equal numbers of trials pseudorandomly ordered). For unimodal trials, participants either only hear a tone checkerboards. For each trial, participants respond by pressing button with their right index finger if they hear any stimuli (Shafto et al., 2014). For each subject, the experiment lasts less than 4 min.
the
The signals pre-processed using low-pass filtering 125 Hz to remove slow drifts in the data, and are resampled to 150 Hz to limit the atom size in the CDL. CDL is computed using alphacsc (Dupré la Tour et al., 2018) with the GreedyCDL method. Twenty atoms of duration 0.5 are extracted from the signals. The extracted atoms' activations are binarized and, similarly as previous experiments, the events time shifted make them correspond to the peak amplitude time in the atom. Then, for every atom, the intensity function is estimated using the E.M- based algorithm with 200 iterations and the "smart start" initialization strategy. Kernels' truncation values are at 0.. Two drivers per atom are considered. The first driver contains timestamps of all bimodal trials (all frequencies combined) with in addition the auditory unimodal stimuli (denoted catch The second driver is similar to the first one, but instead of the auditory unimodal stimuli, it contains the 4 visual unimodal stimuli (denoted catchl).
This experiment is performed for 3 subjects, for each one plot the 5 atoms that have the highest ratio a p. so that we autonatically exhibit atoms that are highly linked to the presented stimuli. Similarly as results presented in Section 4, we plot the spatial and temporal representation of each atom, as well as the two learned intensity functions that we plot "at kernel". Results are presented on Figure 11 Figure A. and Figure A.13.
22

---

Published as a conference paper at ICLR 2022
Atom 1
Atom 2
Atom 4
Intensity
10-1
Spatial
7 10
Atom
6 10-1
0.2
M
Temporal
0.0
5 x 10-1
-0.2
10
0.25
0.5
0.25
0.50
0.25
0.5 0
1
Time (s)
Time (s)
Time (s)
Time (s)
Atom 1
Atom
Atom
Figure A. Spatial and temporal patterns of 3 μ-wavee atoms from somato dataset, alongside with their respective estimated intensity functions. Below, in order to provide further information, are the corresponding brain locations for each atom obtained with dipole fitting.
First, we can observe that, as for experiments on sample and somato datasets, we recover the heartbeat artifact (atom 0 in Figure A.12) well the eye-blink artifact (atoms 0 in Figure A.11 and in Figure A. 13). Because of the paradigm of the experiment, and in particular the instructions given to the subjects that is, press button after every stimulus, the majority of them being visual it is not surprising if the eye blink artifact is slightly linked to the stimuli. Indeed, it is often observed in such experiments where there no designated time blink that the subject "allows" themself to blink after the visual stimulus.
As the majority of the presented stimuli are a combination of visual and auditory, the CDL model struggles to separate the two corresponding neural responses. Hence, that is why it can be observed that most the first atoms reported exhibit mixture between auditory and visual response their topographies (first row). However, one can observe that for some such atoms, DriPP learned intensity able to indicate what is the main contributing stimulus in the apparition of the atom. For instance, for atom 10 in Figure A. I, the auditory stimulus is the main responsible stimulus of the presence of this atom, despite this latter presenting a mixture of both neural responses. similar analysis can be made for atom 1 in Figure A. 12, but this time for the auditory stimulus.
23

---

Published as a conference paper at ICLR 2022
Atom 10
Atom
Atom
Atom
Atom 9
Spatial
0.25
Temporal
0.00
-0.25
0.00
0.25
0.50.0
0.0 0.0
0.5 0.0
0.5 0.0
0.5
40
audivis catcho
Intensity
audivis_catchi
20
0.45 Time (s)
0.45 Time (s)
0.45 Time (s)
0.45 Time (s)
0.45 Time (s)
Spatial and temporal patterns of the 5 atoms from Cam-CAN dataset subject CC620264, ylarr-y female, and their respective estimated intensity functions following an audiovisual stimulus (cue at time = 0 s). Atoms are ordered by their bigger ratio a/μ.
Figure A.11:
Atom
Atom
Atom
Atom 18
Atom
Spatial
Temporal
0.0
-0.5
0.00
0.25
0.50.
0.5 0.0
0.5 0.0
0.0
0.5 0.0
0.5
Ludivis__atchie
Intensity
50
budivis_catch1 catch1
0.45 Time (s)
0.45 Time (s)
0.45 Time (s)
0.45 Time (s)
Time (s)
Figure A.12: Spatial and temporal patterns the 5 atoms from Cam-CAN dataset subject CC520597, 64.25-year-old male, and their respective estimated intensity functions following an audiovisual stimulus (cue at time = 0 s). Atoms are ordered by their bigger ratio a/p.
24

---

Published as a conference paper at ICLR 2022
Atom
Atom
Atom 2
Atom
Atom
Spatial
Temporal
0.25
0.00
-0.25
0.00
0.25
0.50.0
0.5 0.0
0.5
0.5
0.5
audivis_catchie
Intensity
50
aadivis_catchi
25
0
0.45
0.45
0
0.45
0.45
0
0.45
Time (s)
Time (s)
Time (s)
Time (s)
Time (s)
Figure A.13: Spatial and temporal patterns of the 5 atoms from Cam-CAN dataset subject CC723395, 86. 16.08-yea female, and their respective estimated intensity functions following audiovisual stimulus (cue at time = 0 s). Atoms are ordered by their bigger ratio a
25


Reviews:
[{"review_id": "p9BG_z2JT6N", "reviewer": "Reviewer_RDNE", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "1. The model is simple yet flexible. It captures the relative strength of the linkage through the ratio of the weight of the Gaussian and the baseline rate. The mean and standard deviation of the waveform occurrence timing relative to the stimuli timing is intuitive for neuroscientific experiments. \n\n2. The approach seems to work well in practice. Most of the paper is clear and the figures illustrate the modeling approach.\n\n3. The modeling approach avoids a 'tedious' approach of examining the rate of waveform occurrences in windows following the stimuli occurrence versus background cases (which could be done through point process hypothesis sampling with or without binning the waveform occurrences).\n\nConcerns:\n1. At times the manuscript could benefit from more theoretical or formal analysis as phenomenon are noted but not investigated: \" This suggests that the estimators are slightly biased.\" \n\n2. The pathological case at the beginning of step 5 seems speculative, rather than cast in a rigorous case for convergence.  Shouldn't $m \\in [a,b]$ in order for the Gaussian model of jitter for the timing make sense? Furthermore, wouldn't it make sense to restrict $m\\pm 3\\sigma \\in [a,b]$ for the same reason? If Gaussian jitter is not preferred why not consider other parametric families such as uniform or exponential distributions? These questions undermine the global choice of $a,b$ while $m$ and $\\sigma$ are process and driver dependent. Later on page 7, it appears this choice is due to the inter-stimuli interval. Beside the mention of the degeneracy it is not clear why $m$ couldn't be restricted to lie in this range. \n\n3.  The $\\ell_\\infty$ norm definition needs more clarification.  Equation 15 is a norm between functions on the real line, but its left-hand side definition implies the difference between parameter vectors.  It is not exactly a norm between intensity functions on the inside (which it could be using the known driver timings), but rather the maximum difference in intensity at various times since a driver event. This definition makes the y-axis of Figure 3 hard to interpret. \n\n4. ICA overlap (cosine similarity) is mentioned for eye blink. Was this examined for the auditory or visually evoked responses? \n\n5. Equation 7 the $\\pi_{\\mathbb{R}^+}$ function can be assumed to be a non-negative thresholding operating but it is not stated. Later after equation 10 it is stated as a projection operator to the ensemble $E$. Its not clear to me the choice of word ensemble, as compared to set.\n\n6. I'd like to see a comparison of the proposed approach in determining whether or not a waveform is associated with a stimuli, versus using a test of independence (chi-squared) based on how many waveform occurrences fall into the window $[a,b]$ following a driver process. \n\n7. The paper focuses on modeling the timing of waveforms already estimated through a different process from the original time-series. It would be useful to mention whether a joint model would make sense.  I can imagine that some original waveform occurrences could have been missed due to noise obscuring the waveform, but with the known driver process as well as the linking function additional occurrences may be able to be recovered.  Is this feasible?\n\nMinor point:\n\n\"common techniques relies\" -> \"common techniques rely\"\n \n\n", "summary_of_the_paper": "The paper presents a model of the occurrence timing of waveforms in neural potential recordings in terms of a baseline rate and a set of driver processes. The waveforms and their timing are already extracted through a process of convolutional dictionary learning and convolutional sparse coding. The driver processes are known experimental stimuli timings. The driver to waveform model is done in a point process framework where the driver/stimuli increases the likelihood of the waveform in a window [a,b]. Specifically the paper proses a inhomogeneous Poisson process with a baseline rate and a truncated Gaussian distribution for each driver process. \n\n", "summary_of_the_review": "The paper discusses a modeling approach relevant for understanding how waveforms in neural potential time series relate to known stimuli or background processes. The modeling approach seems consistent and performs well on the given data sets, but wider implications and insights as well as formal analysis are limited. Due to some minor concerns and large questions, my current opinion is borderline and I hope the authors can help clarify any issues raised. \n\n++++++++++++++++++++++=\nUpdate after revision:\nI think my minor concerns and questions were answered. I encourage the authors to carefully proofread the revision. \n\nMy concern with not raising the score more are the limited scope/impact of the work in the learning representation community. Nonetheless, I think that meaningful practical utility may be gained in neuroscience analysis from the proposed approach.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "TZrubHS68sh", "reviewer": "Reviewer_ZvAC", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The methodology is introduced precisely. The inspiring application of the proposed DriPP technique wonderfully extracts peripheral (ECG) and brain-related (ERP components), which are crucial for the proper analysis of M/EEG. The illustrative examples are apparent, and provided code works fine. The authors stated that the current approach utilizes some domain knowledge to set hyperparameter values in their limitations. An exploration of data-driven strategies to infer these hyperparameters shall follow. It would also be interesting to show whether, for example, the extracted ECG/heartbeat atom is periodic after extending the analysis window. A similarity of the proposed method to old-school ICA is somehow disappointing. The authors shall explain this part in more detail since the limitations of ICA applications are somewhat known in the community. ", "summary_of_the_paper": "This is an exciting manuscript. It was a delight to read and evaluate it. The manuscript is well written, and the methodology is presented and validated with illustrative cases. Provided code runs without problems, and it is excellently documented. An introduced point process-based technique seems to be adequately fitted to model brain external/reactive and induced/active/imagery experiments. An influence on M/EEG captured spatiotemporal patterns appears to be also straightforward. The introduced technique extracts spatiotemporal patterns through the utilization of convolutional dictionary learning. The proposed parameterization procedures work thoroughly. The authors clearly explain the limitations of the current approach as related to the present preliminary application.", "summary_of_the_review": "To summarize, the authors propose an interesting DriPP technique to decompose complex M/EEG signals with a bit disappointing twist. The resulting components are 95% similar to what ICA could do, and the ICA is scrutinized at the beginning of the submission, which seems confusing and disappointing.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works.", "empirical_novelty_and_significance": "4: The contributions are significant, and do not exist in prior works."}, "extra_text": {}}, {"review_id": "BIdPFDus7Q1", "reviewer": "Reviewer_vbHx", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Strengths\n-\tInclusion of results from synthetic and empirical data to demonstrate the effectiveness of their method. Inclusion of results from evoked and induced responses support the generalizability of the DriPP model.\n-\tOverall, the paper is very well written and illustrated. Authors provide some details to reproduce experiments.\n\nWeakness\n-  A comparison with other related work (e.g., epoch averaging) is needed to support authors claims that their work addresses limitations of these prior work. There appears to be some similarities with ICA, which is given limited consideration in the text. \n-   Results presented are from limited number of subjects and no quantitative measures are provided to assess model performance. It should be noted that authors include results from subjects from another dataset in the supplementary material, but reviewer relied more on results presented in the main text to make assessment. \n-\tIt is not clear why some model choices are made (e.g., Gaussian kernel, number of atoms, \"smart start\" initialization) and how domain knowledge is used to determine some model parameters (e.g., kernel hyperparameter, authors acknowledge this limitation). Selection of kernel supports seems arbitrary (authors suggest based on ISI). \n-\tHow is the number of atoms and the duration of atoms determined for each dataset? (20 for somato and 40 for sample)\n       “For the sample dataset, 40 atoms of duration 1 s each are extracted, and for the somato dataset 20 atoms of duration 0:53 s are estimated.”\n-\tWhat does the latency of atom 6 (Figure 4) correspond to? \n“Regarding atom 6, topography is right lateralized in the occipital region suggesting a visual evoked reponse. This is confirmed by the intensity function estimated that reports a relationship between this atoms and the visual stimuli.”\n-\tModel considers only excitatory effects and not inhibitory effects (authors acknowledge this limitation). \n-\tTo model a PP, only stimulus onset is considered. So potential for DriPP model to be limited in scenarios with continuous stimuli. \n\nAdditional Feedback\n-\tHave a transition at the end of the introduction to state the main contributions of the paper and the overall goal of the experiment.\n-\tFigure 1: Define CDL in caption\n-\tFigure 3: More informative caption. Label the legend.\n-\tFigure 5: Specify that atom 0 is related to eye blink. What of atom 7? Typo: An somatosensory stimulus will induced\n-\tDescribe the N100 response to relate to the stimulus.\n        “This is the N100 response well known in the MEG/EEG literature.”\n-\tAuthors provide justification for the choice of the triggering kernel (Gaussian) in the Discussion. This should be moved to earlier in the text when the kernel is first mentioned.\n-\tFix typos (not an exhaustive list):\n        “Regarding atom 6, topography is right lateralized in the occipital region suggesting a visual evoked reponse. This is confirmed by the intensity function estimated that reports a relationship between this atoms and the visual stimuli.”\n\n", "summary_of_the_paper": "This paper presents a novel model for revealing underlying stimulus-response patterns in MEG/EEG data. The main contributions include a driven temporal point process (DriPP) model that leverages prior work on convolutional dictionary learning (CDL) and point processes to be amenable to MEG/EEG data, an EM algorithm for inferring the DriPP model parameters and application of the model to synthetic data and empirical data. ", "summary_of_the_review": "The proposed DriPP model is broadly applicable to neuroscience applications to isolate underlying stimulus-driven neural patterns and unrelated temporal patterns (background neural activity, artifacts, etc.).  \n\nAuthors need to provide justifications to support model parameter choices. A comparison with other related work (e.g., epoch averaging) is needed to support authors claims that their work address limitations of these prior work. There appears to be some similarities with ICA, which is given limited consideration in the text. \n\n-------------------\nUpdate post-rebuttal: Reviewer read the author responses and the revised manuscript, as well as other reviewer comments and respective responses. Appreciated the author responses and revisions to the issues raised by this reviewer and other reviewers, and most of the main critiques are addressed. Potential to improve/extend DriPP method in future work. There is an applications track (neuroscience) in this conference, and this work is of relevance to the broader learning representations community. Revising score upwards. Authors should proof-read the manuscript for typos (also see other reviewer comments). \n\n", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "43q_wuUcM58", "reviewer": "Reviewer_87RT", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The paper introduces an interesting technique for the characterization of neuronal events that can be acquired on E/MEG signals. Such an approach could help better understanding and characterization these signals, which contain important information but whose analysis highly depend on the clinician and the analysis technique used.\nI have however several main issues with the study\nThe authors do not seem to compare their technique with any approach from the state-of-the-art, and it is therefore quite hard to judge the performance of the technique.\nThe technique requires the use of several pre-processing technique, and most importantly an event detection technique (the authors suggested and used a CDL for this task). Have the authors assessed how the DriPP performance is dependent on the quality of the event detection?\nThe authors do provide a single quantitative assessment of the technique, but results assembled on figure 3 does not allow the reader to understand which sets of parameters to choose, or what would be an acceptable error.\nThe authors suggest a smart initialization procedure, which is used for all experiments but the effect or benefit of this procedure has not been evaluated.\nThe notation used for the description of the DriPP is quite confusing, and in equation 2, I am not sure to distinguish the sets of timepoints (t_i^k) from A_k and T_p. \n", "summary_of_the_paper": "This paper introduces a novel technique for the analysis of electro or magneto encephalogram based on an innovative point process approach. The point process model is estimated using an expectation maximization algorithm and is aimed at characterising events in the signals representing the brain activity. The technique is finally evaluated on synthetic data and real data from openly available databases. ", "summary_of_the_review": "The paper introduces an interesting technique for the characterization of neuronal events, but I do feel that the technique should be evaluated more thoroughly and compared with other approaches.", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "2: Several of the paper’s claims are incorrect or not well-supported.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper presents a novel method for identifying simuli induced
patterns in MEG and EEG signals.  The authors develop a novel statistical
point process model and a fast EM algorithm to learn the parameters.

Discussion of this paper centered around: how to fit hyperparameters,
and similarity and comparison with other algorithms, especially ICA, as
well as the small number of subjects

Comparison to other methods would make the work stronger, as would
adding more datasets but this novel algorithm seems worth publishing.
I recommend acceptance as a poster.
Decision: Accept (Poster)