if (FALSE) {
library(magrittr)
library(dplyr)
library(csv)
library(nlme)
x <- Theoph
# mixed effects model
m1 <- nlme(
conc ~ SSfol(Dose, Time, lKe, lKa, lCl),
data = x,
fixed = lKe + lKa + lCl ~ 1,
random = lKe + lKa + lCl ~ 1
)
# some numeric and categorical properties
names(x) <- tolower(names(x))
x %<>% mutate(arm = ifelse(as.numeric(as.character(subject)) %% 2 == 0, 1, 2))
x %<>% mutate(site = ifelse(as.numeric(as.character(subject)) < 6, 1, 2))
x %<>% mutate(cohort = ifelse(as.numeric(as.character(subject)) %in% c(1:2,6:8), 1,2))
x %<>% mutate(pred = predict(m1,level = 0) %>% signif(4))
x %<>% mutate(ipred = predict(m1) %>% signif(4))
x %<>% mutate(res = residuals(m1) %>% signif(4))
x %<>% mutate(sres = residuals(m1, type = 'pearson') %>% signif(4))
r <- ranef(m1) %>% signif(4)
r$subject <- rownames(r)
x %<>% left_join(r)
# metadata
attr(x$subject,'label') <- 'subject identifier'
attr(x$wt,'label') <- 'subject weight'
attr(x$dose,'label') <- 'theophylline dose'
attr(x$time,'label') <- 'time since dose administration'
attr(x$conc,'label') <- 'theophylline concentration'
attr(x$arm,'label') <- 'trial arm'
attr(x$site,'label') <- 'investigational site'
attr(x$cohort,'label') <- 'recruitment cohort'
attr(x$pred,'label') <- 'population-predicted concentration'
attr(x$ipred,'label') <- 'individual-predicted concentration'
attr(x$res,'label') <- 'residuals'
attr(x$sres,'label') <- 'standardized residuals'
attr(x$lKe,'label') <- 'natural log of elimination rate constant'
attr(x$lKa,'label') <- 'natural log of absorption rate constant'
attr(x$lCl,'label') <- 'natural log of clearance'
attr(x$subject,'guide') <- '....'
attr(x$wt,'guide') <- 'kg'
attr(x$dose,'guide') <- 'mg/kg'
attr(x$time,'guide') <- 'h'
attr(x$conc,'guide') <- 'mg/L'
attr(x$arm,'guide') <- '//1/Arm A//2/Arm B//'
attr(x$site,'guide') <- '//1/Site 1//2/Site 2//'
attr(x$cohort,'guide') <- '//1/Cohort 1//2/Cohort 2//'
attr(x$pred,'guide') <- 'mg/L'
attr(x$ipred,'guide') <- 'mg/L'
attr(x$lKe,'reference') <- 0
attr(x$lKa,'reference') <- 0
attr(x$lCl,'reference') <- 0
attr(x$res,'reference') <- 0
attr(x$sres,'reference') <- '//-1.96//1.96//'
attr(x$subject,'symbol') <- 'ID_i'
attr(x$wt,'symbol') <- 'W_i'
attr(x$dose,'symbol') <- 'A_i'
attr(x$time,'symbol') <- 't_i,j'
attr(x$conc,'symbol') <- 'C_i,j'
attr(x$arm,'symbol') <- 'Arm_i'
attr(x$site,'symbol') <- 'Site_i'
attr(x$cohort,'symbol') <- 'Cohort_i'
attr(x$pred,'symbol') <- 'C_pred_p'
attr(x$ipred,'symbol') <- 'C_pred_i'
attr(x$res,'symbol') <- '\\epsilon'
attr(x$sres,'symbol') <- '\\epsilon_st'
attr(x$lKe,'symbol') <- 'ln(K_e.)'
attr(x$lKa,'symbol') <- 'ln(K_a.)'
attr(x$lCl,'symbol') <- 'ln(Cl_c./F)'
x %>% unpack %>% as.csv('theoph.csv')
}
# \dontshow{
if (FALSE) {
y <- x
y[] <- lapply(y, as.character)
y[] <- lapply(y, as.numeric)
y$arm <- as.factor(y$arm)
y$site <- as.factor(y$site)
y$subject <- as.factor(y$subject)
y %>% metaplot(conc)
y %>% metaplot(site)
y %>% metaplot(wt, arm)
y %>% metaplot(arm, wt)
y %>% metaplot(arm, wt,site)
y %>% metaplot(conc, time)
y %>% metaplot(arm, site)
y %>% metaplot(conc, time, subject)
y %>% metaplot(conc, time, , subject)
y %>% metaplot(conc, time, subject, site)
y %>% metaplot(conc, time, subject, site, arm)
y %>% metaplot(lKe, lKa, lCl)
y %>% scatter(conc, ipred, time)
y %>% scatter(conc, ipred, time, subject)
x %>% metaplot(conc, ipred, time, subject, colors = 'black', points = c(T,F), lines = c(F,T))
y %>% scatter(conc, ipred, time, site, arm)
y %<>% mutate(time = ifelse(time > 15, NA, time))
y %>% scatter(conc, ipred, time, site, arm)
}# }
Run the code above in your browser using DataLab