An Enhanced Segment Particle Swarm Optimization Algorithm for Kinetic Parameters Estimation of the Main Metabolic Model of Escherichia Coli

Abstract

Building a biologic model that describes the behavior of a cell in biologic systems is aimed at understanding the physiology of the cell, predicting the production of enzymes and metabolites, and providing a suitable data that is valid for bio-products. In addition, building a kinetic model requires the estimation of the kinetic parameters, but kinetic parameters estimation in kinetic modeling is a difficult task due to the nonlinearity of the model. As a result, kinetic parameters are mostly reported or estimated from different laboratories in different conditions and time consumption. Hence, based on the aforementioned problems, the optimization algorithm methods played an important role in addressing these problems. In this study, an Enhanced Segment Particle Swarm Optimization algorithm (ESe-PSO) was proposed for kinetic parameters estimation. This method was proposed to increase the exploration and the exploitation of the Segment Particle Swarm Optimization algorithm (Se-PSO). The main metabolic model of E. coli was used as a benchmark which contained 172 kinetic parameters distributed in five pathways. Seven kinetic parameters were well estimated based on the distance minimization between the simulation and the experimental results. The results revealed that the proposed method had the ability to deal with kinetic parameters estimation in terms of time consumption and distance minimization.

1. Introduction

Mathematical modeling is a compelling strategy for predicting and recognizing the biologic behavior of a cell’s system. The mathematical model can create and predict the results of the empirical hypotheses that can be used to examine the process. In this regard, the biologic model studies the chemical metabolism and the pathways to construct a kinetic model that can give a clear understanding of a cell mechanism [1].

Therefore, building and developing a kinetic model requires a large number of parameters such as kinetic parameters, potential of hydrogen (pH) and initial enzymatic and metabolite concentration. In the same vein, the building and the developing of the kinetics model gives a clear picture of cell behavior when it is needed in biotechnology applications [2]. In addition, building and developing a kinetic model is a challenging task because it requires kinetic parameter estimation [3].

Kinetic parameters estimation is the finding of the nearest value that minimizes the distance between the simulated and the real experimental results [4]. Kinetic parameters estimation is a complex task because of the nonlinearity of the kinetic model, and it is usually measured in different conditions and time consumption [5]. Usually, the kinetic parameters stored in the database are insufficient for building accurate kinetic model due to the cell behavior in different measurements stated by different laboratories utilizing various in vitro models and conditions [6].

Normally, the model is built depending on time derivative expressions using Ordinary Differential Equations (ODEs) which describe the changes in a state or quantity of interest over time. Recently, several kinetic models of cell metabolism are based on ODEs which have been developed to detect time-dependent changes in metabolic concentration. For example, models of glycolysis and pentose phosphate (PP) pathways in E. coli [7,8] and the tri-carboxylic acid (TCA) cycle in Dictyostelium discoideum [9] have been constructed. Furthermore, large scale model integration has been performed. For example, the model of [7] has been integrated with the model of (TCA) cycle [8,10] and amino acid biosynthesis [11]. Researchers have studied the glycolysis pathway to clone a new host and observe its changes and effects [12].

Therefore, because of the difficulty in estimating kinetic parameters, many researchers lately have used the metaheuristic optimization algorithms’ methods to estimate the kinetics of the E. coli model and other biologic models. Some of these metaheuristic algorithms were used in order to estimate the kinetic parameters employed by [13,14,15,16,17,18,19,20,21] and some of these algorithms have used experimental data taken from [7,22,23] to investigate their abilities in the kinetic parameters estimation.

In this regards, the biologic kinetic models contain hundreds to thousands of kinetic parameters which have caused serious challenges in parameter estimation due to the high dimension of search space that is required to be explored. Thus, in high dimensional kinetic models (containing hundreds to thousands of kinetic parameters), the procedures are considered to be expensive in terms of computational cost and the performance of the aforementioned algorithms tend to deteriorate, which may lead to low accuracy [6].

The Differential Evolution (DE) algorithm is an example that is widely used and studied in parameter estimation of metabolic model [6,17]. Time consumption is the main weakness of (DE) algorithm; this makes it difficult for the (DE) algorithm to tune its parameters when it involves a large number of processors and multiple local searches. In addition, the Genetic Algorithm (GA) belongs to the class of evolutionary algorithms which almost share the same idea as DE. This algorithm is a popular method that is used in parameter estimation of metabolic model. The main weakness of this algorithm is the high computational time which mostly does not perform effectively when it is compared to DE and Particle Swarm Optimization (PSO) algorithms [24,25,26].

Since 1995, the PSO algorithm has been making great impact in solving many problems in different fields such as [14,27,28,29]. This method was proposed by [30] and it was adopted in simulating flock of birds and a school of fish in order to locate their food. The algorithm particles during their movement in the process of searching share the information with each other and then they search for their destination stochastically and independently. Therefore, the particles’ movements in the search space move towards the optimum solution by taking different paths and directions.

Moreover, the method of Se-PSO was proposed and developed based on the PSO algorithm in order to increase the computational time, and the accuracy of PSO was presented by [31]. This development is based on the recognition of the local and global point problem of PSO [30]. In this regard, the segmentation divides the particles into groups so as to move together towards the optimal solution. This algorithm was adapted for large-scale kinetic parameters of E. coli model [32] and governor–turbine model identification of single area power plant [28,32] was used in a linear inertia weight. Due to the nonlinearity of the model which uses linear inertia weight ($\omega$) by [32], this inertia weight affects the exploration and exploitation of the particles. Consequently, there is need for large exploration at the beginning and small exploitation at the end of the algorithm execution. The purpose for this is to avoid the local optima trap so as to increase the efficiency of the kinetic parameters estimation with a view to minimizing the model distances in reasonable time.

In this study, the method of Se-PSO algorithm was enhanced by adding a damping process which helped the particles to have a wide exploration at the beginning and small exploitation at the end during the algorithm execution. The aim was to investigate the amount of each solution to be found at the search space. Thus, the enhanced Se-PSO was adapted to estimate large-scale kinetic parameters of E. coli and to increase the efficiency of the Se-PSO algorithm in terms of exploration and exploitation.

The other sections of the study are organized as follows: Section 2 describes the problem formulation. Section 3 describes the materials and methods. Section 4 describes the obtained results. Section 5 gives the conclusion of the findings of this study.

2. Problem Formulation

Kinetic parameters estimation of the kinetic model was formulated to find the best set of kinetics that minimizes the model simulation differences with real experimental data. This formulation was set as a global optimization problem. In this regard, optimization algorithm methods can be adapted in order to estimate the kinetic parameters and minimize the differences between the model simulation and real experimental result. The formulation of this problem mathematically was set as a fitness function in Equation (1) [13]:

$$\mathit{f}\mathit{i}\mathit{t}\mathit{n}\mathit{e}\mathit{s}\mathit{s}={\displaystyle \sum }_{\mathit{i}=1}^{\mathit{R}}\left|\frac{{\mathit{y}}_{\mathit{s},\mathit{i}}-{\mathit{y}}_{\mathit{e},\mathit{i}}}{{\mathit{y}}_{\mathit{e},\mathit{i}}}\right|$$

(1)

where $R$ is the number of metabolites $i$ experimental data; so $R=12$, $y_{s,i}$ is the estimated metabolites concentration $s$ of metabolite $i$, $y_{e,i}$ is the real metabolites experimental data $e$ of metabolite $i$.

3. Materials and Methods

3.1. Kinetic Model Structure

The dynamic model of the main metabolic of E. coli formulated by [10] was used as a benchmark. This model, which consists of (glycolysis, pentose phosphate, TCA cycle, gluconeogenesis and glyoxylate) pathways in addition to acetate formation with phosphotransferase system, is described in Figure 1.

This model has 23 metabolites, 28 enzymatic reactions, 10 co-factors (e.g., Adenosine triphosphate (ATP), Coenzyme A (COA), Nicotinamide adenine dinucleotide phosphate (NADPH)) and 172 kinetic parameters. The rate at which the concentration of the metabolite in the considered model changes is given by Equation (2) as follows:

$$\frac{d{C}_i}{dt}={{\displaystyle \sum }}_{j=1}{S}_{i,j}{v}_j-\mu C_i$$

(2)

where $C_i$ is the concentration of metabolite $i$, $v_i$ is the rate of reaction $j$, $S_{i,j}$ is the stoichiometric coefficient of the metabolite $i$ in the reaction $j$ and $\mu$ is the specific growth rate. Thus, the term $\mu C_i$ represents the dilution effect due to biomass growth. The model mass balance equations and the kinetic rate equations of the considered model are stated in Appendix A

Table A1. Mass balance equation.

