Identifying the Oscillatory Mechanism of the Glucose Oxidase

Subscriber access provided by University of Sussex Library

Article

Identifying the Oscillatory Mechanism of the Glucose Oxidase-Catalase Coupled Enzyme System František Muzika, Radovan Jurašek, Lenka Schreiberova, Vuk Radojkovic, and Igor Schreiber J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b08564 • Publication Date (Web): 12 Sep 2017 Downloaded from http://pubs.acs.org on September 23, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Identifying the Oscillatory Mechanism of the Glucose Oxidase-Catalase Coupled Enzyme System František Muzika, Radovan Jurašek, Lenka Schreiberová, Vuk Radojković, Igor Schreiber* Department of Chemical Engineering, University of Chemistry and Technology, Prague, Technická 5, 166 28 Praha 6, Czech Republic

ABSTRACT: We provide experimental evidence of periodic and aperiodic oscillations in an enzymatic system of glucose oxidase – catalase in a continuous-flow stirred reactor coupled by membrane with a continuous-flow reservoir supplied with hydrogen peroxide. To describe such dynamics we formulate a detailed mechanism based on partial results in the literature. Finally, we introduce a novel method for estimation of unknown kinetic parameters. The method is based on matching experimental data at an oscillatory instability with stoichiometric constraints of the mechanism formulated by applying stability theory of reaction networks. This approach has been used to estimate rate coefficients in the catalase part of the mechanism. Remarkably, model simulations show good agreement with the observed oscillatory dynamics, including apparently chaotic intermittent behavior. Our method can be applied to any reaction system with an experimentally observable dynamical instability.

ACS Paragon Plus Environment

1


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 29

INTRODUCTION The enzymes glucose oxidase (GOX) and catalase (CAT) are naturally present in many living organisms, where they play quite different roles. Catalase is ubiquitous in oxygen exposed organisms across all kingdoms where it has a decisive role in preventing oxidative damage of cells by hydrogen peroxide, which it rapidly converts to water and oxygen. On the other hand, glucose-oxidase naturally occurs mostly in fungi and insects, where it serves to produce hydrogen peroxide, locally used as an antibacterial agent. The enzymes can be obtained from Aspergillus niger and various other resources. Apart from many uses of each enzyme separately, by sharing oxygen and hydrogen peroxide among their reactants and products, they are attractive as a coupled pair for various purposes. For instance, the GOX-CAT system can be used for food preservation,1,2 for suppression of growth of undesired organisms,3,4 or it can induce hypoxia5 and influence cancer cells6. It has been shown that the catalase alone7 and later the coupled GOX-CAT system display oscillations8 in a batch membrane reactor, representing in-vitro oscillatory enzyme reactions in addition to other systems such as glycolysis9, peroxidase-oxidase reaction10 and the KaiABC circadian system11. Unlike the well-established reaction mechanism of glucose-oxidase12,13 the mechanism of catalase is known to a limited extent14–19. In this work, we present three mutually supporting aspects of the GOX-CAT system: (i) the experimentally observed oscillatory behavior in a flow reactor-reservoir system, (ii) a mechanism constructed from pieces found in literature14–19 and (iii) a new method for estimation of unknown kinetic parameters based on stoichiometry of the constructed mechanism and our experimental data. The resulting model allows us to simulate major dynamical features of the


2

Page 3 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


experiments rather reliably, including complex intermittent dynamics, which we attribute to deterministic chaos.

EXPERIMENTAL SECTION Chemicals. glucose oxidase (EC 1.1.3.4) type X-S (activity of used charges: 163 400 and 128 200 units/g solid) from Aspergillus niger (Sigma-Aldrich, UK), catalase (EC 1.11.1.6) aqueous suspension (activity of used charges: 42 600 and 35 100 units/mg protein) from bovine liver (Sigma-Aldrich USA), phosphate buffer monobasic (Sigma-Aldrich, USA), phosphate buffer dibasic (Sigma-Aldrich Germany), D-glucose monohydrate (Lach-Ner, Czech Republic), hydrogen peroxide 30% (Penta Czech Republic), deionized water with conductivity < 1 µScm−1 . Glucose oxidase, catalase and glucose are dissolved in phosphate buffer solution (100 mM) set to pH = 5.8.

