Анализ ЭКГ

Основные задачи и этапы непрерывного анализа ЭКГ

При оценке состояния сердца и сердечно-сосудистой системы ключевую роль играет анализ ЭКГ, представляющей собой запись наблюдаемых на поверхности тела проекций объемных электрических процессов, происходящих в сердце. Электрокардиограмма несет информацию как о текущем состоянии сердечно-сосудистой системы, так и о патологических изменениях в самом сердце12 3.

К настоящему времени проблема автоматической обработки и анализа электрокардиосигнала сложилась в отдельное направление, разветвленное на множество более частных задач, связанных с различными применениями и аспектами исследования сердечно-сосудистой системы. Одной из таких задач является длительный непрерывный контроль ЭКС. представляющий собой одну из важнейших функций аппаратуры для кардиологического наблюдения, к которой относятся приборы и системы для прикроватного и амбулаторного мониторного наблюдения ЭКГ, телеметрического наблюдения за состоянием организма, автономного контроля сердечной деятельности, функциональных исследований сердечно-сосудистой системы и решения ряда других задач.

Общим свойством систем кардиологического наблюдения является необходимость автоматической обработки продолжительных непрерывных реализаций электрофизиологических сигналов, к наиболее важным среди которых относится электрокардиосигнал. При этом обработка должна выполняться либо в режиме реального времени, либо в ускоренном режиме (анализ накопленных данных в холтеровских системах). В системах кардиологического наблюдения решаются следующие основные задачи:

  • текущий контроль частоты сердечных сокращений (ЧСС), служащий одним из основных гемодинамических показателей сердечнососудистой системы;
  • распознавание нарушений сердечного ритма (аритмий), которые возникают вследствие нарушений в автоматизме и проводящей системе сердца;
  • анализ ишемических изменений ЭКГ (в частности, изменений ST-сегмента кардиоцикла), служащий индикатором нарушений коронарного кровообращения;
  • анализ вариабельности сердечного ритма, позволяющий получать оценки регуляторных функций вегетативной нервной системы.

Автоматический анализ ЭКГ в системах кардиологического наблюдения характеризуется рядом специфических особенностей, к которым, в первую очередь, относятся следующие:

  • сложность и разнообразие сигнала, что проявляется в наличии его индивидуальных особенностей у каждого пациента, связанных как с типом конституции человека, так и с присутствием тех или иных патологий в проводящей системе сердца;
  • изменчивость сигнала во времени, которая может быть обусловлена как изменениями положения тела, гак и физиологическими феноменами — изменением общего состояния пациента, развитием патологий в сердечно-сосудистой системе и появлением нарушений ритма сердца;
  • наличие различных видов помех, неизбежно возникающих в процессе длительного наблюдения сигнала;
  • необходимость получения результатов анализа в реальном масштабе времени (или в ускоренном режиме, как в случае анализа холтеровских записей), что накладывает ограничения на вычислительную сложность используемых алгоритмов с точки зрения загрузки процессора и используемого объема оперативной памяти;
  • использование ограниченного числа отведений ЭКГ (в некоторых применениях — только одного отведения), что затрудняет получение точных оценок диагностических показателей, ввиду недостаточно полной информации, содержащейся во входном сигнале.

Как правило, процесс автоматического анализа ЭКС в системах кардиологического наблюдения включает в себя следующие основные этапы:

  • Предварительная обработка сигнала:
    • устранение или ослабление той части помех, которая поддается подавлению;
    • оценка уровня помех в сигнале с целью исключения фрагментов, дальнейший анализ которых невозможен вследствие их высокой зашумленности;
    • выделение информативных составляющих сигнала с целью обеспечения наилучших условий для последующих этапов обработки.
  • Обнаружение QRS-комплекса ЭКГ и измерение RR-интервала: обнаружение QRS-комплекса на фоне других элементов кардиоцикла ЭКГ;
    • определение границ комплекса и условной опорной точки, служащей для измерения RR-интервала.
  • Распознавание нарушений ритма, непосредственно угрожающих жизни пациента.
  • Анализ морфологии QRS-комплекса:
    • кластеризация QRS-комплексов по сходству их форм;
    • идентификация сформированных кластеров с точки зрения их принадлежности к фоновому ритму или к какой-либо из категорий патологии.
  • Распознавание нарушений сердечного ритма:
    • распознавание аритмических событий;
    • формирование заключений о характере ритма сердца.
  • Анализ ишемических изменений ЭКГ:
    • измерение смещения ST-сегмента относительно изоэлектрической линии;
    •  анализ типа формы ST сегмента.
  • Анализ вариабельности сердечного ритма:
    • формирование последовательности NN-интервалов (т.е. интервалов между смежными QRS-комплексами фонового ритма);
    • расчет оценок временных и частотных показателей ВСР.

Характеристики электрокардиосигнала и помех

Электрокардиограмма представляет собой запись электрических потенциалов, снимаемых с поверхности тела человека и отражающих электрические процессы, происходящие в сердце на протяжении всех фаз сердечного цикла. На рис. 1 приведен вид типичного цикла ЭКГ здорового человека в норме. На электрокардиограмме различают следующие основные компоненты, соответствующие различным фазам сердечного цикла:

  • Р-зубец (P-волна), соответствует фазе деполяризации (сокращения) предсердий;
  • QRS-комплекс (желудочковый комплекс), соответствует фазе деполяризации (сокращения) желудочков;
  • ST-сегмент, соответствует периоду между окончанием деполяризации и началом рсполяризации желудочков;
  • Т-зубец (Т-волна), соответствует фазе реполяризации желудочков.

Измерение амплитудных и временных параметров перечисленных компонентов кардиоцикла ЭКГ (расшифровка ЭКГ) лежит в основе методов автоматической интерпретации ЭКГ, задачей которой является выявление патологических изменений сердечной мышцы и проводящей системы сердца.

Типичный цикл ЭКГ здорового человека в норме
Рис. 1. Типичный цикл ЭКГ здорового человека в норме

В то же время, в приборах и системах кардиологического наблюдения наибольшее внимание уделяется анализу сердечного ритма и его нарушений или отслеживанию динамических изменений в электрокардиограмме в ходе ее продолжительного наблюдения. Это предполагает непрерывную регистрацию и анализ электрокардиосигнала на протяжении большого числа сердечных циклов. При этом основой автоматического анализа служит измерение интервалов времени между последовательными сокращениями желудочков сердца (RR-интервалов) и анализ характера формы (морфологии) основных волн кардиоцикла.

Нарушения сердечного ритма (аритмии) вызываются патологическими изменениями автоматизма и проводимости сердца.

Еще одной важнейшей задачей анализа электрокардиосигнала при кардиологическом наблюдении является выявление ишемических изменений ЭКГ, которые проявляются в трансформациях формы ST-сегмента (участка кардиоцикла между окончанием S-зубца и началом Т-зубца).

Во многих системах и приборах кардиологического наблюдения обеспечивается возможность анализа вариабельности сердечного ритма (ВСР). Анализ ВСР заключается в оценке статистических и частотных характеристик сигнала, образуемого последовательностью величин смежных RR-интервалов между желудочковыми комплексами, относящимися к фоновому ритму (NN-интервалов). Предполагается, что вычисляемые при этом показатели несут информацию о регуляторных функциях вегетативной нервной системы. .

В условиях продолжительного непрерывного съема электрокардиосигнала, характерных для систем кардиологического наблюдения, неизбежно появление в сигнале помех (артефактов), затрудняющих анализ сигнала и измерение его параметров. Различают следующие основные виды помех:

  • Сетевая наводка (рис. 2, а) представляет собой синусоидальный сигнал с частотой 50 Гц (или 60 Гц в США и Канаде).
  • Контактный шум электродов (рис. 2, б) — это преходящая помеха, вызываемая потерей контакта между электродом и кожей. Контактный шум электродов имеет вид случайным образом возникающего быстрого смещения (скачка) базовой линии.
  • Артефакты, связанные с движением (рис. 2, в) — кратковременное изменение базовой линии, вызванное изменением сопротивления электрод кожа, связанным с движением электрода.
  • Артефакты от мышечных сокращений (миографические шумы, рис. 2, г). Движения мышц пациента вызывают случайные широкополосные колебания, имеющие потенциал порядка одного милливольта.
  • Смещение базовой линии и модуляция амплитуды ЭКГ, вызываемые дыханием, приблизительно могут быть представлены как добавление к ЭКГ синусоидальной составляющей с частотой дыхания.
  • Аппаратный шум, генерируемый электронными устройствами съема и усиления сигналов. Вызывается насыщением входного усилителя. При этом информация об ЭКГ полностью теряется.
  • Электрохирургический шум может быть представлен синусоидой большой амплитуды с частотами приблизительно между 100 кГц и 1 МГц. При этом информация об ЭКГ полностью теряется.

Многообразие помех и их непредсказуемость существенно затрудняют автоматический анализ электрокардиосигнала. В связи с этим алгоритмы обработки ЭКС должны как включать процедуры подавления помех, так и предусматривать специальные меры, обеспечивающие их устойчивость к шумам.

Примеры некоторых видов помех в ЭКГ
Рис. 2. Примеры некоторых видов помех в ЭКГ

Предварительная цифровая фильтрация ЭКГ

На этапе предварительной обработки ЭКС решаются следующие основные задачи:

  • устранение из сигнала той части помех, которая поддается подавлению;
  • выявление фрагментов сигнала, где его автоматический анализ затруднен из-за чрезмерно высокого уровня помех;
  • выделение информативных составляющих сигнала с целью обеспечения наилучших условий для последующих этапов обработки.

Наиболее часто для подавления помех используются методы цифровой фильтрации. Однако применение цифровых фильтров затрудняется тем, что помехи могут иметь спектральный состав, в значительной степени перекрывающийся со спектром полезного сигнала. Поскольку спектр случайной помехи обычно носит широкополосный характер, цифровая фильтрация способна эффективно удалять лишь те составляющие помехи, спектр которых лежит за пределами полосы частот, соответствующей электрокардиосигналу.

Обычно предварительная фильтрация ЭКС включает следующие основные этапы4:

  • устранение сетевой наводки;
  • анализ уровня помех;
  • устранение высокочастотных помех;
  • устранение низкочастотных помех.

Сетевая наводка 50 или 60 Гц носит характер детерминированного сигнала и наиболее легко поддается устранению. Для подавления сетевой наводки обычно применяют режекторные фильтры, частотная характеристика которых имеет провал на частоте помехи. При этом наиболее целесообразно использование адаптивной цифровой фильтрации, которая способна подавлять помеху с частотой 50 Гц, не затрагивая в то же время полезные составляющие сигнала на частотах, близких к 50 Гц.

Основная задача этапа анализа уровня помех — выявление участков сигнала, на которых дальнейший анализ невозможен из-за чрезмерно высокого уровня шумов. Распознавание помех выполняется как по простейшим критериям (выход сигнала за пределы динамического диапазона АЦП, высокоамплитудные импульсы или скачки сигнала), так и путем текущего контроля дисперсии сигнала или относительного содержания в нем высокочастотных составляющих.

