Математическая модель компонентов, обуславливающих вариативность сердечного ритма Текст научной статьи по специальности «Математика»

Романова Г.В., Чертенкова О.С. УДК 57.087

МАТЕМАТИЧЕСКАЯ МОДЕЛЬ КОМПОНЕНТОВ, ОБУСЛАВЛИВАЮЩИХ ВАРИАТИВНОСТЬ СЕРДЕЧНОГО РИТМА

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

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

Исходной концепцией для построения модели явились исследования по волнам кровяного

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

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

вого узла (1ИЯ) с её нестабильностью, наглядно проявляющейся при высоких уровнях рабочей нагрузки.

Выход модели, т.е. сердечный ритм (ИЯ), преобразуется в минутный объём крови, выбрасываемый сердцем с коэффициентом пропорциональности КИЯ. Передаточная функция сосудистой системы описывает колебательное поведение кровяного давления (Р). Воздействие дыхательных движений (РР) на кровяное давление учитывается при суммировании в компараторе 3, сосудистые воздействия (ЛБ) - в компараторе 4 (см. рис.1). Частота афферентных разрядов (БАР) с барроре-цепторов преобразуется в ЦНС в частоту эфферентных разрядов (УС) с коэффициентом пропорциональности КУ8. Нисходящее воздействие ЦНС на синусовый узел осуществляется через волокна вагуса, реализуя быструю регуляцию АД, представленную в модели посредством блоков синусового узла.

Существенным отличием этой модели от моделей [2, 3] заключается в замене дифференциальных уравнений и передаточных функций дискретной математической моделью, более соответствующей импульсному характеру изучаемого процесса и применении соответствующих современным данным исходных значений ряда физиологических переменных.

Модель содержит три нелинейных элемента (рис.2), четыре линейных, описанных передаточными функциями.

1) РО РР

2) РАР А Р5

3) Г«^ УС

Рис. 2. Нелинейные элементы модели: 1) элемент системы дыхания, пропускающий только отрицательную полуволну дыхательного цикла; 2) элемент барорецепторного контроля, минимальная уставка которого соответствует 40 мм рт. ст.; 3) элемент синиаурикулярного узла (АН)

Математическое описание нелинейных элементов дано в формулах:

¡0, если РР > 0 [КО х РР, если РР < 0 , где КО = 1 - коэффициент усиления;

¡0, если р3 < ркр (ркр = 40) [кгР - Ркр), если Р3 > Ркр ;

0, если УС < 0 УС

1) РО =

2) гаг =

3) ГУСГ =

если

ус > 0 •

_А + В х УС

Значения А = 1,74; В = 0,96 приведены в работе [4] .

Четыре передаточных звена включают в себя две пары однотипных звеньев вида:

Ку1,2 Х Т1,2 Х $ ,

К,,^) = -

К3,4 ($) =

1 + т12 X $

Ку3,4 Х

$ + 2£3,4 х $ х -3,4 + -3,4

(1)

(2)

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

Чтобы перейти к разностным уравнениям достаточно от передаточных функций, представляющих собой преобразование Лапласа непрерывной функции, перейти к преобразованию Лапласа дискретной функции, которое часто называют 2-преобразованием, а затем от него уже непосредственно перейти к разностному уравнению.

Пусть имеем непрерывную передаточную функцию К(8)

Чтобы от непрерывной функции времени перейти к дискретной функции на входе системы ставят ключ, который на мгновение замыкает входной сигнал Х, а затем, в течение времени Т держит его разомкнутым. Т называют периодом дискретности.

Ключ математически описывается суммой 8 -функций

ХТ (г) = Т £ 8 (г - кТ)

(3)

к = 0

Передаточная функция ключа есть преобразование Лапласа.

По определению преобразования Лапласа:

Кл (5) = Т - пТ = Т £

п=0 0 п=0

-БпТ

от

СИСТЕМНЫМ АНАЛИЗ И МЕЖДИСЦИПЛИНАРНЫМ ПОДХОД В ИССЛЕДОВАНИЯХ

так как

- а)в~

'ё' = е

-аЯ

из свойств 8 -

функции. Преобразование Лапласа для суммы 8 -

ВД гр

-ЯпТ Т

функций Т

1 - е~

. Это преобразование 1 - 1 - ^ 1 - г Р Р

уже называют ^преобразованием.

Найдем передаточную функцию непрерывного преобразования Лапласа совместно с ключом

Х(Б)

Ккл (Я )* К (Я)

непрерывная Y*(S)

приведенная

часть

Она равна свертке передаточной функции ключа ккл(Я) и непрерывной части К(Я), т.е.

1 С+

К*(Я) = — \кш (Я-Л)К(Л)ёЛ = 2П 1

С-

ГТ1 С+

— У-

2п 1

К (Л)

С-

Обозначим

- е

-(Я-Л)Т

ёЛ

ЯТ

г=е ,тогда

-+

К (г)=У т

К(Л)

ёЛ.

2П - 1 - г-1еЛТ

С -

Это формула перехода от непрерывного преобразования к дискретному ^преобразованию

[5].

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

Рассмотрим конкретно уравнения (1) и (2). Уравнение (1) (индексы опущены)

К (Я) = (тх Я )—Т— .

1 + тх Я

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

