Динамический фрактальный анализ вариабельности сердечного ритма

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

¥"-» * /—' Г"» * Г"» *

ника Е , , так и на краю облака Е „ , причем Е „

т! 7 А т2 7 А т2

больше, чем на стержне по его оси.

Показано, что распределение напряженности поля Е* по поверхности земли под облаком убывает к периферии, а распределение напряженности поля Е* по поверхности облака, наоборот, существенно возрастает к периферии. При этом

следует заметить, что полученные выводы определяются соотношениями размеров исследуемой системы и соответствующими допущениями при расчетах.

В результате исследований установлено, что в выделяемой области 0'т(х, у, г) с некоторой внутренней границей Г5 с учетом применения МИПБС при определении влияния краевых эффектов исследуемого объекта (например, У(х,у)) он не должен соприкасаться с внутренней границей Г5, т. к. это приводит к завышению результатов расчета по Г5.

список литературы

1. D'Alessandro, F. Electric field modelling of structures under thunderstorm conditions [Текст] / F. D'Alessandro, J.R. Gumley // Proc. of the 24th International Conf. on Lightning Protection. -Birmingham, Britain, 1998. -Р. 457-462.

2. Ait-Amar, S. A 3-D numerical model of negative lightning leader interception. Applications to the collection volume construction [Текст] / S. Ait-Amar, G. Berger // Proc. of the 27th International Conf. on Lightning Protection. -Avignon, France, 2004. -Р. 357-362.

3. Резинкина, М.М. Расчет трехмерных электрических полей в системах, содержащих тонкие проволоки [Текст] / М.М. Резинкина // Электричество. -2005. -№ 1. -С. 44-49.

4. Потапенко, А.Н. Метод инверсии для численного расчета распределенных систем типа «плоскость-проводник» [Текст] / А.Н. Потапенко, Е.А. Ка-нунникова, Т.А. Потапенко // Научно-технические ведомости СПбГПУ Сер. Информатика. Телекоммуникации. Управление. -2011. -№ 5. -С. 53-57.

5. Потапенко, А.Н. Численное моделирование электрических полей в системах «электрод - поверхность земли» для элементов молниезащит [Текст] / А.Н. Потапенко, Е.А. Канунникова, М.И. Дыльков //

Изв. вузов. Проблемы энергетики. -2008. -№ 11-12. -С. 72-78.

6. Rezinkina, M.M. Software for determinantion of 3D electrical fields distribution in the vicinity of special installations ans systems with lightning rods during thunderstorm [Текст] / M.M. Rezinkina // Proc. of the 24th International Conf. on Lightning Protection. -Birmingham, Britain, 1998. -Р. 924-928.

7. D'Alessandro, F. A 'Collection Volume Method' for the placement of air terminals for the protection of structures against lightning [Текст] / F. D'Alessandro, J.R. Gumley // J. of Electrostatics. -2001. -№ 50. -Р. 279-302.

8. Ait-Amar, S. Attractive Radius of Elevated Building [Текст] / S. Ait-Amar, G. Berger // Proc. of the 28th International Conf. on Lightning Protection. -Kanasawa, Japan, 2006. -Р. 602-607.

9. Потапенко, А.Н. Исследование распределенных элементов систем молниезащит на основе вычислительных экспериментов [Текст] / А.Н. Потапенко, А.И. Штифанов, Т.А. Потапенко // Изв. Самарского научного центра РАН. -2010. -Т 12. -№ 4 (3). -С. 591-595.

10. Мак-Кракен, Д. Численные методы и программирование на Фортране [Текст] / Д. Мак-Кракен, У Дорн. -М.: Мир, 1977. -584 с.

УДК 519.711.3

В.И. Антонов, А.И. Загайнов, Ву ван Куанг

ДИНАМИЧЕСКИМ ФРАКТАЛЬНЫЙ АНАЛИЗ ВАРИАБЕЛЬНОСТИ СЕРДЕЧНОГО РИТМА

В настоящей статье рассматриваются проблемы нелинейного исследования хаотических временных рядов вариабельности сердечного ритма (ВСР) - временных интервалов между последо-

вательными нормальными QRS-комплексами электрокардиограммы (норма-норма или ЫЫ-интервалов). Благодаря рекомендациям научных сообществ США, Европы [3], Японии, Китая и

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

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

Методика динамической оценки фрактальных компонент ВСР

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

х(?) = (а(?),а(? + тД?), ..., а(? + т(и - 1)Д?)).

Время задержки выбирается различными способами [6, 12]. Дальнейшее исследование временного ряда сводится к исследованию свойств построенного в и-мерном пространстве геометрического объекта, своим видом напоминающего фрактал. Основной характеристикой рассматриваемого объекта в случае его фрактальности

