Aerobic life is sustained by a reaction analogous to that of a hydrogen fuel cell. The reduction of oxygen to water, O

_{2} + 4H

^{+} + 4e

^{−} → 2H

_{2}O, is catalyzed by the protein cytochrome

*c* oxidase (CcO). This reaction generates the electric potential of ≈200 mV across the inner mitochondrial (or bacterial) membrane that powers the production of ATP, the fuel of living cells () [

1–

7]. CcO takes up four protons from the negative (N) side of the membrane and four electrons from the positive (P) side of the membrane to reduce one oxygen molecule to form two water molecules (). Part of the ≈2 eV of chemical energy released during this reaction is used to translocate four protons from the N- to the P-side of the membrane against an membrane potential, resulting in net transport of 8 proton charges across the membrane with a thermodynamic efficiency of ≈8×200 meV/2eV = 80 % (). Unlike molecular motors and transporters that undergo large conformational changes, the proton pump function in CcO is achieved without large-scale changes in protein structure [

8], and its molecular mechanism has remained elusive.

Here, we use a frequency-dependent bias voltage to probe the molecular mechanism of the proton pump. We employ a stochastic-kinetic approach, used widely in studies of molecular machines [

9–

16]. The effects of time-varying electric fields are relevant not only physiologically because of the constantly fluctuating membrane potentials in cells [

17,

18], but also because they provide a unique window into the function of molecular machines [

19–

22]. The frequency dependence of the measurable proton and electron currents [

23] allows us to pinpoint key reaction steps in the pumping function of CcO, and to distinguish between competing models that are both consistent with a large body of experiments.

Our calculations are based on a detailed master-equation description of CcO that is consistent with basic physical principles, is built on the known structure of CcO, and reproduces the rates and equilibria of intermediate reaction steps [

24]. The proton and electron conduction pathways of the real enzyme are represented by three charge sites, two for protons and one for electrons. The proton and electron sites, each being either empty or singly occupied, are kinetically connected to the two sides of the membrane and to each other (). Individual transitions satisfy detailed balance between forward and backward rates, consistent with the second law of thermodynamics. The exergonic reaction of oxygen reduction is described by a product formation step that requires simultaneously occupied electron and proton-1 sites (blue arrow in ). Detailed balance is broken during product formation, with reactant and product concentrations held steady and the backward reaction (i.e., product breakup) slowed by the free energy gain from the chemical reaction. As a result, the system is driven out of equilibrium toward a steady state with constant fluxes of electrons, protons and products.

The dynamics of the populations **P**(*t*) = (*P*_{1}(*t*), *P*_{2}(*t*), ···, *P*_{8}(*t*))^{T} of the 2^{3} = 8 microscopic states is governed by a master equation *d***P**(*t*)/*dt* = **KP**(*t*), where **K** = [*k*_{ij}] is the 8 × 8 rate matrix corresponding to the reaction diagram in being the rate coefficient for transitions from state *j* to *i*, and *k*_{ii} = − Σ_{j}_{≠}_{i} k_{ji}. The forward and backward rates (except for the product formation step) are written as *k*_{ji} = *κ*_{ji} exp [−(*G*_{j} − *G*_{i})/2*k*_{B}*T*] to satisfy detailed balance, *k*_{ji}/*k*_{ij} = exp[−(*G*_{j} − *G*_{i})/*k*_{B}*T*]. *G*_{i} is the free energy of state *i*, and *κ*_{ij} = *κ*_{ji} is the intrinsic rate coefficient in the absence of a driving force (i.e., *G*_{i} = *G*_{j}). For simplicity, we assume that the free energy difference is balanced between the forward and backward reactions, resulting in the factor 1/2 in the exponent. The free energy of state *i* contains 1-body and 2-body contributions,

where

is the intrinsic free energy to occupy site

*μ* with all other sites empty,

is an occupancy indicator equal to 1 (0) if site

*μ* is occupied (empty) in state

*i*,

*ε*_{μν} is the electrostatic coupling between sites

*μ* and

*ν*,

*q*_{μ} = ±

*e* is the charge at site

*μ*,

*z*_{μ} is the distance from the N side of the membrane to site

*μ*,

*L* is the membrane width and

*V*_{m} =

*V*_{P} −

*V*_{N} is the membrane potential. Here we assume a constant electric field

*V*_{m}/

*L* inside the membrane.

Product formation (steps V → 0 and VII → II in ) is driven by a free energy gain of Δ*G*_{p} = 0.5 eV, corresponding to ≈1/4 of the energy released by the formation of two water molecules from one oxygen molecule. In the master equation, we assume that this driving force is realized by multiplying the backward rates of product formation by exp(−Δ*G*_{p}/*k*_{B}*T*).

Previous studies [

24,

25] showed that the three-site kinetic models are the simplest description of CcO that can pump protons across the membrane. These models achieve the stoichiometric efficiency seen in experiments [

1,

7] of

*η* ≈ 1 proton pumped per electron consumed against time-independent, opposing membrane potentials,

*V*_{m} > 0. Here we define the pumping efficiency as

*η* =

*J*_{pump}/

*J*_{el} where

*J*_{pump} and

*J*_{el} are the average fluxes of protons pumped and electrons consumed, respectively. Note that

*η* is related to the thermodynamic efficiency by

*V*_{0}(1 +

*η*)/Δ

*G*_{p}.

Our focus here is to study the effects of time-dependent membrane potentials,

*V*_{m}(

*t*) =

*V*_{0} +

*V*_{1} cos(

*ωt*), on the pumping efficiency of CcO.

*V*_{0} is a constant offset voltage,

*V*_{1} is the amplitude of the oscillatory bias voltage, and

*ω* its frequency. Oscillating voltages

*V*_{m}(

*t*) are applied to models constructed previously [

24] to satisfy experimental data on proton and electron affinities, reaction rates, and equilibria, while pumping protons against constant voltages of

*V*_{0} > 100 mV. Here, we consider two representative models () that differ in their pump mechanism, i.e., in the order of the reaction steps in their dominant pump cycles; states I → V → VI → VII → II → III → I for model 1 (pump cycle 5 in Ref. [

24]; labels as in ), and 0 → I → II → VI → VII → V → 0 for model 2 (pump cycle 2 in Ref. [

24]).

| **TABLE I**Parameter values (in units of *k*_{B}*T* and s^{−1}) of the two three-site pump models [24]. |

Oscillating electric potentials are applied on top of the base voltages of

*V*_{0} = 0,

and

, where

and

are the membrane voltages at which the pumping efficiencies are

*η* = 1/2 and 0, respectively, with

*V*_{1} = 0. The resulting master equation is periodic with period 2

*π*/

