Modelling data from EEG electrodes

Hi there,
I am trying to model EEG data and I am struggling in conceptualising/implementing a model that accounts for the recording of multiple channel within participants. I’ll explain better. I have a dataset containing one outcome variable (Lempel-Ziv Complexity - LZC - centered by subtracting the mean as LZC has values only between 0 and 1) and one categorical predictor (task - 9 levels) - I need to estimate the effects of task on LZC. There are 20 participants for whom I have 63 LZC value, one for each EEG channel. Therefore, in total the data contains:
20 participants X 63 electrodes X 9 tasks

I started modelling the data as:

lzc ~ 1 + task + (1|participant) + (1|channel) + (1|participant:channel)

to account for the crossed-nested relationship between electrodes and participants. I think this model should account for the fact that the electrodes are the same across participants, but there might be idiosyncrasies between electrodes and participants. The model runs mostly fine at the moment (just need to fix the Rhat for the Intercept), but I’m suspicious of the results. The posterior distributions are correctly centered to the experimental mean LZC computed for each task. However, the distributions seem to state that there are small but highly probable effects of tasks on LZC. I am worried that this model does not account for the fact that EEG electrodes are correlated within participants.

Family: gaussian 
  Links: mu = identity 
Formula: lzc ~ task + (1 | participant_id) + (1 | channel) + (1 | participant_id:channel) 
   Data: lzc_data (Number of observations: 16380) 
  Draws: 4 chains, each with iter = 4000; warmup = 1000; thin = 1;
         total post-warmup draws = 12000

Multilevel Hyperparameters:
~channel (Number of levels: 63) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     0.01      0.00     0.01     0.01 1.01     1030     1983

~participant_id (Number of levels: 20) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     0.02      0.00     0.01     0.02 1.00     1214     2500

~participant_id:channel (Number of levels: 1260) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     0.01      0.00     0.01     0.01 1.00     1439     3107

Regression Coefficients:
          Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept     0.16      0.00     0.15     0.17 1.01      545     1064
taskldt      -0.00      0.00    -0.00    -0.00 1.00    12163     8843
taskmix      -0.00      0.00    -0.00    -0.00 1.00    11848     8563
taskmr       -0.00      0.00    -0.01    -0.00 1.00    11795     9527
tasknbk_0    -0.00      0.00    -0.00    -0.00 1.00    11984     8373
tasknbk_1    -0.00      0.00    -0.00    -0.00 1.00    11884     8962
tasknbk_2    -0.00      0.00    -0.00    -0.00 1.00    12182     9306
taskrs_1     -0.00      0.00    -0.00    -0.00 1.00    11629     8898
taskrs_2      0.00      0.00     0.00     0.00 1.00    12091     9111

Further Distributional Parameters:
      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma     0.01      0.00     0.01     0.01 1.00    14510     8674

Does anyone know if there are better ways to model EEG data without aggregating across channels? In my case, I do not believe I can average LZC across channels.

Thank you!

Sample of the data

lzc_data.csv (378.2 KB)

1 Like

I’m following now because I’m interested if anyone else can chip in about the correlations in the group-level effects, but I do have a couple of questions/suggestions on building up the model a little bit:

  1. Is the lzc variable ever exactly 0 or 1? If not, instead of subtracting the mean value you could model the raw values with a Beta regression. If I can pass untransformed values to a model I like to do that, and I’ve seen lots of other people who know loads more than me also recommend that.

  2. Write now, your model has no varying slopes, so it assumes that the effect of task on the outcome variable is identical for each participant and channel. This is typically a bad assumption when it comes to humans, so adding a varying slope of task to the random effects will likely widen up your posteriors a bit as the model will have more uncertainty.

1 Like
  • The beta suggestion is a good idea and it was a rushed decision on my side. I actually realized it just after posting about this yesterday, so it’s good to have confirmation. lzc in the dataset is always <1, but theoretically it could go up to just above above it because of the way it is computed (although this is an unlikely edge case). I think a gamma/lognormal distribution might be a good choice in this situation. Good to know that untransformed data is the general preference, I’ll keep the suggestion in mind in the future!
  • I agree on that. My colleagues and I wanted to sort out the electrodes before adding random slopes, so to add that layer of complexity once we have a good backbone for the model.
1 Like

I’ve been lurking on this because, like @sjp, I wanted to see if anyone more knowledgeable than myself would help with the correlation question. To bump this, if nothing else, here is my thinking process:

if you are explicitly interested in these correlations, you could fit a multivariate model:

mvbind(ch1,...,ch63) ~ 1 + task + (1|p|participand_id)

this would give you correlations among the different channels at the expense of not calculating the lzc for an average channel. But I guess that would be a dealbreaker.

So I wonder whether specifying a model like this:

lzc ~ 1 + task + (1|participant_id) + (1|channel) + (0 + channel|participant_id)

would be a legitimate way to account for these correlations.

I did fit this model to a subset of the data you provided; it looks like it gets the job done. But I am not 100% confident that this specification makes sense.

What i would recommend is to think about the generative model for lzc that include the interelectrode variability and intersubject variability and make synthetic data sets that reflect different reasonable cases. Then you can play with different forms of the model on the synthetic data without concern about double dipping. You can use LOOIC or WAIC to select the model form that best captures the known generative model. Once you feel happy that your model is well validated with synthetic data, you can fit your real data.

Is the expectation that the correlation across electrodes is scaled spatially?

I found some potentially relevant work: