Cribado virtual y dinámica molecular en supercomputadoras para la identificación de inhibidores de la metalo-β-lactamasa NDM-1

1. Introducción

La resistencia antimicrobiana es una de las diez principales amenazas para la salud pública global [1]. Entre los mecanismos de resistencia más preocupantes se encuentra la producción de metalo-β-lactamasas (MBL), enzimas dependientes de zinc capaces de hidrolizar carbapenémicos, considerados antibióticos de último recurso. La variante NDM-1 (New Delhi metallo-β-lactamase), descrita en 2009, se ha diseminado rápidamente en enterobacterias y se ha reportado en brotes hospitalarios en Perú, Brasil, Colombia y México [2,3].

A diferencia de las serín-β-lactamasas, para las cuales existen inhibidores clínicos como el avibactam, no se dispone de inhibidores aprobados de MBL. El sitio activo de NDM-1 contiene dos iones de zinc coordinados por His120, His122, His189, Asp124, Cys208 y His250, y un bucle flexible (L3, residuos 65–73) que modula el acceso de los sustratos [4]. Esta flexibilidad complica el diseño racional de inhibidores y hace necesario el uso de simulaciones que capturen la dinámica conformacional de la enzima.

El cribado virtual de grandes bibliotecas químicas seguido de dinámica molecular (DM) es una estrategia establecida para priorizar candidatos antes de la validación experimental [5,6]. Sin embargo, las simulaciones de DM con una duración suficiente para capturar reorganizaciones del bucle L3 demandan recursos computacionales que exceden a una estación de trabajo convencional. En este estudio aprovechamos un clúster de supercomputación para ejecutar un flujo integrado de docking con funciones de puntuación basadas en redes neuronales convolucionales (GNINA) y simulaciones de DM de 200 ns con GROMACS, con el objetivo de identificar compuestos naturales con potencial inhibitorio de NDM-1.

2. Materiales y métodos

2.1. Preparación del receptor

Se utilizó la estructura cristalográfica de NDM-1 en complejo con hidrolizado de meropenem (PDB 4EYL, resolución 1,05 Å) [4]. Los residuos faltantes del bucle L3 se modelaron con MODELLER 10.4 y se seleccionó el modelo con menor puntuación DOPE. Se añadieron hidrógenos a pH 7,4 con PDB2PQR y se conservaron los dos iones de zinc del sitio activo.

2.2. Biblioteca de compuestos y docking

Se descargaron 42 318 compuestos de origen natural de la base de datos COCONUT [7], filtrados por reglas de Lipinski y ausencia de grupos PAINS. Las conformaciones 3D se generaron con RDKit y se minimizaron con el campo de fuerza MMFF94. El docking se realizó con GNINA 1.1 [8], definiendo la caja de búsqueda alrededor de los iones de zinc (20 × 20 × 20 Å), exhaustiveness = 16 y 9 poses por ligando. Los compuestos se clasificaron por el promedio ponderado de CNNaffinity y Vina affinity. El protocolo se validó mediante redocking del ligando cocristalizado (RMSD = 0,87 Å).

2.3. Dinámica molecular en HPC

Los 12 complejos mejor clasificados se simularon con GROMACS 2024.2 [9] empleando el campo de fuerza CHARMM36m para la proteína y CGenFF para los ligandos. Cada sistema se solvató en una caja dodecaédrica con agua TIP3P y se neutralizó con NaCl 0,15 M. Tras la minimización (steepest descent, 50 000 pasos) y la equilibración NVT y NPT (100 ps cada una, con restricciones posicionales), se ejecutaron 200 ns de producción a 310 K y 1 bar (termostato V-rescale, barostato Parrinello-Rahman, paso de 2 fs). Las simulaciones se ejecutaron como jobs SLURM en nodos con 2 × AMD EPYC 7763 y una GPU NVIDIA A100 por job; se usaron checkpoints cada 15 minutos para reanudar los cálculos divididos por el límite de tiempo de la cola. El consumo de recursos se detalla en la Tabla 1.

2.4. Análisis de trayectorias

Las trayectorias se procesaron con gmx trjconv para corregir las condiciones periódicas de contorno. Se calcularon RMSD, RMSF, radio de giro, puentes de hidrógeno (gmx hbond) y distancias ligando–zinc. Los paisajes de energía libre (FEL) se construyeron a partir de las dos primeras componentes principales del movimiento del bucle L3. Las energías de unión se estimaron mediante MM/GBSA con gmx_MMPBSA sobre 2 000 fotogramas de los últimos 100 ns. Los gráficos se generaron en R 4.4 con ggplot2.

3. Resultados

3.1. Cribado virtual

El docking de los 42 318 compuestos se completó en 38 core-hours. La distribución de puntuaciones CNNaffinity mostró una cola de 214 compuestos con valores superiores a 7,5 (pK_d predicho). Tras la inspección visual de las poses y el filtrado por interacción con al menos un ion de zinc, se seleccionaron 12 candidatos para dinámica molecular. Sus puntuaciones se muestran en la Tabla 2 y las poses representativas en la Figura 2.

3.2. Estabilidad de los complejos

