An Approximate Analytical Method for the Evaluation of the Concentrations and Current for Hybrid Enzyme Biosensor

Mathematical modeling of amperometric biosensor with cyclic reaction is discussed. Analytical expressions pertaining to the concentrationofsubstrate,cosubstrate,reducingagentandmedialproductandcurrentforhybridenzymebiosensorareobtainedintermsofThielemoduleandsaturationparameters.Inthispaper,apowerfulanalyticalmethod,calledhomotopyanalysismethod(HAM)isusedtosolvethesystemofnonlineardifferentialequations.Furthermore,inthisworkthenumericalsimulationoftheproblemisalsoreportedusingScilab/Matlabprogram.Ouranalyticalresultsarecomparedwithsimulationresults.Agoodagreementbetweenanalyticalandnumericalresultsisnoted.


Introduction
Biosensor ( Figure 1) is a device that uses specific biochemical reactions mediated by isolated enzymes, immunosystems, tissues, organelles, or whole cells to detect chemical compounds usually by electrical, thermal, or optical signals [1]. They involve a biological (recognition) element and a transduction element. The biological or recognition element may be an antibody, an enzyme, DNA, RNA, a whole cell, or a whole organ or system. The transduction element, wherein the biological event or signal is converted to a measurable signal, may include anyone of the following forms: chemical, electrical, magnetic, mechanical, optical, or thermal.
The biosensor was first described by Clark and Lyons in 1962, when the term enzyme electrode was adopted [2]. The term "biosensor" was introduced by Cammann in 1977 [3]. Since then, research communities from various fields such as physics, chemistry, and material science have come together to develop more sophisticated, reliable, and mature biosensing devices for applications in the fields of medicine, agriculture, biotechnology, as well as in the military for bioterrorism detection and prevention [4]. Biosensors offer the prospects of simplified, virtually nondestructive analysis of turbid biological fluids. Also, biosensors for medical care have demanded the greatest attention for technical development [5].
Amperometric electrodes have been used in the design of biosensors for glucose, aminoacids, and other molecules [6][7][8][9]. In cases of amperometric enzyme biosensors the potential at the electrode is held constant while the current flow is measured. Amperometric biosensors are quite sensitive and more suited for mass production than the potentiometric ones [10,11]. Electropolymerized films offer wide immobilization capabilities and extremely large diversity in the development of biosensors [12]. The development of amperometric biosensors is an area of growing interest in many branches of science [13][14][15][16][17]. Hybrid biosensors are biosensors with more than one biosensitive material-enzyme, tissue microorganism [18] or other agents [19]. Hybrid biosensors are less selective than enzymatic ones, but the choice of the type of biosensor should be made according to several additional factors like availability of biocatalysts, their stability, and sensitivity required [20].
Recently, Rajendran and coworkers [21][22][23] obtained the analytical expressions of substrate and product for amperometric biosensors. Eswari and Rajendran [24,25] obtained the steady-state current at microdisk biosensor and microcylinder biosensor. Recently, Rangelova [26] found the concentration profiles of substrate, cosubstrate, reducing agent, medial product using finite difference technique and Matlab program. To our knowledge, no general simple analytical expressions for the concentrations and current for hybrid biosensors have been reported for all values of the Thiele module and saturation parameters. However, in general, analytical solution of nonlinear differential equations is more interesting and useful than purely numerical solutions as they are amenable to various kinds of manipulation and data analysis. Analytical expressions are usually derived from the basic physical principles and free from numerical dispersions and other truncation errors that often occurred in numerical simulations. For this reason, in this paper we have derived an analytical expression for the concentrations and current for hybrid enzyme biosensor for all values of the normalized parameters using homotopy analysis method (HAM).

