2023 年 20 巻 3 号 論文ID: e200030
Ion currents associated with channel proteins in the presence of membrane potential are ubiquitous in cellular and organelle membranes. When an ion current occurs through a channel protein, Joule heating should occur. However, this Joule heating seems to have been largely overlooked in biology. Here we show theoretical investigation of Joule heating involving channel proteins in biological processes. We used electrochemical potential to derive the Joule’s law for an ion current through an ion transport protein in the presence of membrane potential, and we suggest that heat production and absorption can occur. Simulation of temperature distribution around a single channel protein with the Joule heating revealed that the temperature increase was as small as <10–3 K, although an ensemble of channel proteins was suggested to exhibit a noticeable temperature increase. Thereby, we theoretically investigated the Joule heating of systems containing ensembles of channel proteins. Nerve is known to undergo rapid heat production followed by heat absorption during the action potential, and our simulation of Joule heating for a squid giant axon combined with the Hodgkin-Huxley model successfully reproduced the feature of the heat. Furthermore, we extended the theory of Joule heating to uncoupling protein 1 (UCP1), a solute carrier family transporter, which is important to the non-shivering thermogenesis in brown adipose tissue mitochondria (BATM). Our calculations showed that the Joule heat involving UCP1 was comparable to the literature calorimetry data of BATM. Joule heating of ion transport proteins is likely to be one of important mechanisms of cellular thermogenesis.
Ion currents through channel proteins in the presence of membrane potential occur everywhere in cellular and organelle membranes. Joule heating is the process of heat production by charge current in a material. Here we theoretically investigated the Joule heating involving ion channel and transporter proteins. We show that the Joule’s law holds for ion currents through an ion transport protein, and both heat production and absorption can occur. This Joule heating explains heat production and absorption of nerve during the action potential and gives an insight into the non-shivering thermogenesis mechanism of uncoupling protein 1 in brown adipose tissue mitochondria.
Thermogenesis plays a vital role in the homeostasis and regulation of body temperature in animals, which affects a plethora of molecular processes in the body. One type of thermogenesis occurs in skeletal muscle, which undergoes repetitive contraction to increase the metabolic rate above the basal level, denoted as shivering thermogenesis [1]. The other types of thermogenesis are denoted as non-shivering thermogenesis and include heat production attributable to basal metabolism and one to mitochondria in brown adipose tissue (BAT) [1,2]. The non-shivering thermogenesis in BAT, in particular, is thought to be associated with aspects in mammals such as cold adaptation, fever, hibernation, and obesity [3]. To better understand the relevance of non-shivering thermogenesis, it is desirable to elucidate how the heat is produced, although the details remain elusive.
Ion currents through channel proteins are ubiquitous in cellular and organelle membranes. Across the membranes, there often exists the difference of electric potentials between one side and the other or the membrane potential at, say, a magnitude of ~100 mV. Therefore, when an ion current occurs through a membrane, Joule heat is likely to be generated. However, to our knowledge, it seems that the Joule heating involving ion transport proteins such as channel proteins has been seldom considered significant as a biological thermogenesis mechanism. Herein, take as an example a situation in which a monovalent cation particle is transported across a membrane through a channel. If the magnitude of the membrane potential is 100 mV, then the amount of the Joule heat will be 0.1 eV or 9.6 kJ/mol, and this value amounts to, e.g., as large as 1/3 of the magnitude of the hydrolysis enthalpy of ATP (–30.11 kJ/mol) [4], where most ion pump proteins are driven by ATP. This suggests that the Joule heating associated with ion transport proteins should not be overlooked.
In this study, we report a theoretical investigation of Joule heating in biological thermogenesis. We begin with deriving the Joule’s law for an ion current through an ion transport protein based upon electrochemical potential, and we note that the ion current can bring about both heat production and heat absorption depending on its direction. We next examine temperature distribution around a channel protein which undergoes Joule heating by heat diffusion simulation. Finally, we extend the present theory of Joule heating to model systems composed of ensembles of ion channels and electrochemical potential-driven transporters [5]. We investigate the heat production of nerve with the action potential and the one involving uncoupling protein 1 in BAT mitochondria, to which both heat production and heat absorption are most likely to be relevant.
We consider a model shown in Figure 1, in which aqueous regions 1 and 2 are separated by a lipid membrane and a flow of ion takes place through a channel. Let μ1 and μ2 be the electrochemical potentials (J/mol) of the ion, V1 and V2 be electric potentials (V), and c1 and c2 be ion concentrations (mol/L) in regions 1 and 2, respectively. The electrochemical potentials are expressed as
Schematic diagram of ion transport between regions 1 and 2 with electric potential difference. The regions are separated by the lipid membrane but connected through the channel. The diagram depicts that Δn mol of the cation is transported from region 1 to region 2 through the channel in the sign convention of Δn>0. The variables of μ, V, and c are electrochemical potential, electric potential, and ion concentration, respectively.
(1) |
and
(2) |
where z is the ion valence, μ’º is a standard transformed chemical potential of the ion at a physiological condition, T is temperature (K), R=8.314 (J/K·mol) is the gas constant, and F=96,485 (C/mol) is the Faraday constant [6]. Herein, we suppose that we observe transport of the ion from region 1 to region 2 by an infinitesimally small amount of Δn (mol) for a sufficiently short time of Δt (>0) so that changes of c1 and c2 are negligibly small, where Δn>0 for μ1>μ2, Δn<0 for μ1<μ2, and Δn=0 for μ1=μ2 (for a specific description of ion current, see Eq (6)). Then, the change of the total Gibbs free energy of this model caused by the ion transport ΔGTr is written as,
(3) |
Herein, we have a look at enthalpy to evaluate the heat production of this ion transport process. By the Gibbs-Helmholtz equation [7], we derive the change of the total enthalpy due to the ion transport process ΔHTr as,
(4) |
where we operationally fixed c’s and V’s in the differential with respect to T by approximation, and the factor of Δn·z·F=Δq corresponds to the amount of transported charge (C). To clarify the meaning of Eq (4), as an example, consider a cation (z>0) and the electric potentials of V2<V1. If Δn>0 (i.e., the cation moves from region 1 with a higher electric potential to region 2 with a lower one), then we have ΔHTr<0, meaning exotherm or heat production. Meanwhile, if Δn<0 (i.e., the cation moves from region 2 to region 1), then we have ΔHTr>0, meaning endotherm or heat absorption. Furthermore, if we express the ion current (A) corresponding to the ion transport in the present model as I=Δq/Δt, then we have the change rate of the enthalpy hI (W) for I as
(5) |
Importantly, the most righthand side of Eq (5) represents the Joule’s law. We note that Eqs (4) and (5) can be applied to ion transport of both channels and electrochemical potential-driven transporters [5], because the derivation of these equations was independent on the detailed mechanism of ion transport. In addition, in electrophysiology, the ion current I may often be expressed in the form of,
(6) |
where G is an ion conductance, and VEq is an equilibrium potential or the Nernst potential for the ion between the regions 1 and 2 and is dependent on c1 and c2 [8]. Comparison between Eq (5) and Eq (6) suggests that it is necessary to distinguish V2–V1 and V2–V1–VEq.
Furthermore, ion transport across membrane is also dealt with by primary active transporters such as ion pumps, molecular motors, synthases, and oxidoreduction-driven transporters [5], which involve the coupled reactions and, most likely, Joule heating. In this case, the enthalpy of the coupled reaction should be added to the Joule heat, and thus, the enthalpy and its change rate in Eqs (4) and (5) should be modified as, respectively,
(7) |
and
(8) |
where ΔrH is the reaction enthalpy of the coupled reaction and hr is its change rate.
On the molecular level, ions undergo collision and interaction with the surrounding water molecules in aqueous solution or with atoms in a channel protein at physiological temperature. This leads to Brownian motion of the ions together with the surrounding water molecules or the atoms, which would be responsible for the electric resistance for the ion current. In aqueous solution, ions such as Na+ and Cl– dynamically interact with the surrounding water molecules [9], indicating Brownian motion of the ions with the water molecules in a concerted manner. As for ions in channel proteins, according to a review paper on molecular dynamics (MD) studies [10], the most direct MD simulation of ion flux can be achieved by solving the Langevin equation with the spatially varying free energy and spatially varying diffusion coefficient, which also indicates frequent interaction of the ions with the atoms in the channel pore and thus Brownian motion. Additionally, a MD study of antiamoebin ion channel, for example, reported that the motion of ions inside the channel was satisfactorily described as diffusive and the diffusion coefficients of K+ and Cl– inside the channel were computed to be 0.46·10–5 cm2/s and 0.70·10–5 cm2/s, respectively [11], whereas those in water are 1.95·10–5 cm2/s and 2.03·10–5 cm2/s, respectively [12], which are comparable to those in the channel. As a thought experiment, if a cation under an electric potential gradient moves to a position of lower electric potential by Brownian motion, then the cation gains the kinetic energy, but relaxation would immediately take place so that the energy gain dissipates to the surrounding water molecules in the aqueous solution or the atoms in the channel pore, resulting in exotherm, or heat production. Conversely, if a cation moves to a position of higher electric potential by Brownian motion, then the cation loses the kinetic energy, but relaxation would immediately take place so that the kinetic energy of the surrounding water molecules in the aqueous solution or the atoms in the channel pore is distributed to the cation, resulting in endotherm or heat absorption. For an ion current through an ion transport protein, the driving force of the ion current is not just the electric potential gradient but also the ion concentration gradient (and free energy liberated from a coupled chemical reaction in the case of primary active ion transporters), suggesting that the ion current may oppose or follow the electric potential gradient (e.g., as described by Eq (6)). Accordingly, we suggest that it is possible that both exothermic and endothermic Joule heating can occur in cells.
Simulation of Heat DiffusionTo simulate time-dependent temperature distribution in a space containing a heat source, we used a continuum model with the Fourier’s law to numerically solve a heat diffusion equation given by,
(9) |
where ΔT and Q are a temperature increase (K) with respect to the initial state of t≤0 and heat supply (W/m3), respectively, as functions of time t (s) and position r in a 3-dimension space, and the whole model space is assumed to be in thermal equilibrium for t≤0, i.e., the temperature T(t≤0, r) is uniform, and ΔT(t≤0, r)=0 for any position r. The variables Cp, ρ, and κ are heat capacity (J/K·kg), density (kg/m3), and thermal conductivity (W/m·K), respectively, as functions of position r. To solve the heat diffusion equation, we used a numerical differential equation solver of NDSolve in Mathematica software (Version, 13.1.0.0) with the Open Cascade Technology modules for the finite element method (FEM). The detailed information of the heat diffusion simulation is provided in Supplementary Text S1 and Supplementary Figure S1A, B, and the parameter values are provided in Supplementary Table S1.
Simulation of Joule Heating of a Nerve Fiber Based upon the Hodgkin-Huxley ModelWe used the Hodgkin-Huxley model [13] to compute the ion currents (μA/cm2) of INa, IK, and IL for Na+, K+, and leak ions, respectively, and the membrane potential of V’ (mV) in a nerve cell. In this study, we used the sign convention of inward negative for the ion currents. The membrane potential V’ is the intracellular electric potential with respect to the intracellular resting potential VR, and we used a value of VR=–60 mV according to ref [13]. In this convention, the depolarization and hyperpolarization correspond to positive and negative values of V’, respectively. Note that the actual intracellular potential with respect to the extracellular potential is given by V’+VR. We considered two types of Joule heat (nW/cm2) in the nerve cell: heat due to the ion movement across the cell membrane through channels (Q1) (nW/cm2) and heat due to ion movement along the longitudinal axis of the nerve (Q2) (nW/cm2) according to ref [14] with modification. By Eq (5), we calculated Q1 as,
(10) |
where gi and Vi are the maximum ion conductance (mS/cm2) and the equilibrium potential (mV), respectively, for ion species i, m, h, and n are the channel gates, and t is time (ms). The Joule heat of Q2 was calculated by
(11) |
where a is a radius (cm) of the cross section of the nerve cell, RCell is an intracellular specific resistance (Ω cm), and θ is a propagation velocity (cm/ms) of the action potential (the derivation is explained in Supplementary Text S2 and Supplementary Figure S2). The factor of 103 is inserted into the righthand of Eq (11) for consistency of the unit prefix with Eq (10).
We considered a model including a single channel protein molecule embedded in a lipid membrane with surrounding solvent (Figure 2A) to estimate the temperature distribution. As an exploratory condition, we assumed a membrane potential of –100 mV across the lipid membrane (with respect to the outside) and an ion current of –10 pA in the inward negative convention which brings about Joule heating at 1 pW in the channel region of the protein (Eq (5), Figure 2A, B). We used the simplified configuration including the channel region in the channel protein as a heater element and other thermal conducting components as an approximation (Figure 2C). To be precise, the electric potential around a channel protein and lipid membrane should be more complex than that shown in Figure 2A, C and their structures should also be more complicated as well [10,15]. However, considering that thermal diffusion coefficients of water, lipid membrane, and proteins are in the order of α~10–8 m2/s [12,16–20], the simplified configuration in Figure 2C would be sufficient to crudely examine the temperature distribution around the channel protein, especially for the outside of the channel on the time scale of >>μs, based upon the mean square distance: (6·α·t)1/2≈(6·10–8·10–6)1/2=0.2·10–6 m, where α and t are a thermal diffusion coefficient and time, respectively.
Heat diffusion simulation involving a single channel protein. (A) A schematic diagram of a channel protein embedded in a lipid membrane in the presence of a membrane potential of –100 mV. (B) An analogous circuit for the ion current involving the channel protein. (C) Spatial configuration of the heater element and thermal conducting materials for the heat diffusion simulation. (D) Distribution of temperature increase around the channel 1 μs after the heater element was turned on. M, lipid membrane; Ch, channel region; P, protein region of the channel. (E) Profiles of temperature increase taken from the lines indicated as “E” in panel (D). z refers to the z-coordinates of the lines in panel (D). (F) Profiles of temperature increase taken from the lines indicated as “F” in panel (D). x refers to the x-coordinates of the lines in panel (D). (G) Time trajectories of temperature increase at points indicated as “G,H” in panel (D). (H) Plots of temperature increase at points indicated as “G,H” in panel (D) as a function of the heating power of the channel. ΔT is a temperature increase with respect to the zero-time temperature.
As we set up the model arrangement of the thermal conducting materials and the heater element (Figure 2C), we performed heat diffusion simulation of the model (Eq (9)). By numerical calculations with a partial differential equation solver, we calculated the solution of the time-dependent distribution of temperature increase with respect to the initial state, ΔT, after the heater element of the channel was turned on. The spatial distribution of ΔT at 1 μs showed a maximum at the center of the channel protein (Figure 2D–F), and the time trajectories of ΔT at the channel center and a point 1 nm off the channel showed that they reached plateaus within <<1 ns (Figure 2G). The plateau temperature at the channel center was ΔT=0.25 mK, whereas the temperature 1 nm off the channel ((x, y, z)=(0, 0, 2.25) in Figure 2D) was ΔT=0.016 mK (Figure 2E).
Furthermore, we tested the effects of changing thermal parameters on ΔT. The value of ΔT was linearly proportional to the heating power of the channel (Figure 2H), because the heat diffusion equation (Eq (9)) is a linear system so that ΔT and Q linearly correlate with each other. A value of ΔT=1 K at the channel center could be possible if the heating power was about 4 nW (Figure 2H), although this heating power should be impossible for a single channel protein molecule. Additionally, solvent thermal conductivity had little effect on the ΔT distribution in the range of 0.03–0.6 W/m·K (Supplementary Figure S1C). The thermal conductivity of the channel protein modestly affected ΔT (Supplementary Figure S1D). The data suggest that, in the present model, the rapid heat dissipation substantially counteracts the accumulation of heat in the channel protein and its vicinity of nm size, and therefore, ΔT reaches a steady state at as small as <1 mK.
Previously, Baffou et al [21] presented a relationship among the power of heater element (P), thermal conductivity (κ), the size (diameter) of the heater element (L), and the heater element’s temperature increase (ΔT) with respect to the temperature before heating, expressed as ΔT=P/(κ·L). This formula was derived from a steady-state solution of the heat diffusion equation for a spheric heater element in a thermal conducting medium. Because we used the heat diffusion equation in the present study without any special nm-scaled effects such as phonon and the Kapitza thermal boundary, the temperature increase ΔT calculated by the present heat diffusion model including the Joule heating should be consistent with the Baffou’s formula. In fact, for a single channel, if we assume that P=1 pW, κ=0.6 W/m·K, and L=2.5 nm, then we have ΔT=0.67 mK, which would be compatible with the present simulation result of ΔT=0.25 mK (Figure 2E–G) as being in the same order of magnitude. The small difference between these values would be mainly due to the effect of shapes for the components in the model.
Although temperature increase ΔT for a single channel protein molecule was very small (Figure 2E–G), an ensemble of channel proteins or an ensemble of cells may exhibit a noticeable value of ΔT. As a trial, if we suppose a cell with, e.g., L=10 μm, P=1 nW (corresponding to 1,000 channels of a heating power of 1 pW), and κ=0.6 W/m·K (the value for water), then ΔT is still as small as 0.17 mK. However, we may encounter a noticeable temperature increase of ΔT>1 K, if we consider an ensemble of cells. Suppose that an ensemble contains N highly-packed cells of a heating power P (W/cell) and the size of each cell is L (m), then the size of this ensemble is ~N1/3·L. Thus, we have an estimation of temperature increase as ΔT=N·P/(κ·N1/3·L)=N2/3·P/(κ·L). A BAT, the thermogenesis of which is relevant to Joule heating (see Joule Heat of the H+ Leak in Brown Adipose Tissue Mitochondria), was reported to produce heat at a power of P=0.8 nW/cell [22] and contain N=8·107 brown adipose cells with an average diameter of L=18 μm for cold-acclimated rat [23]. If we use these parameter values and κ=0.6 W/m·K, then the temperature increase was calculated to be ΔT=14 K by the Baffou’s formula [21], which would be comparable to a temperature increase of ΔT=5.1 K in rat brown adipose tissues with norepinephrine administration as measured by 129Xe magnetic resonance in vivo [24]. This suggests compatibility between the Baffou’s formula with the experimentally measured temperature of cellular ensembles. Additionally, irrespective of theoretical predictions of the very small temperature increases for single cells and single channels by the Baffou’s formula [21] and our heat diffusion equation model (ΔT~0.25 mK, Figure 2E–G), a number of temperature imaging studies of cells have reported intracellular temperatures >1 K higher than the ambient temperatures (this was reviewed in ref [25] with critiques), and temperature observations at the single cell level remain in controversy.
Previously, Chen et al [26] investigated permeating cations in a channel by hydrodynamic simulation and demonstrated that the temperature of the cation particle moving through the channel could be 20–70 K higher than the ambient temperature, but the values may seem very different from the present result of 0.25 mK at the channel center (Figure 2E, F). The hydrodynamic model [26] explicitly included single cation particles undergoing hydrodynamic flow through the channel pore to calculate the temperature of the cations themselves, whereas the present study assumed a protein-like material for the channel pore without explicit ion particles to calculate its temperature (Figure 2C). Thus, we suggest that the difference of the detailed arrangements between the hydrodynamic model [26] and our heat diffusion model including the Joule heating would be responsible for the difference of the calculated temperature results. Taking this into consideration, the present study would not rule out the simulation results by Chen et al [26], and we suggest that the temperature distribution at atomic resolution would be beyond the scope of the present study.
Heat from a Nerve CellOn the cellular membrane of a nerve, a large number of channel proteins are arranged so as to orchestrate ion currents and the membrane potential to initiate and propagate the action potential [27]. This means that the ion currents take place in the presence of membrane potential, which is likely to lead to Joule heating (Figure 3A). In fact, previously, a small but noticeable amount of heat has been observed from nerve fibers during the compound action potential with observed temperature changes of ~μK [28]. Thereby, we tried to theoretically investigate the Joule heating involving a nerve at work (Figure 3A).
Simulation of heat production in a squid giant axon based upon the H-H model [13]. (A) A schematic diagram showing Joule heating due to ion currents through channels (Q1 in panel (D)) and Joule heating due to ion movement along the longitudinal axis of the nerve (Q2 in panel (D)). (B) A time trajectory of the membrane potential as calculated by the H-H model. (C) Time trajectories of ion currents as calculated by the H-H model. (D) Time trajectories of heat production, Q1, Q2, and Q1+Q2 by the two mechanisms as depicted in panel (A). (E) A time trajectory of Q1+Q2 at a time resolution of 1 ms.
We used the Hodgkin-Huxley (H-H) model [13] to simulate the ion currents and the membrane potential during the action potential so that the result could be applied to the Joule heating (Eq (5)). We set up the coupled differential equations according to the H-H model with parameter values for a squid giant axon (Supplementary Text S2; Supplementary Tables S2 and S3) and reproduced on a computer the results of the numerical solutions of the membrane potential (Figure 3B) and the ion currents for Na+, K+, and leak ions (Figure 3C). Subsequently, we used the numerical solutions to examine two types of Joule heat production (Supplementary Text S2) [14]. One type of the Joule heat Q1 was attributed to ion currents passing through the channel proteins on the cellular membrane (Eq (10); Figure 3A). The other type of Joule heat Q2 was attributed to an ion current along the longitudinal axis of the nerve fiber driven by the propagation of the membrane potential, which was also computed by using the solution of the H-H model (Eq (11), Figure 3A). Thereby, we show time courses of Q1 and Q2 in a nerve cell in Figure 3D. The heat Q1 exhibited oscillation for t=0–2 ms, which reflected the rapid changes of conductance for Na+ and K+ and the membrane potential. The negative heat in Q1 for t>2 ms corresponds to heat absorption and reflects the outward small current of K+ with the polarized membrane potential (V’+VR<0 in Eq (10)). Additionally, the heat Q2 contributes to the positive heat production throughout.
We also show the time course of Q1+Q2 in Figure 3D. According to early reports of initial heat in non-myelinated nerves, a single impulse of the action potential was associated with a heat production followed by a heat absorption [28]. Having a closer look at Figure 3D, the heat of Q1+Q2 for t≤2 ms may correspond to the heat production and that for t≥2 ms may correspond to the heat absorption. In fact, by numerical integration, the heat production for t≤2 ms amounted to +5.4 nJ/cm2·impulse and the heat absorption for t≥2 ms amounted to –1.1 nJ/cm2·impulse, where cm2 refers to the surface area of the nerve. As a reference, previous studies reported that the amounts of the initial heat production and absorption from rabbit vagus nerves were reported to be +5.02 and –3.35 nJ/cm2·impulse, respectively [29], and those from pike olfactory nerves were +4.51 and –4.99 nJ/cm2·impulse, respectively [30]. Thus, the amounts of heat production and absorption estimated in the present study would be comparable to the experimentally observed values, although the difference between a single squid giant axon (H-H model) and fibers of rabbit vagus and pike olfactory multiple nerves should be borne in mind. In fact, the propagation velocity and the impulse width for a squid giant axon described in the H-H model were 1,880 cm/s and <7 ms [13], respectively, whereas those for rabbit vagus nerve at 5°C were 12 cm/s and >100 ms, respectively [29], and those for pike olfactory nerves at 0°C were 3.0 cm/s and 240 ms, respectively [30]. If the parameter values of the H-H model for these nerve fibers are determined in the future, the heat production and absorption are likely to be described by the present Joule heating theory more accurately.
Although our theory of Joule heating presented in this study covers not just the heat production but also the heat absorption, the time trajectory of Q1+Q2 (Figure 3D) may not look like ones reported previously [28]. For example, we estimated that the time constant of the response of temperature measurement in the experiment of ref [30] was ~41 ms, and the study calculated the deconvoluted time trajectory of heat production at a resolution of 10 ms. As a trial, we rearranged the time trajectory in Figure 3D using a time resolution of 1 ms (Figure 3E), which clearly showed the initial heat production followed by heat absorption. Thus, although the time scales were somewhat different, the present Joule heat simulation is likely to have reproduced, at least, the traits of the experimentally-observed time trajectories of heat in nerve (e.g., Figure 3 in ref [30]), which would underpin the present theory of Joule heating.
It has long been controversial what the origin of the positive and negative heat in nerve fibers during the action potential is. The mainly-considered theories were the condenser theory, the entropy changes of the cell membrane, and gating charge movement [28]. Recently, de Lichtervelde et al [31] exploited the Poisson-Boltzmann equation and the Maxwell equations to precisely combine the movement of ions surrounding the membrane (the condenser model) and the dielectric change of the membrane (the entropy theory), and they predicted that the depolarization and repolarization led to reversible heat production and absorption at 4–7 nJ/cm2·impulse [31], which was denoted as the electrostatic model. The heat predicted by the electrostatic model [31] was comparable to the amounts of heat calculated by the Joule heating of the present study (+5.4 and –1.1 nJ/cm2·impulse for heat production and absorption, respectively; see above). Thus, it is suggested that the mechanisms of the Joule heating and the electrostatic model contribute to the heat production and absorption in the nerve upon the action potential in comparable magnitudes [31].
Additionally, a theoretical calculation of the Joule heating from a nerve fiber based upon the H-H model was previously performed by Nogueira and Conde Garcia [14], but they computed the heat of channel currents in a form of Gi·(V’–Vi)2, where Gi and Vi are an ion conductivity and the equilibrium potential, respectively, for an ion species i. However, in our study, we corrected their theory as explained in the Section Joule Heating Involving Ion Flow and Electric Potential so that the ion current Gi·(V’–Vi) was multiplied by the membrane potential (V’+VR), which the ion particles experience when they move from one side of a membrane to the other side through a channel. In fact, by our correction, the heat production as well as the heat absorption was at least partly explained by the present framework of the Joule heating.
Joule Heat of the H+ Leak in Brown Adipose Tissue MitochondriaIn mammals, BAT is a major site that is involved in non-shivering thermogenesis [3]. The H+ leak through the inner mitochondrial membrane (IMM) mediated by UCP1 in BAT mitochondria (BATM) is thought to play a key role in the thermogenesis mechanism [2,32], although, to our knowledge, it has little been explicitly explained how the H+ leak produces the heat. Because the membrane potential across IMM is maintained around –200 mV (with respect to the intermembrane space (IMS)) [2], the UCP1-mediated H+ leak is likely to bring about Joule heating. Thereby, although UCP1 is an electrochemical potential-driven transporter [5] and a member of the mitochondrial carrier family (SLC25) [33] but not a channel protein, we extended the present theory of Joule heating to the H+ leak using Eq (4). Furthermore, there are a number of mitochondrial carrier family transporters [33] and primary active transporters [5] on IMM that work with the UCP1-mediated H+ leak, and thus, we also applied the theory of Joule heating to them together with UCP1 (Eqs (4) and (7); Supplementary Text S3).
Previously, Locke et al [34] successfully performed the simultaneous measurement of the membrane potential and H+ conductance in BATM prepared from cold-acclimated hamsters. They observed an increase in the H+ conductance and a modest decrease in the magnitude of the membrane potential upon palmitate stimulation, which were attributed to the H+ leak through UCP1 [2,34]. For crude estimation, we read the data of Figure 1 in ref [34] to use a membrane potential of VIMM=–220 mV and a H+ conductance increment of ΔCIMM=1.63 nmol min–1 mg–1 mV–1 for mitochondria stimulated with palmitate in the presence of ATP and pyruvate, where the unit of mg refers to the total mitochondrial protein. By Eq (5), we estimated the amount of the corresponding Joule heat as,
(12) |
where F=96,485 (C mol–1) is the Faraday constant. As a reference, a calorimetry study of BATM prepared from a cold-exposed rat reported that the amount of heat production in the mitochondria was 0.19–0.76 mW/mg protein [35], although the solution conditions were not the same as that in ref [34]. These data suggest that the amount of Joule heat in Eq (12) would be significant in comparison to the heat generated by BATM.
Because the energy for the mitochondrial heat production, including the Joule heat and metabolic reaction heat, fundamentally originates from metabolic substrates, we worked out the contribution of the Joule heat to the total heat by reaction enthalpy analysis. In the case of thermogenesis associated with the activation of UCP1, the fundamental energy sources can be the substrates for the β-oxidation with its precursor reaction (βOx) and the tricarboxylic acid cycle (TCA). Specifically, we considered palmitate and pyruvate as model substrates and calculated the reaction enthalpies of βOx, the pyruvate oxidation (PO), TCA, the electron transfer chain (ETC), ATP synthesis (AS), and the H+ leak through UCP 1 (Figure 4A). In the UCP1-activated state, a reaction pathway of βOx, TCA, and ETC followed by the H+ leak (βOx→TCA→ETC→H+ leak) is thought to be predominant, whereas, in the resting state, a reaction pathway of PO, TCA, and ETC followed by AS (PO→TCA→ETC→AS) is thought to be predominant [2].
Enthalpy analysis of the metabolism and the Joule heating in BAT mitochondria. (A) Metabolic reactions considered in the enthalpy calculations. Proteins colored in blue are assumed to contribute to the endothermic Joule heating, and proteins colored in orange are assumed to contribute to the exothermic Joule heating. FFA, free fatty acid; ACtn, acylcarnitine; CAT, carnitine-acylcarnitine translocase; βOx, β-oxidation; Pyr, pyruvate; Pi, inorganic phosphate; MPyC, mitochondrial pyruvate carrier; ETC, electron transfer chain; I–IV, Complexes I–IV; ANT, adenine nucleotide translocator; ASy, ATP synthase; MPiC, mitochondrial phosphate carrier; UCP1, uncoupling protein 1; TCA, tricarboxylic acid cycle; IMS, intermembrane space; IMM, inner mitochondrial membrane. (B) A bar graph showing the reaction enthalpies including the Joule heat for βOx→TCA→ETC, βOx→TCA→ETC→H+ leak, and βOx→TCA→ETC→AS at an ionic strength of 0.12 M, pH 7.3, pMg 3, and 298.15 K. Note that negative enthalpy means exotherm. (C) A bar graph showing reaction enthalpies including the Joule heat for PO→TCA→ETC, PO→TCA→ETC→H+ leak, and PO→TCA→ETC→AS at an ionic strength of 0.12 M, pH 7.3, pMg 3, and 298.15 K.
As for palmitate, we combined reactions of acyl-CoA synthase, carnitine acyltransfrase, carnitine-acylcarnitine translocase, TCA, and ETC [36] (βOx→TCA→ETC) and we have the overall reaction of
(13) |
where H+(Matrix→IMS) stands for the transport of H+ from the matrix to IMS (Supplementary Text S4). If the reaction pathway βOx→TCA→ETC is coupled with the H+ leak through UCP1 (βOx→TCA→ETC→H+-leak), then we eliminate the transported H+ in Eq (13) (Supplementary Text S4) and we have the overall reaction as
(14) |
Meanwhile, if the reaction pathway βOx→TCA→ETC is coupled with AS (βOx→TCA→ETC→AS), then we also eliminate the transported H+ in Eq (13) (Supplementary Text S4) and it follows that
(15) |
where 8/3 mol H+ and 1 mol H+ were assumed to be transported through ATP synthase and mitochondrial phosphate carrier, respectively, per synthesis of 1 mol ATP [37].
Furthermore, we also devised the reaction pathways for pyruvate. In analogy to palmitate, we considered the reaction pathway for PO, TCA, and ETC (PO→TCA→ETC), the one for PO→TCA→ETC coupled with the H+ leak through UCP1 (PO→TCA→ETC→H+-leak), and the one for PO→TCA→ETC coupled with AS (PO→TCA→ETC→AS). Accordingly, the overall reactions were derived as (Supplementary Text S4), respectively,
(16) |
(17) |
and
(18) |
As we set up the overall reactions, and thereby, the stoichiometries in Eqs (13)–(18), we evaluated their reaction enthalpies ΔrH’º. For this analysis, we used the formula of
(19) |
where ΔfHi’º is a standard transformed formation enthalpy based upon the extended Debye-Hückel theory for a molecular species i at a physiological condition, νi is a stoichiometric number, and the stochiometric number is positive and negative for product and reactant, respectively [38,39]. Note that it is convenient to denote some molecular species in Eq (19) as a species group composed of specific member species (e.g., by CO2 we mean a species group containing specific members of CO2, HCO3–, and CO32–, which are in equilibrium). In addition to the standard transformed reaction enthalpies, we needed to consider the Joule heat due to the transport of H+, ADP, and ATP (Table 1; Supplementary Text S3) by Eqs (4) and (7). For the reaction pathway of βOx→TCA→ETC, we allowed for the transport of 400H+ from the matrix to IMS per 1 mol palmitate (Eq (13)), and its Joule heat was calculated to be –400·VIMM·F=+8490.7 kJ/mol(palmitate), where the factor of –400 includes the direction of H+ transport, F=96,485 C/mol is the Faraday constant, and we used VIMM=–220 mV for the membrane potential of IMM with respect to IMS. Thereby, we added this Joule heat to the reaction enthalpy (Eq (19)) to compute the overall reaction enthalpy including the Joule heat Δr,JH’º(βOx→TCA→ETC) (Table 2). For the enthalpy of the pathway βOx→TCA→ETC→H+-leak including the Joule heat, Δr,JH’º(βOx→TCA→ETC→H+-leak), we added to Δr,JH’º(βOx→TCA→ETC) the exothermic Joule heat due to H+ transported from IMS to the matrix through UCP1 (+400·VIMM·F=–8490.7 kJ/mol(palmitate)) by Eq (7). For the enthalpy of the pathway βOx→TCA→ETC→AS including the Joule heat, Δr,JH’º(βOx→TCA→ETC→AS), we took the sum of Δr,JH’º(βOx→TCA→ETC), ΔrH’º for ATP synthesis, the exothermic Joule heat for the H+ transport from IMS to the matrix through ATP synthase (+400·VIMM·F=–8490.7 kJ/mol (palmitate)), and the exothermic Joule heat for the ATP/ADP antiport in adenine nucleotide translocator (+115.09·VIMM·F=–2443.0 kJ/mol(palmitate); Supplementary Text S3, S4) (Table 2). Similarly, we also performed calculations of Δr,JH’º(PO→TCA→ETC), Δr,JH’º(PO→TCA→ETC→H+-leak), and Δr,JH’º(PO→TCA→ETC→AS) in the same procedure thereof, and we show the results in Table 2. Additionally, as the pH in the matrix is known to be higher than pH 7.0, we calculated Δr,JH’º at pH 7.0, 7.3, and 7.6 as a trial, but the pH changes were rather small (<2%) (Table 2).
Molecular species | ΔfH’º (kJ/mol) | ||
---|---|---|---|
pH 7.0 | pH 7.3 | pH 7.6 | |
Pyruvate | –596.88 | –596.88 | –596.88 |
O2 (aq) | –11.70 | –11.70 | –11.70 |
CO2 (aq) | –692.97 | –692.48 | –692.18 |
H2O | –286.49 | –286.49 | –286.49 |
ATP | –2979.48 | –2979.35 | –2979.28 |
ADP | –1996.92 | –1996.64 | –1996.49 |
Phosphate | –1299.38 | –1299.62 | –1299.10 |
* The data were calculated according to ref [39] (Supplementary Text S5). “aq” means the water-dissolved state.
Δr,JH’º
(pH 7.0/pH 7.3/pH 7.6) |
||
---|---|---|
kJ/mol (palmitate or pyruvate) | kJ/mol (H+) | |
βOx→TCA→ETC† | –5881/–5879/–5878 | –14.70/–14.70/–14.69 |
βOx→TCA→ETC→H+-leak† | –14370/–14370/–14370 | –35.93/–35.92/–35.92 |
βOx→TCA→ETC→AS† | –11530/–11630/–11690 | –28.82/–29.07/–29.22 |
PO→TCA→ETC | –1019/–1018/–1018 | –22.64/–22.63/–22.63 |
PO→TCA→ETC→H+-leak | –1995/–1995/–1995 | –44.34/–44.33/–44.32 |
PO→TCA→ETC→AS | –1680/–1690/–1697 | –37.33/–37.56/–37.72 |
* The values were calculated based upon the data in Table 1 and methods in Supplementary Text S3–S5. The Joule heat for the electron transfer chain, UCP1, adenine nucleotide translocator, and ATP synthase (Figure 4) were calculated at a membrane potential of –220 mV (with respect to IMS).
† The values were calculated with a substrate of palmitate. In the calculations, we used the formation enthalpy of liquid palmitic acid (ΔfHº=–848.4 kJ/mol) [43].
We also show the standard transformed enthalpies of reaction on the basis of mole of H+ in Table 2. By mole of H+, we meant the amounts of H+ introduced into UCP1 or ATP synthase. Accordingly, we derived an enthalpy value in kJ/mol (H+) from a value in kJ/mol (palmitate) by dividing with a factor of 400 (Eq (13)) or from a value in kJ/mol (pyruvate) by dividing with a factor of 45 (Eq (16)) (Table 2; Figure 4B, C). This H+-based enthalpy (Table 2) can be used to compare it to the previous calorimetry data. If the value of Δr,JH’º(βOx→TCA→ETC→H+-leak)=–35.92 kJ/mol (H+) (pH 7.3, Table 2) is applied to the previously-reported result of the H+ leak measurement [34] (see above), then the amount of heat production at VIMM=–220 mV and ΔCIMM=1.63 nmol min–1 mg–1 mV–1 (see above) was estimated as
(20) |
which is also in the same order as the heat measured by calorimetry reported in ref [35].
In the resting state of mitochondria, pyruvate is thought to be the predominant substrate for the metabolism, whereas in the UCP1-activated state, the fatty acid is thought to be the predominant substrate [2]. In the present enthalpy analysis, we assumed that the H+ leak from IMS to the matrix through UCP1 occurred in the presence of VIMM=–220 mV. The amount of this Joule heating associated with UCP1 should correspond to Δr,JH’º(βOx→TCA→ETC→H+-leak)–Δr,JH’º(βOx→TCA→ETC) (Table 2) or, in a more straightforward manner, VIMM·F, which is –21.2 kJ/mol (H+). This Joule heat amounts to as large as ~60% of Δr,JH’º(βOx→TCA→ETC→H+-leak), suggesting that the H+ leak should be a highly heat producing process in comparison to the overall heat. As for the ATP synthesis pathway (PO→TCA→ETC→AS), the enthalpy change corresponding to the ATP synthesis step is Δr,JH’º(PO→TCA→ETC→AS)–Δr,JH’º(PO→TCA→ETC)=–14.9 kJ/mol (H+), which is still exothermic. Although it is known that a reaction of ADP+Pi → ATP+H2O is endothermic (ΔrH’º=+31.11 kJ/mol (ATP) [39]), the H+ transport from IMS to the matrix through the ATP synthase should undergo exothermic Joule heating alongside the ATP synthesis, leading to the negativity of the reaction enthalpy in total. Moreover, although the values of Δr,JH’º(βOx→TCA→ETC→H+-leak) for the UCP1-activated state and Δr,JH’º(PO→TCA→ETC→AS) for the resting state may be similar to each other (Figure 4B, C), we should also note that the H+ conductance increased by 2-fold with the activation of UCP1 [34], whereas the increase of intracellular ATP concentration in the resting state can give rise to down regulation of the metabolism [36]. Accordingly, we suggest that the interplay of the Joule heating associated with UCP1 and its knock-on acceleration of the metabolism such as βOx, TCA, and ETC most likely plays a key role in the thermogenesis of BATM.
In this study, we theoretically investigated the Joule heating of ion currents involving channel proteins and transporter proteins in biological thermogenesis. We derived the Joule’s law for the ion currents through an ion transport protein from electrochemical potentials and pointed out that both heat production and absorption could occur depending on the direction of the ion current. In the simulation study of the temperature distribution around a single channel protein, we applied the Joule’s law to the ion current occurring in the channel protein to numerically solve the heat diffusion equation, although the temperature increase with the Joule heating was as small as <mK. We then moved onto the model systems in which ensembles of channel proteins were included. We used the H-H model to investigate the Joule heating in a nerve and demonstrated that the Joule heating due to the ion currents across the membrane and along the nerve longitudinal axis could reproduce the traits of experimentally-observed heat in nerve during the action potential. Furthermore, we extended the present theory of Joule heating to electrochemical potential-driven transporters and primary active transporters to investigate the non-shivering thermogenesis associated with the H+ leak through UCP1 in BATM. We applied the Joule’s law to the transport of H+, ATP, and ADP across IMM and showed that the amount of Joule heat of this situation was comparable to the calorimetry data in the literature. Our enthalpy analysis suggests that the interplay of the accelerated rate of the metabolism induced by the H+ leak through the activated UCP1 and the Joule heating of the H+ leak would be the key for the non-shivering thermogenesis of BATM.
Although we have investigated only three cases of Joule heating in biology, ion currents through ion transport proteins in the presence of membrane potential are so ubiquitous that many other biological processes are likely to be associated with Joule heating. For example, the transmembrane potential across the sarcoplasmic reticulum (SR) membrane in skeletal muscle was reported to vary between 0–30 mV (inside positive) as measured from fragmented SR [40]. Previously, Curtin and Woledge reported that a significant part of the sum of the work and heat during muscle contraction was not explained by the energy of the ATP splitting and creatine kinase reactions, denoted as “unexplained heat” [41]. The unexplained heat is currently thought to be attributed to the binding of Ca2+ to troponin C for an exothermic binding enthalpy of –(25–32) kJ/mol [42]. However, the amount of the Joule heat for the Ca2+ transport across the SR membrane from the inside to the outside at a membrane potential of 30 mV would be 6 kJ/mol (Ca2+) (Eq (4)), and thus, the Joule heat, which has been overlooked, should also contribute to the heat of muscle contraction to a significant degree.
We suggested that both heat production and absorption could occur depending on the direction of ion current, and our results showed that this reversible nature was compatible with the previous experimental observations which we have cited. Unlike solid-state materials of semiconductors and metals, concentrations of charge carriers, i.e., intracellular ions, dynamically change as manipulated by ion transport proteins including channels and transporters so that the concentration gradient can sometimes surpass the electric potential gradient. Because of this feature, the reversibility of Joule heating may be more evident in cells than in solid-state electricity conducting materials.
The authors declare no conflict of interest.
TN conceived the idea and coordinated the work. TW designed the models, performed the theoretical calculations, and analyzed the data. TW wrote the manuscript, and TN and TW revised the manuscript.
The evidence data generated and programs coded in this study are available from the corresponding authors upon reasonable request.
This work was partly supported by a grant from CREST, JST (JPMJCR15N3 to TN), and grant-in-aid from Japan Society for the Promotion of Science (23115003, 18H03987 and 18H05410 to TN; 19K05226 and 22K04891 to TW). The authors thank Dr. Kazunori Sugiura at SANKEN, Osaka University for his valuable discussion.