Experimental setup. Dynamic behavior was observed by measuring the concentration of dissolved oxygen (DO) in the reaction solution during individual runs of 45 000 seconds with sampling period of 1 second. The calibration was set according to current temperature and atmospheric pressure, further corrections were done automatically by DO-meter. The DO-probe is capable of reading up to 50 mg/l of DO with precision ±1.5 mg/l. The control measurements of pH of stock solutions were performed using HI 5222 pH meter. The combined reactor (volume 11.95 cm3) and reservoir (volume 1.85 cm3) system was designed and manufactured locally. We used a peristaltic pump ISMATEC IPC-8 with TYGON R3607 tubes of internal diameter 0.25 mm. Operating temperature was controlled by thermostated air in the laboratory at 23 °C.


3


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 29

Experiments were performed according to scheme shown in Figure 1. Hydrogen peroxide is pumped by the peristaltic pump (8) into the reservoir (1), where it diffuses through the membrane (5) to the reactor (2). GOX, CAT and glucose stock solutions are pumped directly into the reactor stirred by stirring bar (6) set to 1020 rpm on the magnetic stirrer (11), while the concentration of dissolved oxygen is measured by DO probe (3). Outflowing solutions from the reservoir and reactor are transported by pressure provided by the pump to the waste tank (9).

Figure 1. Experimental setup: 1 – reservoir, 2 – reactor , 3 – DO probe, 4 – thermometer probe, 5 – Spectra/Por Dialysis Membrane MWCO: 1000, 6 – stirring bar, 7 – stock solutions, 8 – peristaltic pump, 9 – waste tank, 10 – HI98196 DO/pH/ORP meter, 11 - stirrer.

Experimental results. We performed a series of experiments using a flow reactorreservoir system operating under conditions partly similar as in Sasaki et al.8 who used a batch system open to air. Our setup combines a continuous stirred tank reactor with a flow-through reservoir coupled via membrane. The enzymes and glucose are introduced into the reactor, whereas hydrogen peroxide enters the reservoir. Oxygen dissolved in stock solutions enters through all inlets at saturated condition. The system is closed to air. Examples of the measured DO time series at specific input parameters and flow rate are shown in Figure 2. The dynamics


4

Page 5 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


can be either irregular (panel A), which we find to be a more frequent case, or regular periodic oscillations (panel B). The interspike periods range from tens to hundreds of seconds. Series of experiments with successively varying inlet concentration of catalase indicate that the oscillations emerge suddenly when the enzyme concentration exceeds certain value. We assume that this feature indicates presence of a subcritical Hopf bifurcation in the corresponding model, which is exploited in our method of estimating kinetic parameters. The wave form of the oscillations is more resembling the earlier measurements with catalase7 (sharp increase, gradual decrease) rather than the GOX-CAT system, but this may be attributed to differing experimental conditions, most importantly, our setup is closed to air which caused significant supersaturation of oxygen (typically 300-500%).

Figure 2. Experimental time series of dissolved oxygen in the reactor. Reactor inlet concentrations: c 0(GOX ) = 1.3 × 10 −6 M , c 0( glu cos e ) = 0.667 × 10 −3 M , (A) c0 (CAT ) = 0.264 × 10 −7 M ,


5


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 29

(B) c0(CAT ) = 0.33 × 10 −7 M ; reservoir inlet concentration: c 0,res ( H 2O2 ) = 0.574 M ; reactor flow rate: k 0 = 9.73 × 10 −5 s −1 , reservoir flow rate: k0 , res = 2.09 × 10 −4 s −1 .

THEORETICAL SECTION Mechanism. In order to set up a model for simulations of the measured oscillatory dynamics, below we propose a comprehensive mechanism of GOX-CAT reaction occurring in a CSTR, which exchanges substrates/products with a flow-through reservoir. The mechanism is presented as a reaction network in Figure 3. The GOX reaction cycle12,13 starts with the native (oxidized) GOX reacting with β-D-glucose to form an enzyme-substrate complex (step 1), which afterwards decomposes to D-glucono-δ-lactone and a reduced form of GOX (step 2). Step 2 may be affected by glucose at higher concentrations,13 but this effect can be neglected in our case. The D-glucono-δ-lactone hydrolyses to gluconic acid, which does not feed back due to buffered pH. The reduced form of GOX reacts with oxygen, a substrate produced by CAT, forming another enzyme-substrate complex, which decomposes to hydrogen peroxide and the oxidized form of GOX (step 3), returning it to the initial point (step 4) of the GOX cycle. Next, we construct a CAT mechanism as a combination of three pathways. The CAT cycle starts with oxidizing the heme group of the native CAT by hydrogen peroxide to form a highvalence iron intermediate called compound 1 (CPD1) (step 5). CPD1 is assumed further to undergo reduction by hydrogen peroxide (step 6) or by superoxide radical (step 9), initiating two pathways. The first and primary pathway begins by the reaction of CPD1 with hydrogen peroxide, forming a CPD1-H2O2 complex17. According to recent theoretical considerations18 the CPD1-H2O2 complex does not directly split to native CAT and oxygen; rather, it reversibly turns


