PM: Plasma Membrane; EM: Energy Minimization; BLAST: Basic Local Alignment Search Tool

PM (plasma membrane)-H⁺ATPase is a plasma membrane bound proton (H⁺) pumping ATPase and is member of p-type (E1-E2 type) ATPases representing a superfamily of cation transporters involved in the transport of cations like H⁺, Na⁺, K⁺, Ca²⁺, Mg²⁺, etc. at the expense of ATP. Being a member of P-ATPase, the reaction catalysed by it proceed with the formation of a phosphorylated intermediate, unlike that of other ATP hydrolyzing enzymes, V₁V₀-ATPase present in tonoplast and the membranes of the intracellular organelles of eukaryotic cells and F₁F₀ ATPase present in thylakoid and mitochondrial membranes. Like PM-H⁺-ATPase, V₁V⁰-ATPase (or V-ATPase) is also primarily a H⁺ pumping enzyme, the two enzymes differ greatly in their structure, besides being located in different membranes of the cells. Moreover, PM- H⁺-ATPase is unique only to the plants and fungi.

The free energy gradient (ΔµH⁺) produced by plasma membrane H⁺-ATPase (PM-H⁺-ATPase) provides driving force for secondary active transport for the maintenance of ionic balance across the membrane and uptake of necessary elements or molecules, both charge and neutral, mediated by various transporters like uni-, sym- and antiport carriers (Sze,1985). These transporters are thus called secondary transporters, as they utilize the energy in the form of proton (H⁺ ion) gradient generated by the PM H⁺- ATPases. Thus, this membrane enzyme has essential house-keeping functions in the life of plants. It is also involved in special functions such as turgor-dependent processes of growth extension and movements of stomatal guard cells and pulvini (Okazaki, 2002). In most plant species PM-H⁺ATPase is encoded by a large gene family. In tomato (Lycopersicon esculentum cv. large cherry) seven genes have been identified and genomic DNA blot analysis indicates a family of at least 10 genes (Oufattole et al., 2000). Extensive search of cDNA and genomic libraries from Nicotiana plumbaginifolia resulted in identification of nine PM-H⁺-ATPase genes (Perez et al., 1992; Moriau et al., 1993; Oufattole et al., 2000), while complete sequencing of Arabidopsis genome has uncovered 12 genes (Palmgren, 2001) and that of Oryza sativa (Japonica) has led to retrieval of 10 OSA (Oryza sativa PM-H⁺ATPase) genes (Arango et al., 2003). Out of the 10 annotated OSA genes in Oryza sativa, six (OSA1, OSA2, OSA3, OSA4, OSA8 and OSA9) have been characterized at cDNA level (Wada et al., 1992; Ookura et al., 1994; Kikuchi et al., 2003). All these isoforms have been cloned from various organs of plants at different stages of development grown under control condition. The existence of multiple isoforms of the enzyme might indicate that some pumps could be redundant. Nonetheless, isoform diversity may also be related to cellular differentiation with individual isoform exhibiting tissue and developmental specific expression and having slightly different biochemical and regulatory properties (Palmgren, 2001; Arango et al., 2003).

The overall shape of proton pumps has been revealed by the electron microscopy. The ion translocation proceeds through a Ping-Pong mechanism, cation species being sequentially transferred in opposite directions across the membrane. The opening and closure of the access channels are tightly coupled to conformational changes within the cytoplasmic domains induced by ATP binding, phosphorylation, subsequent dephosphorylation and phosphate release from the enzyme. Therefore, ions are sequentially taken up at ion binding sites from one side of the lipid bilayer, occluded within the transmembrane domain, and then released on the other side of the membrane, in exchange for the second ion species (Kühlbrandt, 2004; Møller et al., 1996; Apell, 2004; Albers, 1967; Post et al., 1969).

The X-ray crystal structure of the sarcoplasmic reticulum Ca²⁺-ATPase from rabbit skeletal muscle and AHA2 from Arabidopsis thaliana are also available (Toyoshima et al., 2000; Toyoshima et al., 2002; Sørensen et al., 2004; Pedersen et al., 2007) and these structures provide critical insights into the structure and mechanism of P-type ATPases. These structures explained the general architecture of the catalytic subunit of P-type ATPase which divided the cytoplasm part of the protein into three domains: (i) A domain (actuator), (ii) P domain (phosphorylation), (iii) N (nucleotide binding) domain and a transmembrane segment of 10 helices (TM1 to TM10) bearing the two embedded Ca²⁺ binding sites. However, despite the existence of a large no. of isoforms of PM-H⁺- ATPase with great differences in homologies, a comparative modeling of various isoforms has not been done. This is despite the fact that expressions of many isoforms have been reported to be tissue specific. For example, in Arabidopsis thaliana AHA3 expresses in the phloem companion cells (DeWitt et al., 1991)^[20], AHA9 in anthers (Houlné et al., 1994) and AHA6, AHA8, AHA9 and AHA12 in pollen grains. Reports are also there for tissue specific expression of PM-H⁺ATPase isoforms in other plants, like Nicotiana plumbaginifolia (Lefebvre et al., 2005) and Vicia faba (Hentzen et al., 1996). The only model available for PM-H⁺ATPase is for Arabidopsis thaliana AHA2 (Pedersen et al., 2007).

