Abstract
Seemingly unrelated traits often share the same underlying molecular mechanisms, potentially generating a pleiotropic relationship whereby selection shaping one trait can simultaneously compromise another. While such functional trade-offs are expected to influence evolutionary outcomes, their actual relevance in nature is masked by obscure links between genotype, phenotype, and fitness. Here, we describe functional trade-offs that likely govern a key adaptation and coevolutionary dynamics in a predator-prey system. Several garter snake (Thamnophis spp.) populations have evolved resistance to tetrodotoxin (TTX), a potent chemical defense in their prey, toxic newts (Taricha spp.). Snakes achieve TTX resistance through mutations occurring at toxin-binding sites in the pore of snake skeletal muscle voltage-gated sodium channels (NaV1.4). We hypothesized that these mutations impair basic NaV functions, producing molecular trade-offs that should ultimately scale up to compromised organismal performance. We investigate biophysical costs in two snake species with unique and independently evolved mutations that confer TTX resistance. We show electrophysiological evidence that skeletal muscle sodium channels encoded by toxin-resistant alleles are functionally compromised. Furthermore, skeletal muscles from snakes with resistance genotypes exhibit reduced mechanical performance. Lastly, modeling the molecular stability of these sodium channel variants partially explains the electrophysiological and muscle impairments. Ultimately, adaptive genetic changes favoring toxin resistance appear to negatively impact sodium channel function, skeletal muscle strength, and organismal performance. These functional trade-offs at the cellular and organ levels appear to underpin locomotor deficits observed in resistant snakes and may explain variation in the population-level success of toxin-resistant alleles across the landscape, ultimately shaping the trajectory of snake-newt coevolution.
Introduction
Organismal performance is the integration of many intervening phenotypic steps that translate genetic variation into fitness1. The genetic changes required to adapt to novel selective pressures must therefore integrate with preexisting molecular mechanisms, cellular processes, and physiological functions that cooperate to maintain organismal performance1. Consequently, when selection acts on multiple phenotypes controlled by a single gene, mutations capable of producing phenotypes that overcome one selective demand must also largely preserve the encoded protein’s structure so as to maintain all other molecular functions. When very few, if any, substitutions can accommodate this constraint, mutating a protein that underlies multiple key traits can produce phenotypic trade-offs, whereby adaptive change in one trait jeopardizes another (antagonistic pleiotropy) 2–5. Biophysical constraints on protein structure and function can therefore impose limits on adaptive evolution 2,5–8. As such, decomposing organismal fitness into its constituent molecular phenomena is necessary to understand how phenotypic trade-offs might constrain adaptive evolution, ultimately determining the degree to which evolutionary trajectories are predictable across lineages9.
Here, we demonstrate trade-offs at hierarchical levels of biological organization in a classic coevolutionary system involving toxic newt prey (Taricha spp.) and predatory garter snakes (Thamnophis spp.). Newts employ tetrodotoxin (TTX) as an anti-predator defense10–12, yet some populations of garter snakes have evolved resistance to this potent neurotoxin13–16. Resistance is conferred by mutations directly to the TTX-binding site on voltage-gated sodium channels (proteins: NaV; genes: SCNnA)13,17–21. While these structural changes produce effective resistance to TTX, they occur in highly conserved regions of the protein that are essential to channel function22. Therefore, we hypothesized that the highly conserved structure-function relationship of this critical protein constrains adaptive evolution in this system, whereby mutations that increase TTX resistance are deleterious to channel function and, by extension, physiological and whole-animal performance23–25.
Sodium channels are transmembrane proteins expressed in nerves and muscles that carry vital sensory and motor impulses26–28. The 24 transmembrane passes of these proteins form four voltage-sensing domains (see Figure 2a) that dilate a central, sodium-selective pore through which sodium ions enter the cell, carrying a positive charge that depolarizes the membrane29. TTX occludes the NaV outer pore, thereby arresting nervous and mechanical activity, leading to paralysis and death30,31. This poison’s chemical structure is so compatible with the topology of the nearly universally conserved NaV pore that just nanomolar concentrations of TTX are lethal to nearly all animals30,32. Consequently, newts exploit TTX as a broad-spectrum defense to evade predation, secreting the poison from their skin when attacked33,34, and paralyzing or killing virtually all potential predators11.
Despite this, TTX resistance has evolved repeatedly in nature, generally through target-site insensitivity conferred by mutations directly to the outer pore of the sodium channel19,23,35–38. However, mutations to NaV can have dire phenotypic consequences because many NaV functions are tied directly to protein structure, including selectivity to sodium, the rate of sodium influx by electrodiffusion, and sensitivity to changes in membrane potential 29,39–41. As such, sodium channel mutations are particularly hazardous, often resulting in serious congenital disease or stillborn progeny42. Dietary TTX has nevertheless driven the evolution of snakes carrying very effective NaV mutations, though functional variation is limited to just a few pore sites (Fig. 1a,d,g) 13,23,32. Thus, the tight relationship between molecular structure and sodium channel function is thought to constrain NaV evolution32, and underlie the remarkable convergent molecular evolution seen across TTX-resistant taxa23,25.
The requirement that voltage-gated sodium channels must maintain essential electrophysiological performance while simultaneously preventing TTX blockade therefore imposes competing selective demands on a single locus. Thus, resistant sodium channels are expected to suffer from antagonistic pleiotropy23,25,43. Indeed, previous correlational work15,24,44 suggests that some populations of TTX-resistant snakes suffer locomotor deficits. However, the etiology of whole-animal performance costs remains obscure and has only been considered for a single group of TTX-resistant snakes24,43. Work is needed to directly test the impact of naturally occurring NaV mutations, from the biophysical level to higher order physiology. To determine whether antagonistic pleiotropy creates functional trade-offs that may influence evolutionary outcomes, we need to understand how resistance-conferring mutations disrupt the sodium channel structure-function relationship, resulting in a cascade of performance deficits from molecules to phenotypic trade-offs at the organismal level24,44.
To identify the biophysical mechanism of higher order trade-offs, the molecular effects of TTX-resistant channels must be quantified and then connected to relevant physiological systems (i.e. skeletal muscle) to develop a framework that explains trade-offs at the whole-animal level1,24,43,45. Here, we assessed the biophysical costs of resistance-conferring mutations to the skeletal muscle voltage-gated sodium channel (NaV1.4), primarily finding mutant channels produce less bioelectricity (reduced unitary conductance) and are therefore less able to excite skeletal muscle compared to their TTX-sensitive counterpart. We then investigated mechanical deficits that scale up from effects seen at the molecular level, revealing reduced force output and slowed mechanical activity of skeletal muscle in toxin-resistant snakes. Lastly, we examine how these structural mutations alter the interaction energy dynamics of functional residues in the sodium channel, offering insight into the molecular mechanism of this phenotypic trade-off. The evidence of this trade-off is reproducible across two independently evolved resistant snake species with unique resistance alleles and may underlie the heterogeneous landscape-level distribution of toxin resistance16,20,46,47. We conclude that performance trade-offs at the molecular level are the proximate source of higher-order limitations that give rise to constraints on natural selection, and in so doing, comprise a straightforward map from gene sequence to evolutionary consequence.
Results
Resistant channels produce less current than their toxin-sensitive counterpart
We tested the hypothesis that mutations in the highly conserved sequence of the sodium channel pore (Fig. 1a, d, g; Fig. 2a) underlie organ- and higher-level performance deficits. We characterized two suites of mutations in skeletal muscle voltage-gated sodium channels (NaV1.4) associated with the greatest levels of organismal TTX resistance seen in two species of Thamnophis23. We assessed the electrophysiology of the conserved, TTX-sensitive channel (NaV1.4+, Rattus norvegicus) and two TTX-resistant channel mutants: NaV1.4EPN (D1241E, A1245P, and D1532N) and NaV1.4LVNV (I1520L, I1525V, D1532N, and G1533V) (Fig. 2b). Mutations at these positions on the Rattus channel backbone recreate pore loop (p-loop) substitutions found at analogous sites in wild, TTX-resistant snake NaV1.4 Domains III and IV of Th. atratus (D1277E, A1281P, and D1568N) and Domain IV of Th. sirtalis (I1556L, I1561V, D1568N, and G1569V) (Fig. 1a). As expected, we found that these snake mutations confer strong TTX resistance to the rat channel relative to the conserved protein (Fig. 2c). In search of a possible source of trade-offs, we began by assessing the activation and inactivation parameters of sensitive and resistant channels.
Sodium channels activate in response to depolarizing membranes. Mutations can change the flexibility of the protein, thereby potentially dampening or facilitating the channel’s conformational changes in response to depolarization. Muscle weakness could arise from mutations that cause the sodium channel to be less responsive to depolarization (“hypoexcitability”). The resting membrane potential of skeletal muscle is very negative, varying across phyla from -70mV to nearly as low as -90mV48–50. At such negative resting potentials, healthy sodium channels should remain closed, but ready to open51. When an organism intends to contract a voluntary skeletal muscle, peripheral nerves release acetylcholine onto the neuromuscular junction, depolarizing the postsynaptic membrane through the action of ionotropic Nicotinic Acetylcholine Receptors (nAChR)52. The brief currents raise the local transmembrane voltage. For any sodium channels near enough to be under the electrical influence of this depolarization, the positively-charged residues in the voltage-sensing domains of NaV are repelled extracellularly by this more positive (less negative) internal potential. The conformational change in the voltage-sensing domains applies torsion to the pore-forming domains, dragging them radially outward, thus dilating the channel enough to conduct hydrated sodium ions, producing a current that propagates the depolarizing wavefront53–55. The voltage at which a population of sodium channels opens in near unison defines the “threshold” of an action potential. Channels exhibiting a higher threshold (lower sensitivity to depolarization) will result in muscles that require a greater stimulus to produce a contraction. Thus, the voltage sensitivity of sodium channel activation is a reasonable target for the investigation of trade-offs in toxin resistant snakes.
In whole-cell patch clamp electrophysiology experiments where transmembrane voltage is controlled by an external amplifier, the population of channels in the membrane is held at a constant voltage, conditioning the proteins into nearly uniform molecular configurations. At potentials near or more negative than rest, the voltage-sensing domains are pulled intracellularly, cinching the pore-forming domains, preventing channel patency56–58. As test potentials are increasingly depolarized (i.e. sequentially more positive), greater and greater proportions of the sodium channel population will respond to the activating stimulus, resulting in macroscopic sodium currents that are larger in magnitude. The development of sodium currents from a population of channels conditioned into an equilibrium of nearly uniform conformational states is referred to as “steady-state activation.”
Any mutation near the pore or voltage-sensing domains could feasibly alter the rigidity or flexibility of this orchestrated molecular reconfiguration, thereby requiring greater or lesser depolarizing potentials to achieve the same magnitude of sodium current. However, we found that the voltage dependence of steady-state activation is identical across all sodium channel variants observed in this study (Fig. 2e). That is to say, across a range of membrane potentials from -80mV to +80mV, equivalent proportions of both sensitive and resistant sodium channel populations respond to the same depolarizing stimulus. Therefore, voltage sensitivity of activation does not explain higher order locomotor deficits.
In addition to activating in response to depolarization, healthy sodium channels inactivate in time- and voltage-dependent manners59. Under normal circumstances, this is achieved by a hydrophobic intracellular loop that becomes mobile once the channel is activated, contorting to occlude the channel pore from the intracellular face, preventing further passage of sodium despite voltage-sensing domains remaining in the activated conformation60. The total macroscopic current decays within a few milliseconds as inactivation stochastically sets in for all sodium channels throughout the membrane, known as “fast inactivation”61,62. Slow inactivation takes place over much longer time scales and provides a less convincing etiology for trade-offs in toxin-resistant snakes61. Once inactivated, sodium channels are very unlikely to return to a conducting state until the membrane potential returns to near or below rest, resetting the voltage-sensing domains and removing the hydrophobic inactivation gate. If a mutation were to cause premature inactivation, perhaps by increasing the mobility of the inactivation gate, then less overall current would pass through the membrane, thereby limiting the depolarizing power of the sodium channel and hampering the excitability and contractility of skeletal muscle.
TTX resistance mutations in snakes are found on the extracellular face of the channel pore, separated from the inactivation gate on the intracellular face63. Therefore, pore-domain mutations are not known to affect the voltage-dependence of fast inactivation, which is consistent with our findings here (Fig. 2e)43,64,65. However, holding voltage constant, we do see that one resistant variant inactivates more quickly than others (Supp. Fig. 1g,h). While the triple-point mutant (R.n.NaV1.4EPN) exactly follows the inactivation rate of the ancestral TTX-sensitive channel (R.n.NaV1.4+), the quadruple-point mutant (R.n.NaV1.4LVNV) enters steady-state inactivation significantly faster than the TTX-sensitive channel. This premature onset of fast inactivation in one toxin-resistant channel may underlie higher-order physiological deficits unique to Thamnophis sirtalis (described later). Assessing how quickly inactivation can be removed in these mutant channels, we find that neither the rate (Supp. Fig. 2) nor the voltage dependence of recovery from inactivation (Supp. Fig. 2) are affected by these pore-domain mutations. Thus, these distinct resistance-conferring p-loop mutations appear to exert unique effects on some channel functions while entirely sparing others. However, with identical steady-state activation and faulty inactivation in only one of the two mutants, these findings alone are insufficient to explain all locomotor deficits in both Thamnophis atratus and Thamnophis sirtalis.
Without evidence supporting a kinetic explanation for locomotor deficits (Supp. Fig. 1a-d), then the magnitude of sodium currents, and thus the excitability and contractility of skeletal muscle, must be determined by other fundamental properties of NaV. Across the current-voltage relationship, TTX-resistant channels produce far less macroscopic current than their toxin-sensitive, conserved counterpart (Fig. 2d, 30% less in R.n.NaV1.4EPN and 50% less in R.n.NaV1.4LVNV). However, this global, macroscopic current can be decomposed into three variables: how many channels are active in the membrane (N), the maximum statistical likelihood of finding any of the N channels open (Po,max), and how much sodium current each individual channel conducts when it is actually open (unitary conductance, or gNa)66. The most expedient electrophysiological technique to measure these variables is Non-Stationary Noise Analysis (NSNA, Fig. 2f)66. Assuming all sodium channels encoded by the same allele perform roughly identically, then any variation in macroscopic current not otherwise attributable to thermal, electromagnetic, or vibrational noise, should be quantal: one open sodium channel always produces the same amount of current. Therefore, by taking an average of the inactivating phase of an ensemble of 1000 sodium currents evoked under identical conditions, the statistical variance of the inherent noise around the mean is proportional to the unitary current of the sodium channel. If activation and inactivation determine at what voltages sodium channels open and for how long they remain open, then the unitary current determines the amount of sodium that can flow per unit time and likely underlies the differences in macroscopic currents observed here between ancestral and mutant channels (Fig. 2d).
Through NSNA, we found median single channel conductance, gNa, (Fig. 2g) to be significantly reduced in toxin-resistant mutants (29.3% in R.n.NaV1.4EPN, p = 0.0318; 53.6% in R.n.NaV1.4LVNV, p = 8.84×10-5) without significant changes in channel number, N, (Fig. 2h, Supp. Fig. 1f) or peak open probability, Po,max (Fig. 2i, Supp. Fig. 1.e). Consequently, these mutations appear to reduce the total ionic current by reducing unitary conductance, gNa, as opposed to reduced channel expression or activity. Furthermore, the degree of biophysical dysfunction reported here is very likely to have consequences at higher levels of organization capable of producing the muscular and locomotor deficits observed in previous work24. The electrophysiological data presented here demonstrate that TTX resistance achieved by target-site insensitivity in NaV1.4 is inherently costly by way of reduced sodium channel conductance.
Resistant muscles develop less force than their toxin-sensitive counterparts
Sodium channels propagate the action potential wavefront across the sarcolemma and into the T-tubules, activating voltage-sensing proteins (Dihydropyridine Receptors, DHPR) that physically couple to and activate Ryanodine Receptors (RyR1) in the sarcoplasmic reticulum52,67,68. The opening of RyR1 releases sarcoplasmic calcium stores into the cytosol where Ca2+ can interact with troponin, thereby engaging the contractile machinery69. The greatly reduced conductance of TTX-resistant channels likely limits their ability to depolarize the sarcolemma. We therefore hypothesized that sodium channels generating less conductance will fail to maximally engage this cascade of contractile mechanisms, resulting in weaker and slower muscles.
We contrasted the physiological performance of skeletal muscle from individual Th. atratus and Th. sirtalis carrying ancestral TTX-sensitive pore sequences (Th.spp.NaV1.4+) to their three- and four-point mutant counterparts (TTX-resistant Th.a.NaV1.4EPN and Th.s.NaV1.4LVNV). We first established the concentration-dependent effect of TTX on isolated skeletal muscle twitch force (Fig 3a,e, Supp. Fig. 4a,f). Muscles from both Th. atratus and Th. sirtalis carrying the ancestral pore-domain sequences were paralyzed by low concentrations of TTX (IC50 ∼ 90-150 nM), while those carrying p-loop mutations withstood orders of magnitude greater TTX concentrations (as much as 100-fold) with many individual muscles continuing to contract at concentrations exceeding 70μM. We observed a strong correlation between TTX resistance in skeletal muscle tissue and resistance in the live animal (adjusted R2 = 0.5427, p = 1.4×10-9, df = 46).
To test for a trade-off with increasing TTX resistance, we evaluated the force and speed of toxin-sensitive and toxin-resistant muscles grouped by species and genotype. To avoid potential genetic background-dependent effects, we only made relevant comparisons within species. The mean twitch force produced by skeletal muscles carrying TTX-resistant NaV1.4 was three-fold (Th. atratus) and five-fold (Th. sirtalis) less than TTX-sensitive skeletal muscles (Fig. 3b,f, Supp. Fig. 4c,h). The peak rate of force development (Ḟmax), which is known to be closely related to action potential generation and therefore sodium channel performance67,69, is reduced accordingly (Supp. Fig. 4b,g). Close attention to contraction timing generally shows that deficits in the speed of force development compound throughout the progression of any single twitch (Fig. 3i,j) with variable effects on overall contraction duration (Supp. Fig. 4d,i).
Differential effects of toxin resistance mutations on the progression of force development and contraction timing were observed in each snake species. For toxin-resistant Th. atratus (Fig. 3i) carrying NaV1.4EPN, muscles reached a premature maximum force that quickly resolved back to baseline levels (Fig. 3i, Supp. Fig. 4d,i). However, force development was definitively delayed in Th. sirtalis carrying NaV1.4LVNV (Fig. 3j). These unique p-loop mutations exert distinct effects on contraction timing that scale up to more complex activation as well, such as in tetanus. While peak tetanic force was reduced in both resistant populations (Fig. 3c,g), only resistant Th. sirtalis (NaV1.4LVNV) failed to develop force as quickly as their TTX-sensitive counterparts (Supp. Fig. 4e,j). This slower rate of force generation is consistent with the electrophysiological finding that only NaV1.4LVNV entered steady-state inactivation more quickly than ancestral channels (Supp. Fig. 1g,h). Presumably, as a tetanic contraction progresses, a larger proportion of the channel population inactivates, preventing sustained rapid force development. This likely gave rise to the slower morphology of the mean tetanic contraction in contrast to any of the other species-genotype comparisons.
Comparing peak force (Fmax) to TTX IC50 revealed a clear trade-off between resistance-conferring and ancestral genotypes at the muscular level (Fig. 3d,h). Greater toxin resistance in muscles is accompanied by reduced muscular force in both species, consistent with biophysically compromised channels. Therefore, the pore domain mutations in Th. atratus and Th. sirtalis possess both gain-of-function and loss-of-function features, although it remains unclear which individual substitutions contribute gains in toxin resistance or losses in excitability. As such, we conclude that mutations altering the structure of the sodium channel pore coincide with reduced channel function and concomitant tetrodotoxin resistance. Skeletal muscles bearing channels with these altered structures also show reduced performance despite extraordinary tetrodotoxin resistance. These results sufficiently explain previous findings that some toxin-resistant snakes are slower than their toxin-sensitive counterparts 24.
Resistance mutations perturb the stability of the selectivity filter and S5/S6 helices
Molecular modeling of the sodium channel provides insights into mechanisms that likely contribute to the reduced unitary conductance of toxin-resistant channels. These biophysical findings may also explain why the quadruple-point mutant (Th.s.NaV1.4LVNV) enters into steady-state inactivation faster than either the toxin-sensitive (Th.spp.NaV1.4+) or triple-mutant toxin-resistant (Th.a.NaV1.4EPN) proteins. To understand the molecular stability of these three channel variants, the crystal structure of the human sodium channel was modified to create a homology model of the ancestral snake protein17,70. This ancestral version was then further modified to create homology models of both NaV1.4EPN and NaV1.4LVNV 13. These three models were then run through an optimization routine. We compared the stability (ΔG) of key residues between these three protein models, including the resistance-conferring sites, the selectivity filter, the activation and inactivation gates, and hydrophobic residues at the interface of the pore-forming loops and the S5/S6 transmembrane passes, that cumulatively form the ion-conducting pore71–73.
Calculating the total interaction energy of residues in the ancestral channel and contrasting to mutant toxin-resistant channels demonstrates that the toxin-sensitive channel is the most stable protein, exhibiting a cumulative interaction energy of -393.85 kj/mol across all six residues (sites 1248, 1252, 1527, 1532, 1539, and 1540, using the human numbering system). While the triple-point mutant, NaV1.4EPN, is comparably stable at these sites (-366.84 kj/mol), the quadruple-point mutant, NaV1.4LVNV, is substantially destabilized (-42.43 kj/mol). In fact, NaV1.4+ and NaV1.4EPN demonstrate residue pairings whose combined interaction energies are therefore 6.43x and 6.30x more stable, respectively, than the sum of the interaction energies between the same six sites in NaV1.4LVNV (Supp. Table 12).
When decomposing these total interaction energies into their component interactions, it becomes clear that N1539 strongly destabilizes both NaV1.4EPN and NaV1.4LVNV, accounting for the loss of 122.91 kj/mol stability in NaV1.4EPN and 148.54 kj/mol in NaV1.4LVNV. However, in NaV1.4EPN, E1248 counteracts the instability conferred by N1539. While the conserved residues at positions 1527, 1532 and 1540 each have small, destabilizing interactions with N1539, E1248 interacts favorably with all four residues (I1527, I1532, N1539, and G1540). In fact, the most destabilizing interaction for E1248 is with its nearby neighbor, P1252, also found in resistant Th. atratus channels. So, while D1539N reduces the total favorable interaction energy at site 1539 by 86.4% in NaV1.4EPN, D1248E increases the favorable interaction energy by 353.6%. Therefore, E1248 appears to act as a compensating force against the destabilizing action of N1539.
In NaV1.4LVNV, however, no such countervailing force exists. While D1248 and N1539 appear to have a small, mutually stabilizing interaction (-8.58 kj/mol), the sum of interactions at sites 1527, 1532, 1539 and 1540 produce a net destabilizing effect. The majority of lost stability in NaV1.4LVNV (ΔΔG=+351.41 kj/mol) can be explained by destabilizing interactions between three hydrophobic pairings (L1527-F412, V1532-Y1581, and V1540-G1512) and a hydrophilic interaction (N1539-K1244) (Figure 4). L1527 most likely interferes with the Domain I pore loop residue, F412, due to crowding by the bulky dimethyl group of leucine compared to the single methyl terminus of isoleucine. However, changes in residue volume may not completely explain why V1532 interacts negatively with the aromatic Y1581 of Domain IV S6. Rather, V1532 possesses a stronger attraction to a nearby L1541 in the ascending pore loop than isoleucine does in the toxin-sensitive channel. V1532 is consequently contorted by its interaction with L1541 and interferes with the hydrophilic hydroxyl group on Y1581, thereby yielding an unfavorable interaction. Further, V1540 negatively interacts with G1512 found on an extracellular loop, mostly likely due to steric hindrance compared to the ancestral glycine. Lastly, N1539 demonstrates a deleterious interaction with the positively charged residue, K1244, a key residue that constitutes the Domain III selectivity filter (the “K” of DEKA). This unfavorable interaction between N1539 and the Domain III selectivity filter is present and similar in magnitude in both NaV1.4EPN and NaV1.4LVNV.
In addition to directly altering channel pore topology, mutations in TTX-resistant sodium channels appear to have indirect effects on channel function by destabilizing conserved residues nearby. Specifically, resistance mutations reduce the stability of the hyper-conserved DEKA residues that form the selectivity filter (D406, E761, K1244, and A1536)39. Though not directly mutated, the selectivity filters of the toxin-resistant channels are indeed less stable (-7.16% in NaV1.4EPN and -10.7% in NaV1.4LVNV)(Supp. Table 13). The reduced stability of the ion-conducting pathway may contribute to the reduced conductance of the mutants without negatively affecting sodium selectivity 43,64. We found the toxin-sensitive channel’s median unitary conductance to be 15 pS, which corresponds to conducting around 9400 sodium ions per millisecond (at +10 mV). Presumably, this conductance is facilitated by a stable and rigid pore structure. Future molecular dynamics work will be needed to calculate the flow rate through a less stable pore. Such modeling may provide a full explanation of the reduced unitary conductance in toxin-resistant channels74.
Examining how the network of residue-residue interactions shifts from variant to variant may also help explain, at least from a structural perspective, why steady-state activation is unaffected by pore-domain mutations, and also why NaV1.4LVNV does enter into steady-state inactivation faster than either NaV1.4+ or NaV1.4EPN. Pore-domain mutations on the extracellular face would not be expected to exert direct effects on the activation and inactivation gates near the intracellular face. Indirectly, however, the interior aspect of the pore loops interface with the S5 and S6 helices through hydrophobic interactions, which in turn, directly interact with the activation and inactivation gates (Figure 5). Thus, some aspects of steady-state kinetics can be mechanistically explained by interaction energies between the p-loops and the S5/S6 helices.
The activation gate, a group of hydrophobic amino acids near the intracellular face of the pore, is dilated following the molecular reorganization of the voltage-sensing domains during depolarization. The activation gate in our modified NaV1.4 model is composed of Y450 (Domain I S6), F805 (Domain II S6), F1298 (Domain III S6), and F1601 (Domain IV S6)75,76. Our electrophysiological data show that activation is identical across all three channels. Likewise, the total interaction energy of the activation gate residues is identical in all three molecular models. Steady-state activation is therefore not negatively affected by pore-domain mutations.
Inactivation, however, appears to differ across these channels. Our electrophysiological and myographical data suggest that NaV1.4LVNV inactivates prematurely. At rest, the inactivation gate (I1310, F1311, M1312) is locked in place by interactions with the S5 and S6 helices of Domains III and IV (Figure 5). From the intracellular face, the IFM particle interacts through hydrophobic forces with A1152 (DIII S5) and F1305 (DIII S6). From above, closer to the ion-conducting pathway, the IFM particle interacts, again through hydrophobic forces, with A1481 (DIV S5) and I1594 (DIV S6). As expected, with the pore mutations far away from the inactivation gate, the interaction energies between the IFM particle and the DIII and DIV S5/S6 residues are virtually identical across all channel variants studied here. Thus, differences in tendencies toward inactivation must arise from indirect effects of the pore-domain mutations, likely due to differences in the mobility of S5/S6 domains during the molecular reorganization that initiates inactivation.
The intracellular ends of the DIII and DIV S5 and S6 helices cradle the IFM particle while the extracellular ends of those same transmembrane α-helices form direct hydrophobic interactions with the interior aspect of the reentrant p-loops (Figure 6). For inactivation to initiate, the S5 and S6 helices must slide against one another to liberate the IFM particle, permitting it to contact and occlude the intracellular face of the pore59. Over both Domains III and IV, there are 44 key residues linking the p-loops to the S5/S6 helices. To mobilize the inactivation gate, NaV1.4+ channels must overcome the total interaction energy of those 44 residues (-4044.23 kj/mol) (Supp. Table 14). These residues in the triple-point mutant, NaV1.4EPN, are 76.92 kj/mol less stable than NaV1.4+, spread over both Domains III and IV. However, the p-loop-S5/S6 interactions are least stable in NaV1.4LVNV (+109.06 kj/mol), whose mutations are all concentrated in Domain IV. This decreased interaction energy between the NaV1.4LVNV pore-loops and the S5/S6 domains likely results in hypermobility of the DIV S5/S6 helices, thereby prematurely releasing the IFM particle to engage inactivation.
Discussion
We demonstrate structural design limitations that constrain sodium channels from maximizing both biophysical performance and resistance to tetrodotoxin (TTX). The sites at play in adaptation to TTX serve a structural role in the sodium channel77, such that changes to these residues impose functional costs on the channel’s critical role in excitability. Accordingly, resistant Thamnophis skeletal muscles appear to bear the effects of poorly conducting channels in the form of reduced force and variable speed. Snakes rely on these sodium channels to maximize muscular and organismal performance, essential for foraging, wrestling with prey, escaping predators, hunting, and simply breathing. As such, key physiological and organismal tasks will be compromised by the reduced ionic conductance and myasthenic muscles observed in resistant snakes.
The reduced sodium channel unitary conductance reported here is novel, providing direct evidence of dramatic, traceable consequences to higher orders of organization. Though slow inactivation has been proposed as a biophysical penalty of toxin resistance mutations43, this functional cost is unlikely to underlie the scope and magnitude of skeletal muscle deficits we describe here. Furthermore, our physiological assessment of skeletal muscles in the snake-newt system provides a key link between molecular and organismal variation. Together with insights into the stability of key residues across the sodium channel variants studied here, our work ultimately provides a firm etiology of functional trade-offs in TTX resistant snakes that may ultimately limit phenotypic evolution.
The voltage-gated sodium channels of garter snakes exhibit similar but distinct evolutionary solutions and drawbacks to the steep physiological challenge of dietary TTX. Each solution is forced to preserve the critical function of the sodium channels on which organismal performance and survival depends. Both TTX resistance and muscular performance directly rely on channel pore structure, such that minute topological changes underlie broad phenotypic variation in an apparent zero-sum fashion. The biophysical insights from modeling the interaction energy of three channel variants demonstrates important constraints on protein structure and function, highlighting that even seemingly benign, conservative mutations can have cascading consequences throughout the protein.
The function of the sodium channel depends chiefly on the movement and relative spatial configuration of hydrophobic residues. The mutations seen in these snakes (NaV1.4EPN and NaV1.4LVNV) appear to negatively impact function by interfering with the structural stability of the protein, resulting in hypermobility of these important hydrophobic residues in the pore and, in one case, a premature molecular reorganization that stops current flow. These molecular limitations on the types of mutations that are tolerable to protein function scale up to produce cellular and organ-level phenotypes that are likely suboptimal in many natural landscapes. Though other sodium channel mutations might offer perhaps even greater toxin resistance than NaV1.4EPN or NaV1.4LVNV, such theoretical mutations likely breach the threshold of tolerable trade-offs to persist in nature. Resistance-conferring mutations in Thamnophis persist in just a narrow region of the sodium channel23 compared to the distribution of mutations in TTX-bearing species that must also be resistant to their own defenses35–38. In tetrodotoxic taxa, the persistent exposure to their own poison may require the incorporation of even more costly mutations. Garter snakes, on the other hand, may prey on newts irregularly 78, and as active, widely foraging, semi-aquatic predators79, may demand greater muscular performance, which likely acts as a brake on selection towards extreme TTX resistance.
Regardless, for these extremely resistant snakes, possessing such compromised proteins and impaired physiology may pose dramatic ecological and evolutionary consequences. Resistant snake performance may not meet the demands imposed by predators and prey in some ecological communities, or of harsh abiotic conditions, thereby limiting the geographic scope of resistance-conferring alleles and success of TTX resistance phenotypes. Snakes carrying Th.a.NaV1.4EPN and Th.s.NaV1.4LVNV may thus be confined to only those ecological communities or abiotic environments that can support the costly phenotypes those genotypes produce. Indeed, these TTX-resistant sodium channel genotypes are highly localized and only appear in certain “hotspots” rather than broadly distributed across the landscape as are toxin-sensitive animals18,20,21,46,80.
Antagonistic pleiotropy is thought to drive the convergent molecular evolution of TTX-resistant sodium channels23,25. Here, both snake species (Th. atratus and Th. sirtalis) have evolved the identical replacement D1568N 13, despite being distantly related81. Given that all proteins are universally subject to the chemical and physical forces that govern folding and function5, such biophysical constraints may ultimately underlie the degree to which evolutionary trajectories are predictable across lineages. Consequently, we conclude that the two suites of mutations (NaV1.4EPN and NaV1.4LVNV) that have independently arisen in Th. atratus and Th. sirtalis, respectively, are the outcome of selection constrained by antagonistic pleiotropy to balance two opposing selective demands.
The trade-off we describe suggests that the structure-function relationship of the sodium channel is such that many possible mutations can reduce TTX blockade, but relatively few can do so while simultaneously maintaining sufficient levels of sodium channel function. This in turn must limit the number of mutations that satisfy two opposing selective demands imposed on garter snakes engaged in coevolutionary arms-races with newts. Such antagonistic pleiotropy therefore constitutes a limitation on adaptive evolution that, in the face of natural selection, promotes novel genetic variation in loci otherwise highly constrained by physical forces acting on the encoded proteins. Ultimately, simple biophysical limitations can manifest across levels of organization to shape molecular and phenotypic evolution in natural systems, potentially impacting the fates of populations and even species interactions.
Materials and Methods
Animal sampling and captive care
To measure whole-animal TTX resistance, as well as the properties of skeletal muscles in snakes with different sodium channel genotypes, we field collected individual Thamnophis atratus carrying NaV1.4+ (Th.a.NaV1.4+) and Th. sirtalis carrying NaV1.4+ (Th.s.NaV1.4+) for the ancestral state and their three- and four-point mutant counterparts (Th.a.NaV1.4EPN and Th.s.NaV1.4LVNV). Our final sample included 62 animals from natural populations in California and Oregon, as well as a few lab-born individuals: 19 TTX-sensitive Th. atratus from three sites; 21 TTX-sensitive Th. sirtalis from six sites; 10 TTX-resistant Th. atratus from two sites; and 12 TTX-resistant Th. sirtalis were caught from four sites (Supp. Table. 1).
Field-caught Thamnophis spp. were kept in captivity at either Utah State University (USU) or University of Nevada, Reno (UNR) prior to any work. Snakes were provided water ad libitum and fed either previously frozen fish (tilapia or trout) or mice once or twice weekly with no detectable difference due to diet (data not shown). The snakes were kept at 26 ± 1 °C and on a 12L:12D photoperiod with hide boxes provided. All snakes were housed individually in mouse-style containers or aquaria commensurate with animal size. All cages had one end situated above heat tape to achieve a thermal gradient. All animal care procedures followed USU and UNR IACUC protocols (USU IACUC #1008 and UNR IACUC #00687).
Organismal resistance assay
Animals were encouraged to sprint along a linear racetrack lined with artificial turf. Infrared sensors (every 20 to 50 cm) recorded snake progression along the 4 m track14,82,83. The two shortest time intervals between sensors were averaged and used to calculate the animal’s maximum speed. Snakes that refused to sprint down the track were omitted. All trials were performed at 26 ± 1 °C. Only two individuals (E.D.B., Jr. and C.R.F.) conducted sprint trials to reduce interrater variation.
This assay was repeated following intraperitoneal injection with TTX (MilliporeSigma, Burlington, MA, USA; Tocris, Bristol, UK; Alomone Labs, Jerusalem, Israel; dissolved in citrate buffer pH 6.8; injection volume ≤0.5 mL). All doses were calculated and delivered in mouse units. A mouse unit is the amount of TTX (2.857×10-7 g) required to kill a 20-g mouse in 10 minutes84. Scaling to account for total body mass of the snake to be injected, this mass-adjusted mouse unit (MAMU) contextualizes the snakes’ toxin resistance relative to TTX-sensitive mammalian tissue. Doses are calculated by the equation Dose (MAMU) = (Dose g TTX / m (g)) • (20 g / 2.857×10-7 g TTX), where m is the mass of the individual snake injected. A distribution period of 30 minutes was allowed before racing the individual. Control injections of physiological saline have no effect on snake performance (Brodie and Brodie 1990), and TTX resistance is not affected by repeated short- or long-term exposure83,85. No snake received more than five (5) injections, and all doses were separated by 48 hours to ensure the elimination of residual TTX and to reduce the effect of fatigue83.
Animals were not given standard doses of TTX but rather dosing was informed by the population of origin. We injected TTX starting at 1 MAMU for snakes from known low-resistance populations and 10 MAMU from known high-resistance populations. The dose was then serially increased (3, 5, 10 and sometimes 15 or 25 MAMUs for low-resistance; 25, 50, 100 or 250 MAMUs for known high-resistance). For individuals with exceedingly great resistance, the dose-response protocol was truncated to 250 MAMU to avoid potentially irreversible damage to the animal and preserve costly toxin reserves. The dose of TTX needed to slow an individual to 50% of its baseline sprint speed is reported as the 50% MAMU (analogous to a 50% inhibitory concentration, IC50).
Here, 50% MAMU was calculated by fitting the dose-responses to a sigmoidal curve. A custom R script was designed to handle this fitting and used the equation Speed (m/s) = A2 + (A1-A2)/(1+e((log(Dose (MAMU))-log(x0))/dx)), where A1 and A2 are the upper- and lower-bounds, x0 is the center of the curve (50% MAMU) and dx is the rate of decline with increasing dose. The dose for baseline speed (0 MAMU TTX) was revised to 0.1 MAMU to abide by the limitations of the natural logarithm. A minimum of five (5) dose responses is needed to fulfill this model. Where conditions precluded the administration of more than three (3) doses of TTX, the sigmoidal fit was forced by including sufficiently large doses where speed is assumed to converge to zero m/s.
Genotyping
Regions of SCN4A corresponding to the third and fourth pore domains of NaV1.4 (DIII-PD & DIV-PD) were sequenced for all snakes. Genomic DNA was extracted from tail snips, skeletal muscle, or liver tissue by DNeasy Blood & Tissue Kits (Qiagen Inc., Germantown, MD, USA). Fragments of genomic DNA corresponding to DIII-PD & DIV-PD of SCN4A were amplified by PCR using Thamnophis-specific primers (Supp. Table 2)13,20. Fragments were submitted to the Nevada Genomics Center (UNR) or the DNA Analysis Facility (Yale University) for bidirectional Sanger sequencing using amplification primers, and sequences were resolved on an Applied Biosystems Prism 3730 DNA Analyzer or 3730xl DNA Genetic Analyzer. Nucleotide sequences were aligned and translated to amino acid sequences using a variety of software, including Geneious versions 4.8.5 and 6 (Biomatters Ltd., Auckland, New Zealand; http://www.geneious.com) or Sequencher 4.2 (Gene Codes Corp., Ann Arbor, MI, USA), ClustalW 1.83 (Clustal W and Clustal X version 2.0) and MacClade 4.0886–88. We deposited all sequences in GenBank (accession numbers: Exon 22 ON609732-ON609763; Exon 26 ON609764-ON609824).
Animals from both Th. atratus and Th. sirtalis were scored as NaV1.4+ (TTX-sensitive) for translated SCN4A sequences reading 1276-MDIMYA-1281 in DIII and 1556-ICLFEITTSAGWDG-1569 in DIV. Thamnophis atratus were scored NaV1.4EPN (TTX-resistant) when exhibiting D1277E and A1281P in DIII and D1568N in DIV. Thamnophis sirtalis were scored NaV1.4LVNV (TTX-resistant) when exhibiting I1556L, I1561V, D1568N, and G1569V. Snakes confirmed as heterozygous with respect to SCN4A were omitted from further testing.
Dissection
Immediately prior to myography work, animals were humanely euthanized without the use of chemicals that could affect sodium channel performance. Euthanasia was performed independent of circadian rhythm and without attention to feeding schedules or shedding cycles. All animal procedures followed UNR IACUC protocols (IACUC # 00687). The specimen was submerged at all times in fresh Krebs buffer and supplied with 95%O2/5%CO2. A segment of 2 to 4 cm of the snake body was transferred to a dish supplied with the same solution and gasses for dissection under an upright stereomicroscope. Cleaved down the ventral midline, the segment was skinned and pinned to the dish by the ribs and intercostal muscles. The ventrolateral muscle group known as the iliocostalis muscle was extracted using curved-blade microdissection scissors. All manipulations of the muscle were made to minimize trauma and any undue tension on the muscle prior to suspension, typically in under three minutes.
Skeletal muscle performance and resistance assay
Myography was performed on a Contraction System Chamber 800A which was operated by Dynamic Muscle Control DMC v4 software (Aurora Scientific Inc., Aurora, Ontario, Canada). Tensile strength was recorded by an Aurora Scientific 300C-LR Dual Mode Lever System. Before suspending the iliocostalis muscle, the force transducer was calibrated using a 40 g load. A 0.15 mm nylon suture was driven through the muscle group and secured by modified Timber Hitch knots. The chamber was prefilled with Krebs buffer, continuously perfused with 95%O2/5%CO2, and maintained at 25 °C by water-jacketed circulation (Lauda Alpha A6, Germany). The buffer was allowed to equilibrate for five minutes prior to submerging the muscle. Electric field stimulation (EFS) platinum electrodes ran parallel along the length of the muscle which was situated equidistant from both electrodes.
A baseline tension of 0.1 g was applied to the suspended muscle. The EFS was optimized by applying a 500 μs pulse of direct current in increasing magnitude from 1 to 300 mA. The EFS stimulus was set to 1.5x the smallest magnitude current that provided the greatest recruitment of the muscle. At this optimum electrical stimulus, the length-tension relationship was optimized by increasing the baseline tension until the force output peaked. A 1 min rest period was allowed between the end of the optimization routine and the beginning of the performance trials. All trials were recorded at 10 kHz.
Transient Muscular Force Assessment
Four consecutive EFS pulses, one every five seconds, were applied to the muscle under the optimized conditions as described above. These traces are the basis of comparison for muscular force generation between genotypes and TTX resistance phenotypes within Thamnophis spp. The resulting transient contractions were analyzed for peak force magnitude (N/g of tissue), contraction duration (from center of the stimulus to 50% relaxation), and latency (time from center of the stimulus to the peak force generated). The first derivative of each trace was analyzed for its peak rates of force development and relaxation (kN/g s-1). Of the train of four pulses, comparisons were only made between contractions of the same order.
Tetanic Muscular Force Assessment
Under the same optimized EFS conditions as the transient force assay, the muscle was stimulated at 1000 Hz for 2 seconds. The resulting tetanic contraction was analyzed for peak force magnitude (N/g of tissue), contraction duration (from stimulus onset to 50% relaxation), and latency (time from stimulus onset to the peak force generated). The first derivative of each trace was analyzed for its peak rates of force development and relaxation (kN/g s-1). The tetanus routine is so stressful that it is terminal. A new muscle was dissected and used following this protocol. Transient force, rheobase, and tetanic force protocols were run within a period of 15 minutes with 1 minute of rest between protocols.
Skeletal muscle TTX resistance assay
The transient muscular force assessment and EFS optimization were repeated with a new, freshly dissected muscle. The optimal transient was stimulated and recorded. Known quantities of TTX dissolved in citrate buffer pH 6.8 were added to a bath of known volume. The doses began as low as 1 nM and increased until the transient force peak magnitude decayed. In many cases, [TTX] was increased until the absolute abolishment of muscular activity was achieved, though this was not possible for every experiment. The volume of TTX in citrate buffer added to the bath was kept below 0.01% of the total bath volume to minimize the effects of diluting the Krebs buffer. However, in cases when exceedingly great [TTX] was needed to inhibit the contraction and stock concentrations of TTX dilutions were too low to keep the volume added below 0.01% of the bath volume, a separate experiment where an identical volume of citrate buffer (without TTX) was added to verify the absence of dilution-related effects.
Dose response curves were generated from the peak transient contraction force magnitude. The responses were fit to a sigmoidal curve of the equation Force (N/g) = A2 + (A1-A2)/(1+e((log([TTX] (nM))- log(x0))/dx)), where A1 and A2 are the upper- and lower-bounds, x0 is the center of the curve (toxin concentration at half-maximal inhibition, IC50) and dx is the rate of decline in force with increasing TTX concentration. Comparisons of IC50 (nM) and dx (N/g nM-1) were restricted to animals within the same species, grouping animals by genotype. To abide by the limitations of the natural logarithm, the dose at baseline transient force magnitude (0 nM TTX) was revised to 0.1 nM. If the contraction recorded at 0 nM TTX was unusable or unavailable, the dose responses were normalized to the maximum force produced throughout the routine (frequently occurring within the first 10 nM TTX). Dose response data could not be fit with fewer than four concentration-force points collected. However, when only three points were available, a fourth point was fabricated where the force would visually converge to zero N/g (i.e., 100 µM, which abolishes contractility in even the most resistant snake muscles).
Blinding and randomization
No animals had confirmed genotypes at the time of the whole-animal performance and resistance assays. For Th. atratus, all animals were processed blindly and genotyped post hoc to prevent rater bias. For Th. sirtalis, genotypes were scored prior to the experiment, and an investigator (JSR) made a randomized and blinded list of snakes to process. This list was handed off to the investigator performing the myography experiments (REdC) unaware of any animal’s genotype and whole-animal phenotype until all animals on the list were processed.
Assessment of NaV1.4 function across genotypes
Mammalian Cell Culture
Adherent Human Embryonic Kidney 293 cells (HEK293, authenticated from vendor, ATCC® CRL1573TM, Manassas, VA, USA) were maintained in Dulbecco’s modified Eagle’s medium (DMEM, ThermoFisher Scientific, Waltham, MA, USA) supplemented with 10% v/v heat-inactivated fetal bovine serum (Atlanta Biologicals, Flowery Branch, GA, USA) and 100 U/mL penicillin/streptomycin (ThermoFisher Scientific, Waltham, MA, USA). The cells tested negative for mycoplasma (validated by ATCC and by the authors using Hoechst 33258 staining). Cells were not assessed for cross-contamination by other mammalian cell lines, as our electrophysiological methods are robust to any consequent differences. Cells were plated on nonpyrogenic culture vessels (Sarstedt, Nümbrecht, Germany) and incubated at 37°C supplied with 13.0% CO2 for pH 7.3 at 100% humidity. The age of the cell populations ranged from approximately 7 to 20 passages. All passages were carried out with TrypLE Express trypsin reagent (ThermoFisher Scientific, Waltham, MA, USA). HEK293 is a commonly misidentified cell line (ICLAC) but was used here only as a representative eukaryotic expression system. Cells were only used for transient lipofection and expression of sodium ion channels and subsequent patch clamp electrophysiology.
Plasmid construction
An oocyte-expression plasmid containing a Rattus norvegicus sodium channel clone (rSCN4A) was kindly provided by Mohammed Chahine (Université Laval, Québec, Canada). A rat clone was deemed acceptable for preliminary research provided that the pore domains are entirely conserved between rats and snakes. The insert was extracted by EcoRI restriction enzyme digestion and inserted into a bicistronic mammalian expression vector, pIRES-hrGFP-2a (Agilent, Santa Clara Co., CA, USA). The new vector was then transformed into 5-alpha competent high-efficiency E. coli (strain C2987I, New England Biolabs, Ipswich, MA, USA) and amplified. The pIRES-hrGFP-2a-rSCN4A plasmid was then isolated by mini-prep (Qiagen Inc., Germantown, MD, USA).
Site-directed mutagenesis
The mutations made to pIRES-hrGFP-2a-rSCN4A are identical in relative position and amino acid substitution to those in Thamnophis SCN4A and in some cases converge with tetrodotoxic Taricha and many tetrodotoxin-defended fishes of the Tetraodontidae. We specifically focused on SCN4A mutations and excluded peripheral nerve-type sodium channels (SCN8A/NaV1.6 and SCN9A/NaV1.7) because all garter snakes carry the same mutations in those channels24,89.
The NaV1.4EPN channel variant naturally occurring in Th. atratus was generated by site-directed mutagenesis (QuikChange Lightning II, Agilent, Santa Clara Co., CA, USA; now Stratagene, La Jolla, CA, USA). Using Agilent’s online primer design tool, primers were designed to implement the following changes: D1277E and A1281P in DIII and D1568N in DIV. In the pIRES-hrGFP-2a-SCN4A clone, these amino acid substitutions were achieved by t4640g, g4650c, and g5511a18. Mutagenic primers were synthesized, and their sequences can be found in Supp. Table 2.
Another variant was generated by Mutagenex (Suwanee, GA, USA) to include the mutations found in TTX-resistant Th. sirtalis: NaV1.4LVNV, I1556L, I1561V, D1568N, and G1569V in DIV17. In the pIRES-hrGFP-2a-SCN4A clone, these amino acid substitutions were achieved by a5475c, a5490g, g5511a, and g5515t. Mutagenic primer designs for this mutant were not reported, but the final sequence was confirmed. All three plasmid preps were then amplified by the same service (plasmid.com, Fargo, ND, USA) to obtain comparable DNA quality. Mutagenesis verification for all constructs has been submitted to GenBank (access: ON609825-ON609830).
Transfection
Cells were transiently transfected pursuant to the Lipofectamine 3000 protocol (ThermoFisher Scientific, Waltham, MA, USA). HEK293 cells were passed onto a 35 mm dish 24 h in advance to be 70% confluent at the time of transfection. Cells were transfected with 1 µg of pIRES-hrGFP-2a-rSCN4A or its mutants and allowed to translate protein for 48 to 72 h. This protocol was optimized on the wild-type clone and produced indistinguishable fluorescence reporter activity across all variants tested. Transfections took place within a 4-week period, and the plasmid variants were transfected in random order throughout the experimental period to avoid time-related effects on expression.
Electrophysiology
Clampex of the pCLAMP 10.7 software suite and an Axopatch 200A amplifier with CV202A headstage (Molecular Devices, Sunnyvale, CA, USA) were used to design and carry out voltage-clamp protocols and record INa. An online 10 kHz Bessel filter was applied to the signal, and the resultant INa was digitized at 50 kHz (Axon Digidata 1440A, Molecular Devices, Sunnyvale, CA, USA). Pipettes were fashioned from (1.5 mm OD, 1.1 mm ID) borosilicate capillary glass (BF150-110-7.5, Sutter Instruments Company, Novato, CA, USA). All patch pipettes and sharp electrodes were pulled on a P-87 Flaming/Brown Micropipette Puller (Sutter Instruments), and patch pipettes were polished on a Narishige MF-830 microforge. The pipette electrode resistance varied from 2-5 MΩ with 2.8 MΩ median resistance. The range of cell capacitance was measured to be between 5 and 14 pF. Series resistance was compensated to the edge of oscillation, and tau was set to 3 μs for all cells.
Cells were enzymatically dissociated and plated at 10% confluency. Successfully transfected single HEK293 cells were identified by green fluorescence during irradiation with 488 nm light filtered from a mercury lamp. Prior to making membrane contact, a correction was made for pipette capacitance and a liquid junction potential between internal and external solutions (-4 mV). The whole-cell patch clamp configuration was used to record macroscopic Na+ currents (INa). Following the attainment of a gigaOhm seal, a five-minute rest period was allowed for dialysis of the internal solution and further stabilization. The ambient temperature was 23 ± 2 °C. Series resistance was compensated to ≥90%, τ = 1-5 μs.
The membrane surface area (pF) was estimated by integrating the capacitive current generated by a brief +10 mV pulse. The voltage dependence of INa (I-V) was assessed by 10 ms test pulses from -80 to +80 mV in 10 mV increments. The test pulses followed a 500 ms conditioning period at -100 mV to remove inactivation. Following the test pulse, the membrane potential was returned to -100 mV to provide a symmetric but opposite capacitive current for reference. Currents were normalized to the membrane surface area to provide current density (pA/pF).
Channel kinetics
Steady-state activation (SSA) curves were constructed from the voltage dependence of GNa. This follows GNa= INa,max/ (Em – Erev), where E is the membrane potential during the test pulse and Erev is the empirical reversal potential. Steady-state inactivation (SSI) curves were constructed by plotting INa,max(pA/pF) against each prepulse potential. The resulting points were fit to the Boltzmann equation G/GNa,max or INa,max = 1/[1+exp(E1/2-Em)/kv], where E1/2 is the voltage that half-maximally activates or inactivates an NaV1.4 population and kv is the slope factor. To pool curves from multiple cells, these curves were normalized to the upper bound and forced to approach 0% (for SSI) and 100% (SSA). Window currents are a byproduct of overlapping SSA and SSI behaviors whereby small populations of sodium channels tend to linger in the open state. Differences here could be indicative of hazardous biophysical trade-offs. The E1/2 and kv parameters from the SSA and SSI curves can predict the NaV1.4 open probability during the window period. This is done using the equation Popen,window = [1/(1+e((E1/2,SSA – Em)/kv,SSA)) • [1/(1+e((Em-E1/2,SSI)/kv,SSI)).
Assessment of recovery from inactivation
A paired-pulse protocol was applied to the membrane to assess the time and voltage dependence of recovery from inactivation. The population was first equilibrated (500 ms conditioning prepulse at -100, - 80, or -60 mV) and activated (10 ms test pulse to +10 mV, at peak SSA for all mutants), thereby engaging inactivation. The membrane was then returned to the prepulse potential before a second test pulse (10 ms, +10 mV), with the duration between test pulses serially increasing until the peak current was completely rescued. The ratio (INa,max)test-2/(INa,max)test-1 was plotted against the interpulse interval. The time constant of the recovery from inactivation (τrfi) was calculated by a single exponential fit (OriginPro 2017). Comparisons were made between NaV1.4+, NaV1.4EPN, and NaV1.4LVNV for each prepulse potential.
Onset of fast inactivation
The membrane was conditioned at -100 mV for 500 ms to remove inactivation before each sweep of this protocol. An initial prepulse potential at -30 mV for 0.2 ms preceded a test pulse at 0 mV for 10 ms. All alleles of NaV1.4 have been shown to completely inactivate for sufficiently long exposures to this depolarized potential. The time course of the onset of this inactivation was assessed by repeating this command 100 times, each time increasing the prepulse potential by 0.5 ms. The INa,max of the current developed following the 0.2 ms pulse at -30 mV was taken as the maximum. The maxima of all sweeps were reported as a fraction of this first peak. The fraction of available channels was then plotted against the duration of the inactivating pulse. The time constant of the onset of fast inactivation (τofi) was calculated from a single exponential function and is the basis of comparison between alleles.
Fluctuation Analysis
The cell membrane was conditioned at -100 mV for 20 ms prior to a brief 10 ms pulse to +10 mV. NaV1.4+, NaV1.4EPN, and NaV1.4LVNV are known to have saturated steady-state activation at this test voltage. This protocol was repeated 1000 times, and the subsequent mean current (INa(t), pA/pF) was plotted against the variance over time (σINa(t)2, pA/pF2). The data were fit by a parabola of the equation: σINa(t)2 =[INa(t)]i – [INa(t)]2/N, where t spans from the time INa,max occurs to the end of the test pulse (∼8.5-9.5 ms), i refers to the predicted unitary current of NaV1.4 and N refers to the predicted number of such channels in the tested membrane.
Resistance of NaV1.4 Alleles
The cell membrane was conditioned at -100 mV for 500 ms prior to a 10 ms test pulse to +10 mV. TTX was applied by perfusing aliquots of extracellular solution to which known volumes of TTX diluted in citrate buffer pH 6.8 were added. Added volumes were kept below 0.01% of the aliquot volume. The cell was allowed to equilibrate for two (2) minutes at the intended [TTX] (in nM: 0, 100, 200). Each current trace was then normalized to INa,max (pA/pF) of the 0 nM trace.
Data analysis
Currents were analyzed with pCLAMP software 10.7 (Molecular Devices, Sunnyvale, CA, USA). Important metrics were recorded in Microsoft Excel (Redmond, WA, USA). Statistical tests and tests for approximate normality were performed in R (R Core Team 2017). Graphs were constructed in OriginPro 2017 (OriginLab Corp, Northampton, MA, USA).
Sodium Channel Molecular Modeling and Interaction Energy Estimation
To offer molecular insights into the differences between skeletal muscle voltage-gated sodium channels with different pore-domain residues, we turned to published crystal structures. The structure of the human voltage-gated sodium channel NaV1.4 in complex with β1 was used as a starting point (pdb 6agf, UniProt P35499). The sequence of the α subunit was then aligned (Blastp) to a consensus sequence of Thamnophis spp. skeletal muscle sodium channels from toxin-sensitive Thamnophis sirtalis, Th. atratus, and Th. elegans. The alignment revealed 356 differences between the Homo and Thamnophis sodium channels. In UCSF Chimera, all 356 point mutations were made to the human structure to generate a homology model that approximates the snake channel structure (ignoring differences in the C and N termini)90. This “serpentized” homology model therefore uses the human residue positioning numbers. Therefore, to recreate the mutations found in TTX-resistant snake NaV1.4 Domains III and IV of Th. atratus (D1277E, A1281P, and D1568N), the toxin-sensitive model was mutated at D1248E, A1252P, and D1539N. To recreate the Domain IV mutations of Th. sirtalis (I1556L, I1561V, D1568N, and G1569V), the toxin-sensitive model was mutated at I1527L, I1532V, D1539N, and G1540V. All three of these models were then run through UCSF Chimera’s standard Minimize Structure routine to produce three optimized, Thamnophis homology models for and NaV1.4+, NaV1.4EPN, and NaV1.4LVNV. Each of these models was then uploaded to the Interaction Energy Matrix (IEM) Web Application (https://ip-78-128-251-188.flt.cloud.muni.cz/energy/)71–73. For all six positions in question (residues 1248, 1252, 1527, 1532, 1539, and 1540), the interaction energy of each mutation was calculated by the IEM program for each minimized model, accounting for both Coulombic and Lennard-Jones modeled Van der Waals forces. Furthermore, the total interaction energies of several key residues linking the pore-domain and the S5/S6 domains of the sodium channel, that in turn link to the inactivation gate, were also calculated. The ΔΔG was calculated here using the ancestral, toxin-sensitive model as the reference. Residues that carry a more negative interaction energy (i.e. more stable) than the ancestral residue produce negative values while substitutions that destabilize the interaction energy relative to the ancestral protein produce positive values. Here, key residues were identified due to their known functional significance or by physical proximity when two atoms on adjacent residues lie within 5Å of one another.
Solutions
All experiments on skeletal muscle tissue were performed in Krebs buffer of the following composition (in mM): 145 Na+, 127.7 Cl⁻, 25 HCO3⁻, 5.5 glucose, 5.4 K+, 1.2 Mg2+, 1.8 Ca2+, 1.2 H2PO4⁻, 1.2 SO 2⁻ and supplemented with 95% O and 5% CO to achieve pH ∼7.35. For whole-cell patch clamp recordings of HEK293 cells transfected with pIRES-hrGFP-2a-SCN4A, the pipette solution contained (in mM): 10 Na+, ≥120 Cs+, 110 aspartate-, 32 Cl⁻, 12 tetraethylammonium+, 5 EGTA4⁻, 5 HEPES, 3 Mg2+, and 1.5 ATP4⁻, adjusted to pH 7.2 with CsOH. The extracellular solution was formulated (in mM): 145 Na+, 160.6 Cl⁻, 10 tetraethylammonium+, 5.5 glucose, 5 HEPES, 1.8 Ca2+, and 1 Mg2+, adjusted to pH 7.4 with CsOH (final [Cs+] ∼ 1.6 mM). The measured liquid junction potential between the external and internal patch solutions was -4.4 mV (n = 18, data not presented), for which the pipette offset was adjusted to +4.4 mV prior to cell contact.
Statistics
Mean data are presented alongside the standard error of the mean (sem). However, all statistical tests were performed using nonparametric analysis of variance in R 3.5.1 (R Core Team 2017) using Kruskal-Wallis ANOVA and Dunn’s multiple comparison post hoc test. Nonparametric methods were employed for all tests given that not all data were normally distributed. Statistical significance was determined by comparing the control and mutant groups using the Kruskal-Wallis test with α = 0.05. Preliminary experiments used to secure funding suggested that the effect size between resistant and sensitive animals would be sufficiently large to forego an effect size calculation or other power analysis. All test statistics and observation counts (n) can be found in Supplemental Tables 3-11 in the data repository below. All p values presented were adjusted to account for multiple tests by the conservative Bonferroni correction method. Center lines on box-and-whisker plots represent medians bounded by third and first quartiles; upper caps represent the 95th percentile, while lower caps represent the 5th percentile. Error bars represent the standard error.
Data Availability Statement
Myography, whole-animal assays, and electrophysiology data that support the results can be found at https://github.com/rdelcarlo/Snake_Biophysical_Trade-off.git, in addition to tabular reports on statistical tests performed to form the conclusions presented here. All data sets and/or analyses generated in this study are available from the corresponding author upon reasonable request. All custom codes used to process these data can be found at the listed github.com.
Competing Interest Statement
The authors indicate no conflicts of interest.
Acknowledgements
We thank CAF&W for permits to C.R.F, M.T.J.H. and E.D.B. III, and we thank USU and UNR IACUCs for approval of live animal protocols. We thank Mike Edgehouse, Vicki Thill, Erica Ely, and Kevin Wiseman for aid with field collections, and Jeff Wilcox, John Bailey, Cathy Koehler, Pete Steel, and Margot Rawlins for access to field sites. We are grateful to Gabrielle Blaustein, Vicki Thill, Amber Durfee, Kenzie Wasley, Taylor Disbrow, Sage Kruleski and Aubrey Smith for aid in live animal care at UNR. We thank Mohammed Chahine for providing the wild-type rat NaV1.4 clone. Aspects of this work would not have been possible without the help of Christopher Lingle who provided key aid on non-stationary noise analysis. We thank Kimberly Stanek for producing crystal structure models of the sodium channel and Devon Picklum for snake photographs. Finally, we thank UNR Evol Doers for thoughtful discussions and useful feedback on this manuscript, especially Marjorie Matocq and Jamie Voyles. This work was supported by an NSF grant to C.R.F. and N.L. (IOS-1355221) and NIH grants to N.L. (R01 HL091238 and P20GM130459) and by an NIH grant to S.A. (R01 HL161122).
References
- 1Integrated PhenotypeIntegrative and Comparative Biology 52:64–76https://doi.org/10.1093/icb/ics043
- 2Darwinian evolution can follow only very few mutational paths to fitter proteinsScience 312:111–114https://doi.org/10.1126/science.1123539
- 3Causes of molecular convergence and parallelism in protein evolutionNature Reviews Genetics 17:239–250https://doi.org/10.1038/nrg.2016.11
- 4Constraints on phenotypic evolutionAmerican Naturalist 140:S85–S107https://doi.org/10.1086/285398
- 5Missense meanderings in sequence space: A biophysical view of protein evolutionNature Reviews Genetics 6:678–687https://doi.org/10.1038/nrg1672
- 6The spandrels of San Marco and the Panglossian paradigm: a critique of the adaptationist programmeProc R Soc Lond B Biol Sci 205:581–598
- 7Quantitative genetics and developmental constraints on evolution by selectionJ Theor Biol 110:155–171https://doi.org/10.1016/s0022-5193(84)80050-8
- 8Evolutionary constraint and ecological consequencesEvolution 64:1865–1884https://doi.org/10.1111/j.1558-5646.2010.00960.x
- 9Causes and evolutionary significance of genetic convergenceTrends in Genetics 26:400–405https://doi.org/10.1016/j.tig.2010.06.005
- 10Tarichatoxin--Tetrodotoxin: A Potent NeurotoxinScience 144:1100–1110
- 11Investigations on skin toxin of adult rough-skinned newt Taricha granulosaCopeia
- 12The occurrence of tetrodotoxin (tarichatoxin) in amphibia and the distribution of the toxin in the organs of newts (Taricha)Toxicon 3:195–203https://doi.org/10.1016/0041-0101(66)90021-3
- 13The evolutionary origins of beneficial alleles during the repeated adaptation of garter snakes to deadly preyProceedings of the National Academy of Sciences of the United States of America 106:13415–13420https://doi.org/10.1073/pnas.0901224106
- 14Recovery of garter snakes (Thamnophis sirtalis) from the effects of tetrodotoxinJournal of Herpetology 36:95–98https://doi.org/10.1670/0022-1511(2002)036[0095:rogsts]2.0.co;2
- 15Parallel arms races between garter snakes and newts involving tetrodotoxin as the phenotypic interface of coevolutionJournal of Chemical Ecology 31:343–356https://doi.org/10.1007/s10886-005-1345-x
- 16Phenotypic outcomes of predator-prey coevolution are predicted by landscape variation in climate and community compositionFunctional Ecology 37:2170–2180https://doi.org/10.1111/1365-2435.14360
- 17Evolutionary diversification of TTX-resistant sodium channels in a predator-prey interactionNature 434:759–763https://doi.org/10.1038/nature03444
- 18Genetic architecture of a feeding adaptation: garter snake (Thamnophis) resistance to tetrodotoxin bearing preyProceedings of the Royal Society B-Biological Sciences 277:3317–3325https://doi.org/10.1098/rspb.2010.0748
- 19Historical Contingency in a Multigene Family Facilitates Adaptive Evolution of Toxin ResistanceCurrent Biology 26:1616–1621https://doi.org/10.1016/j.cub.2016.04.056
- 20Convergent adaptation to dangerous prey proceeds through the same first-step mutation in the garter snake Thamnophis sirtalisEvolution 71:1504–1518https://doi.org/10.1111/evo.13244
- 21The geographic mosaic of arms race coevolution is closely matched to prey population structureEvol Lett 4:317–332https://doi.org/10.1002/evl3.184
- 22The tetrodotoxin binding site is within the outer vestibule of the sodium channelMar Drugs 8:219–234https://doi.org/10.3390/md8020219
- 23Constraint shapes convergence in tetrodotoxin-resistant sodium channels of snakesProceedings of the National Academy of Sciences of the United States of America 109:4556–4561https://doi.org/10.1073/pnas.1113468109
- 24III Large-effect mutations generate trade-off between predatory and locomotor ability during arms race coevolution with deadly preyEvolution Letters 2:406–416https://doi.org/10.1002/evl3.76
- 25Predictably Convergent Evolution of Sodium Channels in the Arms Race between Predators and PreyBrain Behavior and Evolution 86:48–57https://doi.org/10.1159/000435905
- 26Ion Channels of Excitable MembranesSinauer Associates
- 27Evidence for electrical transmission in nerve: Part IJ Physiol 90:183–210
- 28Evidence for electrical transmission in nerve: Part IIJ Physiol 90:211–232
- 29Structure and function of voltage-gated sodium channels at atomic resolutionExperimental Physiology 99:35–51https://doi.org/10.1113/expphysiol.2013.071969
- 30Acute Oral Toxicity of Tetrodotoxin in Mice: Determination of Lethal Dose 50 (LD50) and No Observed Adverse Effect Level (NOAEL)Toxins 9https://doi.org/10.3390/toxins9030075
- 31Tetrodotoxin: a brief historyProc Jpn Acad Ser B Phys Biol Sci 84:147–154https://doi.org/10.2183/pjab.84.147
- 32Adaptive evolution of voltage-gated sodium channels: The first 800 million yearsProceedings of the National Academy of Sciences of the United States of America 109:10619–10625https://doi.org/10.1073/pnas.1201884109
- 33Tetrodotoxin levels of the rough-skin newt, Taricha granulosa, increase in long-term captivityToxicon 40:1149–1153https://doi.org/10.1016/s0041-0101(02)00115-0
- 34Variations in tetrodotoxin levels in populations of Taricha granulosa are expressed in the morphology of their cutaneous glandsScientific Reports 9https://doi.org/10.1038/s41598-019-54765-z
- 35Adaptive evolution of tetrodotoxin resistance in animalsTrends in Genetics 22:621–626https://doi.org/10.1016/j.tig.2006.08.010
- 36Toxin-resistant sodium channels: parallel adaptive evolution across a complete gene familyMol Biol Evol 25:1016–1024https://doi.org/10.1093/molbev/msn025
- 37Convergent and parallel evolution in a voltage-gated sodium channel underlies TTX-resistance in the Greater Blue-ringed Octopus: Hapalochlaena lunulataToxicon 170:77–84https://doi.org/10.1016/j.toxicon.2019.09.013
- 38Gene Conversion Facilitates the Adaptive Evolution of Self-Resistance in Highly Toxic NewtsMolecular Biology and Evolution 38:4077–4094https://doi.org/10.1093/molbev/msab182
- 39On the structural basis for ionic selectivity among Na+, K+, and Ca2+ in the voltage-gated sodium channelBiophysical Journal 71:3110–3125
- 40Control of ion flux and selectivity by negatively charged residues in the outer mouth of rat sodium channelsJournal of Physiology-London 491:51–59
- 41Structure of the sodium channel pore revealed by serial cysteine mutagenesisProceedings of the National Academy of Sciences of the United States of America 93:300–304https://doi.org/10.1073/pnas.93.1.300
- 42Inherited disorders of voltage-gated sodium channelsJ Clin Invest 115:1990–1999https://doi.org/10.1172/JCI25505
- 43Biophysical costs associated with tetrodotoxin resistance in the sodium channel pore of the garter snake, Thamnophis sirtalisJ Comp Physiol A Neuroethol Sens Neural Behav Physiol 197:33–43https://doi.org/10.1007/s00359-010-0582-9
- 44Costs of exploiting poisonous prey: Evolutionary trade-offs in a predator-prey arms raceEvolution 53:626–631
- 45.Na Channels from Phyla to Function Vol. 78 Current Topics in MembraneElsevier Academic Press Inc :87–113
- 46The evolutionary response of predators to dangerous prey: Hotspots and coldspots in the geographic mosaic of coevolution between garter snakes and newtsEvolution 56:2067–2082
- 47The geographic mosaic in parallel: Matching patterns of newt tetrodotoxin levels and snake resistance in multiple predator-prey pairsJ Anim Ecol 89:1645–1657https://doi.org/10.1111/1365-2656.13212
- 48The relation of membrane changes ot contraction in twitch muscle fibresJ Physiol 201:589–611https://doi.org/10.1113/jphysiol.1969.sp008774
- 49The resting membrane potential of frog sartorius muscleJ Physiol 297:1–8https://doi.org/10.1113/jphysiol.1979.sp013024
- 50The membrane potential of rat diaphragm muscle fibres and the effect of denervationJ Physiol 255:651–667https://doi.org/10.1113/jphysiol.1976.sp011301
- 51Resting-State Structure and Gating Mechanism of a Voltage-Gated Sodium ChannelCell 178https://doi.org/10.1016/j.cell.2019.06.031
- 52Excitation-contraction coupling in skeletal muscle: recent progress and unanswered questionsBiophys Rev 12:143–153https://doi.org/10.1007/s12551-020-00610-x
- 53Molecular-properties of voltage-sensitive sodium-channelsAnnual Review of Biochemistry 55:953–985https://doi.org/10.1146/annurev.bi.55.070186.004513
- 54Molecular-model of the action-potential sodium-channelProceedings of the National Academy of Sciences of the United States of America 83:508–512https://doi.org/10.1073/pnas.83.2.508
- 55Structural basis for gating charge movement in the voltage sensor of a sodium channelProceedings of the National Academy of Sciences of the United States of America 109:E93–E102https://doi.org/10.1073/pnas.1118434109
- 56Modeling the pore structure of voltage-gated sodium channels in closed, open, and fast-inactivated conformation reveals details of site 1 toxin and local anesthetic bindingJournal of Molecular Modeling 12:813–822https://doi.org/10.1007/s00894-005-0066-y
- 57Evolutionarily conserved intracellular gate of voltage-dependent sodium channelsNature Communications 5https://doi.org/10.1038/ncomms4420
- 58Structures of closed and open states of a voltage-gated sodium channelProceedings of the National Academy of Sciences of the United States of America 114:E3051–E3060https://doi.org/10.1073/pnas.1700761114
- 59Molecular dissection of multiphase inactivation of the bacterial sodium channel NaVAbJournal of General Physiology 151:174–185https://doi.org/10.1085/jgp.201711884
- 60Single-channel analysis of inactivation-defective rat skeletal muscle sodium channels containing the F1304Q mutationBiophys J 71:1285–1294https://doi.org/10.1016/S0006-3495(96)79329-3
- 61Fast and slow activation kinetics of voltage-gated sodium channels in molluscan neuronsJournal of Neurophysiology 77:2373–2384
- 62Rapid and slow voltage-dependent conformational changes in segment IVS6 of voltage-gated Na+ channelsBiophysical Journal 78:2943–2958
- 63Mechanisms of sodium channel inactivationCurrent Opinion in Neurobiology 13:284–290https://doi.org/10.1016/s0959-4388(03)00065-5
- 64Selectivity filter residues contribute unequally to pore stabilization in voltage-gated sodium channelsBiochemistry 44:13874–13882https://doi.org/10.1021/bi0511944
- 65The external pore loop interacts with S6 and S3-S4 linker in domain 4 to assume an essential role in gating control and anticonvulsant action in the Na+ channelJournal of General Physiology 134:95–113https://doi.org/10.1085/jgp.200810158
- 66Nonstationary noise-analysis of membrane currentsBiophysical Journal 33:A265–A265
- 67Voltage sensor of excitation-contraction coupling in skeletal musclePhysiol Rev 71:849–908https://doi.org/10.1152/physrev.1991.71.3.849
- 68Dihydropyridine receptors as voltage sensors for a depolarization-evoked, IP3R-mediated, slow calcium signal in skeletal muscle cellsJ Gen Physiol 121:3–16https://doi.org/10.1085/jgp.20028671
- 69Local control model of excitation-contraction coupling in skeletal muscleJ Gen Physiol 110:415–440https://doi.org/10.1085/jgp.110.4.415
- 70Structure of the human voltage-gated sodium channel Na(v)1.4 in complex with beta 1Science 362https://doi.org/10.1126/science.aau2486
- 71Identifying stabilizing key residues in proteins using interresidue interaction energy matrixProteins-Structure Function and Bioinformatics 72:402–413https://doi.org/10.1002/prot.21938
- 72Optimal Definition of Inter-Residual Contact in Globular Proteins Based on Pairwise Interaction Energy Calculations, Its Robustness, and ApplicationsJournal of Physical Chemistry B 116:12651–12660https://doi.org/10.1021/jp303088n
- 73Amino Acid Interaction (INTAA) web serverNucleic Acids Research 45:W388–W392https://doi.org/10.1093/nar/gkx352
- 74Molecular Dynamics Simulations of Ion Permeation in Human Voltage-Gated Sodium ChannelsJournal of Chemical Theory and Computation 19:2953–2972https://doi.org/10.1021/acs.jctc.2c00990
- 75Voltage-Gated Sodium Channels: Structure, Function, Pharmacology, and Clinical IndicationsJournal of Medicinal Chemistry 58:7093–7118https://doi.org/10.1021/jm501981g
- 76Discovery of a selective NaV1.7 inhibitor from centipede venom with analgesic efficacy exceeding morphine in rodent pain modelsProc Natl Acad Sci U S A 110:17534–17539https://doi.org/10.1073/pnas.1306285110
- 77Possible Roles of Exceptionally Conserved Residues around the Selectivity Filters of Sodium and Calcium ChannelsJournal of Biological Chemistry 286:2998–3006https://doi.org/10.1074/jbc.M110.175406
- 78Stable isotope analysis suggests that tetrodotoxin-resistant Common Gartersnakes (Thamnophis sirtalis) rarely feed on newts in the wildCanadian Journal of Zoology 99:331–338https://doi.org/10.1139/cjz-2020-0215
- 79.The Garter Snakes: Evolution and EcologyUniversity of Oklahoma Press
- 80The road not taken: Evolution of tetrodotoxin resistance in the Sierra garter snake (Thamnophis couchii) by a path less travelledMolecular Ecology 31:3827–3843https://doi.org/10.1111/mec.16538
- 81Phylogenomic analyses resolve relationships among garter snakes (Thamnophis: Natricinae: Colubridae) and elucidate biogeographic history and morphological evolutionMolecular Phylogenetics and Evolution 167https://doi.org/10.1016/j.ympev.2021.107374
- 82Resistance of neonates and field-collected garter snakes (Thamnophis spp.) to tetrodotoxinJ Chem Ecol 30:143–154https://doi.org/10.1023/b:joec.0000013187.79068.d2
- 83Tetrodotoxin resistance in garter snakes: an evolutionary response of predators to dangerous preyEvolution 44:651–659https://doi.org/10.1111/j.1558-5646.1990.tb05945.x
- 84Tarichatoxin: Isolation and purificationScience 140:295–296https://doi.org/10.1126/science.140.3564.295
- 85Repeated injections of TTX do not affect TTX resistance or growth in the garter snake Thamnophis sirtalisCopeia :531–535
- 86Clustal W and clustal X version 2.0Bioinformatics 23:2947–2948https://doi.org/10.1093/bioinformatics/btm404
- 87Phylogenetics - MacClade 4: Analysis of phylogeny and character evolution, version 4.06American Biology Teacher 66:511–512https://doi.org/10.1662/0002-7685(2004)066[0511:ctr]2.0.co;2
- 88Geneious Basic: An integrated and extendable desktop software platform for the organization and analysis of sequence dataBioinformatics 28:1647–1649https://doi.org/10.1093/bioinformatics/bts199
- 89Parallel Evolution of Tetrodotoxin Resistance in Three Voltage-Gated Sodium Channel Genes in the Garter Snake Thamnophis sirtalisMolecular Biology and Evolution 31:2836–2846https://doi.org/10.1093/molbev/msu237
- 90UCSF chimera - A visualization system for exploratory research and analysisJournal of Computational Chemistry 25:1605–1612https://doi.org/10.1002/jcc.20084
Article and author information
Author information
Version history
- Sent for peer review:
- Preprint posted:
- Reviewed Preprint version 1:
Copyright
© 2024, del Carlo et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics
- views
- 273
- downloads
- 20
- citations
- 0
Views, downloads and citations are aggregated across all versions of this paper published by eLife.