Актуальной задачей на сегодняшний день является изучение особенностей проведения электрического возбуждения по сердечной мышце. Нами было проведено сравнение скорости распространения волны возбуждения по трехмерному волокну модели левого желудочка сердца и по одномерному волокну такой же длины. Разработано программное обеспечение, позволяющее проводить расчеты для волокон различной длины и различного пространственного расположения в толще стенки модели левого желудочка. Полученные результаты могут быть полезны в исследованиях физиологии сердечной мышцы.
позволяет использовать различные электрофизиологические модели кардиомиоцитов для расчета охвата левого желудочка волной возбуждения. Для моделирования электрических явлений, происходящих в сердце, мы использовали электрофизиологическую модель кардиомиоцита Алиева – Панфилова, качественно описывающую потенциал действия кардиомиоцита [3].
Миокард как проводящая среда имеет анизотропные свойства, которые заключаются в разнице в скоростях распространения волны возбуждения вдоль и поперек волокон. Кроме того, имеет место вращательная анизотропия миокарда, вызванная трансмуральным вращением волокон в стенке левого желудочка, что также оказывает влияние на скорость проведения электрического возбуждения [4] (рис. 1). Рис. 1: Осесимметричная модель левого желудочка (слева) и трансмуральное вращение волокон в толще его стенки (справа) [5] 3
Описание программного комплекса Вычислительный эксперимент включал следующие этапы: 1. Выбор зоны начальной стимуляции (область миокарда, откуда происходит запуск волны возбуждения).
2. Проведение расчета на модели левого желудочка сердца со встроенной электрофизиологической моделью Алиева – Панфилова для выбранной зоны стимуляции.
3. Выбор волокна на спиральной поверхности модели левого желудочка и нахождение времени прихода волны возбуждения в каждый из узлов волокна.
4. Решение одномерной задачи для волокна той же длины, что и выбранное волокно на спиральной поверхности, с аналогичными параметрами модели Алиева – Панфилова.
Результаты проведенного вычислительного эксперимента представлялись в виде графика скоростей прохождения волны возбуждения по одномерному и трехмерному волокнам. Кроме того, для каждого эксперимента строилась спиральная поверхность с наложенным на нее волокном, для оценки его расположения в толще стенки левого желудочка.
Программное обеспечение (табл. 1) разработано на языке python 2.7, с использованием библиотек matplotlib (построение двумерных графиков) [6], mayavi (построение и работа с трехмерной графикой) [7] и pyside (создание графического интерфейса) [8]. Для проведения одного вычислительного эксперимента требуется в среднем 1 – 3 минуты в зависимости от производительности процессора и установленного модельного времени эксперимента (при фиксированном числе узлов разбиения волокна и шаге по времени). Модуль one_dim_solver:py three_dim_processing:py interpolation:py
gui:py main:py
Назначение Решатель системы дифференциальных уравнений в частных производных модели Алиева–Панфилова для случая одномерных волокон Обработка результатов расчета на модели левого желудочка. Построение волокна на спиральной поверхности и подготовка данных для трехмерной визуализации Реализация метода трилинейной интерполяции Построение графического интерфейса пользователя Модуль запуска программы, осуществляющий вызов остальных модулей и обеспечивающий обмен данными между ними 3.1 Решение задачи распространения волны возбуждения по одномерному волокну Модель Алиева – Панфилова представляет собой систему из двух дифференциальных уравнений в частных производных: (u; v) = 0 + = (u; v)( v ku(u a 1));
1v где фазовая переменная u описывает трансмембранный потенциал (разность потенциалов между внешней и внутренней сторонами клеточной мембраны), фазовая переменная v – мембранный ток восстановления; D1 – коэффициент диффузии вдоль направления волокна; k, a, 0, 1, 2 – параметры модели клетки.
Мышечное волокно миокарда может быть представлено в виде последовательности узлов x0; x1; : : : ; xn. Будем считать, что в начальный момент времени стимулируется только первый из n узлов волокна: Граничные условия второго рода моделируют электрическую изолированность концов волокна: Решая систему (1) – (3), мы получим картину распространения волны возбуждения по одномерному волокну с постоянной скоростью (из-за постоянства коэффициента диффузии). Первое уравнение представляет собой уравнение теплопроводности с функцией внутренних источников. Для численного решения уравнения теплопроводности была выбрана безусловно устойчивая неявная разностная схема. Для численного решения второго уравнения системы (ОДУ) мы воспользовались явным методом Эйлера.
По ходу решения мы также фиксировали момент времени, когда потенциал превышал заранее установленный порог (находили время прихода волны возбуждения). В поставленных вычислительных экспериментах этот порог был равен 0.5 (в модели Алиева – Панфилова безразмерный потенциал u 2 [0; 1]). 3.2
Обработка результатов расчета на трехмерной модели левого желудочка Построение модели левого желудочка осуществляется в специальной системе координат ( , , ), где – положение точки в толще стенки левого желудочка ( = 0 – эндокард, = 1 – эпикард), – аналог географической широты ( = 0 – самая верхняя плоская часть модели левого желудочка, =2 – верхушка, самая узкая часть модели левого желудочка), – аналог географической долготы (от 0 до 2 вокруг оси вращения левого желудочка). Точки в специальной системе координат достаточно просто переводятся в цилиндрическую и декартову системы координат [2]. Модель левого желудочка имеет следующие геометрические параметры (на рис. 2 представлено меридиональное сечение модели): высота левого желудочка Z, толщина стенки на основании d0, толщина стенки на экваторе de, ¾широта¿ экватора на эпикарде l0, ¾широта¿ экватора на эндокарде l1, толщина стенки на верхушке h, радиус желудочка на эпикарде основания r0, радиус желудочка на эпикарде экватора re.
Рис. 2: Меридиональное сечение модели левого желудочка [2] Система дифференциальных уравнений в частных производных модели Алиева–Панфилова для трехмерного случая имеет вид: (u; v) = 0 +
1v = div(D grad u) ku(u a)(u 1) uv; где D – матрица диффузии (с элементами Di;j = D2 i;j + (D1 D2)wiwj, где D1 – коэффициент диффузии вдоль волокон, D2 – коэффициент диффузии поперек волокон, i;j – символ Кронекера, w~ = (w0; w1; w2) – единичный вектор направления волокна).
Для расчета на трехмерной модели левого желудочка используется ранее разработанный программный комплекс [9], обеспечивающий построение геометрической модели и расчет волны возбуждения на трехмерной сетке. В результате данного расчета мы получаем данные о времени прихода волны возбуждения в каждый из узлов сетки левого желудочка в форме: координаты узла ( i; j; k) и время прихода волны в узел Tijk. Задача заключается в построении волокна на спиральной поверхности и нахождении времени прихода волны в точки волокна, исходя из полученных данных. Согласно используемой модели левого желудочка [2], прообразом волокна на спиральной поверхности является хорда Y = const полукруга радиуса P = 1 (рис. 3).
Двигаясь равномерно по полярному углу от одного из концов хорды к другому, можно разбить хорду на конечное число узлов (отрезки разбиения равны не будут): = 1 n
0 ; Рис. 3: Горизонтальные хорды на полукруге. 0 и 1 – полярные углы соответственно правого и левого концов хорды [2] i = 0 + i
; i = 0; 1; 2; : : : ; n;
Y Pi = sin( i) ; i = 0; 1; 2; : : : ; n;
( ) = ; (P ) = (1 P ) : 2 где – шаг разбиения по полярному углу, n – число узлов разбиения. Узлы хорды однозначно отображаются в узлы волокна на спиральной поверхности:
По массиву координат узлов волокна мы вычисляем длину волокна, а также время прихода волны возбуждения в эти узлы, воспользовавшись известными значениями времени прихода волны Tijk в узлы сетки модели левого желудочка (методом трилинейной интерполяции). 3.3 Графический интерфейс пользователя
На рис. 4 представлен графический интерфейс программного комплекса. В левой половине – панель управления, состоящая из блоков: ¾Fiber parameters¿ – задание числа узлов разбиения волокон и выбор волокна на спиральной поверхности; ¾Conditions¿ – загрузка файла, хранящего значения времен прихода волны возбуждения в узлы модели левого желудочка (Time of arrival), и конфигурационного файла, хранящего геометрические параметры модели левого желудочка, параметры модели Алиева – Панфилова и другие параметры вычислительного эксперимента (Initial parameters); ¾Control bar¿ – запуск или остановка вычислительного эксперимента; ¾Status screen¿ – окно, отображающее ход процесса вычисления. В правой половине – окно с графиком сравнения (для вкладки ¾Comparison graph¿).
Вкладка ¾Spiral surface¿ переключает на трехмерное изображение спиральной поверхности с волокном (рис. 5). 4
Пример расчета Для вычислительного эксперимента были выбраны следующие параметры геометрической модели левого желудочка: высота левого желудочка Z = 100 мм, толщина стенки на основании d0 = 10 мм, толщина стенки на экваторе de = 10 мм, ¾широта¿ экватора на эпикарде l0 = 0.2, ¾широта¿ экватора на эндокарде l1 = 0.2, толщина стенки на верхушке h = 8 мм, радиус желудочка на эпикарде основания r0 = 30 мм, радиус желудочка на эпикарде экватора re = 33 мм. Границы расчетной области: 2 [0:1; 0:9]; 2 [0; =2]; 2 [0; 2 ]: Параметры электрофизиологической модели Алиева – Панфилова были равны: k = 8, a = 0:03, 1 = 0:12, 2 = 0:3, = 0:01. Коэффициент диффузии вдоль волокон D1 = 0:05, а поперек волокон – в девять раз меньше: D2 = D1=9. Такое соотношение связано с различием скорости распространения волны возбуждения вдоль и поперек волокон в миокарде по экспериментальным данным (поперек волокон в три раза меньше, чем вдоль волокон) [10]. Рис. 4: Графический интерфейс программы с графиком времени прихода волны как функции расстояния по волокну от первого стимулированного узла 5
Заключение Разработанная методика проведения вычислительного эксперимента и соответствующее программное обеспечение позволяют отслеживать скорость распространения волны возбуждения по выбранному волокну модели левого желудочка сердца, а также сравнивать полученные результаты с результатами на изолированном волокне такой же длины, что довольно сложно осуществить в условиях реального эксперимента. В частности, нами было показано наличие существенного опережения скорости распространения волны возбуждения по трехмерному волокну по сравнению с одномерным, что можно объяснить наличием окружающей трехмерное волокно тканью миокарда, оказывающей значительное влияние на скорость распространения возбуждения по волокну. Полученные результаты могут представлять интерес в вопросах, связанных с изучением электрических явлений, происходящих в миокарде во время сердечных сокращений. Благодарности Список литературы
Авторы выражают благодарность своим коллегам из Института иммунологии и физиологии УрО РАН: Л.Б. Кацнельсону, Н.А. Викуловой, А.Д. Хохловой, О.Э. Соловьевой. Работа поддержана Российским научным фондом №14-35-00005. [1] N.A. Trayanova. Whole Heart Modeling: Applications to Cardiac Electrophysiology and Electromechanics.
Circulation research, 108(1):113-128, 2011. [2] S.F. Pravdin. Nonaxisymmetric mathematical model of the cardiac left ventricle anatomy. Russian journal of biomechanics, 17, (62):84–105, 2013. [3] R.R. Aliev and A.V. Panfilov. A simple two-variable model of cardiac excitation. Chaos, Solitons and
Fractals, 7(3):293–301, 1996. [4] S.F. Pravdin, H. Dierckx, L.B. Katsnelson, O. Solovyova, V.S. Markhasin, A.V. Panfilov. Electrical wave propagation in an anisotropic model of the left ventricle based on analytical description of cardiac architecture. PLOS ONE, 9(5), 2014. [5] S.F. Pravdin, V.I. Berdyshev, A.V. Panfilov, L.B. Katsnelson, O. Solovyova, V.S. Markhasin. Mathematical model of the anatomy and fibre orientation field of the left ventricle of the heart. BioMedical Engineering OnLine, 54(12):1–21, 2013. [6] Matplotlib. http://matplotlib.org/ [7] The MayaVi Data Visualizer. http://mayavi.sourceforge.net/ [8] PySide. https://wiki.qt.io/PySide [9] S.F. Pravdin. A method of solving reaction-diffusion problem on a non-symmetrical model of the cardiac left ventricle, Proceedings of the 47th International Youth School-conference ¾Modern Problems in Mathematics and its Applications¿: 284–296, Yekaterinburg, Russia, 2016. [10] Kanai A., Salama G. Optical mapping reveals that repolarization spreads anisotropically and is guided by fiber orientation in guinea pig hearts. Circ Res, 77(4):784-802, 1995.
Software complex for the comparison properties of myocardial fibers in 1D and 3D of electrophysiological models Timur V. Nezlobinsky2, Sergei F. Pravdin1;2 1 – Krasovskii Institute of Mathematics and Mechanics (Yekaterinburg, Russia) 2 – Ural Federal University (Yekaterinburg, Russia) Keywords: left ventricle, cardiac model, myocardial fibers, rotational anisotropy.
Nowadays, investigations of the heart fibers and their role in myocardium conduction in relation to the electrical stimulation are actual in physiology. We compared the propagation speed of the excitation wave in threedimensional fibers in mathematical model of the left ventricle with the propagation speed in one-dimensional fiber. We developed software that allows us to conduct experiments for fibers with different lengths and spatial locations. Results can be useful in investigations of the cardiac physiology.