#  Identifying the minimal sets of distance restraints for FRET-assisted
protein structural modeling

Zhuoyi Liu  Department of Mechanical Engineering and Materials Science, Yale
University, New Haven, Connecticut, 06520, USA  Integrated Graduate Program in
Physical and Engineering Biology, Yale University, New Haven, Connecticut,
06520, USA  Alex T. Grigas  Integrated Graduate Program in Physical and
Engineering Biology, Yale University, New Haven, Connecticut, 06520, USA
Graduate Program in Computational Biology and Bioinformatics, Yale University,
New Haven, Connecticut, 06520, USA  Jacob Sumner  Integrated Graduate Program
in Physical and Engineering Biology, Yale University, New Haven, Connecticut,
06520, USA  Graduate Program in Computational Biology and Bioinformatics, Yale
University, New Haven, Connecticut, 06520, USA  Edward Knab  Department of
Chemistry, Yale University, New Haven, Connecticut, 06520, USA  Caitlin M.
Davis  Department of Chemistry, Yale University, New Haven, Connecticut,
06520, USA  Corey S. O’Hern  Department of Mechanical Engineering and
Materials Science, Yale University, New Haven, Connecticut, 06520, USA
Integrated Graduate Program in Physical and Engineering Biology, Yale
University, New Haven, Connecticut, 06520, USA  Graduate Program in
Computational Biology and Bioinformatics, Yale University, New Haven,
Connecticut, 06520, USA  Department of Physics, Yale University, New Haven,
Connecticut, 06520, USA  Department of Applied Physics, Yale University, New
Haven, Connecticut, 06520, USA

Manuscript Pages: 15

Total Manuscript Figures: 5

Total Manuscript Tables: 1

Supporting Information Pages: 8

Total Supporting Information Figures/Tables: 6/0

Abstract: Proteins naturally occur in crowded cellular environments and
interact with other proteins, nucleic acids, and organelles. Since most
previous experimental protein structure determination techniques require that
proteins occur in idealized, non-physiological environments, the effects of
realistic cellular environments on protein structure are largely unexplored.
Recently, Förster resonance energy transfer (FRET) has been shown to be an
effective experimental method for investigating protein structure in vivo.
Inter-residue distances measured in vivo can be incorporated as restraints in
molecular dynamics (MD) simulations to model protein structural dynamics in
vivo. Since most FRET studies only obtain inter-residue separations for a
small number of amino acid pairs, it is important to determine the minimum
number of restraints in the MD simulations that are required to achieve a
given root-mean-square deviation (RMSD) from the experimental structural
ensemble. Further, what is the optimal method for selecting these inter-
residue restraints? Here, we implement several methods for selecting the most
important FRET pairs and determine the number of pairs $N_{r}$ that are needed
to induce conformational changes in proteins between two experimentally
determined structures. We find that enforcing only a small fraction of
restraints, $N_{r}/N\lesssim 0.08$, where $N$ is the number of amino acids,
can induce the conformational changes. These results establish the efficacy of
FRET-assisted MD simulations for atomic scale structural modeling of proteins
in vivo.

Significance: Determining protein structure in vivo is essential for
understanding protein function. Most protein structures have been studied in
non-physiological conditions using x-ray crystallography, NMR spectroscopy,
and cryo-electron microscopy. Thus, we do not know whether the cellular
environment significantly affects protein structure. We emphasize the benefits
of FRET-assisted molecular dynamics simulations in characterizing protein
structure in vivo at the atomic scale. We identify the minimum number of FRET
pairs that can induce conformational changes in several proteins, including
one that has been characterized using in-cell NMR.

Keywords: in vivo protein structure $|$ FRET experiments $|$ molecular
dynamics simulations $|$ cellular crowding

##  1 Introduction

Knowing the three-dimensional structure of proteins enables us to understand
the biophysical mechanisms that control protein function, protein-protein
interactions, and cell signaling. Nearly all protein structures that have been
determined experimentally to date have been characterized under idealized,
non-physiological conditions. For example, proteins have been crystallized
into non-native, solid phases for x-ray scattering experiments 1, and proteins
have been dissolved into dilute, non-physiological buffers for NMR
spectroscopy or cryo-electron microscopy 2, 3. However, proteins carry out
their functions in cellular environments that are significantly different from
these in vitro conditions. The cellular environment is crowded with a non-
solvent packing fraction of $0.3$-$0.4$ that includes nucleic acids,
carbohydrates, lipids, organelles, and other components 4, 5, 6. This
environment impacts the physical properties of proteins, including the
protein’s radius of gyration, melting temperature, and rotational diffusion
coefficient 7, 8, 9, 10. Currently, there are over $200$,$000$ in vitro
protein structures deposited in the Protein Data Bank (PDB) 11, as well as
more than $10^{6}$ computational models predicted by AlphaFold and RosettaFold
12, 13. To date, in vivo structures have been obtained for only three proteins
(TTHA1718, GB1, and ubiquitin) using in-cell NMR 14, 15, 16. In-cell NMR is
not widely used since it is difficult to distinguish the isotopic labeled
target protein from its environment, which leads to a low signal-to-noise
ratio and sensitivity 17, 18, 19. Another experimental method for solving
protein structure in vivo is cryo-electron tomography 20, however, this
technique requires proteins to be confined to a thickness less than $500$ nm,
which limits this method to only a subset of cell types and membrane-
associated proteins 21, 22. In addition, there have been numerous
computational studies of proteins in vivo. All-atom molecular dynamics (MD)
simulations have investigated protein structure in the presence of nearly all
components of the cell cytoplasm 23, 24. However, current force fields have
been calibrated to in vitro protein structures and we do not have accurate
potentials for interactions between proteins and nucleic acids, ribosomes, and
organelles 25, 26. Thus, we do not yet have a quantitative, atomistic-level
understanding of protein structure in cells.

Förster resonance energy transfer (FRET) can be used to determine the
separations between donor and acceptor chromophores attached to a pair of
amino acids in a given protein by measuring their energy transfer efficiency.
The distribution of separations between the two amino acids can then be
deduced by calculating the configuration space volumes sampled by the donor
and acceptor chromophores 27, 28, 29. FRET pair labeling is highly specific
and offers high accuracy for the inter-residue separations with uncertainties
less than 2-4 Å 30, 31. Additionally, FRET-labeled proteins can be directly
expressed or injected into cells, enabling single molecule measurements of
protein structure, dynamics, and stability in vivo 32, 7, 33. However, most
FRET experimental studies only obtain inter-residue separations for a small
number of amino acid pairs in a given protein 7, 10, 33. To investigate the
atomic scale structure of proteins in vivo, inter-residue separations obtained
from FRET experiments can be incorporated into molecular dynamics (MD)
simulations. Current all-atom MD simulations of proteins using in vitro
solution conditions have been shown to sample in vitro protein structures
obtained from x-ray crystallography and NMR spectroscopy 34, 35. By including
a sufficient number of inter-residue restraints from FRET studies of proteins
in vivo into the MD force fields that were developed for in vitro solution
conditions, it may be possible to sample the in vivo structural ensemble of
proteins.

Several key questions must be answered before FRET-assisted structural
modeling of protein structure in vivo can be employed. In particular, given
the large number of possible FRET pairs, what is the minimum number of amino
acid pairs for which we need distance restraints to achieve a given accuracy
for the root-mean-square deviations (RMSD) between the Cα positions in the
restrained simulations and those in the experimental structures? However, we
do not yet have access to high-resolution in vivo protein structures. Thus, we
will first develop the methodology for carrying out restrained MD simulations
for protein structural modeling using proteins that are found in multiple
conformational states in vitro. The initial conformational state will be
metastable in the MD force field without restraints, and we will induce a
conformational change in the protein by incorporating restraints between amino
acid pairs that are satisfied in the target state. To our knowledge, there
have only been a few studies aimed at identifying the most important
restraints for efficiently moving to the conformational ensemble of the target
state 36. For example, currently we do not know the minimal number of
restraints and which restraints are necessary to achieve a given RMSD in the
Cα positions from the target state, and how random selection of given number
of restraints compares to other methods for selecting restraints.

