Models for the binding channel of wild type and mutant transthyretin with glabridin

Abstract

Transthyretin (TTR) is a protein whose aggregation and deposition can lead to human amyloid diseases. V30A TTR is a TTR variant observed in some familial amyloidotic polyneuropathy (FAP) patients. Here we report that V30A TTR exhibits a significantly lower glabridin (Glab) binding energy than exhibited by wild-type TTR (WT TTR) in vitro. To compare changes in structural conformation induced by Glab binding to WT TTR and the V30A mutant, molecular dynamics (MD) simulations were used. Our results indicate that additional high-occupancy hydrogen bonds were observed at the binding interface between the two dimers in V30A TTR, while stabilisation hydrophobic interactions between residues in the mutant AB loop decreased. These results suggest that AC and BD contacts between the two dimers are positioned closer together in the mutant. Our results should provide useful clues to guide future TTR studies.

---

1. Introduction

Transthyretin (TTR) is one of 30 amyloidogenic proteins that can form amyloid fibrils through denaturation and misfolding.¹ TTR is mainly produced in the liver, choroid plexus, and retina.² It is one of several transporters of thyroxine (T4) and vitamin A/retinol binding protein (RBP) complexes in human plasma, as well as a primary transporter of T4 in cerebrospinal fluid.^3–5 TTR amyloid fibrils are found in patients afflicted with familial amyloidotic polyneuropathy (FAP), a fatal disease for which liver transplantation is the only available therapy.^6,7 Such fibrils are also present in patients with familial amyloid cardiomyopathy (FAC)^6,8 and senile systemic amyloidosis (SSA), a disease that affects approximately 25% of people over the age of 80.^6,9 To date, more than 100 TTR mutations have been implicated in the induction of amyloidoses such as FAP and FAC.^10–15

The 3D structure of human TTR was first determined by X-ray crystallography in 1974 by Blake et al.¹⁶ To date, more than 100 human TTRs have been solved and uploaded to the Protein Data Bank (PDB). Their 3D structures are highly similar; TTR is a homotetramer that includes four identical 14 kDa β-sandwich monomers, each comprised of 127 amino acid residues. Each monomer consists of eight strands, designated A–H, divided into two sheets (an inner sheet, DAGH, and outer sheet, BCEF) and a nine-residue α-helix. The monomers dimerize and subsequently two dimers form a tetramer through hydrogen bonding and hydrophobic interactions between their GH and AB loops (see Fig. S1†).^6,17

Glabridin (Glab; Fig. S2†) is a prenylated isoflavane originally isolated from the roots of Glycyrrhiza glabra L. (Fabaceae; licourice) and has been used for centuries as a food additive (sweetener) in cosmetics and in traditional medicine.⁶ Glab is a phytoestrogen that exhibits numerous biological properties, including antioxidant, anti-inflammatory, neuroprotective, antitumourigenic and skin-whitening activities.¹⁸ Several small molecules, including prenylated compounds, such as diflunisal, tafamidis and natural flavonoids, are known to bind to the T4-binding site and inhibit TTR amyloid fibril formation.^19–21 In 2014, the crystal structure of TTR bound to a prenylated compound, Glab (PDB ID 4N86), was reported⁶ and revealed that Glab is a potent amyloid fibril formation inhibitor with an inhibitory activity equivalent to that of diflunisal.⁶

Previous experimental and theoretical studies have focused on mechanism of formation of amyloid fibrils from TTR. For example, Zanotti et al.²² demonstrated conformational differences between I84S and I84A mutants by comparing their crystal structures obtained under various pH conditions. Moreover, the C2 space group of an L55P mutant has been shown to induce altered TTR intermolecular packing interactions.²³ Hou et al.¹⁷ demonstrated, using MD simulations, that the central channel of the L55P TTR mutant undergoes an opening-closing fluctuation, which results in a substantially lower binding energy upon binding with T4. Kelly et al.²⁴ revealed that the pathogenic mutation caused the destabilization of the dimer–dimer interactions using NMR and hydrogen-exchange. Lei et al.²⁵ performed two MD simulations to study the differences in the binding channel shapes of TTR upon binding to two inhibitors, flufenamic acid (FLU) and N-phenyl phenoxazine (BPD). These studies may explain the binding mode difference for some ligands (for example, FLU and BDP) for the WT and the mutant TTR, respectively. However, we do not know what is the binding mode of Glab? Are there any differences for other ligands binding? How can induce changes in structural conformation by Glab binding to WT TTR and the V30A mutant, an uncommon TTR mutant in FAP patients?

