Роль структурных особенностей каталитического действия силикатеина-А в процессе биосилификации по данным компьютерного моделирования

В работе использовались молекулярные модели силикатеина-α, структура которых разрешена экспериментально методом ренгеноструктурного анализа, это: модель (pdb 6ZQ3) — силикатеина-α получена из организма T. аurantium (Görlich et al., 2020) и модель химерного белка «4SER» катепсина L человека (pdb 2VHS), приближенная к структуре силикатеина-α и обладающая максимальной кремний полимеризующиеся активностью (Fairched et al., 2009). И модель пре-силикатеина-α предсказанная методом ab-initio по первичной последовательности из T. аurantium (код NCBI: AAC23951.1). Модель 2VHS состоит из четырех гомологичных доменов, в работе использовался домен А.
Лигандами выступали кремниевые кислоты, являющиеся природными субстратами для силикатеина-α с различной длиной цепочки Sin=1-6 и степенью ветвления.
1.2. Конформационные особенности силикатеина-α для начального этапа олигомеризации Si(OH)4
1.2.1. Подготовка молекулярной модели структуры силикатеина-α. Подготовка молекулярной модели силикатеина–α (pdb 2VHS) проходила в программе Gromacs-3.3.1 (Leach, 2001) и включала: очистку координат от примесей, центрирование, растворение полученной структуры в кубе воды размером 77,86*67,01*63,86 Å с потенциалом Amber99 (Zwart et al., 2007). Минимизация и МД проводилась в программе Gromacs, согласно протоколу (Yang et al., 2007), в течение 3нс, при температуре 300К, с учетом модели растворителя PCM (Fennell and Gezelter, 2006) в двух повторах. Структура модели фермента, оценивалась по изменению следующих параметров: среднеквадратичное отклонение, поверхность доступная растворителю, а также подвижности аминокислот активного центра. Полученные данные анализировали с помощью программного пакета VMD (Humphrey et al. 1996). Использовался разработанный подход, который был апробирован для моделирования взаимодействия различных макромолекул с лигандами, в том числе и для изучения взаимодействия линолевой кислоты с ДНК в водном растворе (Тарасов с соавт., 2012).
1.2.2. Подготовка модели субстрата — Si(OH)4. Молекула Si(OH)4 была построена в программе ChemCraft (Chemcraft, 2019). Оптимизация геометрии молекулы проводилась в программе PCGAMESS/Firefly (Gordon et al., 2005), методом теории функционала плотности (DFT) с функционалом X3LYP (Xu and Goddard, 2004) в валетно-расщепленном базисе 6- 31G(d) с последующим пересчетом структуры в базисе 6-311++G(2d, p).
1.2.3. Подготовка 8 конформационных изомеров модели каталитической триады силикатеина-α с субстратом Si(OH)4. Для квантово-химических расчетов взаимодействия Si(OH)4 с АЦ силикатеина–α (pdb 2VHS) были использованы координаты АО остатков,
входящих в состав активного центра (His163, Ser25, Gln19) (АЦ), полученные на последнем шаге МД. Для имитации соединения аминокислот АЦ с полипептидной цепочкой была проведена фиксация атомов С и N, входящих в состав пептидной связи по трем пространственным координатам (X, Y, Z). Стационарные состояния каталитической триады рассчитывались методом теории функционала плотности (DFT) в функционале X3LYP (Xu and Goddard, 2004) в валентно-расщепленном базисе 6-311+G(d,p)/6-31G(d), с учетом модели растворителя PCM. Расчеты производились в программе PC-GAMESS/Firefly (Gordon et al., 2005). Для расчета частичных атомных зарядов использовалась зарядовая схема Малликена.
1.2.4. Построении поверхности потенциальной энергии переноса H+ вдоль донорно акцепторной пары атомов. Вдоль вектора перемещения H+ между донорно- акцепторной парой каталитическая триада-Si(OH)4, была рассчитаны энергии системы. Расчеты проводились в программе PCGAMESS/Firefly (Gordon et al., 2005) методом DFT, в функционале X3LYP (Xu and Goddard, 2004) в валентно-расщепленном базисе 6- 311++G(2d,p), в газовой фазе и с использованием модели поляризованного континуума (PCM). Энергетический барьер прямой и обратной реакций перемещения H+ рассчитан относительно начальной конформации комплекса. Молекулярные структуры комплексов были подготовлены с использованием языка программирования Python (Python, 2021) и методика успешно опробована для построения оптимальных траекторий передачи водорода в позиционно управляемом механосинтезе (Tarasov et al., 2013).
1.3. Поиск областей связывания силикатеина-α с избранными олигомерами кремниевых кислот
1.3.1. Подготовка молекулярных моделей силикатеина-α. Для подготовки модели структуры силикатеина-α была проведена МД в пакете NAMD- multicore CUDA (Phillips et al., 2005) с использованием силового поля CHARMM36 (Brooks et al., 2009). Размер периодической ячейки составил 92 * 87 * 74 Å. Подготовка к МД в NAMD (Phillips et al., 2005) включала следующие этапы: минимизацию (2 нс), нагрев (до 310К, 0,3нс) и эквилибрацию (2нс). МД составила 25нс с шагом итерации 1фс. Координаты структуры силикатеина-α последнего шага МД были оценены на наличие запрещенных углов ψ и φ по карте Рамачандрана (UCSF Chimera, версия 1.15), рассчитаны среднеквадратичное отклонение (RMSD) всей структуры, среднеквадратичное отклонение аминокислот активного центра (RMSF), а так же проведена оценка изменения расстояний между выбранными парами атомов. RMSD, RMSF и динамика расстояний рассчитаны в VMD, визуализированы с использованием Python, модуль matplotlib .
Подготовка модели структуры пре-силикатеина-α заключалась в построении пространственных моделей методом ab-initio по первичной последовательности (код AAC23951.1) пре-силикатеина-α на веб-сервиса Robetta (Kim et al., 2004) и QUARK. Было пересказаны 20 моделей структуры. Оптимальная модель пре-силикатеина-α была выбрана на основании лучшей QMEAN-оценки и отсутствию запрещенных аминокислот на карте Рамачандрана.
1.3.2. Подготовка моделей природных субстратов силикатеина-α. Начальная моделируемая молекулярная система представлена равноудаленными молекулами Si(OH)4, расположенными в периодической ячейке с размерами: 58,44 х 58,44 х 58,44Å, заполненная молекулами воды (TIP3P) с плотностью 1 г/см3.
Исследовались две системы:
система45 – модель с концентрацией ортокремниевой кислоты 0,45 г/см3, содержащая 576 мономеров Si(OH)4 и 5220 молекул H2O, всего атомов 20844.
система56 – модель с концентрацией ортокремниевой кислоты 0,56 г/см3, содержащая 729 мономеров Si(OH)4 и 3026 молекул H2O, всего атомов 15638.
Молекула Si(OH)4 была предварительно подготовлена, согласно пункту (1.2.2).
Эквилибрация системы проводилась при температуре 300К в течение 10пс. Моделирование выполнялось при температуре 2500К, что является частью стандартной методики (Hirakawa et al., 2017; Stawick et al., 2020) для ускорения процесса олигомеризации кремниевых кислот. Интегрирование уравнений движения производилось с помощью алгоритма Верле с шагом 0,001 пс. Траектория движения атомов записывалась с шагом 0,1пс. Продолжительность МД составила 3 нс, что позволяет исследовать начальные этапы олигомеризации (Rao and Gelb, 2004).
В настоящей работе МД кремниевых кислот в водном растворе проводилась с использованием силового поля FG (Feuston and Garofalini, 1990), реализованного в программном пакете GPAMP (Tarasov et al., 2010).
Для обработки большого массива полученной информации соискателем разработан комплекс программ на языке программирования Python (на часть кода получено свидетельство о государственной регистрации программы для ЭВМ No 2021619675).
Распределения длин ковалентных связей и валентных углов были проверена на нормальность распределения, используя критерий Шапиро-Уилка и K’2-тест Д’Агостино. Для сравнения двух выборок применялся U-критерий Манна-Уитни, для множественного сравнения данных – тест Крускала-Уоллиса (Python, библиотека SciPy). При разработки задействованы библиотеки Python: json, matplotlib, scipy, numpy и statistics, NetworkX.
Избранные олигомеры (14) с количеством атомов Sin = 1-6 были оптимизированы с использованием программного пакета ORCA (Neese, 2017), методом теории функционала плотности (DFT) с гибридным функционалом B3LYP и базисным набором 6-311G.
1.3.3. Молекуляный докинг олигомеров кремниевых кислот в силикатеин-α.
Области взаимодействия силикатеина–α с молекулами кремниевых кислот оценивались по результатам молекулярного докинга с гибким лигандом. Поиск областей связывания фермента и лиганда проходил в два этапа:
Глобальный поиск по всей поверхности: лигандом выступала молекула Si(OH)4 – Q01 (подготовка раздел 1.2.1). На основании проведенного исследования были выделены три потенциально значимые области.
Локальный поиск по трем областям связывания был проведен с 14-тью лигандами (подготовка раздел 1.3.2), с радиусом поиска от выбранной точки — 30Å. Для каждой области методом Монте-Карло было сгенерировано 500 положений. Для визуализации областей связывания были использованы 20 состояний с наименьшей энергией связывания. Для проведения молекулярного докинга использовались сервисы BSP-SLIM (Deluca et al., 2015) и ROSIE (Lee and Zhang, 2012; Lyskov et al., 2013) .
2 РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ
2.1. Конформационные особенности силикатеина–α для начального этапа олигомеризации Si(OH)4
2.1.1. Квантово-химические расчеты молекулы Si(OH)4. В оптимизированной молекуле субстрата значение длины связи между О-Н и Si-O составили 0,97Å и 1,65Å, соответственно, значения углов разделяются на 3 типа и составляют: 106°, 116,5° и -113,5°, что соответствует как экспериментальными, так и расчетным данным (Fueston and Garofalini, 1990; Rimsza et al., 2017).
2.1.2. Оценка подвижности силикатеина–α в ходе МД. В ходе МД RMSD мутантного силикатеина-α выходит на плато с 0,6нс и колебания составили не более 1,5Å. Наблюдается незначительное колебание АО АЦ в ходе МД. Изменение площади треугольника, образованного между АО АЦ: τN His, O — радикала Ser (O-Ser) и О — радикала Gln (O-Gln), то средняя площадь треугольника составляет 6,5Å2, при длинах τN His — O-Ser в 5Å; O-Ser — O-Gln в 5,1Å и τN His — O-Gln в 9,6Å.
2.1.3. Характеристика стационарных состояний каталитической триады силикатеина–α с и без Si(OH)4. Были сконструированы и оптимизированы молекулярные комплексы каталитическая триада – Si(OH)4, отличающиеся степенью протонирования азотов имидазольного кольца His и его вращения относительно связи Сα–Сβ (рис. 1). Результаты суммированы в таблице 1. Наименьшей энергией связывания обладает полностью депротонированная форма (состояние 1) с энергией связывания в -32,82 ккал/моль. Для этого варианта характерны 4 водородные связи, между каталитическим центром и субстратом. В ходе оптимизации H+ от молекулы Si(OH)4 переходит на πN, что в дальнейшем может привести к конденсации молекул Si(OH)4 по основному типу (Pereira et al., 1998). Для состояний 5 и 3 наблюдается образование пяти водородных связей между каталитической триадой и субстратом, а также незначительное удлинение ковалентных связей (Si(OH)4 – – O- H – – Nn) на 0,063Å и 0,072Å, соответственно. Состояние 6 интересно тем, что поворот имидазольного кольца приводит к формированию пяти водородных связей между каталитической триадой и субстратом, энергия связывания для этого состояния максимальна
среди остальных и составляет -12,32 ккал/моль.
Рис. 1. Молекулярная модель силикатеина-α и его АЦ (His163, Ser25, Gln19) с помещенным в него субстратом Si(OH)4, все атомы пронумерованы.
Для состояний 2 и 8 комплекса фермента с субстратом образуется цепочка водородных связей, приводящая к небольшому удлинению ковалентных связей. Так в состоянии 2

11
образуются водородные связи His πN-H18 – – O2 Si(OH)4 и Si(OH)4 H6 – -O9 Gln, длина которых 1,96Å и 1,82Å, соответственно. А в состоянии 8 водородные связи выстраиваются следующем образом: His τN-H16 – -O3 Si(OH)4 и Si(OH)4 H7 – -O9 Gln, длина которых 1,65Å и 1,61Å, соответственно. Такое расположение приводит к удлинению ковалентной связи τN- H16 на 0,05Å. И потенциально может быть исходной точкой для образования молекулы воды, удерживаемой какое-то время атомом кремния молекулы Si(OH)4.
Состояния 5 и 6 интересны наличием коротких водородных связей между His τN – – H- O-Si(OH)3, равных 1,64Å и 1,81Å, соответственно. Такое расположение может быть исходной позицией переноса протона от молекулы Si(OH)4 на τN His.
Для состояния 3 наблюдается последовательность водородных связей, в которую вовлечена гидроксильная группа Ser: His πN – -H13 Ser и O6 Ser – -H7 Si(OH)4, длина которых 1,75Å и 1,82Å, соответственно. По данной стартовой конфигурации комплекса можно предположить последовательный переход H+ гидроксила Ser на πN His, что в свою очередь может привести к переносу Н+ с Si(OH)4 на Ser. В результате чего кислород может совершить нуклеофильную атаку на атом кремния другой молекулы Si(OH)4.
Таблица 1. Сводная таблица по 8 состояниям комплекса каталитическая триада силикатеина- α и Si(OH)4 .
Парциальные заряды (АЦ)
Парциальные заряды (АЦ и Si(OH)4)
No Модель изменяющейся
Протонац ия
нет πN
τN
τN, πN
нет
πN
τN
τN, πN
Е, Kp ккал/мол
ь
-32,82 2,681*1024
-26,69 7,339*1019
-28,15 8,961*1020
-25,94 2,029*1019
-26,78 8,564*1019
-12,32 2,11*1012
-16,53 1,479*109
-21,92 2,067*1016
структуры
1
3
4 5
7 8
πN τN -0,522 -0,502
-0,635 -0,427 -0,476 -0,569
-0,525 -0,50
-0,611 -0,480
-0,482 -0,527 -0,613 -0,559
πN τN -0,624 -0,459
-0,631 -0,434 -0,561 -0,561
-0,54 -0,512 -0,606 -0,488
-0,572 -0,595 -0,605 -0,575
Состояние 8 характеризуется наличием 6 водородных связей между каталитической триадой и субстратом, по 2 на каждый аминокислотный остаток, с энергией связывания в – 21,92 ккал/моль. Состояние 4 не является стационарным, поскольку в ходе оптимизации конформация состояния 4, вращаясь вокруг связи Сα—Сβ, переходит в конформацию, близкую состоянию 8 с энергией связывания в -25,94 ккал/моль.
В целом, протонированная форма азота имидазольного кольца His имеет больший (по модулю) парциальный заряд. В случае равного протонирования атомов πN и τN, πN является более электроотрицательным атомом. Данное соотношение характерно как для свободного
АЦ, так и в комплексе с Si(OH)4 (таблица 1). Потенциальными состояниями, в которых возможен начальный этап активации молекулы ортокремниевой кислоты, могут быть состояния 2, 3 или 6. Наибольшая разница между двумя электроотрицательно заряженными атомами азота имидазольного кольца His наблюдается для частично протонированных состояний 2 и 6, в обоих случаях протонирован πN имидазольного кольца His.
2.1.4. ППЭ переноса H+ в пределах электроотрицательной пары атомов АО каталитической триады силикатеина–α и молекулы Si(OH)4. По данным оптимизации комплекса каталитической триады с молекулой Si(OH)4 можно выделить конфигурации АЦ, которые потенциально могут запускать активацию молекулы кремниевой кислоты и процесс конденсации. Полученные конфигурации можно разделить на группы по аминокислотным остаткам, которые могут инициировать активацию Si(OH)4. Первая – группа серина, характерна для состояния 3. Вторая – группа гистидина. Последнюю можно разделить на две подгруппы: протонированный азот His — потенциальный донор H+ (состояния 2 и 7) и депротонированный азот His — потенциальный акцептора H+ (состояния 5 и 6). Для каждого из вариантов построена ППЭ переноса Н+ и рассчитаны энергетические барьеры в пределах электроотрицательных атомов, образующих водородную связь.
В геометрии состояния 2 (протонированный πN имидазольного кольца, πN расположен кнаружи) после оптимизации расстояние между электроотрицательными атомами His πN и O Si(OH)4 составило 2,79Å. ППЭ перемещения Н+ (рисунок 2а, черный цвет) была получена путем генерации серии позиций вдоль вектора His πN→O. При пересчете системы с учетом PCM на ППЭ присутствуют два минимума, разделенные энергетическим барьером. Первый минимум расположен в области 1,1Å от πN His, второй на расстоянии 0,98Å от O–Si(OH)4. Энергетический барьер перемещения Н+ от πN His к O–Si(OH)4 составляет 14,2 ккал/моль, в то время как энергетический барьер обратной реакции низок и составляет 1,7 ккал/моль. Общая энергия системы при расположении перемещаемого Н+ вблизи O–Si(OH)4 составляет на 12,5 ккал/моль больше, чем расположение H+ вблизи πN His. Второй минимум энергии системы вблизи O–Si(OH)4 позволяет предположить, что существует определенная вероятность нахождения там H+. Низкий барьер обратной реакции позволит H+ с легкостью вернуться к πN His.
а)
б)
Рис. 2. ППЭ переноса H+ между электроотрицательной парой атомов His–N и O–Si(OH)4. Энергии представлены относительно позиции с минимальной энергией.
Для состояния 7 (протонированный τN имидазольного кольца His, τN расположен кнаружи) расстояние между электроотрицательными атомами His τN и O–Si(OH)4 составило
2,9Å. ППЭ перемещения H+ (рис. 2а, серый цвет) была получена путем генерации 11 позиций расположения H+ (His NoH16) вдоль вектора τN→O. В системе присутствует один общий минимум, положения H+, на расстоянии в 0,96Å от τN His. Минимум вблизи O2–Si(OH)4 не образуется ни при расчете в вакууме, ни с использованием расчетов в PCM.
Для состояний 5 и 6 построена ППЭ переноса Н+ от О–Si(OH)4 на τN имидазольного кольца His (рис. 2б). ППЭ состояния 5 при пересчете с учетом мнимого растворителя имеет 2 минимума, соединенных энергетическим барьером прямой реакции 4,36 ккал/моль и обратной 3,23 ккал/моль. Разница в энергиях между двумя локальными минимумами составляет 1,13 ккал/моль, что может быть отнесено к погрешности метода. ППЭ состояния 6 при пересчете с учетом мнимого растворителя имеет 2 минимума энергии: первый – на расстоянии 0,98Å от О–Si(OH)4, второй – на расстоянии 1,15Å, от τN–His. Энергетический барьер прямой реакции переноса H+ составляет 7,82 ккал/моль, энергетический барьер обратной – 0,3 ккал/моль. Низкий барьер обратной реакции позволит H+ с легкостью вернуться к О—Si(OH)4. На основании рассмотренных профилей ППЭ переноса H+ потенциальными состояниями, в которых возможен начальный этап активации молекулы Si(OH)4, могут быть состояния 2 и 6, хотя оба характеризуются низким энергетических барьером для обратной реакции. Состояние 1 и 5 могут потенциально являться начальной стадией активации кремниевых кислот, однако при нейтральном значении pH это маловероятно (Zhang et al., 2020). Тем не менее такой вариант может быть возможным, если рассматривать эти формы как короткоживущие интермедиаты, которые работают как части челночных систем переноса H+ подобно классическому механизму работы сериновых протеаз (Blow et al., 1969; Fodor et al., 2006). Следует также отметить, что на ППЭ переноса Н+ могут оказывать влияние АО, не включенные в модель рассматриваемого АЦ. Так присутствие Asp187 может понизить кислотность гистидина (Zhou et al., 1999; Schröder et al., 2012).
2.2. Области связывания силикатеина–α с кремниевыми кислотами.
2.2.1. Олигомерное разнообразия кремниевых кислот в водном растворе.
Природными субстратами для фермента силикатеина-α являются кремниевые кислоты, растворенные в водном растворе. Вопрос олигомерного состава и изомерного разнообразия до сих пор остается открытым (Sjöberg et al., 1996; Cho et al., 2006; Borba et al., 2017). В работе проводился поиск олигомерного разнообразия природных субстратов силикатеина-α, для этого методом МД в реализованном реакционном силовом поле FG (Fueston and Garofalini, 1990). Анализ изменения доли классов олигомеров Sin=1-8 в МД показывает, что в первые 0,5нс динамики происходит резкое увеличение количества небольших олигомеров. В дальнейшем их концентрация начинает плавно уменьшаться в течение всей траектории в обеих системах. В обеих системах доля димеров в ходе МД остается достаточно высокой на всем протяжении МД. Для системы45 доля димеров колеблется около∼0,6, в для системы56 около∼0,525. В системе45 в ходе МД наблюдается уменьшение доли тримеров и медленный рост количества олигомеров Sin=4−9. В системе56 после резкого образования олигомеров Sin=3−5 наблюдается медленный их спад и рост числа структур с количеством атомов кремния от 6 до 16. Большая доля кремниевых кислот приходится на небольшие олигомеры до 9 (система45) и до 16 атомов кремния (система56). За весь период МД кремниевых кислот в водном растворе обнаружены структуры с длиной цепочки до 20 (система45) и 27 атомов кремния (система56).
Общий профиль изменения долей рассмотренных олигомеров кремниевых кислот сопоставим с результатами спектроскопии ЯМР 29Si и спектроскопии УФ-комбинационного рассеяния (Depla et al., 2011), а так же с результатами моделирования олигомеризации кремниевых кислот при нейтральном значении pH (Rao and Gelb, 2004). Авторы работ (Rao and Gelb, 2004; Depla et al., 2011) высказывают предположение, что именно димеры на этом временном промежутке играют ключевую роль в процессах олигомеризации.
Форма олигомеров сильно зависит от условий среды, в которой происходит эксперимент. Возрастание доли циклических форм наблюдается при повышении давления до 0,9 ГПа, присутствия ионов аммония и их производных (Chen et al., 2019; Zhang et al., 2011), катионов тетраалкиламмония (ТАА) и тетраметиламмония (TMA) (Knight et al., 1989), а также использования насыщенных щелочных растворов (Depla et al., 2011; Pilgrim et al., 2018). Предпочтительное образование линейных и разветвленных форм наблюдается при нейтральном и сильно кислом pH и низком парциальном давлении (Lazaro et al., 2017).
Было проведено изучение разнообразия геометрических изомеров образующихся в ходе МД олигомеров и подсчет наиболее высокочастотных вариантов. Доминирующее место занимают линейные и разветвленные структуры, среди двадцати наиболее часто встречаемых изомеров отсутствуют циклические формы (рис. 3). Картина распределения первых 9 олигомеров сопоставима для обеих систем. По результатам МД циклические формы в системе присутствуют в низком процентном соотношении, не превышающие 0,01%, для Sin = 3-6, для олигомеров с большим числом атомов не зафиксированы (рис. 4).
Рис. 3. Двадцать наиболее высокочастотных графовых форм изомеров кремниевых кислот, для системы45 (а) и системы56 (б), за 1 принята доля линейного тримера (3L).
Рис. 4. Доля циклических фрагментов среди
В результате сформирована выборка из 14 олигомеров кремниевых кислот (лигандов) с количеством атомов Sin= 1-6 (рис. 10), которые были оптимизированы и использованы для проведения молекулярного докинга в двух моделях силикатеина-α. Для кремниевых кислот с количеством атомов Sin = 1-8 проведено полное исследование изомерного разнообразия образующихся пространственных изомеров кремниевых кислот, результаты отражены в (Изотова с соавт., 2020). В проведенном исследовании

олигомеров кремниевых кислот, Sin = 3-6. обнаружены ранее не идентифицированные олигомеры.
Полученное частотное разнообразие сопоставимо как с данными экспериментов (Cho et al., 2006; Borba et al., 2017; Chen et al., 2019), так и с модельными расчетами (Catlow et al., 1998; Belton et al., 2012; Pilgrim et al., 2018).
2.2.2. Оценка качества и подвижности молекулярной моделей силикатеина-α и пре-силикатеина-α. Для оценки качества модели использовали метод МД. Продолжительность МД силикатеина-α составила 25нс. Изменение среднеквадратичного отклонения (RMSD) всего белка с течением времени выходит на плато на 10нс. Медиана RMSD составила 4,46Å (3,77 — 5,17Å) для диапазона 10-25нс (рис. 6а). Карта Рамачандрана для последнего шага МД (рис. 5а) показывает, что большинство двугранных углов пептидного остова (ψ и φ) лежат в разрешенных зонах: His163 в области β- листов, Ser25 и Gln19 в области правой α-спирали. Динамика изменения среднеквадратичной флуктуации АО (RMSF) активного центра (рис. 5б) на интервале 10-25нс показывает, что Gln19 и Ser25 менее подвижны на всем интервале МД с медианой 2,72Å (2 – 3,4Å) и 2,78Å (2,2 – 3,5Å), в то время как медиана RMSF His163 4,3Å (3,5 – 5,2Å). Динамика изменений расстояний между функциональными группами АО АЦ (рис. 6б) демонстрирует удаление Gln19 от Ser25, с ~4 Å до ~10Å с одновременным сближение Сα атом His163 с О–Ser25 от ~6Å до ~3Å, что способствует формированию устойчивой водородной связи His163- πN и НО-Ser25.
Конформационная подвижность His163 АЦ силикатеина-α может играть значимую роль в инициации активации кремниевых кислот с различной длиной цепочки (Szőri- Dorogházi et al., 2012). Была проведена оценка изменения расстояний между О–Ser25 и тремя атомами His163, это: Сα, τN и πN (рис. 6а). На всей протяженности МД πN-His163 остается ориентированным к О–Ser25, несмотря на частичное вращение плоскости имидазольного кольца His163. Колебание величины водородной связи His163-πN—НО-Ser25 находится в диапазоне 2-4Å (рис. 6а), что согласуется с экспериментом – 3,6Å (Görlich et al., 2020).
а) б)
Рис. 5. Динамика изменения а) RMSD белка и б) RMSF отдельных АО в течение 25нс МД.
В целом полученная структура силикатеина-α является хорошо подготовленной для проведения молекулярного докинга, аминокислоты АЦ уравновешены и располагаются в допустимых зонах двугранных углов пептидного остова (ψ и φ). АО, которые
предположительно осуществляют катализ, с RMSF в пределах 3,5Å.

