Operational stability of biosensors means “retention of activity of a protein or enzyme when in use” (Gibson, 1999).

It corresponds mainly with the same notion for dynamic systems from stability theory. Traditionally, enzyme–substrate interaction is simulated with the help of the Michaelis–Menten model, which is a nonlinear dynamic system. On the other hand, previous studies have paid little attention to the qualitative behavior of the model from the viewpoint of its stability. Partially, operational stability deals with the asymptotic nature of the stability notion, whereas biochemists conduct experiments during “finite time.”

In turn, the stability theory of dynamic systems, so-called Lyapunov stability, offers powerful tools for enhancing biochemical reactions modeling and saving the qualitative behavior of the systems. Moreover, each model is based on a series of assumptions about biochemical interactions, which allows us to check the validity of phenomena that cannot be verified experimentally, for example, mass action law with the distributed or discrete delays, which will be presented later.

Lactate characteristics: Lactate is an anion of lactic acid and is the final metabolite of the anaerobic breakdown of glucose. It is formed from pyruvate during the processes of glycolysis in the absence of oxygen (Figure 1) (Rabinowitz and Enerbäck, 2020) and is an important substance used in medicine as a marker of hypoxia and a number of other disorders, including diabetes and liver and kidney disorders (Pundir et al., 2016).

Lactate is also actively used in the food industry as an emulsifier, thickener, and acidity regulator. Various salts of lactic acid in the international classification of food additives are numbered E325, E326, E327, E328, and E329. Magnesium lactate is sometimes classified as an antioxidant (Standards, 2021). In addition, lactate can act as an indicator of the activity of bacteria during the fermentation process (Rawoof et al., 2020). It can indicate the freshness and quality of some products—wine (Gamella et al., 2010), juices (Trifirò et al., 1997), etc.

The concentration of lactate in products increases as they spoil due to changes in organoleptic properties. Due to the wide use of lactate as a marker of many processes in medicine and industry, accurate and effective methods of its diagnosis are necessary.

Lactate measurement methods: Traditionally, lactate is measured by colorimetry (Suman et al., 2005; Pundir et al., 2016), spectrophotometry (Pundir et al., 2016), fluorometry (Xue and Yeung, 1994), liquid chromatography (Biagi et al., 2012), and nuclear magnetic resonance (Hu et al., 2012; Park et al., 2015). These methods are certainly effective but have several disadvantages, such as high cost, length, and complexity of the analysis. Some of these methods also require complex preliminary sample preparation.

Biosensors offer effective and sensitive detection methods that can be used in medical institutions to measure the level of lactate in blood, sweat, and other biological fluids. In addition, new methods of lactate detection can be useful for enterprises that must control the processes of manufacturing food products and pharmaceuticals (Rattu et al., 2020).

Classification of lactate biosensors from the viewpoint of a bioselective element: A number of biosensor developments aimed at determining the concentration of lactate are known to be in development. Known biosensors can be divided according to the type of bioselective element: sensors based on lactate oxidase (LOx) and sensors based on lactate dehydrogenase (LHD) (Rassaei et al., 2013). In both cases, the substrate and product of the enzymatic reaction are lactate and pyruvate, respectively.

Basic reactions for lactate biosensors: There is a significant difference in these two reactions. For the LHD reaction, NAD+ is needed as a proton carrier in the reaction of dehydrogenation of lactate to pyruvate (Chatterjee et al., 2023). The lactate oxidase reaction can be much easier because, in it, the role of the proton acceptor is played by oxygen; that is, the use of this enzyme in biosensors does not require additional reagents (Lockridge et al., 1972).

Other variants of biosensors are based on a mixture of LHD and LOx. This configuration of the bioselective element, according to the results obtained by Chaubey et al. (2000), allows measuring lactate at lower concentrations than mono-enzyme biosensors, but this makes the analysis more expensive and is characterized by the complexity of manufacturing.

Classification of lactate biosensors from the viewpoint of the transducer: The main biosensor measurement methods used for lactate analysis are optical—electrochemiluminescence or fluorescence (Pundir and Narwal, 2017) and electrochemical—amperometric, potentiometric, or, less often, conductometric or impedimetric (Rathee et al., 2016). Moreover, various methods of improving the main characteristics of sensors have been applied to known biosensor systems—nanoparticles (Nesakumar et al., 2013; Azzouzi et al., 2015), other nanomaterials (Cui et al., 2007), complex, multi-stage methods of enzyme immobilization (Parra et al., 2006), and multi-enzyme membranes.