Here, the present study is the first to report that the V30A mutant of TTR exhibits significantly lower Glab binding energy than WT TTR in vitro. MD simulations were used to explore the characteristic binding mode of Glab and the molecular mechanism of structural stabilisation during Glab binding to V30A mutant vs. WT TTR.

2. Materials and methods

2.1 Reagents

Glab was purchased from Phytomarker Ltd. (Tianjin, China). The chemical was supplied at >98.0% purity as determined by HPLC.

2.2 rTTR expression and purification

pQE30 vectors bearing WT TTR or V30A TTR cDNA inserts were transformed into M15 competent cells (Qiagen). Proteins were expressed and purified, as previously described.^26,27

2.3 Acid induction of the formation of amyloid fibrils

+++Amyloid formation was detected using turbidity assays, which measure the degree of formation of large insoluble aggregates and amyloid fibrils, as previously described.²⁷ A 90 μL sample of rTTR (8 μM) in 10 mM phosphate buffer (pH 7.0) with 100 mM KCl and 1 mM EDTA was incubated with 10 μL of Glab dilution (stocked in DMSO at a concentration of 2.16 mM, diluted to 14.4, 36, 72, 144, 288 and 576 μM). After 2 h, the samples were diluted with 100 μL of 200 mM sodium acetate buffer containing 100 mM KCl and 1 mM EDTA at the desired pH (pH 4.2 for WT and pH 5.0 for V30A). After further incubation at 37 °C for 72 h, all samples were gently vortexed, and the optical density (OD) at 400 nm was recorded with a Synergy™ 4 Hybrid Microplate Reader (BioTekInstruments, Inc., Winooski, VT, USA). The percentage of amyloid fibril formation was calculated as follows: % fibril formation = 100 × (OD₄₀₀ of TTR at each concentration of Glab − 0% control)/(100% control − 0% control), where 100% control represents the OD₄₀₀ of TTR in the absence of inhibitors and the 0% control represents OD₄₀₀ in the absence of TTR. EC50 values were calculated from the Glab concentrations needed to achieve 50% inhibition of the production of TTR amyloid fibril, as compared with the untreated control groups. The results were plotted using GraphPad Prism 5.0. Data were presented as the means ± standard errors of at least three independent experiments.

2.4 Protocols for MD simulations

The initial structure used for MD simulations employed a published X-ray crystal structure of WT TTR bound to Glab, PDB code 4N86.⁶ The V30A mutant was then modeled using the SWISS-MODEL server.²⁸ All crystallographic water molecules were stored. Each MD simulation was carried out in water at a constant temperature of 300 K and a constant pressure of 1 bar. The protonation states of the ionisable groups were calculated using neutral pH for all simulations. All complex systems were subjected to MD simulations with periodic boundary conditions using the GROMACS 4.5.2 software package with the SPC (simple point charge) water model.²⁹ For the ligand Glab, the Dundee PRODRG-server³⁰ was used to build a GROMACS topology to generate the .itp file containing the molecular information. We added the .itp file to the protein .top file and then performed protein ligand complex simulations using the GROMOS force field force parameter set 53A6.³¹ Na⁺ counter ions were added to retain the neutrality of the entire system. The systems were placed in a cubic box (proteins were placed at least 0.8 nm from the box edge), subjected to periodic boundary conditions, and then solvated using explicit SPC216 modelled water molecules. Before MD simulations were performed, the systems were energy minimised using the steepest descent algorithm to avoid any steric conflicts generated during the initial setup. NVT (canonical ensemble) and NPT (isothermal–isobaric ensemble) equilibration for 500 ps each were performed to allow the system to reach the desired temperature and pressure. Bond lengths and angles were constrained using the P-LINCS algorithm,³² and the geometry of water molecules was constrained using the SETTLE algorithm.³³ A twin-range cutoff of 1.2 nm was used for the van der Waals (vdW) interactions, and long-range electrostatic interactions were determined using the particle mesh Ewald method.³⁴ The equilibration procedure involved solvent thermalisation with the solute atoms fixed for 500 ps at 300 K, and solute atoms were minimised with the solvent coordinates fixed. The complete system was subjected to MD simulations by increasing the temperature from 0 to 300 K in 500 ps increments of 50 K each. Data were produced by sampling for 100 ns. We collected a 100 ns trajectory for each simulation (molecular systems are shown in Table S1†).

2.5 Pathways of the WT and V30A mutant

