Unravelling the Turn‐On Fluorescence Mechanism of a Fluorescein‐Based Probe in GABAA Receptors

Abstract GABAA (γ‐aminobutyric acid type A) receptors are ligand‐gated ion channels mediating fast inhibitory transmission in the mammalian brain. Here we report the molecular and electronic mechanism governing the turn‐on emission of a fluorescein‐based imaging probe able to target the human GABAA receptor. Multiscale calculations evidence a drastic conformational change of the probe from folded in solution to extended upon binding to the receptor. Intramolecular ππ‐stacking interactions present in the folded probe are responsible for quenching fluorescence in solution. In contrast, unfolding within the GABAA receptor changes the nature of the bright excited state triggering emission. Remarkably, this turn‐on effect only manifests for the dianionic prototropic form of the imaging probe, which is found to be the strongest binder to the GABAA receptor. This study is expected to assist the design of new photoactivatable screening tools for allosteric modulators of the GABAA receptor.


Table of contents
Section S1 Selection of the density functional 2 Tables S1, S2 Section S2 Ground state of the fluorescent probe 1 5 Figure S1 Section S3 Emission of the fluorescent probe 1 6 Figure S2 Section S4 Setting up of the systems 7 Table S3 Section S5 MD simulations 8 Figures S3, S4, S5, S6 Section S6 Binding energy analysis 11 Tables S4, S5 Section S7 QM/MM MD simulations 13 Scheme S1, Figure S7 References 16

Section S1 SELECTION OF THE DENSITY FUNCTIONAL
Here we present a benchmark of several density functionals to describe the absorption and emission of 1 against ADC (2).To this aim, we used only the two fluorophores of interest, 2a and 2b.
Absorption.Five density functionals were selected: two hybrid functionals (PBE0, B3LYP), [1,2] two meta-hybrid functionals (M06, M062X), [3] and one long-range corrected hybrid functional (CAM-B3LYP). [4]10][11] Then, the molecules are vertically excited from the optimized S0 geometry to the first excited state (S1) using time-dependent density functional theory (TD-DFT).These results were compared to those obtained on the B3LYP-optimized geometry with ADC(2)/def2-TZVP (our reference method) with implicit water solvation (COSMO), using TURBOMOLE v7.0. [12]Regarding the comparability of PCM and COSMO implicit solvation models: We did not change the dielectric scaling factor (x = 0) used by Gaussian 16 for the integral equation formalism PCM method as it was shown to produce comparable results to the COSMO implicit solvation for ionic solutes (which we have with 2a and 2b). [13]Table S1 collects the energies, oscillator strengths and orbital character (natural transition orbitals, NTO) of the S1 when using TD-DFT and ADC (2).
The five TD-DFT calculations coincide with the results of ADC(2), predicting the same character for the S1 state (  →   * ) for 2a and 2b.In all cases, we see that the S1 state is the brightest state at the S0 geometry.The absorption energies are all overestimated by TD-DFT with LR-PCM, with deviations from 0.49 eV (B3LYP) to 0.64 eV (CAM-B3LYP) in 2a, and from 0.26 eV (B3LYP) to 0.52 eV (CAM-B3LYP) in 2b.The cLR-PCM method overestimates the absorption energies even further.Nevertheless, the trend for all the different functionals is consistent.From these results, we considered the two DFT functionals with the smallest deviations from the ADC(2) results, B3LYP and M06, as most suitable candidates to evaluate the emission process.
Emission.Emission energies and associated orbital characters are investigated with the best performing functionals for absorption, M06 and B3LYP, [1][2][3] as well as CAM-B3LYP, [4] benchmarked against ADC (2).We first optimized the minimum energy geometry on the S1 PES of 2a and 2b with TD-DFT/def2-TZVP and LR-PCM.As for absorption, the emission energies were calculated with the default LR-PCM implicit solvation as well as with the corrected approach cLR-PCM for B3LYP and M06.
As it can be seen in Choice of the level of theory.According to the results obtained above, we consider that B3LYP is an adequate level of theory to perform single point TD-DFT calculations on the OG488 fluorophore, also within the GABAAR.We used the default LR-PCM model, as we do not find appreciable changes when using the corrected LR approach.The results of the static calculations on 2a and 2b are presented in Table 1 of the main manuscript.

