TRCA：提升SSVEP-BCI性能的关键技术解析与Matlab实战

foamz

1023人浏览 · 2025-03-07 08:55:50

foamz · 2025-03-07 08:55:50 发布

文章目录

前言
一、TRCA是什么？
二、TRCA方法讲解
三、matlab实战
总结
参考文献

前言

传统的典型相关分析（CCA）方法在SSVEP信号识别中取得了一定成效，但容易受到自发脑电活动的干扰，且未充分利用相位信息。为此，研究者提出了TRCA方法，通过最大化任务期间神经影像数据的可重复性，提取任务相关成分，增强信号的信噪比（SNR），抑制背景脑电噪声。任务相关成分分析（Task-Related Component Analysis, TRCA）是一种用于脑电信号（EEG）处理的多变量分析方法，主要应用于脑机接口（BCI）领域。其核心目标是从多通道脑电信号中提取与特定任务（如视觉刺激、运动想象等）高度相关的神经活动成分，同时抑制噪声和非任务相关的信号。
在这里插入图片描述

一、TRCA是什么？

在SSVEP-BCI系统中，任务相关成分分析（Task-Related Component Analysis，TRCA）被用于提高SSVEP信号的检测性能。通过增强任务相关成分，TRCA能够提高信号的信噪比，抑制背景脑电活动的干扰，从而提升系统的准确性和稳定性。TRCA 基于以下假设：“同一任务下重复多次的脑电信号中，任务相关成分在时间维度上具有高度一致性。”通过最大化不同试验（trial）之间的协方差，TRCA 可以提取出这些重复出现的稳定成分，而噪声和非任务相关成分因随机性会被抑制。

二、TRCA方法讲解

1.EEG信号的线性表示

TRCA（任务相关成分分析）是一种通过最大化任务期间的可复现性来高效提取任务相关成分的方法。假设观测到的多通道 EEG 信号X(t)由以下线性生成模型表示：
$x_j(t) = a_{1,j} s(t) + a_{2,j} n(t), \quad j = 1,2,\dots,N_c$
其中， $j$ 是通道索引， $a_{1,j}$ 和 $a_{2,j}$ 是混合系数，用于将源信号投影到 EEG 观测信号中。

这里 $S (t)$ 是任务相关信号，n(t)是任务无关信号。

(1)S(t)任务相关信号

$S^{(1)}、S^{(2)}、S^{(3)}$ 之间的协方差是一个正常数，这说明它们在统计意义上具有稳定的线性关系，即它们的相关性在不同时间段或试验条件下保持不变。以下公式可以表示这个关系。
$\text{Cov}(s^{(1)}(t), s^{(2)}(t)) = \text{Cov}(s^{(1)}(t), s^{(3)}(t)) = \text{Cov}(s^{(2)}(t), s^{(3)}(t)) = \text{常数}$

在这里插入图片描述

(2)n(t)任务无关信号

任务无关成分 $n^{(1)}t$ 、 $n^{(2)}t$ 、 $n^{(3)}t$ 之间的协方差为0，并且任务相关信号 $S (t)$ 与任意的 $n (t)$ 之间的协方差也均为0。可以用以下公式表示这个关系

$\text{Cov}(s^{(i)}(t), n^{(j)}(t)) = \text{Cov}(n^{(i)}(t), n^{(j)}(t)) = 0, \quad 1 \leq i, j \leq 3$

2.TRCA 任务相关成分提取方法

(1)问题定义：从观测信号中恢复任务相关成分

对多通道的EEG信号进行加权求和，得到线性的一维信号 $y (t)$ ，表示如下
$\sum_{i=1}^{N} w_i x_i(t) = w^T x(t)$
或者拆解更细致一些，得到
$\sum_{j=1}^{N_c} w_j x_j(t) = \sum_{j=1}^{N_c} \left( w_j a_{1,j} s(t) + w_j a_{2,j} n(t) \right)$

在这里插入图片描述