6

Page 7 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


into another complex called CDP2-like (step 7), which then splits to oxygen and CAT (step 8) thus closing the loop. As a side step, CPD1 undergoes reversible deactivation (step 13)19. The second pathway starts by forming compound 2 (CPD2) and oxygen by reduction of CPD1 with superoxide radical.14,15 CPD2 is further oxidized by hydrogen peroxide to form compound 3 (CPD3) (step 10),17 which reversibly splits to CAT and superoxide radical (step 11)14,15, thereby closing the second loop. Additionally, there is a third possible loop14, in which CPD2 is reduced by superoxide radical to native CAT, while forming another molecule of oxygen and water (step 12). The reaction occurs in a reactor-reservoir system mutually coupled via a membrane to allow for a slow transport of hydrogen peroxide from the reservoir to the reactor along with a slow transport of oxygen and superoxide radical, which is assumed present in a small amount in the H2O2 pool in the reservoir.


7


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 29

Figure 3. Stoichiometric network of the GOX-CAT system in a membrane-coupled CSTRreservoir system. Green arrows – GOX subnetwork; orange arrows – CAT subnetwork; black arrows –inflows, outflows and membrane fluxes.

Method for estimation of kinetic parameters. Once the model is formulated in terms of the stoichiometric network, we need to estimate unknown kinetic parameters necessary for occurrence of a subcritical Hopf bifurcation. This concerns primarily the CAT subnetwork, whereas parameters for the GOX part are known11,12. To that end we use the stoichiometric network theory20 and introduce a novel method of matching experimental results at the condition of emergence of oscillations with the stoichiometric constraints imposed by the assumed mechanism. The entire reactor-reservoir system can be written as dx = Nv (x ), dt

(1)

where N is the (n × m) stoichiometric matrix, x = ( x1 ,L, xn ) is the vector of concentrations of n dynamical species (inert products and buffered/pooled species are not included), v = (v1 ,L , vm ) is the vector of m reaction rates assumed to follow mass action kinetics. N and v accommodate all the inflows, outflows and membrane transport terms treated as pseudoreactions of zeroth or first order. For the species occurring in both the reactor and the reservoir there are separate components in x. Also, different volumes of both subsystems are reflected in N. The rate vector

v is expressed as v = (k1v1 ,L ,k m vm ) = diag (k ) v , where the vector k includes the rate coefficients and v (x) is the vector of reduced rates. The set of all rate vectors v s = v(x s ) at a stationary state xs satisfies Nv s = 0 and thus v s is restricted to the non-negative part of the null


8

Page 9 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


space of N called the flux cone. According to the stoichiometric network theory20 the set of all rate vectors v s can be expressed as vs = Eα, α = (α1 , L , α f ) , α k ≥ 0 ,

(2)

where the columns of E are the (normalized) edges of the flux cone. The edges represent irreducible subnetworks called extreme currents or elementary subnetworks/fluxes. Thus the steady state network v s is obtained as a non-negative linear combination of elementary subnetworks where α k represents a contribution of the k-th subnetwork. For the purpose of parameter determination the rate coefficient vector k and the steady state concentration vector x (for convenience we drop the subscript denoting a steady state) can be split to three parts with: (i) fixed values, known or estimated, (ii) unknown values to be determined, and (iii) indeterminate values implied by solving for the unknowns,

(

)

(

)

k = k ( fv ) , k (uv ) , k (iv ) , x = x ( fv ) , x (uv ) , x (iv ) .

(3)

[

]

(

)