должна являться размерность, вычисленная с учетом вероятности нахождения точек в данной области аттрактора. Наиболее удобна для этих целей корреляционная размерность:

М

1П

D = lim-

2 ln

где М(е) - минимальное количество кубиков со стороной е, полностью покрывающее аттрактор; р. - вероятность посещения г-го кубика фазовой траекторией динамической системы.

При практической реализации корреляционной размерности с помощью метрического определения вероятности выражение (1) переписывается как

]П

D, = lim -1

ln

где C(r) = lim —- У 6(r - p(x, x.)) - корреля-

i, j = 1

ционный интеграл; 9(a) =

- функция

1, a>0 [ö, a < ö

Хевисайда; p - расстояние в и-мерном пространстве.

Для восстановленного по временному ряду аттрактора, состоящего из конечного числа точек, дается следующая оценка корреляционного интеграла:

C

i=ö j=i+1 m(m -1)/2

Алгоритм Грассбергера-Прокаччиа [2] пред-полает, что в случае выполнения соотношения C(r) ~ rDl , то есть ln C(r) ~ D2 ln r , корреляционную размерность можно оценить, получив наклон прямой логарифма корреляционного интеграла. Самый простой способ получения линейной зависимости по последовательности экспериментальных данных - метод наименьших квадратов. В данном случае можно самостоятельно выбирать диапазон значений расстояний r и рассматривать равномерные сетки. Исходя из метода наименьших квадратов, получим следующую систему для оценки корреляционной размерности:

' N N N

У ln Cr )ln r = D2 У ln2 r + ьу ln r

i=1 i=1 i=1

<

N N

У ln C (r) = D2 у ln rt + bN

_ i=1 i=1

где N - количество измерений корреляционного интеграла для различных расстояний ri (вычис-

4

ленных на равномерной сетке); Ь - константа аппроксимирующей корреляционный интеграл прямой у (г) = £>2 г + Ь .

Справедливость приведенного закона ограничена значениями г достаточно малыми по сравнению с размером аттрактора. Очевидно, при увеличении г до размеров аттрактора С(г)^1, а при уменьшении из-за конечности точек на аттракторе С(г)^0, и указанный степенной закон справедлив только в ограниченном диапазоне г (т. н. скейлинговом диапазоне), который может быть использован для определения размерности аттрактора. Этот диапазон необходимо установить на практике для рассматриваемых нами биомедицинских сигналов и менять в зависимости от типа сигнала (типа соответствующего аттрактора сигнала). Разработанные методы определения скейлингового диапазона можно найти в [6].

Следующий шаг - проверка полученной оценки фрактальной размерности. Для этой цели можно вычислить значения размерности восстановленного аттрактора в (п+1)-мерном пространстве вложения по описанному выше алгоритму и оценить значение модуля разности полученных значений \бп2+1 - Ди| = ап.

Если последовательность то можно го-

ворить о сходимости рассматриваемого процесса и, как следствие, фрактальных свойствах самого восстановленного аттрактора ВСР. Однако такие сообщения в научной литературе встречаются редко. Например, в [1] исследована возможность разделения пациентов с гипертрофической кардио-миопатией с помощью корреляционной размерности на две группы: с низким и высоким риском внезапной смерти. Авторы обследовали 16 паци-

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

В проведенных численных экспериментах для записей нормального синусового ритма Таким образом, рассматриваемый процесс находится между случайным и детерминированным хаосом (рис. 1). То есть сам восстановленный аттрактор ВСР не всегда может быть охарактеризован одним параметром фрактальной размерности (показателем скейлинга). Более совершенные методы, состоящие в переходе от одного показателя скейлинга к их спектру, находятся на стадии развития и не имеют достаточных фундаментальных результатов для их использования [11].

Поэтому предложено видоизменение метода фрактального анализа для исследования длинных (~50 000 отсчетов) записей ВСР. Рассмотрено поведение корреляционной размерности в различных частях исходной записи (независимо от ее сходимости или расходимости) [5]. Затем построили полученную кривую в зависимости от времени, сглаживая ее трендом. Установлен следующий ряд закономерностей в поведении тренда:

1) тренд сдвигается практически одинаково в различных частях исходной записи (в случае расходимости в конечномерном фазовом пространстве);

Ш«) = п

Отсутствие насыщения кривой

Детерминированный хаос

Рис. 1. Возможные виды рассматриваемых процессов

2) тренд незначительно меняет свою форму при изменении пространства вложения;