Для устранения высокочастотных помех из электрокардиосигнала (в основном — это высокочастотные составляющие миографической помехи) наиболее часто используются фильтры нижних частот, амплитудно-частотная характеристика которых имеет частоту среза 20…50 Гц. Недостатком таких фильтров является невозможность устранения тех составляющих помехи, спектр которых лежит в пределах полосы пропускания фильтров. Кроме цифровой фильтрации нижних частот предлагаются альтернативные методы, которые в значительной степени избавлены от указанного недостатка и основаны на следующих подходах:

  • согласованная фильтрация;
  • адаптивная фильтрация;
  • использование вейвлет-преобразования;
  • разложение сигнала на главные компоненты;
  • нейронные сети.

Низкочастотная помеха (дрейф изоэлектрической линии) может быть устранена при помощи цифрового фильтра верхних частот с частотой среза 1…3 Гц. Однако это приводит к существенному искажению как формы QRS-комплекса, так и ST-сегмента ЭКГ. Поэтому в случаях, когда одной из функций является анализ ишемических изменений ЭКГ, такой способ устранения помехи неприемлем. Для решения этой задачи наиболее часто прибегают к кусочно-полиноминальной аппроксимации с использованием сплайнов.

На рис. 3 и 4 приведены примеры соответственно фильтров нижних и верхних частот, используемых для предварительной обработки

ЭКГ. Оба фильтра являются КИХ фильтрами и рассчитаны на частоту дискретизации сигнала равную 250 Гц. На этих же рисунках показаны примеры фильтрации фрагментов ЭКГ с помощью данных фильтров.

Слева: амплитудно-частотная и импульсная характеристики ФНЧ; справа: пример фильтрации ЭКГ
Рис. Слева: амплитудно-частотная и импульсная характеристики ФНЧ; справа: пример фильтрации ЭКГ
Слева: амплитудно-частотная и импульсная характеристики ФВЧ; справа: пример фильтрации ЭКГ
Рис. Слева: амплитудно-частотная и импульсная характеристики ФВЧ; справа: пример фильтрации ЭКГ

Как отмечалось выше, наиболее эффективным средством подавления сетевой наводки промышленной частоты 50 Гц является адаптивная цифровая фильтрация56. Это обусловлено тем, что адаптивные фильтры при сравнительной простоте их реализации обеспечивают высокую добротность подавления синусоидальной составляющей заданной частоты. Важным достоинством этих фильтров с точки зрения анализа ЭКС является то, что они обладают способностью настраиваться на сетевую помеху и компенсировать ее, не влияя на близкие по частоте составляющие полезного сигнала. Рассмотрим принцип построения такого фильтра.

Пусть последовательности отсчетов на входе и выходе фильтра (компенсатора помех) обозначены соответственно

x(1),х(2)……x(i), … и у(1),у(2)…..y(i)…..

тогда любая i-я точка выходной последовательности может быть получена из соотношения

y(i) = x(i) — a(i), i =1,2,3….. (3.1)

где a(1),а(2), …. a(i), … — рассчитанные фильтром отсчеты компенсирующего сигнала.

Известно тригонометрическое равенство

sin(θ + δ) = 2sin θ cos δ – sin(θ – δ).

Умножим его на некоторую константу А:

A sin(θ + δ) = 2 cos δ (A sin θ) — 4sin(θ – δ).        (3.2)

Теперь будем считать, что θ текущее значение аргумента синусоидальной функции в точке взятия отсчета a(i), δ — угловое расстояние между соседними отсчетами (по отношению к периоду компенсируемой помехи), а А — амплитуда этой синусоиды. Тогда (3.2) можно переписать в виде

a(i + 1) = 2cosδa(i) — a(i — 1),               (3.3)

где a(i — 1),a(i),a(i + 1) — последовательные отсчеты синусоиды, отстоящие друг от друга на угол δ. Из (3.3) следует, что по двум известным соседним отсчетам синусоиды a(i — 1 ),a(i) всегда можно предсказать значение следующего за ними отсчета a(i + 1). Это обстоятельство и положено в основу алгоритма подавления синусоидальной помехи.

Пусть для момента времени, соответствующего i-му отсчету входного сигнала x(i), известны значение предыдущего отсчета x(i – 1), а также a(i – 1) и a(i), равные ординатам соответствующих точек синусоидальной помехи. Тогда текущий и предшествующий отсчеты выходного сигнала могут быть найдены по формулам

y(i) = x(i) – a(i);

y(i- 1) = х(i – 1) -a(i – 1).

Рассмотрим разность

Δy = y(i)  – y(i- 1) = |x(i) – a(i)| – |х(i – 1) -a(i – 1)|    (3.4)

Оказывается, что чем точнее определены отсчеты компенсирующего сигнала a(i  – 1) и a(i , тем ближе к нулю становится эта разность для каждого значения i, т. е. величина Δу используется в качестве критерия оценки полной выходной мощности фильтра. Поэтому подстройка фильтра под помеху заключается в том, что на каждом шаге фильтрации (для каждого следующего i) оценивается знак выражения (3.4) и в зависимости от того, положителен он или отрицателен, значение a(i) соответственно наращивается или уменьшается на фиксированную величину Δa так, чтобы скорректировать Δу в сторону нуля. Полученное в результате значение a(i) используется для вычисления выходного отсчета по формуле (3.1). Далее с использованием (3.3) вычисляется предполагаемое значение отсчета компенсирующего сигнала a(i + 1), после чего вся процедура повторяется для следующего входного отсчета, и т.д.

Скорость сходимости, устойчивость, ширина полосы подавления и точность настройки фильтра определяются тремя обстоятельствами: шагом адаптации Δa, постоянной соsδ и точностью вычисления промежуточных величин. Выбор константы cos δ выполняется из условия: δ = 360°fсп/fд, где fсп = 50 Гц — частота сетевой помехи, a dд -частота дискретизации сигнала. Например, при fд = 250 Гц получим:

δ = 360°fсп/fд = 360° • 50/250 = 72°,

cos δ = cos 72° = 0,309.

Наибольшее влияние на качество фильтрации оказывает выбор шага адаптации компенсирующего сигнала Δa. На рис. 5 показаны примеры работы фильтра при разных значениях этого параметра. Можно видеть, что при больших значениях Δa фильтр настраивается очень быстро, но заметно искажает полезный сигнал. С уменьшением До искажения сигнала становятся меньше, но растет время настройки.

На рис. 6 показан пример последовательной фильтрации фрагмента электрокардиосигнала с использованием цепочки цифровых фильтров:

  • фильтра нижних частот;
  • фильтра верхних частот;
  • адаптивного фильтра нижних частот.
Устранение сетевой наводки из электрокардиосигнала при разных значениях шага адаптации Де.
Рис. 5. Устранение сетевой наводки из электрокардиосигнала при разных значениях шага адаптации Δa.
Пример последовательной фильтрации фрагмента электрокардиосигнала с использованием цепочки цифровых фильтров
Рис. 6. Пример последовательной фильтрации фрагмента электрокардиосигнала с использованием цепочки цифровых фильтров

Обнаружение опасных аритмий по спектральному описанию электрокардиосигнала

В последнее время все больше внимания уделяется анализу биомедицинских сигналов в частотной области. Особый интерес представляет спектральный анализ ЭКС при нарушениях ритма, приводящих к появлению трепетания и мерцания желудочков сердца. Сигнал при этом приобретает волнообразную форму и имеет довольно характерное спектральное описание. Автоматически эти волны довольно сложно идентифицировать, в то же время задача надежного обнаружения опасных аритмий остается первоочередной для систем кардиологического наблюдения.

Анализ сигналов в частотной области как правило предполагает оценку мощности основных гармонических компонент, а также выявление некоторых особенностей распределения спектральной плотности мощности (СПМ)7. В задачах распознавания биомедицинских сигналов это часто приводит к выбору содержательных признаков, имеющих явно различные для разных классов характеристики (например, наличие узкополосных компонент спектра, выраженность основного пика СПМ и др.)8910.

Исследуя распределения объектов в пространстве выбранных признаков, можно решить задачу классификации сигналов, реализуя линейные решающие функции или конструируя более сложные поверхности решений. Однако часто оптимизация вида решающих правил не приводит к требуемому качеству распознавания сигналов. Это возникает из-за того, что выбранные признаки, отобранные подчас интуитивно или опытным путем, не несут существенной для задач классификации информации. Система признаков может оказаться неполной, а включение в решающие правила некоторых из них может даже снизить эффективность распознающей системы. Действительно, выбор множества признаков, обеспечивающих необходимое качество классификации, представляет собой трудную задачу11. В работах121314  решить ее предлагается путем формирования классификационного описания, рассматривающего полный ряд спектральных коэффициентов как упорядоченный набора параметров. Это приводит к возможности построения в пространстве выбранных признаков эффективной разделяющей гиперплоскости1516.

Пусть исходный сигнал x(t) представляет собой последовательность N отсчетов хn, n = 0, …. (N – 1), заданных с частотой дискретизации fd = 1/h, Гц. Тогда на интервале времени Т = Nh может быть получена выборочная оценка СПМ17:

G\left(f_{k}\right)=\frac{2 h}{N}\left|\sum_{n=0}^{N-1} x_{n} \exp \left(j-\frac{2 \pi k n}{N}\right)\right|^{2}

Далее следует анализ полученной оценки СПМ в области значений fk =k/T =k/Nh, k = 0, 1,2…..N/2, что соответствует ограничению спектра до fc = fd/2- (частоты Найквиста).

Учитывая известный эффект просачивания энергии через боковые лепестки спектрального окна, который возникает в результате усечения исследуемого процесса прямоугольным окном шириной Т, необходимо предварительно провести сглаживание концов исходной реализации. Для этого применяют специальные взвешивающие функции, например косинусоидальную, что позволяет уменьшить величину смещения выборочной оценки спектра.

Полученную оценку СПМ можно использовать для формирования спектральных признаков и построения решающих правил. Однако предварительно необходимо ограничить анализируемую область частот, что уменьшит размерность исходного спектрального описания и упростит процедуру распознавания сигналов. Для этого можно оценить расхождение выборочных описаний сигналов, сформировав обучающие выборки данных для заданных классов (w1 и w2) и вычислив для них оценки СПМ. В качестве меры близости между w1 и w2 можно использовать критерий, основанный на анализе среднего межгруппового расстояния:

\rho^{(k)}\left(\omega_{1}, \omega_{2}\right)=\left(1 / n_{1} n_{2}\right) \cdot \sum_{\mathbf{G}_{q}^{(k)} \in \omega_{1}} \sum_{\mathbf{G}_{j}^{(k)} \in \omega_{2}} \rho\left(\mathbf{G}_{q}^{(k)}, \mathbf{G}_{j}^{(k)}\right)

где w1 и w2 – число объектов, образующих классы w1 и w2; k шаг итерации. Gq(k), Gj(k) — спектральное представление объектов q и j для классов w1 и w2, используемое на к-м шаге (для к частотных составляющих).

Рассмотренный подход использовался при обработке ЭКС для решения задачи распознавания опасного нарушения ритма (фибрилляции желудочков сердца) (класс w1) на фоне других аритмий, нс представляющих опасности для жизни больного (нормальный ритм с одиночными экстрасистолами, групповая экстрасистолия, другие нарушения работы предсердий и желудочков сердца) (класс w2)1819.

В качестве примера на рис. 7 приведены реализации ЭКС и соответствующие выборочные оценки нормированного спектра для желудочковой фибрилляции ЖФ (а, б) и фонового ритма (ФР) (в, г). Фрагменты сигналов представлены записями ЭКГ длительностью 2 с,дискретизованными с частотой 250 Гц (50 ед. АЦП соответствуют 1 мВ).

Реализации ЭКС и соответствующие выборочные оценки нормированного спектра: а), б) для желудочковой фибрилляции ЖФ; в), г) для фонового ритма ФР
Рис. 7. Реализации ЭКС и соответствующие выборочные оценки нормированного спектра: а), б) для желудочковой фибрилляции ЖФ; в), г) для фонового ритма ФР
Зависимость межгруппового расстояния от числа выборок СПМ к: а) интегральная кривая; б) вклад в p(k)( w1, w2) каждого k-ro признака
Рис. 8. Зависимость межгруппового расстояния от числа выборок СПМ к: а) интегральная кривая; б) вклад в p(k)( w1, w2) каждого k-ro признака

