β‐strand interactions at the domain interface critical for the stability of human lens γD‐crystallin

Human age‐onset cataracts are believed to be caused by the aggregation of partially unfolded or covalently damaged lens crystallin proteins; however, the exact molecular mechanism remains largely unknown. We have used microseconds of molecular dynamics simulations with explicit solvent to investigate the unfolding process of human lens γD‐crystallin protein and its isolated domains. A partially unfolded folding intermediate of γD‐crystallin is detected in simulations with its C‐terminal domain (C‐td) folded and N‐terminal domain (N‐td) unstructured, in excellent agreement with biochemical experiments. Our simulations strongly indicate that the stability and the folding mechanism of the N‐td are regulated by the interdomain interactions, consistent with experimental observations. A hydrophobic folding core was identified within the C‐td that is comprised of a and b strands from the Greek key motif 4, the one near the domain interface. Detailed analyses reveal a surprising non‐native surface salt‐bridge between Glu135 and Arg142 located at the end of the ab folded hairpin turn playing a critical role in stabilizing the folding core. On the other hand, an in silico single E135A substitution that disrupts this non‐native Glu135‐Arg142 salt‐bridge causes significant destabilization to the folding core of the isolated C‐td, which, in turn, induces unfolding of the N‐td interface. These findings indicate that certain highly conserved charged residues, that is, Glu135 and Arg142, of γD‐crystallin are crucial for stabilizing its hydrophobic domain interface in native conformation, and disruption of charges on the γD‐crystallin surface might lead to unfolding and subsequent aggregation.     

Sequential unfolding of the hemolysin two‐partner secretion domain from Proteus mirabilis

Protein secretion is a major contributor to Gram‐negative bacterial virulence. Type Vb or two‐partner secretion (TPS) pathways utilize a membrane bound β‐barrel B component (TpsB) to translocate large and predominantly virulent exoproteins (TpsA) through a nucleotide independent mechanism. We focused our studies on a truncated TpsA member termed hemolysin A (HpmA265), a structurally and functionally characterized TPS domain from Proteus mirabilis. Contrary to the expectation that the TPS domain of HpmA265 would denature in a single cooperative transition, we found that the unfolding follows a sequential model with three distinct transitions linking four states. The solvent inaccessible core of HpmA265 can be divided into two different regions. The C‐proximal region contains nonpolar residues and forms a prototypical hydrophobic core as found in globular proteins. The N‐proximal region of the solvent inaccessible core, however, contains polar residues. To understand the contributions of the hydrophobic and polar interiors to overall TPS domain stability, we conducted unfolding studies on HpmA265 and site‐specific mutants of HpmA265. By correlating the effect of individual site‐specific mutations with the sequential unfolding results we were able to divide the HpmA265 TPS domain into polar core, nonpolar core, and C‐terminal subdomains. Moreover, the unfolding studies provide quantitative evidence that the folding free energy for the polar core subdomain is more favorable than for the nonpolar core and C‐terminal subdomains. This study implicates the hydrogen bonds shared among these conserved internal residues as a primary means for stabilizing the N‐proximal polar core subdomain.     

Structure-based optimization of designed Armadillo-repeat proteins

The armadillo domain is a right‐handed super‐helix of repeating units composed of three α‐helices each. Armadillo repeat proteins (ArmRPs) are frequently involved in protein–protein interactions, and because of their modular recognition of extended peptide regions they can serve as templates for the design of artificial peptide binding scaffolds. On the basis of sequential and structural analyses, different consensus‐designed ArmRPs were synthesized and show high thermodynamic stabilities, compared to naturally occurring ArmRPs. We determined the crystal structures of four full‐consensus ArmRPs with three or four identical internal repeats and two different designs for the N‐ and C‐caps. The crystal structures were refined at resolutions ranging from 1.80 to 2.50 Å for the above mentioned designs. A redesign of our initial caps was required to obtain well diffracting crystals. However, the structures with the redesigned caps caused domain swapping events between the N‐caps. To prevent this domain swap, 9 and 6 point mutations were introduced in the N‐ and C‐caps, respectively. Structural and biophysical analysis showed that this subsequent redesign of the N‐cap prevented domain swapping and improved the thermodynamic stability of the proteins. We systematically investigated the best cap combinations. We conclude that designed ArmRPs with optimized caps are intrinsically stable and well‐expressed monomeric proteins and that the high‐resolution structures provide excellent structural templates for the continuation of the design of sequence‐specific modular peptide recognition units based on armadillo repeats.     

