[发明专利]多被试复数fMRI数据移不变CPD分析方法

摘要

一种引入空间源相位稀疏约束的多被试复数fMRI数据移不变CPD分析方法，属于医学信号处理领域。在移不变CPD算法的基础上，利用交替最小二乘法更新被试共享SM成分、共享TC成分和各被试强度；将基于实数数据的被试时延估计扩展到复数数据的时延估计；对共享SM成分采用空间源相位稀疏约束方式再次更新，具体包括两步：首先对共享SM成分进行相位校正；然后对校正后的空间源相位的大相位体素，利用平滑L0范数逼近函数，对共享SM成分更新。本发明能有效提取出多被试复数fMRI数据中更为全面的共享复数脑功能信息，这些信息在今后脑认知和脑疾病研究都有非常好的应用前景。

主权项

1.一种引入空间源相位稀疏约束的多被试复数fMRI数据移不变CPD分析方法，其特征包括以下步骤：第一步：输入多被试复数fMRI数据其中V、J和K分别表示脑内体素数目、扫描次数和被试数目；；第二步：初始化；；设成分个数为N，N为大于0的正整数；；随机初始化共享SM成分共享TC成分和被试强度初始化被试时延为零矩阵；令迭代次数iter＝0，相对误差Δε_iter＝1，计算迭代误差ε_iter：其中，式(1)也是移不变CPD算法模型，τ_k,n表示第k被试第n个成分的时延，令其为整数；b_n(j‑τ_k,n)表示为b_j,n时移了τ_k,n个点，具体地，若τ_k,n>0，第k被试第n个TC成分相对共享TC成分b_n循环左移τ_k,n个点，否则若τ_k,n<0，则循环右移|τ_k,n|个点；第三步：更新共享TC成分B；采用等人提出的移不变CPD算法中共享TC成分B的更新方法(对B进行更新；第四步：对共享SM成分S和X进行降维；对S和X降维成和从而保证第五步的被试时延快速估计；其中为X的1模展开形式，将张量化为第五步：更新被试时延对复数数据的被试时延进行估计；对于第k被试第n'(n'＝1,…,N)个成分的时延τ_k,n'估计，首先进行如下一些定义；令为进行3模展开矩阵的第k行向量，定义向量满足为第n'个成分的信号，其中c_k,n为C的元素；为的第n列向量；的元素表示为(j＝1,…,J；n＝1,…,N)；定义Re{·}和Im{·}分别表示取实部和虚部，并做如下运算：其中，j＝1,…,J，为的第n'个元素；接着将g_RR(j)、g_IR(j)、g_RI(j)和g_II(j)(j＝1,…,J)转换到频域形式和f＝1,…,F；继续做如下运算：其中，上标^*表示取共轭，和分别表示对共享TC成分b_n'的实部和虚部转换到频域形式后的第f个元素；再将和(f＝1,…,F)傅里叶反变换到时域形式φ_RRR(j)、φ_RRI(j)、φ_IRR(j)、φ_IRI(j)、φ_RIR(j)、φ_RII(j)、φ_IIR(j)和φ_III(j)，j＝1,…,J；令(j＝1,…,J)满足的绝对值最大化实际上等价于均方误差最小化，那么第k个被试第n'个成分的时延τ_k,n'通过最大化式(5)来获得，具体实现是：被试时延τ_k,n'为整数；对被试时延的所有元素τ_k,n'，k＝1,…,K，n'＝1,…,N；按照上述方法进行更新；第六步：更新共享SM成分S；采用提出的移不变CPD算法中更新S的方法对S更新；第七步：对共享SM成分S进行相位校正；首先求取联合混合矩阵其元素满足z_j+k(J‑1),n＝c_k,nb_n(j‑τ_k,n)；对每个成分n的共享SM成分s_n，n＝1,…,N，将TC成分a_n替换成联合混合向量z_n，即通过对exp(‑iθ)z_n进行实部能量最大化，获取相位校正旋转角度θ_n，得到相位校正后的共享SM成分0≤θ≤π，exp{·}为指数函数；第八步：对共享SM成分S进行空间源相位稀疏约束更新；将S所有体素相位值从大到小排序，记为第V/3个体素的相位值；对S的V/3大相位体素值，采用L0范数平滑函数进行稀疏约束：其中，s_v,n为S的元素，v＝1,…,V,n＝1,…,N，F_σ(|s_v,n|)＝f_σ(|s_v,n|)，且其中参数σ充分小时，f_σ(|s_v,n|)接近L0范数；σ越大，f_σ(|s_v,n|)越平滑；具体采用最速下降法对S进行更新：其中λ为正定步长，ΔS的每个元素Δs_v,n(v＝1,…,V,n＝1,…,N)满足Δs_v,n＝exp{θ(s_v,n)}f′_σ(|s_v,n|)          (10)在此，θ(s_v,n)为s_v,n(v＝1,…,V；n＝1,…,N)的相位值，f′_σ(|s_v,n|)为f_σ(|s_v,n|)的一阶导数，满足第九步：更新被试强度C；根据提出的移不变CPD算法中更新C的方法对C更新；第十步：计算误差；令iter←iter+1；σ←0.999σ，使σ缓慢变化；根据式(1)，计算本次迭代误差ε_iter，以及相对误差Δε_iter：Δε_iter＝|(ε_iter‑1‑ε_iter)/ε_iter‑1|           (12)第十一步：若ε_iter小于预设误差阈值ε_{iter_min}，跳转到第十四步，否则执行第十二步；第十二步：若Δε_iter小于预设相对误差阈值Δε_{iter_min}，跳转到第十四步，否则执行第十三步；第十三步：若iter大于预设最大迭代次数iter_max，跳转到第十四步，否则执行第三步；第十四步：对共享SM成分S进行相位消噪；对每个成分n的共享SM成分s_n，进行相位消噪；第十五步：输出相位消噪的共享SM成分S、共享TC成分B、被试时延与被试强度C。