На рис.8 приведены графики зависимости межгруппового расстояния от числа выборок СПМ к, включенных в спектральное описание сигналов. На рис. 8, а показана интегральная кривая  p(k)( w1, w2), а на рис. 8, б — зависимость расстояния p1(k)( w1, w2), вычисленного отдельно для каждого спектрального признака.

Зависимости получены по результатам обработки 60 записей ЭКС, специально отобранных из электрокардиографической базы MIT-BIH20.

Анализ полученных зависимостей позволил выделить эффективную область анализа СПМ с верхней частотой, равной fmах = 14 Гц (к = 28), поскольку более высокие частоты, как видно из рисунков, не вносят существенного вклада в оцениваемую меру расхождения классов.

Далее применяется процедура группировки спектральных коэффициентов сглаживанием их значений в границах заданных областей. Известно, что выборочный спектр, вычисляемый для независимых частот fk, подчиняется закону χ2 распределения с двумя степенями свободы. Получаемая для каждой частоты величина не является состоятельной оценкой СПМ, так как ее стандартное отклонение сравнимо по величине со средним значением G(fk) [3.8]. Используя сглаживание по / независимым частотам, перейдем к χ2-распределению с 2l степенями свободы, что приводит к уменьшению величины случайной ошибки в √1 раз. Таким образом, процедура объединения смежных частот способствует уменьшению дисперсии оценки СПМ, т. е. приводит к более устойчивым спектральным признакам.

Формирование набора спектральных признаков можно осуществить одним из предлагаемых способов. Первый способ основан на использовании периодограммной оценки Даньелла, формирующей новый ряд спектральных отсчетов путем усреднения (2Р + 1) соседних частот:

G_{D}\left(f_{i}\right)=\frac{1}{2P+1}\sum_{n=i-P}^{i+P} G(f_{n})

При этом снижается размерность частотного описания сигналов, а также повышается устойчивость спектральных оценок. Учитывая возможное ухудшение разрешения соседних частот, можно использовать перекрывающиеся сегменты, что обеспечит более детальный анализ сигнала в частотной области |3.8|. По в этом случае корреляция соседних оценок станет отличной от нуля, а число оцениваемых признаков практически удвоится. Эффективность выбранного спектрального окна и способа формирования отсчетов СПМ может быть также оценена по критерию максимизации расстояния между заданными классами сигналов (3.5).

Второй способ предусматривает автоматическое формирование эффективных частотных областей путем нахождения оптимальных границ для группировки смежных спектральных коэффициентов. Вначале задается максимально допустимое число интервалов L, в пределах каждого из которых объединяются точечные оценки СПМ. Это значение определяет и число признаков, которые могут быть использованы при построении дискриминантных функций. Затем оцениваются суммарные площади под кривой СИЛА:

G_{L p}=\frac{1}{N h} \sum_{k=m_{p}}^{m_{p+1}} G\left(f_{k}\right), р =1…..L

Перебор вариантов расположения границ областей (mp, mP+1) в каждой из L частотных областей позволяет получить достаточно большой набор векторов {G(L)}, из которых наилучший может быть определен по максимальному значению критерия (3.5). Для сокращения времени анализа может быть использована процедура градиентного поиска экстремума целевой функции. Аналогичные действия выполняются и на последующих шагах обработки при последовательном уменьшении числа интервалов (L, L – 1, L – 2, …). В результате формируется упорядоченный набор векторов {G(L),G(L-1), G(L-2), …} которые сравниваются друг с другом по значениям функции (3.5), и выбирается окончательный вариант спектрального описания G = (С1,…..Gp…..GL).

На рис. 9 представлены зависимости межгруппового расстояния рi( w1, w2) от целочисленного аргумента i, связанного с местоположением границ частотных областей. График построен для двух областей частот (L = 2), и максимум функции  рi( w1, w2) i = 10 приходится на последовательность границ 0,244; 6,832; 14,64 Гц. Все частоты кратны 0,244 Гц, что связано с дополнением нулями последовательности х(n), n = 0…..N – 1 до N = 1024 отсчетов. Аналогичные вычисления для L = 3, 4, 5 привели к следующим результатам:

  • 0,244; 2,44; 6,1; 14,64 Гц при L = 3;
  • 0,244 ; 2,44; 5,124; 10,004; 14,64 Гц при L = 4;
  • 0,244; 2,44; 5,124; 10,004; 12,444; 14,64 Гц при L = 5.
Зависимость межгруппового расстояния рi( w1, w2) от местоположения границ частотных областей x(i)
Рис. 9. Зависимость межгруппового расстояния рi( w1, w2) от местоположения границ частотных областей x(i)

Видно, что при переходе к L = 5 четвертая область (10,004; 14,64 Гц) образовала два новых кластера, остальные области остались прежними. В результате получены наборы признаков, сформированные в виде суммарных значений СПМ по L частотным областям.

Последний этап исследований связан с распознаванием классов ЭКС w1 и w2.  Для построения решающих функций в пространстве спектральных признаков использовался метод нахождения линейного дискриминанта Фишера. Сравнительный анализ полученных спектральных описаний проводился по критерию оптимизации J.

На рис. 10 в виде гистограмм представлены распределения проекций объектов на вектор W (у = WTG(L)) для сглаженных оценок СПМ в области частот 0. 15 Гц. Слева (рис. 10, а) показан результат, полученный с использованием периодограммы Даньелла для L = 15 отсчетов частоты при ширине окна Δf = 0,976 Гц. Справа (рис. 10, б) приведены распределения, полученные при формировании L = 4 частотных областей спектра по критерию максимизации межгрупповых расстояний. Здесь класс w1 представлен реализациями ЖФ, а w2 — сигналами ФР.

Гистограммы распределений проекций объектов на вектор W для классов w1 и w2 использовано окно Даньелла с Δf = 0,976 Гц; б) проведена оптимизация частотных областей
Рис. 10. Гистограммы распределений проекций объектов на вектор W для классов w1 и w2 использовано окно Даньелла с Δf = 0,976 Гц; б) проведена оптимизация частотных областей

Анализ полученных данных показывает, что объекты в обоих случаях образуют достаточно компактные группы, и классификацию сигналов можно эффективно осуществить, используя один из предложенных способов построения пространства спектральных признаков. Важно отметить, что, применяя оптимизацию частотных областей, можно существенно сократить число используемых признаков. В результате задача обнаружения ЖФ была решена путем формирования четырех частотных областей СПМ. Аппроксимация выборочных распределений проекций объектов у на ось решений W (0,433; 0,39; 0,397; 0,703) нормальным законом распределения позволила выбрать порога a= 0,415 и оценить теоретические ошибки распознавания классов w1 и w2 . При этом ошибка классификации даже при условии минимизации риска, связанного с пропуском ЖФ, не превышает 0,5 %. Результаты проверки эффективности обнаружения ЖФ на контрольной выборке данных, включающей более 200 коротких реализаций ЭКГ из базы

MIT-BIH, показали, что чувствительность метода Se = 0,98 и специфичность Sp = 0,9621. Эти результаты отличаются от теоретических оценок ввиду большего разнообразия представленных в контрольной выборке форм ЭКС, но уровень надежности классификации остается достаточно высоким.

В работе22 рассмотрен способ коррекции положения разделяющей функции путем нахождения дополнительного вектора по критерию Фишера J и перехода в двумерное признаковое пространство векторов W1, W2. В ходе экспериментов установлено, что введение дополнительного вектора W2 для спектрального описания ЭКС. использующего периодограммную оценку Даньелла (L = 15, Δf= 0,976 Гц), в 1,5 раза увеличило расстояние между ближайшими соседями рассматриваемых классов w1(ЖФ),   w2(ФР). Такой переход к отображению объектов на плоскости улучшает разделимость классов и позволяет использовать более сложные разделяющие поверхности, аппроксимируя последние кусочно-линейными функциями.

Алгоритм обнаружения QRS-комплекса ЭКГ

QRS-комплекс представляет собой наиболее выраженный элемент кардиоцикла ЭКГ, соответствующий фазе деполяризации (сокращения) желудочков сердца. Обнаружение QRS-комплекса служит основой всех последующих этапов анализа ЭКГ. Решение задачи обнаружения QRS-комплекса в системах длительного непрерывного контроля ЭКГ осложняется разнообразием и изменчивостью форм комплексов во времени, а также неизбежными при длительном наблюдении помехами.

Алгоритмы обнаружения QRS-комплекса решают следующие основные задачи:

  • обнаружение (локализация) QRS-комплекса в ЭКГ на фоне других элементов кардиоцикла;
  • определение условной опорной точки комплекса, служащей для измерения RR-интервала;
  • определение границ QRS-комплекса.

В литературе описываются многочисленные методы обнаружения QRS-комплекса232425, которые могут основываться на следующих подходах:

  • вычисление производных сигнала; линейная цифровая фильтрация;
  • адаптивная цифровая фильтрация;
  • согласованная фильтрация;
  • вейвлет-преобразование;
  • нейронные сети;
  • синтаксические методы.

Рассматриваемый далее алгоритм26 относится к наиболее распространенным решениям. Он основан на выделении QRS-комплексов

при помощи цепочки процедур цифровой фильтрации и последующем анализе полученного сигнала с использованием пороговых правил.

Поскольку большинство временных интервалов кардиоцикла имеют тенденцию к укорочению с ростом частоты сердечных сокращений, то это необходимо учитывать при разработке алгоритма обнаружения желудочковых комплексов, так как частота ритма сердца даже у одного и того же больного может меняться в ходе наблюдения в 2-3 раза. На рис. 11 изображен типичный кардиоцикл и указаны возможные пределы изменения временных интервалов между его основными характерными точками. Для учета возможных вариаций временных и амплитудных характеристик алгоритм обнаружения содержит ряд параметров, адаптирующихся к текущим средним значениям как частоты сердечных сокращений, так и размаха QRS-комплекса.

Типичный кардиоцикл и возможные пределы изменения временных интервалов между его основными характерными точками
Рис. 11. Типичный кардиоцикл и возможные пределы изменения временных интервалов между его основными характерными точками
Амплитудно-частотные характеристики цифровых фильтров, используемых для выделения желудочкового комплекса ЭКГ
Рис. 12. Амплитудно-частотные характеристики цифровых фильтров, используемых для выделения желудочкового комплекса ЭКГ

Алгоритму обнаружения QRS-комплекса ЭКГ обычно предшествует этап выделения QRS-комплекса, включающий дифференцирующий ЦФ, процедуру взятия модуля сигнала и сглаживающий фильтр. Если задаться значением частоты дискретизации fд 250 Гц, то качестве простого дифференцирующего ЦФ может, например, использоваться фильтр, определяемый выражением

