Combining elastic network models and linear response theory as tool to understand the global dynamics in allosteric regulation of HCN channels

Hyperpolarization-activated cyclic nucleotide-gated cation (HCN) channels drive the “funny current,” which controls pace-making in the heart (I_f) and repetitive firing in neurons (I_h) (Robinson and Siegelbaum, 2003; Saponaro and DiFrancesco, 2025). These mixed cation selective channels are present in humans in four isoforms (HCN1–4) and belong from a structural point of view to the superfamily of voltage-gated K⁺ channels. Their monomers are typically made of a transmembrane part of the channel (TMPC) and a cytosolic part of the channel (CPC) (Fig. 1). The TMPC includes six transmembrane domains (TMDs S1–S4) in which TMD S1–S4 form the voltage sensor. In the tetrameric assembly of the functional channel TMDs S5 and S6 form the ion-conducting pore domain (PD) with the selectivity filter (SF) (Fig. 1). The CPC contains the so called HCN domain (HCND) at the N terminus and the cyclic nucleotide-binding domain at the C terminus, connected to the pore via a C-linker formed by “elbow” and “shoulder” (Fig. 1).

Despite their structural resemblance to outward rectifying Kv channels, HCN channels are inward rectifiers. The channel open probability increases with membrane hyperpolarization, providing conductance of a mixed K⁺/Na⁺ current, which is furthermore modulated by the concentration of cytosolic cyclic nucleotides (DiFrancesco and Tortora, 1991).

Binding of the second messenger cAMP to the CNBD increases the voltage-response of the channel and promotes opening at lower, more depolarized potentials. The physiological effect of a rise in the cytosolic concentration of the signaling molecules is an acceleration of pace making.

One challenge in understanding structure/function correlates in these physiologically important channels is how negative voltages can open a channel with the same basic architecture of outward rectifying Kv channels. A second related challenge is to unravel the allosteric nature of gating by both voltage in the voltage-sensing domains (VSDs) and ligand binding in the CNBD.

Major advances in understanding structure/function correlates in HCN channels were made by the acquisition of HCN structures with atomic resolution. While structures of the isolated CNBDs were already resolved some time ago (Zagotta et al., 2003; Xu et al., 2010; Lolicato et al., 2011; Akimoto et al., 2014; Puljung et al., 2014; Saponaro et al., 2014; Goldschen-Ohm et al., 2016), recently also structures of HCN channels comprising the majority of the TMPC and the CPC became available. Currently, structures are known for HCN1 (Lee and MacKinnon, 2017; Lee and MacKinnon, 2019; Kim et al., 2024; Burtscher et al., 2024), HCN3 (Yu et al., 2024), and HCN4 (Saponaro et al., 2021), with and without cAMP in the ligand-binding site. Structures are available not only in the closed pore but also in the open pore state (Saponaro et al., 2021; Burtscher et al., 2024) as well as in states that resemble the voltage-activated state of the channel (Lee and MacKinnon, 2019; Burtscher et al., 2024). A list of the currently available HCN channel structures in the data bank is provided in Table 1.

Based on these high-resolution structural data, many functional features that are specific for HCN channels, like their low unitary conductance (DiFrancesco, 1986; Benndorf et al., 2025), their moderate K⁺/Na⁺ selectivity (Moroni et al., 2000; Lyashchenko and Tibbs, 2008), and their susceptibility to block by Cs⁺ (DiFrancesco, 1982) and ivabradine (Bucchi et al., 2002), could be analyzed and understood at an atomistic level by MD simulations (Saponaro et al., 2021; Saponaro et al., 2024; Bauer et al., 2022; Krumbach et al., 2023; Benndorf et al., 2025).

While these simulations can cover dynamic events in a protein in time frames of nanosecond to maximal microsecond, they are not suitable to tackle the problem of slow structural transitions (Bahar et al., 2010). Such slow conformational changes with dwell times in the range hundreds of microseconds or even milliseconds, however, are dominating voltage and ligand-mediated gating in HCN channels. In addition to advanced time-resolved structural analysis of different conformational states in membrane transport proteins (Lee et al., 2025; Morgan et al., 2025), computational methods are also able to provide insights into such slow conformational changes in proteins (Shrivastava and Bahar, 2006; Hamacher and McCammon, 2006; Eom et al., 2007; Chennubhotla and Bahar, 2007; Bahar et al., 2010; Das et al., 2014; Hoffgaard et al., 2015; Zhao et al., 2024). Available data support the view that the conformational changes and allosteric mechanisms are intrinsically determined by the native fold of proteins, see e.g., Hamacher et al. (2006). The reason for this is that the geometry and topology of a protein under equilibrium conditions restrict its freedom of movement. Due to these inherent structural limitations, the fluctuations in proteins, which range from small local fluctuations to large collective movements of different domains or entire channel subunits in the gating process, are limited.

