Parametrization of crizotinib and ROS1 kinase protein using the Martini3 force field

Autores

Vilacha, J.F. (UNIVERISITY OF GRONINGEN) ; Stevens, J. (UNIVERSITY OF GRONINGEN) ; Marrink, S. (UNIVERSITY OF GRONINGEN)

Resumo

Crizotinib is a kinase inhibitor targeting kinases c-MET and the ALK and, later, the c-ros oncogene1 kinase (ROS1). Unfortunately, it is challenging to experimentally assess the effect of mutations on the proteins' activity and their interaction with kinase inhibitors. In this context, molecular modeling has been previously used to investigate the impact of mutations in the targeted protein and further, drug binding. However, these techniques are often hampered by limited computational resources. The use of coarse-grained force fields is a known model for multiscale simulations and novel advances on the Martini FF enable drug binding studies in microsecond scales requiring moderate computational power. In this work, we provide parametrized models of crizotinib and its target, ROS1.

Palavras chaves

kinases; coarse-grained; Martini

Introdução

3.5. Introdução In lung cancer, kinases are widely explored as drug targets. Mutations in this class of protein are a classical example of oncogene addiction in Non-Small Cell Lung Cancer (NSCLC) and have been successfully inhibited by kinase inhibitors (BHULLAR et al., 2018). Within mutated kinases in NSCLC patients, the incidence of ROS1 rearrangement is proximately 2%. The rearrangement of ROS1 is often through the fusion of its kinase domain with at least 55 possible partners, the most common being the CD74-ROS1 fusion protein (DAVARE et al., 2015). Mutations on the ROS1 kinase domain often correlate with resistance to type I or type II kinase inhibitors and narrow the possible therapeutical choices for patients' treatment (REMON et al., 2021). In this scenario, there is a critical need to develop new molecules that could evolve into potent drugs to tackle the mutant kinase. In medicinal chemistry, computational tools have been widely used but often faced limitations due to the computational cost of simulating a biological system in experimental-like conditions. The computational cost is related to the expense of reproducing all the atoms associated with the biological system of interest (SALO-AHEN et al., 2021). An alternative is to perform such studies using a coarse-grained (CG) force field, which groups different atoms into a representative unit, often called a “bead”. The decrease in the “resolution” of the system accelerates molecular dynamic simulation studies allowing us to aim at bigger scales and time ranges (BARNOUD; MONTICELLI, 2015). Within the field of coarse-grained simulations, the force field Martini has been the top choice for the representation of biomembranes due to its ability to reproduce the behavior of membrane lipids and proteins(MARRINK et al., 2022). The latest updates of this force field, Martini3, improved the interbeds interactions and implemented new classes of beads which, together, lead to the ability of the force field application for ligand-target binding studies (SOUZA et al., 2020). In this scenario, Martini3 raises as a promising tool not only for drug development endeavors but also for drug repurposing studies. This new development facilitates the use of coarse-grained approaches in the study of how small molecules interact with possible targets. Unfortunately, the pool of available small molecules parametrized for the MARTINI3 force field is limited and does not contain any ROS1 inhibitor. In this work, we propose a coarse- grained model for crizotinib based on the Martini3 force field. We also present a CG model for ROS1, a known target for crizotinib in the treatment of lung cancer.

Material e métodos

The atomistic structure and molecular topology for R-crizotinib compatible with the OPLS force field were obtained by feeding the SMILES sequence, obtained from PubChem to the LigParGen server. For the ROS1 model, we used a co-crystal structure with crizotinib (PDB entry 3ZBF). We destitute the crystal structure from the ligand, ions, and water molecules. Missing motifs were modeled using the rotamer library on Modeller (JALILY HASANI; BARAKAT, 2017). For the parametrization studies, drug and target were initially separately simulated using the OPLS force field in cubic boxes with edges distanced 1.0 nm from the substrate. The all-atom simulations were carried out in the presence of TIP3P water. The minimization steps were done using the steep descent method for 5000 cycles followed by 250ps of canonical ensemble equilibration. The equilibrated systems were simulated with GROMACS (version2019.5) during 600ns. For the coarse-grained simulations, we used a combination of the latest available version of Martini3 and GROMACS (version2019.5). The parameterization of crizotinib was performed manually following multiple cycles focused on the optimizations of bounded terms and intramolecular angles. For the optimization of the receptor, we used Martinize to obtain either a MARTINI model alone (model M) or combined with Elastic Network (model M+EN) (PERIOLE et al., 2009). The models M+EN were further optimized by 1) modulations of upper and lower distance cutoff between the backbone beads to be connected through the elastic bond, a harmonical potential between unbounded beads; 2) by addition or removal of selected harmonical bonds from the EN.

Resultado e discussão