Although amperometric (Romero et al., 2010) and potentiometric (Mengarda et al., 2019) biosensors are usually monoenzymatic, that is, based on LHD or Lx, conductometric (Nguyen-Boisse et al., 2013) and impedimetric (Chan et al., 2017) biosensors usually use two-enzyme bioselective elements, for example, based on a mixture of lactate oxidase and peroxidase (Nguyen-Boisse et al., 2013) or a mixture of lactate dehydrogenase and pyruvate oxidase (Chan et al., 2017).

Motivation for lactate biosensor design: Therefore, in connection with the variety of works on the development of biosensor systems for measuring lactate and the prospect of its wide use in various areas of human life, we have concluded that the development of new biosensor methods for determining lactate is an urgent need.

Most existing biosensors are currently not ready for wide implementation and commercialization due to various limitations, such as insufficient sensitivity, selectivity, stability with respect to possible inferents, or too narrow a range of biosensor operation.

According to the model, an enzyme E combines with a substrate S to form an enzyme–substrate complex ES, characterized by a rate constant k ₁. The resulting complex can dissociate into E and S (with a rate constant of k ₋₁) or transform into a product P with a rate constant of k ₂ (Berg et al., 2002).

The speed of the enzyme process is dependent on the ease of formation of the complex of the enzyme with the substrate. For low substrate concentrations, the reaction rate is proportional to the substrate concentration, while at higher concentrations, it tends toward a maximum value and becomes independent of substrate concentration. The general dependence of the rate of an enzyme reaction on substrate concentration is described by an equation called the Michaelis–Menten equation (Eq. 1): v = V max n S K M + n S . (1) Here, v is the rate of reaction (mol/s), V _max is the maximum reaction rate (mol/s), n _S is substrate concentration (mol/dm³), and K _M is the Michaelis–Menten constant (mol/dm³).

The Michaelis–Menten constant from Eq. 1 is an enzyme-specific quantity, dependent on substrate, temperature, and pH and independent of enzyme concentration. This constant is a measure of the affinity of the enzyme for the substrate. The lower the value of the constant, the higher the affinity of the enzyme for the substrate (Berg et al., 2002; Radomska, 2016).

Models of enzyme kinetics are based on compartmental systems, which are dynamic systems characterized by a network of interconnected nodes, each representing a reservoir or compartment where resources are stored (Blanchini et al., 2023). The system’s behavior is governed by the movement of resources, depicted as flows traveling along the edges connecting these compartments. Compartment-based dynamic systems serve as invaluable models across various disciplines, including physiologically based pharmacokinetics (Martsenyuk et al., 2012; Thompson and Beard, 2012), mathematical epidemiology (Reyné et al., 2022; Martsenyuk et al., 2021a; b), enzyme kinetics (Keener and Sneyd, 2009; Craciun et al., 2020; Martsenyuk et al., 2022), demography (Navarro Valencia et al., 2023), and ecology (Krishna et al., 2024). Stability analysis plays a pivotal role in understanding the behavior of these systems under different conditions (Blanchini et al., 2023). In recent years, researchers have made significant strides in advancing stability research methodologies, particularly in the context of compartmental systems with delay (Martsenyuk et al., 2013; Martsenyuk and Gandzyuk, 2013).

One common approach to stability analysis is linearization, which involves approximating nonlinear systems around equilibrium points. This technique has been extensively utilized in physiologically based pharmacokinetic models to assess the stability of drug distribution processes within the body (Martsenyuk et al., 2012). Lyapunov functions represent another powerful tool for stability analysis, offering a rigorous mathematical framework to prove the stability properties of compartment-based systems. In mathematical epidemiology, Lyapunov functions have been employed to establish global stability of disease-free and endemic equilibria in compartmental models of infectious diseases (Martsenyuk et al., 2021a).

Section 2 describes the materials, including the experiment and models used. The methods presented are related to parameter identification and stability research. Section 3 shows the results concerning the parameter identification for models using gamma-distributed delay and with two discrete delays. The qualitative analysis includes the existence and positiveness of the solutions, equilibrium states, marginal stability, and numerical research with the help of Poincaré sections. In Section 4, we discuss the results obtained and open problems.

The objective of the work is to offer the flowchart of the lactate biosensor design, including modeling, parameter identification, and stability analysis.

The enzyme lactate oxidase obtained from Aerococcus viridans with an activity of 100 units (Sigma, United States) was used to create the biosensor. Bovine serum albumin (fraction V) (BSA) and a 25% aqueous solution of glutaraldehyde (GA), glycerol, and KCl were obtained from Sigma-Aldrich (United States). Stock 500 mM sodium L-lactate solution (Sigma-Aldrich, Switzerland) was used as a substrate. HEPES obtained from Sigma-Aldrich (United States) was used to prepare the buffer solution. Other inorganic compounds used in the work were of domestic production and had a purity level of “h.p.” and “p.d.a.”