Moreover, E and α can be split into two parts, E = E (uds ) ,E (uv ) and α = α (uds ) ,α (uv ) , where α (uv ) corresponds to E (uv ) and includes unknown coefficients α k to be determined and α (uds ) ,

E(uds ) represent a preselected unstable dominant subnetwork (uds). The selection of uds is based on stability analysis of the elementary subnetworks20,21 and selecting one (or combination of a few) that represent the decisive instability leading to the emergent dynamics observable in experiments, most notably to an oscillatory instability via a Hopf bifurcation. Then α (uds ) is playing the role of a bifurcation parameter that upon variation lets the system to pass through the bifurcation point. We can now formulate the constraint equations for a set of unknowns by


9


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 29

selecting certain equations from eq 2 given that we can fix: (i) the entire reaction rate (by fixing both the rate coefficient and the relevant concentrations), (ii)

reduced reaction rate (by

measuring or estimating the relevant concentrations), (iii) expression derived from the reaction rate which contains a concentration as a linear unknown (by fixing the rate coefficient and all concentrations except that which is a linear unknown). As a result, we obtain the constraint equations,  E ˆ (uv )  

0 A 0

0  α (uv )   v (fv )      0  k (uv )  =  0  − Eˆ (uds )α (uds ) , B   x (uv )   0 

(4)

where v (fv ) are the fixed rates based on k (fv ) and x (fv ) according to case (i), and matrices A and B are constructed according to cases (ii) and (iii) respectively. The symbol ̂ indicates that only rows of E relevant to the selected reactions are retained in Eˆ . Eqs 4 are expected to have more unknowns than equations due to a typically large number of elementary subnetworks. Thus there is a continuum of solutions and we need a clue how to select a suitable one. A reasonable choice is suggested by the relation between the unstable dominant subnetwork and all the other elementary subnetworks. To ensure that the uds at the Hopf bifurcation is indeed dominant, we may require that the sum of the contributions of the remaining subnetworks is minimal, i.e.,

(

)

!

f α (uv ) , x (uv ) , k (uv ) = ∑ α k(uv ) = min .

(5)

This leads to a linear optimization problem where the unknowns are assumed non-negative and subject to constraints in eqs 4 with the objective function in eq 5. Once the unknown values


10

Page 11 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


α (uv ) , k (uv ) and x (uv ) are determined, v s can be found from eq 2 and subsequently k (iv ) and x (iv )

evaluated.

APPLICATION TO THE GOX-CAT SYSTEM Formulation of the model and parameter determination. The outlined theory assumes an input from experimental measurements and subsequent use of stoichiometric constraints to estimate unknown rate coefficients. More specifically, we use the values of input parameters corresponding to the experimentally observed sudden appearance of oscillations summarized in Table 1. These include the values of flow rates for the reactor and reservoir and inflow concentrations of the enzymes, hydrogen peroxide and oxygen (assumed saturated). Additionally, we assume a small inflow concentration of superoxide radical carried with the hydrogen peroxide. The model is formulated in terms of mass balances for each dynamical chemical species conforming with the format of eq 1. Details are provided in the Supporting information, below we summarize the main features. There are 18 dynamical species, 16 mass action reactions with independent rate coefficients and 3 membrane transport coefficients (permeances) for hydrogen peroxide, oxygen and superoxide radicals (the leak of other species from reactor to reservoir is neglected). Upon decomposing the network into elementary subnetworks and stability analysis we identified an unstable dominant subnetwork involving steps 5, 6, 8, 9, 10 and forward direction of 7 and 11. This uds is complemented by inflow of hydrogen peroxide to the reservoir and its transport across the membrane, as well as outflow of oxygen from the reactor. The fixed values k (fv ) and

x (fv ) needed in the constraint equations (4) are as follows: k (fv ) includes input parameters found


11


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 29

in Table 1 and the rate coefficients of the GOX part taken from literature.13 We also fixed the permeances based on our measurement of the permeance for hydrogen peroxide. The fixed concentrations in x (fv ) include the measured concentration of oxygen and concentrations of certain species fixed at values that are sufficiently small as required by the reaction network theory to observe an instability. Stability analysis20,21 of the unstable dominant subnetwork identifies those species as hydrogen peroxide in the reactor and the following forms of catalase: (uds ) CAT, CPD1, CDP2 and CDP3. We then solve eqs 4,5 successively at stepwise increasing α

