PTM-Psi: A python package to facilitate the computational investigation of post-translational modification on protein structures and their impacts on dynamics and functions.

Post-translational modification (PTM) of a protein occurs after it has been synthesized from its genetic template, and involves chemical modifications of the protein's specific amino acid residues. Despite of the central role played by PTM in regulating molecular interactions, particularly those driven by reversible redox reactions, it remains challenging to interpret PTMs in terms of protein dynamics and function because there are numerous combinatorially enormous means for modifying amino acids in response to changes in the protein environment. In this study, we provide a workflow that allows users to interpret how perturbations caused by PTMs affect a protein's properties, dynamics, and interactions with its binding partners based on inferred or experimentally determined protein structure. This Python-based workflow, called PTM-Psi, integrates several established open-source software packages, thereby enabling the user to infer protein structure from sequence, develop force fields for non-standard amino acids using quantum mechanics, calculate free energy perturbations through molecular dynamics simulations, and score the bound complexes via docking algorithms. Using the S-nitrosylation of several cysteines on the GAP2 protein as an example, we demonstrated the utility of PTM-Psi for interpreting sequence-structure-function relationships derived from thiol redox proteomics data. We demonstrate that the S-nitrosylated cysteine that is exposed to the solvent indirectly affects the catalytic reaction of another buried cysteine over a distance in GAP2 protein through the movement of the two ligands. Our workflow tracks the PTMs on residues that are responsive to changes in the redox environment and lays the foundation for the automation of molecular and systems biology modeling.

NWChem: Recent and Ongoing Developments.

This paper summarizes developments in the NWChem computational chemistry suite since the last major release (NWChem 7.0.0). Specifically, we focus on functionality, along with input blocks, that is accessible in the current stable release (NWChem 7.2.0) and in the "master" development branch, interfaces to quantum computing simulators, interfaces to external libraries, the NWChem github repository, and containerization of NWChem executable images. Some ongoing developments that will be available in the near future are also discussed.

Properties of Gaseous Deprotonated L-Cysteine S-Sulfate Anion [cysS-SO3]-: Intramolecular H-Bond Network, Electron Affinity, Chemically Active Site, and Vibrational Fingerprints.

L-cysteine S-sulfate, Cys-SSO3H, and their derivatives play essential roles in biological chemistry and pharmaceutical synthesis, yet their intrinsic molecular properties have not been studied to date. In this contribution, the deprotonated anion [cysS-SO3]- was introduced in the gas phase by electrospray and characterized by size-selected, cryogenic, negative ion photoelectron spectroscopy. The electron affinity of the [cysS-SO3]• radical was determined to be 4.95 ± 0.10 eV. In combination with theoretical calculations, it was found that the most stable structure of [cysS-SO3]- (S1) is stabilized via three intramolecular hydrogen bonds (HBs); i.e., one O-H⋯⋯N between the -COOH and -NH2 groups, and two N-H⋯⋯O HBs between -NH2 and -SO3, in which the amino group serves as both HB acceptor and donor. In addition, a nearly iso-energetic conformer (S2) with the formation of an O-H⋯⋯N-H⋯⋯O-S chain-type binding motif competes with S1 in the source. The most reactive site of the molecule susceptible for electrophilic attacks is the linkage S atom. Theoretically predicted infrared spectra indicate that O-H and N-H stretching modes are the fingerprint region (2800 to 3600 cm-1) to distinguish different isomers. The obtained information lays out a foundation to better understand the transformation and structure-reactivity correlation of Cys-SSO3H in biologic settings.

Molecular-based analysis of nanoparticle solvation: Classical density functional approach.

Proper statistical mechanics understanding of nanoparticle solvation processes requires an accurate description of the molecular structure of the solvent. Achieving this goal with standard molecular dynamics (MD) simulation methods is challenging due to large length scales. An alternative approach to this problem can be formulated using classical density functional theory (cDFT), where a full configurational description of the positions of all the atoms is replaced by collective atomic site densities in the molecule. Using an example of the negatively charged silica-like system in an aqueous polar environment represented by a two-site water model, we demonstrate that cDFT can reproduce MD data at a fraction of the computational cost. An important implication of this result is the ability to understand how the solvent molecular features may affect the system's properties at the macroscopic scale. A concrete example highlighted in this work is the analysis of nanoparticle interactions with sizes of up to 100 nm in diameter.