Here, we select three proteins that each can take on two, distinct
conformational states to develop the restrained MD simulation methodology: T4
Lysozyme (172L/1L69), Phosphoglycerate Kinase (2XE6/2Y3I), and Adenylate
Kinase (4AKE/1AKE), where the first and second PDBIDs indicate the initial and
target structures, respectively. We initialize the MD simulations with the
crystal structure of the initial state, add a given number of Cα distance
restraints, and calculate the Cα RMSD relative to the target structure. For
each target, we test four methods (normal mode analysis, the largest Cα
separation method, the largest change in pairwise separation method, and
linear discriminant analysis) for selecting the restraints and compare the Cα
RMSD to the target structure to that obtained from random selection of the
restraints. We also vary the number of restraints to determine the minimum
number of restraints needed to achieve a given Cα RMSD from the target. Two of
the methods (normal mode analysis and identifying the largest Cα separations
in the initial structure) do not use information about the target structure,
whereas the other two methods compare the initial and target structures to
identify the most important restraints. We find that for the three proteins
that take on two distinct conformational states in vitro, we can induce the
conformational changes using only a small fraction of restraints,
$N_{r}/N\lesssim 0.02$, where $N$ is the number of amino acids in the protein.
For the proteins that we considered, this result corresponds to $1$-$5$
restraints, which is a number that can readily be achieved in FRET
experiments. In addition, we studied one of the few proteins that has been
characterized using in-cell NMR spectroscopy: the B1 domain of protein G. We
need a slightly larger fraction of restraints, $N_{r}/N\sim 0.08$ (or $5$
restraints), to change the protein conformation from the initial in vitro
structure (2N9K) to the in-cell NMR structure (7QTS). In general, the largest
change in pairwise separation method, which has information about the target
structure, yields the lowest values for the Cα RMSD for a given $N_{r}$. If
the restraint selection method does not have information about the target
structure, the Cα RMSD is still lower than that for random selection. These
results establish the feasibility of FRET-assisted structural modeling of
proteins in vivo using restrained MD simulations.

##  2 Materials and Methods

### Selected Proteins

We identified three proteins from the Protein Data Bank (PDB) that can be
crystallized into two distinct conformational states in vitro and can be used
to develop the restrained MD simulation methodology: T4 Lysozyme (T4L),
Phosphoglycerate Kinase (PGK), and Adenylate Kinase (AK). (See Table 1.) 37,
38, 39, 40, 41, 42. These proteins have been studied extensively in the
context of protein conformational changes 43. While the selected targets are
in vitro structures, recent studies suggest that the differences between the
in vivo and in vitro structures for PGK are largely caused by the hinge motion
of the two domains that occurs between the initial and target in vitro
structures 33, 8. Previous studies have also suggested that the crowded
cellular environment can stabilize the target in vitro structure for AK 44.
For the restrained MD simulations, we require that the initial structure is
stable and the target structure is at least metastable over long time scales
in MD simulations, which ensures that the unrestrained dynamics does not
induce the transition from the initial state to the target. (See Fig. S1 in
Supporting Information.) We also consider one protein that has both an in
vitro and in-cell NMR structure, the B1 domain of protein G (GB1) 45, 16. We
use the first model from the in vitro NMR bundle as the initial structure and
the first model of the in-cell NMR bundle as the target structure. (See Fig.
S6 in Supporting Information.) We find similar results for the restrained and
unrestrained MD simulations when we initialize the system using the other NMR
models.

Protein | initial structure | target structure | $N$ | Cα ${\rm RMSD}$ (Å)
---|---|---|---|---
T4 Lysozyme (T4L) | 172L | 1L69 | 162 | 3.95
Phosphoglycerate Kinase (PGK) | 2XE6 | 2Y3I | 413 | 4.08
Adenylate Kinase (AK) | 4AKE | 1AKE | 214 | 5.96
B1 domain of protein G (GB1) | 2N9K | 7QTS | 57 | 2.86
Table 1: We list the selected proteins including their PDBIDs for the initial
and target structures, number of amino acids $N$, and the Cα RMSD (Eq. 2)
between the initial and target structures. For GB1, the Cα RMSD is calculated
between the first model of the in vitro NMR bundle and the first model of the
in-cell NMR bundle.

### Cα RMSD

To compare the protein structures from the restrained MD simulations (i.e.
structure $S_{i}$) and the target structure (i.e. structure $S_{j}$), we
define the Cα separation vector for the $\beta$th amino acid,

${\vec{\Delta}}(S_{i},S_{j};\beta)=\vec{r}_{i,\beta}-\vec{r}_{j,\beta},$ (1)

where $\vec{r}_{i,\beta}$ is the position of the $C_{\alpha}$ atom on amino
acid $\beta$ in structure $S_{i}$. We define the root-mean-square deviation in
the Cα positions between two structures $S_{i}$ and $S_{j}$ as

${\rm
RMSD}(S_{i},S_{j})=\sqrt{\frac{1}{N}\sum_{\beta=1}^{N}\Delta^{2}(S_{i},S_{j},\beta)}.$
(2)

When comparing two structures $S_{i}$ and $S_{j}$, we typically align them
(i.e. rotate one of them) to achieve the minimum value of the Cα ${\rm
RMSD}(S_{i},S_{j})$ for a given $S_{i}$ and $S_{j}$. We can also calculate the
Cα RMSD between the restrained and target structures averaged over an ensemble
of restrained structures for each amino acid $\beta$:

${\rm
RMSD}(\\{S_{i}\\},S_{j},\beta)=\sqrt{\frac{1}{N_{s}}\sum_{i=1}^{N_{s}}\Delta^{2}(S_{i},S_{j},\beta)},$
(3)

where $N_{s}$ is the number of protein structures in the restrained ensemble
$\\{S_{i}\\}$ and $S_{j}$ represents the target structure.

### Restraint selection methods

Below, we describe the methods that we employ for selecting the Cα distance
restraints. Linear spring restraints will be added to MD simulations of the
initial protein structure, where the rest lengths are obtained from the target
structure, as discussed below. Three of the restraint selection methods,
random selection, normal mode analysis, and identifying the largest Cα
separations in the initial structure, do not use information about the target
structure. Two additional methods, identifying the largest changes in the Cα
pair separations between the initial and target structures and linear
discriminant analysis, compare the initial and target structure to identify
the most informative restraints.

#### a) Random selection

In a protein with $N$ amino acids, there are $N_{p}=N(N-1)/2$ distinct amino
acid pairs. As shown in Fig. 1a, the pair separations between Cα atoms can be
represented using a symmetric distance matrix $R_{\beta\delta}$, where
$\beta$,$\delta=1,\ldots,N$. To establish a baseline for the performance of
the restrained MD simulations, we will first consider random selection of the
restraints. In this approach, we exclude amino acid pairs that are too close
($R_{\beta\delta}<R_{g}$, where $R_{g}$ is the radius of gyration), and pairs
for which at least one amino acid is buried with relative solvent accessible
surface area ${\rm rSASA}<0.1$ since such pairs will preclude FRET
measurements 46, 47, which reduces the pool of restraints to approximately
$N_{p}/3$ amino acid pairs. For each number of restraints $N_{r}\ll N_{p}$ we
consider, we select $100$ sets of $N_{r}$ restraints randomly from the pool of
allowed pairs.

#### b) Normal mode analysis

For this method, we assume that the normal modes of vibration of the initial
structure provide information about how the protein transitions from the
initial to the target structure. We follow the methodology used in other
recent work aimed at selecting the most informative FRET pairs for restrained
MD simulations of proteins 36. First, the method constructs a coarse-grained
elastic network from the atomic positions of the initial protein structure and
calculates the normal modes of the network 48. Then, the initial structure is
displaced by a linear combination of the ten lowest frequency modes with
amplitudes that are inversely proportional to the frequency and have
proportionality constants $-1<\xi<1$ that are chosen randomly. The normal
modes are then recalculated on the displaced structure and the structure is
perturbed again using a linear combination of the ten lowest frequency modes
with amplitudes chosen as before. This process of successive displacements
along the lowest frequency normal modes is continued until $10^{3}$ structures
are obtained 49 and then repeated ten times with independently generated
random mode amplitudes to yield a total of ${\cal N}=10^{4}$ structures. $100$
of the structures for T4L are shown in Fig. 1b.

To identify the amino acid pairs that should be restrained, we seek to
minimize the RMSD of the positions of the Cα atoms (Eq. 2) among the protein
structures in the ensemble, where the weights for each structure in the
ensemble are controlled by the size of the fluctuations in the amino acid pair
separations among structures. Deviations in the pair separations between two
structures $S_{i}$ and $S_{j}$ can be quantified using

