Introducción
Hoy atravesamos tiempos ciertamente difíciles e impredecibles. La situación causada por una única entidad viral, el nuevo coronavirus denominado SARS-CoV-2, ha desencadenado más de 60 millones de infecciones y ha cobrado ya la vida de casi 1 millón y medio de personas a lo largo y ancho de nuestro planeta (https://covid19.who.int/). Este nuevo virus (el séptimo de su tipo en afectar a los humanos) surgió por primera vez en diciembre de 2019 entre un grupo de personas en la ciudad de Wuhan, capital de la provincia de Hubei en China y, luego de su eventual confirmación [1], la Organización Mundial de la Salud (OMS) anunció en febrero de 2020 que la enfermedad causada recibiría el nombre de COVID-19 (Coronavirus Disease 2019). [2] Finalmente, dicha enfermedad fue declarada como pandemia en el mes de marzo [3]. En el contexto particular de Colombia, el COVID-19 está teniendo un impacto sumamente grave, ubicándolo en décimo lugar entre los países con mayor nivel de contagios y letalidad, sumando a la fecha un total de 1,27 millones de nuevo casos y 35.860 muertes. Los departamentos con mayor cantidad de casos de COVID-19 reportados por el Instituto Nacional de Salud (INS) son Bogotá D.C., Antioquia y Valle del Cauca (https://coronaviruscolombia.gov.co/).
La coyuntura de salud pública generada por el COVID-19 ha puesto en jaque a la industria farmacéutica mundial [4] y, por ende, ha instado una profunda urgencia a investigadores de todo el mundo por hallar no solo una vacuna eficaz, sino también un inhibidor contra el virus que logre contener su expansión [5]. Esto viene empujando las fronteras de la informática, la química y la biología, especialmente en el área denominada Diseño de Fármacos Asistido por Computadora (DiFAC). El aumento en el número de trabajos de investigación publicados en revistas indexadas de gran renombre que hacen uso de uno o más métodos de esta disciplina es cada vez más marcado [6]. Entre los métodos derivados de estas tecnologías se destacan algunos como el screening virtual (VS), la dinámica molecular (MD), el diseño por fragmentos (FBDD) y el uso de productos naturales. Este conjunto de técnicas ha demostrado ser exitoso para acelerar la producción de nuevas moléculas, reduciendo los costos relacionados con síntesis/ensayos biológicos y haciendo uso de librerías químicas diversas, novedosas y con un potencial de patentabilidad competitivo que puedan ser llevadas más rápidamente a escalamiento industrial. Con el objetivo de identificar un inhibidor efectivo contra el virus, a pesar del corto tiempo transcurrido, se han realizado importantes avances científicos en investigación aprovechables por el DiFAC, entre los cuales sobresalen dos: la obtención experimental de la estructura de la proteasa principal (Mpro), recientemente liberadas a las bases de datos de acceso mundial, y los estudios que vinculan al producto natural licorina con la inhibición del SARS-CoV-2. La Mpro es una enzima que juega un papel esencial en la replicación y transcripción viral de este virus [7] y, por lo tanto, representa uno de los blancos terapéuticos más prometedores para el diseño de inhibidores (Figura 1). Por su parte, la licorina ha sido vinculada como un alcaloide activo contra SARS-CoV-2 al bloquear directamente su replicación, pero no se ha establecido cuál es el blanco biológico primario sobre el cual interactúa. Sin embargo, y de forma muy interesante, han surgido reportes recientes de otros inhibidores de tipo alcaloide (como los lignoides savinina e hinoquinina) estructuralmente muy similares a la licorina que inhiben la Mpro de otros virus muy relacionados como SARS-CoV-1 . En este estudio, nuestra hipótesis central se basa en que la licorina y compuestos de la misma clase podrían inhibir directamente la Mpro de SARS-CoV-2, aplicando una estrategia combinada in silico con distintas técnicas que puedan aprovechar la información de estructural de los ligandos a evaluar y el enzima blanco Mpro del nuevo coronavirus.
Estrategia metodológica
Generación de librerías extendidas
Las estructuras químicas similares a los fármacos ya conocidos en el mercado (moléculas tipo-fármaco) han sido ampliamente estudiadas y sobreexplotadas durante las últimas décadas, generando una carrera por la búsqueda de nuevos tipos de compuestos [8]. Buscando ampliar la exploración del llamado “espacio químico”, la industria farmacéutica y la academia han comenzado a orientar sus estudios hacia un grupo de compuestos con estructuras bioactivas novedosas de gran potencial y diversidad química: los productos naturales. Estos tienden a poseer una selectividad biológica mayor a objetivos celulares [9] y se pueden encontrar en organismos tan diversos como plantas, invertebrados marinos, bacterias, hongos y otros microorganismos [10]. Gracias a la contribución de la Universidad ICESI, con sede en Cali, Colombia, la licorina (asociadas a la inhibición del SARS-CoV-2) y unas pocas moléculas derivadas se usaron como punto de partida para la generación de una librería extendida de moléculas desde las cuales se pueda llegar a identificar uno o más inhibidores candidatos. Esto se realizó en dos fases distintas. Inicialmente, se obtuvo la representación de las moléculas en formato digital y, a través del servidor web ChemMine Tools (https://chemminetools.ucr.edu/), se llevó a cabo una búsqueda por similitud estructural en la base de datos PubChem usando el coeficiente de Jaccard-Tanimoto, el cual es comúnmente empleado para la comparación y evaluación de diversidad dentro de un conjunto de compuestos. Dicho índice se califica en una escala de 0-1, donde la unidad representa una completa intersección entre los pares de elementos entre dos moléculas comparadas. Finalmente, esta nueva librería generada, dentro de un valor límite de 0.90, se usó como estructura central o “core” para el anclaje de diversos grupos R (Figura 2), obtenidos a partir de un catálogo disponible en el programa Maestro en la Suite Schrödinger v2020-1 [11].
Docking molecular y cálculos de afinidad
Una vez generada la librería extendida de licorina, se buscó entender cuál o cuáles de estos derivados podrían llegar a ser más promisorios llevando a cabo una simulación que permitiese medir la mejor pose y afinidad de estos compuestos en el sitio activo de la enzima de estudio Mpro. Para ello, se realizó un protocolo de acoplamiento o docking molecular de dicha librería usando el programa Glide en la plataforma Schrödinger v2020-1 [11]. La geometría que adopta un ligando cuando interactúa con un receptor es conocida como pose. En el docking molecular se usa un algoritmo que genera varias poses las cuales son posteriormente ubicadas en un ranking que se obtiene usando una función de puntaje [12], [13]. Las poses con mejores puntajes son consideradas como poses que se aproximan o son iguales a la pose experimental.
Con el fin de incluir algunos residuos como flexibles en el proceso, estudiar varias estructuras 3D de la Mpro en complejo con inhibidores ya conocidos. En adición, se implementaron cuatro criterios distintos, basados en el trabajo de Anderson y colaboradores [14]: 1) su presencia en un región de loop o bucle, 2) contactos de van der Waals (vdW) a ligando, 3) posibilidad de movimiento en el sitio activo y 4) evidencia de otros rotámeros conformacionales entre distintos cristales reportados. Posteriormente, se procedió a un proceso híbrido posdocking de refinamiento en la estimación de afinidad de cada compuesto empleando un método híbrido denominado MM/GBSA (Mecánica molecular/Área Superficial Generalizada de Born) [15]. Este tipo de cálculos son similares a los cálculos MMPBSA con la diferencia de que la energía polar de solvatación es calculada con otro tipo de acercamiento [16]. Después de obtener las poses con mejor puntaje de MMGBSA se hizo una inspección visual de las poses obtenidas.
Simulaciones de dinámica molecular y modelamiento de farmacóforos
Para poder justificar el uso de simulaciones de dinámica molecular para el diseño de farmacóforos, es necesario tener indicios de la existencia de una alta flexibilidad del sitio activo de la proteasa incluyendo no solo rotámeros de su cadena lateral (como en el docking), sino también en su backbone o estructura central.
Con este fin, inicialmente se hicieron cálculos de volúmenes y de áreas de cavidades de varios cristales de proteasas en estado libre (apo) y en complejo (holo) utilizando el servidor web CASTp [17], [18], en donde las cavidades de volumen mayor se consideraron como el sitio activo [19]. Adicionalmente, se tomaron 20 inhibidores con actividad inhibitoria conocida y se usaron en simulaciones de dinámica molecular (MD) usando el programa Gromacs v.2020, la cual permitió visualizar el movimiento de los átomos y moléculas de nuestros sistemas proteína-ligando a través del tiempo [20], [21]. Por medio de un estudio de clustering [22], se obtuvieron las poses más representativas de estos compuestos en el sitio activo. Posteriormente, estas poses se usaron para validar el modelo de farmacóforos usando el software comercial LigandScout v4.4 [23], [24].
Estos modelos son abstracciones geométricas que representan el ensamble de características fisicoquímicas compartidas entre compuestos activos de una proteína y que favorecen directamente su interacción con la misma [25]. Si bien pueden estar basados en ligando o estructura de proteína, la estrategia se implementó en este trabajo consiste en usar la información de complejos proteína-ligando, la cual es mucho más robusta y completa [26].
Resultados preliminares
La construcción de las librerías extendidas de licorina y los derivados permitió obtener un conglomerado de un poco más de 280 mil compuestos químicos. En la primera parte, orientada a la búsqueda por similitud, se decidió realizar búsquedas de similitud para un barrido sobre el coeficiente de Jaccard-Tanimoto desde 0.99 hasta 0.90 con el fin de seleccionar librerías suficientemente grandes sin distanciarse demasiado de la molécula inicial. Al finalizar las búsquedas, se exportó la lista de compuestos resultantes en formato SDF e, intuitivamente, se evidenció una tendencia creciente con la disminución en el parámetro de similitud evaluado (Figura 2). Adicionalmente, dada la simplicidad estructural de la licorina y la licorina hidrocloruro, se encuentra que estas aportan cerca del 50 % de la totalidad de las librerías conjuntas encontradas.
Esto sugiere que la simplicidad estructural de la anhidrolicorina podría estar asociada con una elevada cantidad de moléculas a su respectiva librería. No obstante, se observó que es aquella con menor cantidad de moléculas similares. Impedimentos cíclicos, enlaces π y la ausencia de cadenas laterales en el scaffold principal resultan en una evidente restricción a la diversidad dentro del espacio químico.
La presencia de cadenas laterales complejas y extensas (1,2-benzoillicorina) generan asimismo un efecto equivalente. Por otro lado, se estudió el efecto de la enumeración de grupos R en estos derivados iniciales. El análisis de estos compuestos en términos de características fisicoquímicas críticas en el contexto de los fármacos indica que su gran mayoría se mantienen bastante por debajo, como en peso molecular, lipofilicidad y cantidad de aceptores o donadores de puentes de hidrógeno. Lo anterior, sumado a la distribución esencialmente gaussiana que describen, sugiere un gran potencial de contribución a la diversidad y novedad al espacio químico disponible, derivado en buena parte de la librería de grupos R implementada. Dicho potencial permite expandir el campo de estudio alrededor de la licorina y se exploran nuevos aportes que posiblemente contribuyan a la mejoría en la acción inhibitoria del blanco farmacológico de interés.
Después del modelamiento molecular de la proteína Mpro (Figura 3A-B) y la obtención de los residuos a ser considerados como flexibles se procedió a realizar el docking molecular de los compuestos con actividad inhibitoria contra la proteasa de SARS CoV y la proteasa de SARS CoV 2. En este docking molecular se limitó el número de poses generadas por ligando a 20 poses. Luego de obtener el ranking de las mejores poses usando la función de puntaje XP de Glide se optó por hacer re-scoring de las poses con cálculos MMGBSA. Para ello, se seleccionaron las poses iniciales de los inhibidores en el sitio activo de la proteasa (experimentales o predichas) y con estas poses se llevaron a cabo simulaciones de dinámica molecular de la proteasa Mpro en forma apo y 56 simulaciones de dinámica molecular de la proteasa en complejo con varios inhibidores. Todas las simulaciones de dinámica molecular presentadas en este estudio poseen un tiempo de simulación de 0.5 μs, y se analizaron por medio de distintas herramientas como el análisis de componentes principales (PCA) y los cálculos de energía libre MMPBSA. Los resultados provenientes de los cálculos de energía libre MMPBSA se utilizaron para encontrar 29 inhibidores con valores de energía libre de unión proteína ligando menores a 100 kJ/mol. En algunos cristales en el subsitio S2 se podía evidenciar el surgimiento de un bolsillo que no era observado en los demás cristales. La apertura y cierre es este bolsillo de S2 es mediada por Met 49 (Figura 4C). El movimiento de Met 49 se da gracias al cambio en el ángulo del diedro formado por los carbonos beta y gamma de la cadena lateral de la metionina. Este cambio en la cadena lateral de la metionina únicamente es observado cuando un grupo hidrofóbico de un ligando se ubica en el subsitio S2. Por lo tanto, en los cristales apo este movimiento de la metionina de S2 no se observa.
Figura 3.
Protocolo de ampliación en la librería de licorina. La gráfica ilustra el número de compuestos identificados en PubChem en función del valor límite (cut-off) de similitud usado en el índice de Jaccard-Tanimoto.
Figura 4.
(A) Sitio activo y subsitios del sitio activo de la proteasa 3CLpro de SARS CoV-2. (B) Ubicación de la diada catalítica en el sitio activo de 3CLpro de SARS CoV-2 (PDB 6YB7) (C) Movimiento de Met 49 entre los cristales con el bolsillo accesorio cerrado (verde) y bolsillo accesorio abierto (morado).
Finalmente, se usaron las trayectorias de dinámica molecular de estos 20 inhibidores para extraer estructuras representativas por medio de métodos de clustering. Luego, estas estructuras se emplearon para generar 20 modelos de farmacóforos basados en proteína ligando. Estos modelos se superpusieron para poder identificar en qué lugares del sitio activo se encuentran presentes los farmacóforos y como resultado se obtuvo la imagen reportada en la Figura 5. En esta imagen se observa que un gran número de modelos de farmacóforos poseen donadores de puente de hidrogeno en el subsitio S1. Adicionalmente, el modelo correspondiente a la trayectoria del inhibidor 7C8T posee un aceptor de puente de hidrógeno en S1. En el caso del subsitio S1’ se encontraron dos farmacóforos ubicados en la cercanía de Cys 145 y His 41. Los farmacóforos cercanos a Cys 145 tienen la característica de aceptores de puente de hidrogeno y los farmacóforos próximos a His 41 contienen la característica de grupo hidrofóbico. Con respecto a los subsitios S2 y S4 estos incluyen grupos hidrofóbicos, aceptores y donadores de puente de hidrogeno. Finalmente, en S3 se ubican varios farmacóforos con la característica de donador de puente de hidrógeno.
Figura 5.
Modelamiento de farmacóforos para inhibidores de Mpro. (A) Modelos de farmacóforos superpuestos de los 20 complejos con menor valor de energía libre. (B) Pictogramas de los farmacóforos más comunes.
A raíz de los resultados obtenidos tras la superposición de los 20 modelos de farmacóforos (Figura 5), se optó por generar dos modelos de farmacóforos que a su vez contienen las características más repetidas en la superposición: un aceptor y un donador de puente de hidrógeno en S1, un grupo hidrofóbico en S1’ cercano a His 41, un grupo hidrofóbico en S2, un donador y un aceptor de puente de hidrogeno en S3, un donador de puente de hidrogeno y un grupo hidrofóbico en S4 y, finalmente, un donador de puente de hidrogeno ubicado entre S1’, S1, S2, S3 y S4. El hecho de que uno de estos dos modelos contenga la característica de grupo electrofílico es una característica novedosa en el modelamiento por farmacóforos para la proteína de estudio. Por otra parte, la diferencia principal entre los dos modelos de farmacóforos generados es la zona correspondiente al hueco de oxianión, en la cual ocurre el ataque nucleofílico del residuo Cys 145 al sustrato. No obstante, gracias a que los modelos de farmacóforos carecen de flexibilidad puede ocurrir que algunos compuestos activos sean considerados como inactivos. Para solucionar esta carencia en un futuro se pueden crear volúmenes de exclusión flexibles que son generados a partir de trayectorias de dinámica molecular; de esta forma, los modelos de farmacóforos tendrían en cuenta la flexibilidad del sitio activo y mejoraría la capacidad del método para reconocer moléculas activas.
Conclusiones
Los resultados obtenidos a través de este estudio están trazando un camino importante y esperanzador hacia la generación de una nueva molécula con capacidad inhibitoria sobre el nuevo coronavirus SARS-CoV-2. Si bien el análisis de librerías derivadas de la licorina permitirá establecer las mejores rutas de optimizar su armazón estructural o scaffold central, también representa una forma de complementar nuevos modelos computacionales relacionados con la proteína blanco para obtener inhibidores diversos, sintetizables, desarrollables y potentes. Por otra parte, consideramos que los datos computacionales obtenidos pueden incorporarse en una eventual estrategia basada en fragmentos derivados de productos naturales (NPDFs), partiendo de estructuras distintas a la licorina: no solo los candidatos resultantes podrían usarse en el diseño racional de nuevos agentes antivirales, sino que también la estrategia metodológica podría ser aplicada otros blancos proteicos de interés terapéutico.
Finalmente, respecto a nuestra perspectiva a futuro, creemos que el paso siguiente involucrará, sin duda, el desarrollo de estudios de tipo experimental que complementen lo obtenido tras las simulaciones, enfrentando como el principal desafío la necesidad de una alianza bajo el llamado “modelo de la triple hélice” en el que se fortalezca la cooperación entre la academia, el gobierno y la industria farmacéutica. Dicha iniciativa tendrá un valioso impacto en la creación de nuevos conocimiento científicos e innovadores con un papel clave en el desarrollo de medicamentos creados en nuestro país, especialmente en esta apremiante cruzada tecnológica que permita generar confianza durante la recuperación del COVID-19.