Предположительно, в активации кремниевых кислот может участвовать система близко расположенных His, образующих цепь переноса H+. Пространственная структура силикатеина-α имеет область близко расположенных His, это His163, His140, His205. Динамика изменения расстояний между этими His лежит в диапазоне ~11-14Å, что является большим для миграции H+ в пределах одной белковой молекулы, несмотря на то, что для некоторых растворов это расстояние может достигать 10-30Å (Kailaa and Hummera, 2011).
В структурах силикатеина-α (pdb 6ZQ3 и 2VHS) отсутсвуют некоторые АО, поэтому по первичной последовательности была построена и оптимизирована модель пре- силикатеина-α. Предсказанная полная модель пре-силикатеина-α показала хорошую оценку качества, сопоставимую с оценками экспериментальных моделей. Структурное выравнивание предсказанной структуры пре-силикатеина-α согласуется с экспериментально полученными структурами (pdb 6ZQ3 (Görlich et al., 2020) и pdb 2VHS-А (Fairched et al., 2009). Кроме того предсказанная модель является более полной: в ее структуре представлен пропептидный участок и определены координаты атомов всех аминокислотных остатков зрелого белка, в структуре 2VHS-А отсутствует участок Ser174-Gln180, в структуре 6ZQ3 нет участка Ile174- Gln179 (рис. 7).
а) б)
Рис. 6. а) Схема пар расстояний между выбранными атомами; б) – Динамика изменения расстояний избранных пар атомов.
Выравнивание предсказанной структуры пре-силикатеина-α и прокатепсина L (pdb 1CS8), показало хорошее наложение пропептидного участка в районе, содержащем катепсиновый пропептид-ингибиторный домен I29 (рис. 8) (позиции аминокислот 29 – 117).
Рис. 7. Структурное выравнивание предсказанной модели пре-силикатеина-α, pdb 6ZQ3 (Görlich et al., 2020) и pdb 2VHS-А (Fairched et al., 2009).
В целом предсказанная структура пре-силикатеина-α является хорошо подготовленной и
Рис. 8 Структурное выравнивание предсказанной моделипре-силикатеина-αипрокатепсинаL(pdb1CS8).
2.2.3. Поиск областей связывания силикатеина-α с избранными лигандами. Поиск областей связывания силикатеина-α с 14 подготовленными кремниевыми кислотами провился в два этапа. По результатам проведения глобального поиска по всей поверхности структуры пре-силикатеина-α, можно выделить 3 области связывания: Asp68, Thr329 и Ser131, для структуры силикатеина-α были выбраны схожие регионы: Asp6, Thr94 и Ser25 (рис. 9). Для пре-силикатеин-α вблизи Asp329 выделяется одна область взаимодействий, в независимости от типа и размера лиганда. Для областей вблизи Ser131 для лигандов с количеством атомов Sin = 1-4 и Thr329 для лигандов Q12, Q12Q21, Q22Q32, Q23Q32 выявлено по две области. Для силикатеина-α вокруг каждого центра формируется одна область взаимодействия.
Наибольшие значения энергии для каждого вида лигандов представлены на рисунке 10. Максимальная энергия связывания в трех регионах докирования для большинства лигандов приходится на область вблизи активного центра пре-силикатеина-α (Ser131, с
сопоставима как экспериментальными ренгеноструктурного силикатеина-α (pdb 6ZQ3) (Görlich et al., 2020) и химерным белком «4SER» катепсина L человека (pdb 2VHS) (Fairched et al., 2009), так и с прокатепсина L (pdb 1CS8) и может быть использована для проведения молекулярного докинга.
с данными анализа
максимумом -16,1 ккал/моль). Исключение составляют Q25, Q23Q32 и Q12Q24, для которых регионом с максимальной энергией связывания является область вблизи Thr329 (с максимумом -12,55 ккал/моль). Среди групп изомеров максимальными энергиями связывания обладают комплексы фермент-линейный или разветвленный лиганд для областей Ser131 и Thr329 за исключением группы лигандов Sin = 4. Для лигандов с длиной цепочки Sin = 1-4 (по абсолютным значениям) и для Sin = 1-6 по медианам наблюдается прямая корреляция размера лиганда и энергии связывания вблизи трех выделенных областей. Область вблизи Asp68 является наименее предпочтительной, но тем не менее энергия связывания находится в
диапазоне от -3,6 до -9,1 ккал/моль для всех рассмотренных лигандов.
Рис. 9. Области связывания лигандов с а, б) пре-силикатеином-α и с в,г) силикатеином-α
Рис. 10. Максимальные энергии связывания для каждого типа лигандов с пре-силикатеином-α (верхний), силикатеином-α (нижний).
Для силикатеина-α лиганды с количеством атомов Sin = 1-5 лучше связываются в области активного центра (Ser25). Среди изомеров Sin = 4-5 максимальная энергия связывания приходится на линейную и разветвленную форму лигандов с максимумом для этой области в -10,96 ккал/моль. Область связывания вблизи Asp6 интересна тем, что для изомеров из групп Sin = 3-5 максимальной энергией связывания обладают комплексы фермента с циклическими изомерами лигандов с разницей энергии в ~1,4 ккал/моль. Энергия