*ω*. According to Floquet theory, a quasi steady state is established at long times

*t*, and time averages can be replaced by phase averages [

26]. The flux from state

*j* to

*i* can thus be calculated as

, where

*J*_{ij}(

*t*) =

*k*_{ij}(

*t*)

*P*_{j}(

*t*) −

*k*_{ji}(

*t*)

*P*_{i}(

*t*) is the net flux from state

*j* to

*i* at time

*t*. The global electron and proton-pump fluxes,

*J*_{el} and

*J*_{pump}, are then calculated by summing the contributions of contributing elementary fluxes

*J*_{ij} (with

*J*_{el} being exactly twice the rate of forming product water).

To relate the experimentally measurable frequency dependence of the pumping efficiency and product formation rate to the underlying microscopic processes of proton and electron conduction, we combine first- and second-order perturbation theory with an eigenmode analysis of the master equation. The Taylor expansion of the time-dependent rate matrix

**K**(

*t*) in terms of

*V*_{1} is given by

, where

**K**_{0} is the rate matrix of the time-independent system with

*V*_{1} = 0, and

**K**_{1} and

**K**_{2} are the first and second derivatives of

**K** with respect to

*V*_{m}, evaluated at

*V*_{m} =

*V*_{0}, respectively. The master equation is then solved perturbatively for

**P**(

*t*) [

27]. Up to the leading order in

*V*_{1}, the probabilities and the corresponding average fluxes can be written as sums of Lorentzians [

28]:

where

**P**_{0} is the steady-state probability vector,

is the average flux at

*V*_{m} =

*V*_{0}, and

*λ*_{k} and

are the eigenvalues and corresponding right eigenvectors of the time-independent rate matrix

**K**_{0} ordered such that

*λ*_{1} = 0 >

*λ*_{2} > ··· >

*λ*_{8} (all eigenvalues are real in our models; see Ref.

29 for effects of complex eigenvalues). The constants

*e*_{k},

*δ*_{k}, and

*D*_{ij}_{,}_{k} contain terms involving components of

**K**_{0},

**K**_{1} and the left/right eigenvectors of

**K**_{0} [

27].

shows the efficiency *η* and electron flux *J*_{el} of model 1 calculated from first- and second-order perturbation expansions of **P**(*t*) as a function of the amplitude *V*_{1} of the oscillatory bias voltage at a frequency of 10^{3} s^{−1} and without offset, *V*_{0} = 0. The efficiency decreases monotonically to about 90 % as the amplitude increases to *V*_{1} = 100 mV. Interestingly, however, the electron flux, and thus the rate of product formation, increases with the oscillating voltage. This increase reflects the nonlinear dependence of the electron flux on the membrane potential, with gains in the electron flux under low potential outweighing losses under high potential. Based on the excellent agreement with the results of practically exact numerical integration up to *V*_{1} = 60 mV, we will use second-order perturbation theory for the following calculations, unless stated otherwise.