对于 $y (t)$ ，想从 $x (t)$ 的线性组合中，恢复出任务相关成分 $S (t)$ ，需要将 $n (t)$ 信号完全消除掉，然后 $S (t)$ 信号完全保留下来。所以在理想情况下，该问题的解只需要要满足两个条件：

① 任务相关信号 $S (t)$ 的系数为1，公式表示得到 $\sum_{j=1}^{N_c} w_j a_{1,j} = 1$

② 任务无关信号 $n (t)$ 的系数为0，公式表示得到 $\sum_{j=1}^{N_c} w_j a_{2,j} = 0$

只要满足这两个条件，就可以得到最终解： $y (t) = s (t)$

(2)求解方法：跨试次协方差最大化

该问题可以通过跨试次协方差最大化（inter-trial covariance maximization）来求解。
$C_{h_1 h_2} = \text{Cov} \left( y^{(h_1)}(t), y^{(h_2)}(t) \right)$

第 $h$ 次试验的 EEG 信号和估计的任务相关成分分别表示为 $x^{(h)}(t)$ 和 $y^{(h)}(t)$ , 其中 $h =1,2.....N_t$
这里的 $y^{(h)}(t)$ 的时间区间固定为 $t∈[t_h,t_h+T]$ ，此处 𝑇 是每个任务试次的持续时间。

在这里插入图片描述

$y (t)$ 在试次 $h_1$ 和 $h_2$ 之间的协方差，因为 $y (t)$ 是由 $x (t)$ 线性表示的，所以也可以把 $y (t)$ 拆解成用 $x^{(h_1)}_{j_1}(t)$ 来表示，所以该协方差求解可以表示如下：
$C_{h_1 h_2} = \text{Cov} \left( y^{(h_1)}(t), y^{(h_2)}(t) \right) = \sum_{j_1,j_2=1}^{N_c} w_{j_1} w_{j_2} \text{Cov} \left( x^{(h_1)}_{j_1}(t), x^{(h_2)}_{j_2}(t) \right)$

跨试次协方差矩阵表示

上面那个公式只是单一情况，现在对所有情况求和，即可得到跨试次表示
所有可能的试次组合被总结如下：
$\sum_{\substack{h_1, h_2 = 1 \\ h_1 \neq h_2}}^{N_t} C_{h_1 h_2} = \sum_{\substack{h_1, h_2 = 1 \\ h_1 \neq h_2}}^{N_t} \sum_{j_1,j_2=1}^{N_c} w_{j_1} w_{j_2} \text{Cov} \left( x^{(h_1)}_{j_1}(t), x^{(h_2)}_{j_2}(t) \right) = w^T S w$

这里，矩阵 $(S_{j_1 j_2})_{1 \leq j_1, j_2 \leq N_c}$ 会被定义为 $S_{j_1 j_2} = \sum_{\substack{h_1, h_2 = 1 \\ h_1 \neq h_2}}^{N_t} \text{Cov} \left( x^{(h_1)}_{j_1}(t), x^{(h_2)}_{j_2}(t) \right)$

(3)优化约束条件

为了获得有限解，在使用CorrMax之前，需要先对任务相关成分 $y (t)$ 进行标准化处理，标准化之后 $y (t)$ 的方差需要被归一化为1 。
所以用公式表示 $y (t)$ 的方差被约束为：
$\text{Var}(y(t)) = \sum_{j_1,j_2=1}^{N_c} w_{j_1} w_{j_2} \text{Cov} (x_{j_1}(t), x_{j_2}(t)) = w^T Q w = 1$

(4)求解空间滤波器（理解到此即可求解）

该约束优化问题可以求解为：

$\hat{w} = \arg \max_w \frac{w^T S w}{w^T Q w}$

最优的系数向量可以通过矩阵 $Q^{-1} S$ 的特征向量获得。只要求出这个最大特征值对应的特征向量，即可得到最优解。

(5)TRCA 在 SSVEP-BCI 中的应用