Mathematical Formulation.
In enzyme-based catechol biosensor, the cyclic reaction scheme for the substrate, cosubstrate, reducing agent and medial product can be represented as follows [26]: where 1 and 2 are reaction rate constants, 1 is the first product and 2 is the second product (dehydroascorbic acid), is the measured substrate (catechol), is the cosubstrate (oxygen), is the reducing agent ( -ascorbic acid), and is the medial product (1, 2 benzoquinone).
A cyclic reaction between catechol and 1,2-benzoquinone takes place by combining the tyrosinase reaction and the chemical reduction of 1,2 benzoquinone to catechol by Lascorbic acid. This acid is well known as an effective reducing agent, and 1,2 benzoquinone would be reduced to catechol and may drive the reaction in the opposite direction to that of enzymatic oxidation of catechol. The oxygen consumed in the enzymatic reaction is not compensated for by the chemical reduction. Therefore, if L-ascorbic acid does not affect the enzyme activity of tyrosinase, a cyclic reaction should take place and the consumption of the dissolved oxygen will continue until its concentration becomes zero  [27]. The differential equations for this ping-pong mechanism at the steady-state condition are as follows [26]: Here [ ], [ ], [ ], and [ ] (mM) denote the concentration of substrate, cosubstrate, reducing agent and medial product, respectively. , , , and (m 2 /s) are diffusion coefficients, 2 , 3 , and 4 (mmol/(l.s)) are the reaction rate constants, ( m) is the distance coordinate, , , and (mM) are the reaction rate constants, and (mmol/(l.s)) is the maximal rate. The boundary conditions are given by [26] where is the active membrane thickness. The diffusion limiting current can be expressed as follows: where is the number of electrons taking part in the electrochemical reaction on the cathode, is the Faraday's constant, and is the surface of cathode.

Normalized
the nonlinear reaction/diffusion equation (2) take the following normalized form: The transformed boundary conditions are The dimensionless form of the current is given by

Homotopy Analysis Method.
The homotopy analysis method (HAM) [28][29][30][31][32][33][34][35][36] is a general analytic method to get series solutions of various types of nonlinear equations, including ordinary differential equations, partial differential equations, and coupled nonlinear equations. Unlike perturbation methods, the HAM is independent of small/large physical parameters. More importantly, different from all perturbation and traditional nonperturbation methods, the HAM provides us with a simple way to ensure the convergence of solution series. Besides, different from all perturbation and previous nonperturbation methods, the HAM provides us with great freedom to choose proper base functions to approximate a nonlinear problem [29,34].
Liao's book [29] for the homotopy analysis method was first published in 2003. Now, more and more researchers have been successfully applying this method to various nonlinear problems in science and engineering. In this paper we employ HAM to solve the nonlinear differential equations (6), (7), and (9). The basic concept of Homotopy analysis method is given in Appendix A.

Numerical Simulation.
In order to investigate the accuracy of this analytical method with a finite number of terms, the system of differential equations ((6)-(9)) also solved by numerical methods. The function pdepe (finite element method) in Scilab/Matlab software which is a function of solving PDE is used to solve these nonlinear equations [37]. The Scilab/Matlab program is also given in Appendix C.
To validate the results, the convergence studies are carried out. The convergence region of auxiliary parameter h is given in the Appendix D. To show the efficiency of the present analytical method, our results are compared with the numerical solution (Scilab/Matlab program) in Tables 1, 2, 3,  and 4 and Figures 2-6. The average relative errors between our analytical results and numerical results are 0.23%, 0.11%, 0.03%, and 0.01% for the concentrations of substrate, cosubstrate, reducing agent and medial product, respectively.

Effect of the Thiele Module.
The concentration of substrate, co-substrate, reducing agent and medial product depends upon Thiele module and saturation parameters. The Thiele module 2 (= 2 / ) essentially compares the rate of enzyme reaction ( / ) and diffusion in the enzyme layer ( 2 / ). We observe the rise and downfall of concentration profiles in two cases. (i) If Thiele module is small ( 2 < 1), then enzyme kinetics predominates in the biosensor response. The overall kinetics is governed by the total amount of active enzyme. (ii) The response is under diffusion control, if the Thiele module is large ( 2 > 1), which is observed at high catalytic activity and active membrane thickness or at low reaction rate constant ( ) or diffusion coefficient values ( ).
The concentration profiles for the four reactants for some fixed values of parameters are shown in Figure 2. The analytical results are compared with the numerical result, graphically. Upon comparison, it is evident that both results give satisfactory agreement. The concentration of substrate increases and attains its maximum, when Thiele module and saturation parameter increase (refer to Figure 3). Therefore the profile deviates more from the linearity. The concentration profile representing medial product increases as Thiele module and saturation parameter increase (refer to Figure 4). From Figure 5, it is inferred that the increasing value of coordinate decreases the concentration of reducing agent . Figure 6 represents the concentration of co-substrate versus the normalized distance coordinate. It reaches the minimum value zero at electrode interface.