which destabilizes the system through a Hopf bifurcation, and we thus obtain an estimate of the rate coefficients of the CAT part. A more detailed description including parameter values used to compare dynamics of the model with the experiments is found in the Supporting information.

Table 1. Experimental conditions at the onset of oscillations. Parameter value Parameter -5 -1 reservoir flow rate, k 0,res reactor flow rate, k 0 9.73 × 10 s

Value 2.09 × 10 -4 s -1

c0( glu cos e )

0.667 × 10 −3 M

c0 ,res ( H 2O2 )

0.574 M

c0 ( O 2 )

0.266 × 10 −3 M

c0,res ( O 2)

0.266 × 10 −3 M

c0(CAT )

0.166 × 10 −7 M

c 0 ,res (O2 − )

9 × 10 -5 M

c0(GOX )

1.3 × 10 −6 M

Calculations and comparison with experiments. The resulting set of parameters (Table S1 in Supporting information) is used in eq 1 to analyze the dynamics of the model by performing direct simulations as well as using numerical continuation and bifurcation analysis. 22 The results are summarized in Figure 4 as one-parameter bifurcation diagram overlaid with our experimental observations of periodic and aperiodic dynamics.


12

Page 13 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 4. Overlay of experimental results and one-parameter bifurcation diagram with dissolved oxygen plotted against inlet concentration of catalase. Black curve – stationary solutions, red curve – minima and maxima of periodic solutions, full curve – stable, dashed curve – unstable, square – Hopf bifurcation point (HB), triangle – period doubling point (PD). External parameters as in Figure 2. The calculated DO is plotted against the inlet concentration of the catalase for steady state and for the primary periodic solution that bifurcates from the steady state via a subcritical Hopf bifurcation at the left. The periodic solution undergoes a period-doubling bifurcation giving rise to a cascade of periodic orbits with successively doubled period eventually leading to chaotic dynamics. The chaotic region has the familiar pattern with interspersed periodic windows and disappears as the parameter is increased via an inverse period-doubling sequence. Finally, the stabilized primary periodic oscillations shrink to a steady state via another Hopf bifurcation. Major regions of occurrence of periodic and chaotic oscillations found in the model are marked


13


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 29

by shading. Overlaid with the model predictions, the location of our experimental results and their dynamical type are indicated by vertical lines. The experimental results fit rather well within the calculated region of oscillations, the model predicts sudden appearance of oscillations as well as alternation of aperiodic and periodic oscillations. The experimentally observed aperiodic regime is emerging immediately after the oscillatory instability, whereas the model predicts a rather wide range of the inlet catalase concentration where periodic oscillations occur before chaos sets in via period-doubling sequence. Also, as the parameter is further increased, chaos disappears earlier that the experimentally observed aperiodicity. An example of calculated DO oscillations are shown in Figure 5. The inlet catalase concentration was chosen to show a typical chaotic (panel A) and periodic (panel B) patterns. The waveform is in close correspondence with the experimental time series (Figure 2) but the characteristic time is about 25 times larger in the model. This is because in applying the equations 4,5 we did not focus on matching times scales. However, it can be fixed through further systematic refinement of the raw estimate of the kinetic parameters as discussed in the next section. Such refinement and the suggested deterministic nature of the aperiodic experimental time series will be examined in our future work.


14

Page 15 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure

5.

Calculated

DO

in

the

reactor. (A)

c0 (CAT ) = 0.259 × 10 −7 M ,

(B)

c0 (CAT ) = 0.222 × 10 −7 M ; other parameter values as in Table 1 and Table S1 (see Supporting information).

DISCUSSION Even though the agreement is semi-quantitative, the application of the outlined method to the coupled GOX-CAT enzymatic system results in a consistent estimation of DO levels, type of observed dynamics, oscillatory amplitudes and waveform, and apart from time scale provides reasonable correspondence between the model and experimental observations. If the parameters are further refined so that the model is validated by a broad range of experiments, prediction of dynamic behavior of the GOX/CAT pair may prove useful in various applications23-25. However, it needs to be stressed that agreement presented in this work does not imply that the used enzyme mechanisms are without deficiencies. In particular, the source of superoxide radical or an alternative electron donor has been widely discussed in literature.14-17 We simply assume its


15


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 29

