EP3251039A1 - Computerimplementiertes verfahren zur erstellung eines fermentationsmodels - Google Patents

Computerimplementiertes verfahren zur erstellung eines fermentationsmodels

Info

Publication number
EP3251039A1
EP3251039A1 EP16701791.2A EP16701791A EP3251039A1 EP 3251039 A1 EP3251039 A1 EP 3251039A1 EP 16701791 A EP16701791 A EP 16701791A EP 3251039 A1 EP3251039 A1 EP 3251039A1
Authority
EP
European Patent Office
Prior art keywords
reactions
macro
rates
computer
implemented method
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Withdrawn
Application number
EP16701791.2A
Other languages
English (en)
French (fr)
Inventor
Tobias NEYMANN
Lukas Hebing
Sebastian Engell
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Bayer AG
Original Assignee
Bayer AG
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Bayer AG filed Critical Bayer AG
Publication of EP3251039A1 publication Critical patent/EP3251039A1/de
Withdrawn legal-status Critical Current

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B5/00—ICT specially adapted for modelling or simulations in systems biology, e.g. gene-regulatory networks, protein interaction networks or metabolic networks
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B5/00—ICT specially adapted for modelling or simulations in systems biology, e.g. gene-regulatory networks, protein interaction networks or metabolic networks
    • G16B5/30—Dynamic-time models