Influence of Active Membrane
Thickness. The active membrane thickness is one of the important technical parameters and it has a considerable effect on the concentration profiles for the four reactants. Also, the Thiele module is directly proportional to the thickness of the membrane. The concentration profiles depend significantly on the membrane thickness . If the active membrane thickness is large ( > 50 m or thick membrane), the concentration profiles of substrate and medial product increase whereas the concentration of reducing agent and co-substrate decreases. If the membrane thickness is small ( < 50 m or thinner membrane), the concentration profiles of substrate, medial product and reducing agent have uniform values. These uniform values are equal to the concentration of the above reactants at = 0. The concentration profile of co-substrate increases when thickness of the membrane decreases. The increasing value is not significant. Also, the concentration of co-substrate increases when the saturation parameter increases (refer to Figures 3-6).  Figure 5: Normalized concentration profiles of reducing agent for various values of the Thiele module 2 , active membrane thickness , and the normalized parameters are plotted using (14). The key to the graph: stacked line represents (14) and dotted line represents the numerical simulation.

Current
Response. The normalized current versus the diffusion coefficient ratio 1 is calculated at different values of the active membrane thickness . The results obtained for various values of the normalized parameters are depicted in Figure 7(a). The current response increases when the active membrane thickness ( > 50 m) increases. Also, for thinner membrane ( < 50 m) the value of the current is zero. In Figure 7(b), the current response increases, as the saturation parameter increases.

Conclusions
The theoretical model of hybrid amperometric enzyme biosensor with cyclic reaction and biochemical amplification for steady-state condition is discussed. The system of three nonlinear differential equations for ping-pong enzyme kinetics has been solved analytically. Influence of Thiele module and active membrane thickness is investigated. The obtained results have a good agreement with those obtained using numerical method. This analytical result will be useful in sensor design, optimization, and prediction of the electrode response. Using this result, the action of biosensor is analyzed at critical concentration of substrate and enzyme activities.
Theoretical results obtained in this paper can also be used to analyze the effect of different parameters such as active membrane thickness and saturation parameters.

A. Basic Idea of Liao's Homotopy Analysis Method
Consider the following differential equation [38]: where is a nonlinear operator, t denotes an independent variable, and ( ) is an unknown function. For simplicity, we ignore all boundary or initial conditions, which can be treated in a similar way. By means of generalizing the conventional Homotopy method, Liao constructed the so-called zeroorder deformation equation as Applying −1 on both sides of (A.8), we get In this way, it is easy to obtain for ≥ 1, at th order; we have (A.11) When → +∞, we get an accurate approximation of the original equation (A.1). For the convergence of the previous method we refer the reader to Liao [29]. If (A.1) admits unique solution, then this method will produce the unique solution. If (A.1) does not possess unique solution, the HAM will give a solution among many other (possible) solutions.

D. Determining the Validity Region of ℎ
The analytical solution represented by (12), (13), and (15) contains the auxiliary parameter h, which gives the convergence region and rate of approximation for homotopy analysis method. The analytical solution should converge. It should be noted that the auxiliary parameter h controls the convergence and accuracy of the solution series. In order to define region such that the solution series is independent of ℎ, a multiple of ℎ curves are plotted. The region where the distribution of , , and versus h is a horizontal line is known as the convergence region for the corresponding function. The common region among ( ), ( ), and ( ) is known as the overall convergence region. To study the influence of h on the convergence of solution, the h curves of (0.1) are plotted in Figure 8. This figure clearly indicates that the valid region of h is about (−1.5 to − 0.1). Similarly we can find the value of the convergence-control parameter h for different values of constant parameters.