Thermal‐ and urea‐induced unfolding processes of glutathione S‐transferase by molecular dynamics simulation

The Schistosoma juponicum 26 kDa glutathione S‐transferase (sj26GST) consists of the N‐terminal domain (N‐domain), containing three alpha‐helices (named H1‐H3) and four anti‐parallel beta‐strands (S1‐S4), and the C‐terminal domain (C‐domain), comprising five alpha‐helices (named H4‐H8). In present work, molecular dynamics simulations and fluorescence spectroscopic were used to gain insights into the unfolding process of sj26GST. The molecular dynamics simulations on sj26GST subunit both in water and in 8 M urea were carried out at 300 K, 400 K and 500 K, respectively. Spectroscopic measurements were employed to monitor structural changes. Molecular dynamics simulations of sj26GST subunit induced by urea and temperature showed that the initial unfolding step of sj26GST both in water and urea occurred on N‐domain, involving the disruption of helices H2, H3 and strands S3 and S4, whereas H6 was the last region exposed to solution and was the last helix to unfold. Moreover, simulations analyses combining with fluorescence and circular dichroism spectra indicated that N‐domain could not fold independent, suggesting that correct folding of N‐domain depended on its interactions with C‐domain. We further proposed that the folding of GSTs could begin with the hydrophobic collapse of C‐domain whose H4, H5, H6 and H7 could move close to each other and form a hydrophobic core, especially H6 wrapped in the hydrophobic center and beginning spontaneous formation of the helix. S3, S4, H3, and H2 could form in the wake of the interaction between C‐domain and N‐domain. The paper can offer insights into the molecular mechanism of GSTs unfolding. © 2014 Wiley Periodicals, Inc. Biopolymers 103: 247–259, 2015.     

Capping motifs stabilize the leucine‐rich repeat protein PP32 and rigidify adjacent repeats

Capping motifs are found to flank most β‐strand‐containing repeat proteins. To better understand the roles of these capping motifs in organizing structure and stability, we carried out folding and solution NMR studies on the leucine‐rich repeat (LRR) domain of PP32, which is composed of five tandem LRR, capped by α‐helical and β‐hairpin motifs on the N‐ and C‐termini. We were able to purify PP32 constructs lacking either cap and containing destabilizing substitutions. Removing the C‐cap results in complete unfolding of PP32. Removing the N‐cap has a much less severe effect, decreasing stability but retaining much of its secondary structure. In contrast, the dynamics and tertiary structure of the first two repeats are significantly perturbed, based on ¹H‐¹⁵N relaxation studies, chemical shift perturbations, and residual dipolar couplings. However, more distal repeats (3 to C‐cap) retain their native tertiary structure. In this regard, the N‐cap drives the folding of adjacent repeats from what appears to be a molten‐globule‐like state. This interpretation is supported by extensive analysis using core packing substitutions in the full‐length and N‐cap‐truncated PP32. This work highlights the importance of caps to the stability and structural integrity of β‐strand‐containing LRR proteins, and emphasizes the different contributions of the N‐ and C‐terminal caps.     

Predicting repeat protein folding kinetics from an experimentally determined folding energy landscape

The Notch ankyrin domain is a repeat protein whose folding has been characterized through equilibrium and kinetic measurements. In previous work, equilibrium folding free energies of truncated constructs were used to generate an experimentally determined folding energy landscape (Mello and Barrick, Proc Natl Acad Sci USA 2004;101:14102–14107). Here, this folding energy landscape is used to parameterize a kinetic model in which local transition probabilities between partly folded states are based on energy values from the landscape. The landscape‐based model correctly predicts highly diverse experimentally determined folding kinetics of the Notch ankyrin domain and sequence variants. These predictions include monophasic folding and biphasic unfolding, curvature in the unfolding limb of the chevron plot, population of a transient unfolding intermediate, relative folding rates of 19 variants spanning three orders of magnitude, and a change in the folding pathway that results from C‐terminal stabilization. These findings indicate that the folding pathway(s) of the Notch ankyrin domain are thermodynamically selected: the primary determinants of kinetic behavior can be simply deduced from the local stability of individual repeats.     

The Molecular Mechanism of Domain Swapping of the C-Terminal Domain of the SARS-Coronavirus Main Protease