CAVER 3.0 (ref. 35) provides detailed characteristics of the individual transport pathways and their time evolution. The pathways were identified with CAVER 3.0 as the paths in a graph composed of Voronoi vertices and Voronoi edges.³⁶ The lowest cost pathway was selected, and all pathways within a user specified distance were discarded. We used CAVER 3.0 (ref. 35) to compare the WT and V30A mutant TTR channels during 20 ns MD simulations and the RMSD (root-mean-square deviation) plots are shown in Fig. S3.†

2.6 Principal component analysis (PCA)

PCA is a widely used approach for extracting the slow and functional motions of biomolecules from MD trajectories by applying a dimensional reduction method.³⁷ PCA is based on the calculation and diagonalisation of the covariance matrix (C_ij) of the fluctuations of each of the x, y and z coordinates of the C_α atoms from their average motions over 100 ns in the simulations for the two models. Using the displacement vectors Δr_i and Δr_j of atoms i and j, respectively, C_ij was calculated as follows:
where Δr_i (or Δr_j) is the displacement vector corresponding to the ith (or jth) atom of the systems. The eigenvectors of the matrix are also called principal components (PCs), which represent the directions of the concerted motions. The first few PCs describe the slow-motion modes of the system; these modes are related to the functional motions of a biomolecular system.³⁷ The eigenvalues of the matrix indicate the magnitude of the motions along the directions of motion. In this study, PCA was performed to investigate and compare the motion modes of the two systems using GROMACS 4.5.2.

2.7 Molecular mechanical/PBSA calculations

The two MDs were used for two protein–inhibitor complex analysis. The 2000 snapshots isolated from the final 4000 ps MD trajectory for the protein–inhibitor complex were used for the binding free energy calculation using the MM-PBSA method encoded by the Amber 10 software.³⁸ The binding free energy calculation was carried out using the molecular mechanical (MM)/PBSA method.³⁸ The method combines the MM energies with continuum solvent modelling approaches using Amber 10 modules; the Sander module is used to calculate the MM energies, the PBSA module³⁹ is used to calculate the polar solvent energy (the electrostatic contribution) and the Molsurf module⁴⁰ is used to calculate the nonpolar solvent energy.

3. Results and discussion

3.1 Effect of Glab on amyloid fibril formation of TTR proteins

The effects of Glab binding on amyloid fibril formation for WT TTR and the V30A TTR mutant as a function of pH were monitored by measuring the turbidity at 400 nm. The TTR concentration of 3.6 μM was selected based on the average physiologic concentration, whereas Glab was evaluated at concentrations ranging from 0 to 14.4 μM. The optimal pH was selected for each TTR to produce maximal amyloid fibril formation. As shown in Fig. 1a, Glab could suppress amyloid fibril formation of WT TTR at pH 4.2 and of V30A at pH 5.0 in a dose-dependent manner. In addition, Glab stabilised V30A TTR more than WT TTR, as reflected by the lower EC50 of mutant relative to WT TTR (Fig. 1b). However, it was still not known how this mutation could induce a conformation change in TTR. In order to explore the binding mechanics of Glab and explain differences in binding between WT and mutant TTR, further studies are described below.

Fig. 1 Glab stabilizes TTR to inhibit pH-induced amyloid fibril formation. The TTRs under selected pH conditions were incubated with various concentrations of Glab at 37 °C for 72 h. The amyloid fibril formation percentage (a) and EC50 values (b) were calculated from the OD₄₀₀ readings the data were expressed as the mean value ± SD for triplicate experiments. (***indicates p < 0.01).

3.2 Identification of pathways in MD trajectories for WT and V30A mutant TTRs

Channels mediate the transport of ions or molecules across biological membranes.³⁵ Therefore, characterising the central binding channel of TTR is important in designing a highly efficient inhibitor. TTR contains two binding sites (Fig. 2), each of which is divided into three structurally distinct regions (two stretched four-β-stranded antiparallel sheets and a short α-helix, the EF-helix, situated between the E and F strands of TTR).¹⁷ At the channel entrance, Glu54, His56 and Lys15 comprise a charged region. Along the binding site, the methyl groups of Leu17, Thr106, Ala108, Leu110 and Val121 form the hydrophobic pocket (Fig. S4a†).¹⁷

Fig. 2 The quaternary structure of TTR bound to Glab. TTR is represented using a cartoon model, and Glab molecules are represented by a stick model (optimised using Gaussian 09 software with the B3LYP 6/31G* basis set⁴¹). Subunit (A) is cyan, (B) is purple, (C) is orange and (D) is green.