Here, we summarize recent attempts in understanding structure/function correlates of cAMP-mediated gating in HCN channels by means of such coarse-grained approaches. HCN channels are ideal systems for these kinds of studies because high-resolution structures are available for distinct functional states of the protein, like with and without cAMP bound. Also structures with the voltage sensor in the nonactivated and activated conformation are available (Table 1). In this sense, the starting and end points of conformational transitions in HCN channel gating are well-defined by experimental data, and coarse-grained computational approaches can provide insights into the mechanisms and into the mechanics for the transition between these states. We show here that elastic network models (ENMs) like Gaussian network models (GNMs) and similar models (Haliloglu et al., 1997; Tirion, 1996), and in particular their extension to anisotropic network models (ANMs) (Atilgan et al., 2001), were, in combination with linear response theory (LRT) (Ikeguchi et al., 2005), able to reveal the mechanic coupling between domains in these proteins. This analysis provides insight into the dynamics of allosteric changes in these channels upon cAMP binding.

The analysis of ENMs is based on the finding that large dynamic motions in a protein are determined by geometric constrains inherent in the three-dimensional architecture of a polypeptide. These motions, which are dictated simply by the global fold of a protein, are presumably the basis for functional properties such as ligand binding or sensing of physical cues like pressure or voltage. A coarse-grained representation of the contact topology in a protein of interest in an ENM (Bahar et al., 2010) should, together with LRT (Ikeguchi et al., 2005), provide the information of such global motions in relation to function. In combination with their low computational costs both methods are ideal for understanding the mechanical interactions and long-range conformational changes in complex proteins. With these features both approaches are particularly interesting for understanding interactions in HCN channels between the cytosolic CNBD and the distant membrane-resident VSD (Fig. 1).

The native state of a protein is described by the ensemble of its conformational states. These are constantly fluctuating around an equilibrium and all together sample the native fold of a protein (Atilgan et al., 2001; Hamacher and McCammon, 2006; Lezon et al., 2009; Bahar et al., 2010; Mahajan and Sanejouand, 2015). These fluctuations are not only accounting for the passive motions in a protein but also for functionally relevant movement of domains such as those involved in ligand-binding, allosteric regulation, and cooperativity (Lezon et al., 2009; Bahar et al., 2010). In the specific case of HCN channels this could, for example, be complex movements such as voltage-dependent gating, which requires highly collective global motions of the voltage sensor relative to the pore. Such global conformational changes of a protein strongly correlate to its normal modes (Lezon et al., 2009; Bahar et al., 2010). A protein has many modes of which each one represents a different pattern of atomic displacements characterized by a distinct frequency with which it oscillates back and forth (Mahajan and Sanejouand, 2015). Most relevant for the protein’s mechanical functionality are the modes with the lowest frequencies. They are the most interesting for the investigation of protein dynamics since these slow vibrations also generate the largest amplitudes (Atilgan et al., 2001; Hayward and de Groot, 2008; Bahar et al., 2010; Kondrashov, 2011; Sanejouand, 2013). It has been shown that they are responsible for 90–95% of atomic displacements in proteins (Sanejouand, 2013). The relevant normal modes, which are dominating global changes in a protein, can be examined with an ENM (Atilgan et al., 2001; Bahar et al., 2010). The starting point of every ENM is the native fold of a protein obtained from X-ray crystallography, cryo-EM, or NMR. The high-resolution structures with all atomistic coordinates, which are currently available for HCN channels in different functional states, are listed in Table 1. These structures can be used to build coarse-grained models in which each residue is represented by a single node, usually coincident with the C_α atom. All nodes are connected by a spring, if their distance is smaller than a cutoff distance (Hayward and de Groot, 2008; Kondrashov, 2011; Sanejouand, 2013). Therefore, the model does not distinguish between covalent and non-covalent interactions but purely on distance between two nodes.

Fig. 2 shows an example for the transformation of the crystal structure of the CNBD from HCN4 (PDB accession no. 3U11) into an ENM. Each amino acid in the polypeptide is represented by a bead (=node) centered on its C_α atom. All nodes are connected by springs when their distance is smaller than 13 Å, in this case. Backbone connections are displayed in red just for visualization purposes.

To model the interaction between two nodes, e.g., to describe its spring-like character, Newton’s equations of motion have to be solved. Since this is analytically impossible, they need to be approximated, which can be done with a Taylor expansion if the displacements in the vicinity of the equilibrium are small enough. The constant term can be neglected as only differences in the potential energy and not absolute values are of interest. The linear term vanishes too since local minima require the first derivative to be zero. By construction the equilibrium/native state is such a minimum. That leaves the quadratic term as higher derivatives result in even smaller contributions to the potential energy and are thus ignored, too (Sanejouand, 2013). This implies that an ENM is the most general model of biomolecular mechanics—given that the interaction strengths (spring constants) are well chosen.

Generally, there are two classical approaches for using ENMs, which are either GNM or ANMs (Atilgan et al., 2001; Bahar et al., 2010; Sanejouand, 2013). The latter additionally allows to compute the direction of movement.

The interaction potential (1) is the modeling part. In a second step, one has to understand the dynamics induced by it. Typically, this is done by extensive MD simulations. ENMs are, however, explicitly constructed to avoid the need for such time-intensive simulations. Instead, we develop the potential in a Taylor series of polynomials of all coordinates of the amino acids. The first and constant term is uninformative as the physical effects manifest themselves in energy differences, where the first term would cancel itself. The first-order term vanishes by construction as the input structure is assumed to be in the minimal energy state and in a compact domain the first derivative vanishes at extremal points. This leads us to the second-order term that mimics harmonic, spring-like interactions.