связывания лигандов в вблизи Thr94 значительно ниже, чем с другими регионами, но тем не менее максимальные значения располагаются в диапазоне -5,7-8,4 ккал/моль. Полученные значения для лигандов с длиной цепочки Sin=5-6 носят предварительный характер в связи с тем, что рассмотрены лишь крайние структуры в ряду гибкости, это
циклическая и линейная, обладающая максимальной конформационной подвижностью.
В целом связывание лигандов с силикатеином-α и с пре-силикатеином-α с максимальной энергией происходит вблизи АЦ. Лиганды большего размера Sin=5-6 начинают эффективно связываться с областью вблизи Asp6 (для лиганд—силикатеин-α) и с
областью рядом с Thr329 (лиганд—пре-силикатеин-α).
В стабилизации лигандов в областях связывания принимают участие АО,
относящиеся к разных классам по строению радикала. На рисунке 11 (а, б) суммированы данные взаимодействия лигандов с АО пре-силикатеин-α и силикатеин-α. В целом следует отметить, что с увеличением размера лиганда область взаимодействия расширяется. Для пре- силикатеин-α области вблизи активного центра для малых олигомеров разделены на 2 части, которые с увеличением размеров сливаются в одну область с общими АО: Glu132, Ser131 (АО АЦУ), Asp134, Lys130, Ser298 (рис. 11а). Для силикатеин-α регион вблизи активного
центра и Thr94 имеют общие АО, это: Asn66, Cys65 (рис. 11б).
Рис. 11. Поверхность пре-силикатеин-α (а) и силикатеина-α (б), цветами раскрашены области связывания с лигандами с длиной цепочки Sin=1-6.
2.2.4. Белок-белковый докинг силикатеина-α. Силикатеин-α выполняет 2 функции: ферментативную и структурообразующую (Wang et al., 2014; Riesgo et al., 2015). По данным рентгеноструктурного анализа кристаллов силикатеина-α (Fairhead et al., 2008) и комплекса силикатеина-α с силикатами в составе спикулы (Görlich et al., 2020) получены противоречивые данные, и вопросы стыковки функциональных единиц фермента остаются открытыми. В связи с этим становится важным рассмотреть сайты взаимодействия силикатина-α с идентичной субъединицей. Все полученные модели, обладающие максимальной аффинностью, докируются с одной стороны и разделяются на две группы (рисунок 12а, б): первая – силикатеин-α _A-А1 — А2 (две модели); вторая-силикатеин-α _A- А3 — А9 (семь моделей).
АЦ первой группы (силикатеин-α_A-А1—А2) располагаются перпендикулярно друг к другу, в то время как для второй группы силикатеин-α_A-А3—А9 — АЦ располагаются по

