近红外光谱技术可通过检测大脑皮层中的氧合血红蛋白和还原血红蛋白浓度变化信息，提供脑功能活动过程中的大脑皮层血氧代谢信息，从而用于脑功能活动状态检测。与功能性磁共振成像、正电子放射断层扫描、脑电图检测等其它脑功能检测方法相比，近红外光谱技术具有经济、安全、非侵入、使用方便、易实施等优点。但是，在采用近红外光谱技术对脑功能活动进行检测时，人体的心跳、呼吸、人体低频振荡及人体超低频振荡等生理活动会对测量信号产生干扰，也称之为生理干扰。这种生理干扰现象不仅出现在头皮、颅骨和脑脊髓液等外部脑组织中，也同时出现在脑灰质和脑白质等深层脑组织之中，并对脑功能活动检测信号的提取产生严重影响。Near-infrared spectroscopy can provide information on blood oxygen metabolism in the cerebral cortex during brain functional activities by detecting changes in the concentration of oxyhemoglobin and reduced hemoglobin in the cerebral cortex, so as to detect the state of brain functional activity. Compared with other brain function detection methods such as functional magnetic resonance imaging, positron emission tomography, and electroencephalography, near-infrared spectroscopy has the advantages of being economical, safe, non-invasive, convenient to use, and easy to implement. However, when using near-infrared spectroscopy to detect brain function activities, physiological activities such as heartbeat, breathing, low-frequency oscillations and ultra-low-frequency oscillations of the human body will interfere with the measurement signal, which is also called physiological interference. This physiological interference phenomenon not only appears in external brain tissues such as scalp, skull and cerebrospinal fluid, but also in deep brain tissues such as gray matter and white matter, and has a serious impact on the extraction of brain functional activity detection signals.

对于存在生理干扰情况下的近红外脑功能检测，通常是可以采用基于多距测量方法的自适应滤波技术实现对脑功能活动信号的提取。但是，当人体的心跳、呼吸、低频振荡或超低频振荡等生理干扰与脑功能信号的频带严重重叠时，经过该方法处理得到的脑功能活动信号通常存在误差干扰，导致近红外脑功能活动信号的检测精度低，使得该方法在实际应用中受到一定限制。For near-infrared brain function detection in the presence of physiological interference, it is usually possible to use adaptive filtering technology based on multi-distance measurement methods to extract brain function activity signals. However, when the human body's heartbeat, respiration, low-frequency oscillation or ultra-low-frequency oscillation and other physiological interferences seriously overlap with the frequency bands of brain function signals, the brain function activity signals processed by this method usually have error interference, resulting in near-infrared brain function activity signals. The detection accuracy is low, which limits the practical application of this method.

发明内容Contents of the invention

本发明目的是为了解决当生理干扰与脑功能信号频带严重重叠时，采用基于多距测量方法的自适应滤波技术得到的脑功能活动信号存在误差干扰，导致近红外脑功能活动信号的检测精度低的问题，而提出了基于递归最小二乘自适应滤波及最小二乘支持向量机的近红外脑功能信号提取方法。The purpose of the present invention is to solve the problem of error interference in the brain function activity signal obtained by adopting the adaptive filtering technology based on the multi-distance measurement method when the physiological interference and the brain function signal frequency band overlap seriously, resulting in low detection accuracy of the near-infrared brain function activity signal A near-infrared brain function signal extraction method based on recursive least squares adaptive filtering and least squares support vector machine is proposed.

上述的发明目的是通过以下技术方案实现的：Above-mentioned purpose of the invention is achieved through the following technical solutions:

步骤一：在待测脑组织头皮表面放置一个由双波长光源S与检测器D₁和检测器D₂所构成的近红外探头，双波长光源S与检测器D₁之间的直线距离为R₁，双波长光源S与检测器D₂之间的直线距离为R₂，双波长光源S发出两种近红外光的波长分别为λ₁和λ₂，检测器D₁和检测器D₂用于获取大脑安静状态下的漫反射光强和大脑诱发激励状态下的漫反射光强，从而获得两个不同波长的近红外光在不同距离下的光密度变化量的时间信号：和Step 1: Place a near-infrared probe consisting of a dual-wavelength light source S, a detector D₁ and a detector D₂ on the surface of the scalp of the brain tissue to be tested. The linear distance between the dual-wavelength light source S and the detector D₁ is R_1. The linear distance between the dual-wavelength light source S and the detector D₂ is R₂ , the wavelengths of the two near-infrared lights emitted by the dual-wavelength light source S are λ₁ and λ₂ respectively, and the detector D₁ and the detector D₂ are used To obtain the diffuse reflection light intensity in the quiet state of the brain and the diffuse reflection light intensity in the brain-induced excitation state, so as to obtain the time signal of the optical density variation of two different wavelengths of near-infrared light at different distances: and