This is calculated for all interacting atoms/nodes resulting in an NxN matrix consisting of 3 × 3 sub-matrices. Taken together, this results in a 3N × 3N shaped matrix called the Hessian. From this, the singular values and vectors can be computed using a singular value decomposition algorithm (Press et al., 2007). In total, the protein has 3N singular values with the first six values vanishing as they correspond to the three translational and three rotational degrees of freedom (Wilson, 1941). These six values are not of interest since the goal of the analysis is to investigate the conformational fluctuations within the molecule rather than translation and rotation of the whole molecule. Of the remaining 3N-6 singular values, typically the first values, namely the lowest ones, correspond to low-frequency motions. They represent the global structural changes with higher amplitudes, which were shown to be the functionally relevant ones (Hayward and de Groot, 2008; Sanejouand, 2013).

An extension of the ANM, the extended ANM can be employed (Hamacher and McCammon, 2006). Instead of a uniform force constant γ, a distinct term for each contact potential has been included in this network model, which pays respect to the type of amino acid but also to the type of interaction between two nodes. The latter is described by two features, the bond type (covalent vs. non-covalent) and the relative location of interacting nodes (interchain vs. intra-chain). In this way it is possible to weight, for example, covalent interchain interactions or disulfide bonds differently from weak intra-chain interactions. They also distinguish between the difference in hydrophobic force that is experienced by intra-chain interactions versus inter-chain interactions (Hamacher and McCammon, 2006).

The LRT assumes that an external perturbation of a protein will evoke a linear response in its structure. This means that the protein will undergo a conformational change proportional to the strength of the perturbation. Examples for perturbations in an HCN channel are (1) the impact of an electrical field on the voltage sensor domain or (2) the binding of a ligand to a specific binding site. Therefore, it should be possible to predict the conformational transition in HCN channels between the apo- and a holo-form from the analysis of the normal modes. This information will again be contained in the Hessian of an appropriate ANM.

A simple prediction from this model is that binding of a ligand exerts distinct forces on the target protein by interacting with critical residues in the binding site (Fig. 3). Hence, in a computational approach, it should be possible to mimic in HCN channels the impact of cAMP binding to the CNBD on the remaining protein by application of forces to these critical sites in an ENM (Fig. 3). For several proteins, it was already confirmed that these predicted displacements by LRT were, despite the simplicity of this approach, in very good agreement with experimental X-ray and cryo-EM data as well as all-atom MD simulations (Hamacher et al., 2006; Hamacher, 2011; Groß et al., 2017; Groß et al., 2018). For the present review, it is worth mentioning that the accuracy of conformational changes by LRT is even valid for membrane proteins where the membrane was neglected (Das et al., 2014).

In classical LRT, forces are applied to residues in a protein with the intention of mimicking known perturbations to the protein. To contextualize the displacements given by an LRT in an unbiased way, a “null model” of perturbation had been developed (Groß et al., 2017; Groß et al., 2018). Its purpose is to test if the chosen force vector with the directionality of interest is the only one producing particular conformational changes. In this procedure, forces are applied to the residue or residues of interest in all possible directions (Fig. 3). In the case of HCN channels, this could, for example, be the cAMP-binding site. The consequence of these forces on the conformation of the target protein like the conformation of a ligand-binding site are then recorded and clustered according to their directions. Whenever the shifts caused by forces with a different direction than the one of interest cause shifts other than those expected, the LRT result is relevant to the mechanics of the molecule.

Finally, the predicted conformations are compared with experimentally known conformational changes in the protein. Fig. 3 i shows the movement following ligand binding to a hypothetical protein like the CNBD in HCN channels. Fig. 3 ii shows the LRT simulation of ligand binding, producing the same movement. In the null LRT model, shown in Fig. 3 iii, the application of a distinct force, which is represented by the red vector, shifts the ligand-binding domain from an open into a closed conformation, while other forces, represented by green, blue, and light blue vectors, do not. This shows, in turn, that the process of ligand binding will generate a conformational change in the direction of the red vector. In this case, it is possible to evaluate the flow of conformational information over long distances in a protein; in HCN channels, this could be the flow of information between CNBD and VSD. Such a null model approach was first developed for a monomeric protein (Groß et al., 2017). So, to study the conformational changes in case of a tetrameric channel like HCN, it has to be applied to each subunit simultaneously or to individual monomers when examining potential cooperativity (Groß et al., 2018; Kunzmann et al., 2024).

It is important to note that reduced molecular models such as ENMs cannot reproduce the full dynamics of a protein in detail. For example, the assumption of purely harmonic interactions (“springs”) prevents the application of ENMs to large conformational changes, (Miyashita et al., 2003) like, in particular, the “breaking” of contacts. Such a breaking arises, e.g., in large conformational changes so that residues move further apart than the threshold distance upon which in the model no interaction is possible anymore.

First of all, such larger movements might bring the structure into a state where the Taylor expansion to second order (harmonic interactions) is insufficient. Incorporating additional terms, on the other hand, hinders the efficient solution via the singular value decomposition. The respective equations of motions would need to be solved just as traditional MD programs do.

Secondly, breaking of contacts cannot be investigated, as this changes the Kirchhoff and the Hessian matrix, respectively. While some work has be done to merge two or more distinct ENMs—capturing the dynamics of some open or closed state, for example—the overall insight into the transition state remains illusive (Miyashita et al., 2003).

