Моделирование мышечного сокращения

УДК 531:577.353 Е. А. Селицкая

МОДЕЛИРОВАНИЕ МЫШЕЧНОГО СОКРАЩЕНИЯ

Примером высокоспециализированного движения является мышечное сокращение. Для изучения молекулярных механизмов мышечного сокращения предложены математические модели А. Хилла [1], И. В. Дещеревского [2]. Критерием их адекватности является степень совпадения описания макроскопических свойств мышцы в моделях с результатами эксперимента, проведенного над препарированной мышцей, возбуждаемой внешним электрическим током. В работах [3, 4] рассмотрена механическая модель сокращения поперечно-полосатой мышцы при изолированном дельтаобразном возбуждении. На основе модели из работ [5, 6] выполнено моделирование нервного импульса и проведены пробные вычисления мышечного сокращения [7, 8]. В настоящей статье приводится описание модели и численное исследование отдельных этапов сокращения мышцы, дополняющее результаты [7, 8]. В этой модели учитывается связь информационных импульсов, поступающих из центральной нервной системы, с микроскопическими процессами в функциональных элементах мышцы (саркомерах, рис. 1).

Рис. 1. Строение саркомера и схема расположения тонких (актиновых - 1) и толстых (миозиновых - 2) нитей на разных уровнях поперечного сечения саркомера

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

Селицкая Екатерина Александровна — ассистент кафедры механики управляемого движения факультета прикладной математики—процессов управления Санкт-Петербургского государственного университета. Научный руководитель: доктор физико-математических наук, проф. В. С. Новоселов. Количество опубликованных работ: 4. Научные направления: теоретическая механика, математическое моделирование и вычисления в среде Matlab нервной проводимости и сокращения мышцы. E-mail: velikova_e@mail.ru.

© Е. А. Селицкая, 2009

Математическая модель. В модели приняты следующие обозначения: длина сар-комера l(t), где t - время; N(l) - общее число миозиновых головок мышцы, находящихся в перекрывающихся областях между актином и миозином; xi(t) - число миозиновых головок, соединенных с актином; X2(t) - число головок, не соединенных с актином; но присоединивших аденозинтрифосфат (АТФ); x^(t) - число свободных головок; X4 (t) - число активизированных центров, но не связанных с миозиновыми головками; X5(t) - число свободных центров актина; x^(t) - концентрации ионов Ca2+; U(t) - потенциал возбуждения; xpP - равновесная концентрация свободных ионов Ca2+; xr(t) - количество молекул АТФ; xg(t) - количество молекул аденозиндифос-фата (АДФ). Предложенная в [5] математическая модель имеет следующий вид:

ax2x4 — pxi'xi + N-1Nxi,

7x3x7 — ax2x4 + N -1Nx2, pxi'xi — yx3 xi + N-1Nxs, kx5 x6 — ax2x4 + N-1Nx4,

pxi'xi — kx5 x6 + N 1Nx5, vU + ¡j,(xPp — xe),

cxg — ex 1x7 — 7x3 xi — N-1N (x1 + x2),

2^1x7 — cxg.

В этих уравнениях коэффициенты а, в, 7, к, v, ¡л - положительные постоянные, настраиваемые при различных режимах жизнедеятельности. Величина c определяется интенсивностью восстановления АТФ, зависит от работы легких и является коэффициентом настройки. Добавление члена N-1Nxj в правой части уравнений учитывает скорость изменения x*, вызванную изменением N.

Из системы (1) следуют соотношения

x1(t) + x2 (t) + x3(t) = N,

x1

x1(t) + x2(t) + x7(t) + x8(t) = s = const.

Во втором уравнении (2) для упрощения общее число способных к активизации центров полагается также равным N. Последнее уравнение отвечает сохранению числа s молекул аденина или аденозина, входящего в состав АТФ и АДФ.

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

dx1

dt

dx2

dt

dx3

dt

dx4

dt

dx5

dt

dx6

dt

dx7

dt

dx%

dt

Напряжение, уел. ед. 60 г

50

40

30

20

10

0

20

40

60

80 ЮО Время, уел. ед.

Рис. 2. Возбуждение мышцы одиночным импульсом

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

Преобразование уравнений. Введение безразмерных переменных у* =

где г = 1,8, г/е = г = ТУ#, V = Щр-, позволяет привести систему (1), (2) к следующему Х6 Х6

виду:

= ау2у4 - вУ1У7

3у2

~7~ = УУЗУ7 - ау2УА, ат

3у4 хр6

----- = УС------УкУа — СУ.УоУл ,

dy7 с 1 dN

~dr = NV8 ~ УГ ~ 7УЗУГ ~ N~dr + У2 + dy^ с 1 dN