$\chi^{2}(S_{i},S_{j})=\sum_{\\{(\beta,\delta)\\}}\left(\frac{\Delta
R_{\beta\delta}(S_{i},S_{j})}{R_{\beta\delta}^{i}+R^{j}_{\beta\delta}}\right)^{2},$
(4)

where $\Delta
R_{\beta\delta}(S_{i},S_{j})=|R_{\beta\delta}^{i}-R^{j}_{\beta\delta}|$,
$\\{(\beta,\delta)\\}$ is a given set of amino acid pairs, and
$R_{\beta\delta}^{i}$ is the distance between Cα atoms on amino acids $\beta$
and $\delta$ on structure $S_{i}$. To estimate the probability of observing a
mean-square deviation in the pair separations larger than
$\chi^{2}(S_{i},S_{j})$, we calculate

$P(S_{i},S_{j})=\int_{\chi^{2}(S_{i},S_{j})}^{\infty}f(N_{m},\chi^{2})\,d\chi^{2},$
(5)

where $N_{m}$ is the number of amino acid pairs in the set
$\\{(\beta,\delta)\\}$ and $f(N_{m},\chi^{2})$ is the chi-squared distribution
with $N_{m}$ degrees of freedom. To quantify the average Cα RMSD over the
ensemble, we calculate

$\langle{\rm RMSD}\rangle=\frac{1}{\cal N}\sum_{i=1}^{\cal
N}\frac{\sum_{j=1}^{\cal N}P(S_{i},S_{j}){\rm
RMSD}(S_{i},S_{j})}{\sum_{j=1}^{{\cal N}}P(S_{i},S_{j})}.$ (6)

To minimize $\langle{\rm RMSD}\rangle$, the structures for which the selected
amino acid pairs have large deviations in their pair separations (i.e. small
$P(S_{i},S_{j})$) should possess large Cα $\rm RMSD$, and the structures for
which the selected amino acid pairs have small deviations (i.e. large
$P(S_{i},S_{j})$) should possess small Cα $\rm RMSD$. To select the set of
amino acid pairs that minimize $\langle{\rm RMSD}\rangle$, we add one pair to
the set of optimal pairs at a time. We start with identifying the single pair
$(\beta^{*},\delta^{*})$ that minimizes $\langle{\rm RMSD}\rangle$ and use it
as the restraint for $N_{r}=1$. We then consider two possible pairs, but fix
one of the pairs to be $(\beta^{*},\delta^{*})$ and find the new pair
$(\beta^{**},\delta^{**})$ in the set
$\\{(\beta^{*},\delta^{*}),(\beta^{**},\delta^{**})\\}$ that minimizes
$\langle{\rm RMSD}\rangle$. We use these two pairs as the set of restraints
for $N_{r}=2$. This process continues until we have $N_{r}=1,\ldots,N_{\rm
max}$, where $N_{\rm max}/N=0.03$ for T4L, PGK, and AK.

Figure 1: Illustration of several distance restraint selection methods, using
T4L as the example protein. (a) The symmetric distance matrix
$R_{\beta\delta}$, where the color gradient from dark to light indicates
increasing Cα separations between amino acids $\beta$ and $\delta$. (b) For
the normal mode analysis method, we generate an ensemble of $10^{4}$
structures ($100$ structures are shown in green) that have been displaced from
the initial structure along a random superposition of normal modes
corresponding to the ten lowest frequencies. A single amino acid pair that
minimizes the $C_{\alpha}$ RMSD among structures in the ensemble is
highlighted in red and orange in each structure. (c) We identify the amino
acid pairs with the largest Cα separations $R_{\beta\delta}$ in the initial
(green) structure. The pair with the maximum $R_{\beta\delta}$ is highlighted
in red and orange. (d) We can also select restraints by identifying the
largest changes in the pairwise Cα separations between the initial (green) and
target (blue) structures. The Cα atoms of the pair with the maximum $\Delta
R_{\beta\delta}(S_{i},S_{j})$ are shown in red and orange. (e) $20$
conformations from unrestrained MD simulations starting from the initial and
target structures are shown in green and blue, respectively. Using linear
discriminant analysis, we identify the collective variable ${\hat{W}}$ that
maximizes the cross-correlation of the two ensembles.

#### c) Identifying the largest Cα separations in the initial structures

For this method (i.e. the largest Cα separation method) for selecting
important restraints, we identify the amino acid pairs with the largest Cα
separations, or the maximum $R_{\beta\delta}$ over all pairs
$\\{(\beta,\delta)\\}$ in the initial protein structure. (For example, in Fig.
1c, we show the amino acid pair with the largest Cα separation in the x-ray
crystal structure for T4L (172L).) After identifying the pair $(\beta,\delta)$
with the largest Cα separation, we find the pair
$(\beta^{{}^{\prime}},\delta^{{}^{\prime}})$ with the second largest Cα
separation. This process is then continued to find a set of pairs with the
largest Cα separations. However, we seek to identify a minimal set of amino
acid pairs without redundant structural information. Thus, we carried out
unrestrained MD simulations starting from the initial structure and determined
the Pearson correlation between the pairwise Cα separations. In particular, we
include $(\beta^{{}^{\prime}},\delta^{{}^{\prime}})$ in the pool of selected
pairs if the Pearson correlation $\rho$ between $R_{\beta\delta}$ and
$R_{\beta^{{}^{\prime}}\delta^{{}^{\prime}}}$ satisfies $|\rho|<0.9$. We also
require that the Cα atoms of the new pair are not already in the pool of
selected restraints. Thus, we first add the pair $(\beta,\delta)$ with the
largest Cα separation to the pool of restraints ($N_{r}=1$). We add
$(\beta^{\prime},\delta^{\prime})$ with the second largest Cα separation to
the pool of restraints ($N_{r}=2$) as long as it is uncorrelated with
$(\beta,\delta)$ and does not include Cα atom $\beta$ or $\delta$. If so, we
consider the pair $(\beta^{\prime\prime},\delta^{\prime\prime})$ with the next
largest separation. We follow this process until we have sets of restraints
with $N_{r}=1,\ldots,5$.

#### d) Identifying the largest change in pairwise Cα separations between the
initial and target structures

For this method (i.e. the largest change in pairwise separation method), we
assume that the target structure is known. We identify the amino acid pair
with the largest change in pairwise Cα separations between the initial and
target structures $S_{i}$ and $S_{j}$ (i.e. the maximum $\Delta
R_{\beta\delta}(S_{i},S_{j})$). (For example, in Fig. 1d, we show the amino
acid pair with the largest change in Cα separations between the initial (172L)
and target (1L69) x-ray crystal structures for T4L. Note that the amino acids
with the largest change in pairwise separation are not the same as those with
the largest Cα separation.) To identify a set of non-redundant restraints, we
successively implement new restraints in MD simulations for pairs with the
largest deviation in the pair separation from the target. In particular,
assume that pair $(\beta,\delta)$ has the largest $\Delta
R_{\beta\delta}(S_{i},S_{j})$ in the unrestrained MD simulations starting from
the initial structure. We then carry out MD simulations with $R_{\beta\delta}$
restrained to the target value. We identify the pair
$(\beta^{{}^{\prime}},\delta^{{}^{\prime}})$ with the largest $\Delta
R_{\beta^{{}^{\prime}}\delta^{{}^{\prime}}}(S_{i},S_{j})$ and carry out
restrained MD restraints enforcing the target values for both
$R_{\beta\delta}$ and $R_{\beta^{{}^{\prime}}\delta^{{}^{\prime}}}$. We also
require that the Cα atoms of the new pair are not the same as any of the atoms
in the current pool of restraints. We use $(\beta,\delta)$ as the restraint
for $N_{r}=1$ and use $(\beta,\delta)$ with
$(\beta^{{}^{\prime}},\delta^{{}^{\prime}})$ as the set of restraints for
$N_{r}=2$. This process continues until we have $N_{r}=1,\ldots,N_{\rm max}$,
where $N_{\rm max}/N=0.03$ for T4L, PGK, and AK and $0.08$ for GB1.

#### e) Linear discriminant analysis