3) тренд реагирует на внутренние изменения ВСР (сон, физические нагрузки, стрессы, принятие алкоголя и лекарственных препаратов и др.).

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

на одно). Степень постоянства может также являться характеристикой рассматриваемого метода исследований (особенно в случае сходимости в конечномерных пространствах). Однако при таком подходе возникает ряд трудностей с возможностью динамического изменения параметра задержки, которые мы оставили открытыми [6]. Это связано, прежде всего, с отказом реструктуризации формы всего аттрактора в процессе его динамического обновления.

Разработанный программный комплекс

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

Рис. 2. Интерфейс и функциональные возможности комплекса на примере расчета показателей временного ряда, сгенерированного системой нелинейных дифференциальных уравнений

Histogram variability... Spectral analysis...

Autocorrelation function... Information function...

Play... Stop...

2D phase space... 3D phase space... Correlation dimension... Correlation entropy... Derivative of correlation entropy...

Lyapunov's exponent... Hurst exponent...

Playback speed: 1,00 § Window width: 2000 С

Window moving: 1 Delay time: 1

"1 Ï

Space dimension: 3

Hurst constant a: 1,00 [|]

Рис. 3. Функциональные возможности разработанного программного обеспечения

1. Спектральные компонеты ВСР. Мощность в различных частотных диапазонах (LF, HF, VLF, ULF), индекс централизации IC.

2. Гистограмму вариабельности с вычислением характеристик вариационной пульсометрии (Мо, АМо, MxDMn, MxRMn, триангулярный индекс).

3. Автокорреляционную функцию с указанием ее первого нуля.

4. Функцию средней взаимной информации с указанием ее первого локального минимума.

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

1. Графики восстановленного аттрактора в двухмерном фазовом пространстве.

2. Графики трехмерного восстановленного аттрактора в специально реализованном с помощью библиотеки компьютерной графики OpenGL трехмерном фазовом пространстве. Наблюдение за динамическим обновлением аттрактора возможно из любой точки пространства.

3. Корреляционную размерность и ее динамический тренд.

4. Корреляционную (аппроксимированную) энтропию и ее динамический тренд.

5. Старший показатель Ляпунова и его динамический тренд.

6. Экспоненту Херста и ее динамический тренд.

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

фикации меню Parameters. В последней версии меню Parameters содержит следующий список параметров (рис. 3):

скорость «проигрывания» массива (Playback speed);

количество точек ширины окна (Window width);

количество точек, на которое перемещается окно (Window mowing);

параметр задержки (Delay time);

размерность пространства вложения (Space dimension);

константа, необходимая для вычисления показателя Херста (Hurst constant).

Проверка правильности вычислительных процедур комплекса проведена на известных детерминированных моделях [6] (аттракторе Хенона, отображении Икеды, отображении Лози, системе уравнений Лоренца, системе уравнений Ресслера и др.). Сравнение результатов вычислительных процедур на статических выборках проведено со значениями, полученными в программе RQA (авт. C.J. Webber), пакете TISEAN (авт. R. Hegger и соавт.), пакете Фрактан 4.4 (авт. В. Сычев) и др. Установлены значительные преимущества вычисляемых комплексом энтропийных оценок за счет новых алгоритмов их уточнения, разработанных в ходе выполнения проекта.

Динамическая фрактальная обработка нормального кардиоритма

Для этих целей применены записи базы Normal Sinus Rhythm RR Interval Database (nsr2db) сервера [4]. Опишем обработку одной из записей базы. Обследована женщина 64 лет с примерной частотой сигнала 128 Hz. По времени длительность записи составляет 22 ч 33 мин и содержит 106 379 норма-норма интервалов, 68 вентрику-лярных и 13 атриовентрикулярных сокращений. Количество разрывов записи - 375.

Рис. 4. Двухмерный (слева) и трехмерный (справа) восстановленный аттрактор для нормального кардиоритма

и трехмерного аттрактора, а на рис. 5 - область двухмерного аттрактора, содержащая все нормальные сокращения.

Временная задержка находилась с помощью функции средней взаимной информации (рис. 6). Метод ее нахождения с помощью автокорреляционной функции не дает положительных результатов (в данном случае она равнялась 24 230).

Корреляционная размерность в течение суток, как и следовало ожидать [5], делится на два участка - с малыми (~0,5) и большими (~1,5) значениями тренда (рис. 7). Эти участки по времени соответствуют сну и бодрствованию организма. Таким образом, в данном случае рассмотренная нами динамическая характеристика тренда фрактальной размерности ВСР может являться характеристикой систем, регулирующих ритм сердца. Причем эта характеристика универсальна и не зависит от сходимости или расходимости процесса

