Learning Cortico-Muscular Dependence through Orthonormal Decomposition of Density Ratios
Authors: Shihan Ma
Bo Hu
Tianyu Jia
Alexander Kenneth Clarke
Blanka Zicher
Arnault H. Caillet
Dario Farina
Jose C Principe
Venue: NeurIPS 2024 poster (2024)
Keywords: EEG-EMG fusion
statistical dependence
orthonormal decomposition
cortico-muscular connectivity
Primary area: neuroscience_and_cognitive_science

Abstract: The cortico-spinal neural pathway is fundamental for motor control and movement execution, and in humans it is typically studied using concurrent electroencephalography (EEG) and electromyography (EMG) recordings. However, current approaches for capturing high-level and contextual connectivity between these recordings have important limitations. Here, we present a novel application of statistical dependence estimators based on orthonormal decomposition of density ratios to model the relationship between cortical and muscle oscillations. Our method extends from traditional scalar-valued measures by learning eigenvalues, eigenfunctions, and projection spaces of density ratios from realizations of the signal, addressing the interpretability, scalability, and local temporal dependence of cortico-muscular connectivity. We experimentally demonstrate that eigenfunctions learned from cortico-muscular connectivity can accurately classify movements and subjects. Moreover, they reveal channel and temporal dependencies that confirm the activation of specific EEG channels during movement.