The same applies to the application of LRT: scenarios in which ENM and LRT can give insight into the dynamics of proteins are restricted to movements around the native fold. However, this could still provide clues to the precursors of larger conformational changes.

The challenge of understanding allosteric gating in HCN channels is that the VSDs as well as the gates, which determine the voltage-sensitive opening of the channel, are in the TMPC. The ligand-binding domains on the other hand, which modulate the voltage-dependent channel opening, are in the cytosolic domain, e.g., in molecular terms far away from the VSD and the gates (Fig. 1). The key question is how conformational changes in the ligand-binding domain are transmitted to VSD and/or channel gates for a modulation of voltage dependency.

Even before high-resolution structures of HCN channels, including both TMPC and CPC, were available, there was already a bulk of structural information on their isolated CNBDs in holo and apo state (Zagotta et al., 2003; Lolicato et al., 2011; Saponaro et al., 2014). These studies had uncovered the cAMP-binding site as well as the amino acid side chains, which mediate ligand binding in the CNBDs of different HCN channels. Analysis of the cAMP-bound versus -unbound structures had furthermore identified complex rigid body movements of helical domains in the CNBD together with a stabilization of several of the helices in response to binding and unbinding of cAMP (Saponaro et al., 2014). To obtain information on how such ligand-induced local conformational changes in the CNBD are mechanically transmitted to the channel domain in the membrane, Weißgraeber and co-workers (Weißgraeber et al., 2017) performed a coarse-grained molecular analysis on a homology model, comprising the CNBD of HCN4 and the channel domain of the Kv1.2 channel. This was motivated by the idea that such a homology model was sufficiently accurate to provide, in the absence of structural details on the HCN transmembrane part, information on basic mechanical connections between the ligand binding and the channel domain. To this end, an ENM was constructed from the Kv1.2/HCN4 chimera, and the LRT approach (Ikeguchi et al., 2005) was used to mimic the local conformational changes in the CNBD, which are elicited by unbinding of cAMP from its binding site. For this purpose, an external force was applied in each of the four monomers to the six residues in the binding site, which interact with cAMP. This force was tailored such that it caused a distortion of the cAMP-binding site, which mimics the experimentally determined conformation of the unbound state (Saponaro et al., 2014).

Comparative analysis of the ENM with and without this external force revealed interesting short- and long-range movements in the channel protein in response to the local conformational changes in the cAMP-binding site. These conformational changes were not confined to the CNBD but also transmitted via strong rotational movements into the TMPC; notably, only the SF occurred structurally, unaffected by the release of cAMP from its binding site (Weißgraeber et al., 2017).

The data also highlight interesting bidirectional interactions between the C-linker of the CNBD with the cytosolic short linker that connects the transmembrane domains S4 and S5. The most striking observations were torsion like movement of the C-linker away from the TMPC in response to release of cAMP from the binding site (Fig. 4 and Video 1). With the C-linker as a mechanical transmitter between the CNBD and TMPC, a torsion-like movement of the cytosolic domain is than translated into a twist-like motion at the inner gate of the channel (Weißgraeber et al., 2017). Because such an iris-like motion has already been causally linked to the gating in other K⁺ channels (Marchesi et al., 2018), the computational data predict that the same kind of mechanical connection between cAMP-binding site and channel is also part of the allosteric mechanism in HCN channels which couples cAMP- and voltage-dependent gating.

These predictions were made based on a crude homology model and may therefore not hold true for the real structure of an HCN channel. But when the first high-resolution cryo-EM structures of HCN1 in the apo- and holo-conformation (Lee and MacKinnon, 2017) were published, it turned out that the predictions from the ENM/LRT analysis perfectly met the structural findings. The model predicted iris-like movement of the C-linker (Fig. 4 and Video 1) was indeed confirmed by an inspection of the structures of this domain in the presence and absence of cAMP (Fig. 5 C) (Lee and MacKinnon, 2017; James and Zagotta, 2018). The linker of the holo- and the apo-structure exhibit the same rightward iris-like rotation that was predicted from the LRT analysis (Weißgraeber et al., 2017). This good agreement between experimental and computational data underscores the high information content of the coarse-grained approach. It also shows that mechanical connections, which are important for channel gating in HCN channels, are presumably inherent features of the global structure of these proteins and independent of structural details.

The recent availability of high-resolution structures of HCN1 channel with TMPC and CPC (Table 1) paved the way for a more detailed and more realistic analysis of their structure/function correlates. To tackle the mechanism of allosteric gating, Groß and co-workers (Groß et al., 2018) applied the perturbing force in the ENM of the HCN1 channel, the perturbating force of LRT method, not in the cAMP-binding site but further upstream in a single position on the so-called elbow (Zagotta et al., 2003; James and Zagotta, 2018), a critical site in the bend between the A′ and B′ helices of the C-linker. It is well established that the C-linkers are packed in a channel tetramer in an elbow on the shoulder configuration in which the A′ and B′ helix from one subunit interact with the C′ and D′ helix of the adjacent subunit. Hence, the elbow is a key strategic position for imposing a conformational change in the C-linker, e.g., the mechanical element, which connects the CNBD to the TMPC (Zagotta et al., 2003; James and Zagotta, 2018). With this architecture, a realistic perturbation of the elbow by LRT should mimic binding/release of cAMP to/from the binding site; in the domain that connects the CNBD with the membrane-embedded part of the protein, this system should report the respective conformational changes transmitted to the TMPC. A comparison of the structural transition between HCN1 apo- and holo-structure with a computational perturbation of the elbow horizontally and away from the central axis shows the same rotational motion of the CNBD (Fig. 5). This again highlights the very good agreement of computational and experimental data. As a further quality confirmation of the approach, it was also observed that the local application of a force on the elbow, which mimics the transition from the holo-to-the-apo-configuration, also affects the conformation of the cAMP-binding site. As expected for a reciprocal system, the distortion in the downstream direction of the elbow causes an “apo like” conformation of the cAMP-binding site.