其中，t为采样时刻，t＝1，2，…，N，N为正整数；(此处表示t的取值范围是从1到N)；Wherein, t is the sampling moment, t=1, 2, ..., N, N is a positive integer; (the value range of representing t here is from 1 to N);

为双波长光源S发出近红外光的波长为λ₁时，与双波长光源S的直线距离为R₁的检测器D₁获得的光密度变化量的时间信号； When the wavelength of the near-infrared light emitted by the dual-wavelength light source S is λ₁ , the linear distance from the dual-wavelength light source S is the time signal of the optical density variation obtained by the detector D₁ of R₁ ;

为双波长光源S发出近红外光的波长为λ₂时，与双波长光源S的直线距离为R₁的检测器D₁获得的光密度变化量的时间信号； When the wavelength of the near-infrared light emitted by the dual-wavelength light source S is λ₂ , the linear distance from the dual-wavelength light source S is the time signal of the optical density variation obtained by the detector D₁ of R₁ ;

为双波长光源S发出近红外光的波长为λ₁时，与双波长光源S的直线距离为R₂的检测器D₂获得的光密度变化量的时间信号； When the wavelength of the near-infrared light emitted by the dual-wavelength light source S is λ₁ , the linear distance from the dual-wavelength light source S is the time signal of the optical density variation obtained by the detector D₂ of R₂ ;

为双波长光源S发出近红外光的波长为λ₂时，与双波长光源S的直线距离为R₂的检测器D₂获得的光密度变化量的时间信号； When the wavelength of the near-infrared light emitted by the dual-wavelength light source S is λ₂ , the linear distance from the dual-wavelength light source S is the time signal of the optical density variation obtained by the detector D₂ of R₂ ;

步骤二：对步骤一获得的两个不同波长的近红外光在不同距离下的光密度变化量时间信号进行修正郎伯比尔定律计算，从而分别得到检测器D₁和检测器D₂处的氧合血红蛋白浓度变化时间信号和还原血红蛋白浓度变化时间信号：和Step 2: Calculate the time signal of the optical density variation of two different wavelengths of near-infrared light at different distances obtained in step 1 by amending Lambert-Beer's law to obtain the oxygen at detector D₁ and detector D₂ respectively. Combined hemoglobin concentration change time signal and reduced hemoglobin concentration change time signal: and

t为采样时刻，t＝1，2，…，N，N为正整数；t is the sampling moment, t=1, 2, ..., N, N is a positive integer;

为与双波长光源S的直线距离为R₁的检测器D₁测量得到的氧合血红蛋白浓度变化量的时间信号； is the time signal of the change in the concentration of oxyhemoglobin measured by the detector D1 whose straight_- line distance from the dual_- wavelength light source S is R1;

为与双波长光源S的直线距离为R₁的检测器D₁测量得到的还原血红蛋白浓度变化量的时间信号； is the time signal of the reduced hemoglobin concentration change measured by the detector D1 whose straight_- line distance from the dual_- wavelength light source S is R1;

为与双波长光源S的直线距离为R₂的检测器D₂测量得到的氧合血红蛋白浓度变化量的时间信号； is the time signal of the change in the concentration of oxyhemoglobin measured by the detector D2 whose straight_- line distance from the_dual -wavelength light source S is R2;

为与双波长光源S的直线距离为R₂的检测器D₂测量得到的还原血红蛋白浓度变化量的时间信号； is the time signal of the amount of change in the concentration of reduced hemoglobin measured by the detector D₂ whose straight-line distance from the dual-wavelength light source S is R₂ ;

步骤三：利用步骤二得到的和或者和构建自适应滤波函数，从而获得脑功能活动信号，具体表示为：Step 3: Use the result obtained in Step 2 and or and Construct an adaptive filter function to obtain the brain function activity signal, specifically expressed as:

E(t)＝D(t)-X^T(t)W(t)E(t)＝D(t)-X^T (t)W(t)

其中，E(t)为脑功能活动信号，代表经过自适应滤波函数消减生理干扰后获得的氧合血红蛋白浓度变化量的时间信号Δ[HbO₂](t)或还原血红蛋白浓度变化量的时间信号Δ[HHb](t)；Among them, E(t) is the brain function activity signal, which represents the time signal Δ[HbO₂ ](t) of the concentration change of oxygenated hemoglobin obtained after the physiological interference is reduced by the adaptive filter function or the time signal of the change of the concentration of reduced hemoglobin Δ[HHb](t);