The general scheme of amperometric transducers is shown in Figure 2. Platinum disc electrodes were used as amperometric transducers, manufactured in the laboratory of the Department of Biomolecular Electronics of the Institute of Biomolecular Biology and Geosciences using the following technology: a piece of platinum electrode with a diameter of 0.5 mm and a length of 3 mm was placed in a capillary tube with an outer diameter of 3.5 mm, and then the narrowed end of the capillary tube was sealed in a torch flame. The electrical connection between the platinum and the silver wire conductor was made by low-temperature soldering using Wood’s alloy. The open end of the capillary was filled with epoxy resin, with part of the conductor inside the capillary and part outside. A copper contact was soldered to the conductor to connect the transducer to the measuring unit.

Amperometric measurements setup is shown in Figure 3. The PalmSens potentiostat (Palm Instruments BV, the Netherlands) was connected to an auxiliary platinum electrode, a silver chloride (Ag/AgCl) reference electrode, and working electrodes based on platinum disc electrodes. The potentiostat was connected to an 8-channel device (CH-8 multiplexer, Palm Instruments BV, the Netherlands), which allowed it to receive signals simultaneously from several working electrodes or biosensors (up to eight simultaneously). The distance between the auxiliary platinum electrode and all working biosensors during the measurement was the same (approximately 5 mm).

Preparation of bioselective membranes: Bioselective membranes were prepared by immobilizing proteins on the surface of a platinum disc electrode by covalent cross-linking of the enzyme in a bovine serum albumin matrix using glutaraldehyde as a cross-linking agent. The enzyme gel containing 5% lactotoxidase, 5% BSA, and 10% glycerol was mixed with a 1% glutaraldehyde solution in a 1:1 ratio, after which the mixture was applied to the sensitive area of the platinum disc electrode, and the membrane was air-dried for 20 min. After immobilization, the residual glutaraldehyde and unbound membrane components were washed off the membranes in a buffer solution for 5 min with constant stirring, changing the buffer several times.

The measurement was carried out in an open measuring cell with a volume of 2 mL under constant stirring at room temperature. A fixed potential of +0.6 V was applied to the electrodes relative to the chloride silver reference electrode. The working buffer was 25 mM HEPES with pH 7.4. The required concentration of the substrate in the cell was set by adding aliquots of the standard solution, 500 mM sodium L-lactate, to the buffer.

The duration of one response (from the addition of the substrate to the signal output to the baseline) was approximately 4–5 min; between responses, the substrate was washed off the biosensor for 5 min, changing the buffer in the measuring cell several times. All measurements were performed in at least three replicates.

The parameter identification technique mentioned in Section 2.5.1 was used for Model (2). The training data correspond to a set of time series of currents corresponding to given initial substrate concentrations n _S (0) equal to 0.1 mM, 0.5 mM, 1.0 mM, and 2.5 mM sequentially.

The goal is to estimate the parameters Π IR 1 CB = k d , a , m , τ min , Λ m 0 , K ∈ R + 6 . The initial values for the estimation were chosen as Π _init (see Table 1). The lower and upper bounds for Π_IR1CB are shown in Table 1 as Π _lower and Π _upper , respectively.

The most valuable information obtained from Model (2) is the density of the distributed delay distribution (see Figure 4), which will be used in further estimations.

The comparison of predicted and expected values of the currents for the optimal set of the parameters Π _opt is shown in Figure 5. We see that the parameter values, being optimal concerning the cost criterion (8), enable us to find the solution of Model (2) closest to the expected currents for the initial substrate concentration of 1.0 mM. For the smaller initial values of the substrate, we have predicted values smaller than the expected ones, whereas for those bigger than 1.0 mM, the predicted values are larger than the experimental ones. The explanation of such an effect lies in the special kind of stability of the model known as marginal stability, which will be evidenced further.

Model (4) requires the estimating the parameters Π IR 1 C 2 DMM = k 1 , k 2 , k − 1 , τ 1 , τ 2 , Λ m 0 , K ∈ R + 7 . They were obtained as a result of the solution of the optimization problem shown in Eqs 8–10. The training data described in Section 3.1 were used.

The initial values for the estimation were chosen as Π _init (see Table 2). The lower and upper bounds for Π_IR1C2DMM are shown in Table 2 as Π _lower and Π _upper , respectively.

The comparison of predicted and expected values of the currents for the optimal set of the parameters Π _opt for Model (4) is shown in Figure 6. We see similar tendencies of expected and predicted values as we did for Model (2), with gamma-distributed delays.