These findings confirm that the tip of the elbow is indeed a critical position for the mechanical transduction of information from the CNBD to the pore and vice versa. This good agreement between computational data and experimental structures as well as the reciprocal nature of conformational changes in the CNBD/C-linker domain are a strong indication that also the conformational changes in the TMPC are relevant for understanding the impact of the CNBD on voltage-dependent gating.

A close view on the computational data advocates the following model: The membrane embedded TMPC is constrained in its movement particular in the z direction by the lipid bilayer. CNBD and C-linker on the other hand are free to move in all directions. In the transition between apo- and holo-state the CNBD undergoes an anticlockwise rotation with respect to the PD, bringing the CNBD closer to the membrane. This presumably favors direct contacts between the CNBD and residues in the TMPC, including the S4 domain in the VSD (Weißgraeber et al., 2017). Because of the different degrees of freedom the iris-like motion of the CNBD is also transmitted further to the inner gate, which is in direct contact to the C-linker. This process is presumably guided by simple mechanical forces.

Important to mention here is that these predicted complex conformational changes, including the iris-like motion of the CNBD and the change in the distance between the latter domain and TMPC were independently confirmed by high-speed atomic force microscopy studies on the related cAMP gated SthK channel (Marchesi et al., 2018). The latter study reports the same kind of a ligand-triggered iris diaphragm mechanism of C-linker motion that is predicted by the ENM/LRT analysis of HCN structures (Weißgraeber et al., 2017; Groß et al., 2017). We can take this as additional evidence for the basic assumption of this review namely that the three-dimensional geometry and topology in the protein under equilibrium conditions determines the direction and degree of freedom in movements.

It is well established that cAMP binding is insufficient to open HCN channels in the absence of movement of the VSD by negative voltage (Flynn and Zagotta, 2018). Hence, the impact of the iris-like motion of the CNBD following ligand binding on the TMPC bears no direct information on the primary mechanism of HCN channel opening. Nonetheless, we can expect that the induced conformational changes in the TMPC contains information on how ligand binding facilitates voltage-dependent channel opening.

One important information from the LRT data was that all transmembrane domains do not undergo vertical movements but rather distinct tilting motions. This finding predicted a scenario, which was later confirmed by experiments (Lee and MacKinnon, 2019), namely that the S4 domain is not undergoing a rigid-vertical movement in the direction of the electrical field but a tilting motion coupled with a downward movement of the S4 helix, approximately two helical turns relative to the charge transfer center. The same movements, though less pronounced, were later confirmed in experimental structures (Burtscher et al., 2024). Tilting of S4 and the concerted movements of S5 in the opposite direction suggest a scenario in which the ligand-binding status in the CNBD can favor the formation or collapse of a water-filled internal “gating” canal by long distance interactions. Previous work had highlighted the importance of this gating canal between S4 and S5 for voltage-dependent opening of HCN channels (Bell et al., 2004). It was postulated that the shape of this aqueous gating canal would affect the membrane field surrounding the S4 segment and influence in this manner the voltage dependency of the channel (Bell et al., 2004). Hence, the CNBD could have in this way an indirect impact on voltage gating; ligand binding could favor widening of the aqueous gating canal and, in this manner, lower the energy for voltage-dependent channel opening.

A further positive impact of ligand binding on HCN channel opening can be deduced from the effect of the iris type motion of the CNBD at the exit from the channel’s cavity. From structural and experimental data, it had been speculated that HCN channels have an inner gate at the entry from the cytosol into the cavity (Rothberg et al., 2003). This gate was thought to be narrow at depolarized voltages and widened upon hyperpolarization, allowing the passage of ions. A recent cryo-EM structure has confirmed this picture, showing that a small widening of this constriction at the entry of the HCN4 channel is sufficient for transforming the closed channel into a conductive open channel (Saponaro et al., 2021). In MD simulations, this protein exhibits all major functional properties of the HCN4 channel known from experimental studies (Saponaro et al., 2021; Saponaro et al., 2024; Bauer et al., 2022). After determining that cAMP binding and its lateral pushing force on the C-linker elbow also affect conformations in the TMPC, Groß and co-workers (Groß et al., 2018) examined the consequent changes on the radius of the inner gate. This showed that the mechanical motions of the CNBD following ligand binding were indeed favoring a widening of this part of the inner gate. Perturbations in all other directions caused a narrowing of the inner gate. Since cAMP binding by itself is not sufficient to open HCN channels, these data must be interpreted again in the context of a stabilizing effect of cAMP on the open state of the channel. A similar interpretation for the role of cAMP in HCN gating was already proposed based on experimental data (Craven and Zagotta, 2004; Craven et al., 2008; James and Zagotta, 2018).