D(t)表示步骤二中获取的或包含人体生理干扰与诱发响应的脑功能活动信号；D(t) represents the obtained in step 2 or Brain function activity signals including human physiological interference and evoked responses;

X(t)为列向量，表示为X(t)＝[x(t) x(t-1) x(t-2) … x(t-M)]^T，M是自适应滤波器的阶数，为正整数，x(t-M)为信号x(t)的M个单元延迟，T为转置矩阵运算符号；X(t) is a column vector, expressed as X(t)=[x(t) x(t-1) x(t-2) ... x(tM)]^T , M is the order of the adaptive filter, is a positive integer, x(tM) is the M unit delay of the signal x(t), and T is the transpose matrix operation symbol;

x(t－1)为信号x(t)的一个单元延长；x(t－2)为信号x(t)的二个单元延长；x(t)表示步骤二中获取的或主要由人体生理干扰组成，为自适应滤波函数的参考信号；x(t-1) is a unit extension of signal x(t); x(t-2) is a two-unit extension of signal x(t); x(t) represents the obtained in step 2 or Mainly composed of human physiological interference, it is the reference signal of the adaptive filter function;

由于t的取值是1到N的正整数，因此对于任意的t，t-1表示t的前一采样时刻，同理t-2为t的前两个采样时刻，因此此处x(t－1)为信号x(t)的一个单元延迟，x(t－2)为信号x(t)的二个单元延迟。Since the value of t is a positive integer from 1 to N, for any t, t-1 represents the previous sampling moment of t, and similarly t-2 is the first two sampling moments of t, so here x(t -1) is one unit delay of signal x(t), and x(t-2) is two unit delay of signal x(t).

这个符号T是转置矩阵的运算表示符号。即X(t)表示列向量，则X^T(t)为X(t)的转置矩阵，因此X^T(t)为行向量；This symbol T is an operation representation symbol for transposing a matrix. That is, X(t) represents a column vector, then X^T (t) is the transpose matrix of X(t), so X^T (t) is a row vector;

W(t)＝[w₀(t)w₁(t)w₂(t)…w_M(t)]^T为自适应滤波函数的系数向量，w₀(t)、w₁(t)、w₂(t)…w_M(t)为列向量X(t)＝[x(t) x(t-1) x(t-2) … x(t-M)]^T中各对应元素的系数，其中，w₀(t)为x(t)的系数、w₁(t)为x(t-1)的系数、w₂(t)为x(t-2)的系数、w_M(t)为x(t-M)的系数；W(t)=[w₀ (t)w₁ (t)w₂ (t)...w_M (t)]^T is the coefficient vector of the adaptive filter function, w₀ (t), w₁ (t), w₂ (t)...w_M (t) is the coefficient of each corresponding element in the column vector X(t)=[x(t) x(t-1) x(t-2)...x(tM)]^T , Among them, w₀ (t) is the coefficient of x(t), w₁ (t) is the coefficient of x(t-1), w₂ (t) is the coefficient of x(t-2), w_M (t) is the coefficient of x(tM);

步骤四：将步骤三获得的脑功能活动信号E(t)的累积平方误差性能函数J(t)最小化，获得自适应滤波函数的系数向量W(t)的最优化系数向量W^*(t)，累积平方误差性能函数J(t)表示为：Step 4: Minimize the cumulative square error performance function J(t) of the brain function activity signal E(t) obtained in step 3, and obtain the optimal coefficient vector W^* (t) of the coefficient vector W(t) of the adaptive filter function ), the cumulative squared error performance function J(t) is expressed as:

$J J ((t t)) = = {Σ Σ}_{k k = = 00}^{t t} {χ χ}^{t t - - k k} {[[D D. ((k k)) - - {X x}^{T T} ((k k)) W W ((t t))]]}^{22}$

其中，χ为指数加权因子，k＝0，1，2，…，t；k，t均为正整数；Wherein, χ is an exponential weighting factor, k=0,1,2,...,t; k, t are positive integers;

