Novel column generation-based optimization approach for poly-pathway kinetic model applied to CHO cell culture

Abstract

Mathematical modelling can provide precious tools for bioprocess simulation, prediction, control and optimization of mammalian cell-based cultures. In this paper we present a novel method to generate kinetic models of such cultures, rendering complex metabolic networks in a poly-pathway kinetic model. The model is based on subsets of elementary flux modes (EFMs) to generate macro-reactions. Thanks to our column generation-based optimization algorithm, the experimental data are used to identify the EFMs, which are relevant to the data. Here the systematic enumeration of all the EFMs is eliminated and a network including a large number of reactions can be considered. In particular, the poly-pathway model can simulate multiple metabolic behaviors in response to changes in the culture conditions. We apply the method to a network of 126 metabolic reactions describing cultures of antibody-producing Chinese hamster ovary cells, and generate a poly-pathway model that simulates multiple experimental conditions obtained in response to variations in amino acid availability. A good fit between simulated and experimental data is obtained, rendering the variations in the growth, product, and metabolite uptake/secretion rates. The intracellular reaction fluxes simulated by the model are explored, linking variations in metabolic behavior to adaptations of the intracellular metabolism.

Keywords: Amino acid; Chinese hamster ovary cell; Column generation; Elementary flux mode; Kinetic modelling; Metabolic flux analysis; Optimization; Poly-pathway model.

Fig. 1

CG algorithm. Flow chart showing the CG-based identification of an EFM subset using the CG algorithm.

Fig. 2

Metabolic network map. The subscript m denotes mitochondrial metabolites and these metabolites are shown in red. The subscript ext denotes extracellular metabolites. The boundaries between compartments are indicated by dashed lines. Subsystems are indicated by color. The visualization is created with VANTED (Rohn et al., 2012b) (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 3

Determination of the network of macro-reactions by CG. The CG-based algorithm (CGa, yellow), replaces the conventional enumeration of EFMs. The identification of each EFM subset (pink) is carried out independently for each experimental condition. The subsets are then combined into the EFM union.

Fig. 4

Determination of the macro-reactions and kinetics - Step 1. The macro-reactions of the metabolic reaction network are created by an EFM approach: given the metabolic network, the experimental data and the bounds, the ${\mathbf{E}}_{\mathit{union}}$ is generated using the CG algorithm (see Fig. 3). One or several potential kinetic equations with saturation and/or inhibition effects are attributed for each EFM in ${\mathbf{E}}_{\mathit{union}}$, noted as ‘All the potential kinetics equations associated to each macro-reaction’. This generates a large potential model, LargeM, composed of the macro-reactions and several kinetics alternative per macro-reaction. The potential kinetic equations vary in the effect of the metabolites. The saturation parameters ${\mathbf{K}}_{\mathbf{s}}$ of the kinetics are taken from the literature while the inhibition parameters ${\mathbf{K}}_{\mathbf{p}}$ and ${\mathbf{K}}_{\mathbf{r}}$ are obtained by fitting the fluxes estimated by the model to the experimental data in LargeM.

Fig. 5

Identification of the reduced model - Step 2. A cross-validation approach is applied to identify the final model. The data are randomly distributed in a training data set and a testing data set. For each training data set d, the macro-reaction network defined by ${\mathbf{E}}_{\mathit{union}}$ and the kinetics of Step 1 are used to determine a large model LargeM^d. The LargeM^d is reduced into RedM^d. RedM^d is then used on the testing data set to simulate the flux rates and the error between the simulated and measured data is computed. From all the repetition exercises, the RedM^d providing the smallest simulation error for the testing data is selected as final model RedM_final.

Fig. 6

Comparison between EFMs in subsets and union. A comparison between EFMs in the sixteen EFM subsets (y-axis) and in the union ${\mathbf{E}}_{\mathit{union}}$ (x-axis) with green color indicating that the EFM is found in both the EFM subset and the union. In the figure, the 125 unique EFMs in ${\mathbf{E}}_{\mathit{union}}$ have been sorted in descending order according to the number of matched EFM subsets, and occurs along the x-axis in the following order: 4, 8, 16, 18, 26, 27, 1, 3, 22, 30, 31, 2, 9, 17, 20, 28, 5, 15, 32, 53, 10, 29, 47, 36, 13, 19, 44, 65, 6, 11, 24, 33, 42, 50, 55, 97, 100, 25, 34, 41, 43, 46, 54, 56, 62, 66, 67, 69, 85, 89, 91, 7, 12, 14, 21, 23, 35, 37, 38, 39, 40, 45, 48, 49, 51, 52, 57, 58, 59, 60, 61, 63, 64, 68, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 86, 87, 88, 90, 92, 93, 94, 95, 96, 98, 99, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125.

Fig. 7

The x-axis shows the logarithm (base 10) of the cut-off level $w_{\mathit{cut}}$. The left y-axis shows the error (2-norm) indicating the fit between model and data. The right y-axis shows the percentage of equations in LargeM_initial that are kept. $w_{\mathit{cut}}=0.1$ is chosen for the model reduction ($\mathit{log}10(0.1)=-1$).

Fig. 8

Comparison of the simulated fluxes of the training set (yellow triangle) and the test set (green squares) with the experimental data (blue circle) of metabolic rates in pmol/cell, day (y-axis), for biomass, mAb, glucose, lactate, glutamine, ammonia and alanine. For each experimental condition (A0, A200, N0, …, Ctrl), the data and simulated values are shown along the x-axis in ascending order according to time (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 9

Comparison of the simulated fluxes of the training set (yellow triangle) and the test set (green squares) with the experimental data (blue circle) of metabolic rates in pmol/cell, day (y-axis), for glutamate, aspartate, asparagine, cysteine, proline, serine and glycine. For each experimental condition (A0, A200, N0, …, Ctrl), the data and simulated values are shown along the x-axis in ascending order according to time (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 10

Comparison of the simulated fluxes of the training set (yellow triangle) and the test set (green squares) with the experimental data (blue circle) of metabolic rates in pmol/cell, day (y-axis), for isoleucine, leucine, phenylalanine, threonine, tryptophan, tyrosine, valine, methionine and lysine. For each experimental condition (A0, A200, N0, …, Ctrl), the data and simulated values are shown along the x-axis in ascending order according to time (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 11

Flux distributions over the metabolic network in the control medium, subjected to glutamine omission (Q0), asparagine omission (Q0) and serine omission (S0). The thickness of the arrows is proportional to the reaction flux. Dashed lines represent fluxes for which $|v_i|<$ 0.0025. Fluxes that are significantly different compared to the control are colored in the following way: blue indicates an increase, red indicates a decrease, and purple indicates a shift in the direction of the flux. The flux visualizations are created with VANTED (Rohn et al., 2012b) and FluxMap (Rohn et al., 2012a).