The first HCN1 structure obtained by cryo-EM (Lee and Mackinnon, 2017), highlighted the presence of a folded domain comprising three helices (⍺, β, and γ) in their cytosolic N terminus (Fig. 1). The structural organization of this domain, now known as HCND, suggested that it serves for a mechanical coupling between the C-linker and the VSD in these channels. Functional data on voltage- and cAMP-dependent gating, in combination with mutations, proposed that the HCND presumably serves as a sliding crank that transmits the movement between the C-linker and the VSD, once the C-linker is rotating in response to cAMP binding to the CNBD. The β helix of the HCND is on one end connected via salt bridges to the B′ and A′ helices of the C-linker. On the other end, the β helix forms a hydrophobic interaction with the lower end of the S4 domain via a phenylalanine residue on HCND that is inserted in a hydrophobic pocket between S1 and S4. With these connections, the HCND converts the planar iris-like rotational movement of the C-linker into a tilting displacement of the VSD (Porro et al., 2019).

The combination of ENMs and LRTs served as a perfect tool for testing these predictions. For this purpose, the same force, which mimics the aforementioned conformational changes induced by cAMP binding (Groß et al., 2018), was applied to the channel protein. This resulted in a concerted displacement of the C-linker together with the HCND with the effect that both domains rotate with a similar angle to the z-axis. This rotational movement of the C-linker and the HCND was then transmitted to the VSD, where it causes a distinct tilting of the S4 (Porro et al., 2019). The coupled motions of C-linker, HCND and S4 derived from the LRT analysis are illustrated in Fig. 6, A and B.

In the context of HCN gating, the prediction that cAMP binding favors tilting of the S4 domain is most interesting. Notably, tilting with a kink formation of the S4 domain was also proposed in experimental and computational studies as a crucial step in the voltage-dependent activation of HCN channels (Kasimova et al., 2019; Lee and MacKinnon, 2019; Elbahnsi et al., 2023). Fig. 6 C, i–iv, show the S4 domain from cryo-EM structures of HCN1 in the depolarized close and in a presumed hyperpolarized-like pre-open state, obtained by forcing the voltage sensor in the downward conformation with covalent bonds (Lee and Mackinnon, 2019). Within major global structural changes between the two conformations, it occurs that the S4 domain is nearly linear in the depolarized state but kinked in the hyperpolarized state; amino acid Ser272 serves as hinge for a bending of the lower part of the helix away from the pore axis. Fig. 6 C, iii and iv, illustrate the predicted conformations of the same domain from an LRT analysis of HCN1 before and after application of a force on the CNBD, which mimics the effect of cAMP binding (Groß et al., 2018; Porro et al., 2019). An overlay of the S4 domain from the apo-like and the holo-like conformation, in which the two S4 domains are normalized to the linear outward facing part of the α-helix, predicts a similar picture: S4 is approximately linear in the apo form and reveals a kink with S272 as a hinge in the holo-conformation. While both negative voltage and ligand binding seem to favor kink formation in the S4 domain around a Serine hinge, they presumably achieve this in different manners. The force vectors on the S4 domain, which are elicited by release of cAMP from its biding site, are stronger on the upper part of the S4 helix, favoring a bending of this part toward the pore axis (Fig. 6 C). Negative voltage apparently promotes more a bending of the lower part of the helix away from the pore.

Fig. 6 C iv depicts only the predicted tilting/bending movement of S4. Within the global protein, the LRT analysis predicts that this tilting motion also evokes, just like in the experimental structures of HCN1 (Lee and MacKinnon, 2019; Burtscher et al., 2024) and HCN4 (Saponaro et al., 2021), a directional displacement of helices S4, S5, and S6 with a concomitant outward movements of the C-terminal half of S4 and S5 (Groß et al., 2018). Moreover, the LRT prediction implies that this S4–S5 movement is linked to a rotation of S6 and, consequently, a widening of the pore entrance (Groß et al., 2018). The same S6 rotation/pore widening was later experimentally shown both in HCN1 (Burtscher et al., 2024) and in HCN4 (Saponaro et al., 2021) as a part of pore opening.

The finding that negative voltage and cAMP binding seem to favor a tilting movement of the S4 domain (Fig. 6 C) bears a potential explanation for the mechanisms by which ligand binding reduces the energy barrier for voltage-dependent opening of HCN channels. The data suggest that the energy for imposing this crucial conformational change for channel opening may be lower when voltage and ligand binding work in a concerted manner on imposing this tilting of S4 than in the case of voltage alone. This interpretation is consistent with the experimental data according to which cAMP binding cannot open the channel but favors voltage-dependent channel activation (James and Zagotta, 2018).

