事件相关电位的呈现是在特定频带下的时域波形，从这个角度出发的话，其处理就涉及到滤波；此外由于脑电信号易受噪声干扰，所以单独测出的每两次的数据肯定是不同的，甚至会差别很大，即使是在同一条件下，这就需要去除噪声（这里不做详细介绍）、去除基线和叠加平均。
简单的处理流程：

1. 滤波
2. 分段&加标签
3. 去基线
4. 叠加平均

滤波（filt raw）： 滤波的作用有两个，第一，取出我们关注的ERP特征最明显的频段的数据，比如8~13Hz；第二，滤除低频漂移（线性趋势），有的时候观察原始波形，发现会有一个线性上升或者下降，经过低频滤波便可去除，比如过滤掉3Hz以下的数据。
一个典型的MATLAB滤波函数如下：

% 滤波iir
function r=filtraw(raw,fs,lowcut,highcut,order)
sizeraw=size(raw);
udfilter=designfilt('bandpassiir','FilterOrder',order, ...
    'HalfPowerFrequency1',lowcut,'HalfPowerFrequency2',highcut, ...
    'SampleRate',fs);
for i=1:sizeraw(1)
    raw(i,:)=filter(udfilter,raw(i,:));
end
r=raw;
end
1
2
3
4
5
6
7
8
9
10
11

分段&加标签（epoch & label）： 我们一般获取到的数据是连续多个任务的，还涉及到任务之间短暂休息时候的信号，这些信号是没用的，分段&加标签的作用就是把有用的时间段的数据取出来。由于一般我们采用多通道的电极采集信号，所以得到的原始数据一般是二维的，比如得到的数据是16x120000个，可以认为是16通道每个通道120000个数据点，如果采样率是1200，那就是100s的数据。假设任务周期是10s，关注的是每个10s中的第4~7s共3s的数据，那么是不是16x3600x10的数据结构比较好了，如果这10个任务中，有一半是任务1，另外一半是任务2，就需要一个标签文件标注第几个是任务几了（一般采集数据的时候就完成了）；也可以将数据分为两个16x3600x5的，这样就不需要标签文件了，只需要备注好数据名就好。
（由于每个人打标签的方法可能都不同，这里就不放我的代码了，我的代码很乱）
去基线：（remove baseline） 去基线的原因是，我们获取到的脑电信号，有正有负，但是我们不能说这个就是脑电信号的绝对值，脑电信号获取到的都是相对的电压值（所以才会有重参考的概念）；或者可以认为去基线就是为了在画图的时候，任务开始的0时刻之前的数据靠近或者为0，这样才会有任务期间时域波形的相对变化。我们一般取任务开始0时刻之前的100ms~200ms作为基线，每一次任务期间，基线的值可能都不同，所以不去除基线的话，下一步的叠加平均是无效的，因为没有基准，比如有的基线为负，有的基线为正，这样的情况下执行叠加平均，在基线抵消的时候，可能任务期间的波形也会抵消，也就获取不到ERP了。
去基线的具体方式就是，对于单次任务的脑电信号，对基线期间的脑电信号进行叠平均，然后用所有的脑电信号，减去这个平均值，得到去基线之后的信号，新的信号在基线期间是靠近0线的。

% 去除基线
function r=rbase(eeg,fs,tbase)
    sizeeeg=size(eeg);
    base=zeros(sizeeeg(1),1);
    nbase=fs*tbase;
    for i=1:sizeeeg(3)
        for j=1:nbase
            base=base+eeg(:,j,i);
        end
        base=base/nbase;
        eeg(:,:,i)=eeg(:,:,i)-base;
    end
    r=eeg;
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14

叠加平均：（average） 这个就不解释了
代码：

% 叠加平均
function r=average(eeg)
sizeeeg=size(eeg);
sizetemp=size(sizeeeg);
if sizetemp(2)==3
    temp=zeros(sizeeeg(1),sizeeeg(2),1);
    for i=1:sizeeeg(3)
        temp=temp+eeg(:,:,i);
    end
    r=temp/sizeeeg(3);
elseif sizetemp(2)==2
    temp=zeros(sizeeeg(1),1);
    for i=1:sizeeeg(2)
        temp=temp+eeg(:,i);
    end
    r=temp/sizeeeg(2);
end
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

如果你的信号比较好，经过简单的四步，就能得到自己想要的ERP了，如果不行，就需要添加一些去除噪声的步骤了，所以最好的办法就是一开始就获取到好的信号。
题图来源：维基百科
https://en.wikipedia.org/wiki/Event-related_potential