小波时间散射ECG信号的分类
这个例子展示了如何将人类心电图(ECG)信号利用小波时间散射和支持向量机(SVM)分类器。金宝app在小波散射,数据传播通过一系列的小波变换,非线性,平均生产时间序列的低温度差的表示。小波时间收益率散射信号表示对输入信号的变化而不牺牲类辨别力。你必须有小波工具箱™和统计和机器学习工具箱™运行这个示例。在这个例子中使用的数据是公开的生理网。你可以找到一个深度学习在这个例子中这种分类方法问题时间序列分类使用小波分析和深度学习在这个例子和机器学习方法信号使用小波的特性和支持向量机分类金宝app。
关于术语的:在小波散射的上下文中,术语“时间窗”是指样品的数量将采样后得到的输出平滑操作。有关更多信息,请参见时间窗口。
数据描述
这个示例使用心电图数据来自三组,或类,人:人与心律失常、人与充血性心力衰竭,而人与正常窦性节律。162例使用心电图记录生理网从三个数据库:MIT-BIH心律失常数据库[3][5],MIT-BIH正常窦性心律的数据库[3]BIDMC充血性心力衰竭的数据库[2][3]。总共有96个录音与心律失常的人,30录音与充血性心力衰竭的人,36个录音从人与正常窦性节律。目的是训练一个分类器区分心律失常(ARR),充血性心力衰竭(CHF)和正常窦性心律(NSR)。
下载数据
第一步是下载的数据GitHub库。下载数据,点击代码并选择下载ZIP。保存文件physionet_ECG_data-main.zip在一个文件夹,你有写权限。说明对于这个示例假设您已经下载了文件到临时目录,(tempdir在MATLAB中)。修改后续说明压缩和加载数据如果您选择下载一个文件夹中的数据不同tempdir。
该文件physionet_ECG_data-main.zip包含
ECGData.zip
README.md
和ECGData。zip包含
ECGData.mat
Modified_physionet_data.txt
License.txt。
ECGData.mat持有s the data used in this example. The .txt file, Modified_physionet_data.txt, is required by PhysioNet's copying policy and provides the source attributions for the data as well as a description of the pre-processing steps applied to each ECG recording.
加载文件
如果你遵循下载指令在前面的小节中,输入以下命令解压两个存档文件。
解压缩(fullfile (tempdir,“physionet_ECG_data-main.zip”),tempdir)解压缩(fullfile (tempdir,“physionet_ECG_data-main”,“ECGData.zip”),…fullfile (tempdir“ECGData”))
你解压后ECGData。zip文件,数据加载到MATLAB。
负载(fullfile (tempdir“ECGData”,“ECGData.mat”))
ECGData与两个字段是一个结构数组:数据和标签。数据是162 - 65536矩阵每一行是一个心电图记录采样128赫兹。每个ECG时间序列的总持续时间512秒。标签是162 - 1单元阵列的诊断标签,一个用于每一行的数据。这三个诊断类别是:“加勒比海盗”(心律失常),瑞士法郎(充血性心力衰竭),和“NSR”(正常窦性心律)。
创建培训和测试数据
随机将数据分成两组,训练和测试数据集。辅助函数helperRandomSplit执行随机分裂。helperRandomSplit接受培训数据和所需的分离比例ECGData。的helperRandomSplit函数输出两个数据集和一组标签。每一行的trainData和testData是一个心电图信号。的每个元素trainLabels和testLabels包含的类标签对应的行数据矩阵。在这个例子中,我们随机分配70%的数据在每个类的训练集,其余30%是伸出测试(预测)和分配给测试组。
percent_train = 70;[trainData, testData trainLabels testLabels] =…helperRandomSplit (percent_train ECGData);
有113条记录trainData设置和49记录testData。通过设计的训练数据的数据包含69.75% (113/162)。记得ARR类代表59.26%的数据(96/162),瑞士法郎类占18.52%(30/162),和NSR类占22.22% (36/162)。检查每个类的百分比在训练集和测试集。每个符合整个班上百分比百分比的数据集。
Ctrain = countcats(分类(trainLabels)。/元素个数(trainLabels)。* 100
Ctrain =3×159.2920 18.5841 22.1239
ct = countcats(分类(testLabels)。/元素个数(testLabels)。* 100
ct =3×159.1837 18.3673 22.4490
情节样品
情节的几千个样品的四个随机选择的记录ECGData。辅助函数helperPlotRandomRecords这是否。helperPlotRandomRecords接受ECGData和一个随机种子作为输入。最初的种子将在14所以至少一个记录从每个类是策划。你可以执行helperPlotRandomRecords与ECGData作为唯一的输入参数多次你想了解各种ECG波形与每个类有关。你可以找到的源代码和所有辅助函数在支持功能部分的这个例子。金宝app
helperPlotRandomRecords (ECGData 14)
小波时间散射
小波时间散射的关键参数来指定网络的规模不变时,小波变换,和小波的数量每倍频程的每个小波滤波器。在许多应用程序中,两个滤波器的级联能充分实现良好的性能。在这个例子中,我们构造一个小波散射与默认的过滤器银行网络:在第一个过滤器银行每倍频程8小波和1小波每八度第二滤波器组。不变性范围设置为150秒。
N =大小(ECGData.Data 2);sn = waveletScattering (“SignalLength”N“InvarianceScale”,150,“SamplingFrequency”,128);
你可以想象两个滤波器的小波滤波器。
[fb, f, filterparams] = filterbank (sn);图次要情节(211)情节(f, fb {2} .psift) xlim(128[0])网格在标题(1日滤波器组小波滤波器的);次要情节(212)情节(f, fb {3} .psift) xlim(128[0])网格在标题(第二滤波器组小波滤波器的);包含(“赫兹”);
展示的不变性,得到尺度函数的傅里叶反变换和中心在0秒。两个垂直的黑色线条标志着-75年和75年第二边界。此外,情节的实部和虚部coarsest-scale(最低频率)小波从第一个过滤器银行。注意coarsest-scale小波不超过的不变的规模由时间决定支持缩放功能。金宝app这是一个重要的属性散射波的时间。
图;φ= ifftshift(传输线(神奇动物{1}.phift));psiL1 = ifftshift(传输线(神奇动物{2}.psift (:,)));t = (2 ^ 15:2 ^ 15 - 1) * 1/128;scalplt =情节(t,φ);持有在网格在ylim(1.6[-1.5军医的军医]);情节(-75[-75],[-1.5 1.6的军医]的军医,“k——”);情节(75[75],[-1.5 1.6的军医]的军医,“k——”);包含(“秒”);ylabel (“振幅”);wavplt =情节(t,真实(psiL1)图像放大(psiL1)]);传奇([scalplt wavplt (1) wavplt (2)), {“扩展功能”,“Wavelet-Real部分”,“Wavelet-Imaginary部分”});标题({“扩展功能”;“Coarsest-Scale小波第一过滤器银行”})举行从
构建网络散射后,获取训练数据的散射系数矩阵。当您运行featureMatrix与多个信号,每一列都被视为一个单一的信号。
scat_features_train = featureMatrix (sn, trainData ');
的输出featureMatrix在本例中是409 - - 16 - 113。每一页的张量,scat_features_train一个信号的散射变换。小波散射变换非常及时downsampled基于尺度函数的带宽。在这种情况下,这个结果在16个时间窗口的每个409散射路径。
为了获得一个矩阵与支持向量机分类器兼容,重塑multisignal散射变换成一个矩阵中每一列对应一个散射路径和每一行是一个散射时间窗口。在这种情况下,您获得1808行,因为有16次windows的113在训练数据信号。
Nwin =大小(scat_features_train, 2);scat_features_train =排列(scat_features_train [2 3 1]);scat_features_train =重塑(scat_features_train,…大小(scat_features_train 1) *大小(scat_features_train, 2), []);
重复这个过程的测试数据。最初,scat_features_test由- 16 - 409 - 49因为有49心电图波形的SVM分类器训练集。改造后,该特性矩阵是784 - 416。
scat_features_test = featureMatrix (sn, testData ');scat_features_test =排列(scat_features_test [2 3 1]);scat_features_test =重塑(scat_features_test,…大小(scat_features_test 1) *大小(scat_features_test, 2), []);
因为我们每个信号获得16散射窗口,我们需要创建标签匹配窗口的数量。辅助函数createSequenceLabels这是基于windows的数量。
[sequence_labels_train, sequence_labels_test] = createSequenceLabels (Nwin、trainLabels testLabels);
交叉验证
两个分类,分析执行。第一个使用所有的散射数据,适合多层次支持向量机二次内核。总共有2592个散射序列在整个数据集,16为每个162信号。错误率或损失,预计用5倍交叉验证。
scat_features = [scat_features_train;scat_features_test];allLabels_scat = [sequence_labels_train;sequence_labels_test];rng (1);模板= templateSVM (…“KernelFunction”,多项式的,…“PolynomialOrder”2,…“KernelScale”,“汽车”,…“BoxConstraint”,1…“标准化”,真正的);classificationSVM = fitcecoc (…scat_features,…allLabels_scat,…“学习者”模板,…“编码”,“onevsone”,…“类名”,{“加勒比海盗”;瑞士法郎的;“签约”});kfoldmodel = crossval (classificationSVM,“KFold”5);
计算损失和混淆矩阵。显示的准确性。
predLabels = kfoldPredict (kfoldmodel);损失= kfoldLoss (kfoldmodel) * 100;confmatCV = confusionmat (allLabels_scat predLabels)
confmatCV =3×31535 0 1 1 479 0 0 0 576
流(准确性是% 2.2 f %。\ n '100 -损失);
精度为99.92%。
精度是99.88%,这是很好,但实际结果更好的考虑到每个时间窗口是单独分类。有16个独立的每个信号的分类。使用一个简单多数票获得一个类为每个散射表征预测。
类=分类({“加勒比海盗”,瑞士法郎的,“签约”});[ClassVotes, ClassCounts] = helperMajorityVote (predLabels [trainLabels;testLabels)、类);
确定实际的基于模式类的交叉验证的准确性预测每组散射的时间窗口。如果一个给定的模式并不是唯一的,helperMajorityVote返回一个表示分类错误“NoUniqueMode”。这导致一个额外的列混淆矩阵,在这种情况下都是0,因为每个组的独特模式存在散射预测。
CVaccuracy =总和(eq (ClassVotes分类([trainLabels;testLabels]))) / 162 * 100;流(“真正的交叉验证的准确性是% 2.2 f %。\ n”,CVaccuracy);
真正的交叉验证的准确性为100.00%。
MVconfmatCV = confusionmat(分类([trainLabels;testLabels]), ClassVotes);MVconfmatCV
MVconfmatCV =4×496 30 0 0 0 0 0 0 0 0 36 0 0 0 0 0
散射在旨在正确分类的所有信号模型。如果你检查ClassCounts,你会发现2分类错误的时间窗口confmatCV是由于2信号15 16散射的窗户被正确分类。
支持向量机分类
接下来的分析中,我们只适合多层次二次SVM训练数据(70%),然后使用该模型进行预测的30%的数据进行测试。有49个数据记录在测试集。使用多数表决单个散射窗口。
模型= fitcecoc (…scat_features_train,…sequence_labels_train,…“学习者”模板,…“编码”,“onevsone”,…“类名”,{“加勒比海盗”,瑞士法郎的,“签约”});predLabels =预测(模型、scat_features_test);[TestVotes, TestCounts] = helperMajorityVote (predLabels、testLabels、类);testaccuracy =总和(eq (TestVotes分类(testLabels))) /元素个数(testLabels) * 100;流(的测试精度是2.2 f % %。\ n”,testaccuracy);
测试精度为97.96%。
confusionchart(分类(testLabels)、TestVotes)
测试数据集上的分类精度为98%左右。混淆矩阵显示一个记录是瑞士法郎并被错误地归类为加勒比海盗。所有48个其他信号正确分类。
总结
本例中使用小波散射和一个支持向量机分类器对心电图波形进行分类成三个诊断类之一。小波散射器被证明是一个强大的功能,只需要一个用户指定参数的最小集合产生一组强大的特性分类。比较这的例子信号使用小波的特性和支持向量机分类金宝app,这需要大量的专业手工特征用于分类。与小波散射,您只需要指定时间不变性的规模,滤波器的数量(或小波变换),和小波的数量每八度。的组合小波散射变换和支持向量机分类器分类旨在模式和收益率为100% 98%正确分类应用SVM的散射变换时尚未签署的测试集。
引用
然后,J。,Mallat, S。2014. Deep scattering spectrum, IEEE Transactions on Signal Processing, 62, 16, pp. 4114-4128.
拜姆DS,科鲁奇WS Monrad ES,史密斯HS,赖特射频,Lanoue, Gauthier DF, Ransil BJ,格罗斯曼W, Braunwald大肠严重充血性心力衰竭患者的生存与口服药物治疗。美国心脏病学院1986年3月;7 (3):661 - 670。
Goldberger AL、Amaral局域网、玻璃L,豪斯多夫JM,伊万诺夫PCh,马克RG Mietus我,穆迪GB,彭ck,斯坦利。PhysioBank、PhysioToolkit和生理网:组件的一个新的研究资源对于复杂的生理信号。循环。23号卷。101年,2000年6月13日,pp. e215-e220。
http://circ.ahajournals.org/content/101/23/e215.fullMallat, S。,2012年。集团不变的散射。通信在纯和应用数学,65,10,1331 - 1398页。
穆迪GB,马克RG。数据库MIT-BIH心律失常的影响。IEEE Eng在医学和生物20(3):45 - 50(2001年5月- 6月)。(PMID: 11446209)
金宝app支持功能
helperPlotRandomRecords情节四ECG信号随机选择ECGData。
函数helperPlotRandomRecords (ECGData randomSeed)%这个函数仅仅是为了支持XpwWaveletMLExample。金宝app它可能变化百分比或在将来的版本中被删除。如果输入参数个数= = 2 rng (randomSeed)结束M =大小(ECGData.Data, 1);idxsel = randperm (M, 4);为numplot = 1:4次要情节(2,2,numplot)情节(ECGData.Data (idxsel (numplot) 1:3000)) ylabel (“伏”)如果numplot > 2包含(“样本”)结束标题(ECGData.Labels {idxsel (numplot)})结束结束
helperMajorityVote发现的模式预测类标签对每组散射的时间窗口。函数返回的类标签模式和类预测的数量每组散射的时间窗口。如果没有独特的模式,helperMajorityVote返回一个类标签的“错误”,表明散射窗口是一个分类错误。
函数[ClassVotes, ClassCounts] = helperMajorityVote (predLabels、origLabels类)%这个函数支持ECGWaveletTimeS金宝appcatteringExample。它可能变化百分比或在将来的版本中被删除。%进行分类数组如果标签还没有直言predLabels =分类(predLabels);origLabels =分类(origLabels);%预计predLabels和origLabels分类向量npr =元素个数(predLabels);Norig =元素个数(origLabels);Nwin = npr / Norig;predLabels =重塑(predLabels Nwin Norig);ClassCounts = countcats (predLabels);[mxcount, idx] = max (ClassCounts);ClassVotes =类(idx);%检查任何关系的最大值,并确保它们标记为%如果模式出现不止一次错误modecnt = modecount (ClassCounts mxcount);ClassVotes (modecnt > 1) =分类({“错误”});ClassVotes = ClassVotes (:);% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -函数modecnt = modecount (ClassCounts mxcount) modecnt =正(大小(ClassCounts, 2), 1);为数控= 1:尺寸(ClassCounts, 2) modecnt (nc) = histc (ClassCounts(:,数控),mxcount (nc));结束结束% EOF结束
