eFBA_gene: Function: eFBA

Description Usage Arguments Details Value Author(s) See Also Examples

View source: R/eFBA_gene.R

Description

This function performs an expression based flux balance analysis. It takes into account the expression status of genes. First step calculates FB as the max biomass using standard FBA then uses FB=cTx as a new constraint in a MILP with a new objective function Mininmize: Sum(expr(g) != flux(g) where expr(g)=1 if g is expressed , else 0 flux(g)=0 if (all reactions catalyzed by g) have no flux(Threshold Tf and sum of all), 0 otherwise.

Usage

1
2
3
4
5
6
7
8
9
eFBA_gene(model, expressionData,
                  Tf =0.01,#SYBIL_SETTINGS("TOLERANCE"),  
				  pct_objective=100,
 		          lpdir = SYBIL_SETTINGS("OPT_DIRECTION"),
				  solver = SYBIL_SETTINGS("SOLVER"),
			      method = SYBIL_SETTINGS("METHOD"),
				  solverParm=data.frame(CPX_PARAM_EPRHS=1e-6),
				  testgpr=NULL,
              verboseMode = 2, ...)

Arguments

model

An object of class modelorg.

expressionData

a dataframe (gene expression data): gene ID, and expression status (0 :OFF or 1: ON)

pct_objective

Biomass will be garnteed to be at least this value multiplied by max biomass calculated using standard FBA. Values are 0 to 100.

Tf

Threshold in flux default value 1e-5. Fluxes with magnitude less than Tf are considered zeros.

lpdir

Character value, direction of optimisation. Can be set to "min" or "max".
Default: SYBIL_SETTINGS("OPT_DIRECTION").

solver

Single character value. The solver to use. See SYBIL_SETTINGS for possible values.
Default: SYBIL_SETTINGS("SOLVER").

method

Single character value. The optimization algorithm to use. Possible values depend on the setting in solver. See SYBIL_SETTINGS for possible values.
Default: LP_METHOD(SYBIL_SETTINGS).

solverParm

A data frame containing parameters for the specified solver.
Default: SOLVER_CTRL_PARM(SYBIL_SETTINGS).

testgpr

Boolean flag containing flags of whether a reaction is considered or not

verboseMode

An integer value indicating the amount of output to stdout: 0: nothing, 1: status messages, 2: like 1 plus with more details, 3: generates files of the LP problem.
Default: 2.

...

Further arguments passed to optimizeProb. Argument solverParm is a good candidate.

Details

formulates the following LP: the problem: minimize:

| | S | 0 | = b | | —————— | | c^T | 0 | = FB | | —————— | | -1 | -ub | <= 0 ] | | > -ub*yi <= vi <= ub * yi +1 | -ub | <= 0 ] | | —————— lb wt_lb| 0 | ub wt_ub| 1 | | | obj 0 | 2*ei-1 where ei =0 gene i NOT expressed/ 1 otherwise

Value

returns a list containing slots: origfluxes: original fluxes calculted by standard FBA. fluxes: new fluxes that better matches gene expression. gene_stat: state of genes (0 or 1) according to calculated fluxes stat: solution status returned from solver

Author(s)

Abdelmoneim Amer Desouki

See Also

modelorg, optimizeProb, eFBA_rxn

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
## Not run: 
	## The examples here require the package glpkAPI to be
	## installed. If that package is not available, you have to set
	## the argument 'solver' (the default is: solver = "glpkAPI").

	## load the example data set
	data(iAF1260)
        model=iAF1260
	allgenes=allGenes(model)
	##set all genes to OFF
	ex=cbind(LOCUS=allgenes,State=rep(0,length(allgenes)));
	##set state of genes as required
	
	## Optionally optimize only for a set of reactions
	testgpr=rep(FALSE,react_num(model))
	testgpr[react_id(model) %in% c('R_PYK','R_FLDR','R_FEROpp')]=TRUE
	 
	 # read or set the expression status of genes : 0 not expressed, 
	 #    1 Expressed, NA: unknown or ignored
	 ex_on=ex
	 ex_on[ex_on[,1] %in% c('b0123','b1854','b0684','b3924'),2]=1
	 
	 ex_off=ex# all genes are not expressed
	
	 slvr="glpkAPI"
	# run two times, one with genes ON, the other with genes OFF
	 ser_on = eFBA_gene(model, ex_on,solver=slvr,verbose=3,testgpr=testgpr,Tf=0.001);
	 ser_off = eFBA_gene(model, ex_off,solver=slvr,verbose=3,testgpr=testgpr,Tf=0.001)
	
	 #print results
	 cbind(rxn=c('R_PYK','R_FLDR','R_FEROpp'),geneONfluxes=ser_on$rxn[testgpr, "newFlux"],
	     geneOFFfluxes=ser_off$rxn[testgpr, "newFlux"],
	      ruleStateON=ser_on$rxn[testgpr,"expr"],ruleStateOFF=ser_off$rxn[testgpr,"expr"],
	       FBA_flux=ser_off$rxn[testgpr,"origFlux"])
	
	#Count difference between these two solutions
	sum(abs( as.numeric(ser_on$rxn[, "newFlux"])-as.numeric(ser_off$rxn[, "newFlux"]))>0.001)



## End(Not run)

sybilEFBA documentation built on May 29, 2017, 9:35 a.m.