Molecular Dynamics Simulations of Mite Aquaporin DerfAQP1 from the Dust Mite Dermatophagoides farinae (Acariformes: Pyroglyphidae)

Aquaporins are a large family of transmembrane channel proteins that facilitate the passive but highly selective transport of water and other small solutes across biological membranes. House dust mite (Dermatophagoides farinae) is the major source of household immunogens, and we have recently reported six cDNA sequence encoding aquaporins from this mite species. To better understand the structure and role of mite aquaporin, we constructed a tertiary structure for DerfAQP1 by homology modeling from the X-ray structure of malaria aquaporin PfAQP (Protein Data Bank code No. 3C02) and conducted molecular dynamics simulation. The simulation arranged seven water molecules in a single file through the pores of the DerfAQP1. Further, two conserved Asn-Pro-Ala motifs were located on Asn203 and Asn77; residues Arg206, Trp57, Met190, Gly200, and Asp207 constituted an extracellular vestibule of the pore; and residues His75, Val80, Ile65, and Ile182 constituted the cytoplasmic portions. The overall free energy profile for water transport through DerfAQP1 revealed an energy barrier of ~2.5 kcal/mol. These results contribute to the understanding of mite physiology and pathology.


Introduction
House dust mites are a worldwide health problem and are a major risk factor for asthma, rhinitis, dermatitis, and other allergic disorders. More than 15%-20% of the populations in industrialized countries are affected by house dust mite sensitization [1]. The two most important elicitors of mite allergies are house dust mites Dermatophagoides pteronyssinus and D. farinae, which are both members of the Pyroglyphidae family. The prevalence of these two major species of mites depends on the geographical region, with D. pteronyssinus as the most prevalent species in Europe, North and South America, Africa, and Australia. Meanwhile, D. farinae is often misleadingly described as the American house dust mite, although it is dominant only in the northeastern parts of North America [2]. Chemicals have been tested in the laboratory for their ability as acaricides to kill house dust mites, but most have poor control or 100% mortality at very high concentrations, thus posing a threat to the environment and humans.
Aquaporins (AQPs) are a family of small, integral membrane proteins that play critical roles in facilitating water transport across cellular plasma membranes in response to osmotic gradients. AQP performs fundamental tasks in the physiology of virtually every living organism, including higher mammals, other vertebrates, invertebrates, plants, eubacteria, archaebacteria, and other microbes. Proteins are involved in diverse biological processes throughout the natural world, making them attractive targets for the development of novel drug therapies or inhibitors of membrane water permeability, such as mercury, gold, silver, copper, phloretin, tetraethyl ammonium salts, and acetazolamide compounds [3]. Thus, studies of mite aquaporins allow both precise biophysical measurements of their functions and exploration of high-affinity inhibitors for control of these insects.
The AQP primary sequences contain six membranespanning segments (TM1-6) with five connecting loops (A-E). Each AQP monomer consists of an extracellular and cytoplasmic vestibule connected by a central amphipathic pore containing hydrophilic and hydrophobic surfaces. The hydrophilic surface comprises alpha-carbonyl groups from the polypeptide backbone, including NPA motifs, and the hydrophobic residues of the pore prevent the permeation of larger molecules. The extracellular vestibule contains an ar/R constriction region or selectivity filter (SF), which comprises aromatic and arginine residues and is also the narrowest region of the pore [7].
Three-dimensional (3D) structural research provides further insights into the atomic-level mechanisms of water and solute conduction and of proton and/or ion exclusion. However, the AQP structure is classically difficult to characterize because purification and crystallization are impeded in a membrane environment. Further, no current experimental methods can monitor permeation through AQPs on a molecular level with sufficient spatial and temporal resolution. Molecular dynamics simulations (MDS) can go beyond static structures, thus offering an ideal method to understand the working mechanism of water and glycerol transport through AQPs at atomic resolution. On the basis of our reports on human AQP4 [8], PfAQP from the malarial parasite Plasmodium falciparum [9], and Aqy1 from Pichia pastoris [10], the setup for these simulations typically involves a box containing the protein embedded in a simple membrane surrounded with water.
Elucidating the structure of AQPs provides a unique opportunity to understand the mechanisms of these biomolecular systems. By using RNA sequencing and RT-PCR, we identified six AQPs from the common house dust mite Dermatophagoides farinae, namely, DerfAQP1-4, DerfAQP5.01, and DerfAQP5.02 (GenBank Nos. KY231248, KY231249, KY231250, KY231251, KY231252, and KY231253) [11]. To better understand its structure and role in mite physiology and pathology, we constructed the tertiary structure and conducted MDS for DerfAQP1.