CAVER 3.0 was then used to analyse the 400 snapshots from two 20 ns MD simulations of WT and V30A TTR. Although Val30 was located away from the binding site, replacement of Val30 with alanine in the mutant channel resulted in a narrower channel, changing the bottleneck radius from 1.95 Å to 1.65 Å calculated by CAVER 3.0 (Fig. S4b and c,† Table 1) and 1.10 Å to 0.20 Å calculated by CHEXVIS (Table 2). Differences between WT and V30A mutant conformations of the channel possibly contribute to the significantly lower binding energy observed for V30A than for WT and corresponds to a decrease in bottleneck radius in the mutant as well (Tables 1 and 2). This minor but critical change was captured by CHEXVIS visualisations (Fig. 3 and Table 3).⁴² As shown in Fig. 3, hydrophobicity and polarity exhibited a clear change and few polarities were observed in the channel of the mutated structure, which was critical for the stability of the homotetramer channel.

Table 1 The channel bottleneck radii for WT and V30A TTRs, calculated by CAVER 3.0^a

Protein	Avg_BR	Max_BR	Avg_L	Avg_C
a Avg_BR: average bottleneck radius [A]; Max_BR: maximum bottleneck radius [A]; Avg_L: average channel length [A]; Avg_C: average channel curvature.
WT	1.95	2.82	9.68	1.15
V30A	1.65	2.68	8.35	1.23

---

Table 2 The channel bottleneck radii for WT and V30A TTRs, calculated by CHEXVIS

Protein	Length	Bottle-neck	Straightness
WT	18.46	1.10	0.77
V30A	24.10	0.20	0.68

---

Fig. 3 Channel profile visualization of WT TTR and V30A TTR, respectively. Hydrophobicity and electrostatic potential profiles are shown in a split visualization.

Table 3 Summary of hydrogen bonds formed between TTR and V30A mutant monomers/dimers during 100 ns MD simulations

Protein				WT	MUT
A/B	F strand	Val94CO⋯HNGlu89	α helix-Floop	<10	15
	H strand	Ser115B1H⋯OCThr118	H strand	<10	10
C/D	G–Hloop	Thr118NH⋯OCTyr116	H strand	<10	10
	H strand	Tyr114CO⋯HAThr119	H strand	<10	18
	H strand	Thr119BH⋯OGSer115	H strand	<10	26
	H strand	Thr118CO⋯HASer115	H strand	11	52
B/D	H strand	Val122NH⋯OCGly22	A–B loop	16	28

---

3.3 Structural fluctuations induce a conformational change in the Glab-binding channel

Subsequently, we calculated the RMSD of the main chain atoms relative to the corresponding initial structure of free WT TTR, TTR–Glab and V30A TTR–Glab as a function of simulation time. The trend of the RMSDs of the two protein–ligand complexes showed reasonable plateaus during the 100 ns dynamic process (Fig. 4a–d), indicating that statistical convergence was attained in these simulations. Meanwhile, the WT TTR RMSD displayed a more drastic change for unbound A, B, C and D chains than observed for the V30A mutant, suggesting that the homotetramer structure might undergo a larger conformational change while binding Glab to form the WT TTR–Glab complex. Fig. 5a–d display the distributional probability of RMSD values from 100 ns trajectories for the two systems in the A, B, C and D chains, respectively. The mean RMSD values of the WT system four chains were 0.20, 0.21, 0.23 and 0.21 nm, respectively, whereas the mean values for the V30A system chains were 0.18, 0.13, 0.19 and 0.18 nm. Therefore, the mean RMSD values for the four chains of the WT system were all higher than for the V30A system, suggesting that the V30A system is more stable and that RMSD analysis of structural stabilities of MD simulations agree with the experimental results (Fig. 1).

Fig. 4 The RMSDs of main chains for TTR–Glab and V30A TTR–Glab in A, B, C and D chains relative to the corresponding X-ray structures as a function of simulation time.

Fig. 5 The probability distribution of RMSD of main chains for TTR–Glab (black) and V30A TTR–Glab (red) in (a) chain A, (b) chain B, (c) chain C, (d) chain D calculated from 100 ns trajectories.

In consideration of the calculated RMSDs above, root-mean-square fluctuation (RMSF) analyses demonstrated the plausibility of the generalised force fields modelled for Glab. Fig. 6a–d show the RMSF values calculated for 100 ns trajectories for TTR–Glab and V30A TTR–Glab. For the two systems, residues R21 (AB loop), D38 (BC loop), E51 (CD loop), E62 (DE loop) and S100 (FG loop) in the A, B, C and D chains in the V30A TTR mutant exhibited larger fluctuation values than did their WT counterparts. Notably, the V30A mutation correlated with the largest RMSF changes observed for both WT and mutant TTRs studied here, presumably because the AB (residues 19–28) and GH loops (residues 113–114) are located at the dimer–dimer interface. Because the fluctuation at the R21 position was relatively larger for the mutant than the WT TTR, there may be greater flexibility in the region surrounding R21 in the mutant. Thus, the binding of Glab may strengthen dimer–dimer interactions by stabilizing AB loop interaction between dimers. The change of secondary structure of AB loop and CD loop was used (see Fig. S5a and b†), which indicated that there were no drastic changes of secondary structure of AB loop and CD loop during 100 ns.