D(k)表示步骤二中获取采样时刻为k时的或X(k)＝[x(k) x(k-1) x(k-2) … x(k-M)]^T，x(k)表示步骤二中获取采样时刻为k时的或x(k-1)为信号x(k)的一个单元延长；x(k-2)为信号x(k)的二个单元延长；x(k-M)为信号x(k)的M个单元延迟；D(k) means that in step 2 when the sampling time is k or X(k)＝[x(k) x(k-1) x(k-2) … x(kM)]^T , x(k) means the time when the sampling time is k in step 2 or x(k-1) is one unit extension of signal x(k); x(k-2) is two unit extension of signal x(k); x(kM) is M unit delay of signal x(k) ;

对累积平方误差性能函数J(t)相对于W(t)求导，则有：The cumulative square error performance function J(t) is derived relative to W(t), then:

$\frac{J J ((t t))}{\partial \partial W W ((t t))} = = - - 22 {Σ Σ}_{k k = = 00}^{t t} {χ χ}^{t t - - k k} X x ((t t)) [[D D. ((k k)) - - {X x}^{T T} ((k k)) W W ((t t))]]$

令上式为零，则自适应滤波函数最优化系数向量表示为：Let the above formula be zero, then the adaptive filter function optimization coefficient vector is expressed as:

W^*(t)＝R^-1(t)P(t)W^* (t) = R^-1 (t)P(t)

其中，R(t)为X(t)的确定性相关矩阵，P(t)为X(t)与D(t)之间的确定性互相关向量；Among them, R(t) is the deterministic correlation matrix of X(t), and P(t) is the deterministic cross-correlation vector between X(t) and D(t);

步骤五：根据步骤四得到的自适应滤波函数最优化系数向量W^*(t)得到脑功能活动信号E(t)：Step 5: According to the adaptive filter function obtained in step 4, the optimal coefficient vector W^* (t) is obtained to obtain the brain function activity signal E(t):

E(t)＝D(t)-X^T(t)W^*(t)；E(t)＝D(t)-XT(^t )W^* (t);

步骤六：对步骤五得到的脑功能活动信号E(t)进行最小二乘支持向量机回归处理，即采用如下形式的回归函数来表示剔除误差干扰后的脑功能活动信号：Step 6: Perform least squares support vector machine regression processing on the brain function activity signal E(t) obtained in step 5, that is, use the following regression function to represent the brain function activity signal after removing error interference:

其中，E^*(t)表示剔除误差干扰后的脑功能活动信号，V^T为权值向量，为特征向量，t＝1，2，…，N，b为偏置量；Among them, E^* (t) represents the brain function activity signal after removing the error interference, V^T is the weight vector, Is a feature vector, t=1, 2, ..., N, b is a bias;

此时，最小二乘支持向量机回归转化为处理如下的最优化问题：At this point, the least squares support vector machine regression is transformed into the following optimization problem:

$M m i i n no Q Q ((V V,, e e)) = = \frac{11}{22} {V V}^{T T} V V + + \frac{γ γ}{22} {Σ Σ}_{i i = = 11}^{N N} {e e}_{i i}^{22}$

t为正整数 t is a positive integer

其中，γ为规则化参数；e为利用自适应滤波函数所求取得脑功能活动信号与利用最小二乘支持向量机回归函数表示的脑功能活动信号之间的误差向量，表示为e＝[e₁,…,e_N]^T，1≤i≤t；Q表示函数；自变量是V和e；s.t.为使得E(i)满足E(i)为脑功能活动信号；Among them, γ is a regularization parameter; e is the error vector between the brain function activity signal obtained by using the adaptive filter function and the brain function activity signal represented by the least squares support vector machine regression function, expressed as e=[e₁ ,…,e_N ]^T , 1≤i≤t; Q represents a function; independent variables are V and e; st is such that E(i) satisfies E(i) is the brain function activity signal;

步骤七：为求步骤六中的最优化问题，构造如下的拉格朗日函数：Step 7: To solve the optimization problem in step 6, construct the following Lagrangian function:

其中，α_i为拉格朗日乘子，则根据最优化条件得：Among them, α_i is the Lagrangian multiplier, then according to the optimization condition:

对上式消去权值向量V和误差向量e，得如下线性方程组：Eliminating the weight vector V and the error vector e from the above formula, the following linear equations are obtained:

$[\begin{matrix} 00 & 11_{N N}^{T T} \\ 11_{N N} & Ω Ω + + {γ γ}^{- - 11} I I \end{matrix}] [\begin{matrix} b b \\ α α \end{matrix}] = = [\begin{matrix} 00 \\ E E. \end{matrix}]$