y(n) == (1 /4)(x(n) – х(n- I) – х(n – 2) + х(n – 4) + х(n – 5) + х(n – 6)).

АЧХ этого фильтра приведена на рис. 12, а. Как видим, фильтр обладает свойствами дифференциатора в полосе частот от 0 до примерно 30 Гц. Дифференцирование позволяет выделить QRS-комплекс на фоне остальных волн ЭКГ, а также избавиться от низкочастотной помехи. Для удобства дальнейшего анализа сигнал берется по модулю и подвергается дополнительному сглаживанию, что повышает устойчивость обработки. Сглаживание может быть выполнено, например, с помощью такого фильтра:

у(n) = (1 /8)[x(n) – 2x(n – 1) – 2x(n – 2) + 2x(n – 3) + x(n – 4)].

АЧХ этого фильтра приведена на рис. 12, б.

На рис. 13 показаны все этапы выделения QRS-комплекса с использованием цифровых фильтров. На рис. 13, а приведен фрагмент исходной ЭКГ, содержащий три QRS-комплекса. В сигнале присутствуют сетевая помеха, дрейф изоэлектрической линии, а также незначительные высокочастотные шумы. На рис. 13, б показан этот же фрагмент после предварительной обработки фильтрами сетевой наводки и нижних частот, описанными в 3.3. Используемый ФНЧ имеет частоту среза около 30 Гц. Поэтому оставшиеся в сигнале частотные составляющие полностью укладываются в область дифференцирования используемого далее фильтра-дифференциатора. Можно видеть, что сигнал с выхода дифференцирующего ЦФ (рис. 13, в) представляет собой производную от вышеуказанного сигнала: участкам возрастания и убывания соответствуют положительные и отрицательные значения, а точкам экстремума — нулевые. В полученном сигнале амплитуда колебаний в области QRS-комплексов стала существенно выше по отношению к другим элементам ЭКГ, а дрейф изолинии исчез, что и было задачей дифференцирования ЭКГ. На рис. 13, гид показаны соответственно сигналы после процедур взятия модуля и дополнительного сглаживания с целью повышения устойчивости обработки.

Этапы выделения желудочкового комплекса ЭКГ
Рис. 13. Этапы выделения желудочкового комплекса ЭКГ
Иллюстрация принципа работы алгоритма обнаружения желудочкового комплекса ЭКГ
Рис. 14. Иллюстрация принципа работы алгоритма обнаружения желудочкового комплекса ЭКГ

На вход алгоритма обнаружения желудочковых комплексов ЭКГ поступает сигнал после описанных выше процедур выделения. Принцип работы алгоритма иллюстрируется на рис. 14.

Предполагается, что к моменту начала поиска очередного желудочкового комплекса известны:

  • опорная точка последнего обнаруженного желудочкового комплекса (точка в пределах комплекса, относительно которой вычисляется оценка RR-интервала);
  • текущее среднее значение размаха желудочковых комплексов;
  • текущие средние значения RR-интервалов.

На первом этапе работы алгоритма вычисляются значения тех параметров, которые зависят от текущих средних значений размаха желудочковых комплексов и RR-интервалов. Эго обеспечивает адаптивность алгоритма к текущим изменениям характеристик сигнала. Далее алгоритм пропускает промежуток времени от исходной точки до окончания зоны нечувствительности (адаптивно определяемой, как функция от текущего среднего значения RR-интервала), так как пред полагается, что на этом участке не может встретиться очередной комплекс. Начиная с отсчета, соответствующего моменту окончания зоны нечувствительности, выставляется линейно спадающий порог, начальное значение которого определяется с учетом среднего размаха желудочковых комплексов. Наклон и уровень окончания спада порога также являются адаптивно определяемыми параметрами. Такой подход позволяет снизить вероятность ложного срабатывания алгоритма из-за высокоамплитудных Т-зубцов ЭКГ. Обнаружение очередного QRS-комплекса происходит при пересечении сигналом текущего уровня порога.

После того как обнаружен очередной QRS-комплекс, выполняется определение условной опорной точки, относительно которой измеряется значение RR-интервала. Во многих алгоритмах в качестве такой точки используется точка максимума R-зубца. Однако этот способ измерения может давать существенную погрешность в случаях, когда R-зубец имеет расщепленную форму или когда форма QRS-комплекса подвергается трансформациям, связанным с дыханием. Более устойчивые результаты удается получить, используя некоторый интегральный метод (например, площадь под кривой QRS-комплекса). Аналогичная процедура используется и для поиска границ обнаруженного QRS-комплекса27.

Анализ форм ЭКГ и распознавание нарушений ритма

Одной из основных задач систем автоматического непрерывного контроля сердечной деятельности по ЭКГ является распознавание нарушений сердечного ритма (аритмий)2829. На рис. 15 показаны примеры ЭКГ при нормальном ритме сердца и при некоторых его нарушениях:

  • нормальный синусовый ритм (рис. 15, а) – последовательность идентичных кардиоциклов нормальной формы, следующих друг за другом с приблизительно одинаковыми интервалами;
  • суправентрикулярная (предсердная) экстрасистола (рис. 15, б) появление на фоне нормального синусового ритма преждевременного сокращения сердца нормальной формы, вызванного импульсом возбуждения из предсердия; за экстрасистолой, как правило, следует удлиненный интервал, называемый компенсаторной паузой;
  • выпадение желудочкового комплекса (рис. 15, в) — отсутствие QRS-комплекса после Р-зубца на фоне нормального синусового ритма;
  • желудочковая экстрасистола (рис. 15, г) появление на фоне нормального синусового ритма преждевременного сокращения сердца патологической формы, вызванного спонтанным импульсом возбуждения из одного из желудочков; за экстрасистолой, как правило, следует удлиненный интервал, называемый компенсаторной паузой;
  • групповая желудочковая экстрасистола (рис. 15, д) — появление на фоне нормального синусового ритма двух и более последовательных преждевременных сокращений сердца патологической формы, вызванных спонтанными импульсами возбуждения из желудочков.
Примеры ЭКГ при нормальном синусовом ритме (а) и при различных нарушениях ритма (б, в, г, д)
Рис. 15. Примеры ЭКГ при нормальном синусовом ритме (а) и при различных нарушениях ритма (б, в, г, д)

Как можно видеть из приведенных иллюстраций, нарушения ритма сердца проявляют себя на ЭКГ в виде изменения интервалов между смежными желудочковыми комплексами, а также морфологии элементов кардиоцикла, в частности QRS-комплекса. В связи с этим одной из задач анализа ЭКГ является идентификация принадлежности каждого обнаруженного желудочкового комплекса либо к категории комплексов фонового ритма (обычно называемой «норма»), либо к какому-либо виду патологических форм. Решение этой задачи служит основой для ряда этапов последующего анализа сигнала, в частности следующих:

  • распознавание аритмий (основано на одновременном анализе значений RR-интервалов и данных о характере форм желудочковых комплексов);
  • анализ вариабельности сердечного ритма (основан на анализе так называемой последовательности NN-интервалов, т.е. RR-интервалов между смежными QRS-комплексами, относящимися к фоновому ритму);
  • анализ ST-сегмента ЭКГ (выполняется только для QRS-комплексов фонового ритма).

Анализ морфологии QRS-комплексов затрудняется следующими факторами:

  • разнообразие форм ЭКГ у разных людей, а также невозможность выбора каких-либо «образцовых» форм сигнала;
  • изменчивость форм желудочковых комплексов во времени;
  • наличие помех, искажающих форму сигнала;
  • вариабельность формы кардиоцикла ЭКГ, связанная с изменением ЧСС, дыханием, изменением состояния или положения тела пациента.

Как правило, процедура анализа (классификации) форм QRS-комплексов включает следующие основные этапы:

  • измерение информативных признаков QRS-комплексов; формальное разбиение на кластеры по сходству морфологий; идентификация кластеров (присвоение каждому из кластеров
  • условного идентификатора «норма», «патология» или «неопределенность»).

В основе любой процедуры кластеризации QRS-комплексов всегда лежит какой-либо метод вычисления признаков формы, позволяющих оценить сходство или различие морфологий комплексов. Наибольшее распространение находят следующие три вида методов:

  • анализ параметров формы, вычисляемых во временной области; использование взаимной корреляционной функции;
  • анализ спектральных параметров.

При использовании параметров формы, вычисляемых во временной области, для каждого из QRS-комплексов анализируемой ЭКГ рассчитываются четыре признака формы: длительность Т с; размах A, мВ; смещение относительно нулевой линии мВ, и суммарная площадь волн Р, мВ х с (рис. 16). Расчет признаков выполняется по следующим формулам:

Т = К ΔТ;

A=\max _{k=1}^{K} x(k)-\min _{k=1}^{K} x(k)
S=\frac{1}{2}\left[\max _{k=1}^{K} x(k)+\min _{k=1}^{K} x(k)\right]-x(1)

P= \sum_{k=1}^{K} |x(k) - x(1)|

где К — число отсчетов сигнала в пределах QRS-комплекса; ΔТ = 1/fд — интервал дискретизации сигнала при частоте отсчетов fд; x(k) — отсчеты сигнала в пределах данного QRS-комплекса (k= 1,2…..К). При расчете признаков формы значение первого отсчета соответствующего комплекса условно принимается за уровень нулевой линии.

Иллюстрация вычисления признаков формы QRS-комплекса
Рис. 16. Иллюстрация вычисления признаков формы QRS-комплекса

Сопоставление форм двух QRS-комплексов между собой выполняется с помощью критерия близости Dij, рассчитываемого следующим образом:

D1ij = wTCT • |Тi – Tj| + wACA • |Ai – Aj| + wsCs • |Si – Sj| + wpCp • |Pi – Pj|, где признаки формы с индексами i и j относятся к соответствующим комплексам; величины CT, CA, CS и CP — это нормировочные коэффициенты, которые необходимы для приведения всех признаков формы к безразмерным единицам и для уравнивания их масштабов; wT, wA, wS и wP — весовые коэффициенты, служащие для учета значимости соответствующих параметров формы для классификации. Чем ближе между собой значения соответствующих признаков формы двух комплексов, тем меньшей будет величина показателя D1ij.

На рис. 17 показан пример вычисления признаков формы и показателя близости форм D1ij для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (а). На графике ЭКГ под каждым QRS-комплексом через дефис указаны его порядковый номер и условный класс формы. Здесь все QRS-комплексы, кроме комплексов с номерами 7 и 9, относятся к классу «1» (в данном случае он объединяет комплексы фонового ритма, относящиеся к типу «норма»), а комплексы с номерами 7 и 9 к классу «2» (в данном случае — это желудочковые экстрасистолы, относящиеся к типу «патология»). На рисунках б, в и г слева показаны сравниваемые комплексы, а справа — графически представлены вычисленные для них признаки формы (два верхних графика в правом столбце) и разности, рассчитанные для соответствующих признаков (нижний график в правом столбце). Значения всех признаков выражены в нормализованных единицах (н.е.), что соответствует операции нормирования их значений среднеквадратическими отклонениями этих признаков. Приведенная иллюстрация показывает, что для сходных по форме комплексов разности между соответствующими признаками (и следовательно значения D1ij) относительно невелики по сравнению с разностями для комплексов различной формы.

В случае использования корреляционного метода, для двух QRS-комплексов вычисляется взаимная корреляционная функция (ВКФ)30:

