Concentration-Dependent Binding of Small Ligands to Multiple Saturable Sites in Membrane Proteins

Membrane proteins are primary targets for most therapeutic indications in cancer and neurological diseases, binding over 50% of all known small molecule drugs. Understanding how such ligands impact membrane proteins requires knowledge on the molecular structure of ligand binding, a reasoning that has driven relentless efforts in drug discovery and translational research. Binding of small ligands appears however highly complex involving interaction to multiple transmembrane protein sites featuring single or multiple occupancy states. Within this scenario, looking for new developments in the field, we investigate the concentration-dependent binding of ligands to multiple saturable sites in membrane proteins. The study relying on docking and free-energy perturbation provides us with an extensive description of the probability density of protein-ligand states that allows for computation of thermodynamic properties of interest. It also provides one- and three-dimensional spatial descriptions for the ligand density across the protein-membrane system which can be of interest for structural purposes. Illustration and discussion of the results are shown for binding of the general anesthetic sevoflurane against Kv1.2, a mammalian ion channel for which experimental data are available.

computation of μ does not depend upon the choice of concentration, so long as the same thermodynamic state is used for the solution and gas phases, we estimated the excess potential by considering one sevoflurane molecule embedded into a water box of 60 x 60 x 60 Å 3 . W n j * was computed here by taking into considering the whole ligand-channel-membrane system.
All FEP calculations were performed in NAMD 2.9 4 by considering the Charmm-based parameters for sevoflurane as devised by Barber et al. 9 Starting from channel-membrane equilibrated systems containing bound sevoflurane as resolved from docking, forward transformation were carried out by varying the coupling parameter in steps of 0.05 (or for convergence purposes, in steps of 0.025 at final stages of the process). Each transformation then involved a total of 80 windows, each spanning over 32512 steps of simulation. For the purpose of improving statistics, free-energy estimates and associated statistical errors were determined using the simple overlap sampling (SOS) formula 10 based on at least two independent FEP runs. Specifically for ligand-protein calculations, the free-energy change W 1 * for singly-occupied sites was computed as a FEP process that involves ligand coupling to a vacant site. Differently, for doubly-occupied sites, W 2 * was computed as a two-step FEP process involving ligand coupling to a vacant site W 1 * followed by binding of a second ligand at the preoccupied site W 2|1 * . Because is a W 2 * state function, the stepwise approach is equivalent to a single-step process involving simultaneous coupling of two ligands to the protein site that is, W 2 * =W 1 * +W 2|1 * . The colvars module 11 in NAMD 2.9 was used to apply the harmonic restraint potentials when computing these quantities. As described in the main text, the value of W n j * depends on the parameters of the restraint potential adopted in the FEP calculation ie., the reference positions of the ligands in the bound state {R 1 * , ..., R n j * } and the magnitude of force constants {k 1 ,..., k n j } . By minimizing the contribution of the restraint potential to the binding free-energy W n j * , Roux and coworkers in which, ⟨R 1 ⟩,... ,⟨R n j ⟩ and ⟨δ R 1 2 ⟩ ,...,⟨δ R n j 2 ⟩ are respectively the equilibrium average positions for each of the n j bound ligands at site j and their corresponding mean-square fluctuations. Here, these parameters were estimated from the docking configuration space and the resulting force constants, in the range of 1.0 to 10.0 kcal/mol/Å 2 , were considered for computations of the bound state.

Convergence of sampling.
Here, a per-site measure for the ensemble of docking solutions effectively sampled in FEP was determined by quantifying the overlap o( A j , B j ) between the configuration space in both calculations 13 for A j and B j denoting covariance matrices associated respectively to FEP and docking samples at site j and, A j 1/2 and B j 1/ 2 their square roots. Specifically, A j and B j were determined as symmetric 3×3 covariance matrices for centroid positions R j of the ligand at site j and their square roots were solved from the column major eigenvectors { R l , R 2 , R 3 } of the rotation matrix R and the associated eigenvalues is expectedly 1 for identical samplings and 0 for orthogonal spaces. (5) is standard and follows from the coordinate transform (r n )→( R n ,Ω n , I n ) involving the centroid positions R n , orientations Ω n and internal I n degrees of freedom of the n ligands.

Derivation of main text equation (5). Derivation of equation
In this case, the Jacobian of the transformation does not depend on the R n degrees of freedom allowing for cancellation of (Ω n , I n ) contributions. The effective volume [ ∏ i=1 n ( 2π β k i ) 3 2 ] thus results from the 3n -dimensional Gaussian integral appearing in line 4.

Derivation of main text equation (16).
The spatial projection along the z direction of the system shown in equation (16) derives as where, A (z)=Δ x Δ y is the total area of the water-membrane region along the Cartesian x and y directions. gives the probability density and

Coarse
establishes the link between Κ(n) and the standard binding free-energy Δ G o (n) associated to each of the states satisfying supplementary equation (S1).