Appendix A. Supplementary information
Toxicokinetic-toxicodynamic modelling of survival of zebrafish larvae in exposures to metals: model assumptions, parameterization data requirements and predictive power
Yongfei Gao, Jianfeng Feng#,and Lin Zhu
Key laboratory of Pollution process and Environmental Criteria of Ministry of Education and Tianjin Key Laboratory of Environmental Remediation and Pollution Control, College of Environmental Science and Engineering, Nankai University, Tianjin300071, China
#Corresponding author. phone: +86-022-85358315;
E-mail: (Jianfeng Feng)
Number of pages: 7
R Code
The toxicokinetic and toxicodynamic parameters were estimated using the software R 3.1.2. The code for estimating parameters of metal mixtures is provided below.
Example: Code for Cu Toxicokinetic-Toxicodynamic-SD model
# the model equations: #
#------#
model<-function(t,state,parameters){
with(as.list(c(state,parameters)),{# unpack the state variables, parameters
dCt<-Jmax*M*KMBL/(M*KMBL+H*KHBL+1)- Kout*Ct
dHazard<-kk*max(Ct-CT, 0)
list(c(dCt,dHazard)) # the output, packed as a list
})
}
#------#
# the model parameters: #
#------#
parameters<-c(Jmax=15.46, # parameter values
KMBL=2.63e+05,
M=7.87e-07,
KHBL=1.78e+05,
Kout=0.0531,
H=1e-07,
CT=15.86,
kk=7.04e-03
)
#=====Initial=====#
state <-c(
Ct=14.85,
Hazard=0
)
# RUNNING the model: #
#------#
times <-seq(0,94,4) # time: from 0 to 168 hours, steps of 0.5 hours
# integrate the model
require(deSolve) # package with the integration routine:
out <-as.data.frame(ode(state,times,model,parameters))
Hazard<-out$Hazard
S<-exp(-Hazard)
Example: Code for Cu Toxicokinetic-Toxicodynamic-IT model
#------#
# the model equations: #
#------#
model<-function(t,state,parameters){
with(as.list(c(state,parameters)),{# unpack the state variables, parameters
dCt<-Jmax*M*KMBL/(M*KMBL+H*KHBL+1)- Kout*Ct
list(c(dCt)) # the output, packed as a list
})
}
#------#
# the model parameters: #
#------#
parameters<-c(Jmax=15.46, # parameter values
KMBL=2.63e+05,
M=7.87e-07,
KHBL=1.78e+05,
Kout=0.0531,
H=1e-07
)
#=====Initial=====#
state <-c(
Ct=14.85
)
# RUNNING the model: #
#------#
times <-seq(0,100,4) # time: from 0 to 168 hours, steps of 0.5 hours
# integrate the model
require(deSolve) # package with the integration routine:
out <-as.data.frame(ode(state,times,model,parameters))
C<-out$Ct
# IT model
S<-1/(1+(C/44.84)^2.66)
Example: Code for Cu-Cd Toxicokinetic-Toxicodynamic-SD model
#------#
# the model equations: #
#------#
model<-function(t,state,parameters){
with(as.list(c(state,parameters)),{# unpack the state variables, parameters
dCuQ<-Jmax1*Cu*KCuBL/(Cu*KCuBL+Cd*KCdBL+H*KHBL+1)- Kout1*CuQ
dCdQ<-Jmax2*Cd*KCdBL/(Cu*KCuBL+Cd*KCdBL+H*KHBL+1)- Kout2*CdQ
Ct<-1.0*CuQ+0.8*CdQ
dHazard<-kk*max(Ct-threshold, 0)
list(c(dCuQ,dCdQ,Ct,dHazard)) # the output, packed as a list
})
}
#------#
# the model parameters: #
#------#
parameters1<-c(Jmax1=15.46, # parameter values
KCuBL=2.63e+05,
KHBL=1.78e+05,
Kout1=0.0531,
Jmax2=12.02,
KCdBL=4.63e+05,
Kout2=0.0549,
H=1e-07,
Cd=4.45e-07,
Cu=1.57e-07,
threshold=15.86,
kk=0.0072
)
parameters2<-c(Jmax1=15.46, # parameter values
KCuBL=2.63e+05,
KHBL=1.78e+05,
Kout1=0.0531,
Jmax2=12.02,
KCdBL=4.63e+05,
Kout2=0.0549,
H=1e-07,
Cd=4.45e-07,
Cu=1.57e-06,
threshold=15.86,
kk=0.0072
)
parameters3<-c(Jmax1=15.46, # parameter values
KCuBL=2.63e+05,
KHBL=1.78e+05,
Kout1=0.0531,
Jmax2=12.02,
KCdBL=4.63e+05,
Kout2=0.0549,
H=1e-07,
Cd=4.45e-07,
Cu=1.57e-05,
threshold=15.86,
kk=0.0072
)
#=====Initial=====#
state <-c(
CuQ=14.85,
CdQ=0,
Ct=14.85,
Hazard=0
)
#------#
# RUNNING the model: #
#------#
times <-seq(0,24,0.1) # time: from 0 to 168 hours, steps of 0.5 hours
# integrate the model
require(deSolve) # package with the integration routine:
out1 <-as.data.frame(ode(state,times,model,parameters1))
out2 <-as.data.frame(ode(state,times,model,parameters2))
out3 <-as.data.frame(ode(state,times,model,parameters3))
Hazard1<-out1$Hazard
Hazard2<-out2$Hazard
Hazard3<-out3$Hazard
S1<-exp(-Hazard1)
S2<-exp(-Hazard2)
S3<-exp(-Hazard3)
Example: Code for Cu-Cd Toxicokinetic-Toxicodynamic-IT model
#------#
# the model equations: #
#------#
model<-function(t,state,parameters){
with(as.list(c(state,parameters)),{# unpack the state variables, parameters
dCuQ<-Jmax1*Cu*KCuBL/(Cu*KCuBL+Cd*KCdBL+H*KHBL+1)- Kout1*CuQ
dCdQ<-Jmax2*Cd*KCdBL/(Cu*KCuBL+Cd*KCdBL+H*KHBL+1)- Kout2*CdQ
list(c(dCuQ,dCdQ)) # the output, packed as a list
})
}
#------#
# the model parameters: #
#------#
parameters1<-c(Jmax1=15.46, # parameter values
KCuBL=2.63e+05,
KHBL=1.78e+05,
Kout1=0.0531,
Jmax2=12.02,
KCdBL=4.63e+05,
Kout2=0.0549,
H=1e-07,
Cd=4.45e-07,
Cu=1.57e-07
)
parameters2<-c(Jmax1=15.46, # parameter values
KCuBL=2.63e+05,
KHBL=1.78e+05,
Kout1=0.0531,
Jmax2=12.02,
KCdBL=4.63e+05,
Kout2=0.0549,
H=1e-07,
Cd=4.45e-07,
Cu=1.57e-06
)
parameters3<-c(Jmax1=15.46, # parameter values
KCuBL=2.63e+05,
KHBL=1.78e+05,
Kout1=0.0531,
Jmax2=12.02,
KCdBL=4.63e+05,
Kout2=0.0549,
H=1e-07,
Cd=4.45e-07,
Cu=1.57e-05
)
#=====Initial=====#
state <-c(
CuQ=14.85,
CdQ=0
)
#------#
# RUNNING the model: #
#------#
times <-seq(0,24,0.1) # time: from 0 to 168 hours, steps of 0.5 hours
# integrate the model
require(deSolve) # package with the integration routine:
out1 <-as.data.frame(ode(state,times,model,parameters1))
out2 <-as.data.frame(ode(state,times,model,parameters2))
out3 <-as.data.frame(ode(state,times,model,parameters3))
CuQ1<-out1$CuQ
CdQ1<-out1$CdQ
Ct1<-1.0*CuQ1+0.8*CdQ1
CuQ2<-out2$CuQ
CdQ2<-out2$CdQ
Ct2<-1.0*CuQ2+0.8*CdQ2
CuQ3<-out3$CuQ
CdQ3<-out3$CdQ
Ct3<-1.0*CuQ3+0.8*CdQ3
# IT model
S1<-1/(1+(Ct1/91.62)^3.24)
S2<-1/(1+(Ct2/91.62)^3.24)
S3<-1/(1+(Ct3/91.62)^3.24)
S1