УДК 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(t)+x4(t)+x5(t) = N, (2)
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) и слитный тетанус при высокой частоте импульсов (3)
Для оценки адекватности модели уравнений (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 с. (A. V. Hill. First and last experiments in muscle mechanics.)
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 с. (C. R. Bagshaw. Muscle contraction.)
Статья рекомендована к печати проф. В. С. Новоселовым.
Статья принята к печати 28 мая 2009 г.