Comments (6)
andrewhooker
commented on August 17, 2024
Hi Ankur,
This can certainly be implemented! You just need to define different groups with different designs in the trial. See for example https://andrewhooker.github.io/PopED/articles/intro-poped.html. You can post example code here if you like and I can try and help.
Best regards,
Andy
from poped.
Anks2030
commented on August 17, 2024
Hi Andrew,
Thanks for looking into this.
In our trial, we have N=200 in treatment arm until wk-52. Post-52 week we have PK at D1, D2, D3, D4 and D7, D14. And post-52WK the collection is for only N=20.
I dont know how to implement this switch of N=200 to N=20 and have an overall influence of design?
So far I am just assuming 200 for entire duration and have even tried to evaluate the design only till WK-52 (using just N=200)
Please have code for entire study for N=200
#' Define the ODE system
PK.2.comp.sc.ode <- function(Time, State, Pars){
with(as.list(c(State, Pars)), {
CL=CL*(WT/76.8)^(WT_CL)
V1=V1*(WT/76.8)^(WT_V1)
dA1 <- -KA*A1
dA2 <- KA*A1 + A3* Q/V2 -A2*(CL/V1+Q/V1)
dA3 <- A2* Q/V1-A3* Q/V2
return(list(c(dA1, dA2, dA3)))
})
}
#' define the initial conditions and the dosing
ff.PK.2.comp.sc.md.ode <- function(model_switch, xt, parameters, poped.db){
with(as.list(parameters),{
A_ini <- c(A1=0, A2=0, A3=0)
times_xt <- drop(xt)
dose_times = c(seq(from=0,to=336, by=84),seq(from=336,to=8568,by=TAU))
eventdat <- data.frame(var = c("A1"),
time = dose_times,
value = c(DOSE), method = c("add"))
times <- sort(c(times_xt,dose_times))
times <- unique(times) # remove duplicates
out <- ode(A_ini, times, PK.2.comp.sc.ode, parameters, events = list(data = eventdat))#atol=1e-13,rtol=1e-13)
y = (out[, "A2"]/(V1/Favail))+BL
y=y[match(times_xt,out[,"time"])]
y=cbind(y)
return(list(y=y,poped.db=poped.db))
})
}
#' parameter definition function
#' names match parameters in function ff
fg <- function(x,a,bpop,b,bocc){
parameters=c( CL=bpop[1]*exp(b[1]),
V1=bpop[2]*exp(b[2]),
KA=bpop[3],
Q= bpop[4],
V2=bpop[5],
Favail=bpop[6],
WT_CL=bpop[7],
WT_V1=bpop[8],
BL=bpop[9],
DOSE=a[1],
TAU=a[2],
WT= a[3])
return( parameters )
}
-- Residual unexplained variablity (RUV) function
-- Additive + Proportional
feps <- function(model_switch,xt,parameters,epsi,poped.db){
returnArgs <- do.call(poped.db$model$ff_pointer,list(model_switch,xt,parameters,poped.db))
y <- returnArgs[[1]]
poped.db <- returnArgs[[2]]
y = y*(1+epsi[,1])+epsi[,2]
return(list( y= y,poped.db =poped.db ))
}
Original Design
#CL: L/h; V1,V2: L; Q: L/h; V2: L; BL: g/L, KA= h-1
#create poped database original design
poped.db <- create.poped.database(ff_fun="ff.PK.2.comp.sc.md.ode",
fg_fun="fg",
fError_file="feps",
bpop=c(CL=0.00633,V1=8.75,KA=0.02108,Q= 0.0344, V2= 1.8, Favail=0.668,WT_CL=0.796, WT_V1=1.1, BL=4),# TV 0.007475 L/h, L
notfixed_bpop=c(1,1,1,0,1,0,0,0,0),# decides which tested and which are not
d=c(CL=0.282,V1=0.569), # decides the variances of BSV
notfixed_d = c(1,1),
sigma=c(0.00662,0), # decides the variances of RV prop, additive
notfixed_sigma=c(1,0),
m=1, #number of groups
groupsize=200,
#0,7,28,35,56,59,63,70,77,84
xt= c(0,336,671,2016,2688,3360,6048,7392,8568,8592,8616,8640,8664,8736,8904),
minxt=c(0,312,648,1992,2664,3336,6024,7368,8544,8568,8592,8616,8640,8712,8880),
maxxt=c(0,360,696,2040,2712,3384,6072,7416,8592,8616,8640,8664,8688,8760,8928),
bUseGrouped_xt=1,
a=c(DOSE=0.15*76.8,TAU=168, WT=76.8))
#' plot intial design with BSV and RUV in model
#plot_model_prediction(poped.db,IPRED=T,DV=T)
dat_Original <- model_prediction(poped.db)
plot_model_prediction(poped.db,DV=T,
# separate.groups=T,, groupsize_sim=100,PI=T
sample.times=F,alpha.sample.times.DV.points = 0.05,
model_num_points = 1000)+# , PI=T
labs(x = "Time from first dose (Weeks)", y="IgG Conc.(g/L)")+
theme_bw()+ geom_vline(xintercept=336,lty="dashed")+
geom_point(data = dat_Original, aes(Time, PRED ), size=2, alpha=1/2, colour="blue")+#colour="Before Optimization"
ggtitle("Original Design")+coord_cartesian(ylim=c(4,12))+
scale_y_continuous(breaks=sort(c(seq(0,12,1))))+
scale_x_continuous(breaks = seq(0,8928,168), labels=c(seq(1,54,1)))+
theme_bw()
Evaluate Original Design
(dso <- evaluate_design(poped.db))
from poped.
andrewhooker
commented on August 17, 2024
Hi Ankur,
Here is am implementation of what you are trying to do.
First define the model:
library(tidyverse)
library(deSolve)
library(PopED)
packageVersion("PopED") #0.5.0
#' Define the ODE system
PK.2.comp.sc.ode <- function(Time, State, Pars){
with(as.list(c(State, Pars)), {
CL=CL*(WT/76.8)^(WT_CL)
V1=V1*(WT/76.8)^(WT_V1)
dA1 <- -KA*A1
dA2 <- KA*A1 + A3* Q/V2 -A2*(CL/V1+Q/V1)
dA3 <- A2* Q/V1-A3* Q/V2
return(list(c(dA1, dA2, dA3)))
})
}
#' define the initial conditions and the dosing
ff.PK.2.comp.sc.md.ode <- function(model_switch, xt, parameters, poped.db){
with(as.list(parameters),{
A_ini <- c(A1=0, A2=0, A3=0)
times_xt <- drop(xt)
dose_times = c(seq(from=0,to=336, by=84),seq(from=336,to=8568,by=TAU))
eventdat <- data.frame(var = c("A1"),
time = dose_times,
value = c(DOSE), method = c("add"))
times <- sort(c(times_xt,dose_times))
times <- unique(times) # remove duplicates
out <- ode(A_ini, times, PK.2.comp.sc.ode, parameters, events = list(data = eventdat))#atol=1e-13,rtol=1e-13)
y = (out[, "A2"]/(V1/Favail))+BL
y=y[match(times_xt,out[,"time"])]
y=cbind(y)
return(list(y=y,poped.db=poped.db))
})
}
#' parameter definition function
#' names match parameters in function ff
fg <- function(x,a,bpop,b,bocc){
parameters=c( CL=bpop[1]*exp(b[1]),
V1=bpop[2]*exp(b[2]),
KA=bpop[3],
Q= bpop[4],
V2=bpop[5],
Favail=bpop[6],
WT_CL=bpop[7],
WT_V1=bpop[8],
BL=bpop[9],
DOSE=a[1],
TAU=a[2],
WT= a[3])
return( parameters )
}
## -- Residual unexplained variablity (RUV) function
## -- Additive + Proportional
feps <- function(model_switch,xt,parameters,epsi,poped.db){
returnArgs <- do.call(poped.db$model$ff_pointer,list(model_switch,xt,parameters,poped.db))
y <- returnArgs[[1]]
poped.db <- returnArgs[[2]]
y = y*(1+epsi[,1])+epsi[,2]
return(list( y= y,poped.db =poped.db ))
}
Add a design with design with extra observations after 52 wks in 20 patients.
poped_db_1 <- create.poped.database(
ff_fun="ff.PK.2.comp.sc.md.ode",
fg_fun="fg",
fError_file="feps",
bpop=c(CL=0.00633,V1=8.75,KA=0.02108,
Q= 0.0344,
V2= 1.8,
Favail=0.668,WT_CL=0.796, WT_V1=1.1, BL=4),# TV 0.007475 L/h, L
notfixed_bpop=c(1,1,1,
0,
1,
0,0,0,0),# decides which tested and which are not
d=c(CL=0.282,V1=0.569), # decides the variances of BSV
notfixed_d = c(1,1),
sigma=c(0.00662,0), # decides the variances of RV prop, additive
notfixed_sigma=c(1,0),
m=2, #number of groups
groupsize=c(20,180),
#0,7,28,35,56,59,63,70,77,84
# Here we have 20 individuals that have sampling like the rest until wk 52 and
# then more samples post wk52.
# 180 have only sampling until wk 52.
xt= list(c(0,336,671,2016,2688,3360,6048,7392,8568,8592,8616,8640,8664,8736,8904),
c(0,336,671,2016,2688,3360,6048,7392,8568)),
minxt=list(c(0,312,648,1992,2664,3336,6024,7368,8544,8568,8592,8616,8640,8712,8880),
c(0,312,648,1992,2664,3336,6024,7368,8544)),
maxxt=list(c(0,360,696,2040,2712,3384,6072,7416,8592,8616,8640,8664,8688,8760,8928),
c(0,360,696,2040,2712,3384,6072,7416,8592)),
#bUseGrouped_xt=1,
a=list(c(DOSE=0.15*76.8,TAU=168, WT=76.8),
c(DOSE=0.15*76.8,TAU=168, WT=76.8))
)
Evaluate design and plot design.
(p1 <- plot_model_prediction(poped_db_1,separate.groups = T,model_num_points = 1000))
(des1 <- evaluate_design(poped_db_1))
Create a design without the extra observations after 52 weeks:
poped_db_2 <- create.poped.database(
ff_fun="ff.PK.2.comp.sc.md.ode",
fg_fun="fg",
fError_file="feps",
bpop=c(CL=0.00633,V1=8.75,KA=0.02108,
Q= 0.0344,
V2= 1.8,
Favail=0.668,WT_CL=0.796, WT_V1=1.1, BL=4),# TV 0.007475 L/h, L
notfixed_bpop=c(1,1,1,
0,
1,
0,0,0,0),# decides which tested and which are not
d=c(CL=0.282,V1=0.569), # decides the variances of BSV
notfixed_d = c(1,1),
sigma=c(0.00662,0), # decides the variances of RV prop, additive
notfixed_sigma=c(1,0),
m=1, #number of groups
groupsize=c(200),
#0,7,28,35,56,59,63,70,77,84
# Here we have 20 individuals that have sampling like the rest until wk 52 and
# then more samples post wk52.
# 180 have only sampling until wk 52.
xt= c(0,336,671,2016,2688,3360,6048,7392,8568),
minxt=c(0,312,648,1992,2664,3336,6024,7368,8544),
maxxt=c(0,360,696,2040,2712,3384,6072,7416,8592),
#bUseGrouped_xt=1,
a=c(DOSE=0.15*76.8,TAU=168, WT=76.8)
)
Evaluate design and plot design.
(p2 <- plot_model_prediction(poped_db_2,separate.groups = T,model_num_points = 1000))
(des2 <- evaluate_design(poped_db_2))
There is a big difference in parameter RSE between the two designs
tibble::tibble(names(des1$rse),des1$rse,des2$rse)
There is also a big difference in efficiency (52% more individuals in the design without extra observations after 52 wks are needed to achieve same information as the extended design)
efficiency(ofv_init = des2$ofv, ofv_final = des1$ofv, poped_db = poped_db_1)
from poped.
Anks2030
commented on August 17, 2024
Hi Andy,
Thanks for taking time to go over this challenge
I have a question.
I am trying hard to understand if the code is capturing the design appropriately
I am attaching a single design figure (total duration is 1year) in which we have N=180 only till WK-52 and after WK-52 we have ONLY N=20
I trying to understand if this design is applied in below code you shared?
m=2, #number of groups
groupsize=c(20,180),
#0,7,28,35,56,59,63,70,77,84
Here we have 20 individuals that have sampling like the rest until wk 52 and
then more samples post wk52.
180 have only sampling until wk 52.
xt= list(c(0,336,671,2016,2688,3360,6048,7392,8568,8592,8616,8640,8664,8736,8904),
c(0,336,671,2016,2688,3360,6048,7392,8568)),
The output is Group 1 has N=20 only till WK54 and Group 2 is N=180 till WK-52
I wish to see one group which has overall impact of design (N=180 patients only till WK52 and after WK52 I have only 20 subjects for rich PK)
I am not sure if I am interpreting it correctly please advice
Really appreciate you spreading the word for PopED. Very useful tool.
Thanks,
-Ankur
from poped.
andrewhooker
commented on August 17, 2024
Hi
The design I implemented has 180 individuals with measurements between 0 and 52 weeks, and 20 individuals with those measurements PLUS additional measurements after 52 weeks. I assumed that the 20 individuals would be a part of the 200 total that are studies under 52 weeks and then 20 of those individuals are studied for an additional period.
from poped.
Anks2030
commented on August 17, 2024
Hi Andy,
Thank you very much it is clear now. Very useful and I really appreciate for your time and effort.
Do you have the analytical solution for the 2-compartment model.
I will try to do the optimization for this design and seems it is taking forever for ODE's
Thanks,
-Ankur
from poped.
Related Issues (20)
- Manually updating the number of samples as an alternative to opt_samps in R? HOT 1
- Plotting optimized output seems to used wrong times HOT 1
- Where is the optimization parallelized HOT 3
- unused argument in plot_model_prediction() HOT 2
- parallel computation for plot_efficiency_of_windows() HOT 4
- problem with tests run by AppVeyor HOT 2
- Warning "Problems inverting the matrix. Results could be misleading." HOT 2
- Different results when running the optimization R-script in popED in Windows and Mac HOT 2
- Integrating an IV loading dose followed by SC maintenance dose using popEd HOT 1
- Illogical OFV when xts have a lot of significant digits HOT 2
- D-optimal Design for a trial comparing three different drugs
- Gradient function HOT 2
- Handling derived outputs
- `evaluate_fim_map` documentation
- Error in svd(mat) : infinite or missing values in 'x'
- Define a model with lag time HOT 5
- Take an example of the pk.3.comp.olal.ode model
- multiple dose full TMDD model in POPED HOT 3
- Two questions HOT 3
Recommend Projects
-
React
A declarative, efficient, and flexible JavaScript library for building user interfaces.
-
Vue.js
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
-
Typescript
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
-
TensorFlow
An Open Source Machine Learning Framework for Everyone
-
Django
The Web framework for perfectionists with deadlines.
-
Laravel
A PHP framework for web artisans
-
D3
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
-
Recommend Topics
-
javascript
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
-
web
Some thing interesting about web. New door for the world.
-
server
A server is a program made to process requests and deliver data to clients.
-
Machine learning
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
-
Visualization
Some thing interesting about visualization, use data art
-
Game
Some thing interesting about game, make everyone happy.
Recommend Org
-
Facebook
We are working to build community through open source technology. NB: members must have two-factor auth.
-
Microsoft
Open source projects and samples from Microsoft.
-
Google
Google ❤️ Open Source for everyone.
-
Alibaba
Alibaba Open Source for everyone
-
D3
Data-Driven Documents codes.
-
Tencent
China tencent open source team.
from poped.