Site Density Functional Theory and Structural Bioinformatics Analysis of the SARS-CoV Spike Protein and hACE2 Complex.

The entry of the SARS-CoV-2, a causative agent of COVID-19, into human host cells is mediated by the SARS-CoV-2 spike (S) glycoprotein, which critically depends on the formation of complexes involving the spike protein receptor-binding domain (RBD) and the human cellular membrane receptor angiotensin-converting enzyme 2 (hACE2). Using classical site density functional theory (SDFT) and structural bioinformatics methods, we investigate binding and conformational properties of these complexes and study the overlooked role of water-mediated interactions. Analysis of the three-dimensional reference interaction site model (3DRISM) of SDFT indicates that water mediated interactions in the form of additional water bridges strongly increases the binding between SARS-CoV-2 spike protein and hACE2 compared to SARS-CoV-1-hACE2 complex. By analyzing structures of SARS-CoV-2 and SARS-CoV-1, we find that the homotrimer SARS-CoV-2 S receptor-binding domain (RBD) has expanded in size, indicating large conformational change relative to SARS-CoV-1 S protein. Protomer with the up-conformational form of RBD, which binds with hACE2, exhibits stronger intermolecular interactions at the RBD-ACE2 interface, with differential distributions and the inclusion of specific H-bonds in the CoV-2 complex. Further interface analysis has shown that interfacial water promotes and stabilizes the formation of CoV-2/hACE2 complex. This interaction causes a significant structural rigidification of the spike protein, favoring proteolytic processing of the S protein for the fusion of the viral and cellular membrane. Moreover, conformational dynamics simulations of RBD motions in SARS-CoV-2 and SARS-CoV-1 point to the role in modification of the RBD dynamics and their impact on infectivity.

Renormalized site density functional theory for models of ion hydration.

The development of accurate statistical mechanics models of molecular liquid systems is a problem of great practical and fundamental importance. Site-density functional theory (SDFT) is one of the promising directions in this area, but its success hinges upon the ability to efficiently reconcile the co-existence of two distinct intra- and inter-molecular interaction regimes in a molecular liquid. The renormalized formulation of SDFT (RSDFT), which we have recently developed, resolves this problem by introducing an additional potential field variable that decouples two interaction scales and maps the molecular liquid problem onto the effective simple liquid mixture. This work provides a critical assessment of RSDFT for the hydrated ion system-a problem that historically has always been one of the most difficult cases for SDFT applications. Using a two-site model of water, we perform a comprehensive analysis of hydrated alkali metal and halogen ions, including both structural and free energy based characteristics. The results indicate that RSDFT provides a significant improvement over conventional three-dimensional reference interaction site model implementations and may prove useful in coarse grained simulations based on two-site solvent models.

Photoelectron Spectroscopy and Theoretical Study on Monosolvated Cyanate Analogue Clusters ECX-·Sol (ECX- = NCSe-, AsCSe-, and AsCS-; Sol = H2O, CH3CN).

Six monosolvated cyanate analogue clusters ECX-·Sol (ECX- = NCSe-, AsCSe-, and AsCS-; Sol = H2O and CH3CN) were investigated using negative ion photoelectron spectroscopy (NIPES). NIPES experiments show that these clusters possess similar spectra overall compared to their respective isolated ECX- anions but shift to higher electron binding energy with CH3CN solvent, stabilizing the excess electrons slightly more than H2O. For the ECX-·H2O series, vertical detachment energies and their increments relative to the bare species are measured to be 3.700/0.370, 3.085/0.415, and 3.085/0.430 eV for NCSe-, AsCSe- and AsCS-, respectively, while the corresponding values in the ECX-·CH3CN series are 3.835/0.505, 3.145/0.475, and 3.135/0.480 eV. Ab initio electronic structure calculations indicate that the excess charges were located at the terminal N and Se atoms in NCSe- and migrated to the central C atom in AsCSe- and AsCS-. For NCSe-, the solvation is driven by the interactions with the two negatively charged terminal ends, while for AsCSe- and AsCS-, the solvation revolves around the interactions with the central C atom, where all the excess negative charge is concentrated. Two nearly degenerate isomers for NCSe-·H2O are identified, one forming a single strong N···H-O hydrogen bond (HB) and the other featuring a bidentate HB with two hydroxyl H atoms pointing to N and Se ends. In contrast, the negative central C atom in AsCSe-/AsCS- allows the formation of a bifurcated HB with H2O. Similar effects are observed for the acetonitrile case, in which the three H atoms of the methyl group interact with the two negatively charged terminal ends in NCSe-, while preferring to bind to the central negative carbon atom in AsCSe-/AsCS-. The different binding motifs derived in this work may suggest different solvation properties in NCSe- versus AsCSe-/AsCS- with the former anion leading to asymmetric solvation at the N end of the solute, while the latter species creates more "isotropic" solvation around the central C equatorial plane.