Section S2 GROUND STATE OF THE FLUORESCENT PROBE 1
The Cartesian coordinates (see repository: https://phaidra.univie.ac.at/o:1383031) of the fluorescent probe gabazine (Gzn) conjugated with Oregon Green (OG488) [14] (1) were generated as discussed in the Introduction of the main manuscript, with two of its more probable protonation  The ground state geometry for both forms of 1 was optimized using DFT with the B3LYP [15] functional, the def2-SVP [6] basis set, the Grimme D3 dispersion correction [7] and PCM water solvation.The calculations were carried out in the Gaussian 16, v. A.01 suite. [5]

Section S3 EMISSION OF THE FLUORESCENT PROBE 1
The geometries of the first excited state S1 of 1a and 1b were optimized using with TD-B3LYP-D3/def2-SVP/LR-PCM(water), which is the protocol described in Section S1, apart from the basis set size.The choice of the B3LYP functional is justified from the comparison with other DFTfunctionals against ADC(2) as method of reference on 2a and 2b.We note that TD-DFT has been extensively used in the literature to describe the electronic structure of fluorescein-based [37,38] and organic photoprobes. [39,40]The results for the S1 to S0 emission for the whole fluorescent probes 1a and 1b can be seen in Table 1 in the main text.To further validate the DFT results also for 1, we also calculated the emission of 1b with ADC(2)/def2-SVP using the B3LYP-optimized geometry with implicit water solvation (COSMO), as implemented in TURBOMOLE v7.0. [12]We obtain an emission energy of 1.38 eV with an oscillator strength (fosc ) of 0.010, which is comparable to that obtained with B3LYP-D3 (Table 1) of 1.17

Section S4 SETTING UP OF THE SYSTEMS
The Cartesian coordinates of the receptor were obtained from the cryogenic electron microscopy (cryo-EM) solved structure of the human full-length heteromeric GABAAR(α1/β3/γ2) in complex with the anxiolytic alprazolam (Xanax), the transmitter GABA, and the megabody Mb38 (PDB id.

6HUO
). [16] This structure, among other reasons mentioned in the "Results and Discussion" section, was used in this work, because the structure of the same receptor in the presence of another antagonist (bicuculline, see Table S3, PDB id.6HUK), like Gzn, presented a distorted orientation of α1Arg67 by one of the aromatics in bicuculline.Thus, to preserve the former residue as one of the pharmacophoric points for the binding of Gzn at the beginning of our simulation, we selected the structure in the presence of GABA.Alprazolam, GABA, and Mb38 were removed for our simulations.Since the secondary structure of the intracellular domain (ICD) of all subunits is unknown (α1, amino acids 313-391; β3, amino acids 308-423; γ2, amino acids 323-400) the peptide linker of sequence Ser-Gln-Pro-Ala-Arg-Ala-Ala used experimentally by Miller et al. [17] was introduced as replacement for each of the subunits ICD. [18,19][22][23] Care was taken to include the five disulfide bonds (one per subunit) present in the extracellular domain (ECD) of the receptor.The GABAAR was embedded into an explicit 1-palmitoyl-2-oleoyl-snglycero-3-phosphatidylcholine (POPC) lipid membrane (190 lipids on the lower membrane, 200 in the upper membrane) after reorientation with the OPM server. [24]The system was solvated with TIP3P water [25] molecules with a box of size 130x125x195 Å using the server CHARMM-GUI. [26- 30] n addition to 21 neutralizing Cl -ions, a concentration of 0.15 mM of brine (210 Na + and Cl - ions) was added to reproduce the experimental conditions of Sakamoto et al. [14] Before adding the probe, the unbound system was relaxed by means of standard short MD simulations in AMBER18. [31]The placing of the probe inside the GABA binding sites was done with the hypothesis of preserving the two pharmacophoric points of the binding of GABA: the salt-bridges between GABA and residues α1Arg67 and β3Glu155.

Section S5 MD SIMULATIONS
Force field parameters for the fluorescent probe were computed on the DFT(B3LYP-D3/def2-SVP) optimized structure of 1a and 1b using the program Antechamber implemented in AmberTools19. [31]AMBER atom types and RESP charges were used.The receptor was described via the AMBER ff14SB force field [33] and the lipids via the Lipid17 force field. [34]TIP3P water [25] parameters and the ion parameters of Li/Merz were used. [35]Periodic boundary conditions were employed with a cutoff of 10 Å for non-bonded interactions and the Particle Mesh Ewald method for long-range electrostatics.All the MD simulations were run using the AMBER18 software package. [31]e ground state optimized probe forms were neutralized with one or two Na + ions for the 1b and 1a form, respectively, and solvated in a box of TIP3P water [25] with a minimum distance of 40 Å from the solute to the boundaries of the box using tLEaP. [31]Then, the system was minimized in three stages: protons, then solvent molecules and counterions, and finally, the whole system.
Harmonic restraints of 200 kcal mol -1 Å -2 were imposed to the rest of the atoms in each stage.
After that, the system was thermally equilibrated for 20 ps from 100 K to 300 K using the Berendsen thermostat and restraining the solute with 20 kcal mol -1 Å -2 .Finally, and before production, the system was equilibrated in 6 steps from NVT to NPT ensemble for 38 ps.The system was further simulated for 100 ns with the Langevin thermostat with a time step of 0.2 fs.
The simulation of the receptor in the membrane was executed following a similar protocol as for the probe in water.The protocol only deviates in the heating and equilibration process: The system was heated 200 ps from 0 K to 300 K using the Langevin thermostat with harmonic restraint imposed on all the atoms of receptor, probes, and lipids (10 kcal mol -1 Å -2 ); the equilibration was executed in 10 steps for 5 ns in total.Finally, the system was further simulated for 100 ns.
Several parameters of the MD simulation, apart from the ones pictured in the main work, were analyzed.Figure S3 shows the distances d1-d3 along the 100 ns simulation.The choice of the initial receptor structure (PDB id.6HUK vs. 6HUO) was validated in Figure S4. Figure S5 pictures snapshots of 1a and 1b at the GABA binding site B of the GABAAR(α1/β3/γ2) and Figure S6 displays the root-means square deviation along the 100 ns simulation.6HUK [16] (magenta), PDB id.6HUO [16] (purple), and one snapshot from the MD simulation of 1a (green) (B) Plot of the root-means square deviation (RMSD, Å) of the backbone carbon atoms of the receptor along the MD simulation of 1a of 100 ns against PDB id.6HUK (magenta line) and PDB id.6HUO (purple line).The RMSD of 6HUK against 6HUO is shown as a red dot.

Section S6 BINDING ENERGY ANALYSIS
The binding free energies of 1a and 1b in binding sites A and B were calculated using the MM/ISMSA program. [36]It estimates the solvent-corrected binding energies (Table S4) as well as their per-residue decomposition (Table S5).The values reported in the results are averages.The last 15 ns of the MD simulations were considered for the calculation of the binding energies in accordance with their root-means square deviation (RMSD) plots (Figure S6).

Section S7 QM/MM MD SIMULATIONS
As region, the OG488 part of 1a (41 atoms) was treated using the third order semiempirical density-functional tight binding method DFTB3 [41][42][43] implemented in AMBER18. [44]The environment (solvent, protein, and lipid membrane) was treated as point charges in the calculation of the energy of the QM part.First, taking the 100 ns MD simulations of 1a bound to GABAAR and in a box of water molecules, one snapshot per system was randomly selected from the major conformational clusters (58% in site A in the membrane context and 58% in aqueous solution).
The cluster analysis was performed using the k-means algorithm implemented in cpptraj for 3 clusters.Second, QM/MM MD simulations were run for 20 ps in the ground state (electronic state S0) for both systems (see Scheme S1).SHAKE and Ewald were deactivated in the QM region and a cut-off of 15 Å (water)/30 Å (GABAAR) was applied from the center of the QM region for the electrostatic term.Third, from the former trajectories on the ground state, 11 snapshots from the last 10 ps were taken and used as input geometries for a 500 fs Born-Oppenheimer MD simulation on the S1 excited state (11 simulations for each system).The QM region, the OG488 part of 1a, was treated with TD-DFT (TD-BP86-D3/def2-SVP) [45] using the interface AMBER18:ORCA. [46,47]The SlowConv algorithm and the deactivation of the SOSCF algorithm was required to get SCF convergences.8 roots were included in the calculation and the energy of the potential energy surface (PES) of the S1 state was used for the propagation.We also tested the use of B3LYP as DFT-functional for the QM/MM MD propagation on the S1 PES.See Table S6 for the comparison.Although we obtain higher energies with B3LYP compared to BP86, we obtained similar results for oscillator strength (fosc), charge transfer numbers (CT), participation ratio (PR) numbers for the initial/final electron delocalization, and transition character (not shown).
We conclude that the geometry of the molecule is well enough described by BP86.The emission spectrum for 1a was calculated using 5 snapshots from the last 250 fs of each of the 11 S1 excited state Born-Oppenheimer MD simulations in aqueous solution and the protein context (55 snapshots each, see Scheme S1, see repository for xyz-structures: https://phaidra.univie.ac.at/o:1383031).Each snapshot was then used for a TD-DFT single point calculation using Gaussian 16, v. A.01 [5] (B3LYP-D3/def2-SVP) and an electrostatic cut-off of 55 Å (see Section S1 for the discussion about the DFT-functional selection).The QM region is composed of the whole probe in aqueous solution and of only the fluorophore in the protein context.The BP86 functional is sufficient for an efficient S1 propagation, but not for the calculation of the emission spectrum, as it delivers drastic red-shifts.The emission energy is calculated as the vertical excitation energy from ground to first excited state in the optimized S1 geometry.[50] The spectra were convoluted by Lorentzian functions weighted by the oscillator strength with a full width at half maximum of 0.10 eV using TheoDORE [51] (see Figure 5A).To one more time validate the use of the B3LYP functional, a second spectrum of 1a was convoluted based on M06-D3/def2-SVP calculations (30 snapshots).The comparison between B3LYP and M06 is shown in Figure S7.
The spectrum with the M06 functional is shifted to higher energies compared to B3LYP, but the characteristics of the spectrum stay consistent.
states based on a pH value of 7.4 (see FigureS1and Scheme 1) and taking the pharmacophoric points for the binding of GABA into account.The relevant protonatable moieties are the carboxylate of Gzn, the iminium cation of Gzn, the carboxylate of OG488, and the xanthene of OG488.The latter either has a hydroxy group in the monoanionic form (1b) or is fully deprotonated in the dianionic form (1a) (Scheme 1).

Figure S1 .
Figure S1.(A) Chemical structures of three protolytic forms (G1, G2, G3) of 1's gabazine (Gzn).The size of 1 was reduced to Gzn and parts of the linker for the pKa calculations.(2) Plot of the distribution of G1, G2, and G3 against the pH value calculated with ChemSketch.
eV and fosc of 0.001.Moreover, the natural transition orbitals (NTOs) for both methods have the same nπ* transition state character (n COO − ← π OG * ) (shown in Figures 4 for 1a and in S2 for 1b).This result corroborates that the B3LYP functional is sufficiently accurate for the system under study.

Figure S3 .
Figure S3.Plot of distances d1-d3 (Å) along the MD simulation of 100 ns for 1a (orange) and 1b (blue) in site A (A) and in site B (B).

Figure S4 .
Figure S4.(A)Overlayed structures of the aromatic box at the GABA binding site A for PDB id.6HUK[16] (magenta), PDB id.6HUO[16] (purple), and one snapshot from the MD simulation of 1a (green) (B) Plot of the root-means square deviation (RMSD, Å) of the backbone carbon atoms of the receptor along the MD simulation of 1a of 100 ns against PDB id.6HUK (magenta line) and PDB id.6HUO (purple line).The RMSD of 6HUK against 6HUO is shown as a red dot.

Figure S5 .
Figure S5.(A-B) Overlaid snapshots of 1a (A) and 1b (B) at the GABA binding site B of the GABAAR(α1/β3/γ2) along the MD simulation of 100 ns.

Figure S6 .
Figure S6.Plot of the root-means square deviation (RMSD, Å) of the backbone carbon atoms of the receptor (black line) along the MD simulation of 100 ns, as well as RMSD plots for 1a (A) and 1b (B) in site A and in site B for all atoms except hydrogen.

Table S2
* ) and energies by +0.44 eV (similar to the shift seen with the absorption calculations).The errors for the non-emissive 2b are less pronounced: ADC(2) gives a value of 760 nm (1.63 eV), again in reasonable agreement with B3LYP (841 nm, 1.47 eV).We stress that for 2b we found a surprisingly low energy S1 state with ADC(2) (1.63 eV).Out of

Table S4 .
Computed binding energy (kcal mol -1 ) for 1a and 1b at the GABA sites A and B.*

Table S5 .
Computed binding energy (kcal mol -1 ) decomposed per residue for 1a and 1b at the GABA sites A and B.* Energies computed using the last 15 ns of the MD simulations.The initial structures were not relaxed before the calculation of the binding energies.

Table S6 .
First excited electronic state (S1) properties of 5 snapshots on the S1 PES using BP86 and B3LYP as DFT functional in the trajectory of 1a at the GABAAR.*The emission energies were calculated with B3LYP as described in the text.Only the fluorophore is in the QM region.The snapshots were randomly selected.