cran 1.8.1 release

DESCRIPTION

@@ -1,13 +1,14 @@
 Package: BacArena
 Title: Modeling Framework for Cellular Communities in their Environments
-Version: 1.8
+Version: 1.8.1
 Authors@R: c( person("Eugen", "Bauer", , "eugen.bauer@uni.lu", role = c("aut")), 
               person("Johannes", "Zimmermann", , "j.zimmermann@iem.uni-kiel.de", role = c("aut", "cre")))
 Author: Eugen Bauer [aut],
     Johannes Zimmermann [aut, cre]
 Maintainer: Johannes Zimmermann <j.zimmermann@iem.uni-kiel.de>
 Description: Can be used for simulation of organisms living in
-    communities. Each organism is represented individually and genome scale
+    communities (Bauer and Zimmermann (2017) <doi:10.1371/journal.pcbi.1005544>). 
+    Each organism is represented individually and genome scale
     metabolic models determine the uptake and release of compounds. Biological
     processes such as movement, diffusion, chemotaxis and kinetics are available
     along with data analysis techniques.
@@ -28,9 +29,11 @@ Imports:
     reshape2,
     glpkAPI,
     plyr,
-    Rcpp
-Suggests:
+    Rcpp,
     igraph,
+    stringr,
+    R.matlab
+Suggests:
     parallel,
     knitr,
     rmarkdown
@@ -39,5 +42,5 @@ LinkingTo: Rcpp,
     RcppEigen
 License: GPL-3 | file LICENSE
 VignetteBuilder: knitr
-RoxygenNote: 6.1.0
+RoxygenNote: 7.0.2
 LazyData: true

NEWS

@@ -1,4 +1,17 @@
-BacArena v.1.8 (Release data: 2019-02-15)
+BacArena v.1.8.1 (Release date: 2019-12-16)
+================
+Changes:
+- matlab import support
+- fixed several issues when adding organisms
+- corrected default diffusion constant in addSubs
+- igraph is now only suggested
+- complied with new R CRAN guidelines (matrix is of class matrix + array now)
+
+* github reference: https://github.com/euba/BacArena/releases/tag/v1.8.1
+
+
+BacArena v.1.8 (Release date: 2019-02-15)
+==============
 Changes:
 - improved detection of limiting substances (shadow costs) 
 - shadow costs analysis now possible with cplex, too

R/Arena.R

@@ -736,7 +736,6 @@ setMethod("checkPhen", "Arena", function(object, org, cutoff=1e-6, fbasol){
 #' @param cutoff A number giving the cutoff for values of the objective function and fluxes of exchange reactions.
 #' @param fbasol Problem object according to the constraints and then solved with \code{optimizeProb}.
 #' 
-#' @rdname checkPhen_par
 setGeneric("checkPhen_par", function(object, org, cutoff=1e-6, fbasol){standardGeneric("checkPhen_par")})
 #' @export
 #' @rdname checkPhen_par
@@ -2723,7 +2722,7 @@ setMethod("findFeeding3", "Eval", function(object, time, mets, plot=TRUE, cutoff
         idx.flux <- match(interact[k,], colnames(mfluxmat))
         mfluxmat[i, idx.flux]
       })
-      if(class(interact)=="character"){interact = t(as.matrix(interact))}
+      if("character" %in% class(interact)){interact = t(as.matrix(interact))}
       inter = rbind(inter,data.frame(prod=interact[,1],cons=interact[,2],met=rownames(mfluxmat)[i], sim_step=time-1, prod.flux=flux[1,], cons.flux=flux[2,]))
     }
   }

R/Organism.R