C_{x y}(m)=\left\{\begin{array}{ll}{\frac{1}{N \sigma_{x} \sigma_{y}}} & {\sum_{n=0}^{N-m-1} x_{i}(n+m)x_{j}(n), \quad 0 \leqslant m \leqslant M} \\ {\frac{1}{N \sigma_{x} \sigma_{y}}} & {\sum_{n=0}^{N-|m|-1} x_{i}(n) x_{j}(n+m), \quad-M \leqslant m<0}\end{array}\right.

где m — значение сдвига; M — максимальный сдвиг; xi(n) и xj(n) — отсчеты двух сравниваемых комплексов (i-го и j-ro); N — число отсчетов в каждом из образцов сигнала; σi и σj — соответствующие стандартные отклонения.

Максимальный сдвиг M выбирается с учетом возможной ошибки в определении условной опорной точки QRS-комплекса, которая может приводить к неточному совмещению сопоставляемых сигналов. Показатель близости форм двух QRS-комплексов (i-го и j-ro) рассчитывается как максимум из для всех заданных значений сдвига по формуле:

D _{H _{ij}}=\left[1-\max _{m=-M}^{m=M} C_{ij}(m)\right]

благодаря чему обеспечивается то условие, что чем лучше коррелированны два сигнала, тем ближе к нулю будет значение показателя DHij, а случаю наихудшей корреляции (Cij(m)=-1) будет соответствовать значение DHij = 2.

На рис. 18 приведены примеры вычисления ВКФ и показателя близости форм DHij  для QRS-комплексов одинаковых и различных классов форм, взятых из одной и той же реализации ЭКГ. Расчет ВКФ выполнялся для фрагментов сигнала длительностью 200 мс (50 отсчетов). Точками на графиках ВКФ отмечены значения функции в диапазоне величины сдвига m  от -M до +M (|M| = 5), которые использовались для определения соответствующих расстояний DHij  .

В соответствии с подходом, основанным на спектральных параметрах, вычисляется спектральное разложение QRS-комплекса (дискретное преобразование Фурье):

X(k)= \sum_{n=0}^{N-1} x(n)e^{-j \frac{2 \pi}{ N} n k}, k = 0…..N – 1,

где х(k) — отсчеты сигнала, а N — размер фрагмента сигнала, по которому рассчитывается ДПФ. Далее вычисляется амплитудный спектр:

B(k)=\frac{2}{N} \sqrt{\left\{\operatorname{Re}\left[X(k)\right]\right\}^{2}+\left\{\operatorname{Re}\left[X(k)\right]\right\}^{2}}

Пример вычисления признаков формы и показателя близости форм D1ij для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (a)
Рис. 17. Пример вычисления признаков формы и показателя близости форм D1ij для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (a)
Пример вычисления ВКФ и показателя близости форм DHij , для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (а)
Рис. 18. Пример вычисления ВКФ и показателя близости форм DHij , для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (а)
Пример вычисления амплитудных спектров (ЛС) и показателя близости форм D1Hij  для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (а)
Рис. 19. Пример вычисления амплитудных спектров (ЛС) и показателя близости форм D1Hij  для QRS-комплексов одинаковых (б и г) и различных (в) форм, взятых из одной и той же реализации ЭКГ (а)

Компоненты амплитудного спектра могут рассматриваться в качестве признаков формы QRS-комплекса. Так как большая часть энергии ЭКГ в пределах QRS-комплекса сосредоточена в полосе частот от 4 до 50 Гц31, для расчета показателя близости форм используются только спектральные компоненты, относящиеся к этому диапазону частот. Показатель близости форм QRS-комплексов (i-го и j-ro) рассчитывается по формуле32:

D_{1 \mathrm{H}_{ij}}=\frac{1}{\sqrt{\sigma_{i} \sigma_{j}}} \sum_{k=1}^{K}\left|B_{i}(k)-B_{j}(k)\right|

где Вi(к) и Bj(k) — спектральные компоненты с номером к для соответствующих комплексов; К — номер элемента разложения для частотного поддиапазона, включающего частоту 50 Гц; σi и σj соответствующие стандартные отклонения.

На рис. 19 приведены примеры вычисления амплитудных спектров и показателя близости форм D1Hij  для QRS-комплексов одинаковых и различных классов форм, взятых из одной и той же реализации ЭКГ. Расчет амплитудных спектров выполнялся для фрагментов сигнала длительностью 2256 мс (64 отсчета).

Каждый из трех рассмотренных методов обладает определенными преимуществами и недостатками по сравнению с двумя другими. Поэтому более высокое качество классификации может быть достигнуто в случае их совместного применения. В качестве интегрального признака, объединяющего показатели  DIij, DIIij, DIIIij       может быть предложено эвклидово расстояние от начала координат в трехмерном пространстве нормированных значений данных показателей:

D_{i j}=\sqrt{D_{IH_{ij}}^{2}+D_{IIH_{ij}}^{2}+D_{IIIH_{ij}}^{2}}

где DIHij, DIIHij, DIIIHij — значения соответствующих показателей, рассчитанных для некоторых двух комплексов с номерами i и j, нормированные среднеквадратичсскими значениями с целью выравнивания масштабов показателей. Исходя из способа расчета данных показателей, можно предположить, что величина Djj будет тем меньше, чем более сходны между собой сравниваемые комплексы, и наоборот.

На рис. 20 показано распределение значений показателя Dij в пространстве признаков DIHij, DIIHij, DIIIHij. Точками отмечены значения, полученные для QRS-комплексов, относящихся к одному-классу формы, а крестиками к разным классам. Распределение было получено с использованием верифицированной выборки реальных записей ЭКГ. включающей 50 фрагментов сигнала продолжительностью по 60 с. Из рисунка видно, что точки, соответствующие парам комплексов одного класса формы, концентрируются вблизи начала координат, а остальные распределяются на значительном удалении от него. В связи с этим может быть предложено решающее правило в виде порога, определяемого как поверхность сферы с центром в начале координат (см. рис. 20). При этом минимальная ошибка классификация достигается для значения радиуса сферы равного 4,633.

 Распределение значений показателя Dij в пространстве признаков DIHij, DIIHij, DIIIHij. Точками отмечены значения, полученные для QRS-комплексов, относящихся к одному классу формы, а крестиками — к разным классам
Рис. 20. Распределение значений показателя Dij в пространстве признаков DIHij, DIIHij, DIIIHij. Точками отмечены значения, полученные для QRS-комплексов, относящихся к одному классу формы, а крестиками — к разным классам

Конечной целью анализа сердечного ритма в системах кардиологического наблюдения является распознавание аритмий. Распознавание аритмий основывается на информации о временных интервалах между обнаруженными QRS-комплексами (RR-интервалах) и типах форм комплексов, получаемой на предыдущих стадиях анализа сигнала. Непосредственно распознавание нарушений ритма выполняется при помощи комплекса логических правил, описывающих аритмические события в терминах соотношений длительностей RR-интервалов и идентификаторов типа морфологии QRS-комплексов34.

Анализ ишемических изменений ЭКГ

Одной из наиболее важных задач автоматического анализа ЭКГ является выявление так называемых ишемических изменений кардиоцикла ЭКГ, которые отражают наличие недостаточного кровоснабжения сердечной мышцы3536. Ишемические изменения проявляются на ЭКГ в виде трансформации ST-сегмента кардиоцикла (промежутка между окончанием волны S и началом волны Т). В норме ST-сегмент выглядит как горизонтальный участок, расположенный приблизительно на уровне изоэлектрической линии (рис. 21, а). Ишемические изменения могут носить характер как смещения ST-сегмента вверх или вниз по отношению к изоэлектрической линии ЭКГ (значимым обычно считается смещение более чем на ±0,1 мВ), так и изменения формы ST-сегмента (различают горизонтальную, косонисходящую, косовосходящую, вогнутую и выпуклую формы). Примеры кардиоциклов с различными ишемическими изменениями ST-сегмента показаны на рис. 21, б-е.

К задачам алгоритмов анализа ST-сегмента относятся оценка величины смещения ST-сегмента относительно изоэлектрической линии ЭКГ и анализ формы ST-сегмента. При этом анализируются только кардиоциклы, относящиеся к фоновому ритму. Решение задачи анализа ST-сегмента ЭКГ затрудняется следующими факторами:

  • наличием низкочастотной помехи в виде дрейфа изоэлектрической линии ЭКГ;
  • необходимостью предварительного анализа формы QRS-комплекса с целью выбора для анализа ST-сегмента только комплексов фонового ритма;
  • сложностью выявления кратковременных преходящих изменений ST-сегмента;
  • сложностью точного определения границ QRS-комплекса и ST-сегмента при наличии в сигнале миографических помех.

Процедура анализа ST-сегмента, как правило, включает следующие основные этапы:

  • выбор из обнаруженных QRS-комплексов только комплексов, относящихся к фоновому ритму;
  • устранение низкочастотной помехи;
  • измерение смещения ST-сегмента;
  • анализ формы ST-сегмента.

Задача отбора только комплексов фонового ритма решается на этапе классификации форм QRS-комплексов.

Устранение низкочастотной помехи обычно осуществляется с помощью математического аппарата сплайн-интерполяции, так как использование фильтрации верхних частот вносит недопустимо большие искажения в форму самого ST-сегмента. При этом для каждого обнаруженного QRS-комплекса определяется условная точка на участке между волной Р и началом комплекса (на PQ-сегменте), значение сигнала в которой рассматривается как локальная оценка текущего уровня изоэлектрической линии. Найденные таким образом точки используются в качестве узлов сплайн-интерполяции, в результате которой формируется кусочно-непрерывная функция, состоящая из фрагментов полиномов третьей степени и служащая оценкой низкочастотной помехи в сигнале. Вычитая эту функцию из исходного сигнала, получают ЭКГ, очищенную от дрейфа изоэлектрической линии. Описанный процесс иллюстрируется на рис. 22. Па верхнем графике показан фрагмент ЭКГ, на котором наблюдается низкочастотная помеха. Кружками отмечены найденные узловые точки на PQ-сегменте, а жирной линией показан рассчитанный сплайн. Па нижнем графике приведен тот же фрагмент ЭКГ после вычитания сплайна из исходного сигнала. Можно видеть, что дрейф изоэлектрической линии практически полностью устранен.

Примеры различных форм ST-сегмента ЭКГ
Рис. 21. Примеры различных форм ST-сегмента ЭКГ
Иллюстрация процедуры устранения дрейфа изоэлектрической линии с использованием сплайн-ннтерполяции
Рис. 22. Иллюстрация процедуры устранения дрейфа изоэлектрической линии с использованием сплайн-ннтерполяции
Процедура вычисления смешения ST-сегмента ЭКГ
Рис. 23. Процедура вычисления смешения ST-сегмента ЭКГ

После устранения низкочастотной помехи для каждого анализируемого комплекса определяется так называемая точка .1, за которую принимается точка окончания QRS-комплекса. Смешение ST-сегмента относительно изоэлектрической линии измеряется на определенном расстоянии от точки J (наиболее часто — через 60-80 мс после нее). На рис. 23 иллюстрируется процедура вычисления смешения ST-сегмента.

Анализ характера формы ST-сегмента осуществляется с использованием логических правил по положению трех последовательных точек на ST-сегменте (например, по точке J, точке измерения значения смещения и точке начала волны Т).

АНАЛИЗ СИГНАЛОВ ВАРИАБЕЛЬНОСТИ СЕРДЕЧНОГО РИТМА, ЭЛЕКТРОЭНЦЕФАЛОГРАММ И СПИРОГРАММ В МЕДИЦИНСКИХ СИСТЕМАХ

Анализ вариабельности сердечного ритма и артериального давления

Анализ вариабельности сердечного ритма

Анализ ВСР представляет собой метод исследования статистических и частотных свойств сигнала, образуемою последовательностью интервалов времени между смежными сокращениями сердца, относящимися к фоновому ритму. Анализ ВСР позволяет получить количественные индикаторы активности различных отделов вегетативной нервной системы человека, что оказывается полезным для решения многих медицинских диагностических задач3738 39404142.

Интервалы времени между последовательными сокращениями сердца (RR-интервалы) обычно измеряются по ЭКГ как расстояния между вершинами R-зубцов QRS-комплексов. Если за период регистрации ЭКГ в сигнале наблюдались нарушения сердечного ритма или помехи, вызвавшие ошибки в определении RR-интервалов, то участки сигнала, соответствующие этим событиям, должны исключаться из процедуры анализа. В связи с этим при анализе ВСР в качестве входных данных используется не непосредственно последовательность RR-интервалов, а так называемая последовательность NN-интервалов, т.е. интервалов между смежными QRS-комплексами фоновою ритма.

Для анализа ВС.Р наиболее широкое применение находят частотные методы. При этом в качестве основных спектральных показателей используются суммарные мощности сигнала в трех диапазонах частот43: VLF (Very Low Frequency, от 0,003 до 0,04 Гц). LF (Low Frequency, от 0.04 до 0.15 Гц) и HF (High Frequency, от 0,15 до 0.4 Гц), а также ряд показателей, производных от них.

Так как анализируемая последовательность интервалов времени не представляет собой равномерно дискретизованный сигнал, традиционно используемый при частотном анализе, то к ней нельзя непосредственно применить общепринятые методы цифровой обработки. В связи с этим на первом этапе анализа из последовательности NN-интервалов получают равномерно дискретизованный сигнал (рекомендованная частота дискретизации составляет 4 Гц), который рассматривается как описание восстановленной функции управления сердечным ритмом. Именно исследование этой функции, порождающей последовательность регистрируемых RR-интервалов, является предметом анализа ВСР. Для получения равномерно дискретизованного сигнала наиболее часто используется сплайн-интерполяция.

Последовательность NN-интервалов (или RR-интсрвалов) принято графически отображать в виде ритмограммы (тахограммы), на которой значения интервалов изображаются в виде вертикальных линий, а по оси абсцисс откладываются порядковые номера соответствующих интервалов. На рис. 1 приведены примеры ритмограмм для нескольких различных видов сердечного ритма:

  • высокая активность в диапазоне I.F (рис. 1, а),
  • высокая активность в диапазонах LF и HF (рис. 1, б);
  • умеренная активность в диапазонах LF и HF (рис. 1, в);
  • высокая активность в диапазонах VLF, LF и HF (рис. 1, г).

Для получения оценок спектральных показателей вариабельности сердечного ритма используют как непараметрические (классические) методы, основанные на дискретном преобразовании Фурье, так и параметрические методы, в основе которых лежит модель авторегрессии скользящего среднего.

Процедура получения оценок спектральных показателей вариабельности сердечного ритма включает следующие основные этапы:

  1. определение имеющихся в последовательности RR-интервалов нарушений непрерывности, связанных с помехами или аритмиями;
  2. получение последовательности NN-интервалов;
  3. восстановление функции управления сердечным ритмом в виде последовательности равноотстоящих отсчетов;
  4. расчет оценки СПМ;
  5. расчет показателей ВСР.

Если в анализируемом фрагменте сигнала встречаются интервалы, связанные с нарушениями ритма, или интервалы, неверно определенные из-за помех, необходимо предварительно скорректировать сигнал. В противном случае в результирующих спектральных оценках могут появиться завышенные значения высокочастотных составляющих, так как искажения, вызванные указанными причинами, обычно проявляются в виде резких всплесков или впадин в последовательности NN-интервалов. В случае относительно коротких промежутков времени между оставшимися NN-интервалами дефектные значения не используют в качестве узлов интерполяции, но при этом сохраняют непрерывность сигнала во времени. Если же искажается достаточно продолжительный отрезок сигнала, то либо добавляют промежуточные узлы интерполяции, либо отказываются от анализа данной реализации сигнала.

Равномерно дискретизованный сигнал чаще всего получают с помощью математического аппарата сплайн-интерполяции. На рис. 2 иллюстрируется процедура перехода от последовательности NN-интервалов к равномерно дискретизованному сигналу. Верхний график (а) представляет собой фрагмент ритмограммы NN-интервалов, на среднем графике (б) показан процесс получения сплайна, проходящего через узловые точки, соответствующие значениям NN-интервалов, а нижний график (в) демонстрирует расчет по функции сплайна значений равномерно дискретизованного сигнала с частотой отсчетов 4 Гц (или с периодом 0,25 с).

Примеры ритмограмы при различных вилах сердечного ритма
Рис. 1. Примеры ритмограмы при различных вилах сердечного ритма
Иллюстрация перехода от последовательности NN-интервалов (а) к равномерно дискретизованному сигналу с частотой 4 Гц (в) с использованием сплайн-интерполяции (б, в)
Рис. 2. Иллюстрация перехода от последовательности NN-интервалов (а) к равномерно дискретизованному сигналу с частотой 4 Гц (в) с использованием сплайн-интерполяции (б, в)

На рис. 3, а, б соответственно показаны ритмограмма для фрагмента ЭКГ продолжительностью 5 мин и рассчитанный по ней равномерно дискретизованный сигнал.

На рис. 4 показаны полученные для реализации сигнала сердечного ритма, приведенной на рис. 3, оценки СПМ и значения рассчитанных по ним частотных параметров ВСР для двух методов спектрального анализа: метода Уэлча (а) и авторегрессионного (б).

Расчет показателей VLF, LF и HF осуществляется интегрированием СПМ с использованием формулы прямоугольников по соответствующим диапазонам частот:

\begin{aligned} \mathrm{VLF} &=\Delta f \sum_{k} p\left(f_{k}\right), \quad 0,003 \leqslant f_{k} \leqslant 0,04 \mathrm{\Gammau} \\ \mathrm{LF} &=\Delta f \sum_{k} p\left(f_{k}\right), \quad 0,04 \leqslant f_{k} \leqslant 0,15 \mathrm{\Gammau} \\ \mathrm{HF} &=\Delta f \sum_{k} p\left(f_{k}\right), \quad 0,15 \leqslant f_{k} \leqslant 0,4 \Gamma_{\mathrm{ll}} \end{aligned}

где Δf = fд/N (fд = 4 Гц — частота дискретизации, N — число входных значений ДПФ), а p(fk) — значение СПМ с индексом k.

 Ритмограмма фрагмента ЭКГ длительностью 5 минут (а) и интерполирующая кривая, полученная с использованием сплайнов (б)

Рис. 3. Ритмограмма фрагмента ЭКГ длительностью 5 минут (а) и интерполирующая кривая, полученная с использованием сплайнов (б)

Графики оценки СПМ и значения рассчитанных по ним частотных параметров ВСР для двух методов спектрального анализа: метода Уэлча (а) и авторегрессионного (б)
Рис. 4. Графики оценки СПМ и значения рассчитанных по ним частотных параметров ВСР для двух методов спектрального анализа: метода Уэлча (а) и авторегрессионного (б)

Можно видеть, что как вид графиков СПМ, так и значения рассчитанных разными методами параметров несколько различаются между собой, что объясняется случайным характером спектральных оценок, получаемых при использовании любого из методов. Необходимо также обратить внимание на единицы измерения спектральных показателей: СПМ измеряется в мс2/Гц, а мощность — в мс . Это связано с тем, что в качестве исходных данных для расчета спектров использовались длительности интервалов времени, измеряемые в миллисекундах.

Совместный анализ вариабельности сердечного ритма и артериального давления

С появлением технической возможности одновременной регистрации сигналов сердечного ритма и мгновенного кровяного давления широкое распространение получили математические методы, позволяющие выявить и количественно оценить взаимосвязи между этими двумя сигналами с использованием, в частности, взаимной корреляционной функции и взаимной спектральной плотности мощности сигналов4445464748.

На рис. 5 показаны фрагменты сигналов сердечного ритма (RR) и артериального давления (АД), которые были получены с использованием сплайн интерполяции из синхронно записанных последовательностей значений RR-интервалов ЭКГ и мгновенных значений артериального давления. В сигналах хорошо просматриваются сходные по своему характеру высокочастотные колебания.

Фрагменты двух синхронно записанных сигналов сердечного ритма (а) и артериального давления (б)
Рис. 5. Фрагменты двух синхронно записанных сигналов сердечного ритма (а) и артериального давления (б)
Пример графиком оценок СПМ сигналов сердечного ритма (а) и артериального давления (б), а также график ВСПМ этих сигналов (в)
Рис. 6. Пример графиком оценок СПМ сигналов сердечного ритма (а) и артериального давления (б), а также график ВСПМ этих сигналов (в)
Пример графиков оценок фазового спектра сигналов сердечного ритма (а) и артериального давления (б). а также график взаимного фазового спектра этих сигналов (в)
Рис. 7. Пример графиков оценок фазового спектра сигналов сердечного ритма (а) и артериального давления (б). а также график взаимного фазового спектра этих сигналов (в)

На рис. 6 приведен пример графиков оценок СПМ сигналов сердечного ритма и артериального давления, показанных на рис 5, а также график ВСПМ этих сигналов. Значения ВСПМ выражены в мс х мм рт.ст./Гц. На графике отчетливо просматриваются несколько пиков, свидетельствующих о существовании выраженных взаимосвязей между колебаниями сигналов в определенных частотных диапазонах. Видно, что ВСПМ имеет относительно высокие значения только на тех частотах, на которых в спектрах обоих сигналах присутствуют хорошо выраженные составляющие.

На рис. 7 приведен пример графиков оценок фазового спектра сигналов сердечного ритма и артериального давления, показанных на рис 5, а также график взаимного фазового спектра (ВФС) этих сигналов. Значения взаимного фазового спектра выражены в градусах. Из графиков видно, что ВФС принимает устойчивые значения на тех самых участках частотного диапазона, на которых наблюдаются относительно высокие значения ВСПМ В частности, это наиболее хорошо заметно в диапазоне частот от 0,29 до 0,36 Гц, соответствующем дыхательной активности человека. Например, на частоте 0,3 Гц фазовый сдвиг составляет примерно 120 градусов. Так как в данном случае в качестве первого из сигналов испол»»зовался сигнал сердечного ритма, а в качестве второго сигнал мгновенного артериального давления, то положительное значение фазы свидетельствует о запаздывании изменений ритма относительно изменений давления, что соответствует существующим предположениям о физиологических механизмах регуляции сердечного ритма49. Количественно это запаздывание можно оценить следующим образом:

\Delta T\left(f_{k}\right)=\frac{\angle P_{x y}^{\circ}\left(f_{k}\right)}{360^{\circ} f_{k}}=\frac{120^{\circ}}{360^{\circ} 0,3} \approx 1,1 \mathrm{c}

То есть на частоте дыхания колебания сердечного ритма здесь отстают от дыхательных волн примерно на 1 секунду.

Анализ ВСР на основе модели управления водителем сердечного ритма

Методы автоматического анализа ВСР находят применение при исследовании активности различных физиологических систем, управляющих работой сердца. Основным методом оценки влияния таких физиологических факторов, как артериальное давление, терморегуляция и дыхание, является спектральный анализ ВСР. В стандартных методиках принято измерять мощность сигнала в частотных диапазонах: VLF, LF, HF и проводить текущий анализ се изменения во времени. По результатам такого анализа можно получать физиологическую интерпретацию уровней регуляции сердечного ритма пациента, оперативно распознавая вид и степень тяжести его состояния, и даже прогнозировать степень риска внезапной смерти кардиологических больных.

Одним из важнейших этапов обработки рнтмограммы для задач спектрального анализа является непрерывное (или равномерное) представление ВСР. Дискретность сигнала сердечного ритма вызывает неопределенность поведения функции между отсчетами. Устранить подобную неопределенность позволяет достаточно большое число методов, к которым относятся интерполяция кусочно-полиномиальными функциями (нулевого, первого и третьего порядков), представление сигнала в виде счетных событий, интервальной последовательности и другие традиционные методы5051. Их особенность состоит в том, что в основном они восстанавливают функцию, представляющую собой огибающую зарегистрированной последовательности кардиоинтервалов. Существует и другой подход, который основан на использовании модели управления водителем сердечного ритма. В этом случае за исходную функцию, по которой оценивается частотный состав ВСР, принимается функция, управляющая частотой сердечных сокращений (ЧСС), и, как следствие, содержащая в себе наиболее полную информацию об уровнях регуляции вегетативной нервной системы (ВИС).

Изучая взаимосвязь функции автономной управляющей активности и последовательности сердечных событий, была создана математическая модель, называемая интегральным частотно-импульсным модулятором (IPFM integral pulse frequency modulator)52. Эта модель учитывает физиологическую природу и биофизические механизмы функционирования водителя сердечного ритма, что имеет большое значение для распознавания и интерпретации спектральных компонент функции управления ритмом сердца.

Электрофизиологические процессы зарождения импульсов возбуждения в синоатриальном (СЛ) узле и проведения их в систему миокарда находятся под влиянием ряда регуляторных факторов. Нормальное возбуждение сердца начинается в СЛ узле. Когда мембранный потенциал волокон узла достигает порогового значения, генерируется потенциал действия, который активирует другие узловые волокна и способствует распространению возбуждения в сердце. Существует, по меньшей мере, три фактора, определяющие частоту разрядов водителя сердечного ритма. К ним относят:

  1. крутизну диастолической деполяризации;
  2. уровень потенциала покоя;
  3. уровень порогового потенциала.

Влияние их на работу водителя ритма проиллюстрировано на рис. 8.

Изменение частоты появления мембранного потенциала водителя сердечного ритма под влиянием изменения порогового значения — А, уменьшения крутизны диастолической деполяризации — В, увеличения максимального диастолического потенциала — С
Рис. 8. Изменение частоты появления мембранного потенциала водителя сердечного ритма под влиянием изменения порогового значения — А, уменьшения крутизны диастолической деполяризации — В, увеличения максимального диастолического потенциала — С

При этом решающее значение имеет изменение скорости диастолической деполяризации, а изменение порогового потенциала и потенциала покоя оказывает незначительное влияние на ЧСС, которым при разработке модели можно пренебречь. Разные отделы ВНС оказывают противоположное действие на изменение ЧСС. Так, стимуляция симпатических волокон, которые иннервируют сердце, увеличивает частоту разрядов водителя ритма путем увеличения крутизны диастолической деполяризации, в то время как стимуляция парасимпатических волокон дает обратный эффект. Поэтому именно анализ изменения скорости деполяризации лежит в основе исследования различных влияний ВНС на сократительную работу сердца.

Модель управления водителем ритма

Согласно IPFM-модели, существует функция m(t), которая влияет на электрическую активность сердца. Эта управляющая функция представляет собой результат деятельности обеих ветвей BHC. Для оценки вегетативного влияния на работу сердца необходимо соотнести последовательность сердечных событий с этим управляющим сигналом m(t).

На рис. 9 представлена модель управления частотой водителя сердечного ритма посредством изменения крутизны диастолической деполяризации. Данная модель отражает лишь диастолическую фазу мембранного потенциала, гак как значение имеет только время появления потенциала действия.

 IPFM-модель автономного управления частотой водителя сердечного ритма
Рис. 9. IPFM-модель автономного управления частотой водителя сердечного ритма

Как следует из IPFM-модели, диастолическая деполяризация моделируется путем интегрирования входной функции m(t). Как только значение интеграла достигает фиксированного порога Т, соответствующего пороговому значению мембранного потенциала, генерируется импульс, и интегратор устанавливается на значение мембранного потенциала покоя.

Входная функция m(t) и выходная посылка импульсов p(i) математически связаны друг с другом следующим образом. Последовательность p(t) представляется как поток единичных импульсов, расположенных в моменты времени ti, i = 0, 1, 2,   k, в которые наблюдаются желудочковые комплексы электрокардиосигнала:

p(T)=\sum_{i}\delta (t-t_{i})  (4.1)

Моменты времени ti-1,ti определяющие два последовательных сердечных события в (4.1), связаны следующим интегральным соотношением:

\int_{t_{i}-1}^{t_{i}}\left|m(t)+x_{0}\right| d t=T               (4.2)

где Т — постоянная интегрирования, x0 — постоянная величина, обеспечивающая положительное значение интеграла на временном промежутке (ti-1 – ti).

Из выражения (4.2) следует, что данная модель позволяет получить математическую зависимость между временем появления сокращения сердца ti и входной функцией m(t). При этом принимаются следующие допущения и упрощения.

  • A. Соотношения между автономной активностью симпатического и парасимпатического отделов ВИС и крутизной диастолической деполяризации объединены в непрерывную управляющую функцию m(t), задаваемую на входе модели.
  • В. Изменение крутизны диастолической деполяризации является единственным решающим фактором изменения ЧСС.
  • В. В математическом описании отсутствует в явном виде моделирование рефрактерного периода, хотя для максимальной приближения IPFM-модели к своему физиологическому аналогу необходимо его учитывать, поскольку он влияет на формирование длительности кардиоцикла. Эта задача решается путем косвенного учета рефрактерного периода в параметрах модели: пороговом значении компаратора Т и постоянной составляющей x0.

Алгоритм восстановления функции управления m(t)

Восстановление функции m(t) основано на вычислении интегральной функции M(ti)     определяемой следующим образом. Пусть начальное время to = 0. тогда получаем эквивалентное (4.2) выражение:

M(t_{i})=\int_{0}^{t_{i}}\left|m(t)+x_{0}\right| d t=iT (4.3)

где ti, i= 1,2, …,n — моменты времени, связанные с появлением каждого i-го сокращения сердца.

Функция M(ti) имеет дискретный вид, поскольку определена лишь в моменты времени По последовательности M(ti), » = 1,2, …,N можно восстановить непрерывную функцию M(t), применяя линейную или кубическую интерполяцию, например, кусочно-полиномиальными функциями. Последующее дифференцирование M(t) позволяет получить непрерывное описание функции m(t). Производя равномерную дискретизацию восстановленного сигнала, можно решать задачу спектрального анализа сердечного ритма, ожидая объективности оценки спектральной плотности мощности (СПМ) и, следовательно, возможности распознавания функциональных состояний пациента.

На рис. 10 в графическом виде последовательно представлены этапы преобразования ритмограммы в непрерывную функцию m(t). Для наглядности все построения проиллюстрированы на модельном примере.

Первый этап обработки связан с анализом интервальной последовательности RRi, i= 1,2, …,N. Ритмограмма (рис. 10, а) преобразуется в нерегулярную последовательность длительностей RR-интервалов (рис. 10, б), для которой RR(G) = ti — ti-1

Далее осуществляется формирование двумерного массива, содержащего N точек с координатами (ti, M(ti)), i = 1, 2, …. N (рис. 10, в). Значение M(ti) в моменты времени ti, вычисляется как произведение постоянной величины Т на значение целочисленной переменной i.

Следующим шагом (рис. 10, в) является построение непрерывной кривой M(t) по опорным точкам (ti, M(ti)), реализуемое с помощью кубических сплайнов.

Завершающий этап заключается в дифференцировании функции M(t) и нахождении искомой непрерывной функции m(t) (рис. 10, г). Далее осуществляется взятие равномерных отсчетов с заданной частотой (4 Гц) и реализуется спектральный анализ полученного дискретного представления функции m(t).

Блок-схема алгоритма восстановления функции m(t) по ритмограмме приведена на рис. 11. Здесь используются следующие обозначения: RR, — длительность i-го RR-интервала, i = 1, 2…..N; N — размер массива RR-интервалов в отсчетах; t0— начальный момент времени; ti — момент времени появления i-го сокращения сердца; M(t0) — начальный отсчет массива M(ti), i = 1, 2, …. N.

Высокого качества интерполяции по заданным точкам M(ti) можно достичь с помощью глобального интерполяционного сплайна, представленного по базису из кубических В-функций53. Однако организация его вычислений осложняется необходимостью определения большого числа неизвестных, получаемых в процессе решения. Наилучшие результаты с точки зрения вычислительных затрат дает метод интерполяции кубическими сплайнами, в котором определение параметров кубического сплайна сводится к решению системы линейных уравнений, представленной достаточно простой трехдиагональной матрицей [4.10].

Реконструкция функции управления m(t) по ритмограмме
Рис. 10. Реконструкция функции управления m(t) по ритмограмме

а) ритмограмма; б) определение моментов времени появления событий t1 t2...; в)формирование массива Iti, M(ti)), i = 1, 2, .... n, и построение непрерывной функции M(t); г) воспроизведение функции m(t)