In three-dimensional domain swapping, two protein monomers exchange a part of their structures to form an intertwined homodimer, whose subunits resemble the monomer. Several viral proteins domain swap to increase their structural complexity or functional avidity. The main protease (M^pro) of the severe acute respiratory syndrome (SARS) coronavirus proteolyzes viral polyproteins and has been a target for anti-SARS drug design. Domain swapping in the α-helical C-terminal domain of M^pro (M^proC) locks M^pro into a hyperactive octameric form that is hypothesized to promote the early stages of viral replication. However, in the absence of a complete molecular understanding of the mechanism of domain swapping, investigations into the biological relevance of this octameric M^pro have stalled. Isolated M^proC can exist as a monomer or a domain-swapped dimer. Here, we investigate the mechanism of domain swapping of M^proC using coarse-grained structure-based models and molecular dynamics simulations. Our simulations recapitulate several experimental features of M^proC folding. Further, we find that a contact between a tryptophan in the M^proC domain-swapping hinge and an arginine elsewhere forms early during folding, modulates the folding route, and promotes domain swapping to the native structure. An examination of the sequence and the structure of the tryptophan containing hinge loop shows that it has a propensity to form multiple secondary structures and contacts, indicating that it could be stabilized into either the monomer- or dimer-promoting conformations by mutations or ligand binding. Finally, because all residues in the tryptophan loop are identical in SARS-CoV and SARS-CoV-2, mutations that modulate domain swapping may provide insights into the role of octameric M^pro in the early-stage viral replication of both viruses.     

Comparison of the folding mechanism of highly homologous proteins in the lipid‐binding protein family

The folding mechanism of two closely related proteins in the intracellular lipid‐binding protein family, human bile acid‐binding protein (hBABP), and rat bile acid‐binding protein (rBABP) were examined. These proteins are 77% identical (93% similar) in sequence. Both of these single domain proteins fit well to a two‐state model for unfolding by fluorescence and circular dichroism at equilibrium. Three phases were observed during the unfolding of rBABP by fluorescence but only one phase was observed during the unfolding of hBABP, suggesting that at least two kinetic intermediates accumulate during the unfolding of rBABP that are not observed during the unfolding of hBABP. Fluorine NMR was used to examine the equilibrium unfolding behavior of the W49 side chain in 6‐fluorotryptophan‐labeled rBABP and hBABP. The structure of rBABP appears to be more dynamic than that of hBABP in the vicinity of W49 in the absence of denaturant, and urea has a greater effect on this dynamic behavior for rBABP than for hBABP. As such, the folding behavior of highly sequence related proteins in this family can be quite different. These differences imply that moderately sized proteins with high sequence and structural similarity can still populate quite different structures during folding. Proteins 2009. © 2008 Wiley‐Liss, Inc.     

Unfolding study of a trimeric membrane protein AcrB

The folding of a multi‐domain trimeric α‐helical membrane protein, Escherichia coli inner membrane protein AcrB, was investigated. AcrB contains both a transmembrane domain and a large periplasmic domain. Protein unfolding in sodium dodecyl sulfate (SDS) and urea was monitored using the intrinsic fluorescence and circular dichroism spectroscopy. The SDS denaturation curve displayed a sigmoidal profile, which could be fitted with a two‐state unfolding model. To investigate the unfolding of separate domains, a triple mutant was created, in which all three Trp residues in the transmembrane domain were replaced with Phe. The SDS unfolding profile of the mutant was comparable to that of the wild type AcrB, suggesting that the observed signal change was largely originated from the unfolding of the soluble domain. Strengthening of trimer association through the introduction of an inter‐subunit disulfide bond had little effect on the unfolding profile, suggesting that trimer dissociation was not the rate‐limiting step in unfolding monitored by fluorescence emission. Under our experimental condition, AcrB unfolding was not reversible. Furthermore, we experimented with the refolding of a monomeric mutant, AcrB_Δloop, from the SDS unfolded state. The CD spectrum of the refolded AcrB_Δloop superimposed well onto the spectra of the original folded protein, while the fluorescence spectrum was not fully recovered. In summary, our results suggested that the unfolding of the trimeric AcrB started with a local structural rearrangement. While the refolding of secondary structure in individual monomers could be achieved, the re‐association of the trimer might be the limiting factor to obtain folded wild‐type AcrB.     

The C‐terminal calcium‐sensitive disordered motifs regulate isoform‐specific polymerization characteristics of calsequestrin

