Neural decoding from stereotactic EEG: accounting for electrode variability across subjects

Thank you for a very detailed review. Here, we provide our response to your questions. Due to space constraints, for each question we provide only the beginning of your prompt. We apologize for the inconvenience.

Throughout this rebuttal we refer to single-subject as SS and multi-subject as MS.

Weaknesses

1. A main weakness of this paper is the spatial encoding...

Thanks for the constructive criticism. Indeed, our proposed PE only marginally boosts decoding performance, suggesting the attention mechanism in space is the main mechanism by which our model deals with electrode heterogeneity. We will highlight this in our revised manuscript.

2. Furthermore, the ablation study detailed in A.3.1...

This is a valid point. We will highlight this in revisions.

3. Due to the marginal performance boost offered by positional encoding...

Thanks for suggesting this. We compared the performance of our architecture against a variant where 2D attention is performed over time & space. We include the results of this analysis in Table 4 of our General Response (see also Figure 1E of the pdf). The results suggest that our proposed architecture outperforms its 2D attention variant by $\Delta R^{2} = 0.06$.

4. Comparison against other methodology is limited...

Thanks for the constructive criticism. To address this, we compared our model with classic and SOTA SS models. The results, shown in Table 1 of our General Response, indicate that our architecture outperforms all baselines across the board.

