Skip to content

Commit

Permalink
updates 1.6 cran release
Browse files Browse the repository at this point in the history
  • Loading branch information
jotech committed May 23, 2017
1 parent dc53385 commit ddbf7e2
Show file tree
Hide file tree
Showing 22 changed files with 110 additions and 96 deletions.
17 changes: 0 additions & 17 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -51,10 +51,8 @@ export(growExp)
export(growLin)
export(growth)
export(growth_par)
export(lpobj)
export(lsd)
export(lysis)
export(medium)
export(minePheno)
export(move)
export(openArena)
Expand Down Expand Up @@ -89,7 +87,6 @@ export(simEnv_par)
export(simHum)
export(statPheno)
export(stirEnv)
export(type)
export(unit_conversion)
export(usd)
exportClasses(Arena)
Expand Down Expand Up @@ -117,16 +114,13 @@ exportMethods(constrain)
exportMethods(consume)
exportMethods(createGradient)
exportMethods(dat2mat)
exportMethods(deathrate)
exportMethods(diffuse)
exportMethods(diffusePDE)
exportMethods(diffuseR)
exportMethods(diffuse_par)
exportMethods(emptyHood)
exportMethods(evalArena)
exportMethods(extractMed)
exportMethods(fbasol)
exportMethods(feat)
exportMethods(findFeeding)
exportMethods(findFeeding2)
exportMethods(findFeeding3)
Expand All @@ -145,16 +139,8 @@ exportMethods(growExp)
exportMethods(growLin)
exportMethods(growth)
exportMethods(growth_par)
exportMethods(growtype)
exportMethods(kinetics)
exportMethods(lbnd)
exportMethods(lpobj)
exportMethods(lyse)
exportMethods(lysis)
exportMethods(medium)
exportMethods(minePheno)
exportMethods(minweight)
exportMethods(model)
exportMethods(move)
exportMethods(optimizeLP)
exportMethods(plotCurves)
Expand All @@ -172,11 +158,8 @@ exportMethods(simBac_par)
exportMethods(simEnv)
exportMethods(simEnv_par)
exportMethods(simHum)
exportMethods(speed)
exportMethods(statPheno)
exportMethods(stirEnv)
exportMethods(type)
exportMethods(ubnd)
exportMethods(unit_conversion)
import(Matrix)
import(Rcpp)
Expand Down
34 changes: 15 additions & 19 deletions R/Arena.R
Original file line number Diff line number Diff line change
Expand Up @@ -79,6 +79,7 @@ setClass("Arena",
#' @param m A number giving the vertical size of the environment.
#' @param Lx A number giving the horizontal grid size in cm.
#' @param Ly A number giving the vertical grid size in cm.
#' @param seed An integer refering to the random number seed used to be reproducible
#' @param ... Arguments of \code{\link{Arena-class}}
Arena <- function(Lx=NULL, Ly=NULL, n=100, m=100, seed=sample(1:10000,1), ...){
if(is.null(Lx)) Lx <- 0.025/100 * n
Expand Down Expand Up @@ -719,6 +720,7 @@ setMethod("addPhen", "Arena", function(object, org, pvec){
#' @export
#' @rdname unit_conversion
#'
#' @param object An object of class Arena or Eval.
#' @param unit Unit to be converted to fmol/cell
#' @return Conversion factor
setGeneric("unit_conversion", function(object, unit){standardGeneric("unit_conversion")})
Expand Down Expand Up @@ -1776,6 +1778,7 @@ setMethod("evalArena", "Eval", function(object, plot_items='Population', phencol
#' @param retdata A boolean variable indicating if the data used to generate the plots should be returned.
#' @param remove A boolean variable indicating if substances, which don't change in their concentration should be removed from the plot.
#' @param legend Boolean variable indicating if legend should be plotted
#' @param graph True if graphic should be plotted.
#' @return Returns two graphs in one plot: the growth curves and the curves of concentration changes. Optional the data to generate the original plots can be returned.
#' @details The parameter \code{retdata} can be used to access the data used for the returned plots to create own custom plots.
#' @seealso \code{\link{Eval-class}} and \code{\link{Arena-class}}
Expand Down Expand Up @@ -1991,7 +1994,7 @@ setMethod("plotCurves2", "Eval", function(object, legendpos="topright", ignore=c
rownames(mat_nice) <- gsub("\\[e\\]","", gsub("\\(e\\)","", gsub("EX_","",object@mediac[subs_index])))
}
if(num>length(colpal3)) cols <- colpal1[1:num] else cols <- colpal3[1:num]
matplot(t(mat_nice), type='l', col=cols, pch=1, lty=1, lwd=5,
graphics::matplot(t(mat_nice), type='l', col=cols, pch=1, lty=1, lwd=5,
xlab=paste0('time in ', ifelse(object@tstep==1, "", object@tstep), 'h'), ylab='amount of substance in fmol',
main='Strongly changing substances')
if(length(dict) > 0){
Expand Down Expand Up @@ -2037,7 +2040,7 @@ setMethod("plotCurves2", "Eval", function(object, legendpos="topright", ignore=c

len <- dim(mat_with_phen)[1]
if(len>length(colpal3)) cols <- colpal1[1:len] else cols <- colpal3[1:len]
matplot(t(mat_with_phen), type='b', col=cols, pch=1, lty=1, lwd=5,
graphics::matplot(t(mat_with_phen), type='b', col=cols, pch=1, lty=1, lwd=5,
xlab=paste0('time in ', ifelse(object@tstep==1, "", object@tstep), 'h'), ylab='amount of organisms',
main='Growth curve')
legend(legendpos, rownames(mat_with_phen), col=cols, cex=0.7, fill=cols)
Expand Down Expand Up @@ -2076,7 +2079,7 @@ setMethod("plotTotFlux", "Eval", function(object, legendpos="topright", num=20){
mat_nice <- tail(mat[order(mat_var),], num)

if(num>length(colpal3)) cols <- colpal1[1:num] else cols <- colpal3[1:num]
matplot(t(mat_nice), type='l', col=cols, pch=1, lty=1, lwd=3,
graphics::matplot(t(mat_nice), type='l', col=cols, pch=1, lty=1, lwd=3,
xlab='time in h', ylab='reaction activity in mmol/(h * g_DW)',
main='Highly active reactions')
legend(legendpos, rownames(mat_nice), col=cols, cex=0.6, fill=cols)
Expand Down Expand Up @@ -2346,7 +2349,7 @@ setMethod("statPheno", "Eval", function(object, type_nr=1, phenotype_nr, dict=NU
mat_nice <- mat_sub[names(high_corr)[names(high_corr)!="phenotype"],][,t_lb:t_ub]
num <- dim(mat_nice)[1]
if(num>length(colpal3)) cols <- colpal1[1:num] else cols <- colpal3[1:num]
matplot(x=seq(t_lb, t_ub), y=t(mat_nice), type='l', col=cols, pch=1, lty=1, lwd=5,
graphics::matplot(x=seq(t_lb, t_ub), y=t(mat_nice), type='l', col=cols, pch=1, lty=1, lwd=5,
xlab=paste0('time in ', ifelse(object@tstep==1, "", object@tstep), 'h'), ylab='amount of substance in mmol',
main='correlated substances')
legend("topleft", rownames(mat_nice), col=cols, cex=0.4, fill=cols)
Expand Down Expand Up @@ -2659,7 +2662,7 @@ setMethod("statSpec", "Eval", function(object, type_nr=1, dict=NULL,

# plot growth curve with all phenotypes
mat_phen <- do.call(cbind, occ)
matplot(t(mat_phen), type='b', col=cols, pch=1, lty=1, lwd=5,
graphics::matplot(t(mat_phen), type='b', col=cols, pch=1, lty=1, lwd=5,
xlab=paste0('time in ', ifelse(object@tstep==1, "", object@tstep), 'h'), ylab='amount of organisms',
main='Growth curve')
legend(legend_pos, rownames(mat_phen), col=cols, cex=legend_cex, fill=cols)
Expand Down Expand Up @@ -2825,15 +2828,6 @@ setMethod("plotShadowCost", "Eval", function(object, spec_nr=1, sub_nr=10, cutof
#' @param rnd An integer giving the decimal place to which min/max flux should be rounded.
#' @details Returns a list with the minimum and maximum substance usage for each time point.
#' @seealso \code{\link{Eval-class}} and \code{\link{simEnv}}
#' @examples
#' data(Ec_core, envir = environment()) #get Escherichia coli core metabolic model
#' bac <- Bac(Ec_core,deathrate=0.05,
#' minweight=0.05,growtype="exponential") #initialize a bacterium
#' arena <- Arena(n=20,m=20) #initialize the environment
#' arena <- addOrg(arena,bac,amount=10) #add 10 organisms
#' arena <- addSubs(arena,40) #add all possible substances
#' eval <- simEnv(arena,5)
#' fluxlist <- fluxVarSim(eval, 6)
setGeneric("fluxVarSim", function(object, rnd){standardGeneric("fluxVarSim")})
#' @export
#' @rdname fluxVarSim
Expand All @@ -2854,11 +2848,11 @@ setMethod("fluxVarSim", "Eval", function(object, rnd){
mconc = unlist(lapply(arenait@media,function(med,x,y){med@diffmat[x,y]},x=org$x,y=org$y))
mconc = mconc[bact@medium]
bacnum = round((arenait@scale/(bact@cellarea*10^(-8))))
lbs = constrain(bact,names(mconc),-mconc/bacnum,org$biomass,arenait@tstep,arenait@scale,j)
lbs = constrain(bact,names(mconc),-mconc/bacnum,org$biomass,arenait@tstep,arenait@scale,j)[[1]]
fbasl <- optimizeProb(bact@lpobj, react=1:length(lbs), ub=bact@ubnd, lb=lbs)
model = changeBounds(bact@model,names(lbs),lb=lbs)
model = changeBounds(model,model@react_id[which(bact@model@obj_coef==1)],lb=fbasl$obj,ub=fbasl$obj)
nil=capture.output(suppressMessages(fv <- fluxVar(model, bact@medium)))
nil=utils::capture.output(suppressMessages(fv <- fluxVar(model, bact@medium)))

flmat_max[bact@medium,j] = round(minSol(fv,lp_obj),rnd)
flmat_min[bact@medium,j] = round(maxSol(fv,lp_obj),rnd)
Expand Down Expand Up @@ -2950,6 +2944,7 @@ setMethod("findRxnFlux", "Eval", function(object, ex, time, print_reactions=FALS
#'
#' @param object An object of class Eval.
#' @param sub Name of a substance.
#' @param times Time points to be considered.
#' @details Returns a plot with

setGeneric("plotSubDist", function(object, sub, times=NULL){standardGeneric("plotSubDist")})
Expand All @@ -2976,15 +2971,16 @@ setMethod("plotSubDist", "Eval", function(object, sub, times=NULL){
#'
#' @param object An object of class Eval.
#' @param sub Name of a substance.
#' @param times Time points to be considered.
#' @details Returns a plot with

setGeneric("plotSubDist2", function(object, sub, times=NULL){standardGeneric("plotSubDist2")})
#' @export
#' @rdname plotSubDist2
setMethod("plotSubDist2", "Eval", function(object, sub, times=NULL){
if(length(sub) != 1 | !all(sub %in% object@mediac)) stop("Please use exactly one substance.")
if(max(times) > max(seq_along(object@medlist))) stop("Please use another maximum value in the «times» argument. Your input was out of bounds.")
if(min(times) < min(seq_along(object@medlist))) stop("Please use another minimum value in the «times» argument. Your input was out of bounds.")
if(max(times) > max(seq_along(object@medlist))) stop("Please use another maximum value in 'times' argument. Your input was out of bounds.")
if(min(times) < min(seq_along(object@medlist))) stop("Please use another minimum value in 'times argument. Your input was out of bounds.")
if(length(times)==0) times <- seq_along(object@medlist)
all_df <- data.frame()
for(t in times){
Expand All @@ -2993,7 +2989,7 @@ setMethod("plotSubDist2", "Eval", function(object, sub, times=NULL){
df$time = t
all_df <- rbind(all_df, df)
}
q <- ggplot2::ggplot(all_df, ggplot2::aes(x, y)) + ggplot2::geom_tile(ggplot2::aes(fill = value)) +
q <- ggplot2::ggplot(all_df, ggplot2::aes_string("x", "y")) + ggplot2::geom_tile(ggplot2::aes_string(fill = "value")) +
ggplot2::scale_fill_gradient(low = "white", high = "steelblue") +
ggplot2::facet_wrap(~time, labeller = "label_both")+ ggplot2::theme_void() + ggplot2::ggtitle(sub) +
ggplot2::theme(plot.title = ggplot2::element_text(hjust = 0.5))
Expand Down
7 changes: 6 additions & 1 deletion R/BacArena.R
Original file line number Diff line number Diff line change
Expand Up @@ -18,4 +18,9 @@ NULL


#' @import methods
NULL
NULL


.onAttach <- function(libname, pkgname) {
packageStartupMessage("BacArena paper: https://doi.org/10.1371/journal.pcbi.1005544\n Tutorial: https://CRAN.R-project.org/package=BacArena/vignettes/BacArena-Introduction.pdf\n Development and help: https://github.com/euba/bacarena")
}
35 changes: 17 additions & 18 deletions R/Organism.R
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,7 @@
#' @exportClass Organism
#' @import sybil
#' @importFrom stats na.omit
#' @rdname Organism
#'
#' @slot lbnd A numeric vector containing the lower bounds of the model structure.
#' @slot ubnd A numeric vector containing the upper bounds of the model structure.
Expand Down Expand Up @@ -70,7 +71,7 @@ setClass("Organism",

#' Constructor of the S4 class \code{\link{Organism-class}}
#'
#' @export Organism
#' @export
#' @name Organism-constructor
#'
#' @param model Object of class sybil::modelorg containging the genome sclae metabolic model
Expand Down Expand Up @@ -149,50 +150,42 @@ Organism <- function(model, algo="fba", ex="EX_", ex_comp=NA, csuffix="\\[c\\]",
########################################################################################################

setGeneric("lbnd", function(object){standardGeneric("lbnd")})
#' @export
setMethod("lbnd", "Organism", function(object){return(object@lbnd)})
setGeneric("ubnd", function(object){standardGeneric("ubnd")})
#' @export
setMethod("ubnd", "Organism", function(object){return(object@ubnd)})
#' @export
setGeneric("type", function(object){standardGeneric("type")})
#' @export
setMethod("type", "Organism", function(object){return(object@type)})
#' @export
setGeneric("medium", function(object){standardGeneric("medium")})
#' @export
setMethod("medium", "Organism", function(object){return(object@medium)})
#' @export
setGeneric("lpobj", function(object){standardGeneric("lpobj")})
#' @export
setMethod("lpobj", "Organism", function(object){return(object@lpobj)})
setGeneric("fbasol", function(object){standardGeneric("fbasol")})
#' @export
setMethod("fbasol", "Organism", function(object){return(object@fbasol)})
setGeneric("lyse", function(object){standardGeneric("lyse")})
#' @export
setMethod("lyse", "Organism", function(object){return(object@lyse)})
setGeneric("feat", function(object){standardGeneric("feat")})
#' @export
setMethod("feat", "Organism", function(object){return(object@feat)})
setGeneric("deathrate", function(object){standardGeneric("deathrate")})
#' @export
setMethod("deathrate", "Organism", function(object){return(object@deathrate)})
setGeneric("minweight", function(object){standardGeneric("minweight")})
#' @export
setMethod("minweight", "Organism", function(object){return(object@minweight)})
setGeneric("growtype", function(object){standardGeneric("growtype")})
#' @export
setMethod("growtype", "Organism", function(object){return(object@feat)})
setGeneric("kinetics", function(object){standardGeneric("kinetics")})
#' @export
setMethod("kinetics", "Organism", function(object){return(object@kinetics)})
setGeneric("speed", function(object){standardGeneric("speed")})
#' @export
setMethod("speed", "Organism", function(object){return(object@speed)})
setGeneric("model", function(object){standardGeneric("model")})
#' @export
setMethod("model", "Organism", function(object){return(object@model)})
setGeneric("maxweight", function(object){standardGeneric("maxweight")})
setMethod("maxweight", "Organism", function(object){return(object@maxweight)})
setGeneric("cellweight_mean", function(object){standardGeneric("cellweight_mean")})
setMethod("cellweight_mean", "Organism", function(object){return(object@cellweight_mean)})
setGeneric("cellarea", function(object){standardGeneric("cellarea")})
setMethod("cellarea", "Organism", function(object){return(object@cellarea)})
setGeneric("algo", function(object){standardGeneric("algo")})
setMethod("algo", "Organism", function(object){return(object@algo)})


########################################################################################################
###################################### METHODS #########################################################
Expand Down Expand Up @@ -611,6 +604,10 @@ setMethod(show, signature(object="Organism"), function(object){
#' @exportClass Bac
#' @rdname Bac
#'
#' @slot deathrate A numeric value giving the factor by which the growth should be reduced in every iteration (default (E.coli): 0.21 pg)
#' @slot minweight A numeric value giving the growth limit at which the organism dies. (default (E.coli): 0.083 pg)
#' @slot cellarea A numeric value indicating the surface that one organism occupies (default (E.coli): 4.42 mu_m^2)
#' @slot maxweight A numeric value giving the maximal dry weight of single organism (default (E.coli): 1.172 pg)
#' @slot chem A character vector indicating name of substance which is the chemotaxis attractant. Empty character vector if no chemotaxis.
setClass("Bac",
contains="Organism",
Expand Down Expand Up @@ -1031,7 +1028,9 @@ setMethod("cellgrowth", "Human", function(object, population, j, occupyM, fbasol
#' @param arena An object of class Arena defining the environment.
#' @param bacnum integer indicating the number of bacteria individuals per gridcell
#' @param sublb A vector containing the substance concentrations in the current position of the individual of interest.
#' @param cutoff value used to define numeric accuracy.
#' @param pcut A number giving the cutoff value by which value of objective function is considered greater than 0.
#' @param sec_obj character giving the secondary objective for a bi-level LP if wanted.
#' @return Returns the updated enivironment of the \code{arena} parameter with all new positions of individuals on the grid and all new substrate concentrations.
#' @details Human cell individuals undergo the step by step the following procedures: First the individuals are constrained with \code{constrain} to the substrate environment, then flux balance analysis is computed with \code{optimizeLP}, after this the substrate concentrations are updated with \code{consume}, then the cell growth is implemented with \code{cellgrowth}, the potential new phenotypes are added with \code{checkPhen}, finally the conditional function \code{lysis} is performed. Can be used as a wrapper for all important cell functions in a function similar to \code{simEnv}.
#' @seealso \code{\link{Human-class}}, \code{\link{Arena-class}}, \code{\link{simEnv}}, \code{constrain}, \code{optimizeLP}, \code{consume}, \code{cellgrowth}, \code{checkPhen} and \code{lysis}
Expand Down
5 changes: 3 additions & 2 deletions R/Stuff.R
Original file line number Diff line number Diff line change
Expand Up @@ -790,6 +790,7 @@ plotSpecActivity <- function(simlist, subs=list(), var_nr=10, spec_list=NULL, re
#' @rdname plotFVA
#'
#' @param fvares results of FVA results to plot, obtained from function fluxVarSim
#' @param mediac List with substances.
#' @details Returns ggplot objects
plotFVA = function(fvares, mediac){
for(i in 1:length(fvares)){
Expand All @@ -803,8 +804,8 @@ plotFVA = function(fvares, mediac){
rep(i,nrow(fvares[[i]]))))
}
colnames(basedat) = c("met","min","max","time")
fvag = ggplot2::ggplot(basedat, ggplot2::aes(x=time, y=min)) +
geom_ribbon(aes(ymin=min,ymax=max,color=met,fill=met),alpha=0.2) +
fvag = ggplot2::ggplot(basedat, ggplot2::aes_string(x="time", y="min")) +
ggplot2::geom_ribbon(ggplot2::aes_string(ymin="min",ymax="max",color="met",fill="met"),alpha=0.2) +
ggplot2::xlab("Time in h") + ggplot2::ylab("Population flux")
return(fvag)
}
Expand Down
15 changes: 2 additions & 13 deletions R/Substance.R
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,7 @@ setClass("Substance",
id = "character",
difunc = "character",
difspeed = "numeric",
advspeed = "numeric",
diffgeometry = "list",
pde = "character",
boundS = "numeric"
Expand Down Expand Up @@ -76,7 +77,7 @@ Substance <- function(n, m, smax, gridgeometry, difspeed=6.7e-6, advspeed=0, occ
}
}
if(ncol(diffmat)!=n && nrow(diffmat)!=m){
print(paste("arena dimensions :", arena@n, arena@m))
print(paste("arena dimensions :", n, m))
print(paste("diffmat dimension:", ncol(diffmat), nrow(diffmat)))
stop("Dimension of diffmat is invalid")
}
Expand Down Expand Up @@ -125,11 +126,6 @@ setMethod("difspeed", "Substance", function(object){return(object@difspeed)})
#' @param object An object of class Substance.
#' @details The diffusion is implemented by iterating through each cell in the grid and taking the cell with the lowest concentration in the Moore neighbourhood to update the concentration of both by their mean.
#' @seealso \code{\link{Substance-class}} and \code{\link{diffusePDE}}
#' @examples
#' arena <- Arena(n=100, m=100, stir=FALSE, Lx=0.025, Ly=0.025)
#' sub <- Substance(n=20,m=20,smax=40,name='test',difunc='r',
#' gridgeometry=arena@gridgeometry) #initialize test substance
#' diffuseR(sub)
setGeneric("diffuseR", function(object){standardGeneric("diffuseR")})
#' @export
#' @rdname diffuseR
Expand Down Expand Up @@ -183,13 +179,6 @@ setMethod("diffuseR", "Substance", function(object){
#' @param lrw A numeric value needed by solver to estimate array size (by default lwr is estimated in simEnv() by the function estimate_lrw())
#' @details Partial differential equation is solved to model 2d diffusion process in the arena.
#' @seealso \code{\link{Substance-class}} and \code{\link{diffuseR}}
#' @examples
#' arena <- Arena(n=100, m=100, stir=FALSE, Lx=0.025, Ly=0.025)
#' sub <- Substance(n=100,m=100,smax=0,name='test', difspeed=0.1, occupyM=arena@occupyM
#' gridgeometry=arena@gridgeometry) #initialize test substance
#' sub@diffmat[ceiling(100/2),ceiling(100/2)] <- 40
#' diffusePDE(sub, init_mat=matrix(sub@diffmat, nrow=arena@m, ncol=arena@n),
#' gridgeometry=arena@gridgeometry, tstep=arena@tstep)
setGeneric("diffusePDE", function(object, init_mat, gridgeometry, lrw=NULL, tstep){standardGeneric("diffusePDE")})
#' @export
#' @rdname diffusePDE
Expand Down
Loading

0 comments on commit ddbf7e2

Please sign in to comment.