From NWChem to NWChemEx: Evolving with the Computational Chemistry Landscape.

Since the advent of the first computers, chemists have been at the forefront of using computers to understand and solve complex chemical problems. As the hardware and software have evolved, so have the theoretical and computational chemistry methods and algorithms. Parallel computers clearly changed the common computing paradigm in the late 1970s and 80s, and the field has again seen a paradigm shift with the advent of graphical processing units. This review explores the challenges and some of the solutions in transforming software from the terascale to the petascale and now to the upcoming exascale computers. While discussing the field in general, NWChem and its redesign, NWChemEx, will be highlighted as one of the early codesign projects to take advantage of massively parallel computers and emerging software standards to enable large scientific challenges to be tackled.

Distonic radical anion species in cysteine oxidation processes.

Oxidation of cysteine residues constitutes an important regulatory mechanism in the function of biological systems. Much of this behavior is controlled by the specific chemical properties of the thiol side-chain group, where reactions with reactive oxygen species take place. Herein, we investigated the entire cysteine oxidation cycle Cys-SH → Cys-SOnH (n = 1, 2, and 3) using cryogenic negative ion photoelectron spectroscopy and quantum-chemical calculations. The conventional view of the first reversible oxidation step (n = 1) is associated with sulfenate species. Yet our results indicate that an alternative option exists in the form of a novel distonic radical anion, ˙OS-CH2CH(NH2)-COO-, with an unpaired electron on the thiol group and excess negative charge on the carboxylate group. Higher order oxidation states (n = 2 and 3) are thought to be associated with irreversible oxidative damage, and our results show that excess negative charge in those cases migrates to the -SOn- group. Furthermore, these species are stable towards 1e oxidation, as opposed to the n = 1 case that undergoes intra-molecular proton transfer. The molecular level insights reported in this work provide direct spectroscopic evidence of the unique chemical versatility of Cys-sulfenic acid (Cys-SOH) in post-translational modifications of protein systems.

Cryogenic "Iodide-Tagging" Photoelectron Spectroscopy: A Sensitive Probe for Specific Binding Sites of Amino Acids.

This work showcases cryogenic and temperature-dependent "iodide-tagging" photoelectron spectroscopy to probe specific binding sites of amino acids using the glycine-iodide complex (Gly·I-) as a case study. Multiple Gly·I- isomers were generated from ambient electrospray ionization and kinetically isolated in a cryogenic ion trap. These structures were characterized with temperature-dependent "iodide-tagging" negative ion photoelectron spectroscopy (NIPES), where iodide was used as the "messenger" to interpret electronic energetics and structural information of various Gly·I- isomers. Accompanied by theoretical computations and Franck-Condon simulations, a total of five cluster structures have been identified along with their various binding motifs. This work demonstrates that "iodide-tagging" NIPES is a powerful general means for probing specific binding interactions in biological molecules of interest.

Chemical bond effects in classical site density functional theory of inhomogeneous molecular liquids.

Intra-molecular interactions or chemical bonds represent one of the main distinguishing characteristics of molecular fluids. Development of accurate and practical methods to treat these effects is one of the long standing problems in classical site density functional theory (SDFT). One particular instance when these issues become particularly severe is the case of classical interaction potentials with auxiliary sites or dummy atoms. In this situation, current SDFT implementations, such as the three-dimensional reference interaction site model, lead to nonphysical results. We re-examine this issue in this work using our recent reformulation of SDFT (Valiev and Chuev, J. Stat. Mech.: Theory Exp. 2018, 093201). We put forward a simple practical solution to this problem and illustrate its utility for the case of spherical solutes in diatomic liquids.

Photoelectron spectroscopy of solvated dicarboxylate and alkali metal ion clusters, M+[O2C(CH2)2CO2]2- [H2O]n (M = Na, K; n = 1-6).