shows that both the pumping efficiency and the electron flux depend strongly on frequency, even at a small fixed amplitude of

*V*_{1} = 50 mV of the oscillatory bias voltage. The effect of the oscillatory voltage becomes more pronounced as the offset voltage is increased from

*V*_{0} = 0 to

, and

. Remarkably, the two models show distinct frequency dependences at all three offset voltages. At zero offset (

*V*_{0} = 0; top panels), the efficiency in model 2 hardly changes with the frequency, whereas the efficiency drops by about 13 % for model 1 as

*ω* → 0. In contrast, the electron flux (or turnover rate) is insensitive to

*ω* for model 1, whereas it drops by more than 15 % for model 2 at low frequency. At nonzero offset voltages (

and

) the efficiencies show nonmonotonic behavior as a function of

*ω*, whereas the electron fluxes monotonically decrease with increasing

*ω*. In particular, at

model 1 loses the pumping ability for

*ω* 10

^{8} s

^{−1} whereas model 2 pumps protons in the entire frequency range.

For both the efficiency and the electron flux, the two limits of low (

*ω* → 0) and high frequency (

*ω* → ∞) differ from each other and from the value in the case of a time-independent voltage (i.e.,

*V*_{1} = 0). This difference is evident, e.g., in the bottom left panel of , where the efficiency of the pump without oscillatory voltage is

*η* = 0 by construction. In the adiabatic limit,

*ω* → 0, the fluxes can be averaged over one period with time-independent voltages:

, where

is a weight factor arising from the Jacobian associated with changing variables from time

*t* to voltage

*V*. At the other extreme, when

*ω* → ∞, the system dynamics is governed by time-independent effective rate coefficients averaged over a period,

.

The eigenmode analysis can also be used to extract information about the mechanisms of the molecular machine. Perturbation theory,

Eqs. (2) and

(3), allows us to identify the contributions of individual eigenmodes to the average fluxes. We find that in model 1 the coupling to the fifth eigenmode with the eigenvalue −

*λ*_{5} ≈ 10

^{8} s

^{−1} results in the greatest change in the efficiency as a function of the frequency. The corresponding right eigenvector

has two dominant elements at states V and VII of opposite sign. As a consequence, the populations of states V(+0−) and VII(+ + −) oscillate with a ≈180 degree phase shift around their respective steady-state values, according to

Eq. (2). The shift in population from state V to VII is achieved by an internal proton transfer (PT) (V → VI) with a rate coefficient of ≈10

^{8} s

^{−1}, facilitated by the presence of an electron in site 3. This PT is followed by a fast proton uptake (VI → VII, rate ≈10

^{10} s

^{−1}).

In model 2, the second and sixth eigenmodes dominate the frequency-dependence of the efficiency, with |*λ*_{2}| ≈ 10^{5} and |*λ*_{6}| ≈ 10^{7} s^{−1}, respectively. These two eigenmodes shift the population from state I(+00) to VI(0 + −) and VII(+ + −), respectively. These shifts are again achieved by the internal PT (I → II), now in the absence of an electron, but followed by an electron uptake (II → VI) and the protonation of site 1 (VI → VII). Unlike model 1, electron uptake is the fastest reaction with a rate of ≈10^{11} s^{−1}, while the PT and the protonation of site 1 occur on a time scale of ≈10^{6} s^{−1}.

The response of the proton pump currents to oscillating electric fields () is reminiscent of stochastic resonance phenomena observed in many areas of physics and biology, including optical, electronic and magnetic systems, neuronal circuits [

26], and biochemical reaction networks [

29]. Here the output (i.e., the efficiency or electric current) is amplified in the presence of a weak coherent input (i.e., an oscillating voltage) by the assistance of noise inherent in the stochastic systems. Noise in membrane potentials [

17] has been studied for neurons [

18], but little is known for organelles such as mitochondria. Remarkably, oscillatory voltages enhance the pumping efficiency of one of the models essentially over the entire frequency regime, but reduce the efficiency of the other ( left panels). Model 2 thus appears to be better adapted to the fluctuating electric fields in the cellular environment.

The strong effects of oscillating electric fields on a biological proton pump, as found here, are more complex than those in simple bistable systems studied extensively by theory and experiment [

26]. This complex response to oscillating fields reveals details about the microscopic processes of coupled proton and electron transfer events and their contributions to the proton pump function. Different proton pump mechanisms can be identified by the characteristic frequency dependence of their pump and turnover fluxes. Measurements of these fluxes, for instance by using CcO embedded into surface-attached membranes under controlled voltage [

23], will provide important guidance toward a full molecular understanding of the machine that powers all aerobic life. The same formalism used here to characterize a proton pump can be used in studies of other molecular machines, such as molecular motors under oscillatory force load.