presence in the stock solution of hydrogen peroxide, which we qualitatively confirmed by UV spectrophotometry, but there may be other sources, e.g. due to ferric salts that could be present as impurities in the buffer17. On the other hand, our analysis clearly shows, that the presence of the second and/or third pathways involving an electron donor other than hydrogen peroxide is essential for the oscillatory instability (the chosen unstable dominant subnetwork includes the second pathway, in fact, any unstable subnetwork involves either the second or third pathway). Nonetheless, and more importantly, the method of estimating kinetic parameters based on stoichiometry and stability constraints is found to be effective in providing not only qualitative insight into capability of the model to yield a dynamical instability, which has been the traditional realm of the reaction network theories,26,27 but establishes a systematic approach to an estimate of missing kinetic parameters to quantitatively predict that instability. Such an estimate may be initially rough, but provides a starting point for refinements. Also, the more detailed is the input information, the more qualified is the estimate. Hynne et al.28,29 have developed a "global optimization" method based on systematic scanning of the entire reaction rate cone in combination with the concentration space. In principle, there is a unique point of best fit with experimental data, but scanning the entire available space is practically limited to small-sized networks. In addition, any inequalities that would narrow the scanned space must be defined ad hoc26. In contrast, the method proposed here takes any restrictive conditions into account by virtue of formulating the problem as a constrained system. In another work, Hynne et al.30 have used the reaction rate cone to refine certain kinetic parameters in a detailed model of glycolytic oscillations based on their quenching experiments. However, scanning for a Hopf bifurcation was tailored to special features of the glycolytic mechanism rather than based on a general guiding principle as represented by eqs 4,5. It should be mentioned that there are comprehensive


16

Page 17 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


program packages that include an analysis of stoichiometric networks albeit to a limited extent. In particular COPASI31 for analysis of biochemical dynamical systems performs, in addition to many other functions, basic stoichiometric network analysis (decomposition into elementary subnetworks, metabolic flux analysis). The presented method is principally different from standard methods of parameter estimation based on a least-square-fit in that it does not provide a global optimum. Rather, the solution of eq 4,5 at the point of emerging instability matches exactly the available rate coefficients, flows and input/output concentrations with the model but the resulting solution at the point of instability does not necessarily match additional properties, such as the period of oscillations. To refine the raw result, further steps are needed. In our case, we simply modified the results by trial-and-error to obtain agreement with periodic/aperiodic sequence of oscillations observed experimentally and did not focus on matching the time scale. In general, the obtained solution of eq 4,5 at the oscillatory instability can be taken as an initial point for systematic refinement by performing least-square optimization or other suitable procedure such as bifurcation analysis, either directly within the system of eq 4,5 or upon transferring the raw result to eq 1. Our preliminary bifurcation analysis in the parameter plane

of inlet catalase

concentration vs reactor flow rate indicates, that there is a large and parametrically robust region of oscillatory dynamics delimited by Hopf bifurcation in the neighborhood of the parameter set found in Table S1 in Supporting information and a substantial part of that region is corresponding to chaotic dynamics emerging via period-doubling cascade. Moreover, the characteristic time of oscillations can be tuned by rescaling only catalase part and membrane permeances, which is reflecting the dominance of the unstable subnetwork identified in this work that constitutes the core of the oscillator.


17


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 29

CONCLUSIONS The outlined method exploiting stoichiometry and stability holds promise for analysis of other reaction systems with oscillatory dynamics, where mechanisms are available but kinetic data are scarce and difficult to obtain by established methods of reaction kinetics. If alternative mechanisms are available, the method can help in distinguishing among them and revealing their potential for making reliable predictions. Such systems occur in heterogeneous catalysis with important industrial applications, in vitro enzyme reactions as treated here, and oscillatory reactions occurring in living organisms. In particular, applications to biological oscillators are of interest. In this regard, models of biological oscillating systems would need to be formulated in terms of power law kinetics based on (pseudo)elementary steps rather than using rational polynomial functions as is currently prevailing. The extension to mass action type kinetics implies possibly large networks, which, however, does not preclude applications of the type of analysis presented here. ASSOCIATED CONTENT

Supporting information. Details of the reaction network model for the GOX-CAT system and description of the procedure toward estimating the kinetic parameters. AUTHOR INFORMATION

Corresponding Author [email protected]


18

