python数据挖掘分析案例_CCA典型关联分析原理与Python案例

点击上面"脑机接口社区"关注我们

更多技术干货第一时间送达

eec1666a0d98cbc8afc9e5f5e42e158e.gif

Hello,大家好!

Rose今天分享一下CCA的相关原理以及Python应用,CCA在EEG等脑电数据的特征提取中使用很多,很有必要熟悉其原理。

CCA典型相关分析


CCA(canonical correlation analysis)利用综合变量对之间的相关关系来反映两组指标之间的整体相关性的多元统计分析方法。它的基本原理是:为了从总体上把握两组指标之间的相关关系,分别在两组变量中提取有代表性的两个综合变量U1和V1(分别为两个变量组中各变量的线性组合),利用这两个综合变量之间的相关关系来反映两组指标之间的整体相关性。


1936年,Hotelling提出典型相关分析。考虑两组变量的线性组合, 并研究它们之间的相关系数p(u,v).在所有的线性组合中, 找一对相关系数最大的线性组合, 用这个组合的单相关系数来表示两组变量的相关性, 叫做两组变量的典型相关系数, 而这两个线性组合叫做一对典型变量。在两组多变量的情形下, 需要用若干对典型变量才能完全反映出它们之间的相关性。下一步, 再在两组变量的与u1,v1不相关的线性组合中, 找一对相关系数最大的线性组合, 它就是第二对典型变量, 而且p(u2,v2)就是第二个典型相关系数。这样下去, 可以得到若干对典型变量, 从而提取出两组变量间的全部信息。

典型相关分析的实质就是在两组随机变量中选取若干个有代表性的综合指标(变量的线性组合), 用这些指标的相关关系来表示原来的两组变量的相关关系。这在两组变量的相关性分析中, 可以起到合理的简化变量的作用; 当典型相关系数足够大时, 可以像回归分析那样, 由- 组变量的数值预测另一组变量的线性组合的数值。

原理描述


3fe0e74f2819e105da13c7c7f0aba157.png

b1092be7118b30603f6399bc365e87d1.png

bade7c5246f6eb4905df5e59ce9204f4.png

4b1ed0c020c3f81730b51b2784054c60.png

f8d5fa937f76e9992fa1d57eafb8789a.png

8f3f10486beb52855bb5d59410f05683.png

875cb8b1a9f3fc9b5f97eec6c6812285.png

6a652807cf52037eaa18639f12e5022d.png

cecee05d2f565e8d662b92e97c8a7951.png

b6f265886022be7ba23d1c139b189a5e.png

案例实现

# 导入工具包import h5pyimport rccaimport sysimport numpy as npimport cortexzscore = lambda d: (d-d.mean(0))/d.std(0)

第一步:加载数据

请从CRCNS下载数据:http://crcns.org/data-sets/vc/vim-2以下分析假定该数据位于当前目录中名为“ data”的目录中。

data = []vdata = []numSubjects = 3# subjects 是3个受试者列表.subjects = ['S1', 'S2', 'S3']# xfms 是Pycortex中变换名称的列表,该名称用于对齐每个受试者的功能和解剖数据。xfms = ['S1_xfm', 'S2_xfm', 'S3_xfm']dataPath ="./data/VoxelResponses_subject%d.mat"for subj in range(numSubjects):    # 打开数据    f = h5py.File(dataPath % (subj+1),'r')    #  获取数据大小    datasize = (int(f["ei"]["datasize"].value[2]),int(f["ei"]["datasize"].value[1]),int(f["ei"]["datasize"].value[0]))    # 从Pycortex获取皮质面罩    mask = cortex.db.get_mask(subjects[subj], xfms[subj], type = 'thick')    # 获取该受试者的训练数据    data_subj = np.nan_to_num(zscore(np.nan_to_num(f["rt"].value.T)))    data.append(data_subj.reshape((data_subj.shape[0],)+datasize)[:, mask])    # 获取受试者的验证数据    vdata_subj = np.nan_to_num(zscore(np.nan_to_num(f["rv"].value.T)))    vdata.append(vdata_subj.reshape((vdata_subj.shape[0],)+datasize)[:, mask])