Unfortunately, UKF and SEEG-Net are missing from those baselines. We have been unable to get good results using UKF, likely due to the nature of our data, since the trial response times (included in UKF's state) are not stationary, making state estimation unreliable. We also looked into sEEG-Net, but the code is not public and the method requires substantial modification before use, since it operates on the electrode level while our problems requires operation on the subject level. We are still trying to make those models work and hope to include them as baselines in revisions.

Thank you for suggesting we discuss multi-task literature in revisions. We plan to do so, including Crawshaw (2020), in the related work section of our manuscript.

Minor

• Line 210: "long-range dependency" might be more appropriate...

We will amend this in revisions.

• Please provide the test/train/validation splits...

The training/validation/test split is 70/15/15. We will include this in revisions.

Questions

1. Line 230: What is C? Should it be E?

Yes. We will fix the typo in revisions.

• I am assuming that an MRI for each subject in the study was available...

This is a great idea. Unfortunately, for this dataset we only have the MNI coordinates of electrodes and not whole MRIs or other imaging data. Therefore, we are unable to test this out. We look forward to exploring this option in the future.

• I do not understand the point of using m-different 1-D RBFs for modeling the positional encoding...

Your intuition is correct; we used RBFs with different variances to encode position at different spatial scales. This is relevant because studies (Meunier (2009)) show that the brain has a hierarchical modular structure. Therefore, encoding both short and long electrode connections could be useful. The scales $\sigma_{j}$ for each RBF were hyperparameters with values: [1, 2, 4, ..., 64].

• Is the rationale for using 1-D RBFs, i.e., the RBF for each coordinate...

We used 1D RBFs instead of 3D because 1D RBFs: (1) inform the model with direction (as you correctly identified) and (2) introduce less parameters to the model compared to 3D RBFs. Specifically, assuming that sampling an interval in 1D requires m points, sampling a 3D interval, separately in each dimension, at the same scale requires 3m points. Sampling the same 3D interval in 3 dimensions combined requires m$^{3}$ points.

• Lines 205-207: The authors state that the projected positional encodings are added to the latents...

Thank you for this insightful suggestion. Your description of the operation correct. To ensure that $z^{3}$ & $p$ have comparable magnitudes, we computed $||z^{3}||$ and $||p||$ for 512 random samples from our test set. The norms were $||z^{3}|| = 0.790 \pm 0.53$ & $||p|| = 0.013 \pm 0.01$ (mean + sdev), suggesting that $||z^{3}||$ is greater than $||p||$ but not to the extent that the model's output would be unaffected by this operation.

3. Why is Huber loss used for training multi-subject model whereas MSE loss used for training single-subject...

We optimized SS models using MSE Loss because $R^{2}$ is directly proportional to the MSE Loss. Following your suggestion we trained both SS & MS models using both MSE & Huber Loss. The results show that training with Huber Loss boosts performance for all models. We have update the SS model results to reflect this.

Table 5. Comparison of model performance trained with Huber vs MSE Loss

4. Another point of potential unfair comparison...

Thank you for this thoughtful concern. Your intuition is correct. However, the MS and SS models were trained with batch size 1024 and 64, respectively (see line 599 of manuscript). Our dataset contains 3685 trials (175 trials per subject). Therefore, during 1K training epochs MS and SS models underwent 4K and 3K gradient descent steps (GDS), respectively.

To ensure that this did not artificially inflate the performance difference between the MS and SS models, we retrained all SS models for 1.5K epochs (4.5K GDS). The average per-subject test set R$^2$ of the SS models was 0.29 $\pm$ 0.05 (mean $\pm$ sem). The MS model still outperformed those by $\Delta R^{2} = 0.10$.

In practice, we used 1K epochs because we observed that they were enough for the training objective of all model to converge.

Thank you for your thoughtful feedback and question. We were really excited to read that "The proposed framework effectively addresses the heterogeneity across subjects, a significant challenge in sEEG data processing". Here, we provide our response to your questions and the results from new experiments that we ran to address your concerns. If you would like any further clarifications, please let us know.

Weaknesses

The paper discusses the impact of the positional decoding in the ablation study, but this is the least impactful structure compared to other layers (Fig 6). It would be beneficial to explore other position decoding mechanisms besides using MNI locations for future work, perhaps based on different brain atlases.

Thank you for the insightful suggestion! We echo your observation that our proposed spatial positional encoding only improves decoding performance slightly and agree that future work should identify better ways to encode the location of sEEG electrodes in the brain. Using different brain atlases is a great suggestion that we will definitely look into. Unfortunately, for this dataset, we only have the MNI coordinates of the sEEG electrodes and therefore we are not able to test whether encoding position with some other atlases would achieve better decoding results. However, we certainly look forward to exploring such ideas in the future.

The best we could do to adress your feedback was to compare our spatial positional encoding approach against other positional encoding approaches. The results are summarized in Table 2 of our General Response and suggest that our proposed positional encoding, works as well as, or better than other approaches. However, as you pointed out, the explored positional encoding schemes are far from comprehensive, indicating that there is still a lot of margin for improving spatial positional encoding mechanisms. We will make sure to discuss the idea of encoding electrode location using other brain atlases in section 4.5 of our revised manuscript.

Questions

1. With the number of electrodes differing across subjects, how did you handle the different input sizes?

This is a great question. Our models is capable of handling different input sizes due to the ability of convolutional neural networks and transformers to accept inputs of varying lengths. Prior to procesing the latents with our architecture's MLP, we pad the latents of all subjects to a fixed length. For the mathematical details, we refer you to section 3.2 of our manuscript. Here, we try to provide a more intuitive explanation of why our network is capable of handling inputs of different sizes.

The input to our network, $X \in R^{E \times T}$, is a 2D array of electrodes $\times$ timepoints, where $E$ refers to the electrode dimension and $T$ to the time dimension. The tokenizer performs $K$ 1-D convolutions on the input along the time dimension $T$, while parallelizing the computation across the $E$ dimension (electrodes). This operation returns a latent $z \in R^{E \times T \times K}$, where $K$ refers to the number of convolutional kernels. Then, self-attention is performed on $z$ along the time dimension, while computation is parallelized across the $E$ dimension. Self-attention can accept sequences of arbitrary lengths and therefore no matter the number of timepoints of the latent $z$, this operation can be performed successfully. The output of the self-attention has the same size as the input, returning a latent of the form $z \in R^{E \times T \times K}$. Then positional encodings are added to the latents, which does not alter the shape of $z$. Another self-attention operation is performed, this time along the $E$ dimension of $z$, parallelized across the $T$ dimension (timepoints), which again preserves the size of the latent $z \in R^{E \times T \times K}$. At this point, the latents are unrolled to get a latent of the form $z \in R^{E \cdot T \cdot K}$. $T$ and $K$ have the same dimensionality across subjects; $E$ however (the number of electrodes) does not. Therefore, each subject is going to have a latent $z \in R^{E \cdot T \cdot K}$ that has a unique dimensionality. This latent will be projected to a lower dimensional space through a multi-layer perceptron, which requires a fixed input length. Therefore, at this stage the latents $z \in R^{E \cdot T \cdot K}$ are padded with zeros to obtain a fixed length $E_{max} \cdot T \cdot K$ that is common to all subjects. For convenience, $E_{max}$ is set to the maximum number of electrodes that a subject has in our cohort. The MLP projects the latent to $z \in R^{F}$, where $F$ represents a small number of features that will then be projected through the subject-specific MLPs to a single number corresponding to the response time of a subject for a given trial. We hope that this explanation answers your question.

2. Instead of plotting all the sEEG channels, it would be helpful to show the selected channels actually used in training.

This is a great suggestion. Please take a look at Figure 2 of the pdf we submitted with our General Response and let us know if you have any further suggestions. We will appreciate your input on this.

3. What is the (relative) root mean square error for the prediction?

For the multi-subject model, the root mean square error for all predicted response times in the test set, across subjects is RMSE = 0.082 $\pm$ 0.002 (mean $\pm$ sem) sec. For reference, the mean response time in the test set across subjects is $0.41$ $\pm$ 0.005 sec (mean $\pm$ sem). We will make sure to include this in section 4.3 of our revised manuscript.

Thank you for your though provoking feedback. We appreciate reading that "combining data from multiple individuals and training a unified model is a substantial step forward compared to traditional single-subject approaches"! Here, we provide our response to your questions and the results from new experiments that we ran to address your concerns. We are available to provide further clarification, if needed.

1. Scale and diversity of the dataset: The paper mentions using data from 21 subjects, totaling 3600 behavioral trials and 100 hours of sEEG recordings. Is this amount of data sufficient to represent the neural activity patterns of all subjects? Could the diversity of the dataset be insufficient to support the model's generalization capabilities?

This is a very thoughtful concern. It is likely that the amount of data in our dataset is insufficient to represent the neural activity patterns of all humans in a global sense. However, our dataset is larger than that of other recent sEEG datasets (Angrick (2021), Petrosyan (2022), Meng (2021), Wu (2022), Wu (2024) have data from 12 subjects or less), suggesting that it is more likely to represent general neural activity patterns compared to most other sEEG datasets.

In terms of our dataset's diversity: (1) our cohort is composed of 13 females and 8 males, (2) subjects' age ranges from 16 to 57 years old, and (3) electrode number & placement for each subject is unique (solely based on clinical needs). Across subjects electrodes span white and grey matter; cortical, subcortical, and deep structures; and every major structural region of the brain (see Figure 2 of our General Response pdf). Our dataset is very diverse, suggesting that our model, trained on the combined data of all subjects, is likely to generalize.

Thanks for the food for thought, we will include this discussion in our revised manuscript.

2. Effectiveness of spatial positional encoding: The paper proposes a spatial positional encoding method based on MNI coordinates to handle electrode placement variability across subjects. Has the effectiveness of this method been thoroughly validated? Are there comparative experiments showing that this encoding method is superior to other possible encoding methods?

Thanks for pointing out this limitation. In response to this concern, we validated our spatial positional encoding against other positional encodings. The results are summarized in Table 2 of our General Response and suggest that our proposed positional encoding, based on MNI coordinates, works as well as or better than other approaches.

We believe it is also important to note that the MNI-Fourier and MNI-RBF positional encoding schemes outperform other approaches, indicating that informing models with sEEG electrode location boosts their performance. However, the performance gains are small, indicating that there is still a lot of margin for improving spatial positional encoding mechanisms. We will make sure to discuss these findings in section 4.5 of our revised manuscript.

3. Individualized decoder heads: The paper mentions that each subject has a personalized task head for downstream decoding tasks. How effective is this approach when dealing with new subjects? Are there detailed experimental results demonstrating the performance of this method on new subjects?

This is a great question. We echo that the ability of dealing with new subjects is very important to support our model's generalization capabilities.

In section 4.4 of our submitted manuscript we demonstrate that our approach is effective when dealing with new subjects. Specifically, in that section, we take a leave-one-out cross validation approach in which we train 21 models (one for each subject in our cohort) each of which is trained on the combined data of all subjects except one. Then, we use the weights of the pre-trained models as a basis of finetune, upon which we train a new single-subject model with the data of the left-out subject. Across subjects, the models trained on other subjects and finetuned to a new one achieved an average test set $R^{2}$ score of $0.37 \pm 0.06$ (mean $\pm$ sem). Importantly, the performance of those models is superior to that of models trained from scratch and outperforms all baseline models (see Table 1 of our General Response) by $\Delta R^{2} >= 0.10$.

4. Computational complexity and scalability: The paper outlines a complex model architecture involving convolutional tokenization, self-attention in both time and electrode dimensions, and individualized regression heads. What are the computational requirements for training and running this model? Is the approach scalable to larger datasets or real-time applications?

The computational resources used to train our model are provided in lines 617-619 of our submitted manuscript. Single- and multi-subject models train within 5 mins and 1 hour, respectively, on a machine with an AMD EPYC 7502P Processor and one Nvidia A40 GPU, which is well within the resources of most computational research labs. The memory requirement to train the multi-subject model with a batch size of 1024 is $\sim$ 8 GiB, making model training tractable using less powerful GPUs as well.

In terms of scalability, assuming that the batch size remains fixed, irrespective of the dataset size, the aforementioned computational resources should be enough to train the model. Naturally, the time requirement for training would scale linearly with the number available trials. However, even a dataset of $\sim$ 80K trials would train in roughly 24 hrs, which is very manageable considering the dataset size.

To investigate whether our model can be used in real time, we measured the model's inference time on a sever and a laptop (see Table 3 of our General Response). All inference times are on the order of a few msec, which easily allows for integrating our model with real-time systems. We will make sure to add this result to our revised manuscript.

Thank you for your response and for updating your score! We are very happy that our responses and new experiments addressed at least some of your concerns. Please find our answers to your new comments below:

Thank you for the authors' efforts and detailed responses. However, I still have concerns regarding the scale of the dataset. While the authors mention that recent sEEG datasets contain fewer than 12 subjects, this does not sufficiently justify that 23 subjects are adequate. In related fields, datasets often include 50, 100, or even more subjects to ensure robustness and generalizability.

We echo that a larger dataset (50-100 subjects) would provide stronger evidence for the model's generalization capabilities. However, as reviewer o6UF mentioned "gathering this 50-100 sEEG subjects could take years, if not decades, not to mention the additional difficulty of recruiting them for behavioral studies". Importantly, if the scientific community found value in studies containing 12 subjects or less, it is likely that it will find value in our 21-subject study as well.

I acknowledge the challenges of data collection from a medical perspective. However, there are still larger, publicly available sEEG datasets.

This is very insightful. If you could point us towards those datasets, by specifying either the publication with which they are associated or the link to the dataset repository, we would be very thankful and we look forward to working with those datasets in the future.

If collecting sEEG data in epilepsy is too difficult, it might be worth considering using other types of brain signals, which could be more suitable for applying machine learning techniques.

Thank you for your suggestion. While there are certainly many different types of brain signals that can be analyzed using machine learning techniques, possibly with larger datasets, we believe that there is a lot of value in applying machine learning techniques to sEEG datasets because:

1. sEEG is currently the gold standard invasive neural recording modality used in humans. Therefore, improving neural decoding based on sEEG bring it closer to clinical translation compared to other neural recording modalities which might provide a lot of data from animal models but have very rarely been used in humans due to safety or other concerns (such as microelectrode array recordings for example). Edit: We would also like to emphasize that even if another neural recording modality was used, the only clinical population approved for collecting intracranial neural recordings would be epilepsy patients. Therefore, the challenges associated with collecting a large dataset would carry over to any other modality as well.
2. there is a lot of scientific value in showing that machine learning tools can be used on smaller datasets, as in numerous fields (healthcare & medicine, astronomy, environmental science, to name a few) collecting large datasets can be extremely costly and time-consuming. Our work is especially valuable in that sense, since we show that combining many small datasets, despite data heterogeneity, can lead to better machine learning models, compared to training on the smaller datasets individually.

Thanks for engaging in the discussion! We very much appreaciate your insights and truly believe that our work has become stonger based on your feedback and suggestions.

Thank you for your thoughtful feedback and suggestions. We are excited to read that our "method has a value in itself and provides ideas for future research"! Here, we provide our response to your comments/questions and the results from new experiments that we ran to address your concerns. We hope that our answers will clear any lingering concerns. We are also available to provide further clarification, if needed.

Weaknesses

1. No comparison with other methods during within-subject training was done. Authors only provide the results for their own architecture for within-subject experiments. Therefore, they only show that their multi-subject trained model is stronger compared to their single-subject model. However, it might be possible that some of the other existing single-subject State-of-the-Art models will be more effective than authors’ multi-subject model on authors’ data. The importance of these comparisons are magnified by the fact that authors’ dataset will stay private (therefore only the authors of this paper can provide metrics for other models applied to this dataset); and by the reasonably large size of authors’ model (it might be too big to be efficiently trained on a single person data, while other models, which were developed for single-subject tasks, might be lighter and train better on small datasets).

We thank the reviewer for this constructive criticism. To address this, we trained classic and state-of-the-art single-subject models to compare against our model. The results are summarized in Table 1 of our General Response (see also Figure 1 of the attached pdf) and suggest that our architecture outperforms the baselines when trained on either single-subject or multi-subject data. Specifically, our single-subject models outperformed all baselines by $\Delta R^{2} >= 0.03$, while our multi-subject models outperformed all baseline by $\Delta R^{2} >= 0.12$. Importantly, our transfer learned single-subject models also beat the baselines by $\Delta R^{2} >= 0.10$. We hope that these comparisons are sufficient to address your concerns. We plan to include this analysis in section 4 of our revised manuscript.

2. Decoding performance with and without spatial positional encoding seems very similar. For such a difference it is interesting to see if it has a statistical significance.

To identify whether the difference in the decoding performance with and without spatial positional encoding is statistically significant, we performed a Wilcoxon rank-sum test between the groups: (1) $R^{2}$ of all subjects obtained by training the multi-subject model with spatial positional encoding, and (2) $R^{2}$ of all subjects obtained by training the multi-subject model without spatial positional encoding. The test returned a test-statistic = 0.34 and a p-value = 0.73, indicating that there is no significant difference between the performance of the model with and without spatial positional encoding. We plan to include this result in section 4.5 of our revised manuscript.

Questions

1. How many trainable and frozen parameters there were during training and finetuning of the model?

The total number of trainable parameters for the multi-subject model is 797,095. From those parameters 753,394 are shared across subjects and the rest 43,701 parameters are subject-specific (2,081 parameters per subject). When training single-subject models (section 4.2) and multi-session, multi-subject models (section 4.3), all parameters of the model (shared + subject specific) were trained. When transferring the pretrained, multi-subject models to new subjects (section 4.4), for the first 400 epochs the shared parameters were kept frozen (753,394 parameters) while the subject-specific parameters were being trained (2,081 parameters per subject). For the remaining 600 epochs, all parameters were unfrozen and trained. We plan to include those numbers in section A.2 of our revised manuscript.

2. How large is the variance in the response time within each subject?

For convenience, the mean of the variance of the response times across all 21 subjects of our dataset is $\sigma^{2} = 0.011 \pm 0.0017$ (mean $\pm$ sem) sec$^{2}$. The variance of the response times for each subject is: (18, 9, 5, 6, 7, 3, 4, 5, 17, 25, 4, 9, 9, 13, 21, 19, 8, 31, 4, 11, 5) $\times 10^{-3}$ sec$^{2}$.

3. Are all K convolution kernels in Convolution tokenizer the same for all electrodes?

Yes. The K convolutional kernels that our tokenizer uses are the same for all electrodes of all subjects. This design choice increases the effective number of training samples available to the tokenizer, making it less prone to overfitting and increasing the model's robustness. We will make sure to ephasize this in section 3.2.1 of our revised manuscript.

4. How much recording hours of data was used overall and on average per subject (only the electrode-hours are provided)?

The total recording hours of data used during model training were 1.54 hrs. The average per-subject recording time that was used for model training was 4.39 $\pm$ 0.39 (mean $\pm$ sem) mins. We will make sure to add those is section A.2 of our manuscript.

5. What is the statistical significance of the metrics increase while using positional encoding compared to not using it?

The increase in the performance due to spatial positional encoding compared to not using it is not statistically significant. Please refer to our response to your comment under Weaknesses (bullet point No. 2) for a detailed description of how we obtained this result.