@@ -124,6 +124,8 @@ findrBiomass <- function(model, keys=c("biom", "cpd11416")){
 #' @param lyse A boolean variable indicating if the organism should lyse after death.
 #' @param setExInf Enable if all lower bounds of exchange reaction which are set to zero (i.e. no uptake possible!) should be set to -infitity (default: true)
 #' @param setAllExInf Enable if all lower bounds of exchange reaction should be set to -infitity (default: false)
+#' @param coupling_constraints List with coupling parameters.
+#' @param predator Name of organism which can kill this one.
 #' @param ... Arguments of \code{\link{Organism-class}}
 #' @return Object of class Organism
 Organism <- function(model, algo="fba", ex="EX_", ex_comp=NA, csuffix="\\[c\\]", esuffix="\\[e\\]", lyse=FALSE, feat=list(), 
@@ -185,7 +187,7 @@ Organism <- function(model, algo="fba", ex="EX_", ex_comp=NA, csuffix="\\[c\\]",
     lobnd[ which(names(lobnd) %in% medc & lobnd==0) ] <- -1000
   }
   lobnd.ex <- lobnd[match(medc, rxname)]
-  lobnd.ex.med <- median(lobnd.ex[ which( lobnd.ex < 0 & lobnd.ex > -1000 ) ])
+  lobnd.ex.med <- stats::median(lobnd.ex[ which( lobnd.ex < 0 & lobnd.ex > -1000 ) ])
   if( !is.na(lobnd.ex.med) ){
     print(paste0("Median lower bound for non-zero and non-Inf exchanges is:", round(lobnd.ex.med), 6))
     if( lobnd.ex.med > -10 ){
@@ -198,8 +200,8 @@ Organism <- function(model, algo="fba", ex="EX_", ex_comp=NA, csuffix="\\[c\\]",
                                           ifelse( length(sybil::mod_name(model)) > 0, sybil::mod_name(model), 
                                                   ifelse( length(sybil::mod_id(model)) > 0, sybil::mod_id(model), "test" )))
   if(algo=="coupling"){
-    #lpobject <- sybil::sysBiolAlg(model, algorithm="mtfEasyConstraint2", easyConstraint=coupling, pFBAcoeff = 1e-5)
-    lpobject <- sybil::sysBiolAlg(model, algorithm="mtfEasyConstraint", easyConstraint=coupling)
+    #lpobject <- sybil::sysBiolAlg(model, algorithm="mtfEasyConstraint2", easyConstraint=coupling_constraints, pFBAcoeff = 1e-5)
+    lpobject <- sybil::sysBiolAlg(model, algorithm="mtfEasyConstraint", easyConstraint=coupling_constraints)
   }else lpobject <- sybil::sysBiolAlg(model, algorithm=algo)
   fbasol <- sybil::optimizeProb(lpobject, react=1:length(lobnd), ub=upbnd, lb=lobnd)
   names(fbasol$fluxes) = rxname

R/Stuff.R

@@ -927,7 +927,7 @@ findFeeding3rep <- function(simlist, time, mets, plot=TRUE, mfunction="mean"){
       if(length(which(x>0))!=0){interact = rbind(interact,cbind(names(which(x>0)),j))}
     }
     interact = interact[-1,]
-    if(class(interact)=="character"){interact = t(as.matrix(interact))}
+    if("character" %in% class(interact)){interact = t(as.matrix(interact))}
     if(nrow(interact)!=0){inter = rbind(inter,data.frame(prod=interact[,1],cons=interact[,2],met=i))}
   }
   if(any(dim(inter)==0)) {

R/data.R

@@ -24,7 +24,7 @@
 #' 
 #' @format A vector with 22 color hex-codes 
 #' 
-#' @source \url{www.iscc.org/pdf/PC54_1724_001.pdf}
+#' @source \url{http://www.iscc-archive.org/pdf/PC54_1724_001.pdf}
 "colpal3"
 
 

R/matlab.R

@@ -7,11 +7,8 @@
 #' @param file Full path to matlab model file
 #' @details Returns sybil model object (time needed: bacterial model ~ 10s, recon2 ~ 60s)
 readMATmod <- function(file){
-  library(sybil)
-  library(R.matlab)
-  library(stringr)
-
-  print(system.time(data <- readMat(file)))
+
+  print(system.time(data <- R.matlab::readMat(file)))
   dat.mat <- data[[1]]
 
   mod.var <- dimnames(dat.mat)[[1]]
@@ -70,7 +67,7 @@ readMATmod <- function(file){
         mod.genes[[i]] <- ""
         next
       } 
-      j <- as.numeric(unlist(str_extract_all(rule.tmp, "(?<=x\\()[0-9]+?(?=\\))")))
+      j <- as.numeric(unlist(stringr::str_extract_all(rule.tmp, "(?<=x\\()[0-9]+?(?=\\))")))
       mod.GeneMat[i,j] <- 1
       mod.genes[[i]] <- mod.gene_id[j]
     }
@@ -80,9 +77,9 @@ readMATmod <- function(file){
     mod.gprRules <- sapply(sapply(dat.mat[[which(mod.var=="rules")]], unlist), function(entry){
       if( length(entry) == 0 ) "" 
       else {
-        numbers <- as.numeric(unlist(str_extract_all(entry, "[0-9]+")))
+        numbers <- as.numeric(unlist(stringr::str_extract_all(entry, "[0-9]+")))
         dict <- as.character(numbers - min(numbers) + 1); names(dict) <- as.character(numbers)
-        gsub("\\(([0-9]+)\\)","\\[\\1\\]",str_replace_all(entry, dict)) }
+        gsub("\\(([0-9]+)\\)","\\[\\1\\]",stringr::str_replace_all(entry, dict)) }
     })
   }else{ # if 'rules' is not present construct own rules from gpr+genes
     mod.gprRules <- sapply(seq_along(mod.gpr), function(i){
@@ -91,7 +88,7 @@ readMATmod <- function(file){
         genes <- unlist(mod.genes[i])
         dict <- paste0("x[",seq_along(genes),"]")
         names(dict) <- genes
-        gsub("and","&", gsub("or","|",str_replace_all(mod.gpr[i], dict)))
+        gsub("and","&", gsub("or","|",stringr::str_replace_all(mod.gpr[i], dict)))
       }
     })
 
@@ -102,9 +99,9 @@ readMATmod <- function(file){
   mod.ub <- as.vector(dat.mat[[which(mod.var=="ub")]])
 
   # 7) compartments
-  met_comp <- str_extract_all(mod.met_id, "(?<=\\[)[a-z](?=\\])")
+  met_comp <- stringr::str_extract_all(mod.met_id, "(?<=\\[)[a-z](?=\\])")
   if( all(sapply(met_comp, length) == 0 ) ){
-    met_comp <- str_extract_all(mod.met_id, "(?<=_)[a-z][0-9]?(?=$)")
+    met_comp <- stringr::str_extract_all(mod.met_id, "(?<=_)[a-z][0-9]?(?=$)")
   }
   mod.mod_compart <- unique(unlist(met_comp))
   mod.met_comp    <- match(met_comp, mod.mod_compart)
@@ -121,7 +118,7 @@ readMATmod <- function(file){
 
 
   # create new model
-  model <- modelorg(id = mod.id, name = mod.name)
+  model <- sybil::modelorg(id = mod.id, name = mod.name)
   model@mod_desc <- mod.desc
   model@S <- mod.S
   model@lowbnd <- mod.lb
@@ -143,8 +140,8 @@ readMATmod <- function(file){
 
   obj.idx <- which(dat.mat[[which(mod.var == "c")]]!=0)
   if( length(obj.idx) > 0 ){
-    model <- changeObjFunc(model, react = obj.idx)  
-    print(optimizeProb(model))
+    model <- sybil::changeObjFunc(model, react = obj.idx)  
+    print(sybil::optimizeProb(model))
   }
 
   return(model)

examples/dev.R

@@ -13,7 +13,7 @@ devtools::document()
 
 devtools::check() # includes document()
 
-install_local(path="~/uni/bacarena", quick=TRUE, threads=2, quiet=F, force=T)
+install_local(path="~/uni/bacarena", quick=TRUE, threads=2, quiet=F, force=T, upgrade=F)
 #install_github(repo="euba/bacarena", quick=TRUE)
 
 

examples/ecoli.R

@@ -10,11 +10,26 @@ arena <- addOrg(arena,bac,amount=20)
 arena <- addSubs(arena, smax=0.005, mediac="EX_glc(e)", unit="mM")
 arena <- addSubs(arena, smax=1, mediac=c("EX_pi(e)", "EX_h2o(e)",
                                          "EX_o2(e)", "EX_nh4(e)"), unit="mM") 
+arena@removeM[c(1,20),1:20] <- 1
+arena@removeM[1:20,c(1,20)] <- 1
+image(arena@removeM)
 sim   <- simEnv(arena, time=5, with_shadow = T, sec_obj = F)
 
+
+
 getSubHist(sim, "EX_glc(e)")
 
 sim@shadowlist[[3]]
 plotShadowCost(sim, 1)[[2]] # fru, glc -2
 
 plotGrowthCurve(sim.old)
+
+sim <- simEnv(arena, time=2)
+arena2 <- getArena(sim,2)
+arena2@orgdat <- arena2@orgdat[-sample(nrow(arena2@orgdat), round(nrow(arena2@orgdat) / 10)), ] #remove ~10% of individuals randomly 
+simEnv(arena2, time=2)
+
+arena2 <- getArena(sim,2)
+neworgdat <- arena2@orgdat #get the current orgdat
+neworgdat <- neworgdat[-sample(nrow(neworgdat), round(nrow(neworgdat) / 10)), ] #remove ~10% of individuals randomly 
+arena2 <- changeOrg(arena2,neworgdat)