Recently, expression analysis of all the PM-H⁺ATPase isoforms of Arabidopsis thaliana in root tissue using Affymetrix GeneChip revealed maximum expression of AHA1 and AHA2 in comparison to other isoforms (Dinneny et al., 2008). Palmgren et al., (1999), however, have reported no expression of AHA2 in stem tissue. Besides, a comparative study of the activity of the two enzyme isoforms expressed heterologously it has been observed that the specific activity of AHA2 is approximately double than that of AHA1 (Palmgren et al., (1999)). Hence, it was intriguing to know whether there could be difference in the two isoforms of the enzyme at the level of their association with the membrane and substrate affinity/availability. Therefore, in this report we have modeled AHA1 isoform for a putative proton binding site in the transmembarne region by homology modeling method based on the template of the auto-inhibited AHA2.

Most comparative protein modeling tools adopt rigid body assembly approach which constructs the model from a few core regions and loops obtained from related template protein structures. Another approach adopts segment matching which relies on structurally conserved atoms from template and calculate the co-ordinates of other atoms. The MODELLER uses an approach involving distance geometry or optimization techniques to satisfy spatial restraints obtained from the alignment of target sequence with the template structures. MODELLER allows degree of flexibility and automation which is widely used program in comparative modeling and developed a database of protein models (Modbase). In general, the homology modeling method is based on the assumption that the structure of an unknown protein is similar to known structures of some reference proteins. The MODELLER software employs probability density functions (PDFs) as the spatial restraints rather than energy (Sali et al., 1993; 1994; 1995). The main chain conformation of a given residue in the model is described by restraints that depend upon the residue type, the main chain conformation of equivalent residues in the reference proteins and the local sequence similarity. The PDFs used in restraining the model structure are derived from correlations between structural features in a database of families of homologous proteins aligned on the basis of their 3D structures. These functions are used to restrain Cα- Cα distances, main chain N-O distances, main-chain and side-chain dihedral angles, etc. The 3D model of a protein is obtained by optimization of the molecular pdf such that the model violates the input restraints as little as possible. The molecular pdf is derives as a combination of pdfs restraining individual spatial features of the whole molecule. The optimization procedure is a variable target function method that applies the conjugate gradients algorithm to positions of all non-hydrogen atoms. In our case, the sequence identity between the target and template is 94.2%, and hence the fold (secondary structures arrangement in 3D space well conserved), C-alpha trace, and backbone in the model derived is likely to be highly accurate. Besides, since there is high sequence identity, most of the side chains of the model will have similar side sign torsion angles and any bad steric contacts are removed during energy minimization of the modeled structure. Hence, it was decided to build the model of AHA1 using MODELLER programme. The model built was validated using widely followed Procheck and Verify3D programmes.