其中，E＝[E(1),E(2),…,E(N)]，1_N＝[1，…，1]^T，I为N×N的单位矩阵，α＝[α₁,…,α_N]^T为拉格朗日乘子向量，Ω为N×N的对称矩阵，令j＝1,2,…,t，t为正整数，则Ω的第ij项元素表示为：Wherein, E=[E(1),E(2),…,E(N)], 1_N =[1,…,1]^T , I is the unit matrix of N×N, α=[α₁ , …,α_N ]^T is the Lagrange multiplier vector, Ω is a N×N symmetric matrix, let j=1,2,…,t, t is a positive integer, then the ijth element of Ω is expressed as:

其中K(i,j)表示内积核函数，其定义为： Where K(i,j) represents the inner product kernel function, which is defined as:

式中，为t中第i个特征向量，i＝1，2，…，t；为t中第j个特征向量，j＝1，2，…，t；In the formula, is the i-th eigenvector in t, i=1, 2,..., t; is the jth eigenvector in t, j=1, 2,..., t;

步骤八：求解步骤七中的线性方程组，得到拉格朗日乘子向量α和偏置量b的数值解，此时利用最小二乘支持向量机回归函数处理后的脑功能活动信号表示为：Step 8: Solve the linear equations in step 7 to obtain the numerical solution of the Lagrangian multiplier vector α and the offset b. At this time, the brain function activity signal processed by the least squares support vector machine regression function is expressed as :

${E E.}^{* *} ((t t)) = = {Σ Σ}_{i i = = 11}^{N N} {α α}_{i i} K K ((t t,, j j)) + + b b$

K(t,j)表示内积核函数。K(t,j) represents the inner product kernel function.

发明效果Invention effect

本发明的优点在于针对多距离自适应滤波技术实现近红外脑功能活动检测时，当人体生理干扰与脑功能活动信号频带严重重叠情况下，获取的近红外脑功能活动信号中存在误差干扰，导致近红外脑功能活动信号的检测精度低这一问题，将基于多距离自适应滤波技术获取的脑功能活动信号看作是真实近红外脑功能活动信号与误差干扰的线性组合，并利用最小二乘支持向量机回归算法对真实近红外脑功能活动信号进行回归处理，从而实现消除误差干扰影响，实现对真实近红外脑功能活动信号的高精度检测，获得高质量的脑功能活动检测信号。最小二乘支持向量机回归方法是基于结构风险最小化准则的机器学习方法，与传统基于经验风险最小化的学习方法相比，它体现了模型复杂度与学习能力之间的折中思想，可有效抑制过拟合问题并提高回归处理精度，同时该方法采用等式约束将优化问题简化为线性方程组求解，从而消减该方法的计算复杂度，使其适用于近红外脑功能活动信号的处理分析。The advantage of the present invention is that when the multi-distance self-adaptive filtering technology is used to detect near-infrared brain function activity, when the frequency bands of human physiological interference and brain function activity signal overlap seriously, there will be error interference in the acquired near-infrared brain function activity signal, resulting in To solve the problem of low detection accuracy of near-infrared brain function activity signals, the brain function activity signals obtained based on multi-distance adaptive filtering technology are regarded as a linear combination of real near-infrared brain function activity signals and error interference, and the least squares method is used to The support vector machine regression algorithm performs regression processing on real near-infrared brain function activity signals, thereby eliminating the influence of error interference, realizing high-precision detection of real near-infrared brain function activity signals, and obtaining high-quality brain function activity detection signals. The least squares support vector machine regression method is a machine learning method based on the structural risk minimization criterion. Compared with the traditional learning method based on empirical risk minimization, it reflects the idea of compromise between model complexity and learning ability. Effectively suppress the overfitting problem and improve the accuracy of regression processing. At the same time, the method uses equality constraints to simplify the optimization problem into a linear equation solution, thereby reducing the computational complexity of the method and making it suitable for the processing of near-infrared brain function activity signals. analyze.

附图说明Description of drawings

图1为本发明流程图；Fig. 1 is a flowchart of the present invention;

图2为本发明中采用的多距离近红外脑功能活动检测探头与五层脑组织模型结构示意图，其中S代表双波长光源，D₁表示近端检测器，D₂表示远端检测器，R₁表示双波长光源S与近端检测器D₁之间的直线距离，R₂表示双波长光源S与远端检测器D₂之间的直线距离，L₁为头皮，L₂为颅骨，L₃为脑脊髓液，L₄为脑灰质，L₅为脑白质。Fig. 2 is the multi-distance near-infrared brain function activity detection probe that adopts in the present invention and the structural representation of five-layer brain tissue model, wherein S represents dual-wavelength light source, D₁ represents near-end detector, D₂ represents far-end detector, R₁ represents the straight-line distance between the dual-wavelength light source S and the near_- end detector D1, R2 represents the straight_- line distance between the_dual_- wavelength light source S and the far_- end detector D2, L1 is the scalp, L2 is the skull, L₃ is cerebrospinal fluid, L₄ is gray matter, and L₅ is white matter.