For the R crizotinib, a mapping based on the Martini bible was proposed. We combined both small and tiny beads aiming to better reproduce the volume of the all-atoms models. In the CG model, we focused on the intramolecular angles aiming to maintain the stiffness of the molecule around the ethoxy moiety, as would be expected in the all-atoms model due to the chiral center. Despite having good accordance between our model and all atoms’ simulations for the angles involving the chiral center with either the 2,6-dichloro-3-fluorophenyl or the 2-aminopyridine rings, our CG model was not able to reproduce de bimodality of the angle between the aminopyridine ring and the pyrazole ring. This bimodality is also observed when we compare the crystallographic structure of different targets with R- crizotinib. Despite having a perfect overlap in the chiral center and substituents on all three targets, ALK (PDB entry 2XP2), cMET (PDB entry 2WGJ), and ROS1 (PDB entry 3ZBF), we could identify two major conformers. The two conformers differ by the rotation of the pyrazole ring, which can be described by the position of the nitrogen at position 2 of the pyrazole ring according to the aromatic nitrogen from the aminopyridine core. Our model was not capable of reproducing the bimodal distribution of this angle and its corresponding bonds. Another remarkable assessment of our model when compared with the all-atom simulations is the solvent- accessible surface area (SASA). Despite the average values being comparable (CG: 7.2 nm2 AA: 7.5 nm2), the distribution throughout the simulations presents a different behavior. While we have a narrow peak for the all- atom simulations, the CG simulations provide a broader distribution. We hypothesize that such broad distribution occurs due to a combination of the nature of the mapping for the chiral center, as observed in the Conolly Surface, and the loser bonds and weaker that are limited by the simulation stability For the parametrization of the receptor, we used Martinize, a tool capable of automatizing the process of protein parametrization. Martinize allows the user to use the MARTINI3 force field to model the desired protein using predetermined beads and backbone interactions. However, it is well described that MARTINI alone struggles to maintain the tertiary structure of proteins. This artifact can be represented by high values of Root Mean Square Deviation (RMSD) or Radius of Gyration (RoG) when compared with the starting coarse- grained model. To overcome such a bottleneck, we combined the Martini force field with additional harmonical interaction between non-bonded backbone beads. This interaction can be placed based on a cut- off distance between beads. For this study, we settle our initial elastic bond force constant to 500 kJ.mol-1.nm-2 and the lower and upper elastic bond cutoff to 0.5 and 0.9 nm respectively. Our initial simulations showed that, as expected, the MARTINI3 fails to properly reproduce the behavior observed for the ROS1 in the all-atoms simulations. But we observed a much more promising profile for the combined MARTINI3+EN model, with RMSD values being closely related to the all-atoms simulations. However, analysis of the Root Mean Square Fluctuation (RMSF) shows different behavior depending on the region of the protein. As observed, the combination of MARTINI3+EN can lead to over and underestimation of the thermal fluctuation, simultaneously. As the long-term goal of this project is to perform ligand-binding studies, we focused on, primarily, optimizing the regions involved in drug binding. For such, we went back to the co-crystal ROS1-crizotinib and obtained the fingerprint of the interaction. Based on this interactive map, we observed that the G-loop and the regulatory aC-helix in our CG model were presenting a more restrained motion than the same region in the AA simulation. To address this limitation, we omitted the EN interactions connecting beads of both regions. Despite presenting a conspicuous improvement in the overall RMSF of these motifs, we can infer that further improvements are in order.

Conclusões

In this work we propose a model for a ligand (R crizotinib) and its validated target (ROS1) using the MARTINI 3 force field. Our efforts indicate that crizotinib is a hard-to-parametrize molecule, due to the combination of a still end, with the presence of the chiral center, and a motif with rotational freedom. While the parametrization of crizotinib relied mostly on manual labor to determine the intramolecular terms, we could rely on an automatized approach for the parametrization of the receptor. But even this automatized approach still needed human intervention to manually determine which Elastic bonds should be conserved and which should be released to better represent the protein plasticity. This work settles the based for further optimization of such a crucial small molecule. The potency of crizotinib against three distinctive kinases makes us wonder which others could be sensitive to this chiral molecule. Our interest in parametrizing the kinase domain of ROS1 is also justifiable; this kinase has been linked with lung cancer and the raise of mutations is quickly raising as a worrying point for cancer treatment. Thus, our protocols can also be applied to obtain CG models of the mutated ROS1 Furthermore, this works sets the grounds for binding studies using either the ligand or the receptor here described to develop novel molecular entities or even repurposing studies of “old” drugs, already approved by the responsible agencies

Agradecimentos

The authors would like to thank the Rijkuniversiteit Groningen for the grant for J.F Vilacha Ph.D. and the European Research Council for the ERC advance grant for prof Marrink.

Referências

BARNOUD, J.; MONTICELLI, L. Coarse-Grained Force Fields for Molecular Simulations. In: [s.l: s.n.]p. 125–149.
BHULLAR, K. S. et al. Kinase-targeted cancer therapies: Progress, challenges and future directions. Molecular Cancer, v. 17, n. 1, p. 1–20, 2018.
DAVARE, M. A. et al. Structural insight into selectivity and resistance profiles of ROS1 tyrosine kinase inhibitors. Proceedings of the National Academy of Sciences of the United States of America, v. 112, n. 39, p. E5381–E5390, 2015.
JALILY HASANI, H.; BARAKAT, K. Homology modeling: An overview of fundamentals and tools. International Review on Modelling and Simulations, v. 10, n. 2, p. 129–145, 2017.
MARRINK, S. J. et al. Two decades of Martini: Better beads, broader scope. Wiley Interdisciplinary Reviews: Computational Molecular Science, n. April, p. 1–42, 2022.
PERIOLE, X. et al. Combining an elastic network with a coarse-grained molecular force field: Structure, dynamics, and intermolecular recognition. Journal of Chemical Theory and Computation, v. 5, n. 9, p. 2531–2543, 2009.
REMON, J. et al. Current treatment and future challenges in ROS1- and ALK-rearranged advanced non-small cell lung cancer. Cancer Treatment Reviews, v. 95, n. February, p. 102178, 2021. Disponível em: <https://doi.org/10.1016/j.ctrv.2021.102178>.
SALO-AHEN, O. M. H. et al. Molecular Dynamics Simulations in Drug Discoveryand Pharmaceutical Development. p. 1–60, 2021.
SOUZA, P. C. T. et al. Protein–ligand binding with the coarse-grained Martini model. Nature Communications, v. 11, n. 1, p. 1–11, 2020.