Блок-схема алгоритма воспроизведения функции управления водителем сердечного ритма
Рис. 11. Блок-схема алгоритма воспроизведения функции управления водителем сердечного ритма

В этом случае возможно организовать вычисление функции управления m(t) в ходе непрерывного наблюдения за динамикой сердечного ритма.

Результаты экспериментов

Проведен сравнительный анализ различных методов преобразования сигнала сердечного ритма к равномерной временной шкале545556. Исследования включали ряд экспериментов, выполненных на модельных и реальных сигналах.

Эксперименты на модельных сигналах проводились с целью оценки эффективности отображения спектрального состава сигналов с заранее заданными частотными свойствами. Модельные сигналы m1(t), m2(t) и m3(t) представленные с частотой дискретизации 250 Гц, включали одну, две и три гармонических компоненты с частотами f1= 0,10 Гц, f2 =0,03 Гц и f3 = 0,25 Гц, соответствующими диапазонам LF, VLF, и HF. С помощью алгоритма IPFM формировались последовательности единичных импульсов. Затем, используя конкретный способ преобразования ритмограммы, осуществлялся переход к равномерно дискретизованному сигналу и вычислялась оценка СПМ.

В качестве показателей эффективности были выбраны:

  • интенсивность утечки, SL (определяется как относительное значение суммарной мощности спектральных компонентов, находящихся за пределами изолирующих окон, центрированных относительно частот заданных гармонических функций);
  • число ложных спектральных компонент nα с амплитудами более чем α (1, 5 и 10) % от полной мощности сигнала.