具体实施方式detailed description

具体实施方式一：结合图1和图2说明本实施方式，本实施方式的基于最小二乘支持向量机的近红外脑功能信号提取方法，具体是按照以下步骤制备的：Specific embodiment one: this embodiment is described in conjunction with Fig. 1 and Fig. 2, and the near-infrared brain function signal extraction method based on the least squares support vector machine of this embodiment is specifically prepared according to the following steps:

E(t)＝D(t)-X^T(t)W(t)E(t)＝D(t)-X^T (t)W(t)

此处公式的含义是t一共有N个采样时刻，那么对于其中的每一个采样时刻，即上式中的J(t)和W(t)中的t为1至N中的任意一个正整数Q，1《Q《N；此时等式右侧中的k的取值范围就是k＝0，1，2，…，Q；因此k的取值范围是随着当前t的取值不同而变化的，其最大值为当前的t的取值；The meaning of the formula here is that t has a total of N sampling moments, then for each of them, t in J(t) and W(t) in the above formula is any positive integer from 1 to N Q, 1<Q<N; at this moment, the value range of k in the right side of the equation is exactly k=0, 1, 2, ..., Q; therefore the value range of k is different with the value of current t changing, its maximum value is the current value of t;

此处公式中存在D(k)和X^T(k)，其中k的定义是从属于前面对采样时刻t的定义，即k的取值范围最大为当前的t值；t＝1，2，…，N；前面已经对D(t)进行定义如下：D(t)表示步骤二中获取的或即D(k)表示步骤二中获取采样时刻为k时的或前面已经对X(t)和x(t)进行定义如下：There are D(k) and X^T (k) in the formula here, where the definition of k is subordinate to the previous definition of sampling time t, that is, the value range of k is at most the current t value; t=1,2 ,..., N; D(t) has been defined as follows: D(t) represents the or That is, D(k) represents the time when the sampling time obtained in step 2 is k or X(t) and x(t) have been defined as follows:

X(t)为列向量，表示为X(t)＝[x(t)x(t-1)x(t-2)…x(t-M)]^T，M是自适应滤波器的阶数，为正整数，x(t-M)为信号x(t)的M个单元延迟，T为转置矩阵运算符号；x(t)表示步骤二中获取的或X(t) is a column vector, expressed as X(t)=[x(t)x(t-1)x(t-2)...x(tM)]^T , M is the order of the adaptive filter, is a positive integer, x(tM) is the M unit delay of the signal x(t), T is the transpose matrix operation symbol; x(t) represents the obtained in step 2 or

因此x(k)表示步骤二中获取采样时刻为k时的或X(k)＝[x(k) x(k-1) x(k-2) … x(k-M)]^TTherefore, x(k) means that when the sampling time is k in step 2, or X(k)＝[x(k) x(k-1) x(k-2) … x(kM)]^T

W^*(t)＝R^-1(t)P(t)W^* (t) = R^-1 (t)P(t)

前面定义了x(t)为自适应滤波函数的参考信号，X(t)＝[x(t) x(t-1) x(t-2) …x(t-M)]^T是由自适应滤波函数的参考信号所构成的参考信号向量；X^T(t)是X(t)的转置；It is defined that x(t) is the reference signal of the adaptive filter function, X(t)=[x(t) x(t-1) x(t-2) ... x(tM)]^T is determined by the adaptive filter The reference signal vector formed by the reference signal of the function; X^T (t) is the transposition of X (t);

此外，D(t)的含义就是前面所定义的：D(t)表示步骤二中获取的或x(t)的含义就是前面所定义的：x(t)表示步骤二中获取的或In addition, the meaning of D(t) is as defined above: D(t) means that obtained in step 2 or The meaning of x(t) is as defined above: x(t) represents the value obtained in step 2 or

在前面引入D(t)和x(t)的定义是为了后续简化介绍算法，D(t)表示步骤二中获取的或包含人体生理干扰与诱发响应的脑功能活动信号；x(t)表示步骤二中获取的或主要由人体生理干扰组成，为自适应滤波函数的参考信号；The definition of D(t) and x(t) introduced earlier is to simplify the introduction of the algorithm later, and D(t) represents the obtained in step 2 or Contains the brain function activity signal of human physiological interference and evoked response; x(t) represents the signal obtained in step 2 or Mainly composed of human physiological interference, it is the reference signal of the adaptive filter function;