We present results of combined experimental photoelectron spectroscopy and theoretical modeling studies of solvated dicarboxylate species (-O2C(CH2)2CO2-) in complex with Na+ and K+ metal cations. These ternary clusters serve as simple models for the investigation of aqueous ion/solute specific effects that play an important role in biological systems. The experimental characterization of these systems was performed in the presence of up to six solvating waters. In both Na+ and K+ cases, we observe the presence of one major broad band that gradually shifts to higher electron binding energy (EBE) with an increasing number of waters. In the Na+ case further detailed analysis of experimental spectra was performed using ab initio calculations. In particular, we have identified the structures of the lowest energy clusters whose EBE values match well the major band in the experimental spectra. Our results show that evolution of an aqueous solvation shell emphasizes the coordination of the negatively charged carboxylate groups accompanied by simultaneous interaction with metal cations. Calculations also indicate that in the solvation range investigated experimentally (up to 6 waters), Na+ retains direct contact with the dicarboxylate species, i.e. a contact ion-pair (CIP) complex. Preliminary modeling studies show evidence of an alternative solvent separated ion-pair complex once the solvation range approaches 8 waters, however its energy still remains above that of (∼7-8 kcal/ mol-1) the CIP complex. At a higher number of waters (n = 3 for Na+ and n = 5 for K+), the experimental spectra also show the development of a weak low energy band. Its origin cannot be precisely identified. Our calculations in the Na+ case point out the existence of a quaternary complex consisting of Na+, H2O, OH- and a singly protonated dicarboxylate anion (HO2C(CH)2CO2-). Such a complex appears to be stabilized in the solvation range corresponding to the appearance of the low EBE band and does match its peak, even though the energy of such a complex is fairly high compared to the ternary structure.

Formation of (HCOO-)(H2SO4) Anion Clusters: Violation of Gas-Phase Acidity Predictions.

Sulfuric acid is commonly known to be a strong acid and, by all counts, should readily donate its proton to formate, which has much higher proton affinity. This conventional wisdom is challenged in this work, where temperature-dependent negative ion photoelectron spectroscopy and theoretical studies demonstrate the existence of the (HCOO-)(H2SO4) pair at an energy slightly below the conventional (HCOOH)(HSO4-) structure. Analysis of quantum-mechanical calculations indicates that a large proton affinity difference (∼36 kcal/mol), favoring proton transfer to formate, is offset by the gain in intermolecular interaction energy between HCOO- and H2SO4 through the electron delocalization and formation of two strong hydrogen bonds. However, this stabilization comes with a severe entropic penalty, requiring the two species in the precise alignment. As a result, the population of (HCOO-)(H2SO4) drops significantly at higher temperatures, rendering (HCOOH)(HSO4-) to be the dominant species. This phenomenon is consistent with the photoelectron data, which shows depletion in the spectra assigned to (HCOO-)(H2SO4), and has also been verified by ab initio molecular dynamics (AIMD) simulations.

Combined quantum-mechanical molecular mechanics calculations with NWChem and AMBER: Excited state properties of green fluorescent protein chromophore analogue in aqueous solution.

Combined quantum mechanical molecular mechanics (QM/MM) calculations have become a popular methodology for efficient and accurate description of large molecular systems. In this work we introduce our development of a QM/MM framework based on two well-known codes-NWChem and AMBER. As an initial application area we are focused on excited state properties of small molecules in an aqueous phase using an analogue of the green fluorescent protein (GFP) chromophore as a particular test case. Our approach incorporates high level coupled cluster theory for the analysis of excited states providing a reliable theoretical analysis of effects of an aqueous solvation environment on the photochemical properties of the GFP chromophore. Using a systematic approach, which involves comparison of gas phase and aqueous phase results for different protonation states and conformations, we resolve existing uncertainties regarding the theoretical interpretation of experimental data. We observe that the impact of aqueous environment on charged states generally results in blue shifts of the absorption spectra, but the magnitude of the effect is sensitive to both protonation state and conformation and can be rationalized based on charge movement into the area of higher/lower external electrostatic potentials. At neutral pH levels the experimentally observed absorption signal is most likely coming from the phenol protonated form. Our results also show that the high level electron correlated method is essential for a proper description of excited states of GFP. © 2017 Wiley Periodicals, Inc.

Structures and energetics of hydrated deprotonated cis-pinonic acid anion clusters and their atmospheric relevance.