Page 19 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ACKNOWLEDGMENT This work was supported by the project Kontakt II LH14009 from the Ministry of Education, Youth and Sport of the Czech Republic and by the grant 15-17367S from the Czech Science Foundation.

REFERENCES (1) Parpinello, G.P.; Chinnici F.; Versari, A.; Riponi, C. Preliminary study on glucose oxidase–catalase enzyme system to control the browning of apple and pear purees. Lebensm.-Wiss. u.-Technol. 2002, 35, 239–243. (2) Chaudhari, A.; Nitin, N. Role of oxygen scavengers in limiting oxygen permeation into emulsions and improving stability of encapsulated retinol. J. Food Eng. 2015, 157, 7–13. (3) Massa, S.; Petruccioli, M.; Brocchi, G. F.; Altieri, C.; Sinigaglia, M.; Spano, G. Growth inhibition by glucose oxidase system of enterotoxic escherichia coli and salmonela derby: in vitro studies. World J. of Microbiol. Biotechnol. 2001, 17, 287-291. (4) Sobotta, M. C.; Barata, A. G.; Schmidt, U.; Mueller, S.; Millonig, G.; Dick T. P. Exposing cells to H2O2: A quantitative comparison between continuous low-dose and one-time high-dose treatments. Free Radic. Biol. Med. 2013, 60, 325-335. (5) Mueller, S.; Millonig, G.; Waite, G. N. The GOX/CAT system: a novel enzymatic method to independently control hydrogen peroxide and hypoxia in cell culture. Adv. Med. Sci. 2009, 54, 121-135. (6) Askoxylakis, V.; Millonig, G; Wirkner, U.; Schwager, Ch.; Rana, S.; Altmann, A.; Haberkorn, U.; Debus, J.; Mueller, S.; Huber, P.E. Investigation of tumor hypoxia using a two-enzyme system for in vitro generation of oxygen deficiency. Radiat. Oncol. 2011, 6, 1-12.


19


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 29

(7) Hideshima, T.; Inoue, T. Nonlinear oscillatory reaction of catalase induced by gradual entry of substrate. Biophys. Chem., 1997, 63, 81-86. (8) Sasaki, S.; Uchida, T.; Nakano, D.; Hideshima, T. Glucose measurement using an enzymatic oscillatory reaction. Electroanalysis 2004, 16, 1598-1602. (9) Chance, F.; Schoener, B., Elsaesser, S.; Control of the waveform of oscillations of the reduce pyridine nucletide level in a cell-free extract. Proc. Natl. Acad. Sci. U.S.A. 1964, 52, 337-341. (10) Olsen, L. F.; Degn, H.; Oscillatory kinetics of the peroxidase-oxidase reaction in an open system. Experimental and theoretical studies. Biochim. Biophys. Acta. 1978, 523, 321334. (11) Tomita, J.; Nakajima, M.; Kondo, T.; Iwasaki, H. No transcription-translation feedback in circadian rhythm of kaic phosphorylation. Science 2005, 307, 251-254. (12) Gibson, Q. H.; Swoboda, B. E.P.; Massey, V. Kinetics and mechanism of action of glucose oxidase. J. Biol. Chem. 1964, 239, 3927-3934. (13) Duke, F. R.; Weibel, M.; Page, D. S.; Bulgrin, V. G.; Luthy, J. The glucose oxidase mechanism. Enzyme activation by substrate. J. Am. Chem. Soc. 1969, 91, 3904-3909. (14) Kono, Y.; Fridovich, I. Superoxide radical inhibits catalase. J. Biol. Chem., 1982, 257, 5751-5754. (15) Shimizu, N.; Kobayashi, K.; Hayashi, K. The reaction of superoxide radical with catalase – mechanism of the inhibition of catalase by superoxide radical. J. Biol. Chem.

1984, 259, 4414-4418. (16) Davison, A.J.; Kettle, A.J.; Fatur, D.J. Mechanism of the inhibition of catalase by ascorbate – roles of active oxygen species, copper and semidehydroascorbate. J. Biol. Chem. 1986, 261, 1193-1200.


20