Siete de los doce complejos mantuvieron un RMSD del ligando inferior a 2,5 Å durante los 200 ns, mientras que en cinco casos el ligando abandonó el sitio activo antes de los 80 ns (Figura 3). Los valores de RMSF revelaron que el bucle L3 fue la región más flexible en la enzima libre (RMSF máximo 4,8 Å) y que esta flexibilidad se redujo significativamente en presencia de los compuestos NP-0871, NP-1134 y NP-2290 (RMSF máximo 2,1–2,6 Å), lo que sugiere un efecto de estabilización del bucle por parte del ligando.

3.3. Interacciones y energía de unión

Los cinco compuestos con mejor desempeño conservaron la coordinación con al menos uno de los iones de zinc durante más del 85 % de la trayectoria, además de puentes de hidrógeno persistentes con Asp124 y Lys211 y una interacción catión-π con His250. Las energías de unión MM/GBSA variaron entre −38,2 y −47,6 kcal/mol (Tabla 3). El paisaje de energía libre del complejo con NP-0871 (Figura 4) mostró un único mínimo correspondiente a una conformación cerrada del bucle L3, en contraste con los dos mínimos observados en la enzima libre.

4. Discusión

El presente trabajo demuestra la viabilidad de un flujo integrado de cribado virtual y dinámica molecular ejecutado íntegramente en un entorno de supercomputación, con un consumo total de 612 core-hours. Este nivel de recursos es accesible para grupos de investigación de la región a través de alianzas con proveedores de HPC, lo que amplía las posibilidades de investigación en diseño de fármacos en países con infraestructura computacional limitada.

Los cinco compuestos identificados comparten un patrón estructural común: un grupo carboxilato o tiol capaz de coordinar los iones de zinc, junto con un sistema aromático que interactúa con el bucle L3. Este farmacóforo es consistente con inhibidores de MBL descritos previamente, como el ácido dipicolínico y el captopril [10,11], pero con una mayor complejidad estructural que podría conferir selectividad frente a metaloenzimas humanas.

La estabilización del bucle L3 en conformación cerrada, observada en el FEL del complejo con NP-0871, representa un mecanismo de inhibición atractivo, ya que impide el acceso del sustrato al sitio activo. No obstante, los resultados computacionales deben interpretarse con cautela: la precisión de las energías MM/GBSA es limitada y no se consideraron efectos de solvatación explícita en el cálculo de la energía de unión. La validación mediante ensayos enzimáticos y de sinergia con meropenem constituye el siguiente paso natural de este trabajo.

5. Conclusiones

El cribado virtual de 42 318 compuestos naturales con GNINA, combinado con dinámica molecular de 200 ns en GROMACS, permitió identificar cinco candidatos a inhibidores de NDM-1 con interacciones estables en el sitio activo y energías de unión favorables. El flujo computacional desarrollado es reproducible en clústeres académicos y sirve como referencia metodológica para estudios de diseño de fármacos asistido por supercomputación en América Latina.

Agradecimientos

Financiamiento

Conflicto de intereses

Referencias

1. World Health Organization. Global antimicrobial resistance and use surveillance system (GLASS) report 2022. Ginebra: WHO; 2022.
2. Yong D, Toleman MA, Giske CG, et al. Characterization of a new metallo-β-lactamase gene, bla_NDM-1, and a novel erythromycin esterase gene carried on a unique genetic structure in Klebsiella pneumoniae sequence type 14 from India. Antimicrob Agents Chemother. 2009;53(12):5046–5054.
3. Resurrección-Delgado C, Montenegro-Idrogo JJ, Chiappe-Gonzalez A, et al. Klebsiella pneumoniae New Delhi metallo-β-lactamase in a Peruvian national hospital. Rev Peru Med Exp Salud Publica. 2017;34(2):261–267.
4. King DT, Worrall LJ, Gruninger R, Strynadka NCJ. New Delhi metallo-β-lactamase: structural insights into β-lactam recognition and inhibition. J Am Chem Soc. 2012;134(28):11362–11365.
5. Lionta E, Spyrou G, Vassilatis DK, Cournia Z. Structure-based virtual screening for drug discovery: principles, applications and recent advances. Curr Top Med Chem. 2014;14(16):1923–1938.
6. De Vivo M, Masetti M, Bottegoni G, Cavalli A. Role of molecular dynamics and related methods in drug discovery. J Med Chem. 2016;59(9):4035–4061.
7. Sorokina M, Merseburger P, Rajan K, Yirik MA, Steinbeck C. COCONUT online: collection of open natural products database. J Cheminform. 2021;13:2.
8. McNutt AT, Francoeur P, Aggarwal R, et al. GNINA 1.0: molecular docking with deep learning. J Cheminform. 2021;13:43.
9. Abraham MJ, Murtola T, Schulz R, et al. GROMACS: high performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX. 2015;1–2:19–25.
10. Chen AY, Thomas PW, Stewart AC, et al. Dipicolinic acid derivatives as inhibitors of New Delhi metallo-β-lactamase-1. J Med Chem. 2017;60(17):7267–7283.
11. Brem J, van Berkel SS, Zollman D, et al. Structural basis of metallo-β-lactamase inhibition by captopril stereoisomers. Antimicrob Agents Chemother. 2016;60(1):142–150.

Material suplementario