# sysBiolAlg_mtfClass.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/>. #------------------------------------------------------------------------------# # definition of the class sysBiolAlg_mtf # #------------------------------------------------------------------------------# setClass(Class = "sysBiolAlg_mtf", representation( maxobj = "numeric" ), contains = "sysBiolAlg" ) #------------------------------------------------------------------------------# # default constructor # #------------------------------------------------------------------------------# # contructor for class sysBiolAlg_mtf setMethod(f = "initialize", signature = "sysBiolAlg_mtf", definition = function(.Object, model, wtobj = NULL, react = NULL, lb = NULL, ub = NULL, costcoefbw = NULL, costcoeffw = NULL, absMAX = SYBIL_SETTINGS("MAXIMUM"), useNames = SYBIL_SETTINGS("USE_NAMES"), cnames = NULL, rnames = NULL, pname = NULL, scaling = NULL, writeProbToFileName = NULL, ...) { if ( ! missing(model) ) { if (is.null(wtobj)) { tmp <- .generateWT(model, react, lb, ub, ...) wtobj <- tmp[["obj"]] } stopifnot(is(model, "modelorg"), is(wtobj, "numeric"), is(absMAX, "numeric")) # If wtobj is longer than 1, mtf algorithm has to run several # times. In that case, wtobj is not written in the problem # object, it is written separately (maxobj) and used for # each iteration. if (length(wtobj) > 1) { maxobj <- wtobj currmo <- 0 } else { maxobj <- NULL currmo <- wtobj[1] } # the problem: minimize: # # | | | # S | 0 | 0 | = b # | | | # ------------------------- # | | | # 1 | 1 | 0 | >= 0 # | | | # ------------------------- # | | | # -1 | 0 | 1 | >= 0 # | | | # ------------------------- # c_wt | 0 | 0 | >= c^T * v_wt # | | | # lb wt_lb| 0 | 0 | # ub wt_ub|10000 |10000 | # | | | # obj 0 | 1 | 1 | # --------------------------------------------- # problem dimensions # --------------------------------------------- nc <- react_num(model) nr <- met_num(model) nCols <- 3*nc nRows <- nr + 2*nc + 1 # --------------------------------------------- # constraint matrix # --------------------------------------------- # the initial matrix dimensions LHS <- Matrix::Matrix(0, nrow = nRows, ncol = nCols, sparse = TRUE) # rows for the mutant strain LHS[1:nr,1:nc] <- S(model) # location of the mutant strain fi <- c(1:nc) # rows for the delta match matrix diag(LHS[(nr+1) :(nr+nc) ,1 :nc ]) <- 1 diag(LHS[(nr+1) :(nr+nc) ,(nc+1) :(2*nc)]) <- 1 diag(LHS[(nr+nc+1):(nr+2*nc),1 :nc ]) <- -1 diag(LHS[(nr+nc+1):(nr+2*nc),(2*nc+1):(3*nc)]) <- 1 # fix the value of the objective function LHS[(nr+2*nc+1),1:nc] <- obj_coef(model) # --------------------------------------------- # columns # --------------------------------------------- lower <- c(lowbnd(model), rep(0, 2*nc)) upper <- c(uppbnd(model), rep(absMAX, 2*nc)) # --------------------------------------------- # rows # --------------------------------------------- #rlower <- c(rhs(model), rep(0, 2*nc), currmo) #rupper <- c(rhs(model), rep(absMAX, 2*nc + 1)) rlower <- c(rep(0, nr), rep(0, 2*nc), currmo) rupper <- c(rep(0, nr), rep(absMAX, 2*nc + 1)) rtype <- c(rep("E", nr), rep("L", 2*nc + 1)) # --------------------------------------------- # objective function # --------------------------------------------- if (is.null(costcoeffw)) { fw <- rep(1, nc) } else { stopifnot(is(costcoeffw, "numeric"), (length(costcoeffw) == nc)) fw <- costcoeffw } if (is.null(costcoefbw)) { bw <- fw } else { stopifnot(is(costcoefbw, "numeric"), (length(costcoefbw) == nc)) bw <- costcoefbw } cobj <- c(rep(0, nc), bw, fw) # --------------------------------------------- # row and column names for the problem object # --------------------------------------------- if (isTRUE(useNames)) { if (is.null(cnames)) { cn <- c(react_id(model), paste("bw", react_id(model), sep = "_"), paste("fw", 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("bw", 1:nc, sep = "_"), paste("fw", 1:nc, sep = "_"), "obj_wt" ) rowNames <- .makeLPcompatible(rn, prefix = "r") } else { stopifnot(is(rnames, "character"), length(rnames) == nRows) rowNames <- rnames } if (is.null(pname)) { probName <- .makeLPcompatible( paste("MTF", 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 = "mtf", pType = "lp", scaling = scaling, fi = fi, nCols = nCols, nRows = nRows, mat = LHS, ub = upper, lb = lower, obj = cobj, rlb = rlower, rub = rupper, rtype = rtype, lpdir = "min", ctype = NULL, cnames = colNames, rnames = rowNames, pname = probName, algPar = list("wtobj" = wtobj, "costcoefbw" = bw, "costcoeffw" = fw), ...) .Object@maxobj <- as.numeric(maxobj) if (!is.null(writeProbToFileName)) { writeProb(problem(.Object), fname = as.character(writeProbToFileName)) } # # # --------------------------------------------- # # build problem object # # --------------------------------------------- # # lp <- optObj(solver = solver, method = method) # 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 = upper, # lb = lower, # obj = cobj, # rlb = rlower, # rub = rupper, # rtype = rtype, # lpdir = "min" # ) # # if (!is.null(scaling)) { # scaleProb(lp, scaling) # } # # .Object@problem <- lp # .Object@algorithm <- "mtf" # .Object@nr <- as.integer(nRows) # .Object@nc <- as.integer(nCols) # .Object@fldind <- as.integer(fi) # validObject(.Object) } return(.Object) } ) #------------------------------------------------------------------------------# # other methods # #------------------------------------------------------------------------------# setMethod("changeMaxObj", signature(object = "sysBiolAlg_mtf"), function(object, j) { if (!is.null(object@maxobj)) { changeRowsBnds(problem(object), i = nr(object), lb = object@maxobj[j], ub = SYBIL_SETTINGS("MAXIMUM")) } return(invisible(TRUE)) } )