В этом примере показано, как сегментировать человеческую электрокардиограмму (ECG) сигналы с помощью текущих нейронных сетей для глубокого обучения и частотно-временного анализа.
Электрическое действие в человеческом сердце может быть измерено как последовательность амплитуд далеко от базового сигнала. Для одного нормального цикла heartbeat сигнал ECG может быть разделен на следующую морфологию удара [1]:
P волна — маленькое отклонение перед комплексом QRS представление предсердной деполяризации
Комплекс QRS — Само-амплитудный фрагмент heartbeat
T волна — маленькое отклонение после комплекса QRS представление желудочковой реполяризации
Сегментация этих областей форм волны ECG может обеспечить основание для измерений, полезных для оценки полного здоровья человеческого сердца и присутствия отклонений [2]. Вручную аннотирование каждой области сигнала ECG может быть утомительной и длительной задачей. Обработка сигналов и методы глубокого обучения потенциально могут помочь оптимизировать и автоматизировать аннотацию необходимой области.
Этот пример использует сигналы ECG от общедоступной Базы данных QT 3[] 4[]. Данные состоят примерно из 15 минут записей ECG, с частотой дискретизации 250 Гц, измеренных от в общей сложности 105 пациентов. Чтобы получить каждую запись, ревизоры поместили два электрода в другие места на груди пациента, приводящей к двухканальному сигналу. База данных обеспечивает метки области сигнала, сгенерированные автоматизированной экспертной системой [2]. Этот пример стремится использовать решение для глубокого обучения обеспечить метку для каждой выборки сигнала ECG согласно области, где выборка расположена. Этот процесс маркировки необходимых областей через сигнал часто упоминается как сегментация формы волны.
Чтобы обучить глубокую нейронную сеть классифицировать области сигнала, можно использовать сеть Long Short-Term Memory (LSTM). В этом примере показано, как методы предварительной обработки сигнала и частотно-временной анализ могут использоваться, чтобы улучшать производительность сегментации LSTM. В частности, пример использует synchrosqueezed преобразование Фурье, чтобы представлять неустановившееся поведение сигнала ECG.
Каждый канал 105 двухканальных сигналов ECG был помечен независимо автоматизированной экспертной системой и обработан независимо для в общей сложности 210 сигналов ECG, которые хранились вместе с метками области в 210 MAT-файлах. Файлы доступны в следующем местоположении: https://www.mathworks.com/supportfiles/SPT/data/QTDatabaseECGData.zip.
Загрузите файлы данных в свою временную директорию, местоположение которой задано tempdir
MATLAB® команда. Если вы хотите поместить файлы данных в папку, отличающуюся от
tempdir
, измените имя каталога в последующих инструкциях.
% Download the data dataURL = 'https://www.mathworks.com/supportfiles/SPT/data/QTDatabaseECGData.zip'; datasetFolder = fullfile(tempdir,'QTDataset'); zipFile = fullfile(tempdir,'QTDatabaseECGData.zip'); if ~exist(datasetFolder,'dir') websave(zipFile,dataURL); end % Unzip the data unzip(zipFile,tempdir)
unzip
операция создает QTDatabaseECGData
папка в вашей временной директории с 210 MAT-файлами в нем. Каждый файл содержит сигнал ECG в переменной ecgSignal
и таблица области помечает в переменной signalRegionLabels
. Каждый файл также содержит частоту дискретизации сигнала в переменной Fs
. В этом примере все сигналы имеют частоту дискретизации 250 Гц.
Создайте datastore сигнала, чтобы получить доступ к данным в файлах. Этот пример принимает, что набор данных хранился в вашей временной директории под QTDatabaseECGData
папка. Если дело обстоит не так, изменение путь к данным в коде ниже. Задайте имена переменной сигнала, которые вы хотите считать из каждого файла с помощью SignalVariableNames
параметр.
sds = signalDatastore(datasetFolder,'SignalVariableNames',["ecgSignal","signalRegionLabels"])
sds = signalDatastore with properties: Files:{ '/tmp/QTDataset/ecg1.mat'; '/tmp/QTDataset/ecg10.mat'; '/tmp/QTDataset/ecg100.mat' ... and 207 more } AlternateFileSystemRoots: [0×0 string] ReadSize: 1 SignalVariableNames: ["ecgSignal" "signalRegionLabels"]
Datastore возвращает двухэлементный массив ячеек с сигналом ECG, и таблица области помечает каждый раз, когда вы вызываете read
функция. Используйте preview
функция datastore, чтобы видеть, что содержимое первого файла является 225,000 выборки длинный сигнал ECG и таблица, содержащая 3 385 меток области.
data = preview(sds)
data=2×1 cell array
{225000×1 double}
{ 3385×2 table }
Посмотрите на первые несколько строк таблицы меток области и заметьте, что каждая строка содержит предельные индексы области и значение класса области (P, T, или QRS).
head(data{2})
ans=8×2 table
ROILimits Value
__________ _____
83 117 P
130 153 QRS
201 246 T
285 319 P
332 357 QRS
412 457 T
477 507 P
524 547 QRS
Визуализируйте метки для первых 1 000 выборок с помощью displayWaveformLabels
функция помощника.
displayWaveformLabels(data,true,1000)
Обычная процедура классификации машинного обучения следующая:
Разделите базу данных на обучение и тестирующий наборы данных.
Обучите сеть с помощью обучающего набора данных.
Используйте обучивший сеть, чтобы сделать предсказания на наборе данных тестирования.
Сеть обучена с 70% данных и протестирована с остающимися 30%.
Для восстанавливаемых результатов, сброс генератор случайных чисел. Используйте dividerand
функция, чтобы заставить случайные индексы переставлять файлы и subset
функция signalDatastore
разделить данные на хранилища данных обучения и тестирования.
rng default
[trainIdx,~,testIdx] = dividerand(numel(sds.Files),0.7,0,0.3);
trainDs = subset(sds,trainIdx);
testDs = subset(sds,testIdx);
В этой проблеме сегментации вход к сети LSTM является сигналом ECG, и выход является последовательностью или маской меток с той же длиной как входной сигнал. Сетевая задача состоит в том, чтобы пометить каждую выборку сигнала именем области, которой это принадлежит. Поэтому необходимо преобразовать метки области на наборе данных к последовательностям, содержащим одну метку на выборку сигнала. Используйте преобразованный datastore и roi2mask
функция помощника, чтобы преобразовать метки. Предварительно просмотрите преобразованный datastore, чтобы заметить, что он возвращает сигнальный вектор и вектор метки равных длин. Постройте первые 1 000 элементов категориального вектора маски. roi2mask
функция добавляет категорию меток, n/a
, пометить выборки, которые не принадлежат никакой необходимой области.
trainDs = transform(trainDs, @roi2mask); testDs = transform(testDs, @roi2mask); transformedData = preview(trainDs)
transformedData=1×2 cell array
{224993×1 double} {224993×1 categorical}
plot(transformedData{2}(1:1000))
Передача очень длинных входных сигналов в сеть LSTM может привести к ухудшению производительности оценки и чрезмерному использованию памяти. Чтобы избежать этих эффектов, повредите сигналы ECG и их соответствующие маски метки с помощью преобразованного datastore и resizeData
функция помощника. Функция помощника создает как можно больше сегментов с 5000 выборками и отбрасывает остающиеся выборки. Предварительный просмотр выхода преобразованного datastore показывает, что первый сигнал ECG и его маска метки повреждаются в сегменты с 5000 выборками. Обратите внимание на то, что предварительный просмотр преобразованного datastore только показывает первые 8 элементов в противном случае floor(224993/5000)
= 44 массива ячеек элемента, которые закончились бы, если бы мы вызвали datastore read
функция.
trainDs = transform(trainDs,@resizeData); testDs = transform(testDs,@resizeData); preview(trainDs)
ans=8×2 cell array
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
{1×5000 double} {1×5000 categorical}
Во-первых, обучите сеть LSTM с помощью необработанных сигналов ECG от обучающего набора данных.
Задайте сетевую архитектуру перед обучением. Задайте sequenceInputLayer
из размера 1, чтобы принять одномерные временные ряды. Задайте слой LSTM с 'sequence'
режим вывода, чтобы обеспечить классификацию для каждой выборки в сигнале. Используйте 200 скрытых узлов в оптимальной производительности. Задайте fullyConnectedLayer
с выходным размером 4, один для каждого из классов формы волны. Добавьте softmaxLayer
и classificationLayer
выводить предполагаемые метки.
layers = [ ... sequenceInputLayer(1) lstmLayer(200,'OutputMode','sequence') fullyConnectedLayer(4) softmaxLayer classificationLayer];
Выберите опции для учебного процесса, которые гарантируют хорошую производительность сети. Обратитесь к trainingOptions
документация для описания каждого параметра.
options = trainingOptions('adam', ... 'MaxEpochs',10, ... 'MiniBatchSize',50, ... 'InitialLearnRate',0.01, ... 'LearnRateDropPeriod',3, ... 'LearnRateSchedule','piecewise', ... 'GradientThreshold',1, ... 'Plots','training-progress',... 'shuffle','every-epoch',... 'Verbose',0,... 'DispatchInBackground',true);
Поскольку целый обучающий набор данных умещается в памяти, возможно использовать tall
функция datastore, чтобы преобразовать данные параллельно, если Parallel Computing Toolbox™ доступен, и затем собирает его в рабочую область. Обучение нейронной сети является итеративным. В каждой итерации datastore считывает данные из файлов и преобразовывает данные прежде, чем обновить сетевые коэффициенты. Если совпадения данных в память о вашем компьютере, импортируя данные в рабочую область включают более быстрое обучение, потому что данные считаны и преобразованы только однажды. Обратите внимание на то, что, если данные не умещаются в памяти, вы должны, чтобы передать datastore в учебную функцию, и преобразования выполняются в каждую учебную эпоху.
Создайте длинные массивы и для наборов обучающих данных и для наборов тестов. В зависимости от вашей системы количество рабочих в параллельном пуле, который создает MATLAB, может отличаться.
tallTrainSet = tall(trainDs);
Starting parallel pool (parpool) using the 'local' profile ... Connected to the parallel pool (number of workers: 8).
tallTestSet = tall(testDs);
Теперь вызовите gather
функция длинных массивов, чтобы вычислить преобразования по набору данных в целом и получить массивы ячеек с учебными и тестовыми сигналами и метками.
trainData = gather(tallTrainSet);
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: 0% complete Evaluation 0% complete
- Pass 1 of 1: Completed in 28 sec Evaluation completed in 28 sec
trainData(1,:)
ans=1×2 cell array
{1×5000 double} {1×5000 categorical}
testData = gather(tallTestSet);
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 9.8 sec Evaluation completed in 10 sec
Используйте trainNetwork
команда, чтобы обучить сеть LSTM. Из-за большого размера набора данных, этот процесс может занять несколько минут. Если ваша машина имеет графический процессор и Parallel Computing Toolbox™, то MATLAB автоматически использует графический процессор в обучении. В противном случае это использует центральный процессор.
Учебные подграфики точности и потери в фигуре отслеживают процесс обучения через все итерации. Используя необработанные данные сигнала, сеть правильно классифицирует приблизительно 77% выборок как принадлежащий волне P, комплексу QRS, волне T или N/A.
net = trainNetwork(trainData(:,1),trainData(:,2),layers,options);
Классифицируйте данные о тестировании с помощью обученной сети LSTM и classify
команда. Задайте мини-пакетный размер 50, чтобы совпадать с опциями обучения.
predTest = classify(net,testData(:,1),'MiniBatchSize',50);
Матрица беспорядка обеспечивает интуитивное, и информативное означает визуализировать производительность классификации. Используйте confusionchart
команда, чтобы вычислить полную точность классификации для предсказаний данных о тестировании. Для каждого входа преобразуйте массив ячеек категориальных меток к вектору-строке. Задайте нормированное столбцом отображение, чтобы просмотреть результаты как проценты выборок для каждого класса.
confusionchart([predTest{:}],[testData{:,2}],'Normalization','column-normalized');
Используя необработанный сигнал ECG, как введено к сети, только приблизительно 60% выборок T-волны, 40% выборок P-волны и 60% QRS-комплексных выборок были правильны. Чтобы улучшать производительность, примените некоторое знание характеристик сигнала ECG до входа к нейронной сети для глубокого обучения, например, базовая линия, блуждающая вызванный дыхательным движением пациента.
Три морфологии удара занимает различные диапазоны частот. Спектр комплекса QRS обычно имеет центральную частоту приблизительно 10-25 Гц, и ее компоненты лежат ниже 40 Гц. P и T-волны происходят на еще более низких частотах: P-волновые-компоненты ниже 20 Гц, и T-волновые-компоненты ниже 10 Гц [5].
Базовое блуждание является низкочастотным (<0,5 Гц) колебание, вызванное движением дыхания пациента. Это колебание независимо от морфологии удара и не предоставляет значимую информацию [6].
Спроектируйте полосовой фильтр с частотным диапазоном полосы пропускания [0.5, 40] Гц, чтобы удалить блуждание и любой высокочастотный шум. Удаление этих компонентов улучшает обучение LSTM, потому что сеть не изучает несоответствующие функции. Используйте cellfun
на высоких массивах ячейки данных, чтобы отфильтровать набор данных параллельно.
% Bandpass filter design hFilt = designfilt('bandpassiir', 'StopbandFrequency1',0.4215,'PassbandFrequency1', 0.5, ... 'PassbandFrequency2',40,'StopbandFrequency2',53.345,... 'StopbandAttenuation1',60,'PassbandRipple',0.1,'StopbandAttenuation2',60,... 'SampleRate',250,'DesignMethod','ellip'); % Create tall arrays from the transformed datastores and filter the signals tallTrainSet = tall(trainDs); tallTestSet = tall(testDs); filteredTrainSignals = gather(cellfun(@(x)filter(hFilt,x),tallTrainSet(:,1),'UniformOutput',false));
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: 0% complete Evaluation 0% complete
- Pass 1 of 1: Completed in 19 sec Evaluation completed in 19 sec
trainLabels = gather(tallTrainSet(:,2));
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 19 sec Evaluation completed in 19 sec
filteredTestSignals = gather(cellfun(@(x)filter(hFilt,x),tallTestSet(:,1),'UniformOutput',false));
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 8.5 sec Evaluation completed in 8.6 sec
testLabels = gather(tallTestSet(:,2));
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 8.3 sec Evaluation completed in 8.5 sec
Постройте сырые данные и отфильтрованные сигналы для типичного случая.
trainData = gather(tallTrainSet);
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 19 sec Evaluation completed in 19 sec
figure subplot(2,1,1) plot(trainData{95,1}(2001:3000)) title('Raw') grid subplot(2,1,2) plot(filteredTrainSignals{95}(2001:3000)) title('Filtered') grid
Обучите сеть LSTM на отфильтрованных сигналах ECG с помощью той же сетевой архитектуры как прежде.
filteredNet = trainNetwork(filteredTrainSignals,trainLabels,layers,options);
Предварительная обработка сигналов улучшает учебную точность до лучше, чем 80%.
Классифицируйте предварительно обработанные тестовые данные с обновленной сетью LSTM.
predFilteredTest = classify(filteredNet,filteredTestSignals,'MiniBatchSize',50);
Визуализируйте производительность классификации как матрицу беспорядка.
figure confusionchart([predFilteredTest{:}],[testLabels{:}],'Normalization','column-normalized');
Простая предварительная обработка улучшает классификацию T-волн приблизительно на 15% и QRS-комплексную классификацию на 10%. Точность классификации P-волн не изменилась.
Общий подход для успешной классификации данных timeseries должен извлечь функции частоты времени и накормить ими сеть вместо исходных данных. Сеть затем изучает шаблоны через время и частоту одновременно [7].
Synchrosqueezed преобразование Фурье (FSST) вычисляет спектр частоты для каждой выборки сигнала, таким образом, это идеально для проблемы сегментации под рукой, где мы должны обеспечить то же разрешение времени как исходные сигналы. Используйте fsst
функция, чтобы смотреть преобразование одного из учебных сигналов. Задайте окно Кайзера длины 128, чтобы обеспечить соответствующее разрешение частоты.
data = preview(trainDs);
figure
fsst(data{1,1},250,kaiser(128),'yaxis')
Вычислите FSST каждого сигнала в обучающем наборе данных по частотному диапазону интереса, [0.5, 40] Гц. Обработайте действительные и мнимые части FSST как отдельные функции и подайте оба компонента в сеть. Кроме того, стандартизируйте учебные функции путем вычитания среднего значения и деления на стандартное отклонение. Используйте преобразованный datastore, extractFSSTFeatures
функция помощника и tall
функция, чтобы обработать данные параллельно.
fsstTrainDs = transform(trainDs,@(x)extractFSSTFeatures(x,250)); fsstTallTrainSet = tall(fsstTrainDs); fsstTrainData = gather(fsstTallTrainSet);
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: 0% complete Evaluation 0% complete
- Pass 1 of 1: Completed in 2 min 35 sec Evaluation completed in 2 min 35 sec
Повторите эту процедуру для данных о тестировании.
fsstTTestDs = transform(testDs,@(x)extractFSSTFeatures(x,250)); fsstTallTestSet = tall(fsstTTestDs); fsstTestData = gather(fsstTallTestSet);
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 1 min 23 sec Evaluation completed in 1 min 23 sec
Измените архитектуру LSTM так, чтобы сеть приняла спектр частоты для каждой выборки вместо одного значения. Смотрите размер FSST, чтобы видеть количество частот.
size(fsstTrainData{1,1})
ans = 1×2
40 5000
Задайте sequenceInputLayer
из 40 входных функций. Сохраните остальную часть сетевых параметров неизменной.
layers = [ ... sequenceInputLayer(40) lstmLayer(200,'OutputMode','sequence') fullyConnectedLayer(4) softmaxLayer classificationLayer];
Обучите обновленную сеть LSTM с преобразованным набором данных.
fsstNet = trainNetwork(fsstTrainData(:,1),fsstTrainData(:,2),layers,options);
Использование функций частоты времени улучшает учебную точность, которая теперь превышает 90%.
Используя обновленную сеть LSTM и извлеченные функции FSST, классифицируйте данные о тестировании.
predFsstTest = classify(fsstNet,fsstTestData(:,1),'MiniBatchSize',50);
Визуализируйте производительность классификации как матрицу беспорядка.
confusionchart([predFsstTest{:}],[fsstTestData{:,2}],'Normalization','column-normalized');
Используя частоту времени представление улучшает классификацию T-волн приблизительно на 20%, классификацию P-волн приблизительно на 40% и QRS-комплексную классификацию на 30%, когда по сравнению с необработанными данными заканчивается.
Используйте displayWaveformLabels
функция помощника, чтобы сравнить сетевое предсказание с основной истиной помечает для одного сигнала ECG.
testData = gather(tall(testDs));
Evaluating tall expression using the Parallel Pool 'local': - Pass 1 of 1: Completed in 8.3 sec Evaluation completed in 8.5 sec
figure subplot(2,1,1) displayWaveformLabels({testData{1,1}(3000:4000),testData{1,2}(3000:4000)},false) title('Ground Truth') subplot(2,1,2) displayWaveformLabels({testData{1,1}(3000:4000),predFsstTest{1}(3000:4000)},false) title('Predicted')
Этот пример показал, как предварительная обработка сигнала и частотно-временной анализ могут улучшать производительность сегментации формы волны LSTM. Фильтрация полосы пропускания и основанный на Фурье synchrosqueezing приводят к среднему улучшению через все выходные классы от 55% приблизительно до 85%.
[1] Макшарри, Патрик Э., и др. "Динамическая модель для генерации синтетических сигналов электрокардиограммы". IEEE® Transactions на Биоинженерии. Издание 50, № 3, 2003, стр 289–294.
[2] Laguna, Пабло, Рэймон Джейне и Пере Каминаль. "Автоматическое обнаружение контуров волны в мультиведущих сигналах ECG: Валидация с базой данных CSE". Компьютеры и Биомедицинское Исследование. Издание 27, № 1, 1994, стр 45–60.
[3] Голдбергер, Ари Л., Луис А. Н. Амараль, Леон Гласс, Джеффри М. Гаусдорф, Plamen Ch. Иванов, Роджер Г. Марк, Джозеф Э. Митус, Джордж Б. Муди, Чанг-Канг Пенг и Х. Юджин Стэнли. "PhysioBank, PhysioToolkit и PhysioNet: Компоненты Нового Ресурса Исследования для Комплексных Физиологических Сигналов". Циркуляция. Издание 101, № 23, 2000, стр e215–e220. [Циркуляция Электронные Страницы; http://circ.ahajournals.org/content/101/23/e215.full].
[4] Laguna, Пабло, Роджер Г. Марк, Ари Л. Голдбергер и Джордж Б. Муди. "База данных для Оценки Алгоритмов для Измерения QT и Других Интервалов Формы волны в ECG". Компьютеры в Кардиологии. Vol.24, 1997, стр 673–676.
[5] Sörnmo, Леиф и Пабло Лагуна. "Электрокардиограмма (ECG) обработка сигналов". Вайли Энкиклопедия Биоинженерии, 2006.
[6] Колер, B-U., Карстен Хенниг и Райнхольд Орглмайстер. "Принципы обнаружения программного обеспечения QRS". Разработка IEEE в Журнале Медицины и Биологии. Издание 21, № 1, 2002, стр 42–57.
[7] Salamon, Джастин и Хуан Пабло Белло. "Глубокие сверточные нейронные сети и увеличение данных для экологической звуковой классификации". Буквы Обработки сигналов IEEE. Издание 24, № 3, 2017, стр 279–283.