For the linear discriminant analysis method of selecting restraints, we also
assume that the target structure is known. We seek to identify the pairwise
separation between Cα atoms that can serve as a collective variable enabling
the protein to move from the initial to the target structure 50. The
inspiration for this method is the linear discriminant analysis method that
was used to identify the collective variables for small molecule
conformational changes 51. Here, we apply the method in the context of large
conformational changes in proteins.

Linear discriminant analysis requires a distribution of protein structures,
not a single structure. Thus, we first carried out short $20$ ns unrestrained
MD simulations starting from both the initial and target structures to sample
an ensemble of structures near the initial and target structures. In Fig. 1e,
we show these two distributions of structures for T4L as an example. In this
case, we aim to find the collective variable that tracks the hinge closure
motion of the two domains. Each protein conformation can be described by the
set of all pairwise separations between Cα atoms,

${\vec{X}}=(R_{\beta\delta},R_{\beta^{{}^{\prime}}\delta^{{}^{\prime}}},...)^{T},$
(7)

where the length of the vector is the number of distinct amino acid pairs
$N_{p}$ minus the ones that are too close together ($R_{\beta\delta}<R_{g}$)
and the ones where at least one amino acid is buried. We can describe the
distribution of protein conformations sampled near the initial and target
structures as $\rho_{A}({\vec{X}}_{A})$ and $\rho_{A}({\vec{X}}_{B})$,
respectively. We aim to find the direction ${\hat{W}}$ in the space of Cα
separations such that its projection onto the cross-correlation matrix
(between initial and target distributions) is maximized, while its projection
onto the covariance matrix for the initial or target distributions is
minimized. To determine ${\hat{W}}$, we first calculate the mean separations
$\langle{\vec{X}}_{A}\rangle$ and $\langle{\vec{X}}_{B}\rangle$ from the
distributions of the initial and target structures. Second, we define the
(cross-correlation) scattering matrix

$S_{b}=(\langle{\vec{X}}_{A}\rangle-\langle{\vec{X}}_{B}\rangle)(\langle{\vec{X}}_{A}\rangle+\langle{\vec{X}}_{B}\rangle)^{T},$
(8)

which quantifies the covariance of the initial distribution with the target
distribution, and define the projection of ${\hat{W}}$ onto $S_{b}$ as
${\hat{W}}^{T}S_{b}{\hat{W}}$. We also quantify the covariance of the initial
and target distributions, separately:

$\Sigma_{A,B}=(X_{A,B}-\langle X_{A,B}\rangle)(X_{A,B}-\langle
X_{A,B}\rangle)^{T}.$ (9)

We can then calculate the (direct correlation) scattering matrix,

$S_{w}=\Sigma_{A}+\Sigma_{B},$ (10)

with projection ${\hat{W}}^{T}S_{w}{\hat{W}}$. To find the direction
${\hat{W}}$ onto which the projections of the two distributions
$\rho_{A}({\vec{X}}_{A})$ and $\rho_{B}({\vec{X}}_{B})$ are best separated, we
can maximize the Rayleigh ratio:

$\mathcal{R}(\hat{W})=\frac{{\hat{W}}^{T}S_{b}{\hat{W}}}{{\hat{W}}^{T}S_{w}{\hat{W}}},$
(11)

which is equivalent to solving for the eigenvector ${\hat{W}}_{\lambda}$
associated with the largest eigenvalue $\lambda$ of $S_{w}^{-1}S_{b}$:

$S_{w}^{-1}S_{b}{\hat{W}}_{\lambda}=\lambda{\hat{W}}_{\lambda}.$ (12)

We identify the most important pairwise distances as those that have the
largest weights in ${\hat{W}}_{\lambda}$. For $N_{r}=1$, we choose the amino
acid pair $(\beta,\delta)$ with the largest weight. For $N_{r}=2$, we choose
the two amino acid pairs $(\beta,\delta)$ and
$(\beta^{{}^{\prime}},\delta^{{}^{\prime}})$ with the two largest weights, and
so on.

### Restrained MD simulations

To assess the performance of the FRET pair selection methods, we carry out
all-atom MD simulations starting from the initial structure and incorporating
the Cα-Cα separations of the selected pairs from the target structure as
linear spring restraints. We then monitor the RMSD of the Cα atom positions
(Eq. 2) from the target structure as a function of time. Unrestrained and
restrained MD simulations were performed using the AMBER99SB-ILDN force field
52, 53 in the GROMACS software package 54. The MD simulations were carried out
in a periodic dodecahedron-shaped box that is sufficiently large such that the
protein surface is at least $20$Å from the box edges. The simulation box was
solvated with water molecules modeled using TIP3P at neutral pH and $0.15$M
NaCl 55, 56. Short-range van der Waals and screened Coulomb interactions were
truncated at $1.2$ nm, while long-ranged electrostatic interactions were
tabulated using the Particle Mesh Ewald summation method. The LINCS algorithm
was used to constrain the bond lengths. We performed two energy minimization
runs to first relax the protein and then the water molecules and the protein
together using the steepest decent algorithm until the maximum net force
magnitude on an atom is smaller than $500$ $\text{kJmol}^{-1}\text{nm}^{-1}$.
We perform NVT simulations of the system for $20$ ns at temperature $T=300$K
using a velocity rescaling thermostat for sampling the canonical ensemble 57.
The equations of motion for the atomic coordinates and velocities are
integrated using a leapfrog algorithm with a $1$ fs time step.

For the restrained simulations, we employ a linear spring potential for each
restrained amino acid pair $(\beta,\delta)$:

$E_{r}(R_{\beta\delta})=\frac{k_{r}}{2}(R_{\beta\delta}-R^{0}_{\beta\delta})^{2},$
(13)

where $R^{0}_{\beta\delta}$ is the Cα-Cα separation for the $(\beta,\delta)$
pair in the target structure and the spring constant $k_{r}=5000$
$\text{kJmol}^{-1}\text{nm}^{-2}$ is chosen so that the averaged root-mean-
square deviation between $R_{\beta\delta}$ and $R^{0}_{\beta\delta}$ is
$<0.2$Å. (See Fig. S2 in Supporting Information.)

For the unrestrained MD simulations for each protein, we ran $N_{s}=100$
simulations for 20 ns starting from the initial structure, but with different
initial velocities for each atom randomly selected from a Maxwell-Boltzmann
distribution at $T=300$K. Using these simulations, we can calculate the
average Cα RMSD between the initial structures and the target structure, which
serves as a baseline for comparison to the results for the restrained MD
simulations. For the random restraint selection method at each $N_{r}$, we
randomly select $N_{s}=100$ sets of $N_{r}$ amino acid pairs and run one
restrained MD simulation (for 20 ns) for each set of restraint pairs. For all
other restraint selection methods, we carry out $N_{s}=50$ simulations for
each set of $N_{r}$ restraints, but with different random initial velocities.
The number of samples $N_{s}$ is chosen so that average Cα $\rm RMSD$
converges at large $N_{s}$ as shown in Fig. S3 in the Supporting Information.
For each $N_{r}$ and for each restraint selection method, we average the Cα
$\rm RMSD$ over all $N_{s}$ simulations. For the restrained MD simulations of
the in vitro protein structures, the Cα RMSD converges within the first $5$
ns, as shown in Fig. S4 in the Supporting Information. Therefore, we only use
the data generated from $5$ to $20$ ns for the calculations of the Cα $\rm
RMSD$. For the in vivo target GB1, we extended the MD simulations to $50$ ns
and use the time points from $20$ to $50$ ns to calculate the Cα RMSD. (See
Fig. S6c in Supporting Information.)

To determine the ability of the restraints to move the protein conformation
from the initial to the target structure, we need to measure a reference Cα
RMSD (i.e. RMSDr) that quantifies thermal fluctuations around the target
structure. For the proteins T4L, PGK, and AK, we determine RMSDr by running
unrestrained MD simulations starting from the target structure. We set RMSDr =
${\rm RMSD}(S_{i},S_{j})$, where $S_{i}$ is the central structure from the
unrestrained MD simulations starting from the target structure, and $S_{j}$ is
the target structure. We find the central structure by identifying the
structure that has the largest number of neighboring structures in the
ensemble using a cutoff Cα RMSD of 1Å. If the average Cα RMSD obtained from
the restrained MD simulations is below RMSDr, the restraints were successful
in moving the initial structure to the target structure since the differences
are comparable to thermal fluctuations observed near the target structure. For
the in-cell NMR structure GB1, we set RMSD${}_{r}=N_{m}^{-1}\sum_{i>j}{\rm
RMSD}(S_{i},S_{j})$, where $S_{i}$ and $S_{j}$ are distinct models in the in-
cell NMR bundle and $N_{m}$ is the number of distinct pairs of models.

