Integrated Experimental and Theoretical Studies on an Electrochemical Immunosensor

Abstract

Electrochemical immunosensors (EIs) integrate biorecognition molecules (e.g., antibodies) with redox enzymes (e.g., horseradish peroxidase) to combine the advantages of immunoassays (high sensitivity and selectivity) with those of electrochemical biosensors (quantitative electrical signal). However, the complex network of mass-transfer, catalysis, and electrochemical reaction steps that produce the electrical signal makes the design and optimization of EI systems challenging. This paper presents an integrated experimental and modeling framework to address this challenge. The framework includes (1) a mechanistic mathematical model that describes the rate of key mass-transfer and reaction steps; (2) a statistical-design-of-experiments study to optimize operating conditions and validate the mechanistic model; and (3) a novel dimensional analysis to assess the degree to which individual mass-transfer and reaction steps limit the EI’s signal amplitude and sensitivity. The validated mechanistic model was able to predict the effect of four independent variables (working electrode overpotential, pH, and concentrations of catechol and hydrogen peroxide) on the EI’s signal magnitude. The model was then used to calculate dimensionless groups, including Damkohler numbers, novel current-control coefficients, and sensitivity-control coefficients that indicated the extent to which the individual mass-transfer or reaction steps limited the EI’s signal amplitude and sensitivity.

1. Introduction

Electrochemical biosensors are analytical devices that detect analytes by transforming a biochemical reaction into a quantitative, electrical signal. This class of biosensors has proven valuable in research, quality control, food safety, medical diagnosis, and monitoring of therapeutic efficacy [1]. Miniaturized amperometric biosensors that use redox enzymes to generate an electric current in response to voltage applied at a working electrode have been successfully commercialized; personalized blood glucose meters used by diabetics represented 85% of the total biosensor market in 2008 [2]. By 2013, the worldwide market for glucose-monitoring biosensor systems was estimated to be billions of dollars per year, with screen-printed-electrode (SPE) arrays that served as single-use biosensor “strips” representing two thirds of that market [3]. The disposable, redox-enzyme-based biosensor market is being further expanded via the commercialization of glucose-monitoring systems for animals [4].

Immunoassays based on the exceptionally high binding selectivity and affinity of biological recognition molecules (predominantly antibodies, but also aptamers [5]) have a projected global market that is expected to reach $8 billion in 2022 [6]. Immunoassays typically involve a “sandwich” molecular architecture, in which immobilized capture antibodies first bind target-analyte molecules present in the liquid sample, and then secondary antibodies labeled with reporter molecules that generate an optical signal also bind the analyte molecules. The resulting molecular “sandwiches” consist of an analyte molecule held between primary and secondary antibody molecules. To estimate the analyte concentration, the surface concentration of bound reporter molecules is measured by the intensity of the optical signal they generate, and a calibration curve is used to convert the reporter molecule’s concentration into the analyte concentration [7]. Commonly used reporter molecules for immunoassays include redox enzymes whose products can be measured optically, such as horseradish peroxidase (HRP). HRP offers multiple advantages as a reporter. It is robust, has a relatively small molecular size, is inexpensive, is readily bound to antibodies in an active form, has a high turnover rate, and can oxidize a wide range of substrates to yield optically active products [8,9].

Whereas virtually all commercial immunoassay systems involve optical detection, the benefits of integrating electrochemical biosensors and immunoassays have been recognized [10]. Such hybrid electrochemical immunosensors (EI) have the potential to combine the advantages of immunoassays (extremely high sensitivity and selectivity) with those of electrochemical biosensors (reproducible, quantitative, continuous electrical output). The electrical output is achieved by forming a sandwich molecular architecture on the working electrode, and the reporter molecule triggers an electrical signal. Redox enzymes are commonly used as EI reporters, because some of their reaction products can be either oxidized or reduced at the working electrode, resulting in an electric current that serves as the EI’s output. This approach offers exceptional versatility, because an EI biosensor could be developed for virtually any analyte for which antibodies can be developed. Also, inexpensive, disposable SPE arrays designed to be read by portable meters similar to glucose meters EI could be mass-produced. The resulting EI platform would enable an extremely wide range of molecular and cellular analytes to be accurately measured with high sensitivity and selectivity, ease of use, low cost, and portability [11,12,13,14,15]. Prototype EI systems have been developed for healthcare applications. Sanchez-Tirado et al. fabricated an EI to measure cytokines used as markers of inflammation [16]. Tallapragada et al. developed an EI for human epidermal growth factor receptor 2 (HER2) that had a detection limit of 4 ng/mL [17]. Dempsey et al. described a disposable, printed lateral flow EI for human cardiac troponin T (cTnT) [18]. The reporter used in all of these studies, HRP, generated an oxidized product that was electrochemically reduced at the working electrode, resulting in a continuous amperometric output.