Calsequestrin (CASQ) exists as two distinct isoforms CASQ1 and CASQ2 in all vertebrates. Although the isoforms exhibit unique functional characteristic, the structural basis for the same is yet to be fully defined. Interestingly, the C‐terminal region of the two isoforms exhibit significant differences both in length and amino acid composition; forming Dn‐motif and DEXn‐motif in CASQ1 and CASQ2, respectively. Here, we investigated if the unique C‐terminal motifs possess Ca²⁺‐sensitivity and affect protein function. Sequence analysis shows that both the Dn‐ and DEXn‐motifs are intrinsically disordered regions (IDRs) of the protein, a feature that is conserved from fish to man. Using purified synthetic peptides, we show that these motifs undergo distinctive Ca²⁺‐mediated folding suggesting that these disordered motifs are Ca²⁺‐sensitivity. We generated chimeric proteins by swapping the C‐terminal portions between CASQ1 and CASQ2. Our studies show that the C‐terminal portions do not play significant role in protein folding. An interesting finding of the current study is that the switching of the C‐terminal portion completely reverses the polymerization kinetics. Collectively, these data suggest that these Ca²⁺‐sensitivity IDRs located at the back‐to‐back dimer interface influence isoform‐specific Ca²⁺‐dependent polymerization properties of CASQ. © 2014 Wiley Periodicals, Inc. Biopolymers 103: 15–22, 2015.     

Temperature‐induced partially unfolded state of hUBF HMG Box‐5: Conformational and dynamic investigations of the Box‐5 thermal intermediate ensemble

Human upstream binding factor (hUBF) HMG Box‐5 is a highly conserved protein domain, containing 84 amino acids and belonging to the family of the nonspecific DNA‐binding HMG boxes. Its native structure adopts a twisted L shape, which consists of three α‐helices and two hydrophobic cores: the major wing and the minor wing. In this article, we report a reversible three‐state thermal unfolding equilibrium of hUBF HMG Box‐5, which is investigated by differential scanning calorimetry (DSC), circular dichroism spectroscopy, fluorescence spectroscopy, and NMR spectroscopy. DSC data show that Box‐5 unfolds reversibly in two separate stages. Spectroscopic analyses suggest that different structural elements exhibit noncooperative transitions during the unfolding process and that the major form of the Box‐5 thermal intermediate ensemble at 55°C shows partially unfolded characteristics. Compared with previous thermal stability studies of other boxes, it appears that Box‐5 possesses a more stable major wing and two well separated subdomains. NMR chemical shift index and sequential ¹H^N_i‐¹H^N_i+1 NOE analyses indicate that helices 1 and 2 are native‐like in the thermal intermediate ensemble, while helix 3 is partially unfolded. Detailed NMR relaxation dynamics are compared between the native state and the intermediate ensemble. Our results implicate a fluid helix‐turn‐helix folding model of Box‐5, where helices 1 and 2 potentially form the helix 1‐turn‐helix 2 motif in the intermediate, while helix 3 is consolidated only as two hydrophobic cores form to stabilize the native structure. Proteins 2009. © 2009 Wiley‐Liss, Inc.     

X‐ray crystal structure of the N‐terminal region of Moloney murine leukemia virus integrase and its implications for viral DNA recognition

The retroviral integrase (IN) carries out the integration of a dsDNA copy of the viral genome into the host DNA, an essential step for viral replication. All IN proteins have three general domains, the N‐terminal domain (NTD), the catalytic core domain, and the C‐terminal domain. The NTD includes an HHCC zinc finger‐like motif, which is conserved in all retroviral IN proteins. Two crystal structures of Moloney murine leukemia virus (M‐MuLV) IN N‐terminal region (NTR) constructs that both include an N‐terminal extension domain (NED, residues 1–44) and an HHCC zinc‐finger NTD (residues 45–105), in two crystal forms are reported. The structures of IN NTR constructs encoding residues 1–105 (NTR_1–105) and 8–105 (NTR_8–105) were determined at 2.7 and 2.15 Å resolution, respectively and belong to different space groups. While both crystal forms have similar protomer structures, NTR_1–105 packs as a dimer and NTR_8–105 packs as a tetramer in the asymmetric unit. The structure of the NED consists of three anti‐parallel β‐strands and an α‐helix, similar to the NED of prototype foamy virus (PFV) IN. These three β‐strands form an extended β‐sheet with another β‐strand in the HHCC Zn²⁺ binding domain, which is a unique structural feature for the M‐MuLV IN. The HHCC Zn²⁺ binding domain structure is similar to that in HIV and PFV INs, with variations within the loop regions. Differences between the PFV and MLV IN NEDs localize at regions identified to interact with the PFV LTR and are compared with established biochemical and virological data for M‐MuLV. Proteins 2017; 85:647–656. © 2016 Wiley Periodicals, Inc.     