одну сторону. АО, принимающие участие в формировании комплекса рецептор-лиганд, располагаются в одной области молекулы силикатеина-α (рис. 13).
Рис. 12. Взаимное расположение двух субъединиц силикатеина-α (рецептор и лиганд), где лиганды: бирюзовая — первая группа, салатовая — вторая, рецептор: сине-крассно-белая.
Рис. 13. АО областей связывания рецептора (силикатеин-α) и лигандов (силикатеин-α), где а) — поверхность рецептора, б) — поверхность лигандов. Условные обозначения: голубой цвет — АО первой группы, зеленый цвет — второй группы, разовый цвет — обе.
Аминокислотные остатки, принимающие участие во взаимодействии рецептора с лигандами в первой группе (рис. 13а): Glu75, Phe112, Ser160, Leu161; во второй группе: Asp71, Tyr72, Lys117, Ser158, Leu161. Со стороны молекулы лиганда АО, принимающие участие во взаимодействии для первой группы (рис. 13б): Tyr72, Glu75, Phe112, Ser159, Ser160; для второй группы: Tyr72, Glu75, Phe112, Lys117, Ser159, Ser160, Arg206. Выделенные АО частично совпадают с данными оценки поверхности силикатеина-α на степень гидрофобности, доступности для связывания с молекулой ортокремниевой кислоты (Murr and Daniel, 2005) и включают области: 72–80 АО и 151-159 АО. АО на поверхности рецептора и лиганда располагаются в общей области связывания (рис. 13). Такое расположение позволяет сделать вывод, что стыковка субъединиц силикатеина-α замыкается, формируя димер. При этом активные центры каждого силикатеина для первой группы располагаются перпендикулярно плоскости стыковки, а для второй группы располагаются по одну сторону.
Выявлено отсутствие перекрываний областей АЦ, центров потенциального связывания с природными субстратами (кремниевыми кислотами) и областью связывания со второй молекулой силикатеина-α (рис. 12б). В структуре пре-силикатеина-α выявлены