Fig. 6 Plots of the RMSFs compared with thermal B-factors per residue in the A, B, C and D chains for the WT TTR (black) and V30A mutant TTR (red).

Large thermal B-factor values of residues in the BC loop region (residues 35–41) were observed for crystal structures of both WT TTR and V30M TTR.⁴³ In the present study, D38 (in the BC loop) correlated with a B-factor peak in MD simulations and exhibited larger fluctuations in mutant than WT TTR (Fig. 6e–h). The results suggest that the mutant residue A30 shifts the BC loop position toward the monomer–monomer interface, bringing the adjacent outer sheets in closer proximity to one another.

The backbone dihedral fluctuations (BDFs) of residues R21, D38, E51, E62, S100, N98, L17 and L110 discussed above are shown in Fig. 7a–d. The large fluctuations in dihedral angles involved mainly D38 in the BC loop and N98 and S100 in the FG loop. This result is consistent with the RMSF and B-factor results. The large changes observed in dihedral angles of the main chain between WT and mutant TTRs predominantly involve D38, N98 and S100 residues. Stabilisation of dimer–dimer interactions involving these residues is disrupted in V30A TTR.

Fig. 7 The BDFs of WT and V30A TTR in the A, B, C and D chains. Black and red lines represent WT and V30A TTR, respectively.

We calculated the single residue Solvent-Accessible Surface Area (SASA) of the total surface area of the channel along the course of our simulations to reveal the influence of the single-point mutation on the hydrophobicity of the central channel. As shown in Fig. 8a–d, residues R21, D38, D39, E62, E63, D99, N98, S100, L17, L110, K15, A108 and T119 in the V30A mutant exhibited pronounced changes relative to WT for each of these quantities. This result indicates that the shape change of the channel of V30A TTR may undergo substantial fluctuation. Notably, R21 is located in the AB loop (within the dimer–dimer interface in the AC contact), and may enhance the dimer–dimer interaction of V30A TTR. The hydrogen bond networks in the monomer–monomer and dimer–dimer interfaces are also presented in the study. The MD simulation results show that the hydrogen bonds in TTR greatly contribute to better stabilisation of WT TTR vs. V30A TTR, as shown in Table 3. Among these hydrogen bonding residues, Y114 was found in the GH loop and the dimer–dimer interface. Therefore, V30 is located in an intrinsically unstable site in TTR, and mutation of this reside to alanine leads to substantial local and global structural changes.

Fig. 8 Time-dependent SASA of the residues for WT TTR and V30A TTR in A, B, C and D chains. Black and red columns present WT TTR and V30A TTR, respectively.

Essential dynamics analysis was applied to the MD trajectories to isolate low-frequency motions. In principle, this method allows for the extrapolation of motions in the direction of the selected eigenvectors. In this study, this analysis was used to describe the low-frequency persistent motions observed in the molecular trajectories. Two principal measurements of motion, PC-1 and PC-2, were obtained (Fig. S6†). Fig. S6a–d† show the lowest relative free energy change between WT TTR and WT TTR–Glab. Notably, V30A–Glab TTR exhibited a more compact D chain conformation than did WT TTR–Glab. This result indicates that V30A TTR conformation was more dramatically affected by the binding of one Glab ligand than was WT TTR.