Definitions

  • the invention relates to a computer-implemented method for creating a model of a bioreaction - esp. Fermentation or whole-cell catalysis - using an organism.
  • Organism in the sense of the application are cultures of plant or animal cells such as mammalian cells, yeasts, bacteria, algae, etc., which are used in bioreactions.
  • background knowledge of the organism means knowledge about the biochemical reactions of the organism - specific and unspecific reactions - and in particular the cell-internal reactions, or macro-reactions describing the organism-specific metabolic networks (also called SN or metabolic networks) consisting of substrates, metabolites (also called nodes of the metabolic network), products as well as the biochemical reactions between them. These biochemical reactions are through their:
  • the measurement data is mainly used for qualitative monitoring of the process.
  • the following section presents a selection of technical issues that require dynamic process models to solve.
  • a technical use of the process knowledge in the sense of the application provides the model-based state estimation of a process in a bioreactor.
  • Methods such as the extended Kalman filter allow a continuous estimation of process variables over which discontinuous measurements are available [Welch G, Bishop G. 1995. Chapel Hill, NC, USA: University of North Carolina at Chapel Hill The course of non-measurable quantities can also be calculated from other measurements, provided that a process model is correctly described for the underlying process.
  • a dynamic process model is used to optimize process control in terms of product quantity, product characteristics or formation of by-products or other target variables in a model-based, predictive closed-loop control system.
  • Frahm et al. for a hybridoma cell culture [Fever B, Lane P, Atzert H, Munack A, Hoffmann M, Hass VC, Portner R. 2002. Adaptive, Model-Based Control by the Open Loop Feedback Optimal (OLFO) Controller for the Effective Fed-Batch Cultivation of Hybridoma Cells. Biotechnol. Prog. 18 (5): 1095-1103J.
  • Typical product features within the meaning of the application are, for example, glycosylation patterns of proteins or protein integrity, without being limited thereto. Dynamic models used in the above context do not yet have this property. The present approach enables a simple model-based integration of product features.
  • This method is useful for various organisms or strains of the same organism.
  • the selection of macro-reactions for the process model is done using process knowledge. However, process sections are defined for which a predefined number of macro-reactions are selected separately at random.
  • the method described here provides one of many possible combinations of elementary modes. The number of macro reactions, and thus the model complexity, is fixed and unchangeable.
  • the method gives separate models for each process section.
  • the kinetics of the individual macro reactions are selected taking into account the stoichiometry of the selected macro reactions.
  • the parameters of the kinetics (model parameters) are not adapted to the process data. Instead, the use of the separate process section models maps the changes in the process data.
  • the object has been achieved by a method for creating a model of bioreaction with an organism in a bioreactor as described below.
  • the subject of the application is a computer-implemented method for creating a model of a bioreaction - in particular fermentation or whole-cell catalysis - with an organism comprising the following steps:
  • the EMs are summarized in a matrix K, where the EMs summarize the metabolic pathways from a) into macro reactions.
  • This matrix K contains the stoichiometry and the reversibility properties of all possible macro reactions from the background.
  • the measurement data (also called process knowledge) for bioreaction with the organism are entered.
  • the rates specific to the organism - excretion and uptake rates of one or more input variables and output variables - of the input metabolic pathways are calculated.
  • Growth rates, particularly preferably also mortality rates of the organism, are preferably also calculated.
  • the subsets are displayed graphically.
  • reaction rates of the macro reactions of the subset r (t) are calculated on the basis of the input measured data from c) and / or the rates from d).
  • a first adaptation of the model parameter values for each macro reaction is carried out separately to the calculated reaction rates from f) and a check of the adaptation quality.
  • steps g) and h) are repeated until a predefined quality of adaptation is achieved.
  • the model parameter values are adapted to the measurement data from c).
  • the matrix L, the kinetics from g) and the model parameter values from j) form the model and are output and / or transferred to a process control or process development module.
  • the process control module communicates on-line with a process control system commonly used to control the bioreactor.
  • process development modules are used to off-line optimize the process or design experiments.
  • the bioreaction modeling according to the invention is essentially based on the assumption of representative macro reactions which simplify internal metabolic processes.
  • the selection of reactions requires both biochemical background knowledge and process knowledge.
  • the reactions of the metabolic network, their stoichiometry and reversibility property are entered by the user via a user interface or, ideally, automatically by the selection of an organism and its stored metabolic pathways from a database module in which the background information on the organism is stored.
  • the metabolic network also called stoichiometric network in the state of the art
  • the metabolic network preferably contains reactions from metabolic pathways which are important for the organism, for example reactions of glycolysis.
  • the selection contains external reactions.
  • An external reaction in the sense of the application contains at least one component outside the cell, typically at least one input variable and / or at least one output variable (product, by-product, etc.). More preferably, the metabolic network contains reactions that describe cell growth, e.g. B. in the form of a simplified reaction of internal metabolites to external biomass.
  • Figure 5 and Table 1 in the example describe but are not limited to an applicable metabolic network.
  • each elementary mode is a linear combination of reaction rates from the metabolic pathways - d. H. internal and external responses of the metabolic network, which satisfy both the steady state condition for internal metabolites and the reversibility or irreversibility of reactions, in linear combinations of reactions that take into account the steady state condition of internal metabolites , no internal metabolites can accumulate.
  • a macro reaction in the sense of the application summarizes all reactions that lead from one or more input variables to one or more output variables (n) Each elementary mode thus describes a macro reaction Compared to the method of Leifheit et al the macro reactions are determined on the basis of the entered background knowledge.
  • the elementary modes (EMs) are combined in a matrix E, preferably in a module for matrix engineering, which is configured with a corresponding algorithm.
  • EMs The elementary modes
  • Known algorithms can be used outside the Elementary Modes Matrix.
  • METATOOL is mentioned as an example without being limited to: [Pfeiffer T, Montero F, Schuster, 1999. METATOOL: for studying metabolic networks. Bioinformatics 15 (3): 251-257.]
  • METATOOL generates a first matrix E describing the input internal and external responses.
  • step b) a matrix consisting of possible macro reactions K is generated with the aid of the (external) stoichiometric matrix N p from the matrix (E).
  • the column vectors of the matrix K describe the macro reactions.
  • the row vectors describe the components of the macro reactions (input and output variables).
  • the stoichiometry of the macro reactions is entered.
  • the available measurement data (process knowledge) for bioreaction with the organism are entered.
  • cell number, cell vitality, concentrations of substrates such as carbon sources (eg, glucose), amino acids or O 2 , products and by-products (eg, lactate or CO 2 ), process parameters such as temperature and / or pH, or product characteristics determined.
  • substrates such as carbon sources (eg, glucose), amino acids or O 2 , products and by-products (eg, lactate or CO 2 ), process parameters such as temperature and / or pH, or product characteristics determined.
  • process parameters such as temperature and / or pH, or product characteristics determined.
  • step d) the cell-specific excretion and uptake rates of substrates and (ancillary) products - together specific rates q (t) - called, and optionally the
  • ⁇ ( ⁇ ), M d CO Growth and death rates of the organism ( ⁇ ( ⁇ ), M d CO) are calculated. Prerequisite for the calculation is the interpolation of the vital cell count, the total cell count and media concentrations with the help of a From these temporal changes of the measured variables can be determined.
  • the calculated rates q (t), ⁇ (, give information about the observed dynamic behavior of the organism over time.
  • One or more different methods of interpolation can be used in combination to calculate the above rates.
  • Leifheit et al. the determination of the temporal changes of measured variables - z.
  • B. the total cell count, the vital cell count or from other media concentrations of measurement data using spline-interpolated measurement data [Leifheit, J., Heine, T., Kample, M., & King, R. (2007).
  • Computer-aided semi-automatic modeling of biotechnological processes (Semiautomatic Modeling of Biotechnical Processes). at-automation technology, 55 (5), 21 1 -218]. This method is hereby incorporated into the application by reference.
  • the growth rate of the organism ⁇ ( ⁇ ) can be calculated from spline-interpolated values of the total cell number X t (t) and the living cell number X v (t) as well as the temporal values calculated from them
  • D (t) is the dilution rate
  • the rate of death ⁇ ⁇ (£) can be determined from the course of X v (t) and the
  • the calculation of the specific rates of another component tq j (t) can be obtained from spline-interpolated values of the living cell number X v (t) and the concentration of the component Q (t), as well as the variation of the temporal change ⁇ (t) from 8 ⁇ 1 ⁇ 6 - ⁇ 6 ⁇ 1 ⁇ 6 ⁇ values of Q (t) can be determined using formula 4:
  • the measurement data from step c) are prepared before the first interpolation as follows:
  • the measurement data are shifted (in the application called shifts).
  • the amount AQ (t), by which the concentration measurement is shifted, can be calculated according to formula 5:
  • This treatment (shifts) of the measured data prevents a sudden change in the calculated specific rates when turning on or off a feed (feed peak), in particular in a fed-batch process.
  • Figure 1 shows the processing / shifting of the measured data in the sense of this application.
  • the processed data then becomes
  • the lysis is included in the calculation on the basis of a lysis factor K t (formula 9). This can z. B. be assumed to be constant over the course of the process.
  • the processed data are also used for the calculation of the death rate ⁇ ⁇ (t).
  • ICell hi own When the composition of this macromolecule is estimated, e.g. Based on its C-mol content, its specific rate can be varied from [g] to [C-mol], so that the specific rate has the unit [ ⁇ e ⁇ ].
  • the specific rates q (t) form one of the bases for the further step e) of the method, namely the selection of the relevant macro reactions.
  • a subset (L) of the EMs is selected on the basis of the data with which the specific rates q (t) from d) and / or the measurement data from c) can be well mapped according to a mathematical quality criterion.
  • the number of EMs in the subset (L) should be as small as possible in order to ensure the lowest possible complexity of the process model.
  • the subset L should however ensure a good description of the process knowledge.
  • the selection of EMs reduces the solution space compared to the original EMs set (K) from a), but still contains the determined physiologically important area of the cells.
  • Figure 3 shows a representation of the solution space, where the original set of EMs (K) is reduced to a subset (L).
  • step e the calculated specific rates q (t) and the measurement data from c) are usually transferred to a module for selecting the relevant macro reactions, which is configured with corresponding algorithms.
  • step e) i) a data-independent and / or a data-dependent prereduction of the matrix K takes place in any order:
  • the data-independent prereduction is preferably carried out by a geometric reduction.
  • all cosine similarities to all other modes are calculated for a randomly selected EM.
  • the most similar EM will be removed from the set. This procedure is repeated until a predefined number of EMs has been reached.
  • the desired number is usually defined for the procedure in advance.
  • the volume of the solution space can be used. Surprisingly, it was found that a significant reduction in the number of macroreactions while maintaining 90 to 98%, preferably 92 to 95% of the clamped volume, compared to the original volume is possible.
  • the data-dependent prereduction can be calculated by comparing yield coefficients of the EMs (Y EM ) to the yield coefficients calculated from the specific rates q (t) of d)
  • the yield coefficient of the kth EM (i) ' k ) is given by formula 10 Division of the corresponding stoichiometric coefficients of the external metabolites i and j determined. For the kth EM these are the matrix entries K ik and K Jik .
  • the yield coefficient Y j (t) gives the ratio of two measured or d) calculated cell-specific rates (i (t), q j (t)) to each other according to formula 1 1:
  • an upper and a lower limit can be determined for each possible combination of two external components t and j.
  • the lowest yield coefficient of two external metabolites i and j ⁇ TM ( ⁇ ) can be used as the lower bound and the largest value of ⁇ TM ( ⁇ ) as the upper bound, but other limits are possible.
  • EMs whose yield coefficients Y TM are above the upper limit or below the lower limit are removed from the matrix K. If the yield coefficient of an EMs Y TM can not be determined, this remains in the matrix K.
  • step e) ii. a subset of macro reactions is selected with an algorithm: For the selection, a quality criterion with which it can be quantified how well the specific rates q (t) from d) and / or the measurement data from c) are mapped with a subset (L) and an algorithm for selecting the subset.
  • the value for SSR q should be as small as possible.
  • the vector r (tj) has to be determined beforehand for each considered time with the aid of a non-negative-least-squares algorithm such that the following applies:
  • the advantage of this method is that the calculations according to formula 12 - 14 can also be carried out for very large subsets with many EMs.
  • a significant disadvantage is that the computed specific rates q (t) are required for this calculation. Since these are obtained from interpolated measured values, they are present with great uncertainty regarding their true values. Measurement inaccuracies may under certain circumstances have a strong effect on the calculated specific rates c / (t). The quality criterion SSR q can therefore only be determined under great uncertainty.
  • this method also yields an estimated course of the reaction rates r (t) of the subset L as a result of the minimization according to formulas 13 and 14.
  • Leighty, R. et al. describes another method in which the measured values (concentration measurements) are directly approximated by a linear estimate of volumetric reaction rates over time.
  • the course of the reactions can be estimated quickly, assuming that it proceeds linearly between interpolation sites [Leighty, RW, & Antoniewicz, MR (2011), Dynamic metabolic flux analysis (DMFA ): a framework for determining fluxes at metabolic non-steady state. Metabolic engineering, 13 (6), 745-755].
  • the method can also be used for irreversible reactions - such as the elementary modes. If the dimensions of the macro reactions and the measured values do not match, the dimension of the measured values can be adapted to the macro reactions via suitable correlations.
  • This combination of the linear estimation according to Leighty et al. the enhancements to this claim are hereafter referred to as "linear estimation of reaction rates of selected macro reactions".
  • the non-negative least squares solver (NNLS) from Lawson et al is used to solve the linear optimization problem [Lawson, CL and RJ Hanson, Solving Least Squares Problems, Prentice-Hall, 1974, Chapter 23, p. 161.]. This makes it possible to check the quality of larger subsets with the method.
  • the maximum number of macro reactions can also be significantly greater than the number of existing measurements divided by the number of nodes.
  • the method according to the invention of the "linear estimation of reaction rates of selected macro reactions with NNLS" can additionally be used as a further data-dependent method for prereduction of the EMs in step e) i)
  • a very large set K of macro reactions can be used on the one hand the value for SSR C and on the other the course of the reaction rates r (t) EMs with small values of the corresponding rate r (t) are removed from the matrix K. This procedure is repeated until a predefined number of EMs is reached or the value of SSR C exceeds a specified limit.
  • Algorithms for selecting the subset are z. ß. from Provost et al. and Soons et al. known [Provost A. 2006. Metabolism design of dynamic bioreaction models. Faculty of Sciences Appliquees, Universite Catholique de Louvain, Louvain-la-Neuve, Louvain-la-Neuve; Soons, Z. 1. T.A., Ferreira, E.C., Rocha, 1. (2010). Selection of Elementary Modes for Bioprocess Control. 1 Ith International Symposium on Computer Applications in Biotechnology, Leuven, Belgium, July 7-9, 2010, 156-161 J.
  • an evolutionary, in particular a genetic, algorithm is used to select the relevant macro-reactions, ie to select the EM subsets L.
  • a genetic algorithm is z. ß. from Baker et al. [Syed Murtuza Baker, Kai Schallau, Björn H. Junker. 2010. Comparison of different algorithms for simultaneous estimation of multiple parameters in kinetic metabolic models. J. Integrative Bioinformatics:
  • a genetic algorithm can be used, in whose objective function for various combinations of EMs the respective value SSR C is calculated with the method "linear estimation of reaction rates of selected macro reactions.”
  • the matrix L contains the necessary macro reactions (step iii).
  • step iii) the validity of the EM subset L is checked graphically.
  • the flux map from step d) can be used as a projection of the EM subset L.
  • Figure 4 shows the Flux Map with the projection of a subset of six EMs. If the EM subset L is valid, the measurement data remains within the EM subset L. This representation allows a quick graphical check of the validity of the selection.
  • reaction rates of the macro reactions of the subset L are calculated with the specific rates q (t) from d) and / or the measured data from c).
  • the calculation of r (t) can be based on the specific rates q (t) as described in e) according to Soons et al. [Soons, Z. 1.
  • step g) of the method the kinetics of the macro reactions are designed.
  • the determined kinetics should quantify the dynamic influences of the process state on the respective reaction rates r k :
  • the kinetics result in the model parameters to be determined.
  • step g) The generic kinetics are determined in step g) i. designed from the stoichiometry of macro reactions. For substrates of the macro reaction, a limitation of the monodype is assumed.
  • the monod constants K mki and the Hill coefficients ⁇ ⁇ represent the parameters of the equation whose first values are entered manually.
  • the monod constants K mk i are set to one tenth of the respective maximum measured concentrations and the hill coefficients ⁇ ⁇ to the value 1.
  • the determination of generic kinetics from the reaction stoichiometries is described by Provost or by Gao et al. [Provost A. 2006. Metabolic design of dynamic bioreaction models. Faculty of the Sciences Applique, Universite Catholique de Louvain, Louvain-la-Neuve, p.
  • step g) ii. the influencing variables are determined on the reaction rates r (t) determined in f). All variables that describe the process state (ie also bioreaction conditions such as the pH, the reactor temperature, partial pressures that can not be derived from the stoichiometry of the macro reaction) are considered.
  • the influencing variables can be determined manually, for example using a statistical method such as partial least squares. For this purpose, the correlation between the process state (which is summarized in a matrix) and the reaction rates r (t) from f) is determined.
  • iii. the g in g) ii. determined influences quantified and the kinetics of i. extended by corresponding terms.
  • K I ki denotes the inhibition constant and represents another model parameter whose first value is input manually and is usually set to one-tenth of the respective maximum measured concentrations.
  • the model parameter values p of the kinetics are adapted to the reaction rates of the macro reactions r (t) determined in f): m ninn ⁇ (r k () - r k ) (formula 19)
  • model parameter estimation A numerical solution of one or more differential equations according to the formulas 2 to 4 in this step can be dispensed with; the model parameter values can be adjusted in independent groups with usually 3 to 10 parameters separately for each macro reaction k.
  • the adaptation is done by a common method such as the Gauss-Newton method [Bates DM, Watts DG. 1988. Nonlinear regression analysis and its applications. New York: Wiley. xiv, 365.].
  • This model parameter value estimation which is separate for each macro reaction, is particularly advantageous for steps i) and j), since it can be carried out quickly and also provides improved starting values for adapting the model parameter values to measurement data from c) in step j).
  • the goodness of fit is calculated, for example, with the sum of the squared residuals SSR r according to formula 20: N t
  • the check of the quality of fit is done by a graphical comparison of f k ⁇ and r k .
  • step i) the kinetics of the macro reactions selected in g) are checked for their quality of fit.
  • the basis is the value SSR r calculated in step h), which quantifies the quality of fit of the model parameter value estimation.
  • step j) the adaptation of the parameter values of the kinetics from g) to the measured data from c) can be carried out according to a method customary for adaptations.
  • the starting values from step h) are preferably used for this adaptation.
  • the model parameter value adjustment takes place with the inclusion of o. G. Differential equations (formulas 2 to 4), z. B. using the Gauss-Newton method [Bates DM, Watts DG. 1988. Nonlinear regression analysis and its applications. New York: Wiley. xiv, 365. J or using a multiple-shoot algorithm [Peifer M, Timmer J. 2007. Parameter estimation in ordinary differential equations for biochemical processes using the method of multiple shooting. Systems Biology, IET 1 (2): 78-88.J.
  • product features can be integrated into the model. Most preferably, this may be introduced for product characteristics that depend on the concentration of by-products or intermediates. Concentrations of by-products that are external components of the metabolic network entered in a) are already integrated into the model and can be calculated. If necessary, however, other by-products or intermediates may be grouped together in one or more separate metabolic networks. This is advantageous if the expected excretion or uptake rates are in different orders of magnitude or certain metabolic processes are to be considered in different degrees of detail.
  • steps a) to j) can be used to generate a separate model for the calculation of the product features, which also describes the course of the process of the external components of the separate metabolic network with a set of macro reactions with their own kinetics.
  • By-products or intermediates that are not outside the organism but whose intracellular accumulation affects one or more product characteristics may be externalized in step (a) and (b) in the calculation of the EMs and the formulation of macro-reactions, ie classified as external components, become.
  • the involvement of Product features that are dependent on intracellular or out-of-cell concentrations may then be achieved through the additional integration of quantitative or qualitative relationships between concentrations and product characteristics.
  • the model provided by the method according to the invention can be used for process control or planning of the process control as well as investigation of the process in the reactor.
  • the matrix with the calculated EMs E was obtained in step a).
  • the matrix N p contains the stoichiometry, ie the stoichiometric coefficients, of the external metabolite.
  • Possible macro reactions of the stoichiometric network were summarized in the matrix K with formula 21:
  • the measurement data of the process were taken from Baughman et al. which reports various measures of a fermentation of hybridoma cells over the course of a batch process (see Figure 6). In: Computers & Chemical Engineering (2010) 34 (2), pp. 210-222.]. The measurement data was entered in the procedure.
  • the C-mol-related formation rate of the antibody can be estimated.
  • M mAb c _ mol 22.45 - ⁇ .
  • the rate of formation of the antibody then resulted from the formula:
  • step e an EM subset of macro reactions was generated, with which the data record was reproduced in the best possible way.
  • the matrix K from step b) was needed. Since the number of more than 300,000 macro-reactions would have resulted in too many possible combinations, a data-dependent pre-reduction was performed first:
  • the rates q (t) determined in step d) were used to calculate the yield coefficients Y m for all combinations of two external metabolites.
  • the lower limit of a yield coefficient Y tj was chosen such that 99% of the determined yield coefficients Y j ⁇ t) are above this value.
  • the upper limit was chosen so that 99% of the determined Yield coefficients Y ⁇ it) are below this value.
  • some determined limits and the proportion of EMs whose yield coefficients Y TM within these limits are given in Table 2. Overall, the number of EMs could be reduced to about 3000.
  • the subset was selected using a genetic algorithm.
  • the linear optimization problem addressed in the "Linear Estimation of Reaction Rates of Selected Macro Reactions" was solved.
  • the final sum of the least squares of the linear optimization problem calculated here was also the value of the objective function for the respective selection of macro reactions.
  • the optimization was performed repeatedly with a different number of macro reactions in L.
  • the count represents a trade-off between model complexity and rendering accuracy.
  • either the selection of the subset L may be repeated for a varying number of macro-reactions, or a penalty for the number of responses will be used directly
  • Target function of the genetic algorithm can be added.
  • several optimizations were performed with a predefined number of macro reactions (10, 7, 5, 4, and 3). The smallest sum of squares found with the genetic algorithm is plotted against the number of macro reactions in Figure 9. It turned out that in this case fewer than seven macro reactions are too few to represent the course of the process sufficiently well.
  • the selected macro reactions are given in Table 3.
  • Table 3 Selected subset of macro reactions (I). Non-underlined components are not included in the model as there are no measurements for this.
  • the reaction rates over time were determined.
  • the measured values shown in Figure 10 were approximated by an estimate of the reaction rates r (t)
  • the result of the method is a piecewise linear progression of the individual (volumetric) reaction rates Division with the interpolated progression of the vital cell number X v (t), the cell-specific reaction rates r (t) of the macro reactions shown in Table 3 were obtained and the reaction rates r (t) obtained are shown in Figure 10.
  • r k max is the maximum reaction rate
  • Ni the number of limitations taken into account
  • Q the concentration of the component t, K mki the associated monod constants and n; represent the Hill parameter for the reaction order. Their values are adjusted in steps h) and j).
  • the lysis rate Ki was assumed to be constant over the process.
  • the parameters Ciac.cr critical lactate concentration
  • ⁇ , ⁇ maximum rate of death introduced by apoptosis and lysis
  • K d Lac Monitoring parameter for describing the influence of lactate concentration
  • K t lysis rate
  • Target function has many local optima. If you start a deterministic optimization algorithm, such. If, for example, the Levenberg-Marquardt algorithm at the starting values of the parameters known from step h), the chances of success are greatly increased.
  • the adjusted process flow is shown in Figure 12.
  • the adjusted parameters are shown in Table 5. Table 5: Parameters of kinetics as well as apoptosis and lysis
  • the model consisting of the matrix L, the kinetics from Table 4 and the kinetics of apoptosis with the associated parameter values from Table 5 was output.
  • ⁇ (Index k) Denotes the kth element of a vector
  • Figure 1 shows the shift of measured data: The actual course of a measured quantity (C t (t)) is shown, which changes abruptly when the dilution rate changes (D (t)). The shifted history ( is (t)) comes about only by changes caused by the cell.
  • Figure 2 shows the flux map of two specific rates q x and q 2 .
  • the contour lines indicate the frequency with which the respective combination of rates occurs in the measured data.
  • Figure 3 shows a three-dimensional representation of the solution space, which is spanned by a positive linear combination of EMs.
  • the solution space of the entire set is shown, in gray, that of a subset.
  • Figure 4 shows the flux map of two specific rates q i and q 2 .
  • the 2-dimensional projections of the macro reactions of a set I are shown.
  • Figure 5 shows a schematic representation of the metabolic network of Niu et al.
  • the limitation of the cell is shown as a box.
  • the cell-internal demarcation of the mitochondrion is indicated by a dashed line.
  • External components are marked with the index "xt"
  • the arrows and dotted arrows indicate reactions.
  • FIG. 6 shows the measurement data of a fermentation with hybridoma cells from Baughman et al.
  • the total cell count (total cells) is calculated here from the sum of the living cells (vital cells) and dead cells (dead cells).
  • the abbreviations GLC, GLN, ASP, ASN, LAC, ALA and PRO denote the substrate glucose and the amino acids glutamine, aspartic acid, asparagine, alanine and proline as well as the metabolite lactate.
  • the abbreviation MAB designates the product of monoclonal antibodies and BM the biomass.
  • Figure 7 shows the growth and death rates as well as the cell-specific uptake and release rates
  • Figure 8 shows the concentrations approximated by the "linear estimation of reaction rates of selected macro reactions" with the selected reaction set, and the total number of cells (X t ) and the antibody concentration (MAB) were converted to C-mol.
  • Figure 9 shows the smallest calculated sum of error squares ("minimum error") plotted against the number of macro reactions in the subset (n R ).
  • FIG. 10 shows the reaction rates of the macro reactions r (t) determined by the method according to the invention "linear estimation of reaction rates of selected macro reactions".
  • Figure 1 1 shows the reaction rates of macro reactions r (t) (continuous) determined by the method according to the invention "linear estimation of reaction rates of selected macro reactions"
  • Figure 12 shows a comparison of the measured concentrations C m (i) (points) and the simulated process flow (i) (solid line). The concentrations are given in [mM]. Exceptions are the vital and total number of cells (X v / X t in [10 9 cells and the concentration of antibody (mAb in [10 -4 mM]).

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Theoretical Computer Science (AREA)
  • Physiology (AREA)
  • Molecular Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Biotechnology (AREA)
  • Evolutionary Biology (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Biophysics (AREA)
  • Data Mining & Analysis (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Algebra (AREA)
  • Computational Mathematics (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Apparatus Associated With Microorganisms And Enzymes (AREA)

Abstract

Die Anmeldung betrifft ein computerimplementiertes Verfahren zur Erstellung eines Models einer Bioreaktion - Fermentations- oder Ganzzellkatalyseprozess - mit einem Organismus auf Basis von Messdaten.

Description

Computerimplementiertes Verfahren zur Erstellung eines Fermentationsmodels
Die Erfindung betrifft ein Computerimplementiertes Verfahren zur Erstellung eines Models einer Bioreaktion - insb. Fermentation oder Ganzzellkatalyse - unter Verwendung eines Organismus.
Organismus im Sinne der Anmeldung sind Kulturen aus pflanzlichen oder tierischen Zellen wie Säugetierzellen, Hefen, Bakterien, Algen, usw., die in Bioreaktionen eingesetzt werden.
Die sensorische Überwachung eines Fermentationsprozesses und Analyse von Proben aus einem Prozess z. B. mit Hilfe der Quality-by-Design Analytik-Automatisierungsplattform BaychroMAT® der Bayer Technology Services GmbH liefert verschiedene Informationen über den Zustand des Prozesses in dem Bioreaktor in Echtzeit. Typischerweise werden Zellenanzahl, Zellvitalität, Konzentrationen von Substraten, wie Kohlenstoffquellen (z. B. Glukose), Aminosäuren oder O2, Produkten und Nebenprodukten (z. B. Laktat oder CO2,), Prozessparameter wie Temperatur und/oder pH-Wert oder Produktmerkmale ermittelt. Diese Daten können durch berechnete Daten und/oder Extrapolationen z. B. aus dem Stand der Technik ergänzt werden. Zusammen bilden diese Daten die Messdaten oder das Prozesswissen im Sinne der Anmeldung.
Hintergrundwissen zu dem Organismus bedeutet im Sinne der Anmeldung Wissen über die biochemischen Reaktionen des Organismus - spezifische und unspezifische Reaktionen - und insbesondere die zellinternen Reaktionen, bzw. Makroreaktionen zur Beschreibung der organismusspezifischen Stoffwechselnetzwerke (auch SN oder metabolische Netzwerke genannt), die aus Substraten, Metaboliten (auch Knoten des metabolischen Netzwerks genannt), Produkten sowie den biochemischen Reaktionen zwischen ihnen bestehen. Diese biochemischen Reaktionen sind durch ihre:
(a) Stöchiometrie
(b) Reversibilität (unter biologischen Bedingungen) ,
(c) Einbindung in ein stöchiometrisches Netzwerk,
definiert.
Bisher werden die Messdaten hauptsächlich für die qualitative Überwachung des Prozesses eingesetzt. Im Folgenden wird eine Auswahl von technischen Fragestellungen dargestellt, für deren Lösung dynamische Prozessmodelle erforderlich sind. Eine technische Verwendung des Prozesswissens im Sinne der Anmeldung bietet die modellbasierte Zustandsschätzung eines Prozesses in einem Bioreaktor. Methoden wie das„extended Kaiman Filter" ermöglichen eine kontinuierliche Abschätzung von Prozessgrößen, über die diskontinuierlich Messungen vorliegen [Welch G, Bishop G. 1995. An Introduction to the Kaiman Filter. Chapel Hill, NC, USA: University of North Carolina at Chapel Hill.]. Auch der Verlauf nicht messbarer Größen kann aus anderen Messungen berechnet werden. Voraussetzung hierfür ist ein den zugrunde liegenden Prozess richtig beschreibendes Prozessmodell.
Eine weitere Anwendung ist die modellbasierte, optimale Prozessführung. Hierbei wird ein dynamisches Prozessmodell verwendet, um in einem modellbasierten, prädiktiven geschlossenen Regelkreis die Prozessführung hinsichtlich Produktmenge, Produktmerkmalen oder Bildung von Nebenprodukten oder anderen Zielgrößen zu optimieren. Dieses wird beispielsweise von Frahm et al. für eine Hybridoma Zellkultur gezeigt [Frahm B, Lane P, Atzert H, Munack A, Hoffmann M, Hass VC, Portner R. 2002. Adaptive, Model-Based Control by the Open-Loop-Feedback-Optimal (OLFO) Controller for the Effective Fed-Batch Cultivation of Hybridoma Cells. Biotechnol. Prog. 18(5): 1095- 1103J.
Für beide genannten technischen Anwendungen ist es wichtig, dass das erstellte Prozessmodel eine möglichst geringe Komplexität, also eine begrenzte Anzahl an Zustandsvariablen und/oder Gleichungen, bei gleichzeitig guter Genauigkeit der Wiedergabe des Prozesses aufweist.
Neben den genannten Anwendungen zur Prozessführung können dynamische Prozessmodelle auch während der Prozessentwicklung genutzt werden um Experimente mit optimalem Informationsgewinn zu planen. Dieses Vorgehen wird als modellbasierte Versuchsplanung bezeichnet [Franceschini G, Macchietto S. 2008. Model-based design of experiments for parameter precision: State of the art. Chemical Engineering Science 63(19):4846^1872]. Neben den oben genannten Voraussetzungen für die Komplexität des Modells ist es für diese technische Anwendung erforderlich, dass ein dynamisches Prozessmodell bereits während der Entwicklungsphase vorhanden ist. Dieses sollte möglichst schnell aus bereits vorhandenem Prozesswissen erstellt werden können, um den Zeitaufwand der Prozessentwicklung gering zu halten.
Es bestand daher der Bedarf, ein Verfahren bereit zu stellen, dass die Erstellung eines dynamischen Prozessmodells unter der Verwendung von Hintergrundwissen und Prozesswissen ermöglicht. Um dieses Modell zum Beispiel für eine Zustandsschätzung, eine optimale Prozessführung oder zur modellbasierten Versuchsplanung nutzen zu können muss hierzu die Komplexität des Modelles gering sein. Abhängigkeiten, also Einflüsse der Prozessgrößen oder des Prozesszustandes auf das Prozessverhalten, sollen innerhalb des Design Space hinreichend genau quantifiziert sein. Alle vorliegenden Informationen über den Prozesszustand sollen hierfür genutzt werden. Die modellbasierte Beschreibung von Produktmerkmalen soll bei Bedarf in das Modell integrierbar sein. Als Design Space wird der Bereich bezeichnet indem Prozesswissen vorliegt. Das Verfahren soll auf oben aufgeführte Bioreaktionen anwendbar sein und die Entwicklungszeit solcher dynamischen Modelle wesentlich verkürzen. Bisherige Ansätze zur Entwicklung von dynamischen Modellen benötigen Monate bis Jahre bis zur Fertigstellung eines Prozessmodells. Der vorliegende Ansatz reduziert erfahrungsgemäß die Entwicklungszeit auf wenige Wochen.
Typische Produktmerkmale im Sinne der Anmeldung sind beispielsweise Glykosylierungsmuster von Proteinen oder die Proteinintegrität, ohne sich darauf zu begrenzen. Dynamische Modelle, die in dem oben genannten Kontext verwendet werden, weisen diese Eigenschaft bislang nicht auf. Der vorliegende Ansatz ermöglicht eine einfache modellbasierte Integration von Produktmerkmalen.
Die modellbasierte Prozessführung von Fermentationen wird von Frahm et. al am Beispiel einer Hybridoma Zellkultur gezeigt (Frahm B, Lane P, Atzert H, Munack A, Hoffmann M, Hass VC, Portner R. 2002. Adaptive, Model-Based Control by the Open-Loop-Feedback-Üptimal (OLFO) Controller for the Effective Fed-Batch Cultivation of Hybridoma Cells. Biotechnol. Prog. 18(5): 1095-1103). Grundlegende Prozessgrößen werden hier modellbasiert gesteuert. Eine Integration von Produktmerkmalen findet hier nicht statt. Das mathematische Modell der Zelle wurde für diesen spezifischen Prozess entworfen und lässt sich nur mit großem Aufwand auf Prozesse mit dem gleichen oder anderen Organismen oder Stämmen des gleichen Organismus übertragen. Hintergrundwissen in Form von zellinternen Reaktionen wird im Modell nicht explizit berücksichtigt. Eine Integration weiterer Messgrößen in das Modell und somit eine vollständige Nutzung der Informationen über den Prozesszustand kann hier nur mit sehr großem Aufwand erfolgen. Der Ansatz stellt somit eine individuelle Lösung dar, die weder auf andere Prozesse übertragbar ist, noch die volle Nutzung der ermittelten Daten ermöglicht. Die genannte Methode löst aufgrund der zu erwartenden Entwicklungszeit des Modells und der aufwändigen Übertragbarkeit der Lösung auf andere Prozesse mit den gleichen oder mit anderen Organismen das oben genannte technische Problem nicht.
Eine weitergehende Modellierung, die auch Produktmerkmale wie die Glykosylierung mit einschließt, findet sich in den Veröffentlichungen von Kontoravdi et. al. Das Modell, das den Hauptstoffwechsel beschreibt, bezieht kein Hintergrundwissen in Form von zellinternen Reaktionen mit ein und lässt sich auch nicht auf andere Prozesse mit dem gleichen oder anderen Organismen übertragen. Eine Integration weiterer Messgrößen in das Modell kann hier nicht erfolgen [Kontoravdi C, Asprey SP, Pistikopoulos EN, Mantalaris A. 2007. Development of a dynamic model of monoclonal antibody production and glycosylation for product quality monitoring. Computers & Chemical Engineering 31(5-6):392-400.]. Diese Methode ermöglicht auch keine vollständige Nutzung der Informationen über den Prozesszustand, erfordert eine lange Entwicklungszeit des Modells und ist auf andere Organismen oder Stämme nicht übertragbar. Diese Methode stellt somit keine Lösung des technischen Problems dar.
Die Modellierungen der Glykosylierung unter Einbindung des Nucleotidzuckerstoffwechsels von Jedrzejewski et al. und Jimenez et al. beziehen Hintergrundwissen in Form von Bilanzgleichungen interner Stoffwechselintermediate mit ein [Jedrzejewski PM, del Val, Ioscani Jimenez, Constantinou A, Dell A, Haslam SM, Polizzi KM, Kontoravdi C. 2014. Towards Controlling the Glycoform: A Model Framework Linking Extracellular Metabolites to Antibody Glycosylation. International journal of molecular sciences 15(3):4492-4522.; Jimenez del Val, Ioscani, Nagy JM, Kontoravdi C. 2011. A dynamic mathematical model for monoclonal antibody N-linked glycosylation and nucleotide sugar donor transport within a maturing Golgi apparatus. Biotechnology progress 27(6): 1730-1743 J. Bei einer Verwendung dieses Modells für die Prozessführung sind die Komplexität des gesamten Modells und die mangelnde Beobachtbarkeit zellinterner Stoffwechselintermediate aber nachteilig. Zudem lässt das Modell des Hauptstoffwechsels keine Übertragung auf andere Prozesse oder die vollständige Nutzung der Informationen über den Prozesszustand zu. Diese Methode stellt somit keine Lösung des technischen Problems dar.
Eine flexible Modellgenerierung für Bioprozesse wird von Leifheit et al. adressiert [Leifheit J, Heine T, Kawohl M, King R. 2007. Rechnergestützte halbautomatische Modellierung biotechnologischer Prozesse (Semiautomatic Modeling of Biotechnical Processes). at - Automatisierungstechnik 55(5)]. Die Modellgenerierung erfolgt unter Zuhilfenahme von Prozesswissen, allerdings ohne Hintergrundwissen. Die Prozedur ist für verschiedene Prozesse mit dem gleichen oder anderen Organismen verwendbar. Basis sind hier Makroreaktionen, die der Anwender selbst vorgibt. Deren genaue Stöchiometrien werden in dem Verfahren ermittelt. Die Methode wird für eine geringe Anzahl von Zustands- bzw. Messgrößen beschrieben. Eine Integration weiterer Zustands- bzw. Messgrößen wäre mit einem signifikanten Anstieg der Komplexität des Verfahrens verbunden. Bei einer Verwendung von umfassenden Datengrundlagen, wie sie beispielsweise durch die BaychroMAT®- Plattform bereitgestellt wird, wäre diese Methode nicht mehr durchführbar. Eine Integration von Produktmerkmalen ermöglicht die Methode nicht. Sie stellt somit keine Lösung des oben genannten technischen Problems dar. Die Verwendung des Hintergrundwissens in Form von Makroreaktionen, die als Elementary Modes (EM) aus den bekannten metabolischen (stöchiometrischen) Netzwerken eines Organismus gewonnen werden, wird von Provost beschrieben [Provost A. 2006. Metabolie design of dynamic bioreaction models. Faculte des Sciences Appliquees, Universite catholique de Louvain, Louvain-la-Neuve, Louvain-la-Neuve, S. 81 ff., S.107 ff. S. 118 ff,]. Diese Methode ist für verschiedene Organismen oder Stämme des gleichen Organismus verwendbar. Die Auswahl der Makroreaktionen für das Prozessmodell erfolgt unter Verwendung von Prozesswissen. Jedoch werden Prozessabschnitte definiert für die eine vordefinierte Anzahl von Makroreaktionen separat zufällig ausgewählt werden. Das beschriebene Verfahren liefert hier eine von vielen möglichen Kombinationen von Elementary Modes. Die Anzahl der Makroreaktionen, und somit die Modellkomplexität, ist festgelegt und nicht veränderbar. Das Verfahren ergibt separate Modelle für jeden Prozessabschnitt. Die Auswahl der Kinetiken der einzelnen Makroreaktionen erfolgt unter Berücksichtigung der Stöchiometrie der gewählten Makroreaktionen. Die Parameter der Kinetiken (Modellparameter) werden jedoch nicht an die Prozessdaten angepasst. Stattdessen bildet die Verwendung der separaten Prozessabschnittsmodelle die Änderungen in den Prozessdaten ab. Die zufällige Auswahl der Reaktionen kann zwar auch auf Basis einer umfassenden Datengrundlage erfolgen, allerdings können durch den beschriebenen Ansatz zur Auswahl der Kinetiken und die ausgewählten Kinetiken den Verlauf des Prozesses bzw. das Verhalten des Organismus im Prozess nicht abbilden. Die Verwendung mehrerer Prozessabschnittsmodelle führt zudem zu einer unnötigen Komplexitätssteigerung des Prozessmodells. Die Abhängigkeiten, also Einflüsse von Prozessgrößen oder des Prozesszustandes auf das Prozessverhalten, werden mit dieser Methode nicht quantifiziert. Auch findet hier eine Integration von Produktmerkmalen nicht statt. Diese Methode stellt somit keine Lösung des o. g. technischen Problems dar.
Es bestand daher der Bedarf ein Verfahren bereit zu stellen, das die schnelle und effiziente Bereitstellung eines Modells auf Basis von Prozesswissen und -Messdaten und die Optimierung des Produktumsatzes und der kritischen Produktmerkmale unter Berücksichtigung von Hintergrundwissen ermöglicht und die oben genannten Nachteile nicht aufweist.
Die Aufgabe wurde durch ein Verfahren zur Erstellung eines Models einer Bioreaktion mit einem Organismus in einem Bioreaktor wie folgend beschrieben gelöst. Gegenstand der Anmeldung ist ein computerimplementiertes Verfahren zur Erstellung eines Models einer Bioreaktion - insb. Fermentation oder Ganzzellkatalyse - mit einem Organismus, das folgende Schritte umfasst:
a. Ausgewählte Stoffwechselwege des Organismus, deren Stöchiometrie- sowie Reversibilitätseigenschaften werden als Hintergrundwissen in das Verfahren eingegeben. Mit anderen Worten werden ein oder mehrere metabolische Netzwerke des Organismus in das Verfahren eingegeben. Elementary Modes (EMs) werden aus dieser Eingabe berechnet.
b. Die EMs werden in einer Matrix K zusammengefasst, wobei die EMs die Stoffwechselwege aus a) in Makroreaktionen zusammenfassen. Diese Matrix K enthält hiermit die Stöchiometrie und die Reversibilitätseigenschaften aller möglichen Makroreaktionen aus dem Hintergrundwissen.
c. Die Messdaten (auch Prozesswissen genannt) zur Bioreaktion mit dem Organismus werden eingegeben.
d. Mit Hilfe einer Interpolationsmethode werden auf Basis der eingegebenen Messdaten aus c) die für den Organismus spezifischen Raten - Ausscheide- und Aufnahmeraten von einer oder mehreren Eingangsgrößen und Ausgangsgrößen - der eingegebenen Stoffwechselwege berechnet. Bevorzugt werden auch Wachstumsraten, besonders bevorzugt auch Absterberaten des Organismus berechnet.
e. Relevante Makroreaktionen werden in Form eines Subsets der Elementary Modes aus b) ausgewählt durch
i. Datenunabhängige und/oder datenabhängige Vorreduktion der Anzahl der EMs aus b). ii. Auswahl des Subsets aus der Vorreduktion aus e) i. mit den Messdaten aus c) und / oder einer oder mehrerer Raten aus d), bevorzugt mit dem Messdaten aus c) mittels eines Algorithmus nach einem mathematischen Gütekriterium und Zusammenfassung des Subsets in einer Matrix L.
iii. Optional erfolgt eine Darstellung des Subsets grafisch.
f. Mit Hilfe einer Interpolationsmethode werden auf Basis der eingegebenen Messdaten aus c) und / oder der Raten aus d), die Reaktionsraten der Makroreaktionen des Subsets r(t) berechnet.
g. Kinetiken der Makroreaktionen des Subsets aus e) ii. werden mit folgenden Zwischenschritten entworfen; dadurch werden die Modellparameter definiert.
i. Aus der Stöchiometrie der Makroreaktionen werden generische Kinetiken entworfen. ii. Einflussgrößen auf die Makroreaktionen werden aus den Reaktionsraten aus f) ermittelt. iii. Die generischen Kinetiken aus g) i. werden um Terme erweitert, die die in g) ii. ermittelten Einflussgrößen quantifizieren.
h. Optional erfolgt für die Kinetiken aus g) eine erste Anpassung der Modellparameterwerte für jede Makroreaktion separat an die berechneten Reaktionsraten aus f) und eine Überprüfung der Anpassungsgüte.
i. Optional werden die Schritte g) und h) wiederholt bis eine vordefinierte Anpassungsgüte erreicht wird.
j. Die Modellparameterwerte werden an die Messdaten aus c) angepasst.
k. Die Matrix L, die Kinetiken aus g) und die Modellparameterwerte aus j) bilden das Modell und werden ausgegeben und / oder in ein Prozessführungs- oder Prozessentwicklungsmodul übertragen.
Typischerweise kommuniziert das Prozessführungsmodul on-line mit einem Prozessleitsystem, das üblicherweise zur Steuerung des Bioreaktors verwendet wird.
Typischerweise werden Prozessentwicklungsmodule zur off-line Optimierung des Prozesses oder zur Planung von Experimenten verwendet.
Die erfindungsgemäße Modellierung der Bioreaktion baut im Wesentlichen auf der Annahme repräsentativer Makroreaktionen auf, die interne Stoffwechselvorgänge vereinfacht darstellen. Für die Auswahl der Reaktionen wird sowohl biochemisches Hintergrundwissen, als auch Prozesswissen benötigt.
Im ersten Schritt des Verfahrens erfolgt die Eingabe der Reaktionen des metabolischen Netzwerkes, deren Stöchiometrie und Reversibilitätseigenschaft durch den Nutzer über eine Nutzeroberfläche oder idealerweise automatisch durch die Selektion eines Organismus und seiner hinterlegten Stoffwechselwege aus einem Datenbankmodul, in dem das Hintergrundwissen zum Organismus gespeichert ist. Das metabolische Netzwerk (im Stand der Technik auch stöchiometrische Netzwerk benannt) und die Eigenschaften seiner einzelnen Reaktionen stellen das Hintergrundwissen des Organismus dar. Bevorzugt enthält das metabolische Netzwerk Reaktionen aus für den Organismus wichtigen Stoffwechselwegen, beispielsweise Reaktionen der Glykolyse. Besonders bevorzugt enthält die Auswahl externe Reaktionen. Eine externe Reaktion im Sinne der Anmeldung enthält mindestens eine Komponente außerhalb der Zelle, typischerweise mindestens eine Eingangsgröße und/oder mindestens eine Ausgangsgröße (Produkt, Nebenprodukt, usw.). Besonders bevorzugt enthält das metabolische Netzwerk Reaktionen, die das Zellwachstum beschreiben, z. B. in Form einer vereinfachten Reaktion von internen Metaboliten zur externen Biomasse. Abb. 5 und Tabelle 1 im Beispiel beschreiben ein anwendbares metabolisches Netzwerk, ohne sich darauf zu begrenzen.
Dann werden aus den eingegebenen, in einem oder mehreren stöchiometrischen Netzwerken zusammengefassten, Stoffwechselwegen des Organismus Elementary Modes berechnet. Jeder Elementary Mode ist eine Linearkombination von Reaktionsraten aus den Stoffwechselwegen - d. h. internen und externen Reaktionen des metabolischen Netzwerkes, bei der sowohl die„steady State" Bedingung für interne Metabolite erfüllt, als auch die Reversibilität bzw. Irreversibilität von Reaktionen berücksichtigt sind. Bei Linearkombinationen von Reaktionen, die die „steady State" Bedingung für interne Metabolite berücksichtigen, können sich keine internen Metabolite akkumulieren.
Eine interne Reaktion im Sinne der Anmeldung erfolgt ausschließlich innerhalb der Zelle.
Durch eine Externalisierung einer internen Komponente, d. h. durch eine Klassifizierung einer eigentlich internen Komponente als Ein- oder Ausgangsgröße, ist es möglich, die mit der externalisierten, internen Komponenten verbundene interne Reaktion als externe Reaktion zu modellieren und somit die„steady-state" Bedingung für interne Metabolite in diesem Fall zu umgehen. Eine Makroreaktion im Sinne der Anmeldung fasst alle Reaktionen zusammen, die von einer oder mehreren Eingangsgrößen zu einer oder mehreren Ausgangsgröße(n) führen. Jeder Elementary Mode beschreibt somit eine Makroreaktion. Im Vergleich zu der Methode von Leifheit. et al werden im Sinne der Anmeldung die Makroreaktionen auf Basis des eingegebenen Hintergrundwissen ermittelt.
Die Elementary Modes (EMs) werden in einer Matrix E zusammengefasst, vorzugsweise in einem Modul zum Matrixauibau, das mit einem entsprechenden Algorithmus konfiguriert ist. Bekannte Algorithmen können zum Au&au der Elementary Modes Matrix verwendet werden. METATOOL wird als Beispiel genannt, ohne sich darauf zu begrenzen: [Pfeiffer T, Montero F, Schuster S. 1999. METATOOL: for studying metabolic networks. Bioinformatics 15(3):251 -257.]
Mit METATOOL wird eine erste Matrix E, die die eingegebenen internen und externen Reaktionen beschreibt, erzeugt.
In Schritt b) wird mit Hilfe der (externen) stöchiometrischen Matrix Np aus der Matrix (E) eine Matrix bestehend aus möglichen Makroreaktionen K erzeugt.
K = NV - E (Formel l ) Die Transformation der Matrix E nach K ist aus Provost bekannt [Provost A. 2006. Metabolie design of dynamic bioreaction models. Faculte des Sciences Appliquees, Universite catholique de Louvain, Louvain-la-Neuve, Louvain-la-Neuve, S. 81].
Die Spaltenvektoren der Matrix K beschreiben die Makroreaktionen. Die Zeilenvektoren beschreiben die Komponenten der Makroreaktionen (Eingangs- und Ausgangsgrößen). In die Matrix K wird die Stöchiometrie der Makroreaktionen eingetragen.
Jede im Sinne des metabolischen Netzwerks mögliche Reaktionsrate kann als positive Linearkombination dieser Makroreaktionen dargestellt werden.
Die Verwendung von EMs als Basis eines Prozessmodells ist im Stand der Technik bekannt (z. B. aus Provost A. 2006. Metabolie design of dynamic bioreaction models. Faculte des Sciences Appliquees, Universite catholique de Louvain, Louvain-la-Neuve, Louvain-la-Neuve S. 87, S. 118 ff. und Gao, J. et al. (2007). Dynamic metabolic modeling for a MAB bioprocess. Biotechnology pmgress, 23(1), 168- 181).
In einem weiteren Schritt c) werden die verfügbaren Messdaten (Prozesswissen) zur Bioreaktion mit dem Organismus eingegeben. Typischerweise werden Zellenanzahl, Zellvitalität, Konzentrationen von Substraten, wie Kohlenstoffquellen (z. B. Glukose), Aminosäuren oder O2, Produkten und Nebenprodukten (z. B. Laktat oder CO2), Prozessparameter wie Temperatur und/oder pH-Wert oder Produktmerkmale ermittelt. Diese Eingabe kann manuell durch den Nutzer oder automatisch erfolgen, wie z. B. durch Selektion aus einem Datenbankmodul zur Speicherung von Messdaten und Transfer der selektierten Daten in ein Datenanalysemodul, das mit dem Datenbankmodul verbunden ist.
Aus diesen Messdaten werden im Schritt d) die zellspezifischen Ausscheide- und Aufnahmeraten von Substraten und (Neben)produkten - zusammen spezifische Raten q (t) - genannt, sowie optional die
Wachstums- und Absterberaten des Organismus (μ(ί), MdCO) berechnet. Voraussetzung für die Berechnung ist die Interpolation der vitalen Zellzahl, der Gesamtzellzahl und Medienkonzentrationen mit Hilfe einer Aus diesen können zeitliche Änderungen der Messgrößen ermittelt werden. Die berechneten Raten q(t), μ( , geben Auskunft über das zu beobachtende dynamische Verhalten des Organismus über der Zeit.
Für die Berechnung der o. g. Raten kann eine oder mehrere unterschiedliche Interpolationsmethode(n) in Kombination verwendet werden. Exemplarisch beschreibt Leifheit et al. die Ermittlung der zeitlichen Änderungen von Messgrößen - z. B. die der Gesamtzellzahl, die der vitalen Zellzahl oder die von anderen Medienkonzentrationen aus Messdaten mit Hilfe von Spline-interpolierten Messdaten [Leifheit, J., Heine, T., Kawohl, M., & King, R. (2007). Rechnergestützte halbautomatische Modellierung biotechnologischer Prozesse (Semiautomatic Modeling of Biotechnical Processes). at-Automatisierungstechnik, 55(5), 21 1 -218]. Diese Methode wird hiermit per Referenz in die Anmeldung integriert.
Die o. g. Raten q (t), μ(ί), dC werden aus diesen zeitlichen Änderungen berechnet:
Beispielsweise kann die Wachstumsrate des Organismus μ(ΐ) aus Spline-interpolierten Werten der Gesamtzellzahl Xt (t) und der Lebendzellzahl Xv (t) sowie der daraus berechenbaren zeitlichen
lX
Änderung der Gesamtzellzahl (t) mit Formel 2 berechnet werden:
^ (C) = -D(t) ■ Xt (t) + μ(ΐ) ■ Xv(t) (Formel 2)
Wobei D (t) die Verdünnungsrate ist.
Die Absterberate μ^ (£) kann bei bekanntem Verlauf von μ(ί) aus dem Verlauf von Xv (t) und dem
iX
Verlauf der zeitlichen Änderung der vitalen Zellzahl (t) mit Formel 3 berechnet werden:
^ (t) = -D t) ■ Xv (t) + (μ( - μα (ΐ))■ Xv (t) (Formel 3)
Die Berechnung der spezifischen Raten einer anderen Komponente t qj(t) kann aus Spline- interpolierten Werten der Lebendzellzahl Xv (t) und der Konzentration der Komponente Q (t), sowie dem Verlauf der zeitlichen Änderung ^ (t), der aus 8ρ1ίη6-ίηΐ6 θ1ίειΐ6η Werten von Q(t) ermittelt werden kann, mit Formel 4 erfolgen:
^ (t) = D (t) · (cUn - Q( ) + qi (t) ■ Xv(t) (Formel 4)
In einer bevorzugten Ausführungsform des Verfahrens werden die Messdaten aus Schritt c) vor der ersten Interpolation folgendermaßen aufbereitet: Um alle nicht von den Zellen verursachten Konzentrationsänderungen zu berücksichtigen und aus den Messdaten einen stetigen Verlauf der Konzentrationsänderungen zu erhalten, werden die Messdaten verschoben (in der Anmeldung shiften genannt). Der Betrag AQ(t), um den die Konzentrationsmessung verschoben wird, kann nach Formel 5 berechnet werden:
AQ( = J (T) · (Q(T) - CUN (T)) dr (Formel 5)
Wobei D ( ) die Verdünnungsrate ist. Der geshiftete Konzentrationsverlauf Ci s (t) ergibt sich dann nach Formel 6: Q,s( = Q( - AQ(C) (Formel 6)
Die Differentialgleichung, die den Verlauf der geshifteten Konzentration Ci s (t) angibt, ergibt sich somit aus den Formeln 4 und 6 zu:
dC,„ d (Ci (t) -ACi(t)) , . T, , . ._ ,
~ä = dt = <?i (0 · χν W (Formel 7)
Diese Aufbereitung (shiften) der Messdaten verhindert eine sprunghafte Änderung der berechneten spezifischen Raten beim An- oder Abstellen eines Feeds (Feedpeak) insbesondere in einem Fed-Batch Prozess.
Abbildung 1 zeigt das Aufbereitung/Shiften der Messdaten im Sinne dieser Anmeldung.
In einer besonderen Ausführungsform des Verfahrens werden die aufbereiteten Daten dann zur
dX I
Berechnung eines Gradienten der Gesamtzellzahl — -\ (t) mit der Methode von Leifheit et al.
dt l s
verwendet. Dieser wird mit einer Spline-Interpolation gemäß Differentialgleichung 8 angenähert: ^| (t) = μ(ί)■ Xv (t (Formel 8)
Besonders bevorzugt wird die Lysis in der Berechnung anhand eines Lysis-Faktors Kt, einbezogen (Formel 9). Dieser kann z. B. als konstant über dem Prozessverlauf angenommen werden.
fl y
= μ(ί) · Xv (t) - Kt ■ (Xt (t) - Xv (t)) (Formel 9)
Das Absinken der geshifteten Gesamtzellzahl XtiS (t) kann so durch Lysis erklärt werden, wodurch negative Werte für die Wachstumsrate μ(ΐ) vermieden werden können.
Bevorzugt werden die aufbereiteten Daten auch für die Berechnung der Absterberate μά (t) verwendet.
Bevorzugt werden die möglichen Kombinationen spezifischer Raten q (t) in einem Flux-Map
Diagramm wie beispielhaft in Abbildung 2 dargestellt. Diese Darstellung ermöglicht eine gute Übersicht über die berechneten spezifischen Raten q(t). Die Höhenlinien geben hier an, welche
Bereiche physiologisch wichtig sind.
Wenn die spezifischen Raten q(t) und optional die weiteren Raten ( und μ^ (ί) unterschiedliche
Größenordnungen und Einheiten besitzen, werden diese üblicherweise mit Hilfe von Vereinfachungen zu einem spezifischen Ratenvektor q (t) mit gleichen Einheiten zusammengefasst. Beispielsweise wird
3 die spezifische Rate eines Makromoleküls, das in Gramm [g] gemessen wird, die Einheit
ICell hi besitzen. Wenn die Zusammensetzung dieses Makromoleküls geschätzt wird, z. B. basierend auf seinem C-mol-Gehalt, lässt sich seine spezifische Rate von [g] auf [C— mol] verändert darstellen, so dass die spezifische Rate die Einheit [^e^] besitzt.
Die spezifischen Raten q(t) bilden eine der Grundlagen für den weiteren Schritt e) des Verfahrens, nämlich die Auswahl der relevanten Makroreaktionen.
Im Schritt e) wird ein Subset (L) der EMs auf Basis der Daten ausgewählt, mit dem die spezifischen Raten q(t) aus d) und / oder die Messdaten aus c) nach einem mathematischen Gütekriterium gut abgebildet werden können. Die Anzahl der EMs in dem Subset (L) soll möglichst klein sein, um eine möglichst geringe Komplexität des Prozessmodells zu gewährleisten. Das Subset L soll aber eine gute Beschreibung des Prozesswissens sicherstellen.
Die Selektion von EMs verkleinert den Lösungsraum im Vergleich zum original EMs-Set (K) aus a), aber enthält weiterhin den ermittelten physiologisch wichtigen Bereich der Zellen.
Abbildung 3 zeigt eine Darstellung des Lösungsraums, wobei das originale Set von EMs (K) zu einem Subset (L) reduziert wird.
Für den Schritt e) werden die berechneten spezifischen Raten q (t) sowie die Messdaten aus c) üblicherweise in ein Modul zur Auswahl der relevanten Makroreaktionen transferiert, das mit entsprechenden Algorithmen konfiguriert ist.
Im Schritt e) i) erfolgt eine datenunabhängige und/oder eine datenabhängige Vorreduktion der Matrix K in beliebiger Reihenfolge:
Die datenunabhängige Vorreduktion erfolgt vorzugsweise durch eine geometrische Reduktion. Dabei werden für einen zufällig ausgewählten EM alle Kosinus-Ähnlichkeiten zu allen anderen Modes berechnet. Der EM mit der höchsten Ähnlichkeit wird aus dem Set entfernt. Diese Vorgehensweise wird wiederholt, bis eine vordefinierte Anzahl an EMs erreicht ist. Die gewünschte Anzahl wird üblicherweise für das Verfahren im Vorfeld definiert. Als Kontrollvariable kann das Volumen des Lösungsraumes verwendet werden. Überraschenderweise wurde festgestellt, dass eine deutliche Reduktion der Anzahl der Makroreaktionen unter Beibehaltung von 90 bis 98 %, bevorzugt 92 bis 95% des aufgespannten Volumens, im Vergleich zum Originalvolumen möglich ist.
Die datenabhängige Vorreduktion kann durch den Vergleich von Ausbeutekoeffizienten der EMs (YEM) zu den Ausbeutekoeffizienten, die aus den spezifischen Raten q (t) aus d) berechnet werden
(Ym) erfolgen. Der Ausbeutekoeffizient des k-ten EMs (i ) 'k) wird nach Formel 10 durch die Division der entsprechenden stöchiometrischen Koeffizienten der externen Metaboliten i und j ermittelt. Für den k-ten EM sind dies die Matrixeinträge Ki k und KJik.
yEM.k = ^ (Formd
Ist der stöchiometrische Koeffizient KJik = 0, kann der Ausbeutekoeffizient nicht ermittelt werden. Der Ausbeutekoeffizient Yj (t) gibt das Verhältnis zweier gemessener oder nach d) berechneter zell- spezifischen Raten ( i (t), qj (t)) zueinander nach Formel 1 1 an:
*w (0 = H (Formel l l)
Aus den Ausbeutekoeffizienten Ym können für jede mögliche Kombination zweier externer Komponenten t und j eine obere und eine untere Grenze bestimmt werden. Beispielsweise kann der kleinste Ausbeutekoeffizient zweier externer Metabolite i und j Υ™(ΐ) als untere Grenze und der größte Wert von Υ™(ΐ) als obere Grenze verwendet werden, andere Grenzen sind jedoch auch möglich. EMs, deren Ausbeutekoeffizienten Y™ oberhalb der oberen Grenze oder unterhalb der unteren Grenze liegen, werden aus der Matrix K entfernt. Wenn der Ausbeutekoeffizient eines EMs Y™ nicht ermittelt werden kann, verbleibt dieser in der Matrix K. Bevorzugt kann auch die auf Seite 15 beschriebene, erfindungsgemäße Methode „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen mit NNLS" zur datenabhängigen Vorreduktion verwendet werden. Eine Erläuterung der Methode im Kontext der Vorreduktion findet sich dort. Das Heranziehen der Prozessdaten in der datenabhängigen Vorreduktion in dem Verfahren hat den Vorteil, dass die Reduktion prozessbezogen und damit effektiver und fokussierter erfolgt.
Im Schritt e) ii. wird ein Subset aus Makroreaktionen mit einem Algorithmus ausgewählt: Für die Auswahl wird ein Gütekriterium, mit dem quantifiziert werden kann, wie gut die spezifischen Raten q (t) aus d) und / oder die Messdaten aus c) mit einem Subset (L) abgebildet werden können, sowie ein Algorithmus zur Auswahl des Subsets benötigt.
Als Gütekriterium für die Abbildung von berechneten spezifischen Raten q(t) mit einem Subset L kann nach Soons et al. die Summe der quadrierten Residuen der spezifischen Raten (SSRq) nach Formel 12 verwendet werden [Soons, Z. I. T. A., Ferreira, E. C, Rocha, I. (2010). Selection of Elementary Modes for Bioprocess Control. 1 Ith International Symposium on Computer Applications in Biotechnology, Leuven, Belgium, July 7-9, 2010, 156-161]. SSR, (Formel 12)
t=l
Der Wert für SSRq soll möglichst klein sein.
Für die Minimierung von SSRq muss zuvor für jeden betrachteten Zeitpunkt der Vektor r(tj) mit Hilfe eines non-negative-Least-squares Algorithmus so bestimmt werden, dass gilt:
(Formel 13)
mit der zusätzlichen Randbedingung:
r(td > 0 (Formel 14)
Vorteil dieser Methode ist, dass die Berechnungen nach Formel 12 - 14 auch für sehr große Subsets mit vielen EMs durchgeführt werden können. Ein signifikanter Nachteil ist, dass für diese Berechnung die berechneten spezifischen Raten q(t) erforderlich sind. Da diese aus interpolierten Messwerten gewonnen werden, liegen diese mit großer Unsicherheit bezüglich ihren wahren Werte vor. Messungenauigkeiten können sich unter Umständen sehr stark auf die berechneten spezifischen Raten c/(t) auswirken. Das Gütekriterium SSRq kann folglich auch nur unter großer Unsicherheit bestimmt werden. Neben der Information über die üüte der Abbildung wird mit dieser Methode auch ein geschätzter Verlauf der Reaktionsraten r (t) des Subsets L als Ergebnis der Minimierung nach Formel 13 und 14 erhalten.
Leighty, R. et al. beschreibt eine weitere Methode in der die Messwerte (Konzentrationsmessungen) direkt durch eine lineare Schätzung von volumetrischen Reaktionsraten über der Zeit approximiert werden. Durch Lösen eines linearen Optimierungsproblems mit einem Linear Least-Squares Solver kann der Verlauf der Reaktionen schnell geschätzt werden, unter der Annahme, dass dieser linear zwischen Stützstellen verläuft [Leighty, R. W., & Antoniewicz, M. R. (2011), Dynamic metabolic flux analysis (DMFA): a framework for determining fluxes at metabolic non-steady State. Metabolic engineering, 13(6), 745-755]. Diese Methode bezieht sich lediglich auf reversible Makroreaktionen (wie den in der Quelle genannten „Free Fluxes"), zudem können Verdünnungseffekte (also Konzentrationsänderungen, die nicht durch die Zellen verursacht werden) nicht berücksichtigt werden. Stimmen die Dimensionen der Makroreaktionen und der Messwerte nicht überein können diese Messwerte nicht zur Schätzung von Reaktionsraten verwendet werden. Dies ist zum Beispiel der Fall wenn das Zellwachstum in Form von Bildung externer Biomasse Teil der Makroreaktionen ist und aus den Messwerten nur die Zelltrockenmasse bekannt ist. In dieser Form ist diese Methode daher nicht auf die Anwendung von irreversiblen Makroreaktionen und Fed-Batch Prozessen geeignet. Verwendet man das Konzept von Leighty et al, das hiermit per Referenz in diese Anmeldung integriert wird, mit den erfindungsgemäß aufbereiteten (geshifteten) Daten, lässt sich diese Methode jetzt auch auf Fed-ßatch Prozesse anwenden. Darüber hinaus kann durch das Hinzufügen einer unteren Grenze für die Reaktionsraten der Makroreaktionen als Randbedingung des linearen Optimierungsproblems die Methode auch für irreversible Reaktionen - wie die Elementary Modes - verwendet werden. Stimmen die Dimensionen der Makroreaktionen und der Messwerte nicht überein kann über geeignete Korrelationen die Dimension der Messwerte an die der Makroreaktionen angepasst werden. Diese Kombination der linearen Schätzung nach Leighty et al. mit den Verbesserungen aus diesem Anspruch wird im Folgenden als „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" bezeichnet.
Es lässt sich so überprüfen, ob sich mit den gewählten Makroreaktionen eines Subsets L des originalen EM-Sets K die gemessenen Daten hinreichend darstellen lassen. Die hier errechnete finale Summe der quadrierten Residuen SSRC nach Formel 15 zwischen den mit der Methode ermittelten geshifteten Konzentrationen Cs(t) und den geshifteten Konzentrationen C;(t) gibt an wie gut die Messdaten mit dem Subset abgebildet werden können.
(Formel 15)
Je kleiner der Wert von SSRC , desto besser ist das Subset L. Diese Methode wird gegenüber der Methode von Soons et al., insbesondere für die Modellierung von Fed-Batch Prozessen bevorzugt, da eine Überprüfung der Güte eines Subsets schnell auch ohne eine eventuell fehlerbehaftete vorherige Bestimmung der spezifischen Raten möglich ist. Durch die Annahme, dass die geschätzten Reaktionsraten linear zwischen Stützstellen verlaufen, wirken sich Messabweichungen sehr wenig auf die Abschätzung der Reaktionsraten aus. Der Nachteil dieser Methode ist, dass die Größe des zu untersuchenden Subsets L durch die Lösung des linearen Optimierungsproblems begrenzt ist. Die maximale Anzahl der Reaktionen in dem Subset ist gleich der Anzahl an vorhandenen Messungen geteilt durch die Anzahl der Stützstellen.
Neben der Information über die Güte der Abbildung wird mit dieser Methode auch ein geschätzter Verlauf der Reaktionsraten des Subsets r (t) erhalten.
In einer bevorzugten Ausführung der erfindungsgemäßen„lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" wird zur Lösung des linearen Optimierungsproblems statt eines Linear Least-Squares Solvers der Non-negative Least-Squares Solver (NNLS) von Lawson et al verwendet [Lawson, C.L. and R.J. Hanson, Solving Least Squares Problems, Prentice-Hall, 1974, Chapter 23, p. 161.]. Dadurch ist es möglich, auch die Güte von größeren Subsets mit der Methode zu überprüfen. Die maximale Anzahl der Makroreaktionen kann hierbei auch signifikant größer als die Anzahl an vorhandenen Messungen geteilt durch die Anzahl der Stützstellen sein. Diese Kombination der„lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" mit der Verwendung des Non-negative Least-Squares Solvers wird im Folgenden als„lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen mit NNLS" bezeichnet.
Die erfindungsgemäße Methode der „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen mit NNLS" kann zusätzlich als weitere datenabhängige Methode zu Vorreduktion der EMs in Schritt e) i) verwendet werden. Hierfür kann hier ein sehr großes Set K von Makroreaktionen verwendet werden. Ergebnis der Methode ist zum einen der Wert für SSRC und zum anderen der Verlauf der Reaktionsraten r(t). EMs mit kleinen Werten der dazugehörigen Rate r(t) werden aus der Matrix K entfernt. Diese Vorgehensweise wird wiederholt, bis eine vordefinierte Anzahl an EMs erreicht ist oder der Wert von SSRC einen festgelegten Grenzwert überschreitet.
Algorithmen zur Auswahl des Subsets sind z. ß. aus Provost et al. und Soons et al. bekannt [Provost A. 2006. Metabolie design of dynamic bioreaction models. Faculte des Sciences Appliquees, Universite catholique de Louvain, Louvain-la-Neuve, Louvain-la-Neuve ; Soons, Z. 1. T. A., Ferreira, E. C, Rocha, 1. (2010). Selection of Elementary Modes for Bioprocess Control. 1 Ith International Symposium on Computer Applications in Biotechnology, Leuven, Belgium, July 7-9, 2010, 156-161 J.
Soons et al. beschreiben die Bildung eines EM-Subsets in einem zweistufigem Optimierungsverfahren. Für verschiedene, zufällig ausgewählte EMs werden jeweils die Werte für SSRq wie oben beschrieben minimiert. Das Set mit dem kleinsten minimierten SSRq Wert wird ausgewählt. Bei einer großen Zahl von EMs ist jedoch die zufällige Auswahl ineffektiv, da die Anzahl der möglichen Kombinationen sehr stark wächst. Beispielsweise gibt es bei der Auswahl von 10 Reaktionen aus einem Set von 20.000 EMs mehr als 2,8 · 1036 Kombinationen. Die Wahrscheinlichkeit, dass die optimale Kombination gefunden wird ist hier sehr gering. Durch die Verwendung des Gütekriteriums SSRq ist diese Methode anfällig gegenüber Messunsicherheiten und Messabweichungen.
Provost beschreibt einen alternativen Algorithmus, bei dem für verschiedene spezifische Werte von q (tj) t = Ι, .,. , η alle möglichen positiven Linearkombinationen von Elementary Modes ermittelt werden bei denen gilt: SSRq (t = tj) = 0. Von diesen zahlreichen möglichen Kombinationen wird anschließend eine Kombination zufällig ausgewählt. Dieses Verfahren verwendet jeweils nur einen Vektor q (tj) und nicht den gesamten Verlauf über der Zeit. Eine Auswahl von EMs für den gesamten
Prozess ist deshalb nicht möglich. Mit der zufälligen Auswahl lässt sich zwar der Vektor q(ti) darstellen, inwiefern der Rest des Prozesses hiermit dargestellt werden kann, wird jedoch nicht ermittelt. Ein weiterer Nachteil dieses Verfahrens ist, das ein Vektor q tt) , der nicht im Lösungsraum aller EMs liegt, nicht verwendet werden kann. Eine Näherungslösung lässt sich nicht ermitteln. Dies ist vor allem bei mit Unsicherheit behafteten Messungen und spezifischen Raten ein großer Nachteil. Durch die Verwendung des Gütekriteriums SSRq ist diese Methode ebenfalls anfällig gegenüben Messunsicherheiten und Messabweichungen.
In einer weiteren und bevorzugten Ausführungsform des Verfahrens wird daher zur Auswahl der relevanten Makroreaktionen, d. h. zur Auswahl des EM-Subsets L ein evolutionärer, insbesondere ein genetischer, Algorithmus verwendet. Ein solcher Algorithmus ist z. ß. aus Baker et al. bekannt [Syed Murtuza Baker, Kai Schallau, Björn H. Junker. 2010. Comparison of different algorithms for simultaneous estimation of multiple parameters in kinetic metabolic models. J. Integrative Bioinformatics:-l- l.J. Besonders bevorzugt kann ein genetischer Algorithmus verwendet werden, in dessen Zielfunktion für verschiedene Kombinationen von EMs der jeweilige Wert SSRC mit der Methode „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" berechnet wird. Alternativ kann auch eine zufällige Auswahl verwendet werden. Nach Abschluss des Schritt ii) beinhaltet die Matrix L die notwendigen Makroreaktionen (Schritt iii).
In einem optionalen Schritt iii) wird die Validität des EM-Subsets L grafisch überprüft. Es kann hier die Flux Map aus Schritt d) als Projektion des EM-Subsets L herangezogen werden. Abbildung 4 zeigt die Flux Map mit der Projektion eines Subsets von sechs EMs. Ist der EM-Subset L valide, befinden sich die Messdaten weiterhin innerhalb des EM-Subsets L. Diese Darstellung ermöglicht eine schnelle grafische Kontrolle der Validität der Auswahl.
In einem weiteren Schritt f) werden mit den spezifischen Raten q(t) aus d) und / oder den Messdaten aus c) die Reaktionsraten der Makroreaktionen des Subsets L berechnet. Die Berechnung von r(t) kann anhand der spezifischen Raten q(t) wie in e) beschrieben nach Soons et al. erfolgen [Soons, Z. 1.
T. A., Ferreira, E. C, Rocha, 1. (2010). Selection of Elementary Modes for Bioprocess Control. 1 Ith International Symposium on Computer Applications in Biotechnology, Leuven, Belgium, July 7-9, 2010, 156-161], bevorzugt erfolgt die Berechnung von r(t) anhand der Messdaten aus c) mit der erfindungsgemäßen„linearen Schätzung von Reaktionsraten ausgewählter Makroreaktionen".
Im Schritt g) des Verfahrens werden die Kinetiken der Makroreaktionen entworfen. Die ermittelten Kinetiken sollen die dynamischen Einflüsse des Prozesszustands auf die jeweiligen Reaktionsraten rk quantifizieren:
rk = f( , pH, T, ... ) (Formel 16)
Aus den Kinetiken ergeben sich die zu bestimmenden Modellparameter.
Die generischen Kinetiken werden in Schritt g) i. aus der Stöchiometrie der Makroreaktionen entworfen. Für Substrate der Makroreaktion wird eine Limitierung vom Monod Typ angenommen.
itierungen j multipliziert:
(Formel 17)
Wobei die Monod-Konstanten Km k i und die Hill-Koeffizienten ηέ die Parameter der Gleichung darstellen, deren erste Werte manuell eingegeben werden. Üblicherweise werden die Monod- Konstanten Km k i auf ein Zehntel der jeweiligen maximalen gemessenen Konzentrationen und die Hill-Koeffizienten ηέ auf den Wert 1 gesetzt. Die Ermittlung von generischen Kinetiken aus den Reaktionsstöchiometrien wird von Provost oder von Gao et al. beschrieben [Provost A. 2006. Metabolic design of dynamic bioreaction models. Faculte des Sciences Appliquees, Universite catholique de Louvain, Louvain-la-Neuve, S. 126 ; [Gao, J., Gorenflo, V. M., Scharer, J. M., & Budman, H. M. (2007). Dynamic metabolic modeling for a MAB bioprocess. Biotechnology progress, 23(1 ), 168-181]. Diese Methoden werden hiermit in die Anmeldung per Referenz integriert. In diesen Methoden werden der Substratlimitierungen des Monod-Typs für die jeweiligen Substrate einer Reaktion verwendet. Obwohl Provost oder Gao dies nicht beschreiben, sind Inhibierungen durch toxische Produkte ebenfalls mit dieser Methode aus der Reaktionsstöchiometrie ableitbar.
Im Schritt g) ii. werden die Einflussgrößen auf die in f) ermittelten Reaktionsraten r(t) ermittelt. Dabei werden alle Größen, die den Prozesszustand beschreiben (also auch Bioreaktionsbedingungen wie z. B. der pH-Wert, die Reaktortemperatur, Partialdrücke, die nicht aus der Stöchiometrie der Makroreaktion ableitbar sind) betrachtet. Die Einflussgrößen können manuell beispielsweise mit Hilfe einer statistischen Methode wie Partial-Least-Squares ermittelt werden. Hierzu wird die Korrelation zwischen dem Prozesszustand (der in einer Matrix zusammengefasst wird) und den Reaktionsraten r(t) aus f) ermittelt. In einem Schritt g) iii. werden dann die in g) ii. ermittelten Einflüsse quantifiziert und die Kinetik aus i. um entsprechende Terme erweitert. Ein in g) ii. gefundener Einfluss von einer Größe des Prozesszustandes auf eine Reaktionsrate kann dann mit einem Term j belegt werden. Der Term j ist eine beliebige Funktion, die in Abhängigkeit des Prozesszustandes einen Wert zwischen 0 und 1 ergibt. Die in g) i. aufgestellte generische Kinetik der Reaktion wird dann mit diesem Term multipliziert. Beispielsweise deutet eine gefundene negative Korrelation zwischen der Konzentration einer Komponente t und der Reaktion k auf eine Beeinflussung der Reaktionsrate k durch die Konzentration vom (Q) hin. Dieses kann beispielsweise mit einer Inhibierungskinetik nach Haidane belegt werden: I ( T) = K KI,k, Ki + 'r m W (Formel 18)
Wobei KI k i die Inhibierungskonstante bezeichnet und einen weiteren Modellparameter darstellt, dessen erster Wert manuell eingegeben wird und üblicherweise auf ein Zehntel der jeweiligen maximalen gemessenen Konzentrationen gesetzt wird.
In einem optionalen Schritt h) werden die Modellparameterwerte p der Kinetiken an die in f) ermittelten Reaktionsraten der Makroreaktionen r (t) angepasst: m niinn ^ (rk ( ) - rk ) (Formel 19)
v ■
- k=l
Dieses wird im Folgenden als Modellparameterwertschätzung bezeichnet. Auf ein numerisches Lösen einer oder mehrerer Differentialgleichungen nach den Formeln 2 bis 4 in diesem Schritt kann verzichtet werden; die Modellparameterwerte können in unabhängigen Gruppen mit üblicherweise 3 bis 10 Parametern separat für jede Makroreaktion k angepasst werden. Die Anpassung erfolgt durch eine übliche Methode wie der Gauß-Newton Methode [Bates DM, Watts DG. 1988. Nonlinear regression analysis and its applications. New York: Wiley. xiv, 365.].
Diese für jede Makroreaktion separate Modellparameterwertschätzung ist für die Schritte i) und j) besonders vorteilhaft, da sie einerseits schnell durchführbar ist und anderseits verbesserte Startwerte für die Anpassung der Modellparameterwerte an Messdaten aus c) in Schritt j) liefert.
Die Anpassungsgüte wird beispielsweise mit der Summe der quadrierten Residuen SSRr nach Formel 20 berechnet: Nt
SSRr = ^ (rk ( ) - rk) (Formel 20)
k= l
Je kleiner der Wert für SSRr, desto besser ist die Anpassung. Alternativ erfolgt die Überprüfung der Anpassungsgüte durch einen graphischen Vergleich von fk~ und rk .
In einem optionalen Schritt i) werden die in g) gewählten Kinetiken der Makroreaktionen auf ihre Anpassungsgüte hin überprüft. Grundlage ist der in Schritt h) berechnete Wert SSRr, der die Anpassungsgüte der Modellparameterwertschätzung quantifiziert. Bei einer nicht zufriedenstellenden Anpassungsgüte können die Schritte g) und h) wiederholt werden, bis eine vordefinierte Anpassungsgüte erreicht wird.
In einem weiteren Schritt j) kann die Anpassung der Parameterwerte der Kinetiken aus g) an die Messdaten aus c) nach einer für Anpassungen üblichen Methode erfolgen. Bevorzugt werden für diese Anpassung die Startwerte aus Schritt h) verwendet. Die Modellparameterwertanpassung erfolgt unter Einbeziehung der o. g. Differentialgleichungen (Formeln 2 bis 4), z. B. mit Hilfe der Gauß-Newton Methode [Bates DM, Watts DG. 1988. Nonlinear regression analysis and its applications. New York: Wiley. xiv, 365. J oder mit einem multiple Shooting Algorithmus [Peifer M, Timmer J. 2007. Parameter estimation in ordinary differential equations for biochemical processes using the method of multiple shooting. Systems Biology, IET 1(2):78-88.J.
Bevorzugt können Produktmerkmale in das Modell integriert werden. Besonders bevorzugt kann dies für Produktmerkmale, die von der Konzentration von Neben- oder Zwischenprodukten abhängen, eingeführt werden. Konzentrationen von Nebenprodukten, die externe Komponenten des in a) eingegebenen metabolischen Netzwerkes sind, sind bereits in das Modell integriert und können berechnet werden. Bei Bedarf lassen sich allerdings auch andere Neben-oder Zwischenprodukte in einem oder mehreren separaten metabolischen Netze zusammenfassen. Dies ist vorteilhaft, wenn die erwarteten Ausscheide- oder Aufnahmeraten sich in unterschiedlichen Größenordnungen bewegen oder bestimmte Stoffwechselvorgänge in unterschiedlichen Detailgraden betrachtet werden sollen. Alternativ zur einem integrierten Modell kann mit den Schritten a) bis j) ein separates Modell zur Berechnung der Produktmerkmale erzeugt werden, dass den Prozessverlauf der externen Komponenten des separaten metabolischen Netzwerkes ebenfalls mit einem Satz von Makroreaktionen mit eigenen Kinetiken beschreibt. Neben- oder Zwischenprodukte, die nicht außerhalb des Organismus vorliegen, aber deren zellinterne Akkumulation sich auf ein- oder mehrere Produktmerkmale auswirkt können in Schritt a) und b) bei der Berechnung der EMs und der Formulierung der Makroreaktionen externalisiert, also als externe Komponenten klassifiziert, werden. Die Einbindung von Produktmerkmalen, die von zellinternen oder zellexternen Konzentrationen abhängig sind, kann dann über die zusätzliche Integration von quantitativen oder qualitativen Zusammenhängen zwischen Konzentrationen und Produktmerkmalen erfolgen.
Weitere Gegenstände sind ein Computerprogram oder eine Software zur Durchführung des erfindungsgemäßen Verfahrens.
Das mit dem erfindungsgemäßen Verfahren bereitgestellte Modell kann zur Prozessführung oder Planung der Prozessführung sowie Untersuchung des Verfahrens im Reaktor verwendet werden.
Eine besondere Ausführungsform des erfindungsgemäßen Verfahrens wird exemplarisch beschrieben, ohne sich darauf zu beschränken. Mit diesem Verfahren wurde ebenfalls exemplarisch ein Modell einer Fermentation mit Hybridoma-Zellen erstellt und dessen Validität wie beschrieben geprüft.
Beispiel: Modellierung einer Hybridoma-Zellkultur
1 Schritt a)
Das Hintergrundwissen in Form eines metabolischen Netzwerkes wurde aus der Veröffentlichung von Niu et al. (Metabolie pathway analysis and reduetion for mammalian cell cultures— Towards macroscopic modeling. Chemical Engineering Science (2013) 102, S. 461^173. Dül: 10.1016/j.ces.2013.07.034.) entnommen. Das hier beschriebene metabolische Netzwerk einer tierischen Zelle enthält 35 Reaktionen, die 37 interne und externe Metaboliten verknüpfen (siehe Abbildung 5, siehe Tabelle 1).
Tabelle 1 : Reaktionen des Metabolischen Netzwerkes nach Niu et al. (Metabolie pathway analysis and reduetion for mammalian cell cultures— Towards macroscopic modeling. Chemical Engineering Science (2013) 102, S. 461-473. DOI: 10.1016/j.ces.2013.07.034.)
1 Glucose -» 1 G6P
1 G6P + 2 NAD -> 2 Pyruvate
1 Pyruvate -> 1 Lactate + 1 NAD
1 Pyruvate -» 1 Pyruvate m
1 NADm + 1 Pyruvate_m -» 1 Acetyl coA_m
1 Acetyl coA_m + 1 NADm + 1 Oxaloacetat_m -» 1 a-ketoglutarate_m 1 -ketoglutarate_m + 1 NADm -> 1 Succinyl CoA_m
1 FADm + 1 Succinyl CoA m -> 1 Fumarate
1 Fumarate -> 1 Malat_m
1 Malat_m + 1 NADm -> 1 Oxaloacetat_m
1 Glutamine -> 1 Glutamate + 1 NH3
1 Glutamate + 1 NADm ->1 a-ketoglutarate_m + 1 NH3
1 Malat_m -> 1 Malate
1 Malate + 1 NAD -> 1 Pyruvate
1 Glutamate + 1 Pyruvate -> 1 a-ketoglutarate_m + 1 Alanine
1 Glutamate + 1 Oxaloacetat m -> 1 α-ketoglutarate m + 1 Aspartate
1 Arginine + 2 NADm -> 1 Glutamate + 3 NH3
1 Asparagine -> 1 Aspartate + 1 NH3
2 Glycine + 1 NADm -» 1 NH3
1 Histidine + 1 NADm ->1 Glutamate + 2 NH3
1 Isoleucine + 2 NADm -> 1 Acetyl coA_m + 1 NH3 + 1 Succinyl CoA_m
1 Leucine + 3 NADm -» 3 Acetyl coA_m
1 Lysine + 6 NADm -» 2 Acetyl coA_m
1 Methionine + 4 NADm -> 1 NH3 + 1 Succinyl CoA_m
1 NADm + 1 Phenalanine -» 1 Tyrosine
1 Serine -» 1 NH3 + 1 Pyruvate
1 NADm + 1 Threonine->l NH3 + 1 Succinyl CoA_m
19 NADm + 1 TRP -> 3 Acetyl coA_m
5 NADm + 1 Tyrosine ->2 Acetyl coA_m + 1 Fumarate
5 NADm + 1 Valine - 1 NH3 + 1 Succinyl CoA_m
1 NADm ->1 NAD
0.5 Oxygen(02) - 1 NADm
1 NADm - 1 FADm
0.0156 Alanine + 0.0082 Arginine + 0.0287 Aspartate + 0.0167 G6P + 0.0245 Glutamine + 0.0039 Glutamate + 0.0196 Glycine + 0.0038 Histidine + 0.0099 Isoleucine + 0.0156 Leucine + 0.0119 Lysine + 0.0039 Methionine + 0.0065 Phenalanine + 0.016 Serine + 0.0094 Threomne + 0.0047 Tyrosine + 0.0113 Valine - 1 X (Biomasse) + 0.0981 NAD
0.01101 Alanine + 0.005033 Arginine + 0.007235 Asparagine + 0.0081787 Aspartate + 0.010381 Glutamine + 0.010695 Glutamate + 0.01447 Glycine + 0.0034602 Histidine + 0.005033 Isoleucine + 0.014155 Leucine + 0.01447 Lysine + 0.002831 1 Methionine + 0.007235 Phenalanine + 0.026738 Serine + 0.016043 Threonine + 0.0084932 Tyrosine + 0.018874 Valine -» 1 IgG (Antibody)
In der Veröffentlichung ist die Reversibilität der Reaktionen nicht explizit angegeben. Stattdessen wurden die Daten zur„Metabolischen Flussanalyse" aus derselben Veröffentlichung ausgewertet und zur Identifizierung der irreversiblen Reaktionen verwendet.
Mit der stöchiometrischen Matrix N, die die Stöchiometrie, d. h. die stöchiometrischen Koeffizienten, der internen Metabolite enthält, und den Informationen über die Reversibilität der Reaktionen wurden mit Hilfe von METATOOL 5.1 (Pfeiffer et al. METATOOL: for studying metabolic networks, Bioinformatics 199915 (3), S. 251-257.) alle Elementary Modes (EMs) des Netzwerkes berechnet. Die Zahl der EMs liegt hier bei über 300.000.
2 Schritt b)
Die Matrix mit den berechneten EMs E wurde in Schritt a) gewonnen. Analog zur Matrix N, enthält die Matrix Np die Stöchiometrie, d. h. die stöchiometrischen Koeffizienten, der externen Metabolite. Mögliche Makroreaktionen des stöchiometrischen Netzwerks wurden in der Matrix K mit Formel 21 zusammengefasst:
K = Np - E ( Formel 21)
3 Schritt c)
Die Messdaten des Prozesses wurden aus Baughman et al. entnommen, der verschiedene Messgrößen einer Fermentation von Hybridoma-Zellen über den Verlauf eines Batch-Prozess (vgl. Abbildung 6) angibt [On the dynamic modeling of mammalian cell metabolism and mAb production. In: Computers & Chemical Engineering (2010) 34 (2), S. 210-222.]. Die Messdaten wurden in das Verfahren eingegeben.
4 Schritt d)
Mit Hilfe von Spline-interpolierten Messwerten aus c) (Cm ) wurden die Wachstums- und Absterberaten sowie die spezifischen Aufnahme- bzw. Ausscheideraten berechnet (vgl. Abbildung 7). Die Lysis wurde mit einem vordefinierten und in das Verfahren eingegebenen, über den Prozessenzeitraum konstanten Lysis-Faktor A' j = 0.1 einbezogen. Ein Shiften der Messdaten war nicht notwendig, da es sich hier um Daten eines Batch-Prozesses ohne weitere Zugaben handelt. Die Daten zeigen dementsprechend einen steigenden Verlauf da alle Konzentrationsänderungen von den Zellen- und nicht von Zugaben verursacht sind.
Zur Berechnung der Raten q werden zusätzliche Informationen herangezogen. So konnte mit Hilfe der
[C—τγιοΐΛ — -— J) und der gesamten Zellzahl ein durchschnittlicher C-mol Gehalt von fc-moi,x = 18.41 errechnet werden. Die C-mol bezogene Wachstumsrate konnte nun mit Formel 22:
mol mol
μ μ fc C-mol,Xv ( Formel 22)
h■ W9cells W9cells
berechnet werden. Analog kann die C-mol bezogene Bildungsrate des Antikörpers abgeschätzt werden. Hierzu wurde die molare Zusammensetzung des Antiköpers zu CHi 5800 31N0 27S0 004 mit einer formellen molaren Masse von MmAb c_mol = 22.45 -^ geschätzt. Hier wird angenommen, dass die molare Zusammensetzung einer durchschnittlichen molaren Zusammensetzung von Proteinen wie von Villadsen et al. [Bioreaction engineering principles (201 1), Kapitel 3, Elemental and Redox Balances, S. 73, Springer Verlag, ISBN : 978-1 -4419-9687-9J angegeben entspricht. Die molare Masse des gesamten Antikörpers wurde zu MmAb = 150.000 -^- geschätzt. Die Bildungsrate des Antikörpers ergab sich dann aus der Formel:
C— mol 10~4mol M, mAb,C-mol
QmAb QmAb w4 ( Formel 23)
h■ 109cells h■ 109cells MmAb
Der zeitliche Verlauf der Raten q (t) konnte dann zur Auswahl der Makroreaktionen herangezogen werden.
5 Schritt e)
Im Schritt e) wurde ein EM-Subset von Makroreaktionen erzeugt, mit dem der Datensatz bestmöglich wiedergegeben wurde. Hierzu wurde die Matrix K aus Schritt b) benötigt. Da die Anzahl von über 300.000 Makroreaktionen zu einer zu großen Anzahl von möglichen Kombinationen geführt hätte, wurde zuerst eine datenabhängige Vorreduktion durchgeführt:
Hierzu wurden die in Schritt d) ermittelten Raten q(t) zur Berechnung der Ausbeutekoeffizienten Ym für alle Kombinationen zweier externer Metabolite verwendet. Die untere Grenze eines Ausbeutekoeffizienten Ytj wurde so gewählt, dass 99% der ermittelten Ausbeutekoeffizienten Yj {t) über diesem Wert liegen. Die obere Grenze wurde so gewählt, dass 99% der ermittelten Ausbeutekoeffizienten Y^it) unter diesem Wert liegen. Exemplarisch sind einige ermittelten Grenzen sowie der Anteil der EMs, dessen Ausbeutekoeffizienten Y™ innerhalb dieser Grenzen liegen in Tabelle 2 angegeben. Insgesamt konnte die Anzahl der EMs so auf ca. 3000 reduziert werden.
Tabelle 2: Externe Metabolite, deren maximale und minimale Ausbeutekoeffizienten sowie der Anteil der EMs, dessen Ausbeutekoeffizienten innerhalb der angegebenen Grenzen liegen
Im Anschluss an die datenabhängige Reduktion wurde anschließend noch eine datenunabhängige Reduktion durchgeführt. Hierbei wurde ein Maximalwert für die Kosinus -Ähnlichkeit zweier EMs von 0,995 definiert. Beginnend mit der ersten Reaktion wurden so alle Makroreaktionen aus der Matrix K entfernt, die diesen Wert überschritten. Es verblieben ca. 500 Makroreaktionen aus der Matrix K (auch reduzierte Matrix K genannt), die weiterhin mehr als 95% des Volumens des durch die ca. 3000 EMs aufgespannten Lösungsraum abdecken.
Vor dem Auswahlprozess erfolgte zudem ein Abgleich der im metabolischen Netzwerk nach Niu et al. angegebenen Komponenten (die den externen Metaboliten des metabolischen Netzwerkes aus a) entsprechen) und den gemessenen Komponentenkonzentrationen aus c). Bis auf Prolin werden alle von Baughman et al. gemessenen Konzentrationen auch im metabolischen Netzwerk nach Niu et al. berücksichtigt. Um die Messung der Prolinkonzentration nutzen zu können wäre entweder ein anderes vereinfachtes Netzwerk, das Prolin als externen Metaboliten enthält, verwendbar oder eine Erweiterung des bestehenden metabolischen Netzwerks möglich.
Komponenten, die zwar in den berechneten Makroreaktionen vorkamen, aber von denen keine Daten verfügbar waren, wurden ebenfalls im Folgenden ignoriert. Die entsprechenden Zeilen der Matrix K wurden dementsprechend aus der Matrix gestrichen. Das Streichen der entsprechenden Zeilen bedeutet nicht, dass diese Ein- bzw. Ausgänge von der Zelle nicht genutzt werden. Sie bestehen weiterhin im metabolischen Netzwerk, nur fehlen hierzu Messungen, die diese bilanzierbar machen können. In diesem Beispiel wurden die Ein- bzw. Ausgänge von Arginine, Glutamate, Glycin, Histidin, Leucin, Lysine, Methionin, Ammonium, Sauerstoff, Phenylalanin, Serin, Threonin, Tryptophan, Tyrosin und Valin vernachlässigt. In der folgenden Schritten des Verfahrens wird nun die reduzierte Matrix K - die die das Hintergrundwissen abbildet - und die Raten q(t) aus d) sowie die Messdaten aus c) - die das
Prozesswissen bilden - dazu verwendet, ein möglichst kleines Subset L der Makroreaktionen aus K zu gewinnen.
Es wurde die erfindungsgemäße „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" als Gütekriterium verwendet.
Analog zu den Raten q(i) wurden die Messwerte der Zellzahl und des Antikörpers hier auf C-mol normiert. Dies ist erforderlich damit die Dimension der Makroreaktionen mit denen der Messwerte übereinstimmt.
Die Auswahl des Subsets erfolgte mit einem genetischen Algorithmus. In der Berechnung der Zielfunktion dieses genetischen Algorithmus wurde hierbei das bei der „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" angesprochene lineare Optimierungsproblem gelöst. Die hier errechnete finale Summe der kleinsten Fehlerquadrate des linearen Optimierungsproblems war zugleich der Wert der Zielfunktion für die jeweilige Auswahl der Makroreaktionen.
Zur Auswahl der Größe des Subsets L aus K wurde die Optimierung mit wiederholt mit einer unterschiedlichen Anzahl an Makroreaktionen in L durchgeführt. Die Anzahl stellt einen Kompromiss zwischen der Modellkomplexität und der Wiedergabegenauigkeit dar. Um zu ermitteln, wie viele Reaktionen ausreichend für die Wiedergabe sind, kann entweder die Auswahl des Subsets L für eine variierende Anzahl von Makroreaktionen wiederholt werden oder ein Strafterm für die Reaktionsanzahl wird direkt der Zielfunktion des genetischen Algorithmus hinzugefügt werden. In diesem Fall wurden mehrere Optimierungen mit einer vordefinierten Anzahl von Makroreaktionen (10, 7, 5, 4 und 3) durchgeführt. Die mit dem genetischen Algorithmus gefundene kleinste Summe der Fehlerquadrate ist in Abbildung 9 gegenüber der Zahl der Makroreaktionen aufgetragen. Es zeigte sich, dass in diesem Fall weniger als sieben Makroreaktionen zu wenig sind um den Prozessverlauf hinreichend gut darzustellen. Die gewählten Makroreaktionen sind in Tabelle 3 angegeben.
Tabelle 3: Gewähltes Subset der Makroreaktionen (I). Nicht unterstrichene Komponenten werden in dem Modell nicht berücksichtigt, da hierzu keine Messungen vorliegen.
0.474 Alanine + 0.474 Methionine
-* 0.158 Asparag ine + 0.316 Aspartate + 0.632 Glycine
+ 0.158 Tryptophan + 0.00789 Arginine + 0.0304 Asparagine + 0.0161 Glucose
+ 0.0236 Glutamine + 0.00375 Glutamate + 0.00366 Histidine + 0.00953 Isoleucine + 0.015 Leucine + 0.112 Methionine
+ 0.00626 Phenalanine + 0.0154 Serine + 0.0109 Valine
-> 0.963 X (Biomass) + 0.00276 Aspartate + 0.24 Glycine
+ 0.0208 Tryptophan ine + 0.147 Glutamate
-> 0.295 Aspartate + 0.885 Glycine + 0.147 Lactate ne + 0.113 Asparagine + 0.0603 Glucose + 0.0225 Glutamine
+ 0.0824 Histidine + 0.00909 Isoleucine + 0.00597 Phenalanine + 0.0216 Tryptophan + 0.00431 Tyrosine + 0.0104 Valine
-> 0.918 X (Biomass) + 0.061 Alanine + 0.0865 Aspartate
+ 0.343 Glycine + 0.0631 Methionine e + 0.412 Aspartate + 0.00991 Glucose + 0.0145 Glutamine
+ 0.554 Glycine + 0.00226 Histidine + 0.00588 Isoleucine
+ 0.00926 Leucine + 0.00706 Lysine + 0.0649 Phenalanine + 0.0095 Serine + 0.00671 Valine
-* 0.594 X (Biomass) + 0.049 Alanine + 0.395 Asparagine
+ 0.0503 Threonine + 0.0388 Tryptophan
0.0077 Ar ginine + 0.179 Aspartate + 0.0157 Glucose + 0.104 Glutamine
+ 0.216 Glycine + 0.00357 Histidine + 0.00929 Isoleucine
+ 0.0146 Leucine + 0.0112 Lysine + 0.03Q Tyrosine + 0.0106 Valine
-* 0.939 X (Biomass) + 0.0624 Alanine + 0.152 Asparag ine
+ 0.0183 Tryptophan
0.0342 Arginine + 0.211 Aspartate + 0.00762 Glucose + 0.0195 Glutamine
+ 0.244 Glycine + 0.00452 Histidine + 0.0546 Isoleucine
+ 0.0185 Leucine + 0.0171 Lysine + 0.00406 Methionine
+ 0.0178 Tyrosine + 0.0203 Valine
-> 0.457 X (Biomass) + 0.804 IgG (Antibody) + 0.185 Asparagine
+ 0.0153 Tryptophan
In den gezeigten Makroreaktionen sind alle externen Metaboliten des metabolischen Netzwerkes aus a) angegeben. Teil des Modells sind jedoch nur die unterstrichenen externen Metabolite, da nur für diese Messdaten aus c) vorliegen.
6 Schritt f)
Für das gewählte Set der Makroreaktionen wurden die Reaktionsraten über der Zeit ermittelt. In diesem Beispiel wurden mit der erfindungsgemäßen Methode „lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" die in Abbildung 10 gezeigten Messwerte durch eine Schätzung der Reaktionsraten r(t) approximiert. Das Ergebnis der Methode ist ein stückweise linearer Verlauf der einzelnen (volumetrischen) Reaktionsraten. Durch Division mit dem interpoliertem Verlauf der vitalen Zellzahl Xv(t) wurden die zellspezifischen Reaktionsraten r(t) der in Tabelle 3 gezeigten Makroreaktionen gewonnen. Die so erhaltenen Reaktionsraten r(t) sind in Abbildung 10 abgebildet.
7 Schritt g)
Für alle in Tabelle 3 gezeigten Makroreaktionen wurden generische Kinetiken gemäß Formel 24 angenommen: (Formel 24)
In diesem Fall wurden sie durch eine Monod-Kinetik realisiert D. h. das bei jeder Reaktion k für jedes Substrat t eine Limitierung gemäß Formel 25 :
(Formel 25)
eingeführt wurde. Hierbei ist rk max die maximale Reaktionsgeschwindigkeit, Ni die Anzahl der berücksichtigten Limitierungen, Q die Konzentration der Komponente t, Km k i die dazugehörigen Monod-Konstanten und n; der Hill-Parameter für die Reaktionsordnung darstellen. Deren Werte werden in den Schritten h) und j) angepasst.
Weitere Terme ergeben sich aus der Analyse der Reaktionsraten r(t) aus f). In diesem Beispiel wurden neben Substrathmitierungen auch Inhibierungen gemäß Formel 26 berücksichtigt.
(Formel 26)
Auch für diese Limitierung mussten die Werte der Parameter Kl k ii und ηέ angepasst werden. Die verwendeten kinetischen Terme der Reaktionen sind in Tabelle 4 angegeben.
Tabelle 4: Kinetische Terme der gewählten Makroreaktionen aus L
[Glc] (t) \ ( [Gln] (t) \ ( [Asn] (t) r2 (t = r2,r.
Km,Gic,2 + [Glc] (t ) KmiGlni2 + [Gln] (t) J KmiAsrii2 + [Asn] (t)
[Ala] (t)
Km,Ala,2 + [Ala\ (t) j
[Glc] (t) \ ( [Asn] (t) \ ( KliLaCi3 r3 (t = r3,;
Km,Gic,3 + [Glc] (t I KmiGlni3 + [Asn] (t) l KIiLaCi3 + [Lac] (t
8 Schritt h)
Für jede Reaktionsrate konnte mit der in Tabelle 4 angegebenen Kinetik und den interpolierten Werten der in der Kinetik berücksichtigten Konzentrationen Cm (t) der Verlauf der Reaktionsrate algebraisch berechnet werden.
Die Anpassung der Parameter dieser Kinetiken erfolgte separat für jede Reaktion t an die in Schritt f) ermittelte Reaktionsrate ri(t). Die Zielfunktion für eine Optimierung der in Reaktion t vorkommenden Parameter lautete in diesem Beispiel:
(Formel 27)
Die so angepassten Verläufe aller berechneten rk pk, C_int (t) J sind zusammen mit den entsprechenden rk (t) in Abbildung 11 dargestellt. Die Verläufe der ersteren sind gestrichelt, die der letzteren durchgehend dargestellt. Es ist erkennbar, dass der Verlauf qualitativ übereinstimmt. Dies bedeutet, dass mit den gewählten Kinetiken auch die Dynamik des Prozesses zufriedenstellend wiedergegeben werden kann. Diese Information ist in diesem Modellierungsschritt sehr hilfreich, da bei einer nicht zufriedenstellenden Wiedergabe die schnell durchzuführenden Schritte g) (Auswahl anderer Kinetiken) und h) (Parameterwertschätzung) wiederholt werden können, bis der gewünschte Anpassungsgrad erreicht ist. Schritt i) war hier also nicht erforderlich. 9 Schritt j)
Die weiterführende Anpassung der Modelparameterwerte p erfolgte mit den Messdaten aus c). Hierfür wurden alle Parameter zugleich optimiert. Zudem wurden die bisher nicht betrachteten Vorgänge Apoptose und Lysis mit einbezogen. Diese werden in den Differentialgleichungen, die die Entwicklung der vitalen- und gesamten Zellzahl beschreiben, benötigt:
dX1
μά)■ Xv (Formel 28)
dt
1 (Formel 29)
Die gewählte Kinetik für die Beschreibung der Apoptose lautete:
([Lac] (t) - CLaCiCr)
d( = μα, ,max , [Lac] > C • lhac,cr (Formel 30)
Kd,Lac + ([^ac]( — ^Lac,cr)
Md( = 0 , [Lac] < CLaCi ,cr (Formel 31 )
Die Lysisrate Ki wurde als konstant über dem Prozess angenommen. Neben den Parametern der Reaktionsraten wurden in diesem Schritt die durch die Apoptose und Lysis eingeführten Parameter Ciac.cr (kritische Laktat Konzentration), μα,τηαχ (maximale Absterberate), Kd Lac (Monod-Parameter zur Beschreibung des Einflusses der Laktat-konzentration auf die Absterberate) und Kt (Lysisrate) bestimmt. In dem Beispiel wurde ausgehend von den Startwerten des Datensatzes durch numerisches Lösen des ÜDE-Systems der Verlauf der geschätzten Konzentrationen C(t) ermittelt. Die Differenz zwischen den gemessenen Konzentrationen C_m (t) und den geschätzten Konzentrationen C (t) wurde dabei mit üblichen Methoden mit folgender Zielfunktion minimiert:
(Formel 32)
Mit insgesamt 33 Parametern p ist diese Optimierung i. d. R schwierig durchzuführen, da die
Zielfunktion viele lokale Optima besitzt. Startet man einen deterministischen Optimierungsalgorithmus, wie z. B. den Levenberg-Marquardt Algorithmus an den aus Schritt h) bekannten Startwerten der Parameter ist die Erfolgsaussicht hingegen stark erhöht. Der angepasste Prozessverlauf ist in Abbildung 12 abgebildet. Die angepassten Parameter sind in Tabelle 5 abgebildet. Tabelle 5: Parameter der Kinetiken sowie der Apoptose und Lysis
10 Schritt k)
Das Modell, bestehend aus der Matrix L, den Kinetiken aus Tabelle 4 sowie den Kinetiken der Apoptose mit den mit den dazugehörigen Parameterwerten aus Tabelle 5 wurde ausgegeben.
Symbolverzeichnis
(Unterstrich) Bezeichnet einen Vektor
i (Index ) Bezeichnet das i-te Element eines Vektors
^ (Index k) Bezeichnet das k-te Element eines Vektors
[ ] Bezeichnet die Konzentration der in der Klammer stehenden Komponente
C Konzentration AC Konzentrationsdifferenz
Interpolierte Konzentration
c Geschätzte Konzentration (z.B. durch das Lösen einer Differentialgleichung)
Cs Geshiftete Konzentration
r Kritische Konzentration
Cm Gemessene Konzentration
D Verdünnungsrate
q Ermittelte zellspezifischen Ausscheide- und Aufnahmerate
q Ermittelte zellspezifischen Ausscheide- und Aufnahmerate, die von einer beliebigen Einheit auf — : ——\ umgerechnet wurde
° LZeit Zellzahll °
r Ermittelte Reaktionsrate
f Geschätzte Reaktionsrate (z.B. durch das Berechnen einer Reaktionskinetik) f Limitierung einer Kinetik
Parameter einer Reaktionskinetik
N Stöchiometrische Matrix
Np Externe stöchiometrische Matrix
K Matrix, die Makroreaktionen enthält
E Matrix, die alle Elementary Modes enthält
xt Gesamtzellzahl
Xv Vitale Zellzahl
μ Wachstumsrate
Absterberate
μ Wachstumsrate, die von einer beliebigen Einheit auf
\Sto f fmenqe~\ . , .
umgerechnet wurde
LZ eit-Z eilzahl 1 a
Kd Lysis-rate
Ki Parameter einer Inhibierungs-Limitierung
KM Parameter einer Substratlimitierung
n Hill-Parameter einer Inhibierungs- oder Substratlimitierung
L Subset der Makroreaktionen, das für das Modell verwendet wird
V Modellparameter
S Substrat SSRq Summe der quadrierten Residuen der spezifischen Aufnahme- oder Abgaberaten
SSRC Summe der quadrierten Residuen der Konzentration
SSRr Summe der quadrierten Residuen der Reaktionsraten
Beschreibung der Abbildungen:
Abbildung 1 zeigt das Shiften von Messdaten: Dargestellt ist der tatsächliche Verlauf einer gemessenen Größe (Ct (t)), der sich bei Änderungen der Verdünnungsrate (D(t)) sprunghaft ändert. Der geshiftete Verlauf ( i s(t)) kommt nur durch von der Zelle verursachte Änderungen zu Stande.
Abbildung 2 zeigt die Flux Map zweier spezifischer Raten qx und q2. Die Höhenlinien geben die Häufigkeit an, mit der die jeweilige Kombination der Raten in den gemessenen Daten vorkommt.
Abbildung 3 zeigt eine dreidimensionale Darstellung des Lösungsraumes, der durch eine positive Linearkombination von EMs aufgespannt wird. In schwarz ist der Lösungsraum des gesamten Sets, in grau der eines Subsets dargestellt.
Abbildung 4 zeigt die Flux Map zweier spezifischer Raten qi und q2. Als Vektoren sind die 2- dimensionalen Projektionen der Makroreaktionen eines Sets I dargestellt.
Abbildung 5 zeigt eine schematische Darstellung des Metabolischen Netzwerkes aus Niu et al. Hierbei ist die Begrenzung der Zelle als Kasten gezeigt. Die zellinterne Abgrenzung des Mitochondriums ist mit einer gestrichelten Linie gekennzeichnet. Externe Komponenten sind mit dem Index „xt" gekennzeichnet. Die Pfeile und gepunkteten Pfeile kennzeichnen Reaktionen.
Abbildung 6 zeigt die Messdaten einer Fermentation mit Hybridoma-Zellen aus Baughman et al. Die Gesamtzellzahl (total Cells) ist hierbei aus der Summe der lebenden Zellen (vital Cells) und toten Zellen (dead Cells) berechnet. Die Abkürzungen GLC, GLN, ASP, ASN, LAC, ALA und PRO bezeichnen das Substrat Glucose und die Aminosäuren Glutamin, Asparaginsäure, Asparagin, Alanin und Prolin sowie das Stoffwechselprodukt Laktat. Die Abkürzung MAB bezeichnet das Produkt monoklonaler Antikörper und BM die Biomasse. Abbildung 7 zeigt die Wachstums- und Absterberaten sowie die zellspezifischen Aufnahme- und
[ mM 1
h io9ceiis \ an8e8erjen- Die Rate
1« n-
Abbildung 8 zeigt die mit der„linearen Schätzung von Reaktionsraten ausgewählter Makroreaktionen" mit dem gewählten Reaktionsset approximierten Konzentrationen. Die Gesamtzellzahl (Xt) sowie die Antikörperkonzentration (MAB) wurden dafür auf C-mol umgerechnet.
Abbildung 9 zeigt die kleinste ermittelte Summe der Fehlerquadrate ("Minimum error") aufgetragen über der Zahl der Makroreaktionen im Subset (nR).
Abbildung 10 zeigt die mit der erfindungsgemäßen Methode "lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" ermittelten Reaktionsraten der Makroreaktionen r(t).
Abbildung 1 1 zeigt die mit der erfindungsgemäßen Methode "lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen" ermittelten Reaktionsraten der Makroreaktionen r(t) (durchgehende
Linie) zusammen mit den algebraisch berechneten Reaktionsraten r (gestrichelte Linie)
Abbildung 12 zeigt einen Vergleich der gemessenen Konzentrationen Cm(i) (Punkte) und dem simuliertem Prozessverlauf (i) (durchgehende Linie). Die Konzentrationen sind in [mM] angegeben. Ausnahmen sind die vitale und gesamte Zellzahl (Xv/Xt in [lO9 cells und die Konzentration des Antikörpers (mAb in [10-4 mM]).

Claims

Ansprüche:
1. Computerimplementiertes Verfahren zur Erstellung eines Models einer Bioreaktion mit einem Organismus, das folgende Schritte umfasst:
a. Ausgewählte Stoffwechselwege des Organismus, deren Stöchiometrie- sowie Reversibilitätseigenschaften werden als Hintergrundwissen in das Verfahren eingegeben und Elementary Modes werden aus dieser Eingabe berechnet.
b. Die Elementary Modes aus a) werden in einer Matrix K zusammengefasst, wobei die Elementary Modes die Stoffwechselwege aus a) in Makroreaktionen zusammenfassen und die Matrix K die Stöchiometrie und die Reversibilitätseigenschaften aller Makroreaktionen enthält.
c. Die Messdaten zur Bioreaktion mit dem Organismus werden eingegeben. d. Mit Hilfe einer Interpolationsmethode werden auf Basis der eingegebenen Messdaten aus c) die für den Organismus spezifischen Raten - Ausscheide- und Aufnahmeraten von einer oder mehreren Eingangsgrößen und Ausgangsgrößen - der eingegebenen Stoffwechselwege berechnet.
e. Relevante Makroreaktionen werden in Form eines Subsets der Elementary Modes aus a) ausgewählt durch
i. datenunabhängige und / oder datenabhängige Vorreduktion der Anzahl der Elementary Modes aus a),
ii. Auswahl des Subsets aus der Vorreduktion aus e) i. mit den Messdaten aus c) und / oder einer oder mehrerer Raten aus d) mittels eines Algorithmus nach einem mathematischen Gütekriterium und Zusammenfassung des Subsets in einer Matrix L.,
iii. Optional erfolgt eine Darstellung des Subsets grafisch.
f. Mit Hilfe einer Interpolationsmethode werden auf Basis der eigegebenen Messdaten aus c) und / oder der Raten aus d), die Reaktionsraten der Makroreaktionen des Subsets r(t) berechnet.
g. Kinetiken der Makroreaktionen des Subsets aus e) ii. werden mit folgenden Zwischenschritten entworfen; dadurch werden die Modellparameter definiert.
i. Aus der Stöchiometrie der Makroreaktionen werden generische Kinetiken entworfen. ii. Einflussgrößen auf die Makroreaktionen werden aus den Reaktionsraten aus f) ermittelt.
iii. Die generischen Kinetiken aus g) i. werden um Terme erweitert, die die in g) ii. ermittelten Einflussgrößen quantifizieren.
h. Optional erfolgt aus den Kinetiken aus g) eine erste Anpassung der Modellparameterwerte für jede Makroreaktion separat an die berechneten Reaktionsraten aus f) für jede Makroreaktion separat.
i. Optional werden die Schritte g) und h) wiederholt bis eine vordefinierte Anpassungsgüte erreicht wird.
j. Die Modellparameterwerte werden an die Messdaten aus c) angepasst.
k. Die Matrix L, die Kinetiken aus g) und die Modellparameterwerte aus j) bilden das
Modell und werden ausgegeben und / oder in ein Prozessführungs- oder
Prozessentwicklungsmodul übertragen.
2. Computer-implementiertes Verfahren nach Anspruch 1, wobei im Schritt d) auch Wachstumsraten, besonders bevorzugt auch Absterberaten des Organismus berechnet werden.
3. Computer-implementiertes Verfahren nach einem der Ansprüche 1 oder 2, wobei im Schritt g) eine individuelle Anpassung der Kinetiken basierend auf einer Analyse der Reaktionsraten aus f) erfolgt.
4. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 3, wobei im Schritt h) die Anpassung der Parameterwerte der Kinetiken aus g) durch Kombination von mehreren Anpassungsmethoden erfolgt.
5. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 4, wobei im Schritt e) ii. zur Auswahl des Subsets von Makroreaktionen eine lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen durchgeführt wird.
6. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 5, wobei im Schritt e) ii. zur Auswahl des Subsets von Makroreaktionen eine lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen in Kombination mit einem evolutionären Algorithmus durchgeführt wird.
7. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 6, wobei die Messdaten vor Anwendung der Interpolationsmethode im Schritt d) geshiftet werden, um die Beschreibung eines konstanten Verbrauchs ohne Feedpeaks zu erreichen.
8. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 7, wobei im Schritt f) eine lineare Schätzung von Reaktionsraten ausgewählter Makroreaktionen durchgeführt wird.
9. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 8, wobei im Schritt e) i. eine datenabhängige Vorreduktion durchgeführt wird und für diese die Methode der linearen Schätzung von Reaktionsraten ausgewählter Makroreaktionen mit N LS verwendet wird.
10. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 9, wobei im Schritt e) iii. die Validität der Auswahl des Subsets von Makroreaktionen mit Hilfe einer Flux-Map geprüft wird.
11. Computer-implementiertes Verfahren nach einem der Ansprüche 1 bis 10, wobei im Schritt e) ii. die Auswahl aus der Vorreduktion aus e) i. mit den Messdaten aus c) erfolgt.
12. Computerprogram zur Durchführung der Verfahrensschritte nach einem der Ansprüche 1 bis 10.
13. Software zur Durchführung der Verfahrensschritte nach einem der Ansprüche 1 bis 11.
EP16701791.2A 2015-01-29 2016-01-28 Computerimplementiertes verfahren zur erstellung eines fermentationsmodels Withdrawn EP3251039A1 (de)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
EP15153052.4A EP3051449A1 (de) 2015-01-29 2015-01-29 Computerimplementiertes Verfahren zur Erstellung eines Fermentationsmodels
PCT/EP2016/051753 WO2016120361A1 (de) 2015-01-29 2016-01-28 Computerimplementiertes verfahren zur erstellung eines fermentationsmodels