The beauty of coarse-grained approaches with their low requirement for computational resources is that they allow an easy testing of model predictions. The above-mentioned data imply that C-linker, HCND, and VSD form a mechanical continuum, and that conformational changes in the C-linker can be transmitted in this manner via the HCND to the VSD. If this is true, the same flow of information should also occur in the opposite direction, e.g., from the VSD to the C-linker and presumably to the CNBD. A test of this prediction confirmed this assumption in that the critical rotational movement in the protein is elicited irrespectively on whether force is applied to the C-linker, the HCND, or the S4 helix (Fig. 6 D i–iii). Collectively, these data confirm the central position of the HCND as a mechanical connection, which couples in a bevel gear manner the CNBD with the VSD. A recent computational study that examined by elastic network analysis the mechanical connectivity between HCND and VSD in channels suggests that the flexibility of this connection changed during evolution. A high degree of flexibility is apparently crucial for an allosteric modulation of voltage gating by cAMP with the effect that the HCN isoforms 1, 2, and 4 are sensitive to this second messenger, while HCN3 and spHCN are not (Alvarez-Villagómez and Balleza, 2025, Preprint).

This view of HCN channels as a distinct mechanical continuum in which individual domains can communicate via long-range interactions connections also has implications for the development of drugs. By exploiting the knowledge of long-range interactions in these proteins, it should even be possible to modulate the impact of cAMP binding by small molecules, which bind on the extracellular side to HCN channel in sites that are mechanical affected by ligand binding.

The structures of tetrameric HCN channels exhibit an apparent fourfold symmetry in the TMPC as well as in their cytosolic CNBDs (Lee and MacKinnon, 2017; Saponaro et al., 2021). Despite this structural symmetry, functional data suggest an explicit asymmetric behavior in both cAMP binding and channel gating. For example, gating of HCN2 concatemers, which carry between 0 and 4 mutations in the cAMP-binding site, is best explained by a cooperativity between the ligand-binding sites (Liu et al., 1998; Ulens and Siegelbaum, 2003). The data advocate a scenario in which cAMP binding favors the formation of CNBD dimers and that ligand-mediated full channel opening is only occurring after a pairing of these dimers. More evidence for cooperativity in cAMP binding in HCN channel gating was derived from patch-clamp fluorometric recordings, a technique that monitors binding of a fluorescent cAMP analog to their binding site while simultaneously recording channel activation (Kusch et al., 2012). These data suggest that cAMP binds HCN2 with different modes of cooperativity depending on whether the channel is in a voltage-dependent closed or open state (Kusch et al., 2012; Thon et al., 2015). The most complex and physiological relevant mode of allostery occurs in the voltage-activated channel in which binding of a cAMP molecule to one binding site modulates the affinity of the remaining empty sites. Experimental data for this voltage-activated channel advocate a complex allosteric model for sequential binding with a positive-negative-positive cooperativity. It predicts that the affinity for binding of a second and fourth cAMP ligand is enhanced by a preceding binding step while binding of the third ligand molecule is penalized (Kusch et al., 2012).

The general importance of cooperativity in ligand binding to HCN channels in their closed state was later questioned by single-molecule studies in which the binding dynamics of cAMP to individual detergent-solubilized HCN1 and HNC2 channels was directly monitored at single protein level (White et al., 2021). In this study, the ligand molecules were found to bind independently to all four subunits, thus questioning a molecular basis for cooperative cAMP-mediated HCN channel gating. In contrast to that, subsequent measurements also with single-molecule resolution on HCN channels either in lipid nanodiscs (Idikuda et al., 2024, Preprint) or in isolated membrane patches (Kuschke et al., 2024) reported in contrast a distinct positive cooperativity. The affinity of each empty binding site occurred to increase as a function of a preceding binding step. Taken together, these data underpin that HCN channel proteins have an inherent tendency for a cooperative ligand binding; this allostery seems to be affected by membrane voltage and by constrains of the lipid bilayer on the channel protein (Idikuda et al., 2024, Preprint).

To further examine the question of cooperativity in cAMP-mediated gating of HCN channels with an independent computational approach, the ENM analysis of HCN in combination with LRT described before was employed (Kunzmann et al., 2024). The special twist of this study was that the force on the cAMP-binding site, which mimics the conformational change in cAMP binding, was not applied to all four subunits but only to individual binding sites in the tetramer (Fig. 7 A). It was reasoned that a local impact of one subunit should have no consequences on the remaining subunits in the case that ligand binding is noncooperative. If, however, the protein exhibits cooperativity between the ligand-binding sites, one should detect an impact of one subunit on the others.

The systematic analysis of the LRT data on the cryo-EM structures of HCN4 and HCN1 strongly supports a cooperativity in cAMP binding in these channels (Kunzmann et al., 2024). Stepwise binding generated in these channels, which were not in the voltage-activated state, essentially the same sequence of positive, negative, and positive cooperativity in the channel, which was observed in patch-clamp fluorometric recordings of the voltage-activated HCN2 channel (Fig. 7, B and C). The computational data further imply a scenario in which the initial cAMP-binding step increases the affinity for binding in the subunit diagonally opposite to the first binding monomer while imposing a negative effect on the adjacent monomers (Fig. 7 B). Unlike expected from the experimental data (Idikuda et al., 2024, Preprint), this distinct pattern of cooperativity seemed to be unaffected by the presence or absence of a membrane around a protein (Kunzmann et al., 2024). Such a potential constraining impact of a membrane round a protein can be considered in the vibrational analysis of the protein by an explicit ANM for a bilayer (Lezon and Bahar, 2012; Lomize et al., 2022; Kunzmann et al., 2024).