Page 21 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(17) Lardinois, O. M.; Mestdagh, M. M.; Rouxghet, P. E. Reversible inhibition and irreversible inactivation of catalase in presence of hydrogen peroxide. Biochim. Biophys. Acta 1996, 1295, 222-238. (18) Alfonso-Prieto, M.; Biarnés, X.; Vidossich, P.; Rovira, C. The molecular mechanism of the catalase reaction. J. Am. Chem. Soc. 2009, 131, 11751–11761. (19) Alfonso-Prieto, M.; Biarnés, X. ; Rovira, C. The reaction mechanisms of heme catalases: an atomistic view by ab initio molecular dynamics. Arch. Biochem. Biophys. 2012, 525, 121–130. (20) Clarke, B. L. Stability of complex reaction networks. Adv. Chem. Phys. 1980, 43, 1-216. (21) Ross, J.; Schreiber, I.; Vlad, M. O. Determination of complex reaction mechanisms: analysis of chemical, biological and genetic networks. Oxford University Press: Oxford, UK, 2006. (22) Kohout, M.; Schreiber, I.; Kubicek, M. A computational tool for nonlinear dynamical and bifurcation analysis of chemical engineering problems. Comput. Chem. Eng., 2002, 26, 517-527. (23) Gonzales-Ramos, M.; de Frutos, S.; Griera, M.; Luengo, A.; Olmos, G.; RodriguezPuyol, D.; Calleros, L.; Rodriguez-Puyol, M. Integrin-linked kinase mediates the hydrogen peroxide-dependent transforming growth factor-β1 up-regulation. Free Radic. Biol. Med. 2013, 61, 416-427. (24) Trifonov, A.; Herkendell, K.; Tel-Vered, R.; Yehezkeli, O.; Woerner, M.; Willner, I. Enzyme-capped

relay-functionalized

mesoporous

carbon

nanoparticles:

eﬀective

bioelectrocatalytic matrices for sensing and biofuel cell applications. ACS Nano 2013, 7, 11358-11368. (25) Ma, X.; Jannasch, A.; Albrecht, U.-R.; Hahn, K.; Miguel-Lopez, A.; Schaﬀer, E.; Sanchez, S. Enzyme-powered hollow mesoporous Janus nanomotors. Nano Lett. 2015, 15, 7043−7050.


21


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 29

(26) Craciun, G.; Tang, Y.Z.; Feinberg, M. Understanding bistability in complex enzymedriven reaction networks. Proc. Natl. Acad. Sci. U.S.A., 2006, 103, 8697-8702. (27) Hadac, O.; Muzika, F.; Nevoral, V.; Pribyl, M.; Schreiber, I. Minimal oscillating subnetwork in the Huang-Ferrell model of the MAPK cascade. PLoS ONE, 2017, 12, e0178457. (28) Hynne, F.;

Sorensen, P.G.;

Moller, T.

Complete optimization of models of the

Belousov-Zhabotinsky reaction at a Hopf-bifurcation. J. Chem. Phys., 1993, 98, 219230. (29) Hynne, F.; Sorensen, P.G.; Moller, T. Current and eigenvector analyses of chemicalreaction networks at Hopf bifurcations. J. Chem. Phys., 1993, 98, 211-218. (30) Hynne, R.; Dano, S.; Sorensen, P.G. Full-scale model of glycolysis in Saccharomyces cerevisiae. Biophys. Chem. 2001, 94, 121-163. (31) Hoops, S.; Sahle, S.; Gauges, R.; Lee, C.; Pahle, J.; Simus, N.; Singhal, M.; Xu, L.; Mendes, P.; Kummer, U. COPASI – a COmplex PAthway SImulator. Bioinformatics,

2006, 22, 3067-3074. (32) Press, W.H.; Teukolsky, S.A.; Vetterling, W.T.l; Flannery, B.P. Numerical recipes in Fortran: the art of scientific computing. Cambridge University Press: Cambridge, UK, 1992.

TOC GRAPHIC


22

Page 23 of 29


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 ACS Paragon Plus Environment

23


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

82x28mm (300 x 300 DPI)


Page 24 of 29

Page 25 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


153x87mm (300 x 300 DPI)



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

171x85mm (300 x 300 DPI)


Page 26 of 29

Page 27 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


152x71mm (300 x 300 DPI)



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

154x85mm (300 x 300 DPI)


Page 28 of 29

Page 29 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


51x43mm (300 x 300 DPI)


Identifying the Oscillatory Mechanism of the Glucose Oxidase

Recommend Documents