Nothing
# sysBiolAlg_roomClass.R
# FBA and friends with R.
#
# Copyright (C) 2010-2014 Gabriel Gelius-Dietrich, Dpt. for Bioinformatics,
# Institute for Informatics, Heinrich-Heine-University, Duesseldorf, Germany.
# All right reserved.
# Email: geliudie@uni-duesseldorf.de
#
# This file is part of sybil.
#
# Sybil is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# Sybil is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with sybil. If not, see <http://www.gnu.org/licenses/>.
# Original code written by Marc Andre Daxer during his bachelor thesis:
# "Analysis of Gene Defects with Mixed Integer Linear Programming" (2011
# at the Heinrich-Heine-University Duesseldorf, Dpt. for Bioinformatics).
# He wrote the package sybilROOM.
#------------------------------------------------------------------------------#
# definition of the class sysBiolAlg_room #
#------------------------------------------------------------------------------#
setClass(Class = "sysBiolAlg_room",
representation(
wu = "numeric",
wl = "numeric",
fnc = "integer",
fnr = "integer"
),
contains = "sysBiolAlg"
)
#------------------------------------------------------------------------------#
# default constructor #
#------------------------------------------------------------------------------#
# contructor for class sysBiolAlg_room
setMethod(f = "initialize",
signature = "sysBiolAlg_room",
definition = function(.Object,
model,
wtflux,
delta = 0.03,
epsilon = 0.001,
LPvariant = FALSE,
absMAX = SYBIL_SETTINGS("MAXIMUM"),
useNames = SYBIL_SETTINGS("USE_NAMES"),
cnames = NULL,
rnames = NULL,
pname = NULL,
scaling = NULL,
writeProbToFileName = NULL, ...) {
if ( ! missing(model) ) {
if (missing(wtflux)) {
tmp <- .generateWT(model, ...)
wtflux <- tmp$fluxes[tmp$fldind]
}
stopifnot(is(model, "modelorg"),
is(wtflux, "numeric"),
is(delta, "numeric"),
is(epsilon, "numeric"),
is(LPvariant, "logical"),
is(absMAX, "numeric"))
stopifnot(length(wtflux) == react_num(model))
# the problem: minimize
#
# | |
# S | 0 | = 0
# | |
# -------------------------
# 1 | |
# 1 | -vwl | >= wl
# 1 | |
# -------------------------
# 1 | |
# 1 | -vwu | <= wu
# 1 | |
# -------------------------
#
# lb | 0
# ub | 1
# ctype: C | B
# obj: 0 | 1
# ---------------------------------------------
# problem dimensions
# ---------------------------------------------
nc <- react_num(model)
nr <- met_num(model)
nCols <- (2 * nc)
nRows <- (nr + 2 * nc)
# ---------------------------------------------
# constraint matrix
# ---------------------------------------------
# ROOM-Boundaries
# if we use the formulation as linear program, delta and
# epsilon are set to zero (see Shlomi et al.)
#if (isTRUE(LPvariant)) {
# wu <- wtflux
# wl <- wtflux
#}
#else {
wu <- wtflux + delta * abs(wtflux) + epsilon
wl <- wtflux - delta * abs(wtflux) - epsilon
#}
vwu <- uppbnd(model) - wu
vwl <- lowbnd(model) - wl
# the initial matrix dimensions
LHS <- Matrix::Matrix(0,
nrow = nRows,
ncol = nCols,
sparse = TRUE)
# rows for the flux variables
LHS[1:nr,1:nc] <- S(model)
# location of the wild type strain
fi <- c(1:nc)
# constraint-matrix for ROOM
diag(LHS[(nr+1) :(nr+nc),1 :nc ]) <- 1
diag(LHS[(nr+nc+1) :nRows ,1 :nc ]) <- 1
diag(LHS[(nr+1) :(nr+nc),(nc+1):nCols]) <- vwl * -1
diag(LHS[(nr+nc+1) :nRows ,(nc+1):nCols]) <- vwu * -1
# ---------------------------------------------
# columns
# ---------------------------------------------
clower <- c(lowbnd(model), rep(0, nc))
cupper <- c(uppbnd(model), rep(1, nc))
if (isTRUE(LPvariant)) {
ctype <- NULL
pt <- "lp"
}
else {
ctype <- c(rep("C", nc), rep("B", nc))
pt <- "mip"
}
# ---------------------------------------------
# rows
# ---------------------------------------------
#rlower <- c(rhs(model), wl, wu)
#rupper <- c(rhs(model), rep(absMAX, nc), wu)
rlower <- c(rep(0, nr), wl, wu)
rupper <- c(rep(0, nr), rep(absMAX, nc), wu)
rtype <- c(rep("E", nr), rep("L", nc), rep("U", nc))
# ---------------------------------------------
# objective function
# ---------------------------------------------
cobj <- c(rep(0, nc), rep(1, nc))
# ---------------------------------------------
# row and column names for the problem object
# ---------------------------------------------
if (isTRUE(useNames)) {
if (is.null(cnames)) {
cn <- c(react_id(model),
paste("oo", react_id(model), sep = "_")
)
colNames <- .makeLPcompatible(cn,
prefix = "x")
}
else {
stopifnot(is(cnames, "character"),
length(cnames) == nCols)
colNames <- cnames
}
if (is.null(rnames)) {
rn <- c(met_id(model),
paste("wl", react_id(model), sep = "_"),
paste("wu", react_id(model), sep = "_")
)
rowNames <- .makeLPcompatible(rn,
prefix = "r")
}
else {
stopifnot(is(rnames, "character"),
length(rnames) == nRows)
rowNames <- rnames
}
if (is.null(pname)) {
probName <- .makeLPcompatible(
paste("ROOM", toupper(pt), mod_id(model), sep = "_"),
prefix = "")
}
else {
stopifnot(is(pname, "character"),
length(pname) == 1)
probName <- pname
}
}
else {
colNames <- NULL
rowNames <- NULL
probName <- NULL
}
# ---------------------------------------------
# build problem object
# ---------------------------------------------
.Object <- callNextMethod(.Object,
sbalg = "room",
pType = pt,
scaling = scaling,
fi = fi,
nCols = nCols,
nRows = nRows,
mat = LHS,
ub = cupper,
lb = clower,
obj = cobj,
rlb = rlower,
rub = rupper,
rtype = rtype,
ctype = ctype,
lpdir = "min",
cnames = colNames,
rnames = rowNames,
pname = probName,
algPar = list("wtflux" = wtflux,
"delta" = delta,
"epsilon" = epsilon),
...)
.Object@wu <- as.numeric(wu)
.Object@wl <- as.numeric(wl)
.Object@fnr <- as.integer(nr)
.Object@fnc <- as.integer(nc)
if (!is.null(writeProbToFileName)) {
writeProb(problem(.Object),
fname = as.character(writeProbToFileName))
}
# # ---------------------------------------------
# # build problem object
# # ---------------------------------------------
#
# lp <- optObj(solver = solver, method = method, pType = "mip")
# lp <- initProb(lp, nrows = nRows, ncols = nCols)
#
# # ---------------------------------------------
# # set control parameters
# # ---------------------------------------------
#
# if (!any(is.na(solverParm))) {
# setSolverParm(lp, solverParm)
# }
#
#
# loadLPprob(lp,
# nCols = nCols,
# nRows = nRows,
# mat = LHS,
# ub = cupper,
# lb = clower,
# obj = cobj,
# rlb = rlower,
# rub = rupper,
# rtype = rtype,
# ctype = ctype,
# lpdir = "min"
# )
#
# if (!is.null(scaling)) {
# scaleProb(lp, scaling)
# }
#
# .Object@problem <- lp
# .Object@algorithm <- "room"
# .Object@nr <- as.integer(nRows)
# .Object@nc <- as.integer(nCols)
# .Object@fldind <- as.integer(fi)
# .Object@wu <- as.numeric(wu)
# .Object@wl <- as.numeric(wl)
# .Object@fnr <- as.integer(nr)
# .Object@fnc <- as.integer(nc)
# validObject(.Object)
}
return(.Object)
}
)
#------------------------------------------------------------------------------#
# other methods #
#------------------------------------------------------------------------------#
setMethod("applyChanges", signature(object = "sysBiolAlg_room"),
function(object, del, obj, ld,
react = NULL,
lb = NULL,
ub = NULL,
obj_coef = NULL,
fldind = TRUE,
lpdir = NULL) {
if (!isTRUE(fldind)) {
warning("argument ", sQuote("fldind"), " is currently not used for ROOM")
}
fi <- fldind(object)[react]
if (any(is.na(fi))) {
stop("argument ", sQuote("react"), " must contain reactions only")
}
tmp_val <- list("fi" = react, "lb" = NULL, "ub" = NULL)
wu <- object@wu
wl <- object@wl
lpmod <- problem(object)
if (isTRUE(del)) {
# store default lower and upper bounds
tmp_val[["lb"]] <- getColsLowBnds(lpmod, fi)
tmp_val[["ub"]] <- getColsUppBnds(lpmod, fi)
# change bounds of fluxes in react
check <- changeColsBnds(lpmod, fi, lb, ub)
# change constraint matrix and objective coefficients
vwu <- (ub - wu[fi]) * -1
vwl <- (lb - wl[fi]) * -1
ri <- react + object@fnr
ci <- react + object@fnc
for (i in seq(along = react)) {
changeMatrixRow(lpmod, ri[i], c(react[i], ci[i]), c(1, vwl[i]))
changeMatrixRow(lpmod,
ri[i]+object@fnc,
c(react[i], ci[i]),
c(1, vwu[i]))
}
#changeObjCoefs(lpmod, ci, rep(0, length(react)))
tmp_val[["ri"]] <- ri
tmp_val[["ci"]] <- ci
}
return(tmp_val)
}
)
#------------------------------------------------------------------------------#
setMethod("resetChanges", signature(object = "sysBiolAlg_room"),
function(object, old_val) {
fi <- fldind(object)
wu <- object@wu
wl <- object@wl
lpmod <- problem(object)
if ( (!is.null(old_val[["lb"]])) || (!is.null(old_val[["ub"]])) ) {
check <- changeColsBnds(lpmod,
fi[old_val[["fi"]]],
old_val[["lb"]], old_val[["ub"]])
vwu <- (old_val[["ub"]] - wu[fi[old_val[["react"]]]]) * -1
vwl <- (old_val[["lb"]] - wl[fi[old_val[["react"]]]]) * -1
for (i in seq(along = old_val[["fi"]])) {
changeMatrixRow(lpmod,
old_val[["ri"]][i],
c(old_val[["fi"]][i], old_val[["ci"]][i]),
c(1, vwl[i]))
changeMatrixRow(lpmod,
old_val[["ri"]][i]+object@fnc,
c(old_val[["fi"]][i], old_val[["ci"]][i]),
c(1, vwu[i]))
}
#changeObjCoefs(lpmod, old_val[["ci"]], rep(1, length(react)))
}
return(invisible(TRUE))
}
)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.