Используя эти показатели, необходимо было определить метод непрерывного представления последовательности импульсов, наиболее точно воспроизводящий частотный состав исходного сигнала.

Сравнительный анализ рассмотренного метода и традиционно используемых в практических разработках способов преобразования сигнала сердечного ритма показал, что метод спектрального анализа, основанный на реконструкции функции управления, обеспечивает более высокую эффективность распознавания частотных составляющих m(t) в HF-, LF-, VLF-областях. Значения Sc при использовании данного метода составили 0,76 %, 0,69 % и 2,46 % для сигналов m1(t), m2(t) и m3(t) соответственно, а показатель nα= 6; 4; 3 для сигнала m3(t). Близкие результаты показал метод, основанный на кубической интерполяции неравномерной последовательности кардиоциклов, но значения Sc оказались выше и составили, соответственно, 0,92 %. 1,15% и 4,63 %.

В ходе экспериментов, выполненных на реальных записях ЭКГ, установлено, что переход от спектрального анализа функций, описывающих ритмограмму, к анализу функции управления ЧСС улучшает условия распознавания состояний, связанных с разным уровнем вегетативной регуляции. В качестве характеристики баланса симпатических, соответствующих LF-диапазону, и парасимпатических, соответствующих HF-диапазону, влияний, использовался показатель k = PLF/PHF. где PLF и PHF — значения суммарной мощности в диапазоне низких LF и высоких HF частот. Анализ более 60 реализаций ЭКГ длительностью 5 мин, полученных в клиниках Санкт-Петербурга и представляющих два класса функциональных состояний ССС, включал оценку показателя к двумя способами: с использованием традиционного метода построения интервальной последовательности (равномерной последовательности) и на основе рассмотренного метода воспроизведения функции m(t).

