## ----setup, include = TRUE, echo = FALSE-------------------------------------- knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ## ----globalode, eval = FALSE-------------------------------------------------- # #define the global model dynamics # global.dynamics <- function(t,x,parms){ # with(as.list(c(x, parms)), { # # dA <- rho_A*A + mu_BA*B + mu_CA*C - (mu_AB + mu_AC)*C # dB <- rho_B*B + mu_AB*A + mu_CB*C - (mu_BA + mu_BC)*B # dC <- rho_C*C + mu_AC*A + mu_BC*B - (mu_CA + mu_CB)*C # # return(list(c(dA,dB,dC))) # } # ) # } # ## ----odesim, eval = FALSE----------------------------------------------------- # #define simulation parameter set # pars <- c(rho_A = 0, rho_B = 0, rho_C = 0.1, # mu_AB = 0.2, mu_AC = 0, mu_BA = 0, # mu_CA = 0, mu_BC = 0.05, mu_CB = 0) # # #set initial values for cells # init.vals <- c(A = 100, B = 0, C = 0) # # #simulate data # library(deSolve) # # sim.data <- lsoda(y = init.vals, # times = 0:10, # func = global.dynamics, # parms = pars) # ## ----costode, eval = FALSE---------------------------------------------------- # #cost function for famos # cost.function <- function(parms, binary, data, inits){ # #simulate the data with the current parameter set # fit.data <- lsoda(y = inits, # times = 0:10, # func = global.dynamics, # parms = parms) # #calculate the aic # ls2 <- sum((sim.data[,-1] - fit.data[,-1])^2) # out.aic <- ls2 + 2*sum(binary == 1) # # return(out.aic) # } # ## ----odefamos, eval = FALSE--------------------------------------------------- # # library(FAMoS) # # start.vals <- c(rho_A = 0.1, rho_B = 0.1, rho_C = 0.1, # mu_AB = 0.1, mu_AC = 0.1, mu_BA = 0.1, # mu_CA = 0.1, mu_BC = 0.1, mu_CB = 0.1) # # famos.fit <- famos(init.par = start.vals, # fit.fn = cost.function, # homedir = tempdir(), # init.model.type = c("rho_A"), # data = sim.data, # inits = init.vals) # # print(famos.fit) # ## ----bwt, eval = FALSE-------------------------------------------------------- # library(MASS) # # attach(birthwt) # race <- factor(race, labels = c("white","black","other")) # ptd <- factor(ptl > 0) # levels(ftv)[-c(1:2)] <- "2+" # bwt <- data.frame(low = factor(low), age, lwt, race, # smoke = (smoke>0), ptd, ht = (ht>0), ui =(ui>0), ftv) # detach();rm(race,ptd,ftv) # ## ----glmfittingfunction, eval = FALSE----------------------------------------- # fit_func <- function(parms, data, binary){ # #First transform the parameter names into a formula. # #The to-be-fitted parameters are identified using the binary vector # fitted.pars <- names(parms[which(binary == 1)]) # glm_formula <- as.formula( # paste0("low ~ ", paste0(fitted.pars, collapse = "+")) # ) # #fit the logistic model using glm # out <- summary( # glm(glm_formula, # family = binomial(link = logit), # data = data # ) # ) # # #prepare output parameters, in this case only the names of the fitted # #parameters are relevant. As FAMoS also expects values with these, we will just # #return the value 1 for all parameters # out.par <- rep(1, length(fitted.pars)) # names(out.par) <- fitted.pars # # #return a list containing the first the selection criterion value and second the # #vector of the fitted parameters # return(list(SC = out$aic, params = out.par)) # } ## ----glminits, eval = FALSE--------------------------------------------------- # #first we define the available parameters. As the values are not important, we # #will just set them equal to 1. # inits <- c(age = 1, lwt = 1, race = 1, smoke = 1, ptd = 1, ht = 1, ui = 1, ftv = 1) # #Now, we calculate all possible interactions and name the corresponding vector # #entries accordingly # combinations <- combn(names(inits),2) # for(i in 1:ncol(combinations)){ # inits <- c(inits,1) # names(inits)[i + 8] <- paste0(combinations[1,i],":",combinations[2,i]) # } ## ----glmfamos, eval = FALSE--------------------------------------------------- # library(FAMoS) # famos.glm <- famos(init.par = inits, # fit.fn = fit_func, # init.model.type = names(inits[1:8]), # data = bwt, # use.optim = FALSE, # optim.runs = 1)