районы, которые могут регулировать его созревание в зависимости от условий, при чем кремниевые кислоты потенциально могут участвовать в этом процессе.
ЗАКЛЮЧЕНИЕ
Рассмотрены 8 стартовых конфигураций каталитической триады силикатеина-α,
отличающихся степенью протонированности имидазольного кольца His и поворотом его относительно связи Сα–Сβ. Полученные комплексы каталитическая триада (His163, Ser25, Gln19) – субстрат (Si(OH)4) отличаются энергией связывания, а также количеством и длинной водородных связей. В целом протонированная форма азота имидазольного кольца His имеет больший (по модулю) парциальный заряд. При равном протонировании атомов πN и τN, πN является более электроотрицательным атомом. Данное соотношение характерно как для свободного активного центра, так и в комплексе с Si(OH)4.
ППЭ для состояний 7 и 8 имеют только один минимум расположения H+, и данные конфигурации не приводят к его перемещению вдоль вектора водородной пары. Максимальными энергиями связывания обладают состояния, в которых имидазольное кольцо His полностью депротонировано (-32,82 ккал/моль состояние 1 и -26,78 ккал/моль состояние 5). Оптимизированный комплекс (состояние 1) является оптимальным для активации молекулы кремниевой кислоты, поскольку H+ переходит с О силиката на πN имидазольного кольца. Для комплекса состояния 5 имеется одна укороченная водородная связь между τN His163—H -(Si(OH)4), с расстоянием в 1,64Å. ППЭ перемещения H+ вдоль вектора этой связи имеет 2 минимума, вблизи каждого электроотрицательного атома (N и O). Энергетические барьеры прямой и обратной реакции сопоставимы и составляют 4,36 ккал/моль и 3,23 ккал/моль, соответственно. Для состояний 2 и 6 имеются 2 локальных минимума, разделенных энергетическим барьером. Состояние 6 имеет энергетический барьер перемещения H+ от Si(OH)4 на τN His163 в 7,82 ккал/моль и энергетический барьер обратной реакции в 0,3 ккал/моль. Для состояния 2 перемещение Н+ от πN His163 на Si(OH)4 может привести к инициации реакции конденсации кремниевых кислот по SN2 механизму (Brinker et al., 1989; Trinh et al., 2006). Данная ППЭ характеризуется высоким энергетическим барьером прямой реакции переноса H+ в 14,2ккал/моль и низким для обратной 1,7ккал/моль. Энергия активации для молекулы Si(OH)4 по SN2 механизму приходится на затраты 11,3 ккал/моль, в то время как энергия активации по на боковому механизму несколько выше и соответствует 17,3 ккал/Моль (Brinker et al., 1989; Trinh et al., 2006).
При сравнении полученных результатов следует учитывать, что погрешность вычисления энергии системы методом DFT в базисном наборе Х3LYP может достигать 4-5 ккал/моль (Kozuch and Martin, 2011). Поэтому энергетический барьер в 14,2 ккал/моль может быть сопоставим с энергией активации по SN2 механизму. Кроме того полностью протонированные и депротонированные состояния азотов имидазольного кольца His для водных растворов с нейтральным значением pH маловероятны (Zhang et al., 2011). Это приводит к снижению вероятности инициации активации кремниевой кислоты АЦ в следующих состояниях: 1, 4, 5 и 8. С учетом того, что при pH = 6,6 наступает динамическое равновесие между таутомерными формами имидазольного кольца His, где 80% приходится на протонирование πN и 20% на — τN, и то, что ППЭ перемещения H+ для состояний 2 и 6 имеют два энергетических минимума, позволяют предположить, что именно протонирование πN имидазольного кольца является предпочтительным для инициации
реакции активации молекулы субстрата. В то время как вращения имидазольного кольца His относительно связи Сα –Сβ, приводит к увеличению барьеров прямой и обратной реакции примерно в 2 раза.
Для проведения поиска областей связывания силикатеина-α с кремниевыми кислотами был сформирован набор наиболее высокочастотных и некоторых интересных форм олигомеров кремниевых кислот. Для этого была сформирована система, состоящая из равноудаленных молекул Si(OH)4, растворенная в воде (TIP3P), и реализовано силовое поле FG. Исследование зависимости олигомерного разнообразия и долей олигомеров от концентрации Si(OH)4 в растворе для двух систем, отличающихся концентрацией мономеров показаны общие тенденции доминирования линейных и разветвленных форм олигомеров исходства распределения высокочастотных олигомеров. Для обеих систем на начальном этапе моделирования (менее 0,5 нс) характерен рост количества олигомерных структур с содержанием атомов (Sin=3−5). В диапазоне 0,5–3нс наблюдается медленный рост количества олигомеров с более длинной цепочкой, а именно (Sin=6−8). Установлено, что на протяжении всего эксперимента сохраняется высокое значение доли димеров. За весь период МД до ∼98% атомов Si имели степень координации IV. Наличие в изученных системах атомов Si со степенью координации V может говорить об образовании промежуточных комплексов в процессе полимеризации кремниевых кислот.
В рамках исследования определены траектории формирования циклического тримера, бициклического тетрамера и олигомера типа Q12Q22Q41, для которых получены поэтапные схемы образования и предложен механизм их образования, подразумевающий сценарий их формирования из димерных структур. Отметим, что структура и представленность олигомера типа Q12Q22Q41 получены впервые в рамках МД-моделирования и на данный момент такая структура не идентифицирована при помощи экспериментальных методов.
Полученные результаты позволяют заключить, что реакционное силовое поле FG на рассмотренном временном интервале достаточно корректно описывает геометрию и процесс олигомеризации с точки зрения структуры олигомеров кремниевых кислот и их содержания в водном растворе. Данное силовое поле может быть использовано для дальнейшего изучения олигомерного состава и молекулярного механизма образования кремниевых структур.
Избранные олигомеры были оптимизированы и использованы для проведения молекулярного докинга с целью поиска областей связывания силикатеина-α с силикатами. Лиганды с длиной цепи Sin1-6 с силикатеином-α имеют три региона связывания вблизи АЦ: рядом с Asp6, Ser25 и Thr94 с максимальной аффинностью в -16,1 ккал/моль для пре- силикатеина-а и -10,96 ккал/моль. Максимальной аффинностью к силикатеину среди лигандов обладают линейные и разветвленные формы. Полярные и положительно заряженные АО пропептидного фрагмента (позиции 1..116) участвуют во взаимодействии с лигандами во всех 3-х выявленных областях. При докировании малых лигандов Sin=1-4 в структуру пре-силикатеина-α область связывания вблизи Ser131 разделяется на две группы, которые объединяются в одну область с увеличением размера лиганда. Лиганды Q25, Q23Q32 и Q12Q24 имеют максимальную аффинность к региону вблизи Thr329, циклические изомеры Sin=3-5 имеют максимальную аффинность к району Asp6. При образовании димеров силикатеина-α области стыковки не перекрываются с регионами адгезии кремниевых кислот.
Таким образом, проведенная оценка разнообразия кремниевых кислот в водном растворе определила спектр потенциальных природных субстратов силикатеина-α. Изучение
биофизических особенностей конфигурации аминокислотных остатков активного центра силикатеина-α и областей его взаимодействия с силикатами позволило предложить наиболее вероятные пути активации кремниевых кислот на начальном этапе биосилификации.
ВЫВОДЫ
1) Предпочтительными конфигурациями каталитической триады силикатеина-α, которые
могут инициировать реакцию конденсации Si(OH)4, являются частично протонированные
формы имидазольного кольца His по πN (состояния 2 и 6).
2) Реакционное силовое поле Feuston и Garofalini позволяет моделировать разнообразие
природных лигандов силикатеина-α, при этом длина ковалентной связи Si-O, валентные углы Si-O-Si и O-Si-O сопоставимы с экспериментальными данными; в ~98% случаев Si в составе лигандов имеет VI координационное число, что поддерживает тетраэдрическую форму силикатов; время жизни силанольных связей составило ~110пс, а силоксановых связей ~25пс.
3)Полученное графовое изомерное разнообразие кремниевых кислот со степенью олигомеризации Sin = 3-5 полностью согласуется с экспериментальными данным, кроме этого идентифицированы не выявленные в экспериментах структуры: для Sin = 5 – 1 структура, для Sin = 6 – 23 структуры, для Sin = 7 – 56 структур, среди Sin = 8 – 108 структур.
4) Кремниевые кислоты с длиной цепи Sin1-6 с силикатеином-α имеют три региона связывания вблизи АЦ: рядом с Asp6, Ser25 и Thr94 с максимальной аффинностью в -16,1 ккал/моль для пре-силикатеина-а и -10,96 ккал/моль для силикатеина-α. Максимальной аффинностью к силикатеину среди изомеров лигандов обладают линейные и разветвленные формы.