System Preparation.
Templates of the DerfAQP1 (GenBank Accession No. ATL75744.1) protein were analyzed for homology modeling in the Protein Data Bank (PDB) [12] by using Protein BLAST. The appropriate templates were selected based on high score, low e value, and maximum sequence identity. The tertiary structure of DerfAQP1 was constructed with Modeller 9.0 (https://salilab .org/modeller). PROCHECK [13] was employed to determine the stereochemical quality of the DerfAQP1 protein structure by analyzing the residue-by-residue geometry and overall structural geometry. VERIFY 3D [14] was used to detect the compatibility of an atomic model (3D) of DerfAQP1 protein with its own amino acid sequence (1D). ERRAT [15] was used to analyze the statistics of nonbonded interactions between different atomic types. QMEAN [16] was used to derive both global (the entire structure) and local (per residue) absolute quality estimates based on a single model. The GMQE [17] score was expressed as a number between 0 and 1, which reflects the expected accuracy of a model built with that alignment and template, and the coverage of the target wherein higher numbers indicate higher reliability. The interactive visualization and analysis of molecular structures was completed using UCSF Chimera 1.13.1 [18].
We used the automatic protein structure file (PSF) builder plugin of the VMD software to generate the PSF and PDB files for the monomeric DerfAQP1 and crystallographic water molecules and ions. The protonation states of all residues were selected to conform to the physiological pH range of 7.4 according to the pKa values of their side chains using PSF builder. To solvate the protein, we used the SOLVATE program [19] to fill any empty space inside the pores and to surround the protein with water. By using the VMD software and our previous scripts [20], the allatom model system for the DerfAQP1 monomer was formed approximately along the z-axis [21] at the center of an 80 × 80 × 100 Å 3 water box. We then added 52 Na + and 55 Clions to neutralize and salinize the system to 150 mM, with a total of 59,914 atoms.

DerfAQP1
Simulations. Simulations were performed using the Chemistry at Harvard Macromolecular Mechanics 36 (CHARMM36) force field [21,22], TIP3P water model [23], and molecular dynamics program NAMD version 2.9 software [24] with periodic boundary conditions at a constant temperature of 298 K and constant pressure of 1 atm (NPT ensemble) to conform to standard laboratory conditions. The vDW interactions were calculated with a cutoff distance of 1.0 nm, switching distance of 0.9 nm, and pair-list distance of 1.2 nm. The Langevin dynamics and 2 BioMed Research International piston methods [25] were used to maintain constant temperature and pressure. Full electrostatics were employed with the particle mesh Ewald method [26]. Time steps of 2.0 and 1.0 fs were used for long-range electrostatic and for short-range nonbonded and bonded interactions, respectively. The energy of the bulk water molecules of the system was minimized, while protein and ions were fixed. Next, we fixed alpha-carbons on the central parts of the helices ð−10 Å < z < 10 ÅÞ during the equilibrium molecular dynamics run to preserve the backbone structure [20].

Computation of Water Density in Each Interval through
DerfAQP1. To compute the density of water ρðz i Þ at every interval ½z i , z i + Δz along the channel, we first obtained a time average of the amount of water molecules within that interval using the following equation: where t f = 20 ns is the duration of the simulation and Δt = 1 ps is the time interval between each frame of the trajectory. We then divided by the length of the interval using the following equation: Thus, we obtained the number of water molecules per unit Angstrom.

Free Energy Computation.
As stated by Wambo et al. [27], considering that the simulation is performed in the isobaric−isothermal ensemble, the density function of water along the channel can be considered proportional to the probability distribution function of the ensemble: where ΔGðz i Þ is the Gibbs free energy at point z i . If ΔGðz A Þ ≡ 0 at a reference point z A , the free energy profile along the channel may then be expressed as follows:  Table 1). The GMQE value of the DerfAQP1 protein is between 0 and 1 and is close to 1, indicating a high reliability of the three protein structures (Table. 1). The QMEAN Z-scores showed a good agreement between the DerfAQP1 protein structure and the template protein structure (Table 1). Accordingly,

BioMed Research International
we determined the reliability of the tertiary structure of DerfAQP1 obtained by homologous modeling. By using the automatic PSF builder and Add Solvation Box plugins in the VMD software, the monomer DerfAQP1 was placed in an 80 × 80 × 100 Å 3 box, and the TIP3P water model and 52 Na + and 55 Clions were added throughout the system. The system consisted of a total of 59,914 atoms ( Figure 1).
According to the changes of the root mean square deviation (RMSD) of the whole protein, protein without hydrogen atoms, and protein side chains (Figure 2), DerfAQP1 monomers exhibited significant fluctuation at the beginning of simulations but reached a conformational steady state after 30 ns with a RMSD of 3.8-4.2 Å. We saved the files for free energy computation.