Publications (1)

Publication Number Publication Date
EP3251039A1 true EP3251039A1 (de) 2017-12-06

Family

ID=52434610

Family Applications (2)

Application Number Title Priority Date Filing Date
EP15153052.4A Withdrawn EP3051449A1 (de) 2015-01-29 2015-01-29 Computerimplementiertes Verfahren zur Erstellung eines Fermentationsmodels
EP16701791.2A Withdrawn EP3251039A1 (de) 2015-01-29 2016-01-28 Computerimplementiertes verfahren zur erstellung eines fermentationsmodels

Family Applications Before (1)

Application Number Title Priority Date Filing Date
EP15153052.4A Withdrawn EP3051449A1 (de) 2015-01-29 2015-01-29 Computerimplementiertes Verfahren zur Erstellung eines Fermentationsmodels

Country Status (15)

Country Link
US (2) US10296708B2 (de)
EP (2) EP3051449A1 (de)
JP (1) JP6816003B2 (de)
KR (1) KR20170109629A (de)
CN (1) CN107408161B (de)
AR (1) AR103564A1 (de)
AU (1) AU2016212059B2 (de)
BR (1) BR112017016198A2 (de)
CA (1) CA2975012C (de)
EA (1) EA035276B1 (de)
HK (1) HK1247340A1 (de)
IL (1) IL253584A0 (de)
SG (2) SG11201706166PA (de)
TW (1) TWI690813B (de)
WO (1) WO2016120361A1 (de)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11959063B2 (en) 2017-09-29 2024-04-16 Unibio A/S Optimization of fermentation processes
US11508459B2 (en) * 2018-01-31 2022-11-22 X Development Llc Modified FBA in a production network
JP7059789B2 (ja) * 2018-05-14 2022-04-26 富士通株式会社 逐次制御プログラム、逐次制御方法および逐次制御装置
EP4026072A1 (de) 2019-09-06 2022-07-13 Bayer Aktiengesellschaft System zur planung, wartung, führung und optimierung eines produktionsprozesses
US20220284269A1 (en) * 2021-03-03 2022-09-08 Lanzatech, Inc. System for control and analysis of gas fermentation processes

