УДК 517.15
МОДЕЛИРОВАНИЕ ХИМИЧЕСКИХ ПРЕВРАЩЕНИЙ В АКТИВНЫХ ЦЕНТРАХ ХОЛИНЭСТЕРАЗ МЕТОДАМИ КВАНТОВОЙ ТЕОРИИ
А.В. Немухин, А.М. Кулакова, С.В. Лущекина, А.Ю. Ермилов, С.Д. Варфоломеев
(кафедра физической химии; e-mail: [email protected])
Обсуждена роль метода квантовой механики - молекулярной механики (КМ/ММ) для моделирования химических преобразований в активных центрах холинэстераз. Различными вариантами метода КМ/ММ исследованы молекулярные механизмы реакции реактивации фосфорилированной по остатку каталитического серина бу-тирилхолинэстеразы.
Ключевые слова: холинэстеразы, фосфорилированный фермент БХЭ, молекулярное моделирование, метод КМ/ММ.
Важность изучения химических преобразований в активных центрах ферментов центральной нервной системы ацетилхолинэстеразы (АХЭ) и бутирилхолинэстеразы (БХЭ) трудно переоценить ввиду значения этих ферментов для здоровья человека. Ацетилхолинэстераза отвечает за быстрый гидролиз нейромедиатора ацетилхо-лина, приводящий к обрыву передачи нервного импульса; ингибирование АХЭ используется в качестве терапии ряда тяжелых заболеваний. Бутирилхолинэстераза содержится в больших количествах в плазме крови и практически во всех тканях и способна брать на себя функции АХЭ при недостаточной экспрессии или инактивации последней. Оба фермента представляют мишени для фосфорорганических отравляющих веществ; по этой причине активно ведутся работы над созданием каталитических ловушек на основе БХЭ.
Каталитический механизм гидролиза субстратов АХЭ и БХЭ, относящихся к сериновым гидролазам, принципиально установлен. После формирования фермент-субстратного комплекса инициируется стадия ацилирования - образование первого тетраэдрического интермедиата с участием боковой цепи каталитического серина и субстрата и последующий разрыв химической связи, приводящий к ацил-ферменту. На стадии де-ацилирования повторяется сценарий образования тетраэдрического интермедиата при нуклеофиль-ной атаке каталитической молекулы воды, высвобождаются продукты реакции гидролиза субстрата и восстанавливается боковая цепь серина.
Визуализация подобных молекулярных превращений с атомарным разрешением, достижи-
мая при моделировании реакций ферментативного катализа методами квантовой механики (молекулярной механики (КМ/ММ)) [1], стала реальным подспорьем при экспериментальном исследовании ферментов [2]. Расчеты структур интермедиатов и переходных состояний, а также энергетического профиля для полного цикла химических преобразований в активном центре АХЭ при гидролизе ацетилхолина методом КМ/ММ впервые были описаны в работе [3]. Позднее [4] эти расчеты были повторены с использованием другой версии теории КМ/ММ, отличающейся существенными деталями квантово-механиче-ской методики, и с другой компьютерной программой. Согласие полученных результатов [3, 4] подтверждает надежность моделирования.
В настоящей работе мы применили данный подход для моделирования важнейшей реакции реактивации БХЭ с фосфорилированным остатком серина. Одна из основных целей при исследовании реактивации БХЭ - подбор такого варианта фермента за счет точечных мутаций, который способен эффективно гидролизовать ковалентно связанный аддукт:
В частности, известно [5, 6], что мутация С1у117И18 позволяет облегчить (по сравнению с процессами с натуральным ферментом) разрыв фосфорно-кислородной связи (-8ег-)0-Р за счет атаки нуклеофильной молекулы воды.
Накопленный опыт моделирования превращений в активных центрах ферментов по теории КМ/ММ показывает, что целесообразно включать в квантовую подсистему значительное (больше 100) число атомов и использовать кван-тово-химические методы теории функционала плотности (DFT) для вычислений энергий и сил, необходимых для построения энергетического профиля поверхности потенциальной энергии и анализа стационарных точек (интермедиатов, переходных состояний) на пути реакции. К сожалению, подобное моделирование требует больших затрат вычислительных ресурсов суперкомпьютерных центров и длительного времени. Поэтому поиск надежных, но менее затратных вариантов расчетов с использованием квантово-хи-мических подходов в рамках теории КМ/ММ представляет важную задачу. В настоящей работе мы рассматриваем полуэмпирический вариант DFT, называемый методом функционала плотности в приближении сильной связи DF-TB [7], для расчетов реакции реактивации БХЭ.
Модельные молекулярные системы были построены по мотивам кристаллографических структур комплексов БХЭ, образующихся при реакции с экотиофатом: PDB ID 1XLW для на-тивного фермента [8] и PDB ID 2XMD для мутанта Gly117His [9]. После добавления атомов водорода в предположении общепринятого состояния протонирования полярных аминокислотных остатков, приводящего к отрицательно заряженным боковым цепям Asp, Glu и положительно заряженным Lys, Arg, системы были сольватированы молекулами воды. Проводились предварительные расчеты методом молекулярной динамики в течение 5 нс с использованием следующих параметров: силовое поле Charmm27, ансамбль NPT, баростат Нозе-Гуве-ра (р = 1 атм.), термостат Ланжевена (Т = 298 К). Для расчетов методом КМ/ММ были взяты все атомы белка и ингибитора, а также молекулы воды, находящиеся на расстоянии не более 4 Â от атомов белка. Полученные модели включали почти 15 тыс. атомов и более 120 атомов в квантовой подсистеме. При разделении системы на квантовую и молекулярно-механическую части учитывалось возможное участие молекулярных групп в реакции или в формировании системы водородных связей с компонентами активного центра. Выбранная квантовая часть включала аминокислотные остатки каталитической триады (Ser198, His438, Glu325), связанный с Ser198 фрагмент диэтилфосфата, каталитическую молекулу воды (Wat), остаток глутаминовой кис-
лоты Glu197, принимающий протон в ходе реакции, остатки оксианионного центра (Glyl 16, Gly117 (или His117), Ala199), фрагменты Ser224, Asn322, Pro285 и две молекулы воды, входящие в систему водородных связей активного центра.
При расчетах методом КМ/ММ сопоставлялись результаты, полученные в двух вариантах вычислений энергий и сил в квантовой подсистеме: 1) DFT с использованием функционала BB1K [10] и базиса 6-31G**; 2) DF-TB. Для описания молекулярно-механической подсистемы в обоих случаях применялось силовое поле AMBER99. В первом варианте расчеты проводились по программе NWChem [11], во втором - по программе AMBER12 [12].
Структура реагента (RG) для фермента, содержащего гистидин в позиции 117, представлена на рис. 1 (все группы, изображенные шарами и стержнями, составляют часть квантовой подсистемы). Аминокислотные остатки Ser198, His438, Glu325 относятся к каталитической триаде БХЭ, а Gly116, His117, Ala199 - к оксиани-онному центру; Glu197 активирует молекулу воды Wat. На стадии образования первого ин-термедиата INT через первое переходное состояние TS1 протон от Wat переходит на Glu197, а гидроксил присоединяется к атому фосфора.
В качестве координаты реакции для стадии RG ^ TS1 ^ INT выбиралось расстояние между атомом кислорода воды Wat и атомом фосфора; при его последовательном уменьшении, сопровождающемся оптимизацией всех остальных геометрических параметров системы по методу КМ/ММ, была локализована конфигурация переходного состояния TS1. Спуск по потенциальной поверхности подтвердил, что найдена седловая точка, разделяющая конфигурации минимальной энергии RG и INT. Расчеты выполнены двумя вариантами с использованием либо приближения DFT(BB1K)/6-31G**, либо метода DF-TB в квантовой части. Обе модельные системы рассматривались с нативным ферментом (wt-BChE) и с мутантом Gly117His.
На рис. 2 представлена схема наиболее важных изменений, характерных для всех вариантов расчета вдоль пути RG ^ TS1 ^ INT. Вычисленные активационные барьеры (в ккал/ моль) указаны в таблице на рис. 2.
Анализ результатов показывает, что геометрические конфигурации реагента, интермедиа-та и переходного состояния, вычисленные с использованием DFT и DF-TB, достаточно близки. Оба метода решения квантовых уравнений (DFT и DF-TB) в рамках подхода КМ/ММ пра-
Рис. 1. Структура реагента (Яв) для мутанта 01у117Н1з фосфорилированной бути-рилхолинэстеразы по результатам оптимизации геометрических параметров методом
КМ/ММ
Рис. 2. Преобразования на первой стадии реакции RG ^ TS1 ^ INT (значения высоты энергетических барьеров в таблице приведены в ккал/моль)
вильно описывают роль мутации Gly117His; в обоих случаях активационные барьеры заметно меньше для мутированного фермента, чем для нативного. Даже отличия в самих величинах барьеров вполне терпимые: 3 и 7 ккал/моль в случае мутанта и нативного фермента соответственно.
Стадия перехода от INT к продуктам реакции PR описывалась с координатой реакции, отвечающей расстоянию от протона аминокислотного остатка His438 до атома кислорода (O) серина. Последовательное уменьшение этого расстояния с оптимизацией всех остальных геометрических параметров системы по методу КМ/ММ позволило локализовать конфигурацию второго переходного состояния TS2. Спуск с седловой точки TS2 вдоль потенциальной поверхности приводит к структуре интермедиата INT (для проверки правильности расчета) и к конфигурации продуктов реакции с регенерированной боковой цепью Ser198. На рис. 3 показаны соответствующие превращения в активном центре. Как и для первой стадии реакции, геометрические
конфигурации локализованных стационарных точек (Т82 и РЯ), вычисленные с использованием методов ЭРТ и ЭР-ТВ, достаточно близки. Однако величины энергий барьеров различаются весьма заметно.
По результатам данной работы можно сформулировать следующие выводы. Молекулярное моделирование, основанное на теории КМ/ ММ, позволяет воспроизвести эффект точечной мутации С1у117И18 в реакции реактивации бу-тирилхолинэстеразы, фосфорилированной по аминокислотному остатку каталитического се-рина. Применение полуэмпирического кван-тово-химического подхода функционала плотности в приближении сильной связи в рамках модели КМ/ММ хорошо передает структуры интермедиатов и переходных состояний вдоль реакционного пути, но приводит к определенным различиям в энергетических диаграммах.
Авторы выражают благодарность сотрудникам суперкомпьютерных центров МГУ имени М.В. Ломоносова и РАН за возможность использовать вычислительные ресурсы.
DFT DF -ТВ 1- i
wt-BChE 19 <1 i i
Glyll7His 18 <1 i i
Рис. 3. Преобразования на второй стадии реакции INT ^ TS2 ^ PR (значения высоты энергетических барьеров в таблице приведены в ккал/моль)
При написании данной статьи использованы работы, поддержанные грантом Российского научного фонда (проект № 14-13-00124).
СПИСОК ЛИТЕРАТУРЫ
1. WarshelA., LevittM. // J. Mol. Biol. 1976. Vol. 103. P. 227.
2. Немухин А.В., Григоренко Б.Л., Лущекина С.В., и др. // Успехи химии. 2012. Vol. 81. С. 1011.
3. Nemukhin A., Lushchekina S., Bochenkova A. et al. // J. Molec. Modeling. 2008. Vol. 14. P. 409.
4. Nemukhin A.V., Grigorenko B.L., Morozov D.I., et al. // Chem. Biol. Interact. 2013. Vol. 203. P. 51.
5. Lockridge O., Blong R.M., Masson P. et al. // Biochemistry. 1997. Vol. 36. P. 786.
6. Yao Y, Liu J., Zhan C.G. // Biochemistry. 2012. Vol. 51. P. 8980.
7. Elstner M., Frauenheim T., Suhai S. // J. Mol. Struct. (Theochem). 2003. Vol. 632. P. 29.
8. Nachon F., Asojo O.A., Borgstahl G.E.O. et al. // Biochemistry. 2005. Vol. 44. P. 1154.
9. Nachon F., Carletti E., Wandhammer M. et al. // Biochem J. 2011. Vol. 434. P. 73.
10. Zhao Y., Lynch B.J., Truhlar D.G. // J. Phys. Chem. A 2004. Vol. 108. P. 2715.
11. ValievM., Bylaska E. J., GovindN. et al. // Comput. Phys. Commun. 2010. Vol. 181. P. 1477.
12. Case D.A., Darden T.A., Cheatham T.E. et al. // AMBER 12. 2012. University of California, San Francisco.
Поступила в редакцию.01.08.15
MODELING CHEMICAL TRANSFORMATIONS AT THE ACTIVE SITES OF CHOLINESTERASES BY QUANTUM-BASED SIMULATIONS
A.V. Nemukhin, A.M. Kulakova, S.V. Lushchekina, A.Yu. Ermilov, S.D. Varfolomeev
(Division of Physical Chemistry)
Significance of the quantum mechanical - molecular mechanical (QM/MM) method in modeling chemical transformations at the active sites of cholinesterases is discussed. Diverse versions of the QM/MM approach are applied to understand molecular mechanisms of the reactivation reaction of the phosphorylated by the catalytic serine residue butyrylcholinesterase.
Key words: cholinesterases, phosphorylated BChE, molecular modeling, QM/MM method.
Сведения об авторах: Немухин Александр Владимирович - профессор кафедры физической химии химического факультета МГУ, докт. хим. наук ([email protected]); Кулакова Анна Михайловна - аспирант химического факультета МГУ ([email protected]); Лущекина Софья Владимировна - ст. науч. сотр. Института биохимической физики РАН, канд. хим. наук ([email protected]); Ермилов Александр Юрьевич - доцент кафедры физический химии химического факультета МГУ, канд. хим. наук ([email protected]); Варфоломеев Сергей Дмитриевич - зав. кафедрой химической энзимологии, профессор химического факультета МГУ, научный руководитель Института биохимической физики РАН, чл.-корр. РАН, докт. хим. наук ([email protected]).