##  3 Results

Does the choice of the amino acid restraints have a significant impact on
whether the initial protein structure can be moved to the target structure? As
an example, we consider two possible choices for single amino acid restraints
in restrained MD simulations of T4L. In Fig. 2, we plot the Cα ${\rm
RMSD}(S_{i},S_{j})$ as a function of time, where $S_{i}$ are the structures
sampled in the restrained MD simulations and $S_{j}$ is the target structure,
for two possible choices for single amino acid restraints. Although the two
restrained MD simulations start from the same initial structure (green) in
Fig. 2a, the Cα ${\rm RMSD}(S_{i},S_{j})$ display different time dependence.
The Cα RMSD from the MD simulations using the
$(\beta^{\prime},\delta^{\prime})$ restraint rapidly decreases below the
reference value, ${\rm RMSD}_{r}$, for the target structure. When we visualize
the final frame from this restrained MD simulation in Fig. 2c, we find that
the two subdomains of T4L are now closer together and there is strong
alignment between the final frame (green) and target structure (blue). In
contrast, the Cα RMSD does not decrease below RMSDr for the MD simulations
with the $(\beta,\delta)$ restraint, even though the restraint is well-
satisfied during the simulations. This example emphasizes the importance of
the restraint selection method. Below, we describe the performance of four
restraint selection methods in moving an initial structure toward a target
structure. We compare the performance of restraint selection methods that have
information about the target structure and those that do not, and benchmark
their performance to random restraint selection. In addition, we test the
methods on moving three proteins from one in vitro structure to another and
one protein from an in vitro structure to an in-cell structure.

Figure 2: (a) The top and bottom green images indicate the same starting
structures for restrained MD simulations that enforce two different single
restraints: $R_{\beta\delta}=R^{0}_{\beta\delta}$ (red) and
$R_{\beta^{\prime},\delta^{\prime}}=R^{0}_{\beta^{\prime},\delta^{\prime}}$
(orange). The starting structures are aligned with the target structures shown
in blue. (b) ${\rm RMSD}(S_{i},S_{j})$ as a function of time during restrained
MD simulations for the restraints in (a), where the $S_{i}$ are structures
sampled at each time during the restrained MD simulations and $S_{j}$ is the
target structure. The black dashed line indicates ${\rm RMSD}_{r}$ of the
target. (c) The final frames at $20$ ns from the restrained MD simulations
with $R_{\beta\delta}=R^{0}_{\beta\delta}$ (top, red) and
$R_{\beta^{\prime},\delta^{\prime}}=R^{0}_{\beta^{\prime},\delta^{\prime}}$
(bottom, orange). The structures from the final frames of the restrained MD
simulations (green) are aligned with the target structure (blue).

#### Performance of restraint selection methods

To quantify the performance of the four restraint selection methods in moving
the initial structure to the target structure, in Fig. 3, we plot $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ averaged over structures $S_{i}$ from
the restrained MD simulations and number of samples $N_{s}$, where $S_{j}$ is
the target structure, as a function of the number of restraints $N_{r}$ for
several restraint selection methods and three proteins. In Fig. 3a, we show
the results for T4L. When $N_{r}=0$ (i.e. unrestrained MD simulations),
$\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}\sim 4$Å and the
structures sampled in the MD simulations remain far from the target structure.
When randomly selecting restraints, $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ decreases slowly with increasing
$N_{r}$ and does not reach ${\rm RMSD}_{r}$ even after five restraints have
been included. In contrast, if we use structural information about the target
to select the restraints (i.e. using linear discriminant analysis or
identifying the largest change in pairwise Cα separations between the initial
and target structures), $\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$
rapidly decreases to ${\rm RMSD}_{r}$ after adding only one restraint. The
normal mode analysis and largest Cα separation methods, which do not use
information about the target, achieve intermediate results; $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ decreases to ${\rm RMSD}_{r}$ within
$N_{r}=2$-$3$ restraints. Thus, both of these methods (normal mode analysis
and the largest Cα separation method) are more efficient than random selection
for moving the initial structure to the target for T4L.

Figure 3: The Cα RMSD, $\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$,
averaged over structures $S_{i}$ from the restrained MD simulations and number
of samples $N_{s}$, where $S_{j}$ is the target structure, is plotted versus
the number of restraints $N_{r}$ for (a) T4L, (b) PGK, and (c) AK. We show
results for several restraint selection methods: random selection (blue
circles), normal mode analysis (orange squares), the largest Cα separation
method (green crosses), the largest change in pairwise separation method (red
upward triangles), and linear discriminant analysis (purple stars). The
horizontal black dashed line indicates ${\rm RMSD}_{r}$ for each target.
Snapshots of the initial structures are shown for each protein in the upper
right corner of each panel.

To investigate the performance of restraint selection methods in moving an
initial structure to the target structure for a wider range of proteins, we
performed unrestrained and restrained MD simulations on two larger proteins,
PGK (Fig. 3b) and AK (Fig. 3c), and find similar results for the variation of
$\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ with $N_{r}$. (Note that
for AK, we needed to add a background of intra-domain restraints to
recapitulate the B-factor in the “unrestrained” MD simulations as shown in
Fig. S5 in Supporting Information.) For example, for random restraint
selection, $\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ slowly
decreases with increasing $N_{r}$, and $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}>{\rm RMSD}_{r}$ even for $N_{r}=5$.
For PGK and AK, all restraint selection methods (except linear discriminant
analysis) outperform random restraint selection. Normal mode analysis and the
largest Cα separation method, which do not include information about the
target, achieve lower values of $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ than that for random selection, but
their relative performance depends on the protein. For example, $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ decreases faster with $N_{r}$ for the
largest Cα separation method for PGK, whereas $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ decreases faster for the normal mode
analysis method for AK. As expected, the method of identifying the largest
change in pairwise Cα separations between the initial and target structures is
the best performing method with $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}\sim{\rm RMSD}_{r}$ after only a small
number of restraints ($5$ for PGK and $2$ for AK). The performance of the
restrained MD simulations based on the linear discriminant analysis method of
selecting restraints is comparable to that for random restraint selection even
though it incorporates structural information about the target. This result
likely stems from the fact that a nonlinear method is needed to identify the
collective variables in proteins.

Figure 4: The normalized Cα RMSD, $\overline{{\rm RMSD}}(N_{r})$, is plotted
versus $N_{r}/N$ for T4L (blue circles), PGK (orange triangles), and AK (red
squares) using (a) the largest change in pairwise Cα separation method and (b)
the largest Cα separation method for selecting the restraints.

#### Minimum fraction of restraints

These results demonstrate that the addition of a small number of pairwise
distance restraints in restrained MD simulations can effectively move a
protein from an initial structure to a target structure that is more than $4$
Å away. However, the minimal number of restraints required for $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}\sim{\rm RMSD}_{r}$ depends on the
protein. What controls the minimum number of restraints? To address this
question, we calculate the normalized Cα RMSD:

$\overline{{\rm RMSD}}(N_{r})=\frac{\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}(N_{r})-{\rm RMSD}_{r}}{\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}(0)-{\rm RMSD}_{r}},$ (14)

where $\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}(0)$ is $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ from unrestrained MD simulations with
$N_{r}=0$. Thus, $\overline{{\rm RMSD}}(N_{r})=1$ indicates that the protein
samples the ensemble of initial structures and $\overline{{\rm
RMSD}}(N_{r})=0$ indicates that the protein samples the target ensemble. In
Fig. 4, we show that $\overline{{\rm RMSD}}(N_{r})$ collapses for all three
proteins when plotted versus $N_{r}/N$ and using the optimal restraint
selection method (i.e. the largest change in pairwise separation). In this
case, only a small fraction of restraints, less than 1.5% of the protein size,
is required to move between the two in vitro protein structures. For the
normal mode analysis restraint selection method, which does not consider
information about the target, the collapse of ${\overline{\rm RMSD}}$ with
$N_{r}/N$ is not as complete, but ${\overline{\rm RMSD}}\sim 0$ for
$N_{r}/N\gtrsim 0.02$ for all proteins. Thus, even in the absence of complete
knowledge of the target structure, we can induce changes in protein structure
from an initial in vitro structure to a target in vitro structure using only a
small fraction of amino acid separation restraints.