Pinonic acid, a C10-monocarboxylic acid with a hydrophilic -CO2H group and a hydrophobic hydrocarbon backbone, is a key intermediate oxidation product of α-pinene - an important monoterpene compound in biogenic emission processes that influences the atmosphere. Molecular interaction between cis-pinonic acid and water is essential for understanding its role in the formation and growth of pinene-derived secondary organic aerosols. In this work, we studied the structures, energetics, and optical properties of hydrated clusters of the cis-pinonate anion (cPA-), the deprotonated form of cis-pinonic acid, by negative ion photoelectron spectroscopy and ab initio theoretical calculations. Our results show that cPA- can adopt two different structural configurations - open and folded. In the absence of waters, the open configuration has the lowest energy and provides the best agreement with the experiment. The added waters, which mainly interact with the negatively charged -CO2- group, gradually stabilize the folded configuration and lower its energy difference relative to the most stable open-configured structure. Thermochemical and equilibrium hydrate distribution analyses suggest that the mono- and di-hydrates are likely to exist in humid atmospheric environments with high populations. The detailed molecular description of cPA- hydrated clusters unraveled in this study provides a valuable reference for understanding the initial nucleation process and aerosol formation involving organics containing both hydrophilic and hydrophobic groups, as well as for analyzing the optical properties of those organic aerosols.

Deprotonated Dicarboxylic Acid Homodimers: Hydrogen Bonds and Atmospheric Implications.

Dicarboxylic acids represent an important class of water-soluble organic compounds found in the atmosphere. In this work we are studying properties of dicarboxylic acid homodimer complexes (HO2C(CH2)nCO2(-)[HO2C(CH2)nCO2H], n = 0-12), as potentially important intermediates in aerosol formation processes. Our approach is based on experimental data from negative ion photoelectron spectra of the dimer complexes combined with updated measurements of the corresponding monomer species. These results are analyzed with quantum-mechanical calculations, which provide further information about equilibrium structures, thermochemical parameters associated with the complex formation, and evaporation rates. We find that upon formation of the dimer complexes the electron binding energies increase by 1.3-1.7 eV (30.0-39.2 kcal/mol), indicating increased stability of the dimerized complexes. Calculations indicate that these dimer complexes are characterized by the presence of strong intermolecular hydrogen bonds with high binding energies and are thermodynamically favorable to form with low evaporation rates. Comparison with the previously studied HSO4(-)[HO2C(CH2)2CO2H] complex (J. Phys. Chem. Lett. 2013, 4, 779-785) shows that HO2C(CH2)2CO2(-)[HO2C(CH2)2CO2H] has very similar thermochemical properties. These results imply that dicarboxylic acids not only can contribute to the heterogeneous complexes formation involving sulfuric acid and dicarboxylic acids but also can promote the formation of homogeneous complexes by involving dicarboxylic acids themselves.

Probing the early stages of solvation of cis-pinate dianions by water, acetonitrile, and methanol: a photoelectron spectroscopy and theoretical study.

cis-Pinic acid is one of the most important oxidation products of α-pinene--a key monoterpene compound in biogenic emission processes. Molecular level understanding of its interaction with water in cluster formation is an important and necessary prerequisite for ascertaining its role in the aerosol formation processes. In this work, we studied the structures and energetics of the solvated clusters of cis-pinate (cis-PA(2-)), the doubly deprotonated dicarboxylate of cis-pinic acid, with H2O, CH3OH, and CH3CN by negative ion photoelectron spectroscopy and ab initio theoretical calculations. We found that cis-PA(2-) prefers being solvated alternately on the two -CO2(-) groups with increase of solvent coverage, a well-known solvation pattern that has been observed in microhydrated linear dicarboxylate dianion (DCn(2-)) clusters. Experiments and calculations further reveal an intriguing feature for the existence of the asymmetric type isomers for cis-PA(2-)(H2O)2 and cis-PA(2-)(CH3OH)2, in which both solvent molecules interact with only one of the -CO2(-) groups, a phenomenon that has not been observed in DCn(2-) water clusters and exhibits that the subtle effect of the rigid four-membered carbon ring brought on the cis-PA(2-) solvation. The dominant interactions between cis-PA(2-) and solvent molecules form bidentate O(-)...H-O H-bonds for H2O, O(-)...H-O and O(-)...H-C H-bonds for CH3OH, and tridentate O(-)...H-C H-bonds for CH3CN. The formation of inter-solvent H-bonds between H2O and CH3CN is found to be favorable in mixed solvent clusters, different from that between H2O and CH3OH. These findings have important implications for understanding the mechanism of cluster growth and the formation of atmospheric organic aerosols, as well as for rationalizing the nature of structure-function relationship of proteins containing carboxylate groups in various solvent environments.