在后续的算法介绍中：D(t)表示为时，则相应的x(t)为同理，当D(t)表示为时，则相应的x(t)为必须是同一种浓度的变化量相对应。In the follow-up algorithm introduction: D(t) is expressed as , then the corresponding x(t) is Similarly, when D(t) is expressed as , then the corresponding x(t) is It must correspond to the variation of the same concentration.

此处为减少歧义修改如下：In order to reduce ambiguity, the modification is as follows:

E(t)＝D(t)-X^T(t)W^*(t)；E(t)＝D(t)-XT(^t )W^* (t);

t为正整数 t is a positive integer

对上式消去权值向量V和误差向量e，得如下线性方程组：Eliminating the weight vector V and the error vector e from the above formula, the following linear equations can be obtained:

其中K(i,j)表示内积核函数，其定义为：内积核函数能够将原空间在高维特征空间映射变换后的内积计算，使用原始空间的函数形式进行表述，而无需知道具体的函数映射变换关系，通常对于一组给定的自变量x，共N个自变量 Among them, K(i, j) represents the inner product kernel function, which is defined as: the inner product kernel function can calculate the inner product after the original space is mapped and transformed in the high-dimensional feature space, and express it in the functional form of the original space without knowing The specific function mapping transformation relationship, usually for a given set of independent variables x, a total of N independent variables

其取值范围或是可以表述为x₁、x₂，…，x_NIts value range can be expressed as x₁ , x₂ , ..., x_N

则通常核函数可表示为：Then the kernel function can usually be expressed as:

此处，自变量x在高位特征空间的映射为相应的具体的某一个x_i的映射就是Here, the mapping of the independent variable x in the high-level feature space is The corresponding specific mapping of a certain x_i is

此处x表示所有的自变量x₁、x₂，…，x_N中的任意一个，x_i则只表示特定的自变量x_i； Here x represents any one of all independent variables x₁ , x₂ , ..., x_N , and_xi only represents a specific independent variable_xi ;

由于我们在前面已经定义了t的取值范围是1到N，此处i、j同样是相同的定义1到N；Since we have already defined the value range of t from 1 to N, here i and j are also the same definition 1 to N;

具体实施方式二：本实施方式与具体实施方式一不同的是：所述步骤一中5mm<R₁<10mm；30mm<R₂<40mm；Embodiment 2: This embodiment differs from Embodiment 1 in that: in step 1, 5mm<R₁ <10mm;30mm<R₂ <40mm;

其中，R₁为双波长光源S与检测器D₁之间的直线距离；R₂为双波长光源S与检测器D₂之间的距离。Among them, R1 is the linear distance between the dual_- wavelength light source S and the detector D1_; R2 is the distance between the_dual_- wavelength light source S and the detector D2.

其它步骤及参数与具体实施方式一相同。Other steps and parameters are the same as those in Embodiment 1.

具体实施方式三：本实施方式与具体实施方式一或二不同的是：所述R₁为10mm，R₂为40mm。其它步骤及参数与具体实施方式一或二相同。Embodiment 3: The difference between this embodiment and Embodiment 1 or 2 is that the R₁ is 10 mm, and the R₂ is 40 mm. Other steps and parameters are the same as those in Embodiment 1 or Embodiment 2.

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是：所述步骤一中双波长光源S发出两种近红外光的波长λ₁为770nm，λ₂为850nm。Specific embodiment four: this embodiment is different from one of specific embodiments one to three in that: in the step one, the wavelength λ₁ of two kinds of near-infrared light emitted by the dual-wavelength light source S is 770nm, and λ₂ is 850nm.

其它步骤及参数与具体实施方式一至三之一相同。Other steps and parameters are the same as those in Embodiments 1 to 3.

具体实施方式五：本实施方式与具体实施方式一至四之一不同的是：所述步骤二中和的具体计算公式为：Specific implementation mode five: the difference between this implementation mode and one of the specific implementation modes one to four is: in the step two and The specific calculation formula is:

$Δ Δ {[[{HbO HbO}_{22}]]}^{{R R}_{11}} ((t t)) = = \frac{{ϵ ϵ}_{H h H h b b} (({λ λ}_{11})) {ΔOD ΔOD}_{{λ λ}_{22}}^{{R R}_{11}} ((t t)) - - {ϵ ϵ}_{H h H h b b} (({λ λ}_{22})) {ΔOD ΔOD}_{{λ λ}_{11}}^{{R R}_{11}} ((t t))}{{R R}_{11} (({ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{22})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{11})) - - {ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{11})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{22})))) D D. P P F f}$