Metabolites	Mass Balance Description
cell ($\mathrm{X}$)	$\frac{\mathrm{d}\left[\mathrm{X}\right]}{\mathrm{dt}}=\mathsf{\mu }\left[\mathrm{X}\right]$
extra glucose (${\mathrm{GLC}}^{\mathrm{ex}}$)	$\frac{\mathrm{d}\left[{\mathrm{GLC}}^{\mathrm{ex}}\right]}{\mathrm{dt}}=-{\mathrm{v}}_{\mathrm{PTS}}\left[\mathrm{X}\right]$
glucose-6-phosphate ($\mathrm{G}6\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{G}6\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{PTS}}-{\mathrm{v}}_{\mathrm{PGI}}-{\mathrm{v}}_{\mathrm{G}6\mathrm{PDH}}-\mathsf{\mu }\left[\mathrm{G}6\mathrm{P}\right]$
fructose 6-phospahte ($\mathrm{F}6\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{F}6\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{PGI}}-{\mathrm{v}}_{\mathrm{PFK}}+{\mathrm{v}}_{\mathrm{TKTB}}+{\mathrm{v}}_{\mathrm{TAL}}-\mathsf{\mu }\left[\mathrm{F}6\mathrm{P}\right]$
fructose 1,6-phosphate ($\mathrm{FDP}$)	$\frac{\mathrm{d}\left[\mathrm{FDP}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{PFK}}-{\mathrm{v}}_{\mathrm{ALDO}}-\mathsf{\mu }\left[\mathrm{FDP}\right]$
glyceraldehyde 3-phosphate ($\mathrm{GAP}$)	$\frac{\mathrm{d}\left[\mathrm{GAP}\right]}{\mathrm{dt}}=2{\mathrm{v}}_{\mathrm{ALDO}}-{\mathrm{v}}_{\mathrm{GAPDH}}+{\mathrm{v}}_{\mathrm{TKTA}}+{\mathrm{v}}_{\mathrm{TKTB}}-{\mathrm{v}}_{\mathrm{TAL}}-\mathsf{\mu }\left[\mathrm{GAP}\right]$
phosphoenolpyruvate ($\mathrm{PEP}$)	$\frac{\mathrm{d}\left[\mathrm{PEP}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{GAPDH}}+{\mathrm{v}}_{\mathrm{PCK}}-{\mathrm{v}}_{\mathrm{PTS}}-{\mathrm{v}}_{\mathrm{PYK}}-{\mathrm{v}}_{\mathrm{PPC}}-\mathsf{\mu }\left[\mathrm{PEP}\right]$
pyruvate ($\mathrm{PYR}$)	$\frac{\mathrm{d}\left[\mathrm{PYR}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{PYK}}+{\mathrm{v}}_{\mathrm{PTS}}+{\mathrm{v}}_{\mathrm{MEZ}}-{\mathrm{v}}_{\mathrm{PDH}}-\mathsf{\mu }\left[\mathrm{PYR}\right]$
acetyl-CoA ($\mathrm{AcCoA}$)	$\frac{\mathrm{d}\left[\mathrm{AcCoA}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{PDH}}+{\mathrm{v}}_{\mathrm{ACS}}+{\mathrm{v}}_{\mathrm{CS}}-{\mathrm{v}}_{\mathrm{PTA}}-\mathsf{\mu }\left[\mathrm{AcCoA}\right]$
isocitrate ($\mathrm{ICIT}$)	$\frac{\mathrm{d}\left[\mathrm{ICIT}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{CS}}-{\mathrm{v}}_{\mathrm{ICDH}}-{\mathrm{v}}_{\mathrm{ICL}}-\mathsf{\mu }\left[\mathrm{ICIT}\right]$
2-keto-d-gluconate ($2\mathrm{KG}$)	$\frac{\mathrm{d}\left[2\mathrm{kG}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{ICDH}}-{\mathrm{v}}_{2\mathrm{KGDH}}-\mathsf{\mu }\left[2\mathrm{KG}\right]$
succinate ($\mathrm{SUC})$	$\frac{\mathrm{d}\left[\mathrm{SUC}\right]}{\mathrm{dt}}={\mathrm{v}}_{2\mathrm{KGDH}}+{\mathrm{v}}_{\mathrm{ICL}}-{\mathrm{v}}_{\mathrm{SDH}}-\mathsf{\mu }\left[\mathrm{SUC}\right]$
fumarate ($\mathrm{FUM}$)	$\frac{\mathrm{d}\left[\mathrm{FUM}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{SDH}}-{\mathrm{v}}_{\mathrm{FUM}}-\mathsf{\mu }\left[\mathrm{FUM}\right]$
malate ($\mathrm{MAL}$)	$\frac{\mathrm{d}\left[\mathrm{MAL}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{FUM}}+{\mathrm{v}}_{\mathrm{MS}}-{\mathrm{v}}_{\mathrm{MDH}}-{\mathrm{v}}_{\mathrm{MEZ}}-\mathsf{\mu }\left[\mathrm{MAL}\right]$
oxaloacetate ($\mathrm{OAA}$)	$\frac{\mathrm{d}\left[\mathrm{OAA}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{MDH}}+{\mathrm{v}}_{\mathrm{PPC}}-{\mathrm{v}}_{\mathrm{CS}}-{\mathrm{v}}_{\mathrm{PCK}}-\mathsf{\mu }\left[\mathrm{OAA}\right]$
glyoxylate ($\mathrm{GOX}$)	$\frac{\mathrm{d}\left[\mathrm{GOX}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{ICL}}-{\mathrm{v}}_{\mathrm{MS}}-\mathsf{\mu }\left[\mathrm{GOX}\right]$
acetyl phosphate ($\mathrm{ACP}$)	$\frac{\mathrm{d}\left[\mathrm{ACP}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{PTA}}-{\mathrm{v}}_{\mathrm{ACK}}-\mathsf{\mu }\left[\mathrm{ACP}\right]$
acetate ($\mathrm{ACE}$)	$\frac{\mathrm{d}\left[\mathrm{ACE}\right]}{\mathrm{dt}}=\left({\mathrm{v}}_{\mathrm{ACK}}-{\mathrm{v}}_{\mathrm{ACS}}\right)\left[\mathrm{X}\right]$
6-phosphogluconolactone ($6\mathrm{PG}$)	$\frac{\mathrm{d}\left[6\mathrm{PG}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{G}6\mathrm{PDH}}-{\mathrm{v}}_{6\mathrm{PGDH}}-\mathsf{\mu }\left[6\mathrm{PG}\right]$
ribose 5-phosphate ($\mathrm{Ru}5\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{Ru}5\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{6\mathrm{PGDH}}-{\mathrm{v}}_{\mathrm{RPE}}-{\mathrm{v}}_{\mathrm{RPI}}-\mathsf{\mu }\left[\mathrm{Ru}5\mathrm{P}\right]$
ribulose 5-phosphoenolpyruvate ($\mathrm{R}5\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{R}5\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{RPI}}-{\mathrm{v}}_{\mathrm{TKTA}}-\mathsf{\mu }\left[\mathrm{R}5\mathrm{P}\right]$
xylulose 5-phsophate ($\mathrm{Xu}5\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{Xu}5\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{RPE}}-{\mathrm{v}}_{\mathrm{TKTA}}-{\mathrm{v}}_{\mathrm{TKTB}}-\mathsf{\mu }\left[\mathrm{Xu}5\mathrm{P}\right]$
sedoheptulose 7-phosphate ($\mathrm{S}7\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{S}7\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{TKTA}}-{\mathrm{v}}_{\mathrm{TAL}}-\mathsf{\mu }\left[\mathrm{S}7\mathrm{P}\right]$
erythrose 4-phsophate ($\mathrm{E}4\mathrm{P}$)	$\frac{\mathrm{d}\left[\mathrm{E}4\mathrm{P}\right]}{\mathrm{dt}}={\mathrm{v}}_{\mathrm{TAL}}-{\mathrm{v}}_{\mathrm{TKTB}}-\mathsf{\mu }\left[\mathrm{E}4\mathrm{P}\right]$

and

Table A2. Reaction rate equations.

Reactions	Kinetic Equation
cell growth ($X$)	$\{\begin{array}{c}{\mu }_m\left(1-\frac{\left[X\right]}{X_m}\right)\left(\frac{\left[GL{c}^{ex}\right]}{K_s+\left[GL{c}^{ex}\right]}\right)k_{ATP}{v}_{ATP}(.),\left(\left[GL{c}^{ex}\right]>0\right)\\ \frac{{\mu }_{mA}\left[Ac{e}^{ex}\right]}{K_{sA}+\left[Ac{e}^{ex}\right]} k_{ATP}{v}_{ATP}(.),\left(\left[GL{c}^{ex}\right]\le 1and\left[Ac{e}^{ex}\right]>0\right)\end{array}$
PTS (phosphotransferase system)	$\frac{v_{PTS}^{max}\left[GL{c}^{ex}\right]\frac{\left[PEP\right]}{PYR}}{\left(K_{a1}+K_{a2}\frac{\left[PEP\right]}{\left[PYR\right]}+K_{a3}\left[GL{c}^{ex}\right]+\left[GL{c}^{ex}\right]\frac{\left[PEP\right]}{\left[PYR\right]}\right)\left(1+\frac{{\left[G6P\right]}^{n}G6P}{K_{G6P}}\right)}$
PGI (phosphoglucoseisomerase)	$\frac{v_{PGI}^{max}\left(\left[G6P\right]-\frac{\left[F6P\right]}{K_{eq}}\right)}{K_{G6P}\left(1+\frac{\left[F6P\right]}{K_{F6P}\left(1+\frac{\left[F6P\right]}{K_{6pginh}^{F6P}}\right)}+\frac{\left[6PG\right]}{K_{6pginh}^{G6P}}\right)+G6P}$
PFK (phospho-fructokinase)	$\frac{v_{PFK}^{max}{K}_{ATP}\left[F6P\right]}{K_{(ATP,ADP)}\left(\left[F6P\right]+K_z^{F6P}\frac{K_{b(ADP,AMP)+\frac{\left[PEP\right]}{K_{PEP}}}}{K_{a(ADP,AMP)}}\right)\left(1+\frac{L_{pfk}}{{\left(1+\left[F6P\right]\left(\frac{K_{a(ADP,AMP)}}{K_s^{F6P}\left(K_{b(ADP,AMP)}+\frac{\left[PEP\right]}{K_{PEP}}\right)}\right)\right)}^{n_{PFK}}}\right)}$
Aldo (aldolase)	$\frac{v_{ALDO}^{max}\left(\left[FDP\right]-\frac{\left[DHAP\right]\left[GAP\right]}{K_{eq}}\right)}{\left(K_{FDP}+\left[FDP\right]+\frac{K_{GAP}\left[DAHP\right]}{[K_{eq}{V}_{bif}]}+\frac{K_{DHAP}\left[GAP\right]}{[K_{eq}{V}_{bif}]}+\frac{\left[FDP\right]\left[GAP\right]}{K_{inh}^{PEP}}+\frac{\left[DHAP\right]\left[GAP\right]}{K_{eq}{V}_{bif}}\right)}$
GAPDH (glyceraldehyde 3-phosphate dehydrogenase)	$\frac{v_{GAPDH}^{max}\left(\left[GAP\right]-\frac{\left[PEP\right]\left[NADH\right]}{K_{eq}\left[NAD\right]}\right)}{\left(K_{GAP}\left(1+\frac{\left[PEP\right]}{K_{PGP}}\right)+\left[GAP\right]\right)\left(\frac{K_{NAD}}{NAD}\left(1+\frac{\left[NADH\right]}{K_{NADH}}\right)+1\right)}$
PYK (pyruvate kinase)	$\frac{v_{PYK}^{max}\left[PEP\right]{\left(\frac{PEP}{K_{PEP}}+1\right)}^{n_{pyk}-1}\left[ADP\right]}{K_{PEP}\left(L_{PYK}{\left(\frac{1+\frac{\left[ATP\right]}{K_{ATP}}}{\frac{\left[FDP\right]}{K_{FDP}}+\frac{\left[AMP\right]}{K_{AMP}}+1}\right)}^{n_{pyk}}+{\left(\frac{\left[PEP\right]}{K_{PEP}}+1\right)}^{n_{pyk}}\right)\left(\left[ADP\right]+K_{ADP}\right)}$
Ppc (PEP carboxylase)	$\frac{K_1+K_2\left[AcCOA\right]+K_3\left[FDP\right]+K_4\left[AcCOA\right]\left[FDP\right]}{1+K_5\left[AcCOA\right]+K_6\left[FDP\right]}\left(\frac{\left[PEP\right]}{K_m+\left[PEP\right]}\right)$
G6PDH	$\frac{v_{G6PDH}^{max}\left[G6P\right]\left[NADP\right]}{\left(\left[G6P\right]+K_{g6p}\right)\left(1+\frac{\left[NADPH\right]}{K_{ndph}}\right)\left(K_{nadp}\left(1+\frac{\left[NADPH\right]}{K_{nadph}}\right)+NADP\right)}$
PGDH (phosphogluconate dehydrogenase)	$\frac{v_{PGDH}^{max}\left[6PG\right]\left[NADP\right]}{\left(\left[6PG\right]+K_{6pg}\right)\left(\left[NADP\right]+K_{nadp}\left(1+\frac{\left[NADPH\right]}{K_{nadph}}\right)\left(1+\frac{\left[ATP\right]}{K_{atp}}\right)\right)}$
Rpe (ribulose-phosphate 3-epimerase)	$v_{Rpe}^{max}\left(\left[Ru5P\right]-\frac{\left[R5P\right]}{K_{eq}^{Rpe}}\right)$
Rpi (ribose-5-phosphate isomerase)	$v_{Rpi}^{max}\left(\left[Ru5P\right]-\frac{\left[R5P\right]}{K_{eq}^{Rpi}}\right)$
TktA (transketolase A)	$v_{TKtA}^{max}\left(\left[R5P\right]\left[Xu5P\right]-\frac{\left[S7P\right]\left[GAP\right]}{K_{eq}^{TKtA}}\right)$
TktB (transketolase B)	$v_{TKtB}^{max}\left(\left[Xu5P\right]\left[E4P\right]-\frac{\left[F6P\right]\left[GAP\right]}{K_{eq}^{TKtB}}\right)$
Tal (tyrosine ammonia lyase)	$v_{TaL}^{max}\left(\left[GAP\right]\left[S7P\right]-\frac{\left[E4P\right]\left[F6P\right]}{K_{eq}^{TKtB}}\right)$
PcK (PEP carboxykinase)	$v_{PcK}^{max}\left(\frac{\left[OAA\right]\frac{\left[ATP\right]}{\left[ADP\right]}}{K_m^{OAA}\frac{\left[ATP\right]}{\left[ADP\right]}+\left[OAA\right]\frac{\left[ATP\right]}{\left[ADP\right]}+\frac{K_i^{ATP}{K}_m^{OAA}}{K_i^{ADP}}+\frac{K_i^{ATP}{K}_m^{OAA}}{K_m^{PEP}{K}_i^{ADP}}\left[PEP\right]+\frac{K_i^{ATP}{K}_m^{OAA}}{K_i^{PEP}{K}_l^{ATP}}\frac{\left[ATP\right]\left[PEP\right]}{\left[ADP\right]}+\frac{K_i^{ATP}{K}_m^{OAA}}{K_i^{ADP}{K}_l^{OAA}}\left[OAA\right]}\right)$
PDH (pyruvate dehydrogenase complex)	$\frac{\frac{v_{PDH}^{max}}{\left[NAD\right]}\left(\frac{1}{1+K_i\frac{\left[NADH\right]}{\left[NAD\right]}}\right)\left(\frac{\left[PYR\right]}{K_m^{PYR}}\right)\left(\frac{1}{K_m^{NAD}}\right)\left(\frac{\left[COA\right]}{K_m^{COA}}\right)}{\left(1+\frac{\left[PYR\right]}{K_m^{PYR}}\right)\left(\frac{1}{NAD}+\frac{1}{K_m^{NAD}}+\frac{\left[NADH\right]}{K_m^{NADH}\left[NAD\right]}\right)\left(1+\frac{\left[COA\right]}{K_m^{COA}}+\frac{\left[AcCOA\right]}{K_m^{AcCOA}}\right)}$
Pta (phosphotransacetylase)	$\frac{v_{Pta}^{max}\left(\frac{1}{K_i^{AcCOA}{K}_m^P}\right)\left(\left[AcCoA\right]\left[P\right]-\frac{\left[AcP\right]\left[CoA\right]}{K_{eq}}\right)}{\left(1+\frac{\left[AcCoA\right]}{K_i^{AcCoA}}+\frac{\left[P\right]}{K_i^P}+\frac{\left[ACP\right]}{K_i^{ACP}}+\frac{\left[CoA\right]}{K_i^{CoA}}+\left(\frac{\left[AcCoA\right]\left[P\right]}{K_i^{AcCoA}{K}_m^P}\right)+\left(\frac{\left[AcP\right]\left[CoA\right]}{K_m^{ACP}{K}_i^{CoA}}\right)\right)}$
Ack (acetate kinase)	$\frac{v_{Ack}^{max}\left(\frac{1}{K_m^{ADP}{K}_m^{ACP}}\right)\left(\left[AcP\right]\left[ADP\right]-\frac{\left[ACE\right]\left[ATP\right]}{K_{eq}}\right)}{\left(1+\frac{\left[Acp\right]}{K_m^{AcP}}+\frac{\left[ACE\right]}{K_m^{ACE}}\right)\left(1+\frac{\left[ADP\right]}{K_m^{ADP}}+\frac{\left[ATP\right]}{K_m^{ATP}}\right)}$
Acs (acyl-coenzyme A synthetase)	$\frac{v_{Acs}^{max}\left[ACE\right]\left[NADP\right]}{\left(K_m+\left[ACE\right]\right)\left(K_{eq}+\left[NADP\right]\right)}$
Cs (citrate synthase)	$\frac{v_{CS}^{max}\left[AcCoA\right]\left[OAA\right]}{\left(K_d^{AcCoA}{K}_m^{OAA}+K_m^{AcCoA}\left[OAA\right]\right)+\left(\left[AcCoA\right]K_m^{OAA}\left(1+\frac{\left[NADH\right]}{K_{i1}^{NADH}}\right)\right)+\left(\left[AcCoA\right]\left[OAA\right]\left(1+\frac{\left[NADH\right]}{K_{i2}^{NDAH}}\right)\right)}$
ICDH (isocitrate dehydrogenase)	$\frac{\left[ICDH\right]\frac{K_f}{K_m^{ICIT}}{K}_d^{NADP}\left(\left[ICIT\right]-\frac{\left[NADH\right]\left[2KG\right]}{K_{eq}^{ICDH}\left[NADP\right]}\right)}{\left(\begin{array}{c}\frac{1}{\left[NADP\right]}+\frac{\left[ICIT\right]K_m^{NADP}}{K_m^{ICIT}{K}_d^{NADP}\left[NADP\right]}+\frac{1}{K_d^{NADP}}+\frac{\left[ICIT\right]}{K_m^{ICIT}{K}_d^{NADP}}+\frac{\left[ICIT\right]}{K_d^{ICIT}\left[NADP\right]} \frac{\left[NADPH\right]K_m^{NADP}}{K_m^{ICIT}{K}_m^{NADP}{K}_{einh}^{NADPH}}+\\ \frac{\left[NADPH\right]K_{eknh}^{2KG}}{K_m^{2KG}{K}_{enhe}^{NADPH}\left[NADP\right]}+\frac{\left[2KG\right]K_m^{NADPH}}{K_m^{2KG}{K}_{enhe}^{NADPH}\left[NADP\right]}+\frac{\left[2KG\right]}{K_m^{2KG}}\frac{\left[NADPH\right]}{K_{enhe}^{NADPH}\left[NADP\right]}+\frac{\left[2KG\right]K_m^{NADPH}}{K_m^{2KG}{K}_m^{NADPH}}\frac{\left[NADPH\right]}{K_{ekn}^{NADP}\left[NADP\right]}\end{array}\right)}$
IcL (isocitrate lyase)	$\frac{v_{lcl-f}^{max}\frac{\left[ICIT\right]}{K_m^{ICIT}}}{\left(1+\frac{\left[ICIT\right]}{K_m^{ICIT}}+\frac{\left[SUC\right]}{K_m^{SUC}}+\frac{\left[PEP\right]}{K_m^{PEP}}+\frac{\left[2KG\right]}{K_m^{2KG}}+\frac{1}{K_l}\right)}$
MS (malate synthase)	$\frac{v_{MS}^{max}\frac{\left[GOX\right]}{K_m^{GOX}}\frac{\left[AcCoA\right]}{K_m^{AcCoA}}-v_{MS}^{max}\frac{\left[MAL\right]}{K_m^{MAL}}}{\left(1+\frac{\left[GOX\right]}{K_m^{GOX}}+\frac{\left[MAL\right]}{K_m^{MAL}}+\left(1+\frac{\left[AcCoA\right]}{K_m^{AcCoA}}\right)\right)}$
$\mathrm{aKGDH}$ (akg dehydrogenase)	$\frac{v_{2KGDH}^{max}\left[aKG\right]\left[CoA\right]}{\left(\begin{array}{c}\frac{K_m^{NAD}\left[aKG\right]\left[CoA\right]}{\left[NAD\right]}+K_m^{CoA} \left[aKG\right]+K_m^{2KG}\left[CoA\right]+\left[aKG\right]\left[CoA\right]+\frac{K_m^{2KG}{K}_z\left[aKG\right]\left[SUC\right]\left[NADH\right]}{K_1^{2KG}{K}_1^{SUC}\left[NAD\right]}\\ \frac{K_m^{2KG}{K}_z\left[SUC\right]\left[NADH\right]}{K_1^{SUC}\left[NAD\right]}+\frac{K_m^{NAD}\left[aKG\right]\left[CoA\right]\left[NADH\right]}{K_1^{NADH}\left[NAD\right]}+\frac{K_m^{CoA}\left[aKG\right]\left[SUC\right]}{K_1^{SUC}}\end{array}\right)}$
SDH (succinate dehydrogenase)	$\frac{v_{SDH1}{v}_{SDH2}\left(\left[SUC\right]-\frac{\left[FUM\right]}{K_{eq}}\right)}{K_m^{SUC}{v}_{SDH2}+v_{SDH2}\left[SUC\right]+\frac{v_{SDH1}\left[FUM\right]}{K_{eq}}}$
Fum (fumarase)	$\frac{v_{Fum1}{v}_{Fum2}\left(\left[FUM\right]-\frac{\left[MAL\right]}{k_{Fumeq}}\right)}{K_m^{Fum}{v}_{Fum1}+v_{Fum2}\left[FUM\right]+\frac{V_{Fum1}\left[MAL\right]}{K_{eq}}}$
Mez (malic enzyme)	$\frac{v_{Mez}^{max}\left[MAL\right]\left[NADP\right]}{\left(K_{MAL}+\left[MAL\right]\right)\left(K_{eq}+\left[NADP\right]\right)}$
MDH (malate dehydrogenase)	$\frac{v_{MDH1}{v}_{MDH2}\left(\left[MAL\right]-\frac{\left[OAA\right]}{K_{eq}}\right)}{\left(\begin{array}{c}\frac{K_1^{NAD}{K}_m^{MAL}{v}_{MDH2}}{\left[NAD\right]}+K_m^{MAL} v_{MDH2}+\frac{K_m^{NAD}{v}_{MDH2}\left[MAL\right]}{\left[NAD\right]}+v_{MDH2}\left[MAL\right]+\frac{K_m^{OAA}{v}_{MDH1}\left[NADH\right]}{K_{eq}\left[NAD\right]}+\frac{K_m^{NADH}{v}_{MDH1}\left[OAA\right]}{K_{eq}\left[NAD\right]}+\\ \frac{v_{MDH1}\left[NADH\right]\left[OAA\right]}{K_{eq}\left[NAD\right]}+\frac{v_{MDH1}{K}_m^{OAA}\left[NADH\right]}{K_{eq}{K}_1^{NAD}}+\frac{v_{MDH2}{K}_m^{NAD}\left[MAL\right]\left[OAA\right]}{K_1^{OAA}\left[NAD\right]}+\frac{v_{MDH2}\left[MAL\right]\left[NADH\right]}{K_1^{NAD}}\\ +\frac{v_{MDH1}\left[MAL\right]\left[NADH\right]\left[OAA\right]}{K_{eq}{K}_1^{MAL}\left[NAD\right]}+\frac{v_{MDH2}\left[MAL\right]\left[OAA\right]}{K_{{\rm I}{\rm I}}^{OAA}}+\frac{v_{MDH1}\left[NADH\right]\left[OAA\right]}{K_{{\rm I}{\rm I}}^{NAD}{K}_{eq}}+\frac{K_1^{NAD}{v}_{MDH2}\left[MAL\right]\left[NADH\right]\left[OAA\right]}{K_{{\rm I}{\rm I}}^{NAD}{K}_m^{OAA}{K}_1^{NADH}}\end{array}\right)}$

, respectively.

3.2. Parameter Estimation by (ESe-PSO) Algorithm

In this study, an enhanced (Se-PSO) algorithm was proposed. This method shares the same idea as Se-PSO, but it is different in terms of exploration and exploitation of the inertia weight [27]. Decreasing linear inertia weight ($\omega =0.9to0.4$) [27] was used. The inertia weight controlled the effect which the last iteration speed had on the current speed, and thus allowed the particles to explore wider areas in the beginning and nearby areas in later stages with reduced speed for exploitation. Here, the main plan of ESe-PSO was to add a damping process to the inertia weight so as to support the exploration and exploitation of the particles towards the optimum solution. The development was done by initializing the inertia weight $\omega =1$ and damping parameter called $\omega damp=0.99$. This damp parameter was calculated in the iteration process to increase the controlled convergence and to support the particles for reaching the global solution. In

Table 1. Pseudo-code of Enhanced Segment Particle Swarm Optimization Algorithm (ESe-PSO).

Algorithm ESe-PSO Adoption

1. Begin

2. Initialize B, S, D, ω, and ωdamp

3. Initialize ${\mathrm{v}}_i, X_i, c_1,c_2, r_1, r_2$, number_segment

4. Segment length = initial value/number_segment

5. Adopt the $b_{e1}:b_{en}$ parameters boundaries with respect to D

6. For j = 1 to the number of the segment

7. Determine initial fitness for segment J

8. Assume Best fitness = initial fitness

9. End for

10. For m data account

11. $fitness={\displaystyle \sum }_{i=1}^R\left|\frac{y_{si}-y_{ei}}{y_{ei}}\right|$

12. End for

13. If $fitness>\mathrm{Best}fitness$

14. For each S

15. iter S = 1, S++

16. Updating the velocity $V_i$ towards fitness:

17. $v_i\left(t+1, j\right)=\left(\omega \ast \omega damp\right)v_i\left(t,j\right)+c_1{r}_1\left(p_i\left(t,j\right)-x_i\left(t,j\right)\right)c_2{r}_2\left(G\left(t,j\right)-x\left(t,j\right)\right)$

18. Update the position $X_{\mathrm{i}}j$ towards fitness;

19. $X_i\left(t+1, j\right)=X_i\left(t, j\right)+v_i\left(t+1, j\right)$

20. End for

21. End if

22. If $fitness\le \mathrm{Best}fitness$

23. Print $G_i$ best of each particles

24. Update the $b_{e1}:b_{en}$ based of $G_i$ best of each kinetic parameters

25. If $fitness>\mathrm{Best}fitness$ return step 2 till the iteration is finished or discover high-quality solution

26. End if

27. End if

28. Global point(j) = $G_i$ best

29. Next j

30. ${\mathrm{Optimal}}_{\mathrm{segment}}=\max \left(\mathrm{Global}\mathrm{point}\right)\pm \mathrm{segment}\mathrm{length}/2$

31. Repeat algorithm 1 for the new initial values

below, the (ESe-PSO) algorithm is described:

Table 1 describes the steps of the enhanced segment particle swarm optimization algorithm used to estimate the kinetic parameters. These algorithms were started by initializing the particles that had the highest number of generation $S=30$. The bird number $B=10$ and $b_{e1}:b_{en}$ was the kinetic parameter boundaries $b$ from $e1toen$ parameters. The problem dimension $D=7$, inertia weight $\omega =1$ and damping process $\omega damp=0.99$ were in step 2. Then step 3 initialized the velocity and position of particles with the parameters of PSO which are acceleration coefficients towards the best personal position $p_i$ and the global best position of the entire $G_i,$ respectively, where $c_1\&c_2$ are constants called cognitive and social scaling parameters, respectively, and a random number between 0 and 1 for ($r_1$, $r_2$).

The calculation of the segment length was obtained by utilizing step 4. Each kinetic parameter of the boundaries, both of the upper and lower, regarding the dimension of the problem D was created in step 5. Based on the number of segments, the calculation of the initial fitness of the whole segments was done with the assumption that the initial fitness was the best fitness step for 6, 7 and 8.

Here, the new fitness was calculated based on the fitness equation stated in step 11. Where the fitness was greater than the best fitness set, the iteration updated the velocity and the position towards the new fitness steps of 15, 17 and 19. If the new fitness was smaller than or equal to the best fitness print, the global best particles position updated the kinetic parameter, based on the global best particles position which is supposed to be the new initial kinetic parameters step in 22 and 23. If the fitness was greater than the best fitness, it returned to step 2 until the program highly figured out a solution or the iteration was met in step 25. The global point of segments is set equal to the global best position which was then used to calculate the optimal segment by its Equations (28)–(30). After the optimum segments were identified, they were used as the new initial kinetic parameter; PSO algorithm in [13] was searched among them and ended the program in step 31.

4. Results

The estimation of large-scale kinetic parameters in the kinetic model was a difficult task due to the nonlinearity of the model. For this reason, the kinetic parameters were reported from different laboratories and time consumption [32]. It was stated that during the application of the local sensitivity analysis to [14] model result, each kinetic parameter increased up to 200% by step 0.5 to calculate its sensitivity. It was found out that 7 kinetic parameters were highly affected by the used model’s response for estimating the proposed algorithm. These sensitive kinetic parameters are $v_{max}^{pyk}$,$n_{pk}$,$icdh$,$k_{icdh}^f$,$k_{icdhnap}^d$,$k_{icdhnap}^m$ and $v_{max}^{icl}$ which were involved in this reaction rate of $V_{pyk}, V_{icdh}, and V_{icl}$. In addition, the experimental data of [16] were taken. These experimental data are $Glc$, $G6P$, $F6P$, $FDP$, $PEP$, $PYR$, $6PG$, $Ru5P$, $Xu5P$, $S7P$, $R5P$ and $E4P$.

As a matter of fact, during the ESe-PSO algorithm execution, the segments of the kinetic parameters of [32] were increased by adding 1 segment to each kinetic to increase the possibility of searching the wide-area and finding an accurate solution which is depicted in Appendix A

Table A3. Kinetic parameters segmentation.

Kinetics	Number of Segments
$v_{max}^{pyk}$	2
$n_{pk}$	3
$icdh$	2
$k_{icdh}^f$	2
$k_{icdhnap}^d$	3
$k_{icdhnap}^m$	2
$v_{max}^{icl}$	3

.

The ESe-PSO algorithm parameters such as $c_1 and c_2$ was initialized as mentioned in [32], $\omega =1$ and $\omega damp=0.99$. In this regard, the kinetic boundaries with their upper and lower bounds values were initialized with small increase based on [32] to allow the proposed algorithm search large space. Here, the mechanism for selecting the boundaries was based on the original kinetic parameters values, which assumed that these values are $\pm$ with each kinetics’ small value, so that the optimum result may be around it. These kinetic parameters boundaries are depicted in

Table 2. Kinetic parameters boundaries.

Kinetics	Original	Lower Bound	Upper Bound	Estimated Kinetics
$v_{max}^{pyk}$	1.0850	0.8500	1.5000	1.10000
$n_{pk}$	3.0000	2.3000	3.5000	2.83000
$icdh$	24.421	23.500	24.9000	24.1020
$k_{icdh}^f$	289,800	289,799	289,800	289,799.5
$k_{icdhnadp}^d$	0.0060	0.0030	0.0500	0.00920
$k_{icdhnadp}^m$	0.0170	0.0060	0.0600	0.01920
$v_{max}^{icl}$	3.8315	3.300	4.3000	3.64000

with their estimation as follows:

After the kinetic parameters were initialized, the kinetic parameters were segmented to find the optimal segment for each kinetics. Later on, the updating of position, velocity and the fitness function was done based on the initialization of ESe-PSO parameters such as the inertia weight, damping process and the acceleration coefficient ($c_1 and c_2$). Thus, the optimal segment was determined based on the fitness function. The optimal segments were sent to the PSO algorithm which provided a local best solution so that PSO would search around the optimal segments.

In this regard, the ESe-PSO algorithm estimated the kinetic parameters and minimized the model distance to 28.94%. The model simulation was moved closer to the real experimental data when compared to the considered model; the Se-PSO, GA and DE algorithms are described in

Table 3. Comparison of Enhanced Segment Particle Swarm Optimization algorithm (ESe-PSO) algorithm estimations with Segment Particle Swarm Optimization algorithm (Se-PSO), Differential Evolution (DE) and Genetic Algorithm (GA) algorithms.

Metabolites	Experimental Data	Kadir Data	Se-PSO Estimation	ESe-PSO Estimation	DE Estimation	GA Estimation
$Glc$	0.0556	0.12203	0.1023	0.0923	0.10250	0.11740
$G6P$	3.4800	0.12989	0.3120	1.3281	0.19930	0.14120
$F6P$	0.6000	0.02146	0.0532	0.0872	0.02350	0.02210
$FDP$	0.2720	1.51860	2.5257	1.3335	2.01250	1.97500
$PEP$	2.6700	1.50760	1.8320	2.0321	1.84700	1.98320
$PYR$	2.6700	2.82790	2.9380	2.7256	3.01000	3.00200
$6PG$	0.8080	0.01785	0.0752	0.2043	0.01947	0.01874
$Ru5P$	0.1110	0.02140	0.0628	0.0870	0.02324	0.02280
$Xu5P$	0.1380	0.02652	0.0925	0.0902	0.02847	0.02794
$S7P$	0.2760	0.00473	0.0524	0.0802	0.04827	0.04803
$R5P$	0.3980	0.07639	0.0967	0.0985	0.08425	0.07727
$E4P$	0.0980	0.27837	0.0568	0.0671	0.02835	0.02428

. The GA parameters were set as follows: generation gap = 0.8, crossover = 0.85, mutation = 0.05 [33] and maximum generation = 500. The DE parameters were set as follows: scaling factor $F=0.7$, crossover rate $Cr=0.8$, search strategy was DE/rand/1/bin as stated in [6], population size = 100 and a number of generation = 500 [32].

Table 3 shows how the proposed method moved the model simulation towards the experimental data closely. This simulation increased $PEP$ metabolites and this may be due to the other pathways engaged during the execution of the algorithm such as anaplerotic, TCA, with acetate formation. In addition, this increase may be the $OAA$,$PEP$ and $PYR$ metabolites that were affecting the cell’s growth [10].

The analysis of the model distance minimization showed that 12 out of the 12 datasets were moving towards the real experimental data. These metabolites are $Glc$, $G6P$, $F6P$, $FDP$, $PEP$, $PYR$, $6PG$, $Ru5P$, $Xu5P$, $R5P$, $S7P$ and $E4P;$ thus, the metabolites were still slightly far from the experiment, perhaps because of the model complexity, shortage of data and the portion of some metabolites. Furthermore, the G6P empirical data were 3.48 $\mathrm{mM}$ and the simulation outcome of G6P could be leading to the empirical data rather than the G6P result model under examination when related to the G6P model. This discrepancy may be because of the model resolution of G6P consistency in the presence of $5mM$ $MgS{O}_4$ and $0.48mM$ $NAD{P}^{+}$ after the addition of $0.7Um{l}^{-1}$ of $G6PDH$ as reported by [10].

Moreover, the $FDP$ model simulation was having a better movement towards the experimental data rather than the considered model. This may be due to the nonlinearity of the system or the lumping of $GAP$ and $DHAP$ metabolites in one equation. The changes in $FDP$ may be due to the lumping of $pyk\mathrm{II}$ and $pyk\mathrm{I}$ as stated in [10].

As stated in Table 3 above, the simulation response of the distance model minimization was moving towards the experimental data in 12 metabolites. Moreover, the model under investigation had five pathways than the data of only two pathways where only one metabolite was not well minimized ($PEP$) due to the lack of experimental data. The estimated kinetic parameters were used in the model simulation to compare our model responses to [10] model responses. The result of the comparison is depicted in Figure 2 below:

In order to prove that the result is statistically consistent, the STD, mean and F value were calculated for the considered model. The result of the simulated model is displayed in Table A3 and Equations (3)–(6).

In Appendix A

Table A4. STD and mean of the simulated results.

Metabolites	Experimental Data	The Model under Study Data	Se-PSO Estimation	ESe-PSO Estimation	DE Estimation	GA Estimation
$Glc$	0.0556	0.1220	0.1023	0.0923	0.1025	0.1174
$G6P$	3.48	0.1298	0.3120	1.3281	0.1993	0.1412
$F6P$	0.6	0.0214	0.0532	0.0872	0.0235	0.0221
$FDP$	0.272	1.5186	2.5257	1.3335	2.0125	1.975
$PEP$	2.67	1.5076	1.8320	2.0321	1.8470	1.9832
$PYR$	2.67	2.8279	2.9380	2.7256	3.0100	3.002
$6PG$	0.808	0.0178	0.0752	0.2043	0.0194	0.0187
$Ru5P$	0.111	0.0214	0.0628	0.0870	0.0232	0.0228
$Xu5P$	0.138	0.0265	0.0925	0.0902	0.0284	0.0279
$S7P$	0.276	0.0047	0.0524	0.0802	0.0482	0.0480
$R5P$	0.398	0.0763	0.0967	0.0985	0.0842	0.0772
$E4P$	0.098	0.2783	0.0568	0.0671	0.0283	0.0242
MEAN	0.9647	0.5251	0.6833	0.6855	0.6189	0.6216
STD	1.2271	0.9102	1.0832	0.9320	1.0439	1.0554
Variance	1.5059	0.8286	1.1735	0.8686	1.0898	1.1138
N	12	12	12	12	12	12
Distance model minimization		45.56%	29.16%	28.94%	35.58%	35.55%
F test			1.0497	1.0482	1.2881	1.3165
Critical value of F			0.3678	0.3678	0.3678	0.3678

Table A4 describes the STD, MEAN, Variance and F test of the simulated results.

, the mean and STD were calculated for the proposed algorithm results and the response of the considered model’s result. The proposed algorithm showed 0.6855 mean compared to the considered model (0.5251), DE (0.6189) and GA (0.6216) which indicated that the proposed algorithm result moved to the experimental data mean of 0.9647. In addition, DE algorithms showed a better result than GA in the mean values.

In addition, the proposed algorithm showed 0.9320 STD compared to the considered model’s 0.9198 STD which indicated that the proposed algorithm result moved to the experimental data STD of 1.2271. Thus, the hypothesis of the result in

Table 4. Objective function of each algorithm over 10 runs.

Methods	Mean	STD	Best	Lower	Time-Consuming
ESe-PSO	7.035 ×10−5	2.72 × 10−4	2.3 × 10−6	1.58 × 10−4	19 h
Se-PSO	6.03 × 10−4	3.567 × 10−3	2.24 × 10−5	3.52 × 10−3	19 h
DE	4.9185 × 10−2	2.51213 × 10−1	8.27 × 10−3	1.7643 × 10−1	26 h
GA	1.1476 × 10−1	2.58243 × 10−1	2.579 × 10−3	1.4673 × 10−1	29 h

is calculated as follows:

The hypothesis of Table A4 is calculated and decided below:

$$H_0:ST{D}_E^2\ge ST{D}_D^2$$

(3)

$$H_1:ST{D}_E^2<ST{D}_D^2$$

(4)

$$F_{test}=\frac{ST{D}_E^2}{ST{D}_D^2}=\frac{0.8686}{0.8286}=1.0482$$

(5)

$$F_{1-0.05}=\frac{1}{F_{0.05,n_E,n_D}}=\frac{1}{F_{0.05,11,11}}=\frac{1}{2.7186}=0.3678$$

(6)

where $ST{D}_E^2$ is the standard deviation of the optimized result $E$, $ST{D}_D^2$ is the standard deviation of the considered model $D$ and the $n_E,n_D$ are the number of variables for the optimized and model result. This hypothesis was locking for minimizing the considered model’s errors. This suggests that we accepted $H_0$ and rejected $H_1$ as the accepted result. Thus, the F test showed that the proposed method had moved to the critical F value better than Se-PSO, DE and GA by indicating1.0482, 1.0497, 1.2881 and 1.3156, respectively.

The objective function and time consumption were calculated for the proposed algorithm, DE and GA as depicted in Table 4 below:

Table 4 demonstrates a comparison between the proposed algorithm with DE and GA algorithms. Table indicates that the proposed algorithm in this study produced superior results in terms of objective function values over 10 runs. The mean, STD, best and lowest objective functions were calculated. ESe-PSO algorithm produced the best mean (7.035E–05) and STD (0.000272) when compared to Se-PSO (mean 0.000603; STD 0.003567), DE (mean 0.049185; STD 0.251213) and GA, (mean 0.11476; STD 0.258243). The proposed method achieved the best objective function (2.3E-06) while Se-PSO ((0.00352), DE (0.00827) and GA (0.002579) achieved the best objective function. The consumed time by the proposed algorithm was achieved in 19 h while Se-PSO, DE and GA achieved the result in 19, 26 and 29 h, respectively.

5. Conclusions

The proposed method of ESe-PSO was used in order to estimate large scale kinetic parameters based on a large kinetic model. The estimation was carried by minimizing the model response distances with real experimental data. Seven kinetic parameters were estimated using 12 metabolites and the method optimized the 12 metabolites. These metabolites are $Glc$, $G6P$, $F6P$, $FDP$, $PEP$, $PYR$, $6PG$, $Ru5P$, $Xu5P$, $R5P$ $S7P$ and $E4P$. The metabolites were well minimized with a slightly small error which may be due to the lump of $GAP/DAHP$ metabolites or $tktb$ involvement in the pentose phosphate and glycolysis pathways. However, the proposed method of ESe-PSO algorithm revealed that it had the ability to perform kinetic parameters estimation. It is hereby suggested that further study can be carried out on the evaluation of large scale kinetic parameters segmentation of other algorithms like Firefly and Bat algorithms.