We performed a longer MD simulation of up to 100 ns to examine the dynamic cross-correlation matrices (DCCMs) for both WT and V30A TTR–Glab and thus investigated the effect of one Glab on correlated motions of amino acids in both TTRs (Fig. 9a and b). DCCM is a 3D matrix representation that graphically displays time-correlated information for the atoms within the residues of proteins.⁴⁴ On the DCCM map, each point represents a correlation C_ij of atoms i and j. If C_ij = 1, the fluctuations of atoms i and j are completely correlated (same period and same phase). If C_ij = −1, the fluctuations of atoms i and j are completely anticorrelated (same period and opposite phase). If C_ij = 0, the fluctuations of atoms i and j are not correlated. The DCCMs in Fig. 9a and b show that the global dynamics of WT TTR–Glab and V30A TTR–Glab were highly similar to each other; Glab binding did not change the original overall secondary structures of the proteins. However, major differences existed in some regions (refer to the correlations between WT TTR–Glab in the B chain and V30A TTR–Glab in the D chain in the circled areas in Fig. 9c and d). These differences might reflect how Glab changed the local dynamics and conformational subspace of WT and V30A proteins. For WT TTR, the B chain displayed strong anticorrelated movements upon Glab binding, as indicated by the red colour (Fig. 9a and c). These anticorrelated movements disappeared in the B chain but did not disappear in the D chain of V30A TTR (Fig. 9b and d), indicating a change in the interaction network of the protein. For example, when Glab was bound, the original interaction network existing between protein residues and water molecules broke and a new interaction network formed between the protein residues, which strengthened the dimer–dimer interaction. Changes in these fluctuations in protein dynamics contributed to the entropy of the system and DCCMs demonstrated large-scale correlated motions for both WT TTR and V30A TTR protein systems. The observed difference in DCCMs between WT TTR–Glab and V30A TTR–Glab indicate that their entropic contribution to thermodynamic stability also differed (V30A TTR–Glab was more stable than WT TTR–Glab).

Fig. 9 (a) Cross-correlation matrix of the fluctuations of each of the x, y, and z coordinates of the C_α atoms of the residue from their average during 100 ns MD of TTR–Glab, (b) V30A–Glab, (c) TTR–Glab B chain residues, (d) TTR–Glab D chain residues.

We also calculated the Cα distances between residue 30 (in the B strand) and its native hydrogen bonding partner, Gly47 (in the C strand) (Fig. 10a–d). The average distances between Val30 and Gly47 in the WT TTR–Glab complexes were consistent with structural predictions for V30A TTR A, B and C chains. However, the average C_α distance between residue 30 and Gly47 resulted in large separations in the D chain in V30A TTR, which disrupted the hydrogen bonds in the mutant. This observation further demonstrates that the B strand structure at V30 is intrinsically unstable in WT TTR.

Fig. 10 The C_α distances between residue 30 and its native hydrogen bonding partner, Gly47, in the C strand are shown as a function of time in the A, B, C and D chains.

3.4 MM-PBSA calculations

The binding energies of the five systems were calculated using the g_mmpbsa package and the results are shown in Table 4. The van der Waals energies of WT TTR–Glab (−45.11 kcal mol⁻¹) were lower than those of V30A TTR–Glab (−31.10 kcal mol⁻¹), indicating that the dimers were spaced further apart from each other at the dimer interface in the mutant. Thus, the van der Waals energies between V30A TTR and Glab increased, but the electrostatic energies in the mutant system were much lower than those in the WT system. V30A TTR formed several additional hydrogen bonds at the dimer–dimer interface, resulting in more favourable electrostatic energies. However, the polar residues on the dimer–dimer interface formed many more hydrogen bonds in the V30A TTR system; hence, the polar interactions between the mutant protein and Glab may have subsequently decreased. Overall, the total binding energies between WT TTR and Glab were only slightly higher than those of the V30A TTR–Glab complex. The entropy contribution was not included in the g_mmpbsa calculations; thus, this binding energy will exhibit some deviations compared from the actual binding free energy.^45,46

Table 4 The MM-PBSA score for the two complexes (kcal mol⁻¹)

Energy components	WT	V30A TTR
ΔEvdw	−45.11	−31.10
ΔEele	−5.14	−22.17
ΔEMM = ΔEvdw + ΔEele	−67.28	−36.24
ΔGpolar,sol	48.93	20.53
ΔGnonpolar,sol	−5.93	−4.58
ΔGsol	43	15.95
ΔGpolar,sol + ΔEele	26.76	15.39
ΔGnonpolar,sol + ΔEvdw	−51.04	−35.68
ΔGtotal	−20.28	−24.27

---

Overall, the V30A mutation mainly disrupted the dimer–dimer interaction. High occupancy hydrogen bonds present in the interface between the two dimers, as well as the hydrophobic interaction of residue R21 in the AB loop, were decreased in the mutant, indicating that the two dimers contact each other more closely in WT TTR. Therefore, the V30A mutation decreased the binding energy between TTR and Glab relative to that observed for WT TTR.

4. Conclusions

