dGSMEF: A dynamic genome-scale modelling of E. Coli fermentation...
In Kange2014/dGSMEF: A dynamic genome-scale modelling of E. coli fermentation

Description Usage Arguments Details Value See Also Examples

dGSMEF performs a dynamic flux balance analysis under batch, fed-batch and induction cultures.

dGSMEF(model, substrateRxns, initConcentrations, initBiomass, feedSubstrateRxns,
  feedConcentrations, initVolume = 1, initpH = 7.25, u_fix = 0.15,
  yield_rate = 0.6586, x_ind = 51.6, ngam = 11.36, initRatio = 0.8,
  feedRate_ind = 2.9, timeStep, nSteps, exclUptakeRxns, rcd = FALSE,
  fld = FALSE, verboseMode = 2)

`model`	Sybil model structure (class `modelorg`).
`substrateRxns`	List of exchange reaction names for substrates initially in the media that may change (e.g. not h2o or co2).
`initConcentrations`	Initial concentrations of substrates (in the same structure as substrateRxns).
`initBiomass`	Initial biomass (must be non zero).
`feedSubstrateRxns`	List of exchange reaction names for substrates in the feeding solution. Default: "EX_glc__D(e)".
`feedConcentrations`	Substrate concentrations in the feeding solution (same structure as feedSubstrateRxns). Default: 2775.3 mmol/L, 500 g/L.
`initVolume`	Initial total volume of bioreactor, default: 1.
`initpH`	Initial pH in the bioreactor, default: 7.3.
`u_fix`	Fixed specific growth rate during fed-batch phase before induction.
`yield_rate`	Biomass yield rate on glucose, unit: g DW/g glucose.
`x_ind`	Biomass centration at induction.
`ngam`	Non-growth associated maintenance, to characterize the metabolic burden for the recombinant cell.
`initRatio`	Intial ratio of resources allocated to the recombinant protein synthesis in the beginning time of induction.
`feedRate_ind`	Feeding rate during induction, unit: ml/L.
`timeStep`	Time step size.
`nSteps`	Maximum number of time steps.
`exclUptakeRxns`	List of uptake reactions whose substrate concentrations do not change (Default = 'EX_co2(e)','EX_h2o(e)','EX_h(e)').
`rcd`	Boolean. Save the reduced costs. Default: FALSE.
`fld`	Boolean. Save the resulting flux distribution. Default: FALSE.
`verboseMode`	An integer value indicating the amount of output to stdout: 0: nothing, 1: status messages, 2: like 1 plus a progress indicator, 3: a table containing the reaction id's and the corresponding min max values. Default: 2.

The function is implemented with static optimization approach (SOA) First, the simulation time is divided into small periods which are assumed to be in quasi-steady state; Second, for every time step, an FBA problem is solved and the fluxes are integrated over the time period and extracellular concentrations are calculated, accordingly. Furthermore, the flux optimizaiton is achieved by parsimonius FBA (pFBA), which minimize all fluxes within the solution space of the growth optimum, to minimize the requirement for enzyme expression.

In addition, for the bi-objective optimization of both maximum growth and recombinant protein synthesis rates, a novel proteome resource allocation theory is applied. According to this theory, the growth reduction induced by heterologous protein expression is a simple consequence of proteome allocations. And more heterologous proteins, less growth rate.

returns class optsol_dGSMEF

concentrationMatrix: Matrix of extracellular metabolite concentrations.
excRxnNames: Names of exchange reactions for the EC metabolites.
timeVec: Vector of time points.
biomassVec: Vector of biomass values.
fluxMatrix: Matrix of flux distributions.
redCostMatrix: Matrix of reduced costs for each reaction.
RxnNames: Names of all reactions for the E.Coli metabolites.

modelorg, optsol_dGSMEF.

## Not run: 
library(sybil)
library(sybilSBML)
library(glpkAPI)
library(deSolve)
library(dGSMEF)
model = readSBMLmod("../data/HNC47_without_fbc.xml");
lowbnd(model)[react_id(model) == 'EX_glc__D(e)'] = -9.58; ## for HNC47 with Fc
lowbnd(model)[react_id(model) == 'EX_o2(e)']=-11.5; # for HNC47 with Fc

## Add Fc synthesis reaction into the model
model <- addReact(model,id="RecombinantProtein",
			met=c('ala__L[c]','cys__L[c]','asp__L[c]','glu__L[c]','phe__L[c]',
                   'gly[c]','his__L[c]','ile__L[c]','lys__L[c]','leu__L[c]',
                   'met__L[c]','asn__L[c]','pro__L[c]','gln__L[c]','arg__L[c]',
                   'ser__L[c]','thr__L[c]','val__L[c]','trp__L[c]','tyr__L[c]',
                   'atp[c]','adp[c]','pi[c]'),
			Scoef=c(-c(6,6,9,17,9,11,5,4,17,17,3,11,19,11,6,24,14,23,4,8),
                     -4.306*224,964.544,964.544),
			lb=0,ub=1000,obj=0);

