Heller Phar7383 Wk7 title: “Covariate modeling” date: “2025-03-05” output: html_document — #libraries
rm(list=ls())
library(ggplot2)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(PKNCA)## Warning: package 'PKNCA' was built under R version 4.4.2##
## Attaching package: 'PKNCA'## The following object is masked from 'package:stats':
##
## filterlibrary(knitr)#theme
my_theme<-function(x){theme_bw()+
theme(text = element_text(size=20))+
theme(axis.line.y = element_line(size = 2.0))+
theme(axis.line.x = element_line(size = 2.0))+
theme(axis.ticks = element_line(size = 1.5,colour="black"))+
theme(axis.ticks.length= unit(0.45, "cm"))+
theme(axis.title.y =element_text(vjust=1.2))+
theme(axis.title.x =element_text(vjust=-0.2))+
theme(axis.text=element_text(colour="black"))+
theme(panel.background = element_rect(fill ="white"))}#data import and analysis
unt456<-read.csv("C:\\Heller\\PHAR7383\\week7\\unt456a.csv",stringsAsFactors = F)
unt456sum<-unt456%>%filter(EVID==0)%>%group_by(DOSE,TIME)%>%summarise(cmean=mean(DV),stdev=sd(DV))## `summarise()` has grouped output by 'DOSE'. You can override using the
## `.groups` argument.#Exploratory plots
#Population Plot
ggplot(data=unt456,aes(TIME,DV,group=ID))+
geom_line(size=0.5)+
geom_point(size=1)+
scale_x_continuous(limits = c(0,24),breaks = c(0,1,2,4,6,8,12,24))+
theme_bw()+
my_theme()+
labs(x="Time after dose (hour)",y="Plasma concentration (ng/ml)")+
facet_wrap(vars(DOSE))## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## i Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
## i Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
#average plot with standard deviation (+/-1SD)
ggplot(data=unt456sum,aes(TIME,cmean,group=DOSE))+
geom_line(size=0.5)+
geom_point(size=2)+
scale_x_continuous(limits = c(0,24),breaks = c(0,1,2,4,6,8,12,24))+
geom_errorbar(aes(ymin=cmean-stdev, ymax=cmean+stdev), width=.2)+
theme_bw()+
my_theme()+
labs(x="Time after dose (hour)",y="Plasma concentration (ng/ml)")
#model building` #Run1- Base model- two compartment model with first
order absorption
library(xpose4)## Warning: package 'xpose4' was built under R version 4.4.2## Loading required package: latticerun1<-xpose.data(1,dir="C:\\Heller\\PHAR7383\\week7") ##
## Looking for NONMEM table files.
## Reading C:\Heller\PHAR7383\week7/sdtab1
## Reading C:\Heller\PHAR7383\week7/patab1
## Reading C:\Heller\PHAR7383\week7/catab1
## Reading C:\Heller\PHAR7383\week7/cotab1
## Table files read.
## Reading C:\Heller\PHAR7383\week7/run1.phi
##
## Looking for NONMEM simulation table files.
## No simulated table files read.#check the distribution of covariates (is there enough range in the covariates?) You don't need to check this every run.
#change.xvardef line keeps all covariates needed ready for xpose to plot. You need to run this line for xpose to recognize all our covarites in the dataset.
change.xvardef(run1, "covariates") <- c("WT","CRCL","SEX","CYP2D6","MAR","SHOE")
cov.hist(run1)
#Diagnostic plots
dv.vs.ipred(run1,type="p")
dv.vs.pred(run1,type="p")
cwres.vs.idv(run1,type="p")
cwres.vs.pred(run1,type="p")
ind.plots(run1)
ranpar.hist(run1)
#this code below gives plots of EBEs versus covariates. You need to specify covariates properly in the tables of model file for this code to work properly. ETA2 versus CRCL means clearance relationship versus CRCL and so on. You need to use sensibility of covariates and select few to conduct forward addition and backward elimination.
ranpar.vs.cov(run1)
#Read this instructions carefully #The model files for runs 2-4 are
provided. These will form the step 1 of forward addition. You want to
make sure all the table files are properly numbered and create job files
for each modelfile to run the models. # You need to create model files
for step 2 onwards for adding more than one covarites into the models.
#For backward elimination, take your full model and remove one covarite
and see how much increase in OBJFUN occured and is that significant at
alfa 0.01. #It is time for you to understand how the model files are
being changed to fit different models of covariates. You will need to do
this for the project! # Good luck! #you may see some errors and just
continue with the process. Your covariate model may solve those
problems. But, always mention what errors occured in your presentation.
Run2- Base model+WTonV2
#Diagnostic plots
run2<-xpose.data(2,dir="C:\\Heller\\PHAR7383\\week7") ##
## Looking for NONMEM table files.
## Reading C:\Heller\PHAR7383\week7/sdtab2
## Reading C:\Heller\PHAR7383\week7/patab2
## Reading C:\Heller\PHAR7383\week7/catab2
## Reading C:\Heller\PHAR7383\week7/cotab2
## Table files read.
## Reading C:\Heller\PHAR7383\week7/run2.phi
##
## Looking for NONMEM simulation table files.
## No simulated table files read.dv.vs.ipred(run2,type="p")
dv.vs.pred(run2,type="p")
cwres.vs.idv(run2,type="p")
cwres.vs.pred(run2,type="p")
ind.plots(run2)
ranpar.hist(run2)
ranpar.vs.cov(run2)
Run3- Base model+WTonV3
run3<-xpose.data(3,dir="C:\\Heller\\PHAR7383\\week7") ##
## Looking for NONMEM table files.
## Reading C:\Heller\PHAR7383\week7/sdtab3
## Reading C:\Heller\PHAR7383\week7/patab3
## Reading C:\Heller\PHAR7383\week7/catab3
## Reading C:\Heller\PHAR7383\week7/cotab3
## Table files read.
## Reading C:\Heller\PHAR7383\week7/run3.phi
##
## Looking for NONMEM simulation table files.
## No simulated table files read.dv.vs.ipred(run3,type="p")
dv.vs.pred(run3,type="p")
cwres.vs.idv(run3,type="p")
cwres.vs.pred(run3,type="p")
ind.plots(run3)
ranpar.hist(run3)
ranpar.vs.cov(run3)
Run4- Base model+CRCL on CL (is it justifed to test this relationship either by plots or biological rationale? Comment on this in your presentation of results)
run4corr<-xpose.data(4,dir="C:\\Heller\\PHAR7383\\week7") ##
## Looking for NONMEM table files.
## Reading C:\Heller\PHAR7383\week7/sdtab4
## Reading C:\Heller\PHAR7383\week7/patab4
## Reading C:\Heller\PHAR7383\week7/catab4
## Reading C:\Heller\PHAR7383\week7/cotab4
## Table files read.
## Reading C:\Heller\PHAR7383\week7/run4.phi
##
## Looking for NONMEM simulation table files.
## No simulated table files read.dv.vs.ipred(run4corr,type="p")
dv.vs.pred(run4corr,type="p")
cwres.vs.idv(run4corr,type="p")
cwres.vs.pred(run4corr,type="p")
ind.plots(run4corr)
ranpar.hist(run4corr)
ranpar.vs.cov(run4corr)
#Base model+WTonV2 and WTonV3
run5<-xpose.data(5,dir="C:\\Heller\\PHAR7383\\week7") ##
## Looking for NONMEM table files.
## Reading C:\Heller\PHAR7383\week7/sdtab5
## Reading C:\Heller\PHAR7383\week7/patab5
## Reading C:\Heller\PHAR7383\week7/catab5
## Reading C:\Heller\PHAR7383\week7/cotab5
## Table files read.
## Reading C:\Heller\PHAR7383\week7/run5.phi
##
## Looking for NONMEM simulation table files.
## No simulated table files read.dv.vs.ipred(run5,type="p")
dv.vs.pred(run5,type="p")
cwres.vs.idv(run5,type="p")
cwres.vs.pred(run5,type="p")
ind.plots(run5)
ranpar.hist(run5)
ranpar.vs.cov(run5)
#Base model+WTonV2 and WTonV3 and CRCL on CL
run6<-xpose.data(6,dir="C:\\Heller\\PHAR7383\\week7") ##
## Looking for NONMEM table files.
## Reading C:\Heller\PHAR7383\week7/sdtab6
## Reading C:\Heller\PHAR7383\week7/patab6
## Reading C:\Heller\PHAR7383\week7/catab6
## Reading C:\Heller\PHAR7383\week7/cotab6
## Table files read.
## Reading C:\Heller\PHAR7383\week7/run6.phi
##
## Looking for NONMEM simulation table files.
## No simulated table files read.dv.vs.ipred(run6,type="p")
dv.vs.pred(run6,type="p")
cwres.vs.idv(run6,type="p")
cwres.vs.pred(run6,type="p")
ind.plots(run6)
ranpar.hist(run6)
ranpar.vs.cov(run6)