In the present investigation, the protein sequence of Arabidopsis thaliana AHA1 was retrieved from the NCBI database (NP_179486). BLASTP (Altschul et al., 1990) search was performed, against Brook Heaven Protein Data Bank (PDB) (Berman et al., 2000) with default parameters to find the suitable templates for homology modeling. Based on the maximum identity with high score and lower e-value, crystal structure of the plasma membrane proton pump from Arabidopsis thaliana at 3.6 Å resolution (PDB code: 3b8c) was used as template for homology modeling. The sequence identity was found to be 94 percent (similarity = 97%, E-Value = 0.0). The ClustalW (http://www.ebi.ac.uk/clustalw) (Thompson et al., 1994) program was used for sequence alignment.

3D Structure Generation

The academic version of MODELER (http//:www.salilab.org/modeler) (Sali et al., 1993) was used for 3D structure generation based on the information obtained from sequence alignment. Out of 20 models generated by MODELLER, the one with the best G-score of PROCHECK (Laskoswki et al., 1993) and with the best VERIFY3D (Luthy et al., 1992) profile was subjected to energy minimization. Using the parameters as a distance-dependent dielectric constant ∈ = 1.0 and non bonding cutoff of 14Å, CHARMm (Brooks et al., 1983), force field and CHARM-all-atom charges, initially a steepest descent algorithm was used to remove close van der Waals contacts, followed by conjugate gradient minimization until the energy showed stability in sequential repetition. All hydrogen atoms were included during the calculation. The energy minimization was started with main chain of the core and then all core side chains were subjected to the same. All calculations were performed by using ACCELRYS DS Modeling 2.0 (Accelrys Inc. San Diego, CA 92121, USA) software suite. During these steps, the quality of the initial model was improved. VERIFY3D (a structure evaluation server) were used to check the residue profiles of the obtained three-dimensional models. The weighted root mean square deviation (RMSD) of the modeled protein was calculated using combinatorial extension (CE) algorithm (Shindyalov et al., 1998). In order to assess the stereo-chemical qualities of the three dimensional models PROCHECK analysis was performed.

Transmembrane Topology Predictions

The transmembrane topology of AHA1 was predicted from the amino acid sequence by averaging the results of four different predictive algorithms: DAS (Cserzo et al., 1997), HMMTOP (Tusnady et al., 1998), TMHMM (Krogh et al., 2001) and TMPRED (Hofmann et al., 1993) (Table 1). The accuracy of the prediction were assessed by using the same algorithms on 3b8c protein and the results were compared with the topology defined in the 3.6 Å structure (3b8c) (Table 2).

Table 1: AHA1 topological prediction.

Table 2: 3b8c topological prediction.

Model Building

The structure of the AHA1 protein has not yet determined experimentally (X-ray or NMR), and therefore we built a model following homology modeling protocol. BLASTP search was performed against PDB with default parameters to find the suitable templates for homology modeling. Based on the maximum identity with high score and lower e-value 3b8c (AHA2) was used as the template for homology modeling. Sequence alignment between AHA1 and 3b8c was done using ClustalW program. The sequence – structure alignment used for model building is shown in figure. 1. The alignment is characterized by some insertions and deletions in the loop regions. The first 11 residues and the last 105 residues were deleted because the corresponding residues were not present in the template. Modeling was carried out from 12 to 844 residues followed by a rigorous refinement of the model by means of energy minimization using CHARMm (Brooks et al., 1983) force field. The final stable structure of AHA1 is shown in Figure 2.

Figure 1: Sequence alignment of AHA1 from Arabidopsis thaliana with 3b8c from Arabidopsis thaliana done using CLUSTALW server that was subsequently submitted to MODELLER. The conserved regions are indicated by ‘*’.

Figure 2: Ribbon representation of modeled AHA1. The α-helices and β -sheets are shown as helices and ribbons, respectively. The rest are shown as loops. The figure was prepared by PyMol.

Protein Structure Validation

To validate the model obtained from the homology modeled AHA1 structure, a Ramachandran plot was drawn and the structure was analyzed by PROCHECK (Laskoswki et al., 1993), a well-known protein structure checking program. It was found that the phi/psi angles of 77.4% residues fell in the most favored regions, 17.5% residues lied in the additional allowed regions 4.0% fell in the generously allowed regions. Only 1.1% of the residues lied in the disallowed conformations (Figure 3). The overall PROCHECK G-factor for the homology modeled structure was -0.35. The score indicates that the modeled structure is acceptable as value is greater than the acceptable value of – 0.50. The structural superimposition of C^α trace of the model over the template structure 3b8c (Figure 4) using CE program (http://cl.sdsc.edu/ce.html) gave a root mean square deviation (RMSD) of 0.4 Å (Z-score 8.7) indicating further that the structure of the model was valid.

Figure 3: Ramachandran’s Map of AHA1 protein. The plot calculation on 3D model of AHA1 protein was calculated with the PROCHECK program.

Figure 4: Superimposition of C^α trace of AHA1 (represented in white) and 3b8c (represented in pink color).

The amino acid sequence alignment between the template and the final model of AHA1 was generated using CLUSTALW program. The secondary structures was analyzed and compared by the JOY (Mizuguchi et al., 1998) server (protein sequence-structure representation and analysis) and it was found that the secondary structures of the template and final AHA1 model were highly conserved. This also shows that the final model obtained for AHA1 in this study (Figure. 5) is highly reliable. However, despite difference in the tissue-specific expression of the two proteins, phylogenetic tree constructed for various H⁺ATPase isoforms (Sahu and Shaw, 2009) showed no remarkable difference in their structure, suggesting that the enzyme isoform evolution may not be linked to its tissue-specific expression.

Figure 5: Structure based sequence alignment of template and final structures of the AHA1 using JOY program. The key to the JOY annotation is as follows: lowercase red letter, α-helix; lowercase blue letter, β-strand; lowercase maroon letter, 3₁₀-helix; uppercase letter, solvent-inaccessible residue; lowercase letter, solvent-accessible residue; italic lowercase letter, positive Φ.

The four clearly defined domains in AHA2, i.e. a transmembrane domain with ten helices, M1 to M10, and the three cytosolic domains comprising of N (nucleotide binding; residues 338-488), P (phosphorylation; residues 308- 337 and 489-625) and A (actuator; residues 12-57 and 129- 233) are also present in the model derived for AHA1. The N domain remains inserted into the P domain through a hinge (including the conserved sequence motif DPPR) and with bound nucleotide it can move towards the P domain to assemble the catalytic site, at which Asp329 will become phosphorylated once every pumping cycle. The A domain, is situated on the top of M2 stimulates dephosphorylation of Asp, which protrudes as a pole out of the membrane, and it is further connected to the M1 and M3 transmembrane segments by extended loops. In the transmembrane domain, the M1 helix shows a prominent 90° kink imposed by a proline residue, Pro⁶⁸, which is conserved in type III P-type ATPases.

Two ion-binding sites are present in the SR Ca²⁺-ATPase structure (Toyoshima et al., 2000) but the structure of AHA2 H⁺-ATPase suggests the presence of only one proton-binding site. The structure of AHA1 also suggests the presence of only one proton-binding site. The residue supposed to contribute to the proton-binding site in AHA1is Asp⁶⁸⁴. AHA1 seems to be strictly conserved for the corresponding amino acid residues in the AHA2. In fact, this Asp residue appears to be conserved among the P-type ATPases (http://www.patbase.kvl.dk) across the species as it aligns with Asp⁷³⁰ of PMA1, Asp⁸⁰⁰ of sercala and Asp⁸⁰⁸ of Na/K.

The side chain of Asp⁶⁸⁴ AHA1 (TM6), as well as the backbone carbonyls of the amino acid residues in the transmembrane helix4 (Ile282, Gly283 and Ile 285) could provide the liganding groups for a high-affinity hydronium ion coordination site in the AHA1 similar to that in AHA2. The mutagenesis experiments support for such a model; it has been observed that after removal of the potential protonation group (D684N, D684A, D684V and D684R) of the Asp684AHA2 side chain there occurs no change in E₁P-E₂P enzyme conformational and proton transport although the ATP hydrolysis continues (Buch-Pedersen et al., 2000; Buch-Pedersen et al., 2003). The mutagenesis experiment also led to the suggestion that the two oxygen of the Asp⁶⁸⁴ AHA1 side chain that participate in the H⁺- transporting mechanism can be discriminated with respect to their function . This is because the Asp⁶⁸⁴ can be replaced with asparagines with apparently little effect upon initial binding, and removing both the side chains oxygen completely lowers the rate of initial proton binding markedly. Thus, the side chain of the amino acid at position 684 has been ascribed a role in initial proton binding. As reported for AHA2, the Asp⁶⁸⁴ AHA1 sidechain oxygen seem also to influence the rate of initial proton binding to the E₁ proton-binding site.

Effect of Mutations on Molecular Surface Area and Solvent Accessible Surface Area in AHA1 and AHA2 Protein Structures

Molecular surface area (MSA) and solvent accessible surface area (ASA) were calculated for AHA1 and AHA2 using surface racer program (http://apps.phar.umich.edu/tsodikovlab/index_files/Page756.htm). Table 3 shows the changes in MSA and ASA, both polar and charged, from AHA2 to AHA1 structures due to amino acid mutations. There is a decrease in MSA and ASA from AHA2 (3b8c) to AHA1 due to mutations . However, whether the decrease in the MSA and ASA values could be related to the function of a protein is yet to be addressed. Mutational analysis of the protein could probably provide a suitable answer.

Table 3: Calculated molecular surface area (MSA) and solvent accessible surface area (ASA).

Thus, in this study we built a complete atomic model of Arabidopsis thaliana AHA1 plasma membrane H⁺-ATPase by homology modeling using the 3D structure of AHA2 from Arabidopsis thaliana as a template. The transmembrane topology of the protein was also predicted. The model suggested that a description of the proton translocation pathway through the protein (enzyme) must include a conserved Asp in M6 corresponding to Asp⁶⁸⁴ similar to that in AHA2. This residue is conserved in the P-type H⁺ATPase across the species. The model derived matched to AHA2 in all respect, suggesting that the two proteins should not be very much different in its activity or function despite the difference in their tissue-specific expression, which of course in not the case. The difference in the activity of the two isoforms, nevertheless, could be associated with the changing MSA and ASA or possibly other parameters because of change in the amino acids. Further studies may provide conclusive information in this regard, particularly considering the region which was not included in the modeling.

The study was supported by Department of Biotechnology grant no.BT/BI/04/058/2002 to Distributed Information Sub Center, Institute of Life Sciences, Bhubaneswar, India. Funds contributed by ILS, Bhubaneswar, to BPS and Senior Research Fellowship to BBS awarded by UGC, New Delhi are kindly acknowledged.