Learning Cortico-Muscular Dependence through Orthonormal Decomposition of Density Ratios
Shihan Ma
Bo Hu
Tianyu Jia
Alexander Kenneth Clarke'
Blanka Zicher
Arnault H. Caillet
Dario Farina? José C. Príncipe"
'Equal Contribution
Department of Bioengineering, Imperial College London Department of Electrical and Computer Engineering, University of Florida
Abstract
The cortico-spinal neural pathway is fundamental for motor control and move- ment execution, and in humans it is typically studied using concurrent electroen- cephalography (EEG) and electromyography (EMG) recordings. However, current approaches for capturing high-level and contextual connectivity between these recordings have important limitations. Here, present novel application of statistical dependence estimators based on orthonormal decomposition of density ratios model the relationship between cortical and muscle oscillations. Our method extends from traditional scalar-valued measures by learning eigenvalues, eigenfunctions, and projection spaces of density ratios from realizations the signal, addressing the interpretability. scalability, and local temporal dependence of portcc-muusulla connectivity. We experimentally demonstrate that eigenfunctions learned from cortico-museular connectivity can accurately classify movements and subjects. Moreover, they reveal channel and temporal dependencies that confirm the activation of specific EEG channels during movement. Our code is available at wtt.cooooo/oaooe/ bohpboo/hhha/eee/ cttp:////ooooioooocooooooooao -iia-eeiiinnnner.r...... encoder.
https: //github...t///g//
thub. com
Introduction
The brain communicates with muscles by sending informa- tion the spinal cord. Part of this information is directly transmitted from cortex to spinal motor neurons via the cortico-spinal neural pathway, which is vital for motor control and movement execution. Because motor neurons directly connected to muscles, cortical oscillations tray- eling through the cortico-spinal pathway are coherent with oscillations in muscle electrical activities, both in humans and non-human primates [1-3]. This relationship, known functional cortico-muscular connectivity, critical in neuroscience and is typically studied using concurrent recordings of neural signals such as electroenceehalogra- phy (EEG) and electromyography (EMG) [4, 5] (Fig. I), Practical applications include diagnosing and monitoring of neuromuscular disorders, such as amyotrophic lateral sclerosis [6], stroke [7], and Parkinson's disease [8], well as developing brain-computer interfaces (BCIs) for individuals with motor impairments [9].
informutionn
Motor related
merbori
Cartial
Penpheri
Grtio--nncuarr patlway
dnctt rede
Statitical
BBG
Figure The cortico-muscular pathway allows brain-muscle communication with co- herent cortical and peripheral oscillations. This paper models this connectivity through the statistical dependence between their con- current recordings of EEG and EMG.
38th Conference on Neural Information Processing Systems (NeurIPS 2024). Correspondence to Bo Hu (bohu@uf1. edd) and Shihan Ma (mmasss 1isaasss1200@mmil.ccom

---

Despite several applications, there is still a lack of proper statistical tools to model the relationship between EEG and EMG. The predominant method, Mustullr Coherence (CMC), measures temporal and spectral coherence by computing the normalized cross-spectrum in time intervals [10]. From this analysis, it generally accepted that the EEG's s beta band (13-30 Hz) is linked to steady motor control, while the gamma band (>30 Hz) is associated to motor planning and execution [S, 11-13].
While CMC provides some relevant information on cortico-museular connectivity, there remains lack of generalized and higher-order statistical measures that quantify nonlinear and high-level connectivity. Can high-level contextual information, such as muscle movements and participant identifiers, be directly learned from modeling
connectivity? Our paper explores the potential of using statistical dependence estimators to address this problem.
Statistical dependence estimators typically follow a procedure of defining a measure preferably by probabilistic distributions, deriving a variational bound, and optimizing a variational cost of this bound using function approximator, such as mutual information estimators [14-16] and Kernel Independent Component Analysis (KICA) [17, 18]. These measures are defined for realizations.
However, statistical dependence estimators above have rarely been successfully applied to cortico- muscular analysis, mainly due to three reasons: instability and poor scalability, lack of spatio-temporal resolution, and lack of practical contextual connections. As these estimators typically quantify dependence at the trial level, they overlook the importance of channel and temporal dependence in mousulo-muscular analysis. More importantly, these measures only produce a scalar-valued score, but how this score should be used and its connection to the desired contextual factors are unclear.
This paper successfully applies statistical dependence estimators to EEG-EMG pairs using the concept of orthonormal decomposition of the density ratio. Recently, there has been a shift from scalar-valued measures to decomposing density ratio positive definite function, and learning its eigenvalues, eigenfunctions, and associated projection spaces through neural network optimization with matrix cost funetions such log det and nuclear norm [19-21]. This decomposition, known the Functional Maximal Correlation Algorithm (FMCA), addresses the fundamental issue of relating dependence to contextual information: Eigenvalues define multivariate dependence measure, and eigenfunctions span a feature projection space that captures the contextual factors affecting dependence.
This paper expands on this idea, addressing interpretabilltyy, scalability, and local-Ievel dependence that are missing in existing dependence analyses. Sec. 2.1 explains why eigenfunctions, learned from cortico-museular connectivity, identification. Sec. 2.2 introduces FMCA-T, optimizing a new matrix trace cost for the theory, which demonstrates greater efficiency and stability than the log det cost. Sec. 2.3 shows that while the objective estimates global dependence from trial realizations, localized channel-level and temporal- level dependencies can also be formed in top-down manner, which important in EEG as they indicate channel activations and synchronization of activities. Our framework is illustrated in Fig. 2.
capture contextual factors for motor control and participant
Our main experiment demonstrates that the learned eigenfunctions, without labels, effectively capture factors such as movements and subjects that contribute to high-level cortico-museular connectivity. After training, using EEG's eigenfunction as a feature projector noticeably improves classification accuracy over various baselines. Additionally, channel-level and temporal-level dependencies indicate that specific EEG channels are selectively activated during movements, corroborating neuroscientific findings. Simulated data further confirm that our proposed measure invariant to nonstationary noise, including pink noise and random delays.
2 Methods
2.1 Density ratio decompasition for EEG-EMG signal pairs
Problem formulation. Consider EEG signals X == F.r and EMG signals Y Y 1I.T- Denote as the subject, c as the type of movement, and u other auxiliary contextual factors. These are factors that could potentially affect the statistical dependence between EEG and EMG signals. Each signal is conditioned on these parameters. Denote these factors as == {s, c, u) with distribution P(z). Distributions for EEG and EMG given these conditions are p( X = X|²) respectively. Their joint distribution is given by p(X, Y) p(X pY (1)))p())d Similarly, the marginal distributions are given by p(X) p(X ((X)2))((2 (2)((2)(( and likewise p(Y)) = (2))2.)
(2) and p( Y = Y

---

(e) Optinizing Marris Trace Con
Zr
Ze
(c) Feature Projection Spare
Laspere Navook
22
Chased
sajan
20
Natwoak
Rp
Fox
(X)
Tram: Cost
sc
Yorki 200
Gangeri
Chaned
(0) Channel and Temporal Reseluiion
Narook
Natwook
Rd
6.a
Scx
ZF
Ze
Ziz
Approrimaainn Density Raiio
Clabed
Figure 2: Diagram for learning cortico-museular dependence by decomposing density ratios: (a) Network fa is applied to EEG X1.r and 9., EMG Y 1.r to minimize matrix trace cost. (b) EEG-EMG pairs sampled from a joint distribution, from which a density ratio p(X,, Y) is defined and considered positive definite function. Its linear operator has a spectral decomposition of eigenfunctions {ó, } and eigenvalues (ax)- The networks provably approximate the dominant cigenvalues and cigenfunctions of this dccomposition with
network outputs [f }. and SVD results {Ax Figenvalues here measure multivariate statistical dependence; eigenfunctions are optimal feature projectors. (c) After training. the eigenfunctions, specifically those from EEG, form projection space containing contextual information for motor control and participant identification. (d) To provide channel activation and activity synchronization cortico-musculat analysis, we compute density ratios between channel-leve Ze and temporal-lews features Ze.s against global features Zy to quantify channel-level and temporal-lew dependencies.
Our goal is extract factors s, c, that affect the dependence between modalities, from available sample pairs X and Y, even when 8, c, and u are not given. We propose that this can be achieved by decomposing the density ratio of this probabilistic system.
Decomposition of EEG and EMG density ratios. Following the work FMCA, propose orthonormal decomposition of the density ratio to measure the dependence between EEG and EMG:
Y(X,Y) p:= p(X)p(Y)
Vor (o)ox(X)u
k=1
1 (x)))x (X)
11, if =
(Y) (P() 1, if
0. ifi
0, if
for any i,j 1,2, .... The density ratio p(X, Y) is treated positive definite function associated with a linear operator Lʃ = -( dX for any measurable scalar function f. According to Mercer's theorem, this operator has spectral decomposition with eigenvalues σi, o2 and orthonormal basis functions 61, Ø2. ... and w/, th2 .... In scenarios where X and Y are statistically independent, all eigenvalues are zero. Conversely, larger eigenvalues suggest stronger dependence.
Such eigenfunctions form linear span. Our hypothesis is that this span captures shared contextual factors such as e and s across two modalities, stated as follows.
Lemma 1. Assuming conditional independence given z == {8, c, u), we have p(X, 2|2) p(X.)) p(X (Y)))) Hence, the ratio p(X.Y) == pYX)p(() decomposes as p(X.Y) = p(X(s)n(=) p(Y)=))(=)
p(X)p(z)
p(Y]p(2) -p(z)dz. Assuming z is discrete (e.g.. movement patterns and participant identities 8).
the information of is contained in the span of the basis functions for the density ratio p(X,Y).
Proof. Define the ratios Px(X, = p(X.2)
p())) and pr (Y,z) = p(Yx)
p(Y p(Y)p(z)) Considering the sample
space Z, the sets px (X, and py for are discrete. Under the conditional independence assumption, these ratios satisfy p(X,Y)= E.ez p(=)px (X,
p) (X, (Y,))
This indicates that the sets px (X, =) and py (Y, decompose p(X, Y), similar to the eigenfunctions du and v/a- Since both decompositions represent p( (X,Y), Y), the functions px (X, and pr (Y,)) must lie within the span of these basis functions. Hence, there exist coefficients a=,AA and A2.A such
that px (X, 2) Ex az &x (X) and py (Y, Σx A (Y) for each Thus, learning the
3

---

dependence between X and Y is implicitly learning the dependence between each of them relative the factors 2, even when 2 not observed.
2.2 FMCA-T: Learning decomposition for the matrix trace
When probability densities are unavailable, we approximate eigenvalues and eigenfunctions using learning system with two neural networks and a cost function, typically matrix cost like log det or nuclear norm [19-21]. These costs optimize an aggregation of eigenvalues. The networks learn the dominant eigenvalues and cigenfunctions when optimized.
Aggregation of eigenvalues. To measure the total power the eigenspectrum, define scalar-valued measure using a convex function & R + E with {()) = 0. Assume the eigenvalues are ranked 0i 2 02 2 The truncated total statistical dependence measure of the top K eigenvalues is defined by Te = EE=E((os) ((or). Function {(x) = log (1 x) corresponds to the log det cost.
Prior work: log-determinant cost. Consider two networks, fa :X X R and In : R mapping realizations of X and to K-dimensional outputs, respectively. Assume g for EEG and 9w for EMG. The autocorrelation (ACFs) and cross-correlation functions (CCFs) are defined as:
R PPc
PRG Ra
R Rxlf Rs = = (Y)g) PrG Ex.y RyG
(2)
FMCA minimizes a logdet cost, which reaches the negative value of the total measure Te of E(x) llo((1 x)) when minimized. The cost is defined by:
min w)( = log det RrG log det Rr log det RG, r Σ log(1 ox). 8w
(3)
Normalization trick. After training, normalizations are needed to obtain eigenfunctions. The first step to ensure orthonormality: fo = R ff g_ = R g.. The second step is a singular value decomposition: PFG f (X)g
US V. where S diag (AI, Axx). The third step is to normalize functions such that they are invariant to the linear operator: fo = UT 9 = VTg.. Functions for 9., are the top eigenfunctions of the density ratio, and \1, '2, ... are the top eigenvalues. An approximation of the density ratio given by أ fo T
S g p.
Newly proposed: matrix trace cost. This paper explores alternative convex functions, specifically the simplest case {(x) x, The cost, in the form of a matrix trace, is described below.
Lemma 2. Denote P = PFG. Given neural nets fo and g minimizing the matrix trace
min (f (raceeR;)) 'R-'PR- PT),
8,0
yields (0, - eigenvalues of the density ratio when minimized. We name this algorithm FMCA-T.
= where (9,ww) is the optimal cost, reaching the sum of the top K
Proof. Applying the Schur complement to rz, we obtain rz = log det R PR R
Denoting eigenvalues of matrix \1(.), .x((), the cost becomes rr. = Exx log(1 (xM))
where M R PR pT R Optimizing the sum of eigenvalues instead, use Trace (M) and, based on the trace property Trace (AB) = Trace( BA), derive the trace cost for learning multivariate statistical dependence as rt = Trace( R PR PT).
FMCA-T is more computationally efficient uses only matrix operations of dimension K. Directly optimizing the sum of the cigenvalues is also more stable than optimizing their logarithm.

---

2.3 Channel-level and temporal-level dependencies
Motivations. For EEG X and EMG YIT,, FMCA-T applies two networks, Se and 9w- to minimize the matrix trace cost. Dependence is measured two levels: random-process level, measured by eigenvalues for the overall dataset dependence, and frial level. measured by the density ratio-the higher the ratio, the greater the contribution of this pair of realizations to overall dependence. In cortico-muscular muscular analysis, is vital to understand how individual channels and time steps contribute connectivity, especially EEG signals, as they represent the temporal and spatial dynamics of brain. Hence, we propose localized density ratios to measure temporal-level dependence and channel-level dependence. The core idea is computing density ratios between channel-level and temporal-level features against the global trial-level features.
Channel-level features. We design a specialized network topology to generate features for individual channels and time intervals, ensuring that the internal layers of this network quantify channel-level and temporal-level features, similar to [22-24].
Given 11T = (I), (C)] for channels c= 1, 1..... C, we define temporal network FIN R R that maps single-channel signals to K-dimensional feature space, and a channel network FON R - IR K that maps concatenated channel features to global features:
Z. FrN (X (iT())), c=1,,.,CC ZF = ([ Zz, :. Zcl"),
(5)
where 11 2,, Zc are channel-level features, and ZF is global trial-level features.
Channel-level dependence PC,F (c). The density ratio of Z1. Zz. Zc relative to Zr measures channel-level dependence. Post- bost-training and with fixed parameters, we compute the ACF of the channel features Rc = E c=1 Z Z.Ze TI the ACF of the global features Rs Z Z [... and
the CCF between them Pc,F E [.
Next, the features normalized in the Sec. 2.2: Ze and Zr normalized to Ze R Ze and Zr R Zr for orthonormality. The SVD of R Pc
US V is computed. The outputs are further normalized to Z. U.Z and ZF VTZ to guarantee invariance in the linear operator. Finally, the density ratio can be constructed as Pc. Z. S2
This ratio (C) is function of channel and trial X, implying the dependence between channel and global features. The greater the value, the stronger the activation of the channels, showing which channels contribute the most the cortico-muscular connectivity.
Temporal-level features and dependence. To measure time-domain dependence, we compute density ratios between the internal features of the temporal network FIN and the global features, in two steps: first, computing density ratios between adjacent network layers; second, aggregating these ratios consider all layers.
Step 1: Construct density ratios (T- (TI,,) T2) between adjacent layers. Fix a channel c and feature Z... Consider simple temporal network with S convolution layers with nonlinaer tivation functions: FTN. F(2) .... F(3)
with kernel sizes A1, A2. As, and their outputs c, Ze.2. .... Ze.s. Suppose the time dimensions of these layers are Tı, T2. ... Ts. Fix any layer 8. The r-th element of Ze,., denoted as Ze,((T) is obtained by applying a nonlincar operation to a segment of the previous layer's output:
Ze,((T) = F((-1)
(Zo :T + A.-1)).
(6)
Define the ACF of layer 1: Re,s-1
1 E Σ Ze,o-((T) ZT (T)
To-
Define the ACF of layer s: Re = 7 E Zo.s (T) TT (T)
Define the CCF between them: P
+
+E
Σ, Ad, Z Ze OZ)
-
5

---

Normalization with Re.s-1. Re.e, and P cc.--I.e yields density ratios Ps-1.s,e T2), which quantify the dependence between time Ti and 72 across two layers s and for given trial and channel c. A higher value indicates a stronger dependence between adjacent network layers.
Step 2: Aggregate layer-wise ratios for localized responses (s,c(T). While Ps-1,s,e(71, 72) quanti- fies dependence between two layers, we aggregate these ratios to account for all network layers.
Again, fix the element ce (T) in layer s. We focus its mapping to the next layer + 1. Based on Eq. (6), elements in layer s + 1 that are mapped from (T) by F( are within a window of size Ay (the kernel size). To ensure this window stays within the feature vector's boundary, coordinates I3- = max (0, A, + 1), min (i, T,+1 - 1)]. Feature elements in layer + that are mapped from Z c. (T) fall within these coordinates.
define
We then create a series functions (T)(T) with T E [1, Tal for each layer as the aggregations of the ratios. Starting from (T, = (c) (channel-level density ratio), compute recursively
(T)) = Σ Bottic((t)Pssssiiitii TI TH [1,Tg] ToEL.
That is, starting from the top layer of the network, aggregate the density ratios within the window Ig layer-by-layer until we generate localized measurement for each element of the function (T, at layer s, channel c, and time [1, T. considering all neural network layers.
The final localized responses of the density ratio, (T, obtained in top-down manner, functions of the EEG trial 1 time step т, and channel c, providing both temporal and channel- level resolution. The same analysis applies to EMG signals, differing only in the number of channels.
Experiments
Our experiments have three key findings: (1) Dependence measured by FMCA-T is stable against nonstationary noises and delays in simulated dataset; (2) Learning from unlabeled EEG-EMG pairs extracts movement and subject information from EEG's eigenfunctions; (3) EEG's spatio-tempol dependencies are consistent with ground truth brain activations in simulated dataset and match theoretical evidence in experimental dataset.
3.1 Datasets and baselines
Dataset Wav, We construct a simulated dataset where each data pair {XIrT, Y [iTT } has clean sinusoidal signal X, A (in( and noisy counterpart Y, superimposed with various types noise: stationary white noise E ~ NO0 o²), nonstationary Gaussian noise of x Xe), and nonstationary pink noise S(f) x 1 Random delays are added by padding the start and end of the signal with noise such that Y,- = C is white for < t S Ti and T 72 t S T, and Y, = X,-r is sinusoidal for Ti < t < T - where Ti T2 = is fixed. Since these noises do not change the underlying sinusoid, the signal pairs are statistically independent when conditioned on the sinusoid. Thus, we expect the dependence measure to be unaffected by the noise level.
Dataset 2: EEG-EMG-Fusion. We use public dataset [25] (approved by the Institutional Review Board at Korea University, 1040548- 1--0548-KU--BB -U-IRB--7-188-A4 channel EMG recordings from 25 subjects. The subjects perform three main movements: arm- reaching. hand-grasping, and wrist-twisting. Each main movement contains sub-movements: arm- reaching along six directions, hand-grasping three objects, and wrist-twisting with two motions, and thus 11 movements in total. Subjects perform one sub-movement per trial, and 50 trials are collected per sub-movement. The same recordings repeated for three sessions at one-week intervals. Both EEG and EMG are recorded at 2,500 Hz and downsampled to 1,000 Hz. The dataset cleaned by removing eye-blinking artifacts and baseline wandering, and segmented into 4-second intervals, creating 41,250 paired samples of complete movement cycles.
with paired 60-channel EEG and
Dataset 3: Simulated EEG-EMG Dataset. We simulate 128-channel EEG signals and 7-channel EMG signals for left/right motor and sensory activations from 20 subjects using EEGSourceSim

---

[26]. Motor sources are used to simulate the corresponding EMG signals. Both EEG and EMG are sampled at 1,000 Hz. White Gaussian noise at dB is added to the EEG and EMG signals. FMCA-T is trained on 16 subjects and tested on 4 subjects to compare FMCA-T's spatial-level dependence representation with the ground truth activations, as shown in Fig. 6.
Classification tasks. We conduct three classification experiments: 3-class (three main movements), I1-class (11 sub-movements), and Subj (25 subjects). We also compare inter-subject and cross- subject classifications. Cross-subject means the test set contains only unseen subjects. Inter-subject uses an 80-20 split of trials from all subjects for training and testing, while cross-subject uses 20 subjects for training and 5 for testing with five-fold cross-validation.
Statistical dependence baselines. We compare our proposed dependence measure with established baselines: (1) KICA [17] and HSIC [18], which solve the generalized eigenvalue problem of two kernel Gram matrices, using Σ, A; for HSIC and Σ; log(1 >,) for KICA; (2) MINE [14), which optimizes the Donsker-Varzn Varadhan representation with a three-layer MIP; (3) CC: Pearson correlation coefficient averaged over time; (4) MIR (KNN estimator [27]), which optimizes entropies using k-nearest neighbor distances, and then computes mutual information; (5) Our method uses density ratios for trial-level dependence and eigenvalue aggregations Te for candom-process-ll dependence.
For the EEG-EMG dependence study, we compare with CMC, the correlation coefficient between EEG and EMG spectra of windowed data on the alpha band of channel C4 [25]. We extend CMC by replacing linear correlation with nonlinear measures, computing CMC-KICA and CMC-MIR.
EEG feature projector baselines. After training fa and g., networks for dependence estimation, with parameters fixed, we train a three-layer MLP on EEG's eigenfunctions (X)) for classification. This is compared with baseline EEG classifiers trained from scratch: (1) Supervised: Vanilla Classifier, with the same topology as ours but using a standard log- ligeliieldoood cost; EEG-Nel [22), a specialized network for EEG-based BCI; FFGG-Coffrrmr [28], compact convolutional transformer for EEG decoding and visualization; Deep4 [29], a deep ConvNet for classification using raw FEG; and CSP-RI.DA [30], using common spatial pattern (CSP) for feature extraction and regularized linear discriminant analysis (RLDA), adapted for multi-class classification with majority voting. (2) Self- Supervised: contrastive costs using 1-second windows from the same signal as positive pairs and from different signals as negative pairs, including Barlow Twins [31], SimCLR [32], and VicReg [33]. Experimental and implementation details can refer to the App. C.
3.2 Main results
Robustness of FMCA-T. Fig. 3 shows the robustness of our proposed measure in the Sin dataset, when there are increasing levels of nonstationary noise and delays. Since EEG and EMG signals are often damaged and distorted by environmental noise and the functional coupling occurs with time delays, an effective measure should maintain its robustness against these factors.
(a) Ganesien Noise
(b) Nona ba ticnary Geuazian
(c) Pink Boise
(d) Randon Delay
1.50
PNCA Isssuce
Xeasure
M5asuze
Messure
313. Measuue
Figure 3: Density ratios from FMCA-T are robust to various noise types: (a) stationary white Gaussian noise, (b) nonstationary Gaussian noise, (c) nonstationary pink noise, and (d) random delays. FMCA-T proves the most robust estimations across all noise types and outperforms all linear and nonlinear baselines. Note that as delays increase, estimations using CC produce negative values given the opposite phase between the paired sinusoids.
In Fig. 3, FMCA-T is first trained on noisy data pairs with all four noise types and magnitudes (Sec. 3. Using the trained models, we measure the dependence between a clean sinusoid and its noisy counterpart. A noise level of 1.0 means the noise magnitude matches the sinusoid. The delay level determines the extent to which the clean sinusoid is shifted from its original position. A delay
7

---

level of 1.0 shifts the sinusoid to have no intersection with the original one. FMCA-T consistently shows invariance noise and delays, as the dependence is determined by their frequency but not by the noise and phase shift. MINE fails to converge and produce stable results for this dataset.
Applying FMCA-T to EEG-EMG-Fusion. We confirm primary hypothesis that the projection space defined by EEG eigenfunctions, derived from modeling the statistical dependence between EEG-EMG recordings, captures essential contextual factors like movements and subjects without requiring labels. We visualize the learned eigenfunctions and density ratios in Fig. 4.
Fig. 4(a), Fig. 4(b): eigenfunctions for EEG's eigenfunctions effectively capture relevant con- textual information. After training eigenfunctions using the entire dataset, extract subset of eigenfunctions that belong to specific subject or a movement and apply t-SNE to visualize them.
Fig. 4(a) visualizes the eigenfunctions of all trials for one subject (SUB1). Each trial is color-coded by the type of movement (MOV1~MOV3). Notably, the cigenfunctions form nine clusters, which are verified to correspond to the three movements recorded over three sessions. This demonstrates that the eigenfunctions contain motion-related information. The consistent clustering patterns across all 25 subjects detailed in the App. B.
Fig. 4(b) visualizes the eigenfunctions from a single type of movement (reaching, labeled as MOV1) across ten different subjects (SUB1~SUB10), Each color represents a subject. Distinct clustering patterns are observed, showing that the eigenfunctions also contain subject information which could be useful for participant identification.
Fig. 4(c), 4(d): density ratio p(X, Y). Based on the t-SNE plot for SUB1 trials, we plot the estimated density ratio values between each EEG-EMG pair in Fig. 4(c). In Fig. 4(d), we extract the mean of density ratios for trials in each cluster (C1~~99,, rank them from smallest to largest, and plot them alongside the standard deviation. These figures show that the density ratios remain consistent within each cluster (a movement during a session) while vary across different clusters. We find the highest dependence in reaching, followed by grasping, and the lowest in twisting. Our results are consistent with existing literature that links cortico-musculat connectivity with movement types [34, 35J. The results are consistent across all subjects, detailed in the App. B.
Fig. (h) presents the results of MINE and CMC measures for SUB1 trials. Only MINE produces comparable measure that shows difference across clusters, but with higher variance and instability.
[+) TSNEE SUBL, MDVI -MCV3
(b) FSHE, MOV1, SUB10
(c) FMCA Denaity
Gruup
Resch
140
120
100
Twint
(d) EMCA Neasure
(6) MINE Heasure
CMC Measuce
[o) CMC-MIX Measure
(h) CMC-XICA Heasure
GGGO
G.
Figure 4: Visualizing eigenfunctions and density ratios in EEG-EMG fusion with FMCA-T: (a) t-SNE of EEG's cigenfunctions for a single subject (SUB1) show nine clusters specific to three movements (MOV1~XOV3) across three sessions. (b) t-SNE of EEG's cigenfunctions for reaching movement (MOV1) across 10 subjects (SUB1~SUB10) shows clusters specific to subjects, where each color a subject. (c) Density ratios and (d) their mean and std each cluster (C1~ C9) demonstrate intra-cluster consistency and inter-cluster separability. (e-h) Comparison of baseline measures, where only MINE is comparable but with higher variance and instability.
EEG's eigenfunctions as optimal feature projector. Table validates the claims in Fig. 4 with classification accuracy comparisons. We extract EEG eigenfunctions from the training set after the networks are trained on EEG-EMG pairs. The eigenfunctions are used to train a three-layer MI. for classification. This classifier predicts the class of any EEG test samples using its cigenfunctions
8

---

(output of fo). without requiring EMG samples. As detailed in the Sec. 3.1, scores for 3-class, I I-class, and subj are presented both inter-subject and cross-subject settings.
In inter-subject 3-class classification (80-20 split across trials from 25 subjects), FMCA-T exceeds the supervised baseline (EEGNet) by 7.2%, the classical EEG decoding method (CSP-RLDA) by 1.0%, and self-supervieed methods by over 9.5%. Notably, CSP-RLDA is trained and tested on each individual subject, with the accuracy averaged across 25 subjects, thereby representing upper bound for classical methods. All other methods are trained and tested on the combined subject data. For the 11-class classification task, FMCA-T surpasses all baselines with a classification accuracy of 0.32, significantly higher than the chance level of 0.09. Since CSP-RIDA uses binary classification with majority voting, it is computationally infeasible for 11 subemovements classification.
In the more challenging cross-subject classification (trained on 20 subjects, tested on 5), FMCA-T with trace loss outperforms all baselines by over 10%, achieving an accuracy of 0.54 in the 3-class task. The highest scores from 10,000 iterations are recorded, and experiments are repeated with five-fold validation. The superior performance of FMCA-T in the cross-subject setting suggests that learning EEG-EMG dependence is robust against distribution shifts and nonstationary noise, which is consistent with the observation that sele-super methods outperform supervised ones.
Comparing FMCA-LD (log dct) and FMCA-T (matrix trace), we find that trace cost has greater stability and reduced variance. The sum of eigenvalues, especially during prolonged training, is more stable. While both costs show similar performance at the initial training stages, FMCA-T has notably reduced variance during the convergence stage of training.
Methods
(a) Inter-Subject Ace.
(b) Cross-Subject Ace.
Table 1:
Comparison
3-Class
11-Class
Subj
3-Class
11-Class
of classification accura- cies: supervised, self- supervised, and our EEG- EMG dependence leam- ing. FMCA-T's cigen- functions, trained with trace cost without labels, are optimal feature pro- jectors for EEG. EMG is not required for testing, but only used for train- ing.
Vanilla
0.907=0.02 020
0,220=0.015
0,980=0.010
0.427+0.021
0,110=0.00 005
EBGNet
0.90420.015
0,24620.028
0.98810.077
0.40510.019
0.09520.021
EEG-Conforeer
0.949=0.001
0.26820.001
0,97640.002 002
0.415+0.002
0.105=0.001
Deepd
0.90120.001
0.27410.001
0.94110.001
0.42910.001
0.14020.000
CSP RI
0.985=0.019
0.408+0.018
Siperiii
Barlow Twins
0.893=0.018
0,269=0.012
0.987 008
0.437
0.1 0..00000004
SimCLR
0.890=0.019
0.25720.013
0.979=0.011
0.441 020
0.11720.006
VicReg
0. 0.899=0.016
0,274+0.014
0.980=0.009
0.449*0.010
0.1 0.11000.005
REG-EMG Deversience
FMCA-LD
0.985=0.003
0,257+0.011
0.989-00007
0.509+0.014
0.115:0.003
FMCA-T
0.994=0.002
0.320=000009
0,998=0.004
0.540=0.012
0.121=0.002
Spatio-temporal dependencies real data. We visualize the local density ratio responses of cluster SUB3-C1 (reaching movement) in both spatial (pc)F(c) and temporal domains (T), in Fig. 5. The channel-level dependence is averaged across all trials and displayed in Fig. 5()). We also randomly select nine trials from the same cluster and visualize them in Fig. 5(b). The temporal-level dependence for the first trial T1 in SUB3-C1 shown in Fig. 5(c). Consistent activations are observed in other subjects, details in the App. B.
As illustrated in Fig. the localized density ratio remains stable throughout the 4-second movement, connectivity during stable states [36]. We also find
corroborating the consistency of that in channel-level dependence, the density ratio activates the fronto-central (FC) areas. The sensorimotor area is crucial for movement control, with EEG data from these regions often used to
decode motor intentions. However, motor control also relies on cognitive processes [37]. especially during movement planning, complex tasks, and collaborations [38]. Thus, the region of Brodmann area 6, well acknowledged to playing role in movement planning may contribute differently during various movement tasks [39]. Our findings show that the features extracted from these fronto-central areas play an important role in classification.
Spatio-temporal dependencies simulated data. We implement FMCA-T on paired EEG-EMG samples from the simulated dataset, using 16 subjects for training and 4 subjects for testing. We compare FMCA-T's spatial-level dependence maps for these testing subjects with their ground truth brain activations computed by the motor ROI and forward matrices, shown in Fig. 6. FMCA-T's spatial-level dependencies are highly similar to the ground truth activations, indicating that the learned density ratios effectively captured the real brain activations.

---

[4] Channel Dependance 50B3-C1-AV0
(b) Chaccell Dependance 2053-c1, 11-19
Nemporael Dependence 5083-C1-12
AND
Figure 5: Localized density ratio for real data, Topographies for SUB3, C1 (reaching session) normalized the range [0. 1]. (a) displays the averaged spatial distribution C1 across all trials; (b) presents nine random trials T1-T9 from this session. Dark red indicates prominent activations around channel FC1. (c) shows stable temporal-leee dependence over the 4-second movement.
IMCA-1 Dependance (1)
Ground Truth INCA-T Ground Truth Activations Dependance Activations (L)
PMCA-T Dependence (L)
Ground Truth PNCA-F Activations Dependence Activations
Ground Truth
(R)
(R)
(L)
(R)
(R)
Subject
Subject 19
17
Subject 18
Subject
20
Figure 6: Localized density ratio for simulated data, showing that FMCA-T's distribution is consistent with ground truth brain activations when evaluated on four test subjects. Topographies are normalized to the range [0, 1]. L indicates left brain activities and R indicates right brain activities.
4 Conclusion
This paper introduces a novel approach for estimating cortico-muscull dependence through the orthonormal decomposition of the density ratio. By treating the density ratio as a positive definite function and learning its projection space from EEG and EMG, we unveil the relationship between statistical dependence, contextual factors impacting connectivity, and the apatio-temporal information shared between the brain and muscles. While our method shows promising results, challenges remain. For example, performance drops in cross-subject classification, likely due to the limited dataset of 25 participants. Future work will focus on applying our framework larger datasets and incorporating additional bio-signal modalities to model a broader common space in neural data.
References
[1] S. N. Baker, M. Kilner, E. M. Pinches, and R. N. Lemon. The role of synchrony and oscillations in the motor output. Experimental Brain Research, 128:109-117, 1999.
[2] Stuart N. Baker. Oscillatory interactions between sensorimotor cortex and the periphery. Current Opinion in Neurobiology. 17:649 655, 2007.
[3] Jan Mathijs Schoffelen, Robert Oostenveld, and Pascal Fries. Neuronal coherence mechanism of effective corticospinal interaction. Science, 308: 11--11,,,, 2005.
[4] Cristina Brambilla, leaaa Pirovano, Robert Mihai Mira, Giovanna Rizzo, Alessandro Scano, and Alfonso Mastropictro. Combined use of EMG and EEG techniques for neuromotor assessment in rchabilitative applications: systematic review. Sensors, 21(21). 2021.
[5] P. Grosse, M.J. Cassidy, and P. Brown. TEG-EMG. MEG-ENE and DMG-EMMR frequency analysis: physiological principles and clinical applications. Clinical Neurophysiology 113(10):1523- 1531, 2002.
10

---

[6] Roisin McMackin, Peter Bede, Caroline Ingre, Andrea Malaspina, and Orla Hardiman. Biomarkers in amyotrophic lateral sclerosis: current status and future prospects. Nature Reviews Neurology, 19(1))754 768, 2023.
[7] Zhixian Gao, Shiyang Lv, Xiangying Ran, Yuxi Wang, Mengsheng Xia, Junming Wang, Mengyue Qiu, Yinping Wei, Zhenpeng Shao, Zongya Zhao, al. Influencing factors of corticomuscular coherence in stroke patients. Frontiers in Human Neuroscience, 18: 18:1354332, 2024.
[8] Nahid Zokaci, Andrew J Quinn, Michele T Hu, Masud Husain, Freek van Ede, and Anna Christina Nobre. Reduced cortico-muscular beta coupling Parkinson's 's disease predicts motor impairment. Brain Communications, 2021.
[9] Anirban Chowdhury, Haider Raza, Yogesh Kumar Meena, Ashish Dutta, and Girijesh Prasad. An EEG- EMG correlation-aase brain- compuerr interface for hand orthosis supported neuro- nuuoo-ccabilitation. Journal of Neuroscience Methods, 312:1-11, 2019.
[10] Tatsuya Mima and Mark Hallett. Corticomussularr coherence: review. Journal of Clinical Neurophysiol- ogy, 16(6):501, 1999.
[11] Gent Pfurtscheller and FH Lopes Da Silva. Event-related EEG/MEG synchronization and desynchronization: basic principles. Clinical Neurophysiology, 1 110(1 (01111))142-1857 1999.
[12] Ole Jensen and Laura L Colgin. Cross-frequenn coupling between neuronal oscillations. Trends in Cognitive Sciences, I 1(7):677-269 2007.
[13] Yuhang Xu, Verity M McClelland, Zoran Cvetkovié, and Kerry R Mills. Corticomuscular coherence with time lag with application to delay estimation. IEEE Transactions on Biomedical Engineering, 64(3):588 0000 2016.
[14] Mohamed Ishmael Belghazi, Aristide Baratin, Sai Rajeswar, Sherjil Ozair, Yoshua Bengio, Aaron Courville, and R. Devon Hjelm. MINE: Mutual information neural estimation. In Proceedings of the 35th International Conference on Machine Learning (ICML 2018), volume 80, pages 531-540, 2018.
[15] XuanLong Nguyen, Martin J Wainwright, and Michael Jordan. Estimating divergence functionals and the likelihood ratio by convex risk minimization. IEEE Transactions on Information Theory, 56( 155847-6661 2010.
[16] XuanI Nguyen, Martin Wainwright, and Michael Jordan. On surrogate loss functions and f- divergences. The Annals of Statistics, 37(2)876-904, 2009.
[17] Francis Bach and Michael 1 Jordan. Kernel independent component analysis. Journal of Machine Learning Research, 3(Jul): 2002.
[18] Arthur Gretton, Olivier Bousquet, Alex Smola, and Bernhard Schölkopf. Mcasuring statistical dependence with Hillbert-Schmide norms. International Conference on Algorithmic Learning Theory, pages 63-77. Springer, 2005.
[19] Bo Hu and Jose C. Principe. The normalized cross density functional: A framework to quantify statistical dependence for random processes. arXiv preprint arXiv:22/2.04631, 2024.
[20] Shao-Lun Huang. Gregory W Wornell, and Lizhong Zheng. Gaussian universal features, canonical correlations, and common information. In 2018 IEEE Information Theory Workshop (ITW), pages 1-5. JEEE, 2018
[21] Shao-Lun Huang. Anuran Makur, Gregory W Wornell, and Lizhong Zheng. On universal features for high-dimensional learning and inference. arXiv preprint arXiv:1911.09105. 2019.
[22] Vernon Lawhemn Amelia Solon, Nicholas Waytowich, Stephen Gordon, Chou Hung, Brent Lance. Fegnet: compact convolutional neural network for BEG-hased brain-computer interfaces. Journal of Neural Engincering, 15(5)0056013, 2018
[23] Daoze Zhang. Zhizhang Yuan, Yang Yang, Junru Chen, Jingjing Wang, and Yafeng Brant: Foundation model for intracranial neural signal. Advances in Neural Information Processing Systems, 36, 2024.
[24] Zhizhang Yuan, Daoze Zhang, Junru Chen, Geifei Gu, and Yang Yang. Brant-2: Foundation model for brain signals. arXiv preprint arXiv:2402. 10251, 2024.
[25] Ji-Hoon Jeong, Jeong- Hyun Cho, Kyung- Hwan Shim, Byoung-Hee Kwon, Bycong- Lec, Do-Yeun Lee, Dae-Ilyeok Lee, and Seong-Whin Whan Lee. Multimodal signal dataset for 11 intuitive movement tasks from single upper extremity during multiple recording sessions. GigaScience, 910)ggaa0098 2020.
11

---

[26] Elham Barzegaran, Sebastian Bosse, Peter J Kohler, and Anthony M Norcia, EEGSourceSim: A framework for realistic simulation of EEG scalp data using mri-based forward models and biologically plausible signals and noise. Journal of Neuroxcience Methods, 328:108377, 2019.
[27] Alexander Kraskov, Harald Stöghauer, and Peter Grassherger. Estimating mutual information. Physical Review E. 6966)80 2004.
[28] Yonghao Song, Qingqing Zheng, Bingchuan Liu, and Xiaorong Gao. Eeg conformer: Convolutional transformer for EEG decoding and visualization. IEEE Transactions Neural Systems and Rehabilitation Engineering, 31:710-719, 2022.
[29] Robin Tibor Schirrmeister, Jost Tobias Springenberg, Lukas Dominique Josef Fiederer, Martin Glasstetter, Katharina Eggensperger, Michael Tangermann, Frank Hutter, Wolfram Burgard, and Tonio Ball. Deep learning with convolutional neural networks for EEG decoding and visualization. Human Brain Mapping, 38(1 (::11)9931-5420, 2017.
[30] Kai Keng Ang, Zheng Yang Chin, Chuanchu Wang. Cuntai Guan, and Haihong Zhang. Filter bank common spatial pattern algorithm on BCI competition IV datasets 2a and 2b. Frontiers Neuroscience, 6:2 1:21,02, 2012.
[31] Jure Zhontar, Li Jing, Ishan Misra, Yann LeCun, and Stéphane Deny. Barlow twins: Self-supervised learning via redundancy reduction. In International Conference on Machine Learning, pages 12310-12320. PMLR, 2021.
[32] Ting Chen, Simon Komblith, Mohammad Norouzi, and Geoffrey Hinton. A simple framework for contrastive learning of visual representations. In International Conference on Machine Learning, pages 1597-1607. PMLR, 2020.
[33] Adrien Bardes, Jean Ponce, and Yann LeCun. Vicreg: Vuriance-imarrineecccvvaaanererr regularzzztoo for sule-supe learning. arXiv preprint arXiv:2105.04400 2021.
[34] Cristian D Mneerero-Mended and Andres F Ruiz-Olaye Coherence-baseee connectivity analysis of EEG and EMG signals during reach-to-grasp movement involving two weights. Brain-Computer Interfaces. 9(3):140-154 154, 2022.
[35] Fei Ye, JinSuo Ding, Kai Chen, and Xugang Xi. Investigation of corticomuseular functional coupling during hand movements using vine copula. Brain Sciences, 12(6):754, 2022.
[36] Jinbiao Liu, Yixuan Sheng, and Honghai Liu. Corticomuscular coherence and its applications: a review. Frontiers in Human Neuroscience, 13: 100, 2019.
[37] Jason P Gallivan, Craig S Chapman, Daniel M Wolpert, and Randall Flanagan. Decision Decision-making in sensorimotor control. Nature Reviews Newroscience, 19(9):11-5344 2018.
[38] Tianyu Jia, Jingyao Sun, Ciarán McGeady, Linhong Ji, and Chong Li. Enhancing comin-crmputer interface performance by incorporating brain-to-brand coupling. Cyborg and Bionic Systems, 5:0116, 2024.
[39] Richard P Dum and Peter L Strick. Motor areas in the frontal lobe: the anatomical substrate for the central control of movement. In Motor Cortex in Voluntary Movements, pages 3 48. CRC Press, 2004.
[40] Jingyao Sun, Tianyu Jia, Zhibin Li, Chong Li, and Linhong Ji. Enhancement of EEG-EMG coupling detec- tion using corticomuscular coherence with 20(3):036001, 2023.
optimization. Journal of Neural Engineering,
[41] Igor Susmelj. Mathiis Heller, Philipp Wirth, Jeremy Prescott, and Malte Ebner et al. Lightly. GitHau. Note: https://githubbco htwp://hwwwcooo/iggggtllaiiiiiii 2020.
12

---

A Additional Implementation Details of FMCA-T
Adaptive estimators. Since autocorrelation functions (ACF) and crosscorrelation functions (CCF) are expectations of inner products, we prefer to compute these expectations and the cost gradient using an adaptive filter rather than a batch of data. For any parameter 8, the partial derivative of the cost (T 6. [w) Trace( R PR
PT) has the following form:
urt
RRF R
OP
Tracc R
PRG
PT) Trace R
PT) Rg
00
00
00
(8)
RRG
0P
-Trace R
PRG
R PT) + Trace( R P 00
00
Observe that terms needed for the gradient have two classes: terms RF. Rg, P. as well as their inverse, and terms of the derivatives 0Rr
iterations to estimate the Rr. Rg. P, and then use the estimated values of Rr, Rg. P as the matrices and 8RG äp we use the derivatives of the batch. We 19 88 36 choose smoothing coefficient B = 0.9 across all experiments.
8Ra OP We use smoothing window
36
00
their inverse in the cost. For derivative terms 0Rs
Important parameters. When computing the matrix inverse of ACFs (Rr-1 and R small diagonal matrix scaled by a regularization parameter, eI, added. This constant is important for training stability. We choose
10-5 across all experiments. The output dimension of the network chosen as K = 128, which is also the number of eigenfunctions.
Pseudo-code of the algorithm. We also explain our algorithm with pseudo-code. Algorithm I shows the pseudo-code for FMCA-T with adaptive estimators and computing eigenfunctions, Algorithm shows how to generate channel-level and temporal-lewe level dependencies.
Algorithm 1: FMCA-T w/ Adaptive Estimators.
Alg. Adaptive matrix estimators
Input: Batch estimation A; Tracking matrix A; Iteration i;
Step A + /- À (1-3) A
3
Step A A/(1- - )²)
return Tracking mairix A; Smoothed estimation A
Alg. FMCA-T: Optimize matrix trace cost
Initialize: Neural networks (S 9w 5
Initialize: Tracking matrices Rp, RG. for = 1, do
Sample batch of signals {X,,. Y, Jn }be Compute the ACF and CCF {Ry, Rg, P] with this batch;
Apply adaptive estimators, obtain {Rr. RG. P), and update (Rr. RG, P).
Estimate gradients with smoothed estimator, update networks.
end return 0,w
Alg. Retrieve eigenfunctions and density ratios
Input: Trained networks f a and gw- Dataset {X,,, Y, YnJnE1 Step 1. Obtain outputs of networks. (X,)) (Y,)).
13
Step Re-estimate RF, Rg. P using these outputs.
Step 3. Normalization for orthonormality: Jo R fo. Tw
R Iw
Step 4. Compute SVD: P Ef ())g T(X) = USi V
16
Step Normalization for invariance: fo U fa. 9 V g
18
Step Construct the density ratio: ò -= fo ~ p.
return Figenfinnciions f T eigenvalues S diag .., Ax). density ratio p
13

---

Algorithm 2: Generate channel-level and temporal-level dependencies.
Alg. Generate temporal and channel features:
Input: Trained network f for EEG, consisting of temporal network FrN and channel network FcN: An arbitrary EEG trial X
2 Step 1. Apply FTN to each channel: Ze FFT((X()),, c=
3 Step Apply FON to all channel features: Zp = FON FCN (IZ1. ZclT)
Step Extract internal features of FCN: Ze,1, ... Ze,s from convolution blocks
return Temporal features Ze channel features Ze, global features Zr
Alg. Generate channel-level dependencies
Input: Temporal, channel, global features from the dataset foreach channel do
Compute ACF and CCF for global and channel features, ZF and Ze Use ACF and CCF to compute density ratio PC,F(c) with normalization
end
Alg. Generate temporal-level dependencies:
10
foreach channel e and layer s
Compute ACF and CCF for temporal features between two layers, Ze Z Ze,s+1
Use ACF and CCF compute density ratio ( 72) between two layers. end
Initialize (T, (c): Initialize each (T) have length T..
for 5 - 1, 1,.., do
15
For every ZT( find elements in layer + I that are mapped from it (set I,)
16
Aggregate density ratios (,,c(T)) = E,eL, end
T2),
return Channel dependence PC.F and spatio-temporal dependence (T. that can be applied any arbitrary EEG trial X.
14

---

B Additional Experiments
Visualization of EEG's eigenfunctions. Building on Fig. 4 from the main paper, we further present results for the visualization of EEG's eigenfunctions on all 25 participants. In Fig. 7, each trial is color-coded by the estimated density ratios for all participants (SUB1 SUB25). In Fig. each trial is color-coded by the label of their movements.
be observed that for each participant, the t-SNE projections form individual clusters, each corresponding to one of the three main movements during a session. The density ratios are highly similar within each cluster (intra-cluster) while varying across clusters (inter-cluster). This indicates that using density ratio as a dependence measure effectively captures each movement individually.
FMCA-1 -T Density Ratio, SUB1~SUB25
160
725
142
140
180
160
-80
192
140
240
120
100
800
200
340
1000
120
100
130
:00
800
224
250
300
100
140
130
200
200
1.0
:00
80
100
253
325
200
200
115
273
120
106
125
100
Figure 7: t-SNE visualization of eigenfunctions from trials of 25 participants, color-coded by their corresponding density ratios. Each subplot represents one participant. The density ratios for the three movements (reaching, grasping, and twisting) show consistent values: reaching is the lowest, grasping is in the middle, and twisting is the highest. For cluster labels, refer to Fig. 8.
15

---

TSNE, SUB1~SUB22, MOV1~MMV33
Figure 8: t-SNE visualization of all trials from 25 participants, color-coded by the type of movement: reaching (green), grasping (blue), and twisting (red). Each participant has distinct clusters, with each cluster corresponding to one of the three types of movements observed during the sessions.
From the figure, we observe that the reaching movement consistently has the lowest dependence, twisting highest, grasping falling in between. This pattern matches the differences and com- plexity of different movement tasks [34, 35]. This variability reflects the differences in connectivity between the brain and target muscles, which is the basis of movement recognition and measurement by using CMC [40].
Mean and variance of density ratios in clusters. Same as in Fig. 4 in the main paper, we extract the nine clusters shown in the t-SNE projections (Fig. 8) and compute the mean and standard deviation of the density ratios for each cluster (C1~C9). This process is repeated for all 25 participants (SUB1~SUB25), with results shown Fig. 9. The results further demonstrate the effectiveness of using density ratio dependence measure, as the density ratios are similar within each cluster but vary across different clusters (movements and sessions).
Participant identification. Similar to Fig. the main paper, where we visualize the projected EEG cigenfunctions from 10 subjects during the reaching movement, extend this analysis to include all
16

---

three movements: reaching, grasping, and twisting. The projections shown in Fig. illustrate that participant information is consistently contained in the projection space of the eigenfunctions across all movements, not just limited to reaching.
FMCA Density Ratio Variability,
SUB1~SUB25
101
500
13
130
02 C1 6 0
G 0 01 6 66
12
61 C 6T C8
62
180
32 28 12
190
130
00606000
G00000000
30
191
160
140
120
191
C2
Figure 9: Similar to Fig. 4, we extract the density ratios of cach cluster for the subject, corresponding to one movement and session, rank their means, and plot them along with their standard deviations. can be that for each subject, the density ratio value each cluster different, indicating that our measures provide unique measurement for each movement the subject. This figure extended from Fig. 4 shows that this observation applies to all participants.
TSNE, REACH, GRASP, TWIST,
SUB1~SUBE
8
8
6
6
4
4
2
2
2
0
Figure 10: t-SNE projection of trials from ten participants. Subplots from left to right show the projected EEG eigenfunctions from reaching, grasping, and twisting movements, respectively. In all movements, participant- specific information consistently captured in the clusters of the EEG's eigenfunctions.
17

---

Temporal-level dependence. Extending the analysis in Fig. 5 of the main paper, where the temporal- level dependence was shown for single subject, we randomly select another seven trials from subject SUB3's reaching movement (C1) and visualize the temporal-level dependence in Fig. 11. We also plot the average temporal dependence across all trials in Fig. 11h. We find consistent activations of fronto-central (FC) channel during the s-cono movement in each trial. When we compare each trial's result with the averaged one, still we observe that the activation patterns are highly similar. This indicates that our dependence measure captures information that is consistent and generalizable across subjects.
12:
206
Tier
200
2.00
(a) SUB3, C1, T1
(b) SUB3, C1,
T2
124
148
2014
200
1000
(c) SUB3, C1, T3
SUB3, C1 T4
2000
Tims
(e) SUB3, C1, T5
SUB3, C1 T6
(g) SUB3, C1 T7
(h) SUB3, C1 Averaged
Figure 11: Temporal-level dependence for nine trials from SUB3's reaching movement, confirming consistent activation of the frontal brain region and stable temporal dependence over the movement.
18

---

Channel-level dependence. Similar to Fig. of the main paper, we also quantify the channel-level statistical dependence from other subjects, not just SUB3. We randomly select clusters from subjects and visualize the results in Fig. 12. We find consistent pattern across subjects that the FC channels are activated most strongly. It can also observed that within each cluster (each subplot), channel activations are highly similar across trials. This indicates that our dependence measure robustly captures the cortical-museula. connectivity of the same movement.
(a) SUB6, C1
(b) SUB13, C9
(c) SUB17,
(d) SUB17, C4
(e) SUB6, C1
(I) SUB24, C2
Figure 12:
Channel-level dependence for visualizing activation patterns. We find that there strong activations around the FC area across subjects, not just SUB3. show that multiple clusters from various subjects demonstrate similar activation channels in the FC area. This suggests that these channels are overall the most important to classifying movements and contribute the most to connectivity.
Learning curve comparison. We show the smoothness the training stage FMCA-T, comparing its learning curve with the learning curve of MINE when applied both EEG-EMG-Fusion As shown in Fig. 13, FMCA-T demonstrates superior stability, whereas MINE suffers from greater instability even when smoothing windows are applied to estimate the gradient of the variational cost. in the Donsker- Varadhan representation.
3-class Il 11-class subj
05
Objective Fuaction
Otjerrto Fumeition
Objectvve Furcdiis
06
-30
-30
21-
6200
10030
licratios
Iteratiom
(a) MINE
(b) FMCA-T
(c) Eigenfunction Classifiers
Figure 13: Comparison of learning curves BEG-EMG- dataset. MINE: Variational cost; FMCA- Matrix trace cost; Eigenfunction Classifiers: Training emors. MINE suffers from high noise even when smoothing windows are applied to estimate the gradient of the variational cost. MINE is unable to produce stable and comparable results on Sin Wav. The classifier of eigenfunctions is trained separately from the rest of the networks.
19

---

Temporal-level dependence on San We analyze the temporal activations of the learned dependence measure on SawWav by visualizing its localized density ratios. This performed by computing the density ratios between adjacent layers of feature projectors. The layer-wise density ratios are then aggregated for visualization. We find that the localized density ratios exhibit higher activations at the hills and valleys, and correctly capture the period and phase of the sinusoids when there an increasing delay between the two sinusoids (Fig. 14). As shown Fig. 13, when Gaussian noise with standard deviation equal to 1.0 is added to the clean sinusoidal signal (signal-to-non toisse ratio less than 0 dB), the density ratio can still correctly identify the hills and valleys. This indicates that proposed dependence measure is robust to random noise and delay by filtering out trivial factors like noise while focusing on the primary signals.
1.0
1.0
1.0
0.5
0.3
0.5
fo (X)) 0.0
Vo
0.0
0.0
Normalized localizad dansity ratio
-0.5
-0.5
-1.0
-1.0
200
400
600
800
1000
200
400
600
800
1000
0
200
400
600
800
1000
1.0
1.0
0.5
0.5
0.5
0.0
0.0
0.0
-0.5
-0.5
-0.5
-1.0
-1.0
Ipput
-1.0
200
400
600
800
1000
200
400
600
800
1000
0
200
400
600
800
1000
Time
Time
Delay level 0.
Delay level 0.5
Delay level 1.0
Figure 14: Visualization of localized density ratios for samples from Wan under three delay levels. The first row shows localized density ratios based on the eigenfunctions of the input (clean) signal and the second shows the localized density ratios based on the eigenfunctions of the delayed signal. The density ratio successfully captures the period and phase of the two signals.
0.50
0.2
0.25
0.5
6(XX 0.0
Normalize Sxxmalized localized density ratio
0.00
0.0
00.5
10kg
200
400
800
1000
200
400
600
800
1000
Time
Time
Time
0.2
0.25
0.5
0.0
0.00
0.0
-0.25
0.5
-0.50
200
400
600
800
1000
200
400
600
800
1000
200
400
600
800
1000
Noise level 0.0
Noise level
Noise level
Figure 15: Visualization of localized density ratios for samples from Way under three noisy levels. The first row shows localized density ratios based the eigenfunctions of the input (clean) signal and the second row shows the localized density ratios based on the eigenfunctions of the noisy signal (Gaussian noise). The localized density ratio successfully captures the period and phase of the noisy signals even when the noise level reaches 1.0, where the signal-to-now noino ratio is less than 0 dB.
20

---

C
Implementation Details
This section includes details of data preprocessing, implementation of baselines, network struc- tures, and training configurations. Our code is available at https:/ http:://githbb com/bohu615
aarticomuscular-eggeee eigen-encoder::
C.1 EEG and EMG preprocessing
We performed standard preprocessing procedures for the 60-channel EEG signals, including 1) Down-sampling from 2500 Hz 1000 Hz, Band-pass filtering 1-48 Hz, 3) Removal of signal artefacts with independent component analysis (ICA), and 4) Segmenting EEG at (04s of the onset of the movement cue for each trial. 7-channel EMG signals were preprocessed according to the following procedure: 1) Down-sampling from 2500 Hz to 1000 Hz, 2) High-pass filtering at 5 Hz, 3) Baseline correction, and 4) Segmenting EMG 0-4s of the onset of the movement cue for each trial. The raw dataset (http:// / igada. (lataset/1007550 distributed under a CCO license.
C.2 Baseline implementations
CMC baselines. Cortico-muscular coherence (CMC) measures the linear synchronization between sensorimotor rhythms (present in EEG) and muscular activities (reflected EMG) to analyze brain-musele coupling. Given N trials, for each i-th EEG channel X1:r(n,)) i) and j-th EMG chan- nel (n, signals are filtered using a Butterworth filter and segmented into X i) and Yea-r(n, Coherence is computed the correlation coefficient (cc) between the spectral densities of the windows, eei F X ()) averaged over / to +.
It can be seen that cc is the correlation coefficients between two variables. We replace this linear measure by nonlinear measures, such as measures from KICA and mutual information estimated by KNN (MIR), producing CMC-KICA and CMC-MIR.
MINE implementation [14]. We adapt MINE with the same topology as ours but operating on the joint space of EEG and EMG X' V, producing one-dimensional output. We find that using sigmoid activation function stabilizes training. MINE computes ha X, Y) for samples from the joint distribution and ho X', Y') for samples from respective marginal distributions, and minimizes
the variational cost ming [h X, Y)] - log Efe ho Y') 10-5]. The network is optimized by Adam optimizer with learning of 10-4 ,1 = 0.5, and [3 0.9. MINE's trial-level solution follows the Donsker- Varadhan representation: log (XX, + where y can be any constant.
KICA and HSIC [17, 18]. KICA and HSIC are implemented the following steps. First, individual Gram matrices Rx and Ry are estimated for given Gaussian kernel K[XX,, X,)) =N(X; x with the standard deviation. Then, we construct the normalization matrix Ni.j and normalize the
two Gram matrices Rx NR N and Ry NRy N. For KICA-KGV, matrices A and B are constructed:
A1
0
A
A Rx Ry,
A2 Ry Rx
Ag
B
0 Bg
B
B
Rx + R) + ET), B, Ry + (T) +cl)
Then, solve the generalized eigenvalue problem for KICA Av; =a Bu, where 1,... 2N. This generalized eigenproblem generates 2N eigenvalues that are symmetric over the real line. Only N positive eigenvalues of them are used to compute the measure, obtaining KICA's Kernel Generalized Variance (KGV) measure. For HSIC, we construct matrix C:
C= 7 Bg
(10)
and solve the eigenvalue problem Cr; = U,U; where I,... N. Compute THSIC = Trace(e) C).
HSIC's measure named the Normalized Cross-Covariance Operator (NOCCO). Hyperparameters set as kernel size 0.1 and regularization constant 0.1.
21

---

Self-Supervised baselines: Self-Supervised Learning (SSL) methods are also implemented to com- pare the classification performance, including Barlow Twins [31], SimCLR [32], and VICReg [33]. We mainly use their cost functions. SSL experiments use window size of 1, 000 with windows from the same trial as positive pairs, and windows from different trials as negative pairs. The cost, hyper-parameters, and implementations follow the Lightly package [41]]
EEGNet [22]. EEGNet is convolutional neural network commonly used for EEG signal classifica- tion. After random search hyperparameter optimization, we find optimal performance with settings close to the original paper's 's recommendations. Validation set is split from the training set to enable early stopping regularization. As the original paper, we use the Adam optimizer with a learning rate of 0.0001 and decay factor of 0.1 every 20 epochs. We train and test EEGNet on the same train-test splits as the proposed algorithm for all inter-subject and cross-subject experiments.
CSP-RLDA [30]. Common Spatial Pattern (CSP) has been proven to effectively discriminate two classes of EEG by constructing optimal spatial filters. We use CSP to extract features and Regularized Linear Discriminant Analysis (RLDA) as a classifier. If class sample volumes are unbalanced, the CSP-based classifier may be biased towards the larger sample volume category. Thus particularly for CSP, one participant with one wrist-twisting session with bad EEG quality was discarded from the analysis. We use two pairs of CSP filters and extract four feature dimensions. CSP is trained and tested in both inter-subject and cross-subject settings. To achieve three-class classification with CSP, three classifiers trained for each pair of the three classes, and voting strategy determines the results.
C.3 Network structures
The structure of the temporal network is illustrated in Table 2, which consists of four convolutional blocks and max pooling. Each of these blocks and max pooling are treated layer for computing localized density ratio responses. We apply the temporal network to each channel of the signal, obtaining 1.9. Z2,,,... Zc where the total number of channels is 60 for EEG and C C= for EMG and 9 indicates the ninth layer. The output of the temporal network for each channel is vector with dimension K 128. the paper, we use the localized density ratios of Ze,6 to visualize the temporal resolution.
The channel network is three-layer MLP that takes [Z Z2.9.... c [ (dimension of K x C) as input and also produces an output of dimension K = 128. Each layer uses BN and ReLU with 2, 000 units per layer. The classifier used for eigenfunctions is also three-layer MLP with 500 units per layer.
Layer
In Ch.
Out Ch.
Kernel Size
Padding
Output
Cony, BN, ReLu
32
5
Ze,1
Maxpool
32
32
Za
Conv, BN, ReLu
32
64
11
5
Z
Maxpool
64
64
Ze,4
Conv, BN, ReLu
64
128
5
Za.5
Maxpool
128
128
2
Conv, BN, ReLu
128
256
11
5
Z,,7
Maxpool
256
256
2
Linear BN, ReLu
256 15
1024
Linear BN, ReL
1024
512
Linear, Sigmoid
512
K
Z
Table 2: Architecture of Temporal Network.
C.4 Training configurations
SanWav experiments were conducted on NVIDIA GeForce RTX 3090. REG-EMG-Fusson experiments were conducted on an NVIDIA GeForce A5000. Both Way and EEG-EMGG-Pusion used an Adam optimizer with 81 0.5 and B2 0.9 for network optimization.
22

---

D Limitation
We demonstrate the effectiveness of using FMCA to learn the dependence between EEG and EMG signals and have shown that the learned eigenfunctions embed subject and movement information after optimization. We conducted the experiments on public EEG and EMG dataset, which only contains 11 discrete upper extremity movements from 25 subjects. What we leave in the future is to use the meaningful eigenfunctions for regression tasks, i.e., continuously predicting the kinematics and contraction forces during the movement, Another limitation of the study is that we did not include patients' data due to dearth of such large datasets that collect patients' multi-modal bio-signals. We hope that the promise offered by using our dependence measurement to evaluate cortico-muscul connectivity will further stimulate experimental research in this direction. Last, from a technical point of view, only used convolutional neural networks with concatenated MLP the backbone of our networks in this study. Although we suppose that a CNN model is sufficient for the current scope as the temporal information is processed with CNN and spatial information can be leveraged by the final MLP projection layer, more advanced network structures, such as the attention in the transformer, could be potentially useful when aiming for more complex tasks.
E Broader Impact
This paper proposes to use the statistical dependence between the densities of neural data to evaluate cortico-museular connectivity. Compared with the traditional method (CMC) that computes the linear correlation between EEG and EMG spectra, our measurement shows less variance within movement and subject while is more distinguishable across movements and subjects. This helps exploring the neural information pathways from the brain to the muscle, which could further be used in the field of patient rehabilitation. Moreover, the learned eigenfunetions can be used as "common information" decoders, for example to decode movement and subject, and are more robust to distribution shift when tested on unseen subjects. This could be very useful for developing brain-machine interfaces under inter-subject conditions.
All datasets in this paper are publicly available and are not associated with any privacy or security concerns. Usages of the datasets strictly follow the corresponding licenses.
23

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: In the abstract and introduction, we claim that 1) learning the dependence between EEG and EMG signals by using FMCA produces robust dependence measure of cortico-muscular connectivity, 2) the learned eigenfunctions embed rich information of movement and subject, and 3) the temporal and spatial dependence can be visualised by computing the localized density ratios. The claims are well justified in the Sec. with theoretical support and experimentally demonstrated in the Sec. 3.2, where we show that our dependence measure is more consistent than other dependence measurements and classifiers fed with EEG eigenfunctions outperform other baseline methods. We also show temporal and channel activations that indicate the temporal and spatial distribution of the density ratio in the Sec. 3.2.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: See Conclusion and App. D for discussions on limitations of this study.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?
Answer: [Yes]
Justification: We provide the proof in the Sec. 2.1 to demonstrate the eigenfunctions decomposed from the density ratio form linear span. in the Sec. 2.2, we provide proof to derive the matrix trace cost.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided not)?
Answer: [Yes]
We describe matrix trace cost in detail in the Sec. 2.2, which is the critical com- ponent to reproduce our work. Network structures and training configurations are described in the App. C. Pseudocodes for the algorithms are provided in the App. A. We include an anonymous link (sce App. C) that provides the source codes with all implementation details and implementation of baselines.
Justification:
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: We include an anonymous link (see App. C) that provides the source codes with all implementation details and implementation of baselines. Details of data preparation are provided in the App. C.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Details of experimental settings are provided in the App. C.
24

---

7. Experiment Statistical Significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: We run experiments for five times and report the average value with standard deviation. See Table 1, Figure 3, Figure 4 in the main text and App. B for more details.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: See App. C for information on computer resources.
9. Code Of Ethies
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics https: htttps:/nneurips.oom http://s....o.b..t..//aaaaeaaaaa
Answer: [Yes]
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: See App. E for discussion on broader impacts the work.
11. Safeguards
Question:
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Answer: [NA]
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: We used one public EEG and EMG dataset in this study. We cited the original paper and provided the URL and the license of the asset in the App.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification: We used one toy dataset that generates pairs of sinusoids for preliminary study, We do not regard this toy dataset as new asset. However, we described the details of the dataset in the Sec. 3.1 and provided the codes via the link in the Appendix.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
25

---

Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
We used one public EEG and EMG dataset in this study. We have acknowl- edged the secondary use this public human dataset and included the IRB approval of the original public dataset (approved by the Institutional Review Board Korea University, 1-1--) in the manuscript.
Justification:
26


Reviews:
[{"review_id": "o1dVXdpwaM", "reviewer": "Reviewer_fPCn", "summary": "This paper introduces a novel approach to analyzing cortico-muscular connectivity using statistical dependence measures based on density ratio decomposition. The authors apply a method called Functional Maximal Correlation Algorithm with Trace cost (FMCA-T) to paired EEG and EMG recordings. The key idea is to decompose the density ratio between EEG and EMG signals into eigenvalues and eigenfunctions, which can capture important contextual information that affects the EEG-EMG dependency such as type of movement or subject without having them labeled. They also use the learned eigenfunctions as feature projectors and train a classifier on top for movement type classification tasks.\n\nThe authors test their approach on simulated data (SinWav) and a real EEG-EMG dataset with 25 subjects performing 11 different upper limb movements. They compare FMCA-T against several baseline methods for dependence estimation and classification. They find that the learned eigenfunctions capture factors such as movement type and subject identity. Further, FMCA-T outperforms the baselines, for example by 10% for cross-subject classification of arm-reaching, hand-grasping and  wrist-twisting.", "questions": "\"In scenarios where X and Y are statistically  independent, all eigenvalues are zero\" > doesn't this make the density ratio 0 then? shouldn't the density ratio be 1 if they are independent?", "limitations": "Authors could discuss a bit more under what conditions the assumption of conditional independence may be problematic and when it is fine for EEG/EMG. In terms of what one may expect to see in analyses as performed here under different conditions.\nOf course, evaluation on further datasets would also be helpful for this.", "rating": 7, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "* very sophisticated method with clear motivation\n* original idea cleanly mathematically derived (as far as I can tell)\n* produces good results", "weaknesses": "* classification baselines, could be stronger, e.g. by also using [EEG Conformer](https://pubmed.ncbi.nlm.nih.gov/37015413/)  and [Deep4](https://onlinelibrary.wiley.com/doi/full/10.1002/hbm.23730)\n* text is very dense at times, definitely found some part hard to read, but not sure how much it can be made easier, possibly you could explain some concepts used in 2.2 in more detail in the supplementary", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the comments. \n\n**1. Classification baselines could be stronger.** We have now added EEG-Conformer (available on GitHub) and Deep4 (from the braindecode repository on GitHub) as new baselines in the general response letter. The results show that FMCA-T surpasses the added baselines in almost all tasks, with Deep4 slightly ahead in cross-subject 11-movement classification. \n\n**2. Independence criterion with the eigenvalues.** The reviewer correctly pointed out our mistake. A more precise argument would be: two variables are independent (i.e., $p(X, Y) = p(X)p(Y)$) if and only if the spectrum has a single positive eigenvalue $\\lambda = 1$. \n\nDefine the linear operator $\\psi := \\mathcal{L}\\phi = \\int \\frac{p(X,Y)}{p(X)p(Y)} \\phi(X) \\, d\\mathbb{P}_X$ with $\\phi$ and $\\psi$ as functions of $X$ and $Y$. It can be confirmed that $\\phi=1$ yields $\\psi=1$, making it an eigenfunction of the operator with eigenvalue $\\lambda=1$. Two variables are independent if and only if the only positive eigenvalue is this trivial one. All eigenvalues of the operator are non-negative and bounded by $1$. More than one positive eigenvalue implies dependence. We will correct the error.\n\n**3. Applicability and limitations of conditional independence assumption.**\n\na. We demonstrated the system's ability to capture discrete and easily distinguishable factors like movement types, participants, and sessions. But it struggles with finer-grained factors, such as sub-movements that include arm-reaching in six directions, hand-grasping three objects, and wrist-twisting with two motions. \n\nOne possible explanation is that non-invasive techniques like EEG face more challenges to extract finer movement information compared to invasive techniques. Consequently, the ratio $\\rho(X,z)$ becomes trivial, and finer information may not be easily extracted from the dependent components between EEG and EMG.\n\nb. Future work is needed for scenarios with a continuous conditioned variable $\\mathbf{z}$. In Lemma 1, we proposed the expression $\\frac{p(X,Y)}{p(X)p(Y)} = \\sum_{z\\in \\mathcal{Z}} p(z) \\rho(X, z)\\rho(Y, z) dz = \\sum_{k}\\lambda_k \\phi_k(X)\\psi_k(Y)$, which applies when $\\mathbf{z}$ is discrete. This may not account for continuous variables such as the force of muscle contractions and kinematics (joint angle). \n\n**4. Density of texts.** We appreciate the reviewer's feedback and will improve the formatting.", "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": "", "extra_scores": {}}, {"review_id": "OJ4jYW86Pk", "reviewer": "Reviewer_MC3m", "summary": "The authors apply novel but already existing (https://www.sciencedirect.com/science/article/pii/S0047259X2300074X, https://arxiv.org/pdf/2212.04631) machinery based approach on the orthonormal decomposition of density to decipher the relationship between cortical brain activity and the electromyographic signals during basic hand movements. The work is based on the publicly available dataset and the code is available. The unknown decomposition is modeled by a pair of neural networks concurrently processing EEG and EMG data in order to arrive at the internal representation for each of the modality. The internal representations are then aligned to minimize the rank of the joint covariance matrix. To guide the learning process, the authors propose a somewhat novel loss function equal to the negative trace of the canonical correlation matrix calculated using latent representations. The authors test their approach on the downstream tasks of classifying movement types in both within and across subject designs. They also apply the obtained representations to distinguish between participants based on their EEG data. The authors provide some interpretation to the obtained solution in the form of channel and temporal maps indicating the electrodes and time moments that contribute to the decoding most.", "questions": "1. Having significant experience in the domain of recording and analyzing electrophysiological data I founnd the obtained maps very suspicious. While EEG electrode FC1 can indeed be implicated and be coherent with EMG, I would expect other electrodes such as C3, C5 to have some significant contribution to the EEG derived latents that are maximally aligned with EMG. Instead in addition to FC1 we see the involvement of peripheral electrodes and the frontal electrodes. These electrodes often lose proper contact with the skin and become sensitive to the physical movements due to capacitive effects, when slight body displacements during the actual movement causes significant fluctuation in the electrode-skin capacitance and modulates the signals registered by EEG. The analysis of frequency response of the temporal layer (see W.5) may help to resolve this potential issue.\n\n2. In the dataset used by the authors the reference channel was located in the midline between FC1 and FC2 sensors. Such an arrangement often results in low variance of the signals located close to the reference. The spatial patterns that the authors demonstrate in Figures 5 and 11 show peaks around Fc1 and FC2. The ground electrode was located at the edge of the EEG cap between F1 and F2 and the temporal maps show them as the next best electrodes after Fc1. Could it be that within certain normalization steps the authors explicitly or implicitly  divided the data by the channel variance or multiplied the data by a poorly conditioned and not properly regularized inverse covariance that the role of these electrodes got artificially inflated?\n\n3. The authors show spatial maps for several other *selective* subjects and clusters to illustrate across subject reproducibility. What about the spatial maps corresponding to the other latent channels\\clusters for SUB3? \n\n4. Why did the authors not follow the EEGNet architecture and decided to use non-linearities between the temporal and spatial processing blocks? Avoiding nonlinearities in the front end would improve interpretability and would help to make the presentation more convincing. \n\n5. How do the authors avoid the trivial training result, i.e. that the two networks will simply learn to generate similar EMG and EEG embeddings regardless of the input data?", "limitations": "The authors partly addressed limitations but several items, see Weaknesses section, are left out.", "rating": 5, "confidence": 4, "soundness": 2, "presentation": 3, "contribution": 1, "strengths": "1. The authors applied novel but existing methodology of orthonormal density decomposition to the EEG+EMG dataset for the first time \n2. The authors introduced a novel loss function and showed that it provides better performance in the downstream task of EEG-based classification of movement types\n3. The authors used multi subject dataset \n4. The authors attempted to provide interpretation of the obtained decision rule\n5. The authors present detailed results of their experiments in the appendix", "weaknesses": "1. Several inaccuracies and lack of details in the mathematical expressions:\n1.1 line 100, last expression and additional p(z) is needed in the integral\n1.2 equation 3  - do the eigenvalues need to be normalized? Does the sum exclude the first normalized eigenvalue?\n2. I would argue against the suggested novelty of the proposed loss function as it seems like the loss function that is closely associated with the Canonical correlation analysis (CCA) (equation 4). Generally speaking, the proposed approach boils down to the CCA in the latent variable space with latents computed by means of a CNN.\n3. The authors claim that they “..design a specialized network topology to generate features for individual channels and time intervals, ensuring that the internal layers of this network quantify channel-level and temporal-level features, similar to [22–24]” - however unlike for instance the EEGnet, the authors use non-linearities in the temporal network (prior to the spatial) which in my view prevents the straightforward interpretation of the obtained representations at least using simple correlational measures. \nSee also Q.1 and 2. \n4. The authors did not validate their approach to interpreting the decision rule and obtaining spatial and temporal maps with simulated data. This needs to be done and the simulated data should contain not only the neuronal sources coupled to the simulated EMG but also the sources unrelated to the signal of interest (EMG). The authors then need to demonstrate that their methods infers the proper spatial patterns corresponding to the task-relevant simulated neuronal sources. Ideally,  the obtained maps should be ranked based on their importance for the overall decoding accuracy. If this is not possible within the review cycle, the authors should significantly reduce the proportion of the manuscript dedicated to the physiological plausibility of their solution and instead describe limitations related to potentially non-physiological origin of the extracted features.  \n5. It is disappointing that when interpreting the decision rule the authors did not provide information regarding the EEG frequency domain their network got tuned to during the training.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the detailed insightful feedback. \n\n**1.1 Novelty over CCA.** We acknowledge the reviewer's observation that the final cost is the trace of a normalized canonical correlation matrix between two multivariate neural networks. \n\nBut we emphasize the link between cost optimization and joint density ratio $\\frac{p(X,Y)}{p(X)p(Y)} = \\sum_{k=1}^K \\lambda_k \\widehat{\\phi_k} (X) \\widehat{\\psi_k} (Y)$, which the reviewer might have misunderstood.\n\nVanilla CCA uses linear models to maximize a correlation coefficient. After KICA was introduced (i.e., Kernelized CCA), the problem became minimizing matrix costs using universal approximators, with costs including log determinant, Frobenius norm, and matrix trace. Our trace cost is closer to HSIC. We have thoroughly compared our method with KICA and HSIC as baselines, both of which use kernels while we use neural nets.\n\nOnly recently, studies proved that optimizing matrix costs is equivalent to factorizing joint density ratios, including FMCA and Gaussian universal features. An informal proof: if $\\widehat{\\phi}$ and $\\widehat{\\psi}$ are two orthonormal sets, optimizing the norm of $\\mathbb{E}[\\mathbf{\\widehat{\\phi}}(X)\\mathbf{\\widehat{\\psi}^\\intercal}(Y)] = \\int \\mathbf{\\widehat{\\phi}}(X)\\mathbf{\\widehat{\\psi}}^\\intercal(Y) p(X,Y)dXdY$ decomposes the joint density ratio $\\frac{p(X, Y)}{p(X)p(Y)}$ and finds its eigenfunctions.\n\nThis paper further shows that all such costs, including log determinant, Frobenius norm, and trace, are fundamentally equivalent, differing only in the convex functions applied to eigenvalues. We consider this a major novelty.\n\n**1.2 More discussions about novelties.** Clarifying the question on $\\mathbf{z}$, it is the variable making $\\mathbf{X}$ and $\\mathbf{Y}$ conditionally independent, such as movement types, not the network's latent features that approximate eigenfunctions. Lemma~1 explains why eigenfunction can be used for downstream tasks. Another novelty is analyzing dependence within the network for temporal and spatial activation. \n\n**1.3 Normalization.** Eigenvalues of $\\frac{p(X,Y)}{p(X)p(Y)}$ are well-regularized and don't need normalization: (1) All eigenvalues are positive and bounded by 1; (2) $\\lambda_1=1$ regardless of dependence; (3) independence iif there exists only one non-zero eigenvalue. This is due to the implicit normalization dividing $p(X,Y)$ by $p(X)p(Y)$, corresponding to the cost's normalization $\\mathbf{R}_F^{-\\frac{1}{2}} \\phi$ and $\\mathbf{R}_G^{-\\frac{1}{2}} \\psi$. \n\nThe normalization is costs also enforces orthonormality constraints on features, so outputs won't be constants (trivial solution)\n\n**2.1 Simulated dataset.** We have added additional experiments with simulated dataset, shown in Fig. 3 in the general response. \n\nWe simulated EEG and EMG signals for left/right motor and sensory activations in 20 subjects using EEGSourceSim. Motor sources were used to simulate the corresponding EMG signals. We trained the networks on paired EEG-EMG samples from 16 subjects' data as training samples, and visualized the spatial-level dependence maps for the remaining 4 subjects as test samples. Then we plotted ground truth brain activations calculated from motor ROI and forward matrices.\n\nThe spatial-level dependence is highly similar to ground truth activations, indicating the learned ratio captured real brain activations. \n\n**2.2 Dependence shows major activations around FC areas.** We have added a frequency analysis using event-related desynchronization (ERD) in the beta band (Fig. 2) and complete spatial dependence maps of our method (Fig. 3), both for Subject 3.\n\nThe frequency topology shows that left-central areas (e.g., C3) are commonly activated across sessions and movements while frontal-central (FC) and central-parietal areas have unique patterns for each movement. In comparison, our dependence maps show the most activated areas around FC areas.\n\nWe argue that the difference in activation occurs because we are quantifying ***nonlinear dependence***, not correlations. It has been experimentally proven that the eigenfunctions of this dependence yield high accuracies in downstream classification tasks, which means that the activation maps should show the areas that make each movement most distinguishable. Thus the frontal central and central parietal areas are activated in dependence maps. In frequency topologies, C3 is activated for all movements, but it's unable to be used as an identifier for individual movements. This causes the differences between our maps and frequency maps. More details can be found in the general response.\n\nThe public dataset used has been validated before publication. No activation or deactivation was found on peripheral or frontal electrodes in the frequency analysis-based topographies.\n\n**2.3 Electrode arrangement and normalization steps.** The electrode arrangement is reliable and has been used in commercial equipment. In EEG analysis, if the target electrode is near the reference, a weighted average filter will be employed to enhance the signal. Our study instead had no clear target before analysis. Therefore, we applied an average reference across electrodes, which is commonly used in EEG analysis. During preprocessing, all EEG channels were normalized to the first sample's maximum amplitude. Trials exceeding an absolute amplitude of 5 were discarded. No channel-wise normalization was used that could affect the variance. \n\n**3.1 Choice of nonlinearity differs from EEGNet.** In EEGNet paper, the authors found that nonlinear activations did not improve performance, which is why they chose linear functions, thus making it a matter of preference.\n\nAlthough linearity may help prevent overfitting in cases of long signal durations, our experiments showed no significant overfitting. Using ReLU or Sigmoid may provide a more stable approximation when computing cross-layer ratios, as the ratios should be non-negative and bounded.", "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": "", "extra_scores": {}}, {"review_id": "6ZKw0GUOyc", "reviewer": "Reviewer_J9kE", "summary": "The paper presents a new method to model the relationship between cortical and muscular oscillations using EEG and EMG recordings. Traditional methods like Cortico-Muscular Coherence (CMC) have limitations, so the authors propose using statistical dependence estimators to learn eigenvalues, eigenfunctions, and projection spaces. This approach improves interpretability, scalability, and local temporal dependence. Experimental results show that the method accurately classifies movements and subjects, highlighting specific EEG channel activations during movements, and demonstrates robustness against noise and delays, suggesting its potential for diagnosing neuromuscular disorders and developing brain-computer interfaces.", "questions": "1. The method identifies specific EEG channel activations during movements. Could you provide more detailed explanations or references to how these findings align with existing neuroscientific knowledge? What are the physiological implications of the identified activations, and how do they contribute to our understanding of cortico-muscular connectivity?\n\n2. The author mention potential applications in diagnosing neuromuscular disorders and developing brain-computer interfaces. Can you provide concrete examples or case studies where your method has been or could be applied? What specific benefits or improvements does your method offer over existing approaches in these applications?", "limitations": "Plz go and check weaknesses and questions.", "rating": 5, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "1. The paper combines statistical dependence estimators with neural network optimization techniques. This fusion of methodologies enhances the ability to capture high-level and contextual connectivity between cortical and muscular oscillations.\n\n2. The paper provides a detailed description of the proposed methodology, including the mathematical foundations, algorithmic implementation, and practical considerations. The inclusion of eigenvalues, eigenfunctions, and projection spaces adds depth to the analysis.\n\n3. The authors conduct comprehensive experiments to validate their method. The results demonstrate the method's robustness against nonstationary noise and random delays, confirming its reliability and practical applicability.", "weaknesses": "1. Mathematical and Algorithmic Complexity: The proposed method involves complex mathematical formulations and advanced statistical techniques that may be challenging for a broader audience to grasp. Simplifying some of the mathematical derivations or providing more intuitive explanations and visualizations could make the paper more accessible.\n\n2. Interpretation of Results: While the method highlights specific EEG channel activations during movements, the physiological and neuroscientific significance of these results could be further elaborated. Providing more detailed discussions on how these findings align with or differ from existing neuroscience research would enhance the interpretability and relevance of the results.\n\n3. Scalability: The scalability of the proposed method to larger datasets or longer signal durations is not thoroughly addressed. Discussing the computational complexity and providing benchmarks on how the method performs with varying data sizes would be valuable.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the instructive comments. Please find below our replies to the concerns/questions.\n\n**1. Mathematical and algorithmic complexity.** We thank the reviewer for the suggestion. We will add more explanations of our methodology to the supplementary.\n\n**2. Interpretation of results and physiological implications of the identified activations.** In terms of interpretation, we found that the frontal central areas (FC) are most activated in the spatial-level dependence maps for most subjects. Since we used eigenfunctions that decompose the density ratio for classification, this implies that these FC areas contribute significantly to movement classification. It is reasonable for FC1 to show the most activation in the dependence map. A related previous study also found movement-related cortical potential changes on FC1-FC2 and C2 electrodes [1].\n\nWe have additionally added a frequency analysis using event-related desynchronization (ERD) in the beta band. In the general response, we addressed the differences and consistency between frequency topologies and our dependence ratio map.\n\n**3. Scalability.** We appreciate the reviewer's suggestion. Given the limited time frame and the limited availability of paired EEG-EMG datasets, it is difficult to show scalability on new real-world datasets. Instead, we divided the dataset into smaller sub-datasets and ran experiments with increased dataset sizes to demonstrate how our framework scales up. Fig. 4(a)~4(c) in the attached letter show that: (a) the convergence speed of training errors is almost the same, with smaller datasets converging slightly faster; (b) downstream classification accuracies steadily increase when using more data; (c) dependence decreases when using more data. \n\nWe observed the expected trend in dependence scores - as more data is used, it implies greater uncertainty in the data.\n\nThe computational complexity was not particularly heavy, and we did not see it as an obstacle. We trained the model on an A5000 GPU and did not encounter significant difficulties.\n\n**4. Potential applications.** Corticomuscular coupling is promising for quantitively measuring movement disorders. This application has been validated in stroke populations [2, 3] and in Parkinson’s disease [4]. Therefore, the proposed dependence measure can be used as an effective and clinically relevant neural marker to evaluate movement performance of stroke patients and to indicate Parkinson’s disease pathology. Moreover, as demonstrated in the manuscript, EEG eigenfunctions show excellent performance in movement classifications. Therefore, the proposed method also has great potential in improving brain-computer interfaces. We will add this discussion in the revised manuscript.\n\n*[1] Spring J N, Place N, Borrani F, et al. Movement-related cortical potential amplitude reduction after cycling exercise relates to the extent of neuromuscular fatigue[J]. Frontiers in human neuroscience, 2016, 10: 257.*\n\n*[2] Chen X, Xie P, Zhang Y, et al. Abnormal functional corticomuscular coupling after stroke[J]. NeuroImage: Clinical, 2018, 19: 147-159.*\n\n*[3] Liu J, Wang J, Tan G, et al. Correlation evaluation of functional corticomuscular coupling with abnormal muscle synergy after stroke[J]. IEEE Transactions on Biomedical Engineering, 2021, 68(11): 3261-3272.*\n\n*[4] Zokaei N, Quinn A J, Hu M T, et al. Reduced cortico-muscular beta coupling in Parkinson’s disease predicts motor impairment[J]. Brain communications, 2021, 3(3): fcab179.*", "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": "", "extra_scores": {}}, {"review_id": "bCUJMCQqwZ", "reviewer": "Reviewer_CLW7", "summary": "The paper presents a novel approach called Functional Maximal Correlation Algorithm with Trace cost (FMCA-T) for estimating cortico-muscular dependence by leveraging orthonormal decomposition of density ratios. This method is designed to model the relationship between EEG (electroencephalography) and EMG (electromyography) signals, addressing the challenges of interpretability, scalability, and local temporal dependence in cortico-muscular connectivity. The key contributions include introducing a matrix trace cost optimization for improved stability and efficiency, demonstrating robustness against nonstationary noise and delays, and effectively capturing movement and subject information from EEG features for enhanced classification accuracy. The proposed method outperforms existing baselines, particularly in cross-subject scenarios, and provides insights into channel-level and temporal dependencies, reinforcing its potential applications in brain-computer interface development and neuromuscular disorder diagnostics.", "questions": "1. While the paper discusses the improved stability and efficiency of the FMCA-T method, can the authors provide more detailed about the computational resources required for training and inference?", "limitations": "1.Limited Dataset Size: The cross-subject classification performance drops, potentially due to the limited dataset size of only 25 participants. Larger and more diverse datasets may be needed to validate the method’s robustness comprehensively.\n2. Generalization to Other Modalities: While the method shows promise for EEG and EMG signals, its generalizability to other types of biosignals or broader neural data modalities has not been extensively disscuss.", "rating": 7, "confidence": 4, "soundness": 4, "presentation": 3, "contribution": 4, "strengths": "1. Innovative Method: Introduces the Functional Maximal Correlation Algorithm with Trace cost (FMCA-T), providing a novel approach for estimating cortico-muscular dependence.\n2. Improved Stability and Efficiency: Utilizes matrix trace cost optimization, which is more stable and computationally efficient compared to traditional log-determinant cost methods.\n3. Enhanced Classification Accuracy: Effectively captures movement and subject information from EEG features, significantly improving classification accuracy, especially in cross-subject scenarios\t\n4. Validation on Multiple Datasets: Validated using both simulated and real EEG-EMG datasets, confirming the method’s effectiveness and robustness.\n5. Open Data and Reproducibility: Offers open access to datasets and detailed implementation code, facilitating reproducibility and further research in the field.", "weaknesses": "The provided baselines are relatively few; future work could expand on this.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the comments. Please find below our replies to the concerns/questions.\n\n**1. Limited baselines.**\n\nWe have added EEG-Conformer (available on GitHub) and Deep4 (from the braindecode repository on GitHub) in the attached pdf. All baselines were implemented following the official codes. Results indicate that FMCA-T outperforms the additional baselines on almost all tasks. Deep4 performs slightly better on cross-subject 11-movement classification task. We will update Table 1 in the revised manuscript.\n\n**2. Computational resources for training and inference.**\n\nExperiments on the simulated sinusoidal dataset and the additional simulated EEG and EMG dataset were conducted on an NVIDIA GeForce RTX 3090. Experiments on the experimental EEG and EMG dataset were conducted on an NVIDIA GeForce A5000. All training and inference can be conducted on a single GPU.\n\n**3. Limited dataset size.**\n\nWe agree with the reviewer about needing a larger, diverse dataset for better validation. Although unavailable now, we hope our method will encourage large-scale datasets, which we plan to explore further.\n\n**4. Generalization to other modalities.**\n\nWe thank the reviewer for proposing the possible extension to broader neural data modalities. This is an interesting question for further investigation. We expect our method can still apply to learning the dependence between other modality data, e.g., between EMG and kinematics or between EEG and visual data. We will continue studying such scenarios in future work.\n\n**5. Ethics review.**\n\nThe public EEG and EMG dataset was reviewed and approved by the Institutional Review Board at Korea University (1040548-KU-IRB-17-181-A-2).", "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": "", "extra_scores": {}}]

Author rebuttal:
We appreciate the reviewers' feedback. The following responses address their shared concerns.

We have attached a letter containing the additional results, including the requested classification baselines, a frequency analysis of brain activations, full maps for Subject 3, maps for simulated EEG-EMG data, and model scalability with dataset size.

**1. Additional classification baselines.** We have added EEG-Conformer and Deep4 results as baselines. FMCA-T surpasses these supervised methods in most tasks, except Deep4's slight advantage in cross-subject 11-movement classification. This is shown in Table 1 of the attached letter.

**2. Additional frequency analysis.** To address the reviewer's concern regarding frequency topologies, we have performed a frequency analysis of brain activations using event-related desynchronization (ERD) in the beta band.

Fig. 1 shows results for Subject 3 across all movement types and sessions. We observed that the left central area is commonly activated across all movement types, frontal-central (FC) and central-parietal areas have unique patterns for each movement, and no apparent activation or deactivation was found in peripheral or frontal electrodes.

This is consistent with conventional EEG analysis and validates the dataset.

**3. Dependence map interpretation.** We have added dependence maps for all nine clusters from Subject 3 in Fig. 2, as requested by the reviewer. FC areas are consistently activated in most maps, matching our main paper.

To explain why the activated areas (FC) may differ from frequency topologies (left central), note that we are measuring nonlinear dependence, not correlations. Eigenfunctions for this dependence have been proven to yield high classification accuracies in downstream classification tasks, so we expect the activated areas in the dependence map to contribute significantly to classifying and identifying movements and participants, thus leading to FC activation.

Frequency topologies potentially support this argument. Since C3 is activated across all movement types, it may not be an effective identifier for distinguishing individual movements. Instead, FC has unique activation patterns for each movement.

**4. Additional results on simulated EEG-EMG dataset.** We have added results from a simulated EEG-EMG dataset to validate the spatial-level dependence map. We simulated EEG and EMG signals for left/right motor and sensory activations in 20 subjects using EEGSourceSim. Motor sources were used to simulate the corresponding EMG signals. As shown in Fig. 3 in the letter, the spatial-level density ratio learned from the dataset shows similar activation patterns to the ground truth motor activations computed from motor ROI and forward matrices.

Meta-review: Perhaps the best paper of my lot and the only one I give some consideration as to whether recommend it for oral. After some consideration, I have favoured recommending for poster though because of the heavy mathematical content which is usually better communicated over a poster than in an oral.
Range of rating was only 2 points difference suggesting good agreement between reviewers. Rounded average score was 6 (accept)
Based on the reviews a major strength is the pragmatical good results (R-MC3m and R-fPCn) attributable to the novel cost function (R-MC3m) as well as its practical applicability (R-J9kE). Two of the reviewers timidly mentioned mathematical complexity/sophistication (R-J9kE and R-fPCn), and it is certainly a maths heavy draft to the point that even to person with mathematical inclination like me, feels abusive at times. I truly believe clarity could have been improved by obfuscating less the message with the abuse of maths.
Because the reviews gave more praise than requested further clarifications, rebuttal was relatively succinct. 
It was pointed to the authors that the IRB approval cannot be hidden (the authors did disclosed on the rebuttal though), but otherwise, no further issues with ethics were highlighted.
Decision: Accept (poster)