= 2,%l!n - Nm - N^m На изометрическом этапе сокращения мышцы (N = const) общее число миозино-вых головок связано с потенциалом действия V и постоянными л, с, xpP таким образом, что величины n> W> lv в уравнениях (3) можно рассматривать как «управляющие» элементы системы. Поэтому при изменении значения N и одновременном пропорциональном изменении величин V, л, с, xpP правая часть уравнений остается неизменной. Преобразованная система (3) удобна для изучения процессов включения и выключения возбуждения на последовательных участках времени (рис. 3). При распространении потенциала возбуждения на мембране мышечного волокна в виде синусоиды система выходит из равновесного положения. Вследствие этого процесса усилие в мышце возрастает. При снятии возбуждения V (через 0.8 условных единиц времени) усилие в мышце начинает убывать пропорционально концентрации ионов Ca2+ до равновесного положения.

Напряжение, уел. ед.

Время, уел. ед.

Рис. 3. Возбуждение мышцы из положения равновесия Сплошная линия — ух, пунктирная — уб.

Вычисления, проведенные в среде Matlab, позволяют проследить изменение не только величины усилия yi, но и других переменных системы (3). Ввиду краткости статьи, приведем лишь их начальные значения, которым соответствуют кривые усилия и концентрации ионов Ca2+ на рис. 3:

Переменные системы ........ У2 Уз У4 У5 У7 У8

Значение .................. 0.81949 0.015042 0.00027381 1.8343 0.02716 8.9879

Напряжение, уел. ед.

Время, уел. ед.

Рис. 4- Возбуждение мышц с различными функциональными перекрытиями

1 — у1, N = 1000; 2 — уб; 3 — у\, N = 100.

Напряжение, уел. ед. 0.3 г

0 12 3

Время, уел. ед.

Рис. 5. Сохранение баланса свободных молекул АТФ и АДФ Сплошная линия — у7, пунктирная — ув.

На рис. 4 представлена зависимость усилия для двух вариантов мышечной ткани с различным числом миозиновых головок в саркомере и возбуждения при определенной частоте импульсов. Поведение кривых на рисунке показывает, что при одной и той же частоте потенциала возбуждения, распространяющегося на мембране мышечного волокна, мышцы с различным функциональным числом миозиновых головок реагируют по-разному. Мышца с N = 100 медленнее развивает силу, но при этом удерживает ее. Мышца с N = 1000 быстрее развивает усилие, но не может его удержать. При увеличении частоты импульсов кривая 1 выходит на верхнюю огибающую рис. 4.

На рис. 5 представлена зависимость концентраций свободных молекул АТФ и АДФ при такой же частоте потенциала, что и на рис. 3.

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

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

Литература

1. Хилл А. Механика мышечного сокращения. Старые и новые опыты / пер. с англ. Ю. А. Шаронова; под ред. и с предисл. Г. М. Франка. М.: Мир, 1972. 184 с.

2. Дещеревский В. И. Математические модели мышечного сокращения. М.: Наука, 1977. 160 с.

3. Филиппов Б. В. Усилия и энергетические затраты в поперечно-полосатых мышцах. I // Вестн. С.-Петерб. ун-та. Сер. 1: Математика, механика, астрономия. 2000. Вып. 1. С. 86—91.

4. Филиппов Б. В. Усилия и энергетические затраты в поперечно-полосатых мышцах. II // Вестн. С.-Петерб. ун-та. Сер. 1: Математика, механика, астрономия. 2001. Вып. 1. С. 112—122.

5. Новоселов В. С. Статистические модели нейродинамики. СПб.: С.-Петерб. гос. ун-т, 2004. 64 с.

6. Новоселов В. С., Королев В. С. Модель возбуждения мышцы // Труды IV Междунар. конференции «Идентификация систем и проблемы управления». М., 2005. C. 367—374.

7. Селицкая Е. А. Синаптическая передача нервного возбуждения // Труды XXXVIII Между-нар. науч. конференции «Процессы управления и устойчивость». Санкт-Петербург, 2007 г. / под ред. А. В. Платонова, Н. В. Смирнова. СПб.: Изд-во С.-Петерб. ун-та, 2007. С. 295—299.

8. Селицкая Е. А. Численное моделирование сокращения поперечно-полосатой мышцы // Труды XXXIX Междунар. науч. конференции «Процессы управления и устойчивость». Санкт-Петербург, 2008 г. / под ред. Н. В. Смирнова, Г. Ш. Тамасяна. СПб.: Изд-во С.-Петерб. ун-та, 2008. С. 277-282.

9. Бэгшоу К. Мышечное сокращение / пер. с англ. Н. А. Габеловой. М.: Мир, 1985. 159 с.

Статья рекомендована к печати проф. В. С. Новоселовым.

Статья принята к печати 28 мая 2009 г.