А(п) = [X(п) -X(п - 1)];п = 1,2,...,вд . Второй множитель преобразуем

К

К

1 + тх Я

т\ Я + -

Из таблиц ^преобразования следует

К * (г) = Т Х К Х г Т Х К г

(

Чтобы сделать запись более компактной обозначают еЯТ = г, тогда

Т Т . = -Т.г

тх г - е

V

т

г-е

Входом для звена ъ является А *(г), а выход обозначим У*(ъ). Тогда

Т х К г У *(г)

т -т А*(г)

г - е т

Отсюда

У *( г) хтх

г-е

= Т х К х г хА *( г)

или

У * (г) х г - У * (г) х е т=А *(г) х г х

Т х К

т

Переходя к дискретной функции времени и учитывая теорему о сдвиге получим выражение для разностного уравнения

Т

У (п +1) - е тх У (п) = — х К х А(п +1). т

Сдвинем аргумент на единицу, заменив п + 1 на п, а п на п - 1 и получим окончательно

У (п) = А х У (п -1) + А [X(п) - X(п -1)] (5)

А = е

А =

К х Т

X(п) - X(п -1) = А(п). Найдем ъ-преобразование уравнений (2). Пусть Xвх (') и Xвь1х (') - соответственно функции времени на входе и выходе динамического звена (2). Тогда

Xь (Я)

К х ^

(6)

X,,:(Я) Я2 + 2£х Ях w + w2 В таблицах ъ-преобразований нет функции, точно соответствующей правой части (6). Поэтому преобразуем знаменатель правой части, так чтобы получить полный квадрат. Для этого прибавим и отнимем w2. Тогда

Я2 + 2£Яw + w2 + ^2w2 2 =

= (Я + ^)2 + (^2-#2w2) =

= (Я + ^)2 + w2(1 -#2) = (Я + а)2 +р2,

п=0

Г

7

Г

I

т

где a = Sw и в = w^ 1 -£2 . Окончательно получим

K ( S ) = -

Kw2

(S + a)2 +в2

Этой функции соответствует z-

преобразование

K х w2 x T z х e ~aT x sin BT K * ( z ) =----

v ' n 2 o -aT /JT-1 , -2aT

p z - 2z x e x cos pT + e Обозначим

A8 = 2 x e~aT cos pT = 2 x e^wT cos[ w^ 1 -£2 T ] ;

A, = e-2w ;

-SwT sm[wV 1 -S2T] ,

KTw2 Aio =-î— e

тогда

K *( z) =-

Так как K *( z) = _s«

z A.g x z + A,

Xeux (z)

Xх (z)

то

Xeux (z) _ A10 x z

Xex (z) z - A8 x z + A,

или

(z)[Z2 -A xZ + Al = Xex(z)X Ao XZ ,

Z 2 X Xeux ( Z) - A x Z x ХеЫх ( Z) + A, x Xeux ( Z) = = X ex ( Z ) X Ao X Z

С учетом теоремы о сдвиге получаем разностное уравнение

Xeux (П + 2) - AXeux (n + 1) + AXeux (n) = = AoXex (П + 1) Заменим m = n + 2, тогда Xeux (m) =

AX.« (m -1) - AXeu^x (m - 2) + A10-Xex (m -1) Это окончательное выражение для разностного уравнения, в котором аргумент m можно снова заменить на n.

При Ç = 1 уравнение (2) переходит в

K (S) =

Kw¿

(S + w)

2 '

Ему соответствует K *( z ) =

KTw x z x e ( z - e - wT )2

-wT

Xeux (П) =

= AXeux (П - 1) - A4XeuX (n - 2) + A5 Xex (n - 1)

4 eux -2 wT

; A5 = Kw Te

wT

где A3 = 2e w ; A4 = e

Расчет переходных процессов при T = 0,1 сек и Т = 1/6 сек показал, что надо выбрать Т = 1/6 сек. При таком Т наблюдается устойчивое решение системы при всех условиях и дальнейшее уменьшение Т приводит только к увеличению времени решения задачи.

Можно сделать следующие выводы:

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

2. Переход к разностным уравнениям обеспечивается переходом от преобразования Лапласа непрерывной функции к преобразованию Лапласа дискретной функции.

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

БИБЛИОГРАФИЯ

1. Koepchen, H.P. History of studies and concepts of blood pressure waves / H.P. Koepchen // Mech. Blood Pressure Waves. Berlin et. al. 1984. P. 323.

2. Luckzak, H. Decomposition of heart rate variability under the ergonomics aspects stressor analysis / H. Luczak, U. Philipp, W. Rohmetr // The study of heart rate variability. Oxford, 1980.

3. Miyawaki K. Analysis and simulation of the periodic heart rate fluctuation / K. Miyawaki, T. Ta-kahashi, H. Takamura // Received Nov. 25. 1965. № 709. P. 315.

4. Luczak H. Fractioned heart rate variability. Part I: Analysis in a model of the cardiovascular and car-diorespiratory system / H. Luczak // Ergonomics, 1978. V. 21. № 11. P. 895-911.

5. Anderson D.E. Ambylatory monitoring of respiration: Inhibitory breathing in the natural ehviro-ment / D.E. Anderson, K. Code, J.A.Hagthwite // Psychophysiol, 1992. V. 29. № 5.

и разностное уравнение