However, commercial implementation of EI systems has been hampered by the complexity of the multiple molecular mass-transfer, binding, and reaction steps that give rise to the electrical signal. This complexity complicates efforts to design new EIs that achieve specified performance metrics, including the lower detection limit and sensitivity (defined as change in output per unit change in analyte concentration). Fabrication methods and operating conditions needed to achieve these metrics are expected to vary between EI systems, due to factors including analyte-antibody binding affinities, the concentrations of primary antibodies bound to the electrode, and the kinetics of both the reporter enzyme’s reaction and the electrochemical reaction. These kinetics will, in turn, be influenced by the liquid sample’s properties, including its pH, its concentrations of the analyte, and substrates for the enzymatic reaction. Moreover, the concentrations of redox-active interferents in the sample may limit the working electrode’ voltage.

Development of robust product-design algorithms for new EI systems that meet specified performance metrics would be aided by mechanistic mathematical models that quantitatively describe the rates of the key molecular mass-transfer, binding, and reaction steps. Such models would enable the step(s) that limits performance to be identified and guide strategies to overcome such limitation(s). To date, few mechanistic models of HRP-based EIs have been reported [19,20,21,22,23], and these models have not been sufficiently comprehensive to predict how the output would vary with key independent variables, including the working electrode’s applied voltage (E), the pH, and the concentrations of HRP’s substrates. Such models are needed to help design EIs, identify factors that limit their performance properties, and guide research strategies to optimize EI systems. Mechanistic models would also help support petitions for U.S. Food and Drug Administration (FDA) approval of EI systems for healthcare applications. The FDA requires that stringent accuracy and consistency standards be met by portable glucose monitoring systems while in the hands of lay users [24], and similar requirements would be expected for EIs. Mechanistic models would enable rapid, in-silico hypothesis testing, including “what-if” studies to assess whether non-standard use by lay users could result in dangerously incorrect readings.

This paper addresses the need for such mechanistic models by presenting a novel, integrated experimental and mathematical framework to characterize EI performance, and then applies the framework to optimize performance of a novel EI that can detect a target protein (mouse IgG) at the ng/mL level. The framework includes three components. The first is a detailed mechanistic model that can predict the rates of the individual mass-transfer and reaction steps that give rise to the EI’s amperometric output. The second is a statistical-design-of-experiments approach that generates an empirical, statistical model describing the effects of key independent variables on the EI’s output. This statistical model is used both to optimize the EI system and to help validate the mechanistic model. The third is an integration of dimensional analysis with principles of flux-control theory to quantify the extent to which individual mass-transfer and reaction steps limit the EI’s sensitivity and output current (J). The paper concludes by discussing the utility of the integrated experimental and mathematical framework for future design, optimization, and validation of EI systems.

2. Materials and Methods

2.1. Materials

Thioctic acid, sodium phosphate (monobasic and dibasic), mouse IgG, anti-mouse IgG antibody (ap124), HRP-conjugated-goat anti-mouse IgG (a5278), TWEEN 20, H₂O₂), C, and N-hydroxysulosuccinimide sodium salt (NHS) were obtained from Sigma Aldrich (St. Louis, MO, USA). MES buffered saline packs and 1-ethyl-3-(3-dimethylaminopropyl carbodiimide hydrochloride) (EDC) were purchased from ThermoFisher Scientific (Waltham, MA, USA). Ultrapure water (18.2 MΩ) was produced by a Nanopure-UV four-stage purifier (Barnstead International, Dubuque, IA, USA); the purifier was equipped with a U.V. source and a final 0.2 μm filter. Ultrapure water was used in all aqueous solutions. Screen-printed electrodes were obtained from Conductive Technologies Inc. (York, PA, USA) and Metrohm DropSens (Oviedo, Asturias, Spain) (models DRP-250BT and DRP-110SWCNT).

2.2. Preparation of Immunosensing Layer