Conduction Channel, SF, and NPA Motifs.
Comparison of the radii of DerfAQP1 at the equilibrated state to the crystal structures of bAQP1 (PDB code No. 1J4N), hAQP4 (PDB code No. 3GD8), and PfAQP (PDB code No. 3C02) revealed differences between the DerfAQP1 and PfAQP radii, although the radius of DerAQP1 was similar to that of bAQP1. Through the channel, the narrowest site was located on the SF. The two most conserved NPA motifs formed a canonical "fireman's grip-like" structure ( Figure 3).

Water Arrangement in the Channel.
In the last frame, seven water molecules formed a single file in the channel of the DerfAQP1 monomer (Figure 4(a)). Before the water molecules entered the channel, they stayed in the extracellular vestibule, which contained Arg206, Trp57, Met190, Gly200, and Asp207 (Figure 4(b)). The middle section of the DerfAQP1 channel contained two conserved NPA motifs similar to that in bAQP1 and hAQP4 (Figure 4(c)), resulting in a narrower radius than the channel of PfAQP containing Asn-Leu-Ala and Asn-Pro-Ser [30]. In the channel near the cytoplasm, the water molecules interacted with residues His75, Val80, Ile65, and Ile182 (Figure 4(d)).

Free Energy Computation.
Considering that the single file in the channel of the DerfAQP1 monomer is located between -10 Å and 10 Å ð-10 Å < z < 10 ÅÞ, the free energy profile showed an energy barrier of~2.5 kcal/mol with respect to the SF structure ( Figure 5) [31]. As we described previously [20], the error associated to our free energy calculation originates from the thermal fluctuation energy kBT = 0:6 kcal/mol inherent to molecular dynamics simulations. Since we are subtracting two measurements of energy at different points of the channel, this becomes ±1.2 kcal/mol.

Discussion
The current study investigated the working mechanisms of house dust mite AQP to provide information about the future development of AQP inhibitors. MDS showed that seven water molecules were arranged in a single file through the pores of DerfAQP1. As the water molecules trickled through the pores, their polarization changed around the NPA motifs in Asn203 and Asn77. Residues Arg206, Trp57, Met190, Gly200, and Asp207 constitute the extracellular vestibule of the pore, while residues His75, Val80, Ile65, and Ile182 constitute the cytoplasmic portions. Therefore, these key amino acid residues should be the targets for the development of new inhibitors, which can be designed and tested in silico against the docking energies of AQP4 [32]. In 2004, Law and Sansom performed MDS (totaling almost 40 ns) by using the X-ray structures of Aqp1, a homology model of Aqp1 based on the GlpF structure, and the intermediate resolution structures of Aqp1 derived from the electron microscopy [33]. These results suggest that the homology models based on bacterial homologues could be used to derive functionally relevant information on the structural dynamics of mammalian transport proteins. Here, we performed MDS of DerfAQP1 with homology modeling constructed from the X-ray structure of malaria PfAQP (PDB code No. 3C02) [30]. Subsequently, we placed the monomer DerfAQP1 in an 80 × 80 × 100 Å 3 box and added Gly200 Arg206 Trp57 Figure 4: The water-conducting channel of one monomer of DerfAQP1: (a) the hydrophilic residues are shown in licorice, hydrogen in white, oxygen in red, carbon in yellow, and nitrogen in blue. The licorice ball models represent water, (b) residues around the selectivity filter structure; (c) residues around the NPA structure; (d) residues between the cytoplasm and channels.

BioMed Research
International the TIP3P water model and 52 Na + and 55 Clions throughout the system. The RMSD change in the last few nanoseconds of the MDS is insufficient to demonstrate our success.
All identified AQP structures have a common homotetrameric arrangement, and each AQP monomer folds into a right-handed alpha-barrel architecture. Our monomer AQP model in simulation achieved structural equilibrium after 30 ns according to the RMSD. Moreover, the free energy 1 profile showed an energy barrier of~2.5 kcal/mol, and this condition was observed in the SF structure, which was the narrowest part of the pore. Computational free energy profiles deliver valuable structure-function relationships in regard to selectivity. For example, the free energy profiles for the computational mutants of two channels, Aqp1 and GlpF, where the bipolar field was eliminated by artificially discharging backbone atoms [34], suggested that mutating the energy barrier (i.e., SF structure) may effectively destroy mite AQPs.

Conclusions
We conducted MDS of the mite AQP DerfAQP1 and computed its free-energy profile for the first time. We are conducting similar simulations for other mite AQPs. The analysis of water molecules through the water channel can further elucidate the AQP structure and functional implications. The structural comparison of mite AQPs, such as mutating the key amino acid residues, can also be helpful in the development of inhibitors for mite control.

Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
The authors declare that they have no competing interests.