How cooperative are protein folding and unfolding transitions?

A thermodynamically and kinetically simple picture of protein folding envisages only two states, native (N) and unfolded (U), separated by a single activation free energy barrier, and interconverting by cooperative two‐state transitions. The folding/unfolding transitions of many proteins occur, however, in multiple discrete steps associated with the formation of intermediates, which is indicative of reduced cooperativity. Furthermore, much advancement in experimental and computational approaches has demonstrated entirely non‐cooperative (gradual) transitions via a continuum of states and a multitude of small energetic barriers between the N and U states of some proteins. These findings have been instrumental towards providing a structural rationale for cooperative versus noncooperative transitions, based on the coupling between interaction networks in proteins. The cooperativity inherent in a folding/unfolding reaction appears to be context dependent, and can be tuned via experimental conditions which change the stabilities of N and U. The evolution of cooperativity in protein folding transitions is linked closely to the evolution of function as well as the aggregation propensity of the protein. A large activation energy barrier in a fully cooperative transition can provide the kinetic control required to prevent the accumulation of partially unfolded forms, which may promote aggregation. Nevertheless, increasing evidence for barrier‐less “downhill” folding, as well as for continuous “uphill” unfolding transitions, indicate that gradual non‐cooperative processes may be ubiquitous features on the free energy landscape of protein folding.     

Exploring the folding pathway of green fluorescent protein through disulfide engineering

We have introduced two disulfide crosslinks into the loop regions on opposite ends of the beta barrel in superfolder green fluorescent protein (GFP) in order to better understand the nature of its folding pathway. When the disulfide on the side opposite the N/C‐termini is formed, folding is 2× faster, unfolding is 2000× slower, and the protein is stabilized by 16 kJ/mol. But when the disulfide bond on the side of the termini is formed we see little change in the kinetics and stability. The stabilization upon combining the two crosslinks is approximately additive. When the kinetic effects are broken down into multiple phases, we observe Hammond behavior in the upward shift of the kinetic m‐value of unfolding. We use these results in conjunction with structural analysis to assign folding intermediates to two parallel folding pathways. The data are consistent with a view that the two fastest transition states of folding are "barrel closing" steps. The slower of the two phases passes through an intermediate with the barrel opening occurring between strands 7 and 8, while the faster phase opens between 9 and 4. We conclude that disulfide crosslink‐induced perturbations in kinetics are useful for mapping the protein folding pathway.     

Nonuniform chain collapse during early stages of staphylococcal nuclease folding detected by fluorescence resonance energy transfer and ultrarapid mixing methods

The development of tertiary structure during folding of staphylococcal nuclease (SNase) was studied by time‐resolved fluorescence resonance energy transfer measured using continuous‐ and stopped‐flow techniques. Variants of this two‐domain protein containing intradomain and interdomain fluorescence donor/acceptor pairs (Trp and Cys‐linked fluorophore or quencher) were prepared to probe the intradomain and interdomain structural evolution accompanying SNase folding. The intra‐domain donor/acceptor pairs are within the β‐barrel domain (Trp27/Cys64 and Trp27/Cys97) and the interdomain pair is between the α‐helical domain and the β‐barrel domain (Trp140/Cys64). Time‐resolved energy transfer efficiency accompanying folding and unfolding at different urea concentrations was measured over a time range from 30 μs to ～10 s. Information on average donor/acceptor distances at different stages of the folding process was obtained by using a quantitative kinetic modeling approach. The average distance for the donor/acceptor pairs in the β‐barrel domain decreases to nearly native values whereas that of the interdomain donor/acceptor pairs remains unchanged in the earliest intermediate (<500 μs of refolding). This indicates a rapid nonuniform collapse resulting in an ensemble of heterogeneous conformations in which the central region of the β‐barrel domain is well developed while the C‐terminal α‐helical domain remains disordered. The distance between Trp140 and Cys64 decreases to native values on the 100‐ms time scale, indicating that the α‐helical domain docks onto the preformed β‐barrel at a late stage of the folding. In addition, the unfolded state is found to be more compact under native conditions, suggesting that changes in solvent conditions may induce a nonspecific hydrophobic collapse.     

Folding and association of the human cell cycle regulatory proteins ckshs1 and ckshs2.