Результаты экспериментов показали, что критерий Фишера J, принимающий максимальное значение в условиях наилучшей классификации двух состояний, увеличивает значение в 2 раза при переходе к анализу функции автономного управления.

Footnotes

  1. Дюк В.. Эммануэль В. Информационные технологии в медико-биологических исследованиях. — СПб.: Питер, 2003.
  2. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  3. Рангайян Р. М. Анализ биомедицинских сигналов / Пер. с англ. A. П. Калиниченко; под род. А. II. Немирко. — М.: Физматлит, 2007. — 440 с.
  4. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  5. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  6. Уидроу Б., Стирнз С. Адаптивная обработка сигналов. — М.: Радио и связь, 1989. — 440 с.
  7. Марпл-мл. С. Л. Цифровой спектральный анализ и его приложения / Пер. с англ. — М.: Мир, 1990. — 584 с.
  8. Табак В. Я., Черныш А. М., Немирно А. П.. Манило Л. А. Динамика спектральных характеристик ЭКГ при развитии фибрилляции желудочков сердца // Анестезиология и реаниматология. 1980, № 1. С. 71 -74.
  9. Clayton R. Н., Murray A. Comparison of techniques for time-frequency analysis of the ECG during human ventricular fibrillation / IF.F. Proc, on Sci., Measurement and Technol. 1998, № 145. P. 301-306.
  10. Nemirko A. P., Manilo I.. A., Degtyareva I. N. Recognition of Heart Ventricular Fibrillation on the Basis of Analysis of the Electrocardiogram in the Frequency Domain // Pattern Recognition and Image Analysis. 2001. Vol. 11. № 2. P. 353-355.
  11. Ту Дж., Гонсалес Р. Принципы распознавания образов / Пер. с англ. —М.: Мир, 1978.   411 с.
  12. Манило Л. А. Упорядочение спектральных признаков по эмпирическим оценкам межгруппового расстояния в задачах классификации биосигналов // Изв. вузов России. Радиоэлектроника. 2006. Вып. 3. С. 20-29.
  13. Manilo L. A. Measurement and computerized analysis of signals in biomedical systems И Proc, of the V Intern. Symp. on test and measurement. 2003. Vol. 1. China, June 1-5. P. 867-870.
  14. Manilo L. A., Nemirko A. P. Forming the spectral signs when classifying (he electrocardiosignals within the frequency range // Pattern recognition and image analysis. 2005. Vol. 15, № 4. P. 668-671.
  15. Манило Л. А. Линейный дискриминант Фишера в задачах распознавания биосигналов по частотным свойствам // Докл. XII Всерос. конф. «Математические методы распознавания образов» — ММРО-12, Москва, 20-26 нояб. 2005. М : МАКС Пресс, 2005. С. 371-374.
  16. Манило Л. А. Построение решающих функций в пространстве спектральных признаков для систем кардиологического наблюдения // Изв. СПбГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2006. Вып. I. С. 13-21.
  17. Марпл-мл. С. Л. Цифровой спектральный анализ и его приложения / Пер. с англ. — М.: Мир, 1990. — 584 с.
  18. Манило Л. А. Упорядочение спектральных признаков по эмпирическим оценкам межгруппового расстояния в задачах классификации биосигналов // Изв. вузов России. Радиоэлектроника. 2006. Вып. 3. С. 20-29.
  19.  Manilo L. A., Nemirko A. P. Forming the spectral signs when classifying (he electrocardiosignals within the frequency range // Pattern recognition and image analysis. 2005. Vol. 15, № 4. P. 668-671.
  20. MIT-BIH Arrhythmia Database. Available from MIT-BIH database distribution // Massachusetts   institute   of   technology. URL: http://www.physionet.org/physiobank/database/niitdb; http://ecg.mit.edu.
  21. Manilo L. A., Nemirko A. P. Forming the spectral signs when classifying (he electrocardiosignals within the frequency range // Pattern recognition and image analysis. 2005. Vol. 15, № 4. P. 668-671.
  22. Manilo I.. A., Nemirko A. P. Recognition of biomedical signals based on their spectral description data analysis // Pattern Recognition and Image Analysis. 2016.
  23. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  24. Микрокомпьютерные медицинские системы: Проектирование и применение / Под ред. У. Томпкинса, Дж. Уэбстера; пер. с англ. — М.: Мир, 1983. – 544 с.
  25. Рангайян Р. М. Анализ биомедицинских сигналов / Пер. с англ. A. П. Калиниченко; под род. А. П. Немирко. — М.: Физматлит, 2007. — 440 с.
  26. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  27. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  28. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  29. Рангайян Р. М. Анализ биомедицинских сигналов / Пер. с англ. A. II. Калиниченко; под род. А. II. Немирко. — М.: Физматлит, 2007. — 440 с.
  30. Калиниченко А. Н. Оценка разделяющей способности методов классификации форм ЭКГ // Изв. СП6ГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2006. Вып. 1. С. 21 – 30.
  31. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  32. Калиниченко А. Н. Оценка разделяющей способности методов классификации форм ЭКГ // Изв. СП6ГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2006. Вып. 1. С. 21 30.
  33. Калиниченко А. Н. Оценка разделяющей способности методов классификации форм ЭКГ // Изв. СП6ГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2006. Вып. 1. С. 21 30.
  34. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  35. Кардиомониторы. Аппаратура непрерывного контроля ЭКГ / А. Л. Барановский, А. Н. Калиниченко, Л. А. Манило и др.; под ред. Л. Д. Барановского и Л. П. Немирко. — М.: Радио и связь, 1993. — 248 с.
  36. Рангайян Р. М. Анализ биомедицинских сигналов / Пер. с англ. A. II. Калиниченко; под род. А. II. Немирко. — М.: Физматлит, 2007. — 440 с.
  37. Калиниченко А. II., Юрьева О. Д. Влияние частоты дискретизации ЭКГ на точность вычисления спектральных параметров вариабельности сердечного ритма // Информационно-управляющие системы. 2008, № 2. С. 46 49.
  38. Калиниченко Л. Н., Коляденко М. И. Исследование алгоритмов оценки стационарности сердечного ритма // Изв. СПбГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2006. Вып. 2. С. 101-105.
  39. Калиниченко Л И. О точности и достоверности спектральных методов расчета показателей вариабельности сердечного ритма // Информационно-управляющие системы. 2007, № 6. С. 41-48.
  40. Михайлов В. М. Вариабельность ритма сердца. Опыт практического применения. — Иваново: НейроСофт, 2000. — 200 с.
  41. Рангайян Р. М. Анализ биомедицинских сигналов / Пер. с англ. А. II. Калиниченко; под ред. А. П. Немирно. — М.: Физматлит, 2007. — 440 с.
  42. Heart rate variability. Standards of measurements, physiological interpretation. and clinical use // Circulation. 1996. V. 93 (5). P. 1043-1065.
  43. Heart rate variability. Standards of measurements, physiological interpretation. and clinical use // Circulation. 1996. V. 93 (5). P. 1043-1065.
  44. Богачев М. И., Мамонтов О. В., Конради А. О., Ульяницкий У. Ц. Оценка спонтанного артериального барорефлекса методом совместного анализа показателей кратковременной изменчивости артериального давления и сердечного ритма. // Артериальная гипертензия. 2007. Т. 13, №1.
  45. Калиниченко Л И. О точности и достоверности спектральных методов расчета показателей вариабельности сердечного ритма // Информационно-управляющие системы. 2007, № 6. С. 41-48.
  46. Heart rate variability. Standards of measurements, physiological interpretation. and clinical use // Circulation. 1996. V. 93 (5). P. 1043-1065.
  47. La Rovere M. T., Pinna G. D., Raczak G. Baroreflex sensitivity: measurement and clinical implications // Ann. Noninvasive Elcctrocardioi. 2008 Apr. 13(2). P. 191-207.
  48. Sleight P., La Rovere M. T., Mortara A., Pinna G., Maestri R., Leuzzi S., Bianchini B., Tavazzi L.. Bernardi L. Physiology and pathophysiology of heart rate and blood pressure variability in humans. Is power spectral analysis largely an index of barorcflex gain? // Clinical Science. 1995. V. 88. P. 103-109.
  49. Манило Л. А., Родина Н.И. Новый подход к спектральному анализу вариабельности сердечного ритма // Известия СПбГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2003. Вып. I. С. 16-20.
  50. Guimaraes H. N.. Santos R. A. S. A comparative analysis of preprocessing techniques of cardiac event series for the study of heart rhythm variability using simulated signals// Braz. J. Biol. Res. 1998. Vol. 31(1). P. 421 430.
  51. Heart rate variability. Standards of measurements, physiological interpretation. and clinical use // Circulation. 1996. V. 93 (5). P. 1043-1065.
  52. Hyndman B. U’Mohn R. K. A model of the cardiac pacemaker and its use in decoding the information content of cardiac intervals Ц Automedica. 1975. Vol. 1. P. 239-252.
  53. Павлидис Т. Алгоритмы машинной графики и обработки изображений / Пер. с англ. — М.: Радио и связь, 1986. — 400 с.
  54. Манило Л. А., Родина Н.И. Новый подход к спектральному анализу вариабельности сердечного ритма // Известия СПбГЭТУ «ЛЭТИ». Сер. Биотехнические системы в медицине и экологии. 2003. Вып. I. С. 16-20.
  55. Manilo L. A.. Rodina N. I. Investigation of a Model of the Cardiac Rhythm Pacemaker Control for the Spectral Analysis of a Rhythmogram // Pattern Recognition and Image Analysis. 2001. V. 11, № 2. P. 342 344.
  56. Manilo L.A., Rodina N.I. Estimation of the Frequency Properties of Rhythmograms in Problems of Recognition of Physiological States // Pattern Recognition and Image Analysis. 2003. V. 13, № 2. P. 298-301.