## Define FDM recipe, while we suppose NH4,Ca2+,cl-, fe2+, and cu2+ are abundant enough
substrateRxns = c('EX_glc__D(e)','EX_nh4(e)','EX_pi(e)','EX_so4(e)','EX_mg2(e)',
                  'EX_k(e)','EX_ca2(e)','EX_cl(e)','EX_fe2(e)','EX_mn2(e)',
                  'EX_cu2(e)','EX_zn2(e)','EX_mobd(e)');
initConcentrations = c(55.56,2000,67.6,63.08,2,67.6,1000,1000,2000,0.09171,1000,0.05912,0.01078);

## EX_ni2(e) and EX_cobalt2(e) both are required to support the growth of cell
lowbnd(model)[react_id(model) == 'EX_cobalt2(e)'] = -1000;
lowbnd(model)[react_id(model) == 'EX_fe3(e)'] = 0;
lowbnd(model)[react_id(model) == 'EX_tungs(e)'] = 0;
lowbnd(model)[react_id(model) == 'EX_ni2(e)'] = -1000;
lowbnd(model)[react_id(model) == 'EX_sel(e)'] = 0;
lowbnd(model)[react_id(model) == 'EX_slnt(e)'] = 0;

## 10% YE + 25% glucose for feeding during induction
## here, suppose YE is merely composed of amino acids
AA <- c('EX_ala__L(e)','EX_arg__L(e)','EX_asn__L(e)','EX_asp__L(e)','EX_cys__L(e)','EX_gln__L(e)',
        'EX_glu__L(e)','EX_gly(e)','EX_his__L(e)','EX_ile__L(e)','EX_leu__L(e)','EX_lys__L(e)',
        'EX_met__L(e)','EX_phe__L(e)','EX_pro__L(e)','EX_ser__L(e)','EX_thr__L(e)',
        'EX_trp__L(e)','EX_tyr__L(e)','EX_val__L(e)')

conc_aa_feeding <- 2.5*c(15.018,8.795,0.000,6.655,0.000,0.000,16.459,7.423,0.000,6.955,
                         11.957,5.917,2.207,5.588,3.082,4.451,6.036,0.000,1.175,9.289)

## Define feeding component concentrations,
## where suppose O2 is always ample during the whole fermentation process
feedSubstrateRxns = c('EX_glc__D(e)','EX_nh4(e)','EX_so4(e)','EX_mg2(e)','EX_fe2(e)',
                      'EX_mn2(e)','EX_cu2(e)','EX_zn2(e)','EX_mobd(e)',AA);
feedConcentrations = c(2775.3,0.00924,40.23174,40,0.0719,0.09171,0.00881,0.05912,0.01078);
feedConcentrations = feedConcentrations/2;
feedConcentrations = c(feedConcentrations,conc_aa_feeding);

## This step runs dGSMEF, and will take a long time
Ec_df2 <- dGSMEF(model,substrateRxns = substrateRxns,initConcentrations = initConcentrations,
				initBiomass=0.0142,u_fix=0.15,x_ind= 31.63,
				feedSubstrateRxns = feedSubstrateRxns,
				feedConcentrations = feedConcentrations,
				yield_rate = 0.6586,
				ngam = 6.34,
				initRatio = 0.70,
				feedRate_ind = 2.8,
				timeStep=0.25,nSteps=180,rcd=FALSE,fld=TRUE);

## compare predictions with experimental results during induction phase
locs <- Ec_df2@biomassVec >= 31.63
ind_time = Ec_df2@timeVec[locs][1]
plot(spline(Ec_df2@timeVec[locs]-ind_time,Ec_df2@biomassVec[locs], n = 201, method = "natural"), 
     col = 2,main='Concentrations',xlab='Time(hrs)',ylab="Biomass(g/L)",ylim = c(31,42),
     type="l",lwd=2);

ind_time = 0
points(c(ind_time,ind_time+4,ind_time+8,ind_time+12,ind_time+16,ind_time+20),
       c(31.63,37.84164381,40.34110603,40.00462097,41.14031209,41.68691719),pch=2,col=2)

par(new=T);	

ind_time = Ec_df2@timeVec[locs][1]	
plot(spline(Ec_df2@timeVec[locs]-ind_time, 
            Ec_df2@concentrationMatrix[Ec_df2@excRxnNames %in% "Protein"][locs]*25163.41/1000, 
            n = 201, method = "natural"),axes=F,col = 4,xlab="",ylab="",type="l",lwd=2);

ind_time = 0
points(c(ind_time,ind_time+4,ind_time+8,ind_time+12,ind_time+16,ind_time+20),
       c(0,2.198356194,5.63889397,8.805379035,10.99968791,12.85308281),pch=3,col=4)
axis(4);
mtext("Recombinant protein(g/L)",side=4,cex=0.8);
legend("bottomright", c("Biomass","Recombinant protein"), col=c(2,4), lty=1, 
       pch=c(2,3),text.font = 1.2, cex = 0.6);

## End(Not run)