第二步:定义CCA参数

# 这里设置1e-4和1e2之间的一系列正则化值regs = np.array(np.logspace(-4, 2, 10))# 这里考虑3到10之间的成分数量numCCs = np.arange(3, 11)# 初始化cca模型cca = rcca.CCACrossValidate(numCCs=numCCs, regs=regs)

第三步:对数据训练,分析并保存分析结果

"""说明:由于数据量大,此分析的计算量很大。在笔记本中运行它会花费大量时间,因此建议对其进行并行化和/或在计算机群集上运行它,然后加载结果以进行可视化。"""# 利用cca训练数据cca.train(data)# 利用cca对验证数据进行验证cca.validate(vdata)# 计算方差,解释每个体素中的验证响应cca.compute_ev(vdata)# 保存分析结果cca.save("./data/CCA_results.hdf5")

第四步:可视化分析结果

# 导入可视化工具包%matplotlib inlineimport matplotlib.pyplot as plt# 导入Brewer色彩图以进行可视化from brewer2mpl import qualitativenSubj = len(cca.corrs)nBins = 30bmap = qualitative.Set1[nSubj]f = plt.figure(figsize = (8, 6))ax = f.add_subplot(111)for s in range(nSubj):    # 绘制所有三个对象的所有体素之间的相关性直方图    ax.hist(cca.corrs[s], bins = nBins, color = bmap.mpl_colors[s], histtype="stepfilled", alpha = 0.6)plt.legend(['Subject 1', 'Subject 2', 'Subject 3'], fontsize = 16)ax.set_xlabel('Prediction correlation', fontsize = 20)ax.set_ylabel('Number of voxels', fontsize = 20)ax.set_title("Prediction performance across voxels", fontsize = 20)# p <0.05时的显着性阈值(针对多次比较进行了校正)# 重要性是使用渐近方法计算的(有关详细信息,请参见论文文本)thresh = 0.0893ax.axvline(x = thresh, ymin = 0, ymax = 7000, linewidth = 2, color = 'k')ax.text(thresh+0.05, 5000, 'p<0.05', fontsize = 16)ax.set_xticklabels(0.1*np.arange(-8, 11, 2), fontsize = 16)ax.set_yticklabels(np.arange(0, 10000, 1000), fontsize = 16)

3086562987ea7caf29b55896cfc28cb9.png

该图显示了皮质图上一个对象的跨学科预测结果,即预测的验证响应和实际的验证响应之间的相关性。不重要的相关性(p <0.05,已针对多个比较进行校正)设置为0。

import cortexfrom matplotlib import cmfrom copy import deepcopysubj = 0subjName = "S1"subjTransform = "S1_xfm"corrs = deepcopy(cca.corrs[subj])# 将所有低于显着性阈值的体素设置为0corrs[corrs_ = cortex.quickflat.make_figure(cortex.Volume(corrs, subjName,                                                subjTransform,                                                cmap = cm.PuBuGn_r,                                               vmin = 0., vmax = 1.), with_curvature = True)

dc1f2b8ccd5bb307b2acb8126a650ecd.png

参考:

https://www.cnblogs.com/jerrylead/archive/2011/06/20/2085491.html

https://github.com/gallantlab/pyrcca

文章用于学术交流,不用于商业行为,若有侵权及疑问,请后台留言!

更多阅读

蒲慕明:200亿投资,用于探索脑科学前沿领域

Python处理脑电--Epochs数据可视化

EEGLAB 使用教程 3 -参考电极和重采样

第2期 | 国内脑机接口领域专家教授汇总

精彩长文 | 脑机接口技术的现状与未来!

ICA处理脑电资料汇总

收藏 | 脑电EEG基础与处理汇总

脑机接口BCI学习交流QQ群:903290195

微信群请扫码添加,Rose拉你进群

(请务必填写备注,eg. 姓名+单位+专业/领域/行业)

306f328899144d3f122b266a4464506f.png

长按加群

1492f806d84be4ef64ba89c0246cfc00.png

欢迎点个在看鼓励一下

bf2dfd4c8bb2d76561d66750e0250026.gif
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值