一般来说，从一个 $N_c \times N_c$ 矩阵中可以得到 $N_c$ 个特征值和特征向量。矩阵 $Q^{-1} S$ 的特征值 $\lambda$ 可以按降序排列，它表示相应特征向量 $\hat{w}$ 的成本函数值。因此，这些特征值反映了多个试次之间的任务一致性。
在这里插入图片描述

图 (a) 显示了基于 TRCA 方法的流程图。在基于 SSVEP 的脑机接口（BCI）中，TRCA 可用于设计空间滤波器，以去除头皮记录中的背景 EEG 活动。
针对第 $n$ 个刺激的空间滤波器， $w_n^{(m)} \in \mathbb{R}^{N_c}$ 通过 TRCA 从个体校准数据 $\chi_n^{(m)}$ 中获取，并在应用滤波器组分析后计算得到。在这里，公式 $\text{Var}(y(t)) = w^T Q w = 1$ 中的 $Q$ 矩阵是使用所有训练试次的拼接矩阵 $\chi_n^{(m)}$ 计算得到的。

(6)单试次与平均训练数据的相关性计算

单试次测试数据 $X^{(m)} \in \mathbb{R}^{N_c \times N_s}$ 与跨试次平均的第 $n$ 个视觉刺激训练数据 $\bar{\chi}_n^{(m)} \in \mathbb{R}^{N_c \times N_s}$ 之间的相关系数计算如下：
$r_n^{(m)} = \rho \left( \left( X^{(m)} \right)^T w_n^{(m)}, \bar{\chi}_n^{(m)} w_n^{(m)} \right)$
其中， $p (a, b)$ 表示信号 $a$ 和 $b$ 之间的皮尔逊相关分析。

三、matlab实战

1.数据集

使用清华大学脑机接口团队的公开数据集，benchmark dataset
下载地址： https://bci.med.tsinghua.edu.cn/
关于该数据集的介绍和使用，可以看这篇文章：
清华大学脑机接口研究组Benchmark dataset数据集介绍
在这里插入图片描述

2. test_trca.m

trca的测试代码如下，主要是读取数据文件，然后对数据滤波处理后，送入trca算法。

tic
% 加载所需的数据文件
load('Freq_Phase.mat');  % 载入频率和相位数据
data = load('S1.mat');    % 载入EEG数据
eeg_data = data.data;    % 获取EEG数据矩阵
% 处理数据
stim_data = eeg_data([48, 54, 55, 56, 57, 58, 61, 62, 63], 251:1250, :, :);
low_freq = 5;   % 下边界频率
high_freq = 90; % 上边界频率
fs = 250;       % 采样频率为250Hz
% 获取数据的维度
[N_channels, N_timepoints, N_targets, N_blocks] = size(stim_data);
% 将数据重组为一个 2D 数组
stim_data_reshaped = reshape(stim_data, N_channels, N_timepoints * N_targets * N_blocks);
% 使用 parfor 进行并行处理
parfor channel_idx = 1:N_channels
    stim_data_reshaped(channel_idx, :) = bandpass(stim_data_reshaped(channel_idx, :), [low_freq, high_freq], fs);