To answer the question if this kind of cooperative behavior is a specific feature of HCN channels, the same analysis was performed with structurally related channels with active or inactive cAMP-binding sites. The results showed that also the cyclic nucleotide gated (CNG) channel TAX4 and the Kv channel KAT1 exhibit a similar but much less pronounced pattern of cooperativity between the ligand-binding sites. Interesting to mention in this context is that an independent anisotropic network analysis of the human CNG channel (CNGA3) (Gofman et al., 2014) found the same positive and negative interplay between CNBDs as in the analysis with LRT (Kunzmann et al., 2024). Both studies propose that the pairs of diagonal CNBDs move in an alternating manner toward and away from each other.

Taken together, the data suggest that cooperativity between subunits is an inherent property of channel proteins with the same overall architecture. In the case of HCN channels, this cooperative interaction between ligand-binding sites was presumably optimizing during evolution to allow the fine tuning of cAMP-mediated gating in these channels. Inherent structural constrains in the HCN channel tetramer presumably favor a collective alternating movement of diagonal CNBDs toward and away from each other in the voltage-activated state of the channel. In the next step, it would be interesting to examine how the voltage sensor domain in the active and inactive state affects the geometry and topology of the CNBD.

This review underscores that coarse grained computational methods are, in combination with high-resolution structures of these proteins, able to distill relevant mechanical connections and conformational movements in the complex HCN channels. The combination of ENMs with LRT has predicted a number of conformational movements in these channels in the context of cAMP binding that were either confirmed later by experimental structures or supported by functional data from HCN channels. This underlines that the respective methods are, in spite of their simplicity, suitable for understanding macroscopic, long-range interactions in these channels like the modulation of voltage gating by ligand binding. The good agreement between experimental structural and functional data on one side and model predictions from coarse grained simulations on the other side provide some important information: because ENM/LRT simulations are not considering any detailed chemical interactions in proteins, the geometry and topology of HCN channels must be such that small local fluctuations and large collective movements of different domains are an inherent “material” property in these proteins; under equilibrium conditions, the degree of freedom of movements is to some extent limited and determined by their global architecture. Critical interactions between amino acids in these proteins like salt bridges, hydrogen bonds, or hydrophobic interactions, which were found important in cAMP-mediated gating of HCN channels (e.g., Porro et al., 2019; Porro et al., 2020; Porro et al., 2024; Kondapuram et al., 2022) must be understood as interactions, which enforce and physically stabilize conformational movements, which are already inherent in these proteins.

The systematic analysis of conformational changes in HCN channels following binding of cAMP to its binding site in the CNBD is summarized in Fig. 8, A and B. Ligand binding promotes a sideward movement of the CNBD (arrows a). This generates an iris-like, anticlockwise rotation of the C-linker (arrows b). The rotational movement of the C-linker pulls the HCND and the lower part of the S4 domain along the tangential axis of the rotation (arrow c). The consequence is a tilting of the transmembrane domains, including the S4 domain (arrows d) relative to the z-axis of the lipid bilayer. All these movements are reversed by a loss of cAMP from its binding site.

The conformational changes in HCN channels detected by the combined ENM/LRT analysis have some implications for cAMP-mediated gating, which are summarized in Fig. 8 A: (1) ligand binding in one monomer has positive and negative impacts on binding domains in the remaining monomers, which can explain the complex pattern of cAMP-mediated cooperativity in HCN gating. (2) The C-linker approaches S4, and (3) the rotational movement of the C-linker favors, via the HCND, a tilting of the S4 domain. (4) Tilting of S4 away from S5 promotes opening of a “gating canal,” an aqueous structure which is postulated to modulate the electrical properties of the electrical field surrounding S4. (5) The conformational transitions in the protein promote opening of the inner gate. Important to mention at this point is that cAMP binding is by itself not able to open HCN channels; it only decreases the energy barrier for an opening of the channels by negative voltage. Hence, within the limitations of the coarse-grained method, which is not observing actual movements in the protein, we must interpret these data in a sense that cAMP binding lowers the energy for the real physical movement of the S4 domain and the opening of the inner gate by negative voltage.

Crina M. Nimigean served as editor.

We thank Dr. S. Weißgraeber (TU Darmstadt) for generating the video in Fig. 4. The authors are grateful to anonymous reviewers who provided valuable insight and thereby improved the manuscript.

This work was funded by Deutsche Forschungsgemeinschaft (grants HA5261/6-1 to K. Hamacher and TH558/34-1 to G. Thiel), a Foundation Leducq Research Grant (TNF FANTASY- 19CVD03 to A. Moroni and D. DiFrancesco) and by the Ministero Della Ricerca (MUR) Progetti di Rilevante Interesse Nationale (PRIN) (2022EMA8FA to A. Moroni).

Author contributions: Magnus Behringer: visualization and writing—original draft, review, and editing. Jan Krumbach: visualization and writing—original draft. Alessandro Porro: visualization and writing—original draft. Andrea Saponaro: writing—review and editing. Dario DiFrancesco: writing—review and editing. Anna Moroni: funding acquisition and writing—review and editing. Kay Hamacher: conceptualization, formal analysis, funding acquisition, methodology, project administration, software, supervision, and writing—original draft, review, and editing. Gerhard Thiel: conceptualization and writing—original draft, review, and editing.