The two human proteins ckshs1 and ckshs2 are each 79 amino acids in length and consist of a four-stranded beta-sheet capped at one end by two alpha-helices. They are members of the cks family of essential cell cycle regulatory proteins that can adopt two native states, a monomer and a domain-swapped dimer formed by exchange of a C-terminal beta-strand. ckshs1 and ckshs2 both have marginal thermodynamic stability (the free energies of unfolding at 25 degrees C are 3.0 and 2.5 kcal/mol, respectively) and low kinetic stability (the rates of unfolding in water are approximately 1 s(-1)). Refolding of their denatured states to the monomeric forms of the proteins is slowed by transient oligomerization that is likely to occur via domain swapping. The folding behavior of ckshs1 and ckshs2 is markedly different from that of suc1, the cks protein from Schizosaccharomyces pombe, but the domain swapping propensities are similar. The greater thermodynamic and kinetic stability of suc1 and the population of a folding intermediate are most likely a consequence of its larger size (113 residues). The similarity in the domain swapping propensities, despite the contrast in other biophysical properties, may be attributable to the common double-proline motif in the hinge loop that connects the swapped domain to the rest of the protein. The motif was shown previously for suc1 to control the equilibrium between the monomer and the domain-swapped dimer. Finally, according to our model, the kinetic barrier separating the monomer and the domain-swapped dimer arises because the protein must unfold for beta-strand exchange to occur. Consistent with this, interconversion between the two states is much faster in the human proteins than it is for suc1, reflecting the faster unfolding rates of the former.     

Mechanism of formation of the C‐terminal β‐hairpin of the B3 domain of the immunoglobulin binding protein G from Streptococcus. III. Dynamics of long‐range hydrophobic interactions

A 20‐residue peptide, IG(42–61), derived from the C‐terminal β‐hairpin of the B3 domain of the immunoglobulin binding protein G from Streptoccocus was studied using circular dichroism, nuclear magnetic resonance (NMR) spectroscopy at various temperatures and by differential scanning calorimetry (DSC). Unlike other related peptides studied so far, this peptide displays two heat capacity peaks in DSC measurements (at a scanning rate of 1.5 deg/min at a peptide concentration of 0.07 mM), which suggests a three‐state folding/unfolding process. The results from DSC and NMR measurements suggest the formation of a dynamic network of hydrophobic interactions stabilizing the structure, which resembles a β‐hairpin shape over a wide range of temperatures (283–313 K). Our results show that IG (42–61) possesses a well‐organized three‐dimensional structure stabilized by long‐range hydrophobic interactions (Tyr50 ··· Phe57 and Trp48 ··· Val59) at T = 283 K and (Trp48 ··· Val59) at 305 and 313 K. The mechanism of β‐hairpin folding and unfolding, as well as the influence of peptide length on its conformational properties, are also discussed. Proteins 2010. © 2009 Wiley‐Liss, Inc.     

Constant pH molecular dynamics of proteins in explicit solvent with proton tautomerism

pH is a ubiquitous regulator of biological activity, including protein‐folding, protein‐protein interactions, and enzymatic activity. Existing constant pH molecular dynamics (CPHMD) models that were developed to address questions related to the pH‐dependent properties of proteins are largely based on implicit solvent models. However, implicit solvent models are known to underestimate the desolvation energy of buried charged residues, increasing the error associated with predictions that involve internal ionizable residue that are important in processes like hydrogen transport and electron transfer. Furthermore, discrete water and ions cannot be modeled in implicit solvent, which are important in systems like membrane proteins and ion channels. We report on an explicit solvent constant pH molecular dynamics framework based on multi‐site λ‐dynamics (CPHMD^MSλD). In the CPHMD^MSλD framework, we performed seamless alchemical transitions between protonation and tautomeric states using multi‐site λ‐dynamics, and designed novel biasing potentials to ensure that the physical end‐states are predominantly sampled. We show that explicit solvent CPHMD^MSλD simulations model realistic pH‐dependent properties of proteins such as the Hen‐Egg White Lysozyme (HEWL), binding domain of 2‐oxoglutarate dehydrogenase (BBL) and N‐terminal domain of ribosomal protein L9 (NTL9), and the pK_a predictions are in excellent agreement with experimental values, with a RMSE ranging from 0.72 to 0.84 pK_a units. With the recent development of the explicit solvent CPHMD^MSλD framework for nucleic acids, accurate modeling of pH‐dependent properties of both major class of biomolecules—proteins and nucleic acids is now possible. Proteins 2014; 82:1319–1331. © 2013 Wiley Periodicals, Inc.