TTR is a protein that can participate in development of amyloid diseases in humans. The native state of TTR is a homotetramer. Each monomer is comprised of eight β-strands organised into a β-sandwich. In this study, we first found that V30A TTR, an uncommon TTR mutant observed in FAP patients, exhibited a significantly lower Glab binding energy than did WT TTR in vitro. Serval theoretical methods were used to explore the binding mode of Glab with the WT and V30A mutant, respectively. We calculated the bottleneck radius of the WT and V30A TTR, and a decrease in the mutant demonstrated replacement of Val30 with alanine in the mutant channel resulted in a narrower channel. MD simulations indicated that RMSD, RMSF, B-factor, BDF and DCCMs for both WT and V30A TTR complex elucidated the structural reorganization mechanisms for Glab binding, in which the V30A mutation led to a conformational change upon Glab binding that moved the AC and BD contacts closer together, the lower binding energy for the mutant relative to WT TTR can be rationalized. These results should further our understanding and should guide future studies regarding amyloid disease causation in humans.

Acknowledgements

This work was supported by Major scientific research projects of Jilin Province (20140203025NY), the Fundamental Research Funds for the Jilin Universities (450060481232) and was worked at the High Performance Computing Center of Jilin University.

References

1. T. Yokoyama, M. Mizuguchi, Y. Nabeshima, K. Kusaka, T. Yamada, T. Hosoya, T. Ohhara, K. Kurihara, K. Tomoyori and I. Tanaka, J. Struct. Biol., 2012, 177, 283–290.
2. J. A. Hamilton and M. D. Benson, Cell. Mol. Life Sci., 2001, 58, 1491–1521.
3. L. Bartalena and J. Robbins, Clin. Lab. Med., 1993, 13, 583–598.
4. G. Schreiber and S. J. Richardson, Comp. Biochem. Physiol. B Biochem. Mol. Biol., 1997, 116, 137–160.
5. H. L. Monaco, Biochim. Biophys. Acta. Protein Struct. Mol. Enzymol., 2000, 1482, 65–72.
6. T. Yokoyama, Y. Kosaka and M. Mizuguchi, J. Med. Chem., 2014, 57, 1090–1096.
7. M. M. Reilly and R. H. M. King, Brain Pathol., 1993, 3, 165–176.
8. D. R. Jacobson, R. D. Pastore, R. Yaghoubian, I. Kane, G. Gallo, F. S. Buck and J. N. Buxbaum, N. Engl. J. Med., 1997, 336, 466–473.
9. P. Westermark, K. Sletten, B. Johansson and G. G. Cornwell, Proc. Natl. Acad. Sci. U. S. A., 1990, 87, 2843–2845.
10. J. N. Buxbaum and C. E. Tagoe, Annu. Rev. Med., 2000, 51, 543–569.
11. A. F. Castro-Rodrigues, L. Gales, M. J. Saraiva and A. M. Damas, Acta Crystallogr., 2011, 67, 1035–1044.
12. J. Y. Liu, X. M. Jiang, M. Zhang and Y. J. Guo, Int. J. Neurosci., 2012, 122, 716–718.
13. A. R. Batista, D. Gianni, M. Ventosa, A. V. Coelho, M. R. Almeida, M. Sena-Esteves and M. J. Saraiva, Gene Ther., 2014, 21, 1065.
14. J. P. Liu, Q. Y. Wang, F. Zheng, J. H. Lu, P. Ge, Y. J. Gu and X. G. Sun, Genet. Mol. Res., 2014, 13, 2368–2376.
15. D. Santos, T. Coelho, M. Alves-Ferreira, J. Sequeiros, D. Mendonça, I. Alonso, C. Lemos and A. Sousa, Eur J Hum Genet., 2016, 24, 756–760.
16. C. C. Blake, M. J. Geisow, I. D. Swan, C. Rerat and B. Rerat, J. Mol. Biol., 1974, 88, 1–12.
17. M. Lei, M. Yang and S. Hou, J. Struct. Biol., 2004, 148, 153–168.
18. C. Simmler, G. F. Pauli and S. N. Chen, Fitoterapia, 2013, 90, 160–184.
19. S. L. Adamski-Werner, S. K. Palaninathan, J. C. Sacchettini and J. W. Kelly, J. Med. Chem., 2004, 47, 355–374.
20. C. E. Bulawa, S. Connelly, M. Devit, L. Wang, C. Weigel, J. A. Fleming, J. Packman, E. T. Powers, R. L. Wiseman, T. R. Foss, I. A. Wilson, J. W. Kelly and R. Labaudiniere, Proc. Natl. Acad. Sci. U. S. A., 2012, 109, 9629–9634.
21. D. B. Trivella, C. V. dos Reis, L. M. T. Lima, D. Foguel and I. Polikarpov, J. Struct. Biol., 2012, 180, 143–153.
22. N. Pasquato, R. Berni, C. Folli, B. Alfieri, L. Cendron and G. Zanotti, J. Mol. Biol., 2007, 366, 711–719.
23. M. P. Sebastiao, M. J. Saraiva and A. M. Damas, J. Biol. Chem., 1998, 273, 24715–24722.
24. K. Liu, J. W. Kelly and D. E. Wemmer, J. Mol. Biol., 2002, 320, 821–832.
25. H. Wang, Y. Tang and M. Lei, Arch. Biochem. Biophys., 2007, 466, 85–97.
26. K. Matsubara, M. Mizuguchi and K. Kawano, Protein Expression Purif., 2003, 30, 55–61.
27. F. Zhang, C. Hu, Y. Dong, M. S. Lin, J. Liu, X. Jiang, Y. Ge and Y. Guo, Arch. Biochem. Biophys., 2013, 535, 120–127.
28. M. Biasini, S. Bienert, A. Waterhouse, K. Arnold, G. Studer, T. Schmidt, F. Kiefer, T. G. Cassarino, M. Bertoni, L. Bordoli and T. Schwede, Nucleic Acids Res., 2014, 42, W252–W258.
29. D. Wang, H. Jin, J. Wang, S. Guan, Z. Zhang and W. Han, J. Biomol. Struct. Dyn., 2016, 34, 749–761.
30. A. W. Schuttelkopf and D. M. F. van Aalten, Acta Crystallogr., Sect. D: Biol. Crystallogr., 2004, 60, 1355–1363.
31. H. Jin, Z. Zhou, D. Wang, S. Guan and W. Han, Int. J. Mol. Sci., 2015, 16, 6217–6234.
32. B. Hess, H. Bekker, H. J. C. Berendsen and J. G. E. M. Fraaije, J. Comput. Chem., 1997, 18, 1463–1472.
33. S. Miyamoto and P. A. Kollman, J. Comput. Chem., 1992, 13, 952–962.
34. G. M. Torrie and J. P. Valleau, Chem. Phys. Lett., 1974, 28, 578–581.
35. E. Chovancova, A. Pavelka, P. Benes, O. Strnad, J. Brezovsky, B. Kozlikova, A. Gora, V. Sustr, M. Klvana, P. Medek, L. Biedermannova, J. Sochor and J. Damborsky, PLoS Comput. Biol., 2012, 8, e1002708.
36. F. Aurenhammer, ACM Comput. Surv., 1991, 23, 345–405.
37. M. Laberge and T. Yonetani, Biophys. J., 2008, 94, 2737–2751.
38. D. A. Case, T. A. Darden, T. E. Cheatham III, C. L. Simmerling, J. Wang, R. E. Duke, R. Luo, M. Crowley, W. S. Ross, W. Zhang, K. M. Merz, B. Wang, S. Hayik, A. Roitberg, G. Seabra, I. Kolossváry, K. F. Wong, F. Paesani, J. Vanicek, X. Wu, S. R. Brozell, T. Steinbrecher, H. Gohlke, L. Yang, C. Tan, J. Mongan, V. Hornak, G. Cui, D. H. Mathews, M. G. Seetin, C. Sagui, V. Babin and P. A. Kollman, AMBER 10, University of California, San Francisco, 2008.
39. P. A. Kollman, I. Massova, C. Reyes, B. Kuhn, S. Huo, L. Chong, M. Lee, T. Lee, Y. Duan, W. Wang, O. Donini, P. Cieplak, J. Srinivasan, D. A. Case and T. E. Cheatham, Acc. Chem. Res., 2000, 33, 889–897.
40. D. Sitkoff, K. A. Sharp and B. Honig, J. Phys. Chem., 1994, 98, 1978–1988.
41. J. Z. Vilseck, J. Kostal, J. Tirado-Rives and W. L. Jorgensen, J. Comput. Chem., 2015, 36, 2064–2074.
42. T. B. Masood, S. Sandhya, N. Chandra and V. Natarajan, BMC Bioinf., 2015, 16, 119.
43. C. J. Terry, A. M. Damas, P. Oliveira, M. J. Saraiva, I. L. Alves, P. P. Costa, P. M. Matias, Y. Sakaki and C. C. Blake, EMBO J., 1993, 12, 735–741.
44. X. Y. Wang, C. G. Ji and J. Z. H. Zhang, J. Phys. Chem. B, 2012, 116, 11570–11577.
45. A. Grossfield, Biochim. Biophys. Acta, 2011, 1808, 1868–1878.
46. F. Wang, H. Wan, J. P. Hu and S. Chang, Mol. BioSyst., 2015, 11, 223–231.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c6ra19814g

---