$Δ Δ {[[H h H h b b]]}^{{R R}_{11}} ((t t)) = = \frac{{ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{22})) {ΔOD ΔOD}_{{λ λ}_{11}}^{{R R}_{11}} ((t t)) - - {ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{11})) {ΔOD ΔOD}_{{λ λ}_{22}}^{{R R}_{11}} ((t t))}{{R R}_{11} (({ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{22})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{11})) - - {ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{11})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{22})))) D D. P P F f}$

$Δ Δ {[[{HbO HbO}_{22}]]}^{{R R}_{22}} ((t t)) = = \frac{{ϵ ϵ}_{H h H h b b} (({λ λ}_{11})) {ΔOD ΔOD}_{{λ λ}_{22}}^{{R R}_{22}} ((t t)) - - {ϵ ϵ}_{H h H h b b} (({λ λ}_{22})) {ΔOD ΔOD}_{{λ λ}_{11}}^{{R R}_{22}} ((t t))}{{R R}_{22} (({ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{22})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{11})) - - {ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{11})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{22})))) D D. P P F f}$

$Δ Δ {[[H h H h b b]]}^{{R R}_{22}} ((t t)) = = \frac{{ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{22})) {ΔOD ΔOD}_{{λ λ}_{11}}^{{R R}_{22}} ((t t)) - - {ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{11})) {ΔOD ΔOD}_{{λ λ}_{22}}^{{R R}_{22}} ((t t))}{{R R}_{22} (({ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{22})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{11})) - - {ϵ ϵ}_{{HbO HbO}_{22}} (({λ λ}_{11})) {ϵ ϵ}_{H h H h b b} (({λ λ}_{22})))) D D. P P F f}$

其中，ε_HHb(λ₁)是双波长光源S为波长λ₁时的还原血红蛋白消光系数；ε_HHb(λ₂)是双波长光源S为波长λ₂时的还原血红蛋白消光系数；是双波长光源S为波长λ₁时的氧合血红蛋白消光系数；是双波长光源S为波长λ₂时的氧合血红蛋白消光系数；DPF是差分路径因子。Wherein, ε_HHb (λ₁ ) is the reduced hemoglobin extinction coefficient when the dual-wavelength light source S is the wavelength λ₁ ; ε_HHb (λ₂ ) is the reduced hemoglobin extinction coefficient when the dual-wavelength light source S is the wavelength λ₂ ; is the extinction coefficient of oxyhemoglobin when the dual-wavelength light source S is the wavelength λ₁ ; is the extinction coefficient of oxyhemoglobin when the dual-wavelength light source S is the wavelength λ₂ ; DPF is the differential path factor.

其它步骤及参数与具体实施方式一至四之一相同。Other steps and parameters are the same as in one of the specific embodiments 1 to 4.

具体实施方式六：本实施方式与具体实施方式一至五之一不同的是：所述步骤四中χ＝0.99。Embodiment 6: The difference between this embodiment and one of Embodiments 1 to 5 is that χ=0.99 in step 4.

其它步骤及参数与具体实施方式一至五之一相同。Other steps and parameters are the same as one of the specific embodiments 1 to 5.

具体实施方式七：本实施方式与具体实施方式一至六之一不同的是：所述步骤四中R(t)和P(t)分别表示为：Specific embodiment seven: this embodiment is different from one of specific embodiments one to six in that: R (t) and P (t) in the step four are expressed as:

$R R ((t t)) = = {Σ Σ}_{k k = = 00}^{t t} {χ χ}^{t t - - k k} X x ((k k)) {X x}^{T T} ((k k))$

$P P ((t t)) = = {Σ Σ}_{k k = = 00}^{t t} {χ χ}^{t t - - k k} X x ((k k)) D D. ((k k)) . .$

其它步骤及参数与具体实施方式一至六之一相同。Other steps and parameters are the same as one of the specific embodiments 1 to 6.

具体实施方式八：本实施方式与具体实施方式一至七之一不同的是：所述步骤七中的核函数为径向基核函数。Embodiment 8: The difference between this embodiment and one of Embodiments 1 to 7 is that the kernel function in step 7 is a radial basis kernel function.

其它步骤及参数与具体实施方式一至七之一相同。Other steps and parameters are the same as one of the specific embodiments 1 to 7.