end
% 将数据恢复到原来的四维形式
stim_data = reshape(stim_data_reshaped, N_channels, N_timepoints, N_targets, N_blocks);
eeg = stim_data;
% 留一交叉验证（LOO Cross-Validation）
predicted_labels_all = [];  % 存储所有预测标签
true_labels_all = [];  % 存储所有真实标签
trial_info_all = {};  % 存储每个试次的文字说明
inconsistent_flags = [];  % 标记与真实标签不一致的数据
for loocv_i = 1:N_blocks
    % 提取测试数据
    Testdata = squeeze(eeg(:, :, :, loocv_i));
    Traindata = eeg;
    Traindata(:, :, :, loocv_i) = [];
    
    % 对每个目标频率，计算训练数据的平均值
    for targ_i = 1:N_targets
        aver_Traindata(:, :, targ_i) = squeeze(mean(squeeze(Traindata(:,:,targ_i,:)),3));  % 对每个目标频率进行平均
    end
    
    % 标签赋值
    truelabels = freqs;
    
    % 获取测试数据的试次数量
    N_testTrial = size(Testdata, 3);
    
    % 对每个测试试次进行处理
    outputlabels = zeros(1, N_testTrial);  % 存储预测的标签
    for trial_i = 1:N_testTrial
        coefficience = zeros(1, length(truelabels));
        for targ_j = 1:length(freqs)              
            % 使用训练数据计算空间滤波器 wn
            wn = TRCA(squeeze(Traindata(:,:,targ_j,:)));  
            
            % 计算加权训练数据与测试数据的相关性
            weighted_train = wn' * aver_Traindata(:,:,targ_j);  
            weighted_test = wn' * Testdata(:,:,trial_i);  
            
            % 计算相关系数
            coefficienceMatrix = corrcoef(weighted_test, weighted_train);
            coefficience(targ_j) = coefficienceMatrix(1, 2);
        end
        % 目标检测：选择最大相关系数的频率作为输出标签
        [~, index] = max(coefficience);
        outputlabels(trial_i) = freqs(index);  % 输出预测的标签
    end
    % 计算该次测试的正确分类数量
    trueNum = sum((outputlabels - truelabels) == 0);  % 统计预测标签与真实标签相同的个数
    acc(loocv_i) = trueNum / length(truelabels);  % 计算该次交叉验证的准确率
    
    % 打印当前交叉验证的准确率
    fprintf('第%d试次的交叉验证准确率是: %.4f, 样本数: %d/%d\n', loocv_i, acc(loocv_i), trueNum, N_testTrial)
    
    % 存储预测标签、真实标签、试次信息和一致性标志
    trial_info = repmat({['试次 ' num2str(loocv_i)]}, N_testTrial, 1);  % 每个试次的说明
    predicted_labels_all = [predicted_labels_all; num2cell(outputlabels')];  % 存储预测标签
    true_labels_all = [true_labels_all; num2cell(truelabels')];  % 存储真实标签
    
    % 修改一致性标记：若预测与真实标签一致，则标记为 TRUE；否则标记为 FALSE
    inconsistent_flags = [inconsistent_flags; num2cell(outputlabels' == truelabels')];  % 判断一致的地方
    trial_info_all = [trial_info_all; trial_info];  % 存储试次信息
end
% 将数据存入表格
T = table(trial_info_all, predicted_labels_all(:,1), true_labels_all(:,1), inconsistent_flags(:,1), 'VariableNames', {'试次', '预测标签', '真实标签', '一致性'});
% 输出到 Excel 文件
writetable(T, 'predicted_labels_with_consistency.xlsx');
% 计算总的运行时间
t = toc;
% 打印总的运行时间和平均准确率
fprintf('\n-----------------------------------------\n')
disp(['总运行时间: ', num2str(t), ' s']);
fprintf('6试次数据交叉验证的平均准确率是: %.4f\n', mean(acc));

4.运行结果

在这里插入图片描述

总结

TRCA作为一种有效的任务相关信号提取方法，在SSVEP信号检测中发挥了重要作用。通过最大化任务相关成分的可重复性，TRCA提高了信号的信噪比，增强了SSVEP信号的检测性能。与其他方法相比，TRCA在处理背景脑电噪声方面表现出色，为SSVEP-BCI系统的优化提供了新的思路。

参考文献

M. Nakanishi, Y. Wang, X. Chen, Y. -T. Wang, X. Gao and T. -P. Jung, “Enhancing Detection of SSVEPs for a High-Speed Brain Speller Using Task-Related Component Analysis,” in IEEE Transactions on Biomedical Engineering, vol. 65, no. 1, pp. 104-112, Jan. 2018, doi: 10.1109/TBME.2017.2694818.

Tanaka H, Katura T, Sato H. Task-related component analysis for functional neuroimaging and application to near-infrared spectroscopy data. Neuroimage. 2013 Jan 1;64:308-27. doi: 10.1016/j.neuroimage.2012.08.044. Epub 2012 Aug 24. PMID: 22922468.