Negative ion photoelectron spectra of ISO3-, IS2O3-, and IS2O4- intermediates formed in interfacial reactions of ozone and iodide/sulfite aqueous microdroplets.

Three short-lived, anionic intermediates, ISO3-, IS2O3-, and IS2O4-, are detected during reactions between ozone and aqueous iodine/sulfur oxide microdroplets. These species may play an important role in ozone-driven inorganic aerosol formation; however their chemical properties remain largely unknown. This is the issue addressed in this work using negative ion photoelectron spectroscopy (NIPES) and ab initio modeling. The NIPE spectra reveal that all of the three anionic species are characterized by high adiabatic detachment energies (ADEs) - 4.62 ± 0.10, 4.52 ± 0.10, and 4.60 ± 0.10 eV for ISO3-, IS2O3-, and IS2O4-, respectively. Vibrational progressions with frequencies assigned to the S-O symmetric stretching modes are discernable in the ground state transition features. Density functional theory calculations show the presence of several low-lying isomers involving different bonding scenarios. Further analysis based on high level CCSD(T) calculations reveal that the lowest energy structures are characterized by the formation of I-S and S-S bonds and can be structurally viewed as SO3 linked with I, IS, and ISO for ISO3-, IS2O3-, and IS2O4-, respectively. The calculated ADEs and vertical detachment energies are in excellent agreement with the experimental results, further supporting the identified minimum energy structures. The obtained intrinsic molecular properties of these anionic intermediates and neutral radicals should be useful to help understand their photochemical reactions in the atmosphere.

How Anion Chaotrope Changes the Local Structure of Water: Insights from Photoelectron Spectroscopy and Theoretical Modeling of SCN(-) Water Clusters.

The behavior of charged solute molecules in aqueous solutions is often classified using the concept of kosmotropes ("structure makers") and chaotropes ("structure breakers"). There is a growing consensus that the key to kosmotropic/chaotropic behaviors lies in the local solvent region, but the exact microscopic basis for such differentiation is not well-understood. This issue is examined in this work by analyzing size selective solvation of a well-known chaotrope, a negatively charged SCN(-) molecule. Combining experimental photoelectron spectroscopy measurements with theoretical modeling, we examine evolution of solvation structure up to eight waters. We observe that SCN(-) indeed fits the description of weakly hydrated ion, and its solvation is heavily driven by stabilization of water-water interaction network. However, the impact on water structure is more subtle than that associated with "structure breaker". In particular, we observe that the solvation structure of SCN(-) preserves the "packing" structure of the water network but changes local directionality of hydrogen bonds in the local solvent region. The resulting effect is closer to that of "structure weakener", where solute can be readily accommodated into the native water network, at the cost of compromising its stability due to constraints on hydrogen bonding directionality.

Probing microhydration effect on the electronic structure of the GFP chromophore anion: Photoelectron spectroscopy and theoretical investigations.

The photophysics of the Green Fluorescent Protein (GFP) chromophore is critically dependent on its local structure and on its environment. Despite extensive experimental and computational studies, there remain many open questions regarding the key fundamental variables that govern this process. One outstanding problem is the role of autoionization as a possible relaxation pathway of the excited state under different environmental conditions. This issue is considered in our work through combined experimental and theoretical studies of microsolvated clusters of the deprotonated p-hydroxybenzylidene-2,3-dimethylimidazolinone anion (HBDI(-)), an analog of the GFP chromophore. Through selective generation of microsolvated structures of predetermined size and subsequent analysis of experimental photoelectron spectra by high level ab initio methods, we are able to precisely identify the structure of the system, establish the accuracy of theoretical data, and provide reliable description of auto-ionization process as a function of hydrogen-bonding environment. Our study clearly illustrates the first few water molecules progressively stabilize the excited state of the chromophore anion against the autodetached neutral state, which should be an important trait for crystallographic water molecules in GFPs that has not been fully explored to date.