Figure 5: (a) A ribbon diagram of the initial in vitro structure (green) of
GB1 aligned with its in-cell target structure (blue). The gray arrows indicate
the structural changes from the initial structure to the target structure in
the two top loop regions. (b) $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}(N_{r})$ plotted as a function of
$N_{r}/N$ for restrained MD simulations of GB1 using the largest change in
pairwise Cα separation between the initial and target structures method (red
upward triangles) for selecting restraints. The horizontal black dashed line
indicates ${\rm RMSD}_{r}$ for the target. (c) ${\rm
RMSD}(\\{S_{i}\\},S_{j},\beta)$ (Eq. 3) is plotted versus amino acid index
$\beta$, where $\\{S_{i}\\}$ are structures sampled in the restrained MD
simulations and $S_{j}$ is the target structure. The color from dark blue to
yellow indicates $N_{r}$. The two top loop regions in (a) with
$21\leq\beta\leq 26$ and $47\leq\beta\leq 52$ are highlighted in light gray.

#### Moving from initial in vitro structure to target in-cell structure for
GB1

To investigate whether the proposed methodology can also move an initial in
vitro structure to a target in vivo structure, we also carried out
unrestrained and restrained MD simulations for the B1 domain of protein G
(GB1) 16. We find that inter-bundle $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},S_{j}}\sim 3.0$Å between the in vitro NMR
models $S_{i}$ and the in-cell NMR models $S_{j}$ is rather small, mostly due
to small changes in the two top loop regions shown in Fig. 5a. Despite this,
the Cα RMSD between the in vitro and in-cell NMR bundles is larger than both
of the intra-bundle fluctuations. Can we use similar restraint selection
methods to those above to move from an in vitro structure to an in-cell
structure? In contrast to the in vitro target structures for T4L, PGK, and AK,
the in-cell target structure for GB1 is not a local minimum of the AMBER99SB-
ILDN force field. (See Figs. S8a and b in the Supporting Information.) As a
result, simulating in-cell protein structures is more challenging. In Fig. 5b,
we implement the largest change in pairwise separation method for identifying
restraints and calculate $\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$
as a function of $N_{r}$. We find that the target in-cell structures are
sampled for $N_{r}\geq 5$ in restrained MD simulations, which corresponds to
$N_{r}/N\sim 8$%. Unlike for the in vitro structures, where $\langle{\rm
RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ monotonically decreases as $N_{r}$
increases, $\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ exhibits non-
monotonic behavior with $N_{r}$ for GB1. To investigate this effect, we
calculate the Cα RMSD between the structures from the MD simulations and the
target for each amino acid $\beta$ individually (Eq. 3). As shown in Fig. 5c,
the deviations in the structures sampled in the unrestrained simulations and
target structure occur predominantly in the top two loop regions in Fig. 5a
with $21\leq\beta\leq 26$ and $47\leq\beta\leq 52$. The deviations between the
initial and target structures can be decreased in the loop regions to
$\lesssim 2$ Å for $N_{r}=5$. However, the structural deviations at amino acid
positions adjacent to the loop regions, for example, $17\leq\beta\leq 20$ and
$27\leq\beta\leq 31$, begin to increase as $N_{r}$ increases. These results
suggest that the restraints are competing with the force field when attempting
to move the protein to the in-cell structure. The non-monotonic behavior in
$\langle{\rm RMSD}(S_{i},S_{j})\rangle_{S_{i},N_{s}}$ gives rise to a higher
fraction of restraints that are necessary to move GB1 from the initial in
vitro structure to the in-cell structure.

##  4 Discussion

Numerous studies have emphasized that the in-cell environment can strongly
influence protein structure and interactions 4, 6, 7, 9, 59. However, the
complete atomic structure for proteins in the cellular environment has been
obtained for only three proteins using NMR spectroscopy 14, 15, 16. FRET has
also been employed to determine a small number of separations between amino
acids in several proteins in vivo 7, 10, 33. In contrast, MD simulations,
which have been calibrated to structures in the protein data bank, can provide
the atomic coordinates of proteins in idealized, in vitro conditions. Properly
calibrated force fields that would allow accurate all-atom descriptions of
protein conformations in vivo currently do not exist. For example, GTT WW
domain fails to refold in an all-atom cytoplasm computational model 23. Thus,
incorporating FRET measurements of amino acid separations from proteins in
vivo as restraints in MD simulations may dramatically improve our
understanding of protein conformational dynamics in vivo.

While current all-atom MD simulations allow us to sample in vitro protein
structures that match the x-ray crystal or NMR structures, we want to capture
the all-atom conformational dynamics of in vivo protein structures. To do
this, we can start with an initial in vitro protein structure, and add inter-
residue distances measured by FRET as restraints in the MD simulations.
However, an important question is how do we determine which amino acid pairs
should be measured in the FRET experiments and then incorporated into the
restrained MD simulations? In particular, what is the minimal number of
restraints needed to achieve a given Cα RMSD to the target structure and what
is the optimal method to select these restraints? To answer these questions,
we considered several proteins that can take on two distinct conformations in
vitro and one protein for which we have both the in vitro and in-cell NMR
structures. We employed four methods for selecting the restraints, varied the
number $N_{r}$ and type of restraints that are incorporated into the
restrained MD simulations, and compared the performance (i.e. Cα RMSD relative
to the target structure) of the restraint selection methods to random
selection.

We found that for the three proteins (T4L, PGK and AK) with two distinct
conformational states solved in vitro, it only takes a small fraction of
restraints ($N_{r}/N<2\%$) to induce the conformational changes. The results
emphasize that only a limited amount of information about pairwise distances
between amino acids is needed to induce protein conformational changes from an
initial in vitro structure to a target in vivo structure. The largest change
in pairwise separation method, which uses the target structure information,
gives the lowest values of Cα RMSD at each $N_{r}$. The normalized
$\overline{{\rm RMSD}}(N_{r})$ for the largest change in pairwise separation
method collapses for all three proteins when plotted versus $N_{r}/N$, and
reaches zero at a small fraction of restraints ($N_{r}/N<1.5\%$). The linear
discriminant analysis, which also uses the target structure information,
exhibits inconsistent performance, sometimes comparable to and sometimes
better than random selection. The normal mode analysis and the largest Cα
separation method, which do not need target structure information, give lower
values of Cα RMSD compared to random selection. The performance of these two
methods varies slightly for different proteins. Specifically, for the normal
mode analysis, $\overline{{\rm RMSD}}(N_{r})$ reaches zero when $N_{r}/N>2\%$,
which is significantly lower than the fraction of restraints used in previous
studies of FRET-assisted protein structural modeling 36.

The studies described here provide proof of principle that FRET-assisted MD
simulations will add to our understanding of the atomistic structure of
proteins in vivo since only a small fraction of restraints are required to
lock-in the target structures. However, since the in vivo structures are not
stable minima of the current MD force fields, additional restraints are
necessary to reach a given Cα RMSD and the dependence of the Cα RMSD on
$N_{r}$ is non-monotonic. Therefore, an important future direction is to
develop force fields for which the in vivo structures are potential energy
minima. A first step would be to develop force fields for mimetic systems that
stabilize the in-cell protein structure 33, 8.

The largest change in pairwise separation method for restraint selection gives
the lowest values of Cα RMSD to the target structures. However, this method
requires information about both the initial and target structures, and thus it
cannot be used to predict an in vivo structural ensemble that has not yet been
determined. Given that there are several NMR structures for proteins in
cellular environments available, this selection method can be used in
restrained MD simulations to verify that the simulations can correctly sample
the conformational dynamics for these proteins in vivo. In contrast, if we aim
to predict the in vivo structure of proteins, we have shown that normal mode
analysis can be used to select the FRET-pair labeling positions that will
enable the restrained MD simulations to sample the target in vivo structure.
In this work, which considers large-scale domain motion in proteins, we showed
that the normal mode analysis is a successful restraint selection method.
However, it is not yet clear whether such motion is relevant for proteins in
vivo. Thus, if the normal mode analysis is not efficient in inducing a change
from the initial to the target in vivo structure, the largest Cα separation
method, which gave the second best performance, can also be implemented.

It is important to note that there are several limitations for selecting amino
acid pairs for FRET labeling. For example, the labeled dye molecules should
not perturb the protein structure. Thus, in this study, we did not select
amino acid pairs that occur within the core. FRET-labeling core amino acids
would likely lead to changes in protein structure. In addition, optimal
efficiency for FRET-labeling is achieved when the distance $d$ between donor
and acceptor molecules is comparable to the characteristic Förster distance
$R_{0}$ for each dye pair. For most common FRET chromophores, $50$ Å
$<R_{0}<60$ Å, which is larger than $R_{g}$ for the proteins we considered in
this work. Thus, we excluded amino acid pairs that are separated by $d<R_{g}$.

As mentioned above, we showed in this study that the normal mode analysis
selection method (and largest Cα separation method) can identify a small
number of important Cα-Cα distance restraints that can effectively move an
initial in vitro structure to a target in vivo structure. However, since the
FRET energy transfer efficiency reports on the ensemble of distances between
dye molecules (not Cα-Cα distances), quantitative FRET-assisted protein
structural modeling requires the mapping between dye-dye distances and Cα-Cα
distances 29. Therefore, in future studies, we will identify the mapping by
incorporating atomic-scale modeling of the dye molecules into the restrained
MD simulations.

Acknowledgements: The authors acknowledge support from NIH Grant No.
R35GM151146 (E. K. and C. D.) and NIH Training Grant Nos. T32GM145452 (Z. L.,
A. G., and C. O.) and T15LM007056-37 (J. S.), as well as the High Performance
Computing facilities operated by Yale’s Center for Research Computing.

Conflict of interest statement: The authors declare no conflicts of interests.

Supporting Information: Additional supporting information can be found online
in the Supporting Information section at the end of this article.

## References

  * 1 Smyth M, Martin J (2000) X-ray crystallography. Mol Pathol 53:8.
  * 2 Williamson MP, Havel TF, Wüthrich K (1985) Solution conformation of proteinase inhibitor IIA from bull seminal plasma by 1H nuclear magnetic resonance and distance geometry. J Mol Biol 182:295–315.
  * 3 Carroni M, Saibil HR (2016) Cryo electron microscopy to determine the structure of macromolecular complexes. Methods 95:78–85.
  * 4 Fulton AB (1982) How crowded is the cytoplasm? Cell 30:345–347.
  * 5 Zimmerman SB, Trach SO (1991) Estimation of macromolecule concentrations and excluded volume effects for the cytoplasm of Escherichia coli. J Mol Biol 222:599–620.
  * 6 Ellis RJ (2001) Macromolecular crowding: Obvious but underappreciated. Trends Biochem Sci 26:597–604.
  * 7 Ebbinghaus S, Dhar A, McDonald JD, Gruebele M (2010) Protein folding stability and dynamics imaged in a living cell. Nat Methods 7:319–323.
  * 8 Davis CM, Deutsch J, Gruebele M (2020) An in vitro mimic of in-cell solvation for protein folding studies. Protein Sci 29:1046–1054.
  * 9 Leeb S, Sörensen T, Yang F, Mu X, Oliveberg M, Danielsson J (2020) Diffusive protein interactions in human versus bacterial cells. Curr Res Struct Biol 2:68–78.
  * 10 Wang Y, Sukenik S, Davis CM, Gruebele M (2018) Cell volume controls protein stability and compactness of the unfolded state. J Phys Chem B 122:11762–11770.
  * 11 Berman HM, Westbrook J, Feng Z, Gilliland G, Bhat TN, Weissig H, Shindyalov IN, Bourne PE (2000) The protein data bank. Nucleic Acids Res 28:235–242.
  * 12 Jumper J, Evans R, Pritzel A, Green T, Figurnov M, Ronneberger O, Tunyasuvunakool K, Bates R, Žídek A, Potapenko A, Bridgland A, Meyer C, Kohl SAA, Ballard AJ, Cowie A, Romera-Paredes B, Nikolov S, Jain R, Adler J, Back T, Petersen S, Reiman D, Clancy E, Zielinski M, Steinegger M, Pacholska M, Berghammer T, Bodenstein S, Silver D, Vinyals O, Senior AW, Kavukcuoglu K, Kohli P, Hassabis D (2021) Highly accurate protein structure prediction with Alphafold. Nature 596:583–589.
  * 13 Humphreys IR, Pei J, Baek M, Krishnakumar A, Anishchenko I, Ovchinnikov S, Zhang J, Ness TJ, Banjade S, Bagde SR, Stancheva VG, Li XH, Liu K, Zheng Z, Barrero DJ, Roy U, Kuper J, Fernández IS, Szakal B, Branzei D, Rizo J, Kisker C, Greene EC, Biggins S, Keeney S, Miller EA, Fromme JC, Hendrickson TL, Cong Q, Baker D (2021) Computed structures of core eukaryotic protein complexes. Science 374(6573):eabm4805.
  * 14 Sakakibara D, Sasaki A, Ikeya T, Hamatsu J, Hanashima T, Mishima M, Yoshimasu M, Hayashi N, Mikawa T, Wälchli M, Smith BO, Shirakawa M, Güntert P, Ito Y (2009) Protein structure determination in living cells by in-cell NMR spectroscopy. Nature 458:102–105.
  * 15 Tanaka T, Ikeya T, Kamoshida H, Suemoto Y, Mishima M, Shirakawa M, Güntert P, Ito Y (2019) High-resolution protein 3D structure determination in living eukaryotic cells. Angew Chem Int Ed Engl 131:7362–7366.
  * 16 Gerez JA, Prymaczok NC, Kadavath H, Ghosh D, Bütikofer M, Fleischmann Y, Güntert P, Riek R (2022) Protein structure determination in human cells by in-cell NMR and a reporter system to optimize protein delivery or transexpression. Commun Biol 5:1322.
  * 17 Ikeya T, Güntert P, Ito Y (2019) Protein structure determination in living cells. Int J Mol Sci 20:2442.
  * 18 Serber Z, Dötsch V (2001) In-cell NMR spectroscopy. Biochemistry 40:14317–14323.
  * 19 Luchinat E, Banci L (2023) In-cell NMR: Recent progresses and future challenges. Rend Fis Acc Lincei 34:653–661.
  * 20 Cheng J, Liu T, You X, Zhang F, Sui SF, Wan X, Zhang X (2023) Determining protein structures in cellular lamella at pseudo-atomic resolution by GisSPA. Nat Commun 14:1282.
  * 21 Dunstone MA, de Marco A (2017) Cryo-electron tomography: An ideal method to study membrane-associated proteins. Philos Trans R Soc B Biol Sci 372:20160210.
  * 22 Hylton RK, Swulius MT (2021) Challenges and triumphs in cryo-electron tomography. iScience 24.
  * 23 Rickard M, Zhang Y, Pogorelov T, Gruebele M (2020) Crowding, sticking, and partial folding of GTT WW domain in a small cytoplasm model. J Phys Chem B 124:4732–4740.
  * 24 Stevens JA, Grünewald F, van Tilburg PM, König M, Gilbert BR, Brier TA, Thornburg ZR, Luthey-Schulten Z, Marrink SJ (2023) Molecular dynamics simulation of an entire cell. Front Chem 11:1106495.
  * 25 Tucker MR, Piana S, Tan D, LeVine MV, Shaw DE (2022) Development of force field parameters for the simulation of single- and double-stranded DNA molecules and DNA-protein complexes. J Phys Chem B 126:4442-4457.
  * 26 Love O, Galindo-Murillo R, Zgarbová M, Šponer J, Jurečka P, Cheatham TE (2023) Assessing the current state of Amber force field modifications for DNA-2023 Edition. J Chem Theory Comput 19:4299-4307.
  * 27 Kalinin S, Peulen T, Sindbert S, Rothwell PJ, Berger S, Restle T, Goody RS, Gohlke H, Seidel CA (2012) A toolkit and benchmark study for FRET-restrained high-precision structural modeling. Nat Methods 9:1218–1225.
  * 28 Dimura M, Peulen TO, Hanke CA, Prakash A, Gohlke H, Seidel CA (2016) Quantitative FRET studies and integrative modeling unravel the structure and dynamics of biomolecular systems. Curr Opin Struct Biol 40:163–185.
  * 29 Klose D, Holla A, Gmeiner C, Nettels D, Ritsch I, Bross N, Yulikov M, Allain FHT, Schuler B, Jeschke G (2021) Resolving distance variations by single-molecule FRET and EPR spectroscopy using rotamer libraries. Biophys J 20:4842–4858.
  * 30 Hellenkamp B, Schmid S, Doroshenko O, Opanasyuk O, Kühnemuth R, Rezaei Adariani S, Ambrose B, Aznauryan M, Barth A, Birkedal V, Bowen ME, Chen H, Cordes T, Eilert T, Fijen C, Gebhardt C, Götz M, Gouridis G, Gratton E, Ha T, Hao P, Hanke CA, Hartmann A, Hendrix J, Hildebrandt LL, Hirschfeld V, Hohlbein J, Hua B, Hübner CG, Kallis E, Kapanidis AN, Kim JY, Krainer G, Lamb DC, Lee NK, Lemke EA, Levesque B, Levitus M, McCann JJ, Naredi-Rainer N, Nettels D, Ngo T, Qiu R, Robb NC, Röcker C, Sanabria H, Schlierf M, Schröder T, Schuler B, Seidel H, Streit L, Thurn J, Tinnefeld P, Tyagi S, Vandenberk N, Vera AM, Weninger KR, Wünsch B, Yanez-Orozco IS, Michaelis J, Seidel CAM, Craggs TD, Hugel T (2018) Precision and accuracy of single-molecule FRET measurements—a multi-laboratory benchmark study. Nat Methods 15:669–676.
  * 31 Agam G, Gebhardt C, Popara M, Mächtel R, Folz J, Ambrose B, Chamachi N, Chung SY, Craggs TD, de Boer M, Grohmann D, Ha T, Hartmann A, Hendrix J, Hirschfeld V, Hübner CG, Hugel T, Kammerer D, Kang HS, Kapanidis AN, Krainer G, Kramm K, Lemke EA, Lerner E, Margeat E, Martens K, Michaelis J, Mitra J, Moya Muñoz GG, Quast RB, Robb NC, Sattler M, Schlierf M, Schneider J, Schröder T, Sefer A, Tan PS, Thurn J, Tinnefeld P, van Noort J, Weiss S, Wendler N, Zijlstra N, Barth A, Seidel CAM, Lamb DC, Cordes T 2023. Reliability and accuracy of single-molecule FRET studies for characterization of structural dynamics and distances in proteins. Nat Methods 20:523–535.
  * 32 Feng R, Gruebele M, Davis CM (2019) Quantifying protein dynamics and stability in a living organism. Nat Commun 10:1179.
  * 33 Davis CM, Gruebele M (2018) Non-steric interactions predict the trend and steric interactions the offset of protein stability in cells. ChemPhysChem 19:2290–2294.
  * 34 Lindorff-Larsen K, Piana S, Dror RO, Shaw DE (2011) How fast-folding proteins fold. Science 334:517-20.
  * 35 Lindorff-Larsen K, Maragakis P, Piana S, Eastwood MP, Dror RO, Shaw DE (2012) Systematic validation of protein force fields against experimental data. PLoS One 7:e32131.
  * 36 Dimura M, Peulen TO, Sanabria H, Rodnin D, Hemmen K, Hanke CA, Seidel CA, Gohlke H (2020) Automated and optimally FRET-assisted structural modeling. Nat Commun 11:5394.
  * 37 Zhang XJ, Baase W, Matthews B (1992) Multiple alanine replacements within $\alpha$-helix 126–134 of T4 lysozyme have independent, additive effects on both structure and stability. Protein Sci 1:761–776.
  * 38 Zhang XJ, Wozniak JA, Matthews BW (1995) Protein flexibility and adaptability seen in 25 crystal forms of T4 lysozyme. J Mol Biol 250:527–552.
  * 39 Zerrad L, Merli A, Schröder GF, Varga A, Gráczer É, Pernot P, Round A, Vas M, Bowler MW (2011) A spring-loaded release mechanism regulates domain movement and catalysis in phosphoglycerate kinase. J Biol Chem 286:14040–14048.
  * 40 Lallemand P, Chaloin L, Roy B, Barman T, Bowler MW, Lionne C (2011) Interaction of human 3-phosphoglycerate kinase with its two substrates: Is substrate antagonism a kinetic advantage? J Mol Biol 409:742–757.
  * 41 Müller CW, Schulz GE (1992) Structure of the complex between adenylate kinase from Escherichia coli and the inhibitor Ap5A refined at 1.9 Å resolution: A model for a catalytic transition state. J Mol Biol 224:159–177.
  * 42 Müller C, Schlauderer G, Reinstein J, Schulz GE (1996) Adenylate kinase motions during catalysis: An energetic counterweight balancing substrate binding. Structure 4:147–156.
  * 43 Flores S, Echols N, Milburn D, Hespenheide B, Keating K, Lu J, Wells S, Yu EZ, Thorpe M, Gerstein M (2006) The database of macromolecular motions: New features added at the decade mark. Nucleic Acids Res 34:D296–D301.
  * 44 Li M, Xu W, Zhang JZ, Xia F (2014) Combined effect of confinement and affinity of crowded environment on conformation switching of adenylate kinase. J Mol Model 20:2530.
  * 45 Ikeya T, Hanashima T, Hosoya S, Shimazaki M, Ikeda S, Mishima M, Güntert P, Ito Y (2016) Improved in-cell structure determination of proteins at near-physiological concentration. Sci Rep 6:38312.
  * 46 Grigas AT, Liu Z, Regan L, O’Hern CS (2022) Core packing of well-defined X-ray and NMR structures is the same. Protein Sci 31:e4373.
  * 47 Richards FM (1974) The interpretation of protein structures: Total volume, group volume distributions and packing density. J Mol Biol 82:1–14.
  * 48 Ahmed A, Gohlke H (2006) Multiscale modeling of macromolecular conformational changes combining concepts from rigidity and elastic network theory. Proteins Struct Funct Bioinf 63:1038–1051.
  * 49 Krüger DM, Ahmed A, Gohlke H (2012) NMSim web server: Integrated approach for normal mode-based geometric simulations of biologically relevant conformational transitions in proteins. Nucleic Acids Res 40:W310–W316
  * 50 Sittel F, Stock G (2018) Perspective: Identification of collective variables and metastable states of protein dynamics. J Chem Phys 149.
  * 51 Mendels D, de Pablo JJ (2022) Collective variables for free energy surface tailoring: Understanding and modifying functionality in systems dominated by rare events. J Phys Chem Lett 13:2830–2837.
  * 52 Best RB, Hummer G (2009) Optimized molecular dynamics force fields applied to the helix-coil transition of polypeptides. J Phys Chem B 113:9004–9015.
  * 53 Lindorff-Larsen K, Piana S, Palmo K, Maragakis P, Klepeis JL, Dror RO, Shaw DE (2010) Improved side-chain torsion potentials for the Amber ff99SB protein force field. Proteins Struct Funct Bioinf 78:1950–1958.
  * 54 Abraham MJ, Murtola T, Schulz R, Páll S, Smith JC, Hess B, Lindahl E (2015) GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX 1:19–25.
  * 55 Jorgensen WL (1981) Transferable intermolecular potential functions for water, alcohols, and ethers. Application to liquid water. J Am Chem Soc 103.
  * 56 Mark P, Nilsson L (2001) Structure and dynamics of the TIP3P, SPC, and SPC/E water models at 298 K. J Phys Chem A 105:9954–9960.
  * 57 Bussi G, Donadio D, Parrinello M (2007) Canonical sampling through velocity rescaling. J Chem Phys 126:014101.
  * 58 Hoffmann A, Kane A, Nettels D, Hertzog DE, Baumgärtel P, Lengefeld J, Reichardt G, Horsley DA, Seckler R, Bakajin O, Schuler B (2007) Mapping protein collapse with single-molecule fluorescence and kinetic synchrotron radiation circular dichroism spectroscopy. Proc Natl Acad Sci 104(1):105-10.
  * 59 Speer SL, Zheng W, Jiang X, Chu IT, Guseman AJ, Liu M, Pielak GJ, Li C (2021) The intracellular environment affects protein-protein interactions Proc Natl Acad Sci 118(11):e2019918118.