Рис. 5. Область двухмерного восстановленного аттрактора, содержащая все интервалы без экстрасистол

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

Рис. 6. Временная задержка, найденная с помощью функции средней взаимной информации для ВСР нормального кардиоритма

Рис. 7. Корреляционная размерность нормального кардиоритма в течение суток

в конкретной точке. Ее практическое применение будем развивать в дальнейших исследованиях.

Разработанное в ходе выполнения проекта автоматизированное программное обеспечение позволяет расширить возможности фрактальной

диагностики ВСР. Практическая значимость полученных результатов может использоваться в клинической медицине. Для этих целей планируется создание комплекса холтеровского мони-торирования.

Работа выполнена при поддержке РФФИ (грант № 09-08-01135-а).

список литературы

1. Carvajal, R. Dimensional analysis of HRV in hypertrophic cardiomyopathy patients [Text] / R. Carvajal, J.J. Zebrowski, M. Vallverdi [et al.] // IEEE Engineering In Medicine And Biology. -2002. -№ 21 (4). -P. 71-78.

2. Grassberger, P. Measuring the strangeness of strange attractors [Text] / P. Grassberger, I. Procaccia // Physica D. -1983. -Vol. 9. -P. 189-208.

3. Malik, M. Heart rate variability. Standards of measurement, physiological interpretation and clinical use [Text] / M. Malik, J. T. Bigger, A. J. Camm [et al.] // European Heart J. -1996. -№17. -P. 354-381.

4. PhysioNet. The research resource for complex physiologic signals [Электронный ресурс] / Режим доступа: http:// www.physionet.org

5. Антонов, В.И. Динамический тренд корреляционной размерности как характеристический показатель жизнедеятельности организма [Текст] / В.И. Антонов, А.И. Загайнов, А.Н. Коваленко // Научно-технические ведомости СПбГПУ Сер. Информатика. Телекоммуникации. Управление. -2009. -№. 6 (91). -C. 111-119.

6. Антонов, В.И. Аппаратно-программный комплекс энтропийно-динамического мониторинга состояния кардиоритма [Текст] / В.И. Антонов, А.И. Загайнов, А.Н. Коваленко, Ву ван Куанг // Научно-технические ведомости СПбГПУ. Сер. Информатика. Телекоммуникации. Управление. -2011. -№ 1(115). -C. 143-150.

7. Ахметханов, Р.С. Применение теории фракта-

лов и вейвлет-анализа для выявления особенностей временных рядов при диагностике систем [Текст] / Р.С. Ахметханов // Вестник науч.-техн. развития. -2009.-№ 1 (17). -С. 26-31.

8. Баевский, Р.М. Анализ вариабельности сердечного ритма при использовании различных электрокардиографических систем. Ч. 1 [Текст] / Р.М. Баевский, Г.Г. Иванов, А.П. Гаврилушкин [и др.] // Вестник арит-мологии. -2002. -№ 24. -С. 65-86.

9. Гудков, Г.В. Диагностические возможности определения детерминированного хаоса в структуре вариабельности ритма сердца плода [Текст] / Г.В. Гудков // Вестник муниципального здравоохранения. -2008. -№ 1 (1). -19 с.

10. Колюцкий, А.К. Исследование вариабельности сердечного ритма при анализе аритмий [Текст] / А.К. Колюцкий, Г.Г. Иванов, В.Е. Дворников [и др.] // Вестник Рос. ун-та дружбы народов. Сер. Медицина. -2001. -№ 2. -С. 113-130.

11. Павлов, А.Н. Мультифрактальный анализ сигналов на основе вейвлет-преобразования [Текст] / А.Н. Павлов, В.С. Анищенко // Изв. Саратовского унта. Сер. Физика. -2007. -Т. 7. -Вып. 1. -С. 3-25.

12. Янсон, Н.Б. Моделирование динамических систем по экспериментальным данным [Текст] / Н.Б. Янсон, В.С. Анищенко // Изв. вузов «ПНД». -1995. -Т.3. -№ 3. -С. 112-121.

УДК 621.785

С.Е. Коршиков

МОДЕЛИРОВАНИЕ ПОЛЕЙ УПРУГИх ДЕФОРМАЦИИ В ПРОЦЕССЕ ИНДУКЦИОННОГО НАГРЕВА ВРАщАЮщИхСЯ ЗАГОТОВОК

Индукционные нагревательные установки (ИНУ) проходного и периодического действия широко применяются на практике для индукционного нагрева металлов перед последующей обработкой давлением, поскольку они обладают рядом технико-экономических преимуществ по сравнению с конкурентоспособными технологиями.

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