Computational modeling of the enzymatic activities of biomolecules at different scales: from quantum mechanical reaction studies to systemic understanding of cell behavior

The aim of the thesis is the development of computational models of the enzymatic activity of biomolecules at different scales. Parallel investigations have been carried out at a quantum level to study the reactivity of an enzyme from an electronic point of view, and at a systemic level using simulation techniques to determine the role of enzymes in the network of cellular reactions. Starting from the lowest complexity level, the thesis begins with two computational studies with the aims of describing the molecular mechanism of the catalytic reaction between amavadin and methyl mercaptoacetate and the structural coordination between angiogenin and copper ion, respectively. Since in both cases transition metals are involved, Density Functional Theory (DFT) computations and Quantum Theory of Atoms in Molecules (QTAIM) analysis of the electron density were used. The first study regards the electronic description of the enzymatic mechanism by which molecule amavadin, an unusual octa-coordinated VIV complex isolated from Amanita muscaria mushrooms, catalyzes the oxidation of some thiols to the corresponding disulfides. An experimental work by G.da Silva et al. proposed an inner-sphere mechanism for the reaction (J. Am. Chem. Soc. 1996, 118:7568–7573.) but the electronic mechanism was not identified. In the first step of our investigation, the stereoelectronic features of the V(IV) (inactive) and V(V) (active) states of amavadin were determined by means of DFT. Then, the formation of the VV complex with methyl mercaptoacetate (MMA), which has been chosen as a prototypical substrate, has been characterized both thermodynamically and kinetically. DFT results reveal that protonation of VV amavadin at a carboxylate oxygen not directly involved in the V coordination, favors MMA binding into the first coordination sphere of vanadium, by substitution of the amavadin carboxylate oxygen with that of the substrate and formation of an S–H···O hydrogen bond interaction. The latter interaction can promote SH deprotonation and binding of the thiolate group to vanadium. The kinetic and thermodynamic feasibility of the V(V)–MMA intermediates formation is in agreement, along with electrochemical experimental data, also with the biological role exerted by amavadin. Finally, the presence of an ester functional group as an essential requisite for MMA oxidation has been rationalized. Moreover, the results proposed an indirect evidence on the role of the vanadium (in its d0 active state) based catalyst protonation. The second studies focused on the structural and spectroscopic properties of the complexes formed by Cu2+ and the peptide fragments Ac-PHREN-NH2, which encompasses the putative cell binding domain of angiogenin, as well as its Ac-PHREN-NH2 variant. Analysis of structures, relative energies and EPR parameters has allowed to conclude that the metal coordination environment at pH 8 is formed by a nitrogen atom of His, two deprotonated amide groups, a water molecule and an oxygen atom from the COO- side chain of Glu, in nice agreement with recent experimental results [Dalton Trans, 2010, 39:10678]. Moreover, DFT results allowed to reveal that the Glu sidechain of the Ac-PHREN-NH2 peptide is coordinated in equatorial position, fully disclosing the effects of Cu2+ binding on the structural properties of this key angiogenin portion. We have also investigated the configuration space of the E→Q mutated system Ac-PHRQN-NH2/Cu2+·H2O In this case, computational results led to the conclusion that the H2O molecule is coordinated in equatorial position and the oxygen atom of the carbonyl group of glutamine is weakly coordinated in apical position. Using as a reference the recent experimental results reported by Bonomo and collaborators [Dalton Trans, 2010, 39:10678], we have used classical MM/MD calculations followed by DFT optimizations to explore the configurational space of the Ac-PHREN-NH2/Cu2+·H2O complex. Both of these DFT studies have been carried out in the laboratory of professor Piercarlo Fantucci, Department of Biotechnology and Biosciences, University of Milan-Bicocca. These two works lead to the following publications: “DFT characterization of key intermediates in thiols oxidation catalyzed by amavadin”, Dalton Transaction, 40 (30): 7704-7712, 2011. “Copper coordination to the putative cell binding site of angiogenin”. A DFT investigation., Inorg.Chemistry. Accepted. Acquired some modeling skills at the molecular level, we decided to add another layer of complexity to our investigation and we decided to test if a correlation could exist between the peculiar structure of a protein and a biological effect at a cellular level. The fundamental biological system of the Ras/Sos signaling activation pathway has been chosen for our study. Based on two different previous studies we developed a model that can correlate bistability and the microdomains clustering in this small network. In fact Das et al. (Cell 2009, 136:337-351, 2009) reported that the Sos-dependent Ras activation in lymphocytes, beyond stimulus, causes a digital (on/off) response that implies a positive feedback loop in the regulation mechanism of Sos activity. The deterministic model developed by the authors demonstrates that for low or high levels of Sos there is only one possible state correlated to low or high levels of active Ras, respectively. For intermediate values of Sos, three states of Ras activity are generated (two states are stable and could be simultaneously observed, the third is unstable and slightly perturbated). Bistability occurs because system lies in the lower state until this is no longer possible and then jumps to an high Ras activity correlated state. An interconnected work reproduced the same biological system in a minimal 2D lattice model in order to investigate the relevance of the microdomains clustering in the signaling pathway. It was reported that positive feedback, in the presence of slow diffusion, results in clustering of activated molecules on the plasma membrane and rapid spatial spreading as the front of the cluster of the newly formed membrane bounded Ras-GTP propagates with a constant velocity dependent on the feedback strength (J.Chem.Phys, 2009, 130:245102). Our work, implemented in the software Smoldyn (Phys.Biol, 2000, 1:137-151), presents a new description of both the above mentioned models by utilizing brownian dynamic stochastic simulations and spatial discretization. It describes how the microdomains clustering of the components of the network can induce bistability. Interestingly it has recently been found a correlation between the activation levels of components of networks that show ultrasensitivity (such as the Ras/SOS regulatory network) and the Fermi-Dirac statistics generally used to measure the probability of the distribution of the states of a system of particles (PNAS, 2010, 107:1247-1252). We are carrying out further investigation to test if this kind of description, in which the molecules of the systems can be treated as analogous of fermion gas, could be useful to identify the statistic distribution of the states of a bistable system. This study has being carried out in the laboratory of professor Piercarlo Fantucci, Department of Biotechnology and Biosciences, University of Milan-Bicocca. In the last part of the thesis we were finally able to deal with a whole cellular network. We focused on the heterotrimeric Gq/11 protein signaling network which is known to be critical for evoking varied physiological functions including sensory perceptions, behavioral and mood regulation, regulation of the immune system activity and inflammation. Dysfunctions of the Gq/11 network can lead to diseases like cancers and immune system deficiencies. A detailed network topology of the various pathways that emanate from Gq/11 and that allow signal to flow from the receptor to multiple transcription factors such as AP1, c-fos, c-jun, SRF and CREB has been developed. We have then implemented a multi compartment ordinary differential equation model using the software Virtual Cell to analyze signal flow from the plasma membrane through the cytoplasm to the nucleus. The obtained time courses of activation of key components of the system, such as small GTPases, MAPkinases and transcription factors are in accordance with experimental data provided by the S.Gutkind Lab (Oral and Pharyngeal Cancer Branch, National Institute of Dental and Craniofacial Research, National Institutes of Health, Bethesda, Maryland, USA). The obtained numerical data will soon be implemented in dynamic graphs to identify the key regulatory motifs of the information processing through the Gq/11 protein signaling network that are controlling the temporal and spatial activity of transcription factors such as SRF, c-jun, c-fos, AP1 and CREB. This study has being carried out in the laboratory of the Professor Ravi Iyengar (Department of Pharmacology and Systems Therapeutics) of the Mount Sinai School of Medicine, New York, USA. The present thesis allowed us to study four systems each at a different level of complexity, according to the biological question we wanted to answer. The models we proposed constitute another demonstration of the importance of quantitative modeling in biology. This work also suggests that even if all the approaches we used were mathematically different to each other, they share a common methodology. This implies that an underlying correlation runs through all the scales of modeling.