Family Cites Families (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005506840A (ja) * 2001-10-01 2005-03-10 ディヴァーサ コーポレイション リアルタイム代謝フラックス分析を用いた全細胞操作
US20070048732A1 (en) * 2005-08-30 2007-03-01 Hepahope, Inc. Chemosensitivity tester
KR100718208B1 (ko) * 2006-04-21 2007-05-15 한국과학기술원 Ｃｒｄ 및 ｓｒｄ를 이용한 대사흐름 분석방법
US20080103747A1 (en) * 2006-10-31 2008-05-01 Macharia Maina A Model predictive control of a stillage sub-process in a biofuel production process
US7831318B2 (en) * 2006-10-31 2010-11-09 Rockwell Automation Technologies, Inc. Model predictive control of fermentation temperature in biofuel production
US8571690B2 (en) * 2006-10-31 2013-10-29 Rockwell Automation Technologies, Inc. Nonlinear model predictive control of a biofuel fermentation process
US8571689B2 (en) * 2006-10-31 2013-10-29 Rockwell Automation Technologies, Inc. Model predictive control of fermentation in biofuel production
EP2350299A4 (de) * 2008-10-28 2012-04-25 Univ Rice William M Mikroaerobe kulturen zur umwandlung von glycerin in chemische stoffe
US8507233B1 (en) * 2009-06-30 2013-08-13 Nanobiosym, Inc. NanoBiofuel production using nanoscale control methods and materials
WO2011112260A2 (en) * 2010-03-11 2011-09-15 Pacific Biosciences Of California, Inc. Micromirror arrays having self aligned features
AU2011238099A1 (en) * 2010-04-07 2012-11-29 Novacare Computer based system for predicting treatment outcomes
EP2668281A4 (de) * 2011-01-25 2015-08-26 Cargill Inc Zusammensetzungen und verfahren zur succinatherstellung
WO2012168988A1 (ja) * 2011-06-10 2012-12-13 パナソニック株式会社 抗体を自己組織化膜上に固定する方法
CN103044542B (zh) * 2012-08-17 2015-08-19 常熟理工学院 鲫鱼卵中丝氨酸蛋白酶抑制剂及其基因和应用
US10507235B2 (en) * 2013-03-19 2019-12-17 Globeimmune, Inc. Yeast-based immunotherapy for Chordoma
CN103413066B (zh) * 2013-08-28 2016-12-28 南京工业大学 自吸式反应器放大发酵培养酵母的预测方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
HYUN-SEOB SONG ET AL: "Reduction of a set of elementary modes using yield analysis", BIOTECHNOLOGY AND BIOENGINEERING, vol. 102, no. 2, 22 July 2008 (2008-07-22), pages 554 - 568, XP055618777, ISSN: 0006-3592, DOI: 10.1002/bit.22062 *
LUKAS HEBING ET AL: "Efficient Generation of Models of Fed-Batch Fermentations for Process Design and Control", IFAC-PAPERSONLINE-8TH IFAC SYMPOSIUM ON ADVANCES IN AUTOMOTIVE CONTROL AAC 2016, NORRKÖPING, SWEDEN, 20-23 JUNE 2016, vol. 49, no. 7, 9 August 2016 (2016-08-09), pages 621 - 626, XP055618927, ISSN: 2405-8963, DOI: 10.1016/j.ifacol.2016.07.237 *
See also references of WO2016120361A1 *
ZITA I.T.A. SOONS ET AL: "Selection of Elementary Modes for Bioprocess Control", IFAC THE 2012 IFAC WORKSHOP ON AUTOMATIC CONTROL IN OFFSHORE OIL AND GAS PRODUCTION, vol. 43, no. 6, 31 December 2010 (2010-12-31), Red Hook, NY, pages 156 - 161, XP055618613, ISSN: 1474-6670, ISBN: 978-1-123-47890-7, DOI: 10.3182/20100707-3-BE-2012.0019 *

Also Published As

Publication number Publication date
WO2016120361A1 (de) 2016-08-04
EA035276B1 (ru) 2020-05-22
JP2018511849A (ja) 2018-04-26
SG10202006972VA (en) 2020-08-28
TW201643744A (zh) 2016-12-16
CN107408161B (zh) 2021-02-26
AU2016212059A1 (en) 2017-08-17
IL253584A0 (en) 2017-09-28
US20190228835A1 (en) 2019-07-25
AU2016212059B2 (en) 2021-07-29
AR103564A1 (es) 2017-05-17
CN107408161A (zh) 2017-11-28
KR20170109629A (ko) 2017-09-29
CA2975012A1 (en) 2016-08-04
HK1247340A1 (zh) 2018-09-21
EA201791659A1 (ru) 2018-01-31
CA2975012C (en) 2021-06-15
US10872680B2 (en) 2020-12-22
EP3051449A1 (de) 2016-08-03
JP6816003B2 (ja) 2021-01-20
US20160224721A1 (en) 2016-08-04
SG11201706166PA (en) 2017-08-30
TWI690813B (zh) 2020-04-11
BR112017016198A2 (pt) 2018-04-17
US10296708B2 (en) 2019-05-21

Similar Documents

Publication Publication Date Title
WO2016120361A1 (de) Computerimplementiertes verfahren zur erstellung eines fermentationsmodels
EP3732485A1 (de) Vorhersage des metabolischen zustands einer zellkultur
Cho et al. Organic carbon, influent microbial diversity and temperature strongly influence algal diversity and biomass in raceway ponds treating raw municipal wastewater
EP2990900B1 (de) System und Verfahren zur Steuerung, Erfassung, Regelung und/oder Analyse von biologischen, biochemischen, chemischen und/oder physikalischen Prozessen durch Datenkorrelation
WO2021043712A1 (de) System zur planung, wartung, führung und optimierung eines produktionsprozesses
Kuligowski et al. Application of leather waste fractions and their biochars as organic fertilisers for ryegrass growth: Agri-environmental aspects and plants response modelling
Wolfgramm et al. The research–action interface in sustainable land management in Kyrgyzstan and Tajikistan: challenges and recommendations
Vatsa et al. A sustainable approach to improving agrifood production: getting the balance right between organic soil amendments and chemical fertilizers
Brinc et al. Optimization of process conditions for mammalian fed-batch cell culture in automated micro-bioreactor system using genetic algorithm
Fink et al. Verträge, Vertrauen und Unternehmenserfolg in Automobilclustern
Balakrishnan et al. Application of design of experiments in bioprocessing: process analysis, optimization, and reliability
De la Hoz Siegler Jr Optimization of biomass and lipid production in heterotrophic microalgal cultures
Zhang et al. Reoligotrophication of a High-Nitrogen Reservoir with Phosphorus Removal and Implications for Management
Koch Petrinetze in der Systembiologie
Schäfer et al. Integration of knowledge in inter-and transdisciplinary research projects: use of constellation analysis in a project of sustainable land use management.
Amaro TUBULAR PHOTOBIOREACTORS ON MICROALGAE PLANT
Babor Application of nature-inspired optimization algorithms to improve the production efficiency of small and medium-sized bakeries
Fryxell et al. Algal blooms as a reactive dynamic response to seasonal perturbation in an experimental system
Koffel New insights from niche theory on plant adaptation and ecosystem functioning
EP2513293A1 (de) Methode zur optimierung eines biopharmazeutischen produktionsprozesses
Valecillos et al. The approach of intelligent organizations in implementing new technologies for managing small and medium enterprises (SMEs)
Pieralli The role of soils in production: Aggregation, separability, and yield decomposition in Kenyan agriculture
Mailier et al. A Fast and Systematic Procedure to Develop Dynamic Models of Bioprocesses-Application to Microalgae Cultures
Spencer Final Report Submitted by NIMBioS Sabbatical Visitor: Matthew Spencer School of Environmental Sciences, University of Liverpool, Liverpool, UK Sabbatical dates: September 2012-January 2013
Mehrabi The number of farms globally could halve by the end of the century

Legal Events

Date Code Title Description
STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: THE INTERNATIONAL PUBLICATION HAS BEEN MADE

PUAI Public reference made under article 153(3) epc to a published international application that has entered the european phase

Free format text: ORIGINAL CODE: 0009012

STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: REQUEST FOR EXAMINATION WAS MADE

17P Request for examination filed

Effective date: 20170829

AK Designated contracting states

Kind code of ref document: A1

Designated state(s): AL AT BE BG CH CY CZ DE DK EE ES FI FR GB GR HR HU IE IS IT LI LT LU LV MC MK MT NL NO PL PT RO RS SE SI SK SM TR

AX Request for extension of the european patent

Extension state: BA ME

DAV Request for validation of the european patent (deleted)
DAX Request for extension of the european patent (deleted)
STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: EXAMINATION IS IN PROGRESS

17Q First examination report despatched

Effective date: 20190913

STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: THE APPLICATION IS DEEMED TO BE WITHDRAWN

18D Application deemed to be withdrawn

Effective date: 20200603