The immunosensing layer was prepared by using 1-ethyl-3-(3-dimethylaminopropyl (EDC) and N-hydroxysulfosuccinimide sodium salt (NHS) chemistry to attach the primary (capture) antibodies covalently to carboxylate groups present on the DropSens array’s working electrodes. EDC-NHS chemistry has been widely used to fabricate the immunosensing layers of EIs [25,26,27]. EDC is a zero-length cross-linker that activates carboxylate groups for covalent coupling to primary amines. The addition of NHS with EDC results in an NHS ester intermediate that reacts rapidly with primary amines, thereby increasing the efficiency of the coupling reaction [28]. Cleaned gold SPEs were dipped in 15 mM thioctic acid in ethanol for 1 h. The resulting carboxylated SPEs were washed with ethanol and dried under nitrogen. The carboxyl groups were activated by incubating the SPEs in 100 mM MES buffer containing 5.0 mM EDC and 9.0 mM NHS at pH 4.6 for 1 h at room temperature. Electrodes were then rinsed with MES buffer and dipped in 6 µg/mL goat anti-mouse IgG antibody in 50 mM phosphate buffer at pH 7 for 2 h. The primary-antibody-functionalized SPEs were then washed with phosphate buffer. To block nonspecific binding of the secondary antibody, the SPEs were incubated in 2 % bovine serum albumin (BSA) in phosphate buffer for 1 h at room temperature. The resulting functional SPEs were washed with phosphate buffer at pH 7 and stored in phosphate buffer at 4 °C.

SPEs were each dipped in a standard solution having a known concentration of the target analyte (mouse IgG) in a 2% aqueous BSA solution in 50 mM phosphate buffer at pH 7 for 1 h at room temperature. The SPEs were then washed four times with washing buffer (0.05% TWEEN20 in 50 mM phosphate buffer at pH 7) and incubated in a [1:333] dilution of HRP-conjugated-goat anti-mouse IgG in pH 7, 50 mM phosphate buffer in 2% BSA (Figure 1). After 1 h, the electrodes were rinsed four times with washing buffer and stored in phosphate buffer at 4 °C until the electrochemical measurements were conducted.

2.3. Electrochemical Measurement of EI Signal

The EIs were removed from the refrigerator and allowed to equilibrate at room temperature. Forty µL of a solution (subsequently referred to as the “bulk solution”) containing 50 mM phosphate buffer, 1 mM $H_2{O}_2$, and 8 mM C were added on the SPE. Wire leads from a potentiometer (CHI 660, C.H. Instruments, Austin, TX, USA) were connected to the EI’s working, reference, and auxiliary electrodes, and a reduction potential of −0.2 V relative to an Ag/AgCl reference electrode was applied to the working electrode. After about 1 min, the reduction current (i.e., the EI’s signal (J)) reached a steady-state value, and the current level was recorded as the EI’s output for that set of experimental conditions. Each EI was used once. All electrochemical potentials given in this paper are relative to an Ag/AgCl reference electrode.

2.4. Optimization of EI Operating Conditions and Characterization of EI Performance Properties

A statistical design of experiment (DOE) approach was used for two purposes: (1) to determine the values of key independent variables that optimized the EI’s signal and (2) to obtain an empirical equation that described the effects of the key independent variables on the EI’s signal to help validate the mechanistic model. The independent variables expected to most strongly affect the performance of the EI described above included (1) the working electrode’s E, (2) the bulk solution’s [C], (3) the bulk solution’s [H₂O₂], and (4) the bulk solution’s pH [29].

A two-level half-factorial design with center points and three replicates for each experiment was set up using Minitab^® software (Supplementary Table S1). For each factor, the following three levels, denoted low (−1), center point (0), and high (+1), were chosen: −0.05 V, −0.125 V, and −0.2 V for E; 1.0 mM, 4.5 mM, and 9.0 mM for [C]; 0.5 mM, 1 mM, and 1.5 mM for $\left[H_2{O}_2\right]$; and 6.2, 6.6, and 7.0 for pH, respectively. To avoid electrical noise arising from reduction of redox-active interferents in the bulk solution [30], the lowest E value was set to −0.2 V. To control the rate of C autoxidation [31], 8 mM was selected as the highest [C] value. Experiments were conducted in triplicate for each combination of factors specified by Minitab^® using a constant analyte concentration of 40 ng/mL mouse IgG. Each EI’s signal was calculated as the difference between the J measured first in the absence of analyte and then in the presence of the analyte. All signal data were input to Minitab^®, which provided a statistical analysis of the results. The experimental conditions that Minitab^® indicated were optimal for the EI were used in subsequent experiments to characterize the EI’s performance properties. In these experiments, the EI signal was measured in triplicate for six concentrations of the analyte.

3. Mechanistic Mathematical Model

The mechanistic mathematical model of the EI describes the transport and reaction processes involving catechol ($C$), O-quinone $\left(Q\right)$, and hydrogen peroxide $(H_2{O}_2)$ that generate a current (J) at the EI’s working electrode. Differential mass-balance equations describe the diffusion of these species in the x direction (perpendicular to the electrode), through two layers (Figure 2) that lie between the electrode’s surface at x = 0 and the bulk solution: (1) the immunosensing layer between x = 0 and x = L containing the antibodies and HRP, and (2) a stagnant, aqueous, diffusion layer between x = L and x = L + δ. The HRP-catalyzed conversion of $C$ and $H_2{O}_2$ to $Q$ is assumed to occur uniformly throughout the immunosensing layer, and the electrochemical reduction of Q to C is assumed to occur on the electrode’s surface. The bulk solution is assumed be well-mixed, with the concentrations of all chemical species remaining constant at their initial values [32]. Mass transfer is assumed to follow Fick’s law, with a diffusion coefficient (D) that is assumed to be the same for $Q$,$C$, and $H_2{O}_2$ but to vary between the diffusion layer ($D_{\delta })$ and the immunosensing layer ($D_L$). The HRP concentration and maximum reaction rate constant ($V_{max}$) are assumed to be uniform throughout the immunosensing layer [33].

3.1. Kinetics of Enzymatic and Electrochemical Reactions

The non-linear, ping-pong kinetic mechanism describing HRP oxidation of C in the presence of $H_2{O}_2$ is shown in Reactions (i)–(iii) [33,34,35,36]:

HRP (Fe³⁺) + H₂O₂ → Compound(I) + H₂O

(i)

Compound (I) + C → Compound (II) + Q

(ii)

Compound (II) + C → HRP (Fe³⁺) + Q

(iii)

where compounds (I) and (II) are oxidized intermediates of HRP. The kinetic formula resulting from this mechanism [21,37,38,39,40] is:

$$V=\frac{V_{max}\left[H_2{O}_2\right]\left[C\right]}{K_m^C\left[H_2{O}_2\right]+K_m^{H_2{O}_2}\left[C\right]+\left[H_2{O}_2\right]\left[C\right]}$$

(1)

where $V$ is the reaction rate, $V_{max}$ is the maximum reaction rate constant ($V_{max}=k^{cat }\left[HRP\right])$, $k^{cat }$and [HRP] are turnover number and HRP concentration within the immunosensing layer, respectively; $K_m^C$ and $K_m^{H_2{O}_2}$ are the corresponding Michaelis–Menten constants, and $\left[H_2{O}_2\right]$and $\left[C\right]$ are $H_2{O}_2$ and C concentrations, respectively.

Molecules of Q produced by HRP can be reduced back to C at the surface of the working electrode in a two-electron, two proton reaction shown in Reaction (iv) at a rate described by the Butler–Volmer Equation (2) [41]:

Q + 2e⁻ + 2H⁺ → C

(iv)

$$J=nF{D}_L{\left[\frac{dQ}{dx}\right]}_{x=0}= nF{K}_0{\left[Q\right]}_{x=0} e^{\left(-\frac{\alpha nF\left(E-E_h\right)}{RT}\right)}-nF{K}_0{\left[C\right]}_{x=0} e^{\left(\frac{\left(1-\alpha \right)nF\left(E-E_h\right)}{RT}\right)}$$

(2)

where, $J$ is the electric current density, $n$ is the number of transferred electrons (n = 2 for this reaction), $\alpha$ is the charge transfer coefficient (assumed 0.35), $F$ is the Faraday constant (96,485 C mol⁻¹), $K_0$ is the apparent electron transfer rate constant for Q, $R$ is the universal gas constant (8.314 J K⁻¹ mol⁻¹), $T$ is the absolute temperature (298 K), and $E_h$ is the redox potential for electrochemical reduction of Q to C under the experimental conditions used (0.15 V at pH 6.2). Values of $E_h$ for a given set of experimental conditions were determined as the midpoint potential ($E_{mid})$ between the cathodic peak (for Q reduction) and anodic peak (for C oxidation) of cyclic voltammograms obtained under the same conditions [42]. The calculated value of J was taken to be the current generated by the EI.

The effect of pH on $E_{mid}$ is shown in Equation (3) [43,44], in which m (=2) and n (=2) are the number of protons and electrons involved in the reduction of Q, respectively. This equation indicates that increasing the pH would make $E_{mid}$ more negative and thereby reduce the working electrode’s overpotential, reaction rate, and EI’s signal, according to the Butler–Volmer equation. To simulate the effect of pH on$E_h$, Equation (3) was incorporated in the mechanistic model.

$$E_{mid} ~ \mathrm{const}-2.303\frac{\mathrm{mRT}}{\mathrm{nF}}pH$$

(3)

3.2. Mass Balance Equations

Assuming one-dimensional diffusion in the x-direction, the steady-state, differential, mass balance equations including diffusion and enzymatic reaction for$H_2{O}_2$,$C,$ and Q across the immunosensing layer (0 < x < L) are shown in Equations (4)–(6) [33,45,46,47]:

$$0=D_L\frac{d^2\left[H_2{O}_2\right]}{d{x}^2}-\frac{V_{max}\left[H_2{O}_2\right]\left[C\right]}{K_m^C\left[H_2{O}_2\right]+K_m^{H_2{O}_2}\left[C\right]+\left[H_2{O}_2\right]\left[C\right]}$$

(4)

$$0=D_L\frac{d^2\left[C\right]}{d{x}^2}-\frac{V_{max}\left[H_2{O}_2\right]\left[C\right]}{K_m^C\left[H_2{O}_2\right]+K_m^{H_2{O}_2}\left[C\right]+\left[H_2{O}_2\right]\left[C\right]}$$

(5)

$$0=D_L\frac{d^2\left[Q\right]}{d{x}^2}+\frac{V_{max}\left[H_2{O}_2\right]\left[C\right]}{K_m^C\left[H_2{O}_2\right]+K_m^{H_2{O}_2}\left[C\right]+\left[H_2{O}_2\right]\left[C\right]}$$

(6)

3.3. Boundary Conditions

Previous mathematical models [47,48,49,50] describing electrochemical reduction of Q have assumed the electrochemical driving force (E-E_h) was sufficiently large that [Q] at the electrode’s surface (where x = 0) could be assumed to be approximately zero Equation (7).

$${\left[Q\right]}_{x=0}=0$$

(7)

However, this assumption is likely to be invalid for an EI under some realistic operating conditions. For example, to avoid electrical noise and/or interference by electroactive species in the solution, it may be desirable to use a moderate (E-E_h) value, for which ${\left[Q\right]}_{x=0}$ would not be vanishingly small, and use of Equation (7) would cause significant error in the model’s predictions. For that reason, we used the Butler-Volmer equation (Equation (2)) as a boundary condition at the working electrode surface. This equation is valid over the entire spectrum of positive and negative (E-E_h) values.

Because Q reduction at the electrode generates C in equimolar amounts, the fluxes of Q and C at x = 0 were assumed to be equal in magnitude but opposite in sign (Equation (8)). Control experiments showed that J caused by the reduction of $H_2{O}_2$was close to zero under the experimental conditions (Figure S1, available in supplementary information). Therefore, the flux of $H_2{O}_2$ at x = 0 was assumed to be zero Equation (9).

$$J=nF{D}_L{[\frac{dQ}{dx}]}_{x=0}=-nF{D}_L{[\frac{dC}{dx}]}_{x=0}$$

(8)

$${[\frac{d{H}_2{O}_2}{dx}]}_{x=0}=0$$

(9)

Partitioning kinetics of all reactants between the diffusion layer and the immunosensing layer were assumed to be rapid enough that the interfacial concentrations were assumed to remain at equilibrium. Identical partition coefficients ($k_p$) were assumed for all reacting species Equations (10)–(12).

$${\left[Q\right]}_{x=L-}=k_p{\left[Q\right]}_{x=L+}$$

(10)

$${\left[H_2{O}_2\right]}_{x=L-}=k_p{\left[H_2{O}_2\right]}_{x=L+}$$

(11)

$${\left[C\right]}_{x=L-}=k_p{\left[C\right]}_{x=L+}$$

(12)

The bulk solution (at x = ∞) was assumed to be well mixed and have the concentrations indicated in Equations (13A)–(13C).

$${\left[C\right]}_{x=\infty }=C(\infty )$$

(13A)

$${\left[H_2{O}_2\right]}_{x=\infty }=H_2{O}_2(\infty )$$

(13B)

$${\left[Q\right]}_{x=\infty }=0$$

(13C)

No reaction is assumed to occur in the diffusion layer, so the mass transfer rate of $C, H_2{O}_2$ and Q across this layer is modeled as the product of a mass transfer coefficient (D/δ) and the concentration driving force across the layer. Also, at the interface between the diffusion layer and the immunosensing layer, the diffusive fluxes of $C, H_2{O}_2$ and Q exiting one layer are assumed to be equal to those entering the other layer Equations (14)–(16).

$$D_L{[\frac{\partial Q}{\partial x}]}_{x=L-}=-\frac{D_{\delta }}{ \delta }\left\{{\left[Q\right]}_{x=L+}-0\right\}$$

(14)

$$D_L{[\frac{d{H}_2{O}_2}{dx}]}_{x=L-}=\frac{D_{\delta }}{\delta k_p }\left\{k_p H_2{O}_2(\infty )-{\left[H_2{O}_2\right]}_{x=L-}\right\}$$

(15)

$$D_L{[\frac{dC}{dx}]}_{x=L-}=\frac{D_{\delta }}{\delta k_p }\left\{k_p C(\infty )-{\left[C\right]}_{x=L-}\right\}$$

(16)

The coupled, second-order differential Equations (4)–(6) that described nonlinear kinetics of HRP-catalyzed oxidation of C to Q Equation (1) and electrochemical reduction of C back to Q Equation (2), along with the boundary conditions Equations (8)–(16) were solved numerically using function BVP4C in MATLAB (codes available in supplementary information).

4. Results and Discussion

4.1. EI System’s Properties under Optimal Operating Conditions

Based on the half-factorial experiments with a centerpoint, the experimental conditions that optimized the EI signal were E = −0.2 V, [C] = 8 mM, pH = 6.2, and $\left[H_2{O}_2\right]$ = 1 mM. The subsequent EI characterization experiments, which were conducted under these optimal experimental conditions (Figure 3), indicated that the EI’s limit of detection was 1 ng/mL, its sensitivity was 0.63 nA mL/(ng mm²), and its inter-assay/intra-assay variation was less than 5%.

4.2. Validation of Mechanistic Model

Minitab’s^® statistical analysis of the experimental optimization studies was integrated with the mechanistic mathematical model of EI operation for three purposes: (1) to help validate the mechanistic model, (2) to explain trends seen in the experimental data, and (3) to develop new methods to identify factors that limit an EI system’s signal strength and sensitivity to the target analyte.

Some of the constants used in the mechanistic model (

Table 1. Values of constants used in the mechanistic model. The $K_0$ and [HRP] values were fit to the experimental data obtained using a constant analyte concentration of 40 ng/mL.

Parameter	Value	Units	Source
$k^{cat}$	2.2 × 104	s−1	[46]
[HRP]	0.5	µM	-
$K_m^C$	3.0	mM	[46]
$K_m^{H_2{O}_2}$	0.2	mM	[46]
$K_0$	8.0 × 10−7	cm s−1	-
$E_h$	0.15	V	-
δ	3.0 × 10−3	cm	[55]
L	25	nm	-
$D_{\delta }$	2.2 × 10−5	cm2 s−1	-
$D_L$	2.3 × 10−6	cm2 s−1	-
$k_p$	1.0	-	-

) were obtained from literature data. Others were estimated by fitting the model to the empirical, statistical model that Minitab^® generated from the experimental optimization studies. The statistical model was a best-fit polynomial that expressed the EI’s signal as a function of the four factors. The polynomial had a linear term for each factor and binary, ternary, and quaternary product terms for each combination of factors to simulate interactions between factors (Equation (S1), available in supplementary information).

Values for the kinetic constants of HRP’s kinetic model were obtained from the BRENDA database [51]. The diffusion layer (δ) for the unstirred bulk solution was assumed to remain constant [52] at a value of 200 µm [53,54]. The thickness of immunosensing layer was assumed to be 25 nm [55]. Values of diffusion coefficients in the immunosensing layer and diffusion layer were assumed to be 2.5 × 10⁻⁶ cm² s⁻¹ and 2.2 × 10⁻⁵ cm² s⁻¹, respectively [56]. Values of $K_0$ and [HRP] (5.0 × 10⁻⁷ cm s⁻¹ and 0.5 µM, respectively) were fit to the experimental data obtained using a constant analyte concentration of 40 ng/mL. Because the deposition of HRP molecules in the immunosensing layer results from formation of sandwich molecular architectures, the HRP value is expected to vary with the analyte concentration in the bulk liquid.

To help validate the mechanistic model, trends in the model’s prediction of how each of the four independent variables influenced the EI’s signal were compared to the corresponding experimental data (Figures 5–8). The strength of each independent variable’s effect was quantified as the standardized effect (SE) value [42] in the Pareto chart (Figure 4) generated by Minitab^®. The dotted line marks the minimum SE value for statistical significance at the 95% confidence level (SE = 2.09). These results indicate that all four independent variables significantly affected the signal, with the strength of those effects decreasing in the order E (SE = 11.4) > [C] (SE = 8.9) > pH (SE = 4.6) > [$H_2{O}_2$] (SE = 2.1).

The strong increase in the EI’s signal with E, and thus the magnitude of (E-E_h), is apparent in both the experimental results and the model’s predictions (Figure 5). This trend is attributed to the Butler–Volmer Equation’s (2) exponential dependency of the EI’s amperometric signal on (E-E_h).

The effects of the two HRP substrate concentrations, C and $H_2{O}_2$, predicted by the model are also similar to those observed experimentally (Figure 6 and Figure 7, respectively). The increase in signal with an increase in each substrate’s concentration, is consistent with the ping-pong kinetic model (Equation (1)), which predicts that HRP’s reaction rate would increase as either C or $H_2{O}_2$ increases. However, the SE for C is considerably stronger (SE = 8.9) than that for $H_2{O}_2$ (SE = 2.1), possibly because the $H_2{O}_2$ used in the experiments was much greater than the $K_m^{H_2{O}_2}$ value for HRP.

Both the experimental results and the mechanistic model (Figure 8) indicated a slightly higher EI signal in a mildly acidic bulk solution (pH = 6.2 or 6.6) than a neutral one (pH = 7). This trend is consistent with published reports that HRP oxidized substrates more rapidly in slightly acidic buffer than in neutral buffer [57]. One explanation for this effect is that pH (i.e., proton concentration) affects the thermodynamic driving force for the two-electron, two-proton electrochemical reduction of Q to C at the electrode. The $E_h$ value used in the model was measured as the midpoint potential ($E_{mid}$) of a cyclic voltammograms of an aqueous solution containing C and Q. Equation (3) shows that increasing pH would make $E_{mid}$ more negative, which would reduce the magnitude of (E-E_h)and thereby reduce the EI’s signal [43,44].

4.3. Integration of Dimensional Analysis and Flux Analysis to Determine Rate-Limiting Step

Previous mathematical models developed to describe kinetics of HRP on the electrodes [21,37,38,39,40] focused on the enzyme’s kinetics or were based on an assumption that the J is mass-transfer limited. In contrast, our model explicitly calculates the rates of all key reaction and mass transfer steps, all of which could limit the signal’s magnitude to some extent. Additionally, incorporation of Equations (2) and (3) allows effects of (E-E_h) and pH, respectively, to be predicted, even under conditions in which the commonly used assumption that ${\left[Q\right]}_{x=0}$ = 0 is invalid. Figure 9A shows that ${\left[Q\right]}_{x=0}$ decreases as the magnitude of (E-E_h) and the reduction rate of [Q] increases.

This extension of the model provides a significant improvement in accuracy over models based on the commonly used assumption that ${\left[Q\right]}_{x=0}$= 0. Figure 9B shows the percent error in predicted output current caused by that assumption under the experimentally realistic range of (E-E_h) values between −0.2 and −0.35 V. Notably, the error would have been about 15% for the (E-E_h) value of −0.3 V used by Kohli et al. [47,58], whose model was based on the assumption that that ${\left[Q\right]}_{x=0}$ = 0. A strategy of reducing the magnitude of the overpotential to reduce the modeling error due to that assumption is likely to be counterproductive, because that strategy would also reduce the biosensor’s signal.

The performance properties of an EI are controlled by the dynamics of the underlying transport and reaction steps that give rise to its J. We developed a mathematical framework that leverages dimensional analysis and the mechanistic model’s ability to predict the rates of the underlying steps to quantitatively assess the degree to which individual steps control the magnitude of the EI’s signal and its sensitivity (defined as change in J per unit change in analyte concentration). Examples of the approach are described below.

The dimensionless Damkohler number ($\sigma$) shown in Equation (17) expresses the ratio of the relative rates of enzymatic reaction ($\frac{V_{max}}{K_M}$) and diffusional mass transfer ($\frac{D_L}{L^2})$ of HRP’s substrates within the immunosensing layer [59]. Plugging constants from Table 1 into Equation (17) revealed that $\sigma$ for C and $H_2{O}_2$ were on the order of 10⁻⁵, indicating that the diffusion could provide C and $H_2{O}_2$ to the HRP orders of magnitude faster than the HRP could consume it [60,61]. This result indicates that the EI’s signal is not significantly limited by the diffusion rate within the immunosensing layer.

$${\sigma }^2=\frac{V_{max}{L}^2}{D_L{K}_M}$$

(17)

Flux-control analysis has been used to determine the extent to which the rates of individual enzymatic reactions in a biochemical reaction pathway limit the overall mass flux through that pathway [62]. We used a similar approach to determine the relative degrees to which the enzymatic and electrochemical reaction steps limit the magnitude of EI’s signal. We defined a current-control coefficient ($C_{Vi}^J)$for each reaction step (V_i) as the ratio of the percent change in the EI’s signal to the percent change in V_i while holding all other independent variables constant Equation (18). We used the mechanistic model to calculate an incremental change in J (∆J) resulting from an incremental change (∆$V_i$) in either the enzymatic reaction rate (simulated by changing the [HRP] value) or the electrochemical reaction rate (simulated by changing the (E-E_h) value). The incremental changes (∆J and ∆$V_i)$ were then used in place of the differentials ($\partial J$ and $\partial V_i$) in Equation (18) to calculate the $C_{Vi}^J$ values for both the enzymatic reaction and the electrochemical reaction across the range of (E-E_h) values used in this study.

$$\raisebox{1ex}{$\frac{\partial J}{J}$}\!\left/ \!\raisebox{-1ex}{$\frac{\partial V_i}{V_i}$}\right.=C_{Vi}^J$$

(18)

The $C_{Vi}^J$ values calculated by making incremental changes in [HRP] remained virtually 1.0 across the entire range of (E-E_h), for the HRP value listed in Table 1 (0.5 μm), as well as values ranging from 0.005 µM to 50 µM (results not shown). This result indicates that the EI’s signal is strongly limited by [HRP] over the entire range simulated. Consequently, the EI’s signal has the potential to be linearly correlated with the target analyte’s concentration, depending on the shape of the adsorption isotherm of the immobilized primary antibody for its target analyte.

In contrast, the $C_{Vi}^J$ values for the electrochemical reaction varied significantly across the range of overpotential used in this study (Figure 10) and exhibited a peak at about 1.5 at an (E-E_h) value of about −0.24 V. Although the predicted EI’s signal curve increased monotonically as the magnitude of (E-E_h) increased, the curve exhibited an inflection point at about the same (E-E_h) value the $C_{Vi}^J$ curve peaked. This observation suggests that a transition occurs at this point. For lower (E-E_h) magnitudes, increasing the magnitude strongly increases the EI’s signal. However, for higher (E-E_h) magnitudes, further increases in the (E-E_h) magnitude offer diminishing returns, suggesting that the peak in $C_{Vi}^S$ may mark an optimal operating overpotential in the absence of other overriding considerations, such as the presence of electrochemical interferents. For significantly higher (E-E_h) magnitudes, the J asymptotically approaches a maximum value, and the $C_{Vi}^J$ value approaches 0 (Supplementary Figure S3).

Because [HRP] would be expected to increase with the analyte concentration, the mechanistic model was also used to calculate the EI’s sensitivity (S) to HRP (defined in Equation (19)), as well as sensitivity-control coefficients ($C_{Vi}^S)$ (defined in Equation (20)).

$$S\simeq \frac{dJ}{d\left[HRP\right]}$$

(19)

$$\raisebox{1ex}{$\frac{\partial S}{S}$}\!\left/ \!\raisebox{-1ex}{$\frac{\partial V_i}{V_i}$}\right.=C_{Vi}^S$$

(20)

The S and $C_{Vi}^S$ values were calculated in a manner similar to that used to calculate $C_{Vi}^J$ values. The model was used to calculate incremental ∆J values resulting from incremental ∆[HRP] values. The incremental change values were substituted for differentials in Equations (19) and (20). The resulting S values and $C_{Vi}^S$ values (Figure 11) have shapes similar to the J and $C_{Vi}^J$ curves, respectively, shown in Figure 10. However, the peak in the $C_{Vi}^S$ curve occurs at a slightly different (E-E_h) value (−0.23 V) than the peak in the $C_{Vi}^J$ curve (−0.26 V). If an EI were operated near the peak of the $C_{Vi}^S$ curve, the sensitivity could be adjusted simply by making a relatively small change in the (E-E_h) value. The ability to adjust an EI’s sensitivity may provide users some flexibility to tailor the tradeoff between precision and usable analyte concentration range. Assuming an EI can only be used over a fixed range of amperometric signals, a higher sensitivity would be expected to provide higher precision but a narrower useful concentration range, whereas a lower sensitivity would be expected to provide lower precision but a wider useful concentration range.

Our mechanistic model could readily be adapted to other systems and used for other purposes than the ones shown here. For example, Kergaravat et al. optimized conditions for electrochemically measuring HRP’s oxidation rate for seven redox-active co-substrates [63]. For each co-substrate, they optimized the pH and the concentrations of HRP, $H_2{O}_2$, and the co-substrate. They also reported the measured current as a function of working electrode potential for each co-substrate. Our model could be fit to their experimental data by adjusting parameters, such as HRP’s kinetic constants and the midpoint potential for each co-substrate. In addition, the effect of pH on the electrochemical reduction of the oxidized co-substrate could be modeled using the Butler–Volmer Equation (2). The authors also reported dynamics of the HRP reaction (current vs. time), both to monitor the batch reaction after the co-substrate or $H_2{O}_2$ was added to the reaction mixture, and to monitor the rate of change in current following a change in the working electrode’s potential [63]. These results could be simulated by adding an accumulation term (i.e., a time derivative) for the chemical species being balanced on the left-hand side of Equations (4)–(6). The resulting system of partial differential equations could be solved to generate the types of time-dependent curves described above.

Our model could be adapted to describe other reporter enzymes that can be coupled to redox-active products that could be reversibly oxidized and reduced. For example, the commonly used reporter alkaline phosphatase (AP) [64] can hydrolyze phenylphosphate to phenol, which can then be oxidized to O-quinone by tryrosinase [65] and then measured electrochemically via reduction at an electrode [66,67,68]. The model could also be adapted to other redox-active phenolic or aromatic co-substrates of HRP, including o-phenylenediamine [18], 3,3′,5,5′-tertramethyl benzidine (TMB) [69], and p-aminophenol (PAP) [70]. Each co-substrate’s electrochemical reduction kinetics could be accurately described by substituting appropriate E_h and $K_0$ values into the model.

Another natural extension of the model would be to add equations that describe the equilibrium partitioning of the antigen binding to the immobilized primary antibody and the labeled secondary binding to the antigen-primary-antibody complex to form the sandwich molecular structure. This extension would allow the concentration of reporter enzyme (e.g., HRP) bound to the electrode, and thus allow the current generated by the EI to be predicted as a function of the antigen concentration in the sample. This capability would be useful for designing EI systems that meet desired performance specifications.

The novel dimensionless groups defined in this paper will enable future EI systems to be designed to meet performance specifications. For example, for an EI to have high sensitivity, it should be designed such that the CCC for the labeling enzyme is close to 1, indicating that that enzyme is highly rate-limiting. Also, an EI could be designed so that it’s sensitivity gives the desired balance between the precision and breadth of the useful analyte concentration range. Assuming an EI can only be used over a fixed range of amperometric signal, a higher sensitivity would be expected to provide higher precision but a narrower useful concentration range, whereas a lower sensitivity would be expected to provide lower precision but a wider useful concentration range.

Furthermore, our model could be extended to simulate some nonideal electrochemical effects. For example, when using cyclic voltammetry to determine the E_h value, the voltage difference between the oxidation and reduction peaks (i.e., the peak separation) can vary significantly with experimental conditions. A minimum peak separation of about 60 mV/n would be expected for an ideal system (rapid diffusion and a rapid, reversible redox reaction). However, a larger peak separation would be expected for slower electron-transfer kinetics and/or slower diffusion [71,72]. Thus, variations in peak separation could, in principle, be simulated by adjusting the model’s $D_L$ and $K_0$ values.

5. Conclusions

This study demonstrated the use of a novel, integrated experimental and modeling framework to analyze and optimize the performance of EIs. The experimental component included (1) deposition of an EI interface on the working electrode of miniature SPE arrays; (2) measurement of the performance properties of the resulting EIs for measuring the concentration of a surrogate protein antigen (mouse IgG); (3) use of a response-surface, statistical-design-of-experiments approach to optimize four independent variables: electrode overpotential, pH, and the concentrations of HRP’s two substrates ([C] and $\left[H_2{O}_2\right]$); and (4) development of a statistical model of the experimental data that empirically describes the effect of the four independent variables on the EI’s signal.

The modeling component included (1) development of a detailed, mechanistic model of the EI interface that described the rates of the mass-transfer and reaction steps that gave rise to the EI’s signal; (2) use of the statistical model of the experimental data to help validate the mechanistic model; and (3) integration of dimensional analysis, principles of flux-control analysis, and the mechanistic model’s predictive capabilities to obtain unprecedented insight into which steps control the magnitude of the EI’s signal and its sensitivity to the target analyte.

The EI developed in this study had a limit of detection of 1 ng/mL, and an inter-assay/intra-assay variation of less than 5%. The mechanistic model was able to reproduce experimentally observed effects of the four independent variables on the EI’s signal. Calculation of Damkohler numbers indicated that diffusion of HRP’s substrates in the biocatalytic layer did not limit the EI’s performance at the overpotential of −0.3 V. Calculation of current-control and sensitivity-control coefficients analyses provided new insight into the extent to which the enzymatic and electrochemical reactions limited both the EI’s signal and its sensitivity over the experimentally relevant range of (E-E_h) values.

The novel, integrated experimental and modeling framework presented in this study provides unprecedented capabilities to design, optimize, and validate EIs for diverse applications. Its ability to quickly identify key mass transfer or reaction step(s) that limit(s) could guide strategies to overcome such limitation(s) and thereby reduce time required to develop new commercial EI systems. Also, the predictive power of the mechanistic model could, in principle, enable EIs to be designed a priori to meet specifications and enable rapid, in-silico hypothesis testing that could accelerate FDA approval of EI systems for healthcare applications.