Varying RUV while maintaining IIV #983
Replies: 7 comments 1 reply
-
Hi @thinguyen94. If I'm interpreting your question, you'd like to simulate concentrations with varying levels of RUV, while keeping IIV and covariate terms consistent. Since RUV is tacked on at the end of the simulation, you could always create a separate output for each assay, and use different variance terms for each one. This has the advantage of only having to run a simulation once and guarantees the same subjects are used between the two outputs (in terms of randomly sampled ETAs and baseline characteristics). See example below. I edited your model code to include two more sigma terms and to capture two separate DV terms (one for each assay). ## Load Packages
library(mrgsolve)
library(ggplot2)
library(tidyr)
## Set plot theme
theme_set(theme_classic())
## Model Code
code<- '
$INIT // Initial Conditions for Compartments
CMT1 = 0, // Central Compartment
CMT2 = 0, // Peripheral Compartment
$SET // Set Differential Equation Solver Options
atol = 1e-8, rtol = 1e-8
maxsteps = 100000
$PARAM // Population parameters
POPCL = 0.0493, // Clearance
POPV1 = 0.833, // Central Volume
POPV2 = 0.833, // Peripheral Volume
POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
// Covariate effects
BW_CL = 1.340, // Exponent for Birthweight
PNA_CL = 0.213, // Slope of PNA
IBU_CL = 0.838, // Effect of Ibuprofen Co-administration
CW_V = 0.919 // Exponent for Current Weight
// Default Covariate Values
BW = 1750, // g
CW = 1760, // g
PNA = 2, // days
NSAID = 0 // no co-administration of ibuprofen
$OMEGA // Population Parameter Variability
@labels ETA1
0.0899
// Added variance terms for second assay
$SIGMA // Residual Unexplained Variability
@labels ERRPROP1 ERRADD1 ERRPROP2 ERRADD2
0.0614 0.267 0.02 0.1
$MAIN // Individual Parameter Values
double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) * exp(ETA1);
double V1 = (POPV1 * pow(CW/1760,CW_V));
double V2 = V1;
double Q = (POPQFRAC * CL);
if(NSAID == 1) {
CL = CL * IBU_CL;
Q = Q * IBU_CL;
}
$ODE // Differential Equations
double k10 = CL/V1;
double k12 = Q/V1;
double k21 = Q/V2;
dxdt_CMT1 = - k12*CMT1 + k21*CMT2 - k10*CMT1;
dxdt_CMT2 = k12*CMT1 - k21*CMT2;
$TABLE //Determines Values and Includes in Output
double IPRED = CMT1/V1;
// Capturing DV for both assays below
double DV_ASSAY1 = IPRED*(1 + ERRPROP1) + ERRADD1;
while(DV_ASSAY1 < 0) {
simeps();
DV_ASSAY1 = IPRED*(1 + ERRPROP1) + ERRADD1; // keep concentrations positive while using ERRADD
}
double DV_ASSAY2 = IPRED*(1 + ERRPROP2) + ERRADD2;
while(DV_ASSAY2 < 0) {
simeps();
DV_ASSAY2 = IPRED*(1 + ERRPROP2) + ERRADD2; // keep concentrations positive while using ERRADD
}
$CAPTURE
CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV_ASSAY1 DV_ASSAY2
' Compile and run modelmod_decock <- mcode("AmikacinPopPKModel",code)
# Run Model
res <- mod_decock %>%
data_set(ev(ID = 1:100, amt = 500)) %>%
obsonly()%>%
mrgsim(end = 24) %>%
as.data.frame() %>%
gather(key= "Assay", value = "DV", DV_ASSAY1:DV_ASSAY2) Plot concentration distribution to show effect of Assayggplot(data = filter(res)) +
geom_density(aes(y = DV, fill = Assay), alpha = 0.6, color = NA) +
facet_wrap("time", scales ="free_x", nrow=1, switch = "x") +
labs(x = "Time point (hrs.)", y = "Concentration",
title = "Distribution of Amikacin Concentration") +
scale_fill_manual(values = c( "grey40", "purple2")) +
theme(axis.text.x = element_blank(), axis.ticks.x = element_blank()) +
theme(strip.background =element_rect(fill="white", color = "white")) Of course there’s other ways to accomplish the same thing. If you have many different sets of RUV values of interest, the above method may not be the most convenient. In that case you can check out the |
Beta Was this translation helpful? Give feedback.
-
Hi Mohamed,
Thanks so much for your detailed instructions! Sorry that my question was
unclear. In fact, it was how to sample randomly values between runs with a
fixed set of ERRPROP and ERRADD (as they were characterized in a previous
study). Could you please give me another guide with this? 🙏
Thanks again.
Kind regards,
Thi.
…On Sat, 21 May 2022 at 05:42, Mohamed Ismail ***@***.***> wrote:
Hi @thinguyen94 <https://github.com/thinguyen94>. If I'm interpreting
your question, you'd like to simulate concentrations with varying levels of
RUV, while keeping IIV and covariate terms consistent. Since RUV is tacked
on at the end of the simulation, you could always create a separate output
for each assay, and use different variance terms for each one. This has the
advantage of only having to run a simulation once and guarantees the same
subjects are used between the two outputs (in terms of randomly sampled
ETAs and baseline characteristics).
See example below. I edited your model code to include two more sigma
terms and to capture two separate DV terms (one for each assay).
## Load Packages
library(mrgsolve)
library(ggplot2)
library(tidyr)
## Set plot theme
theme_set(theme_classic())
## Model Code
code<- '
$INIT // Initial Conditions for Compartments
CMT1 = 0, // Central Compartment
CMT2 = 0, // Peripheral Compartment
$SET // Set Differential Equation Solver Options
atol = 1e-8, rtol = 1e-8
maxsteps = 100000
$PARAM // Population parameters
POPCL = 0.0493, // Clearance
POPV1 = 0.833, // Central Volume
POPV2 = 0.833, // Peripheral Volume
POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
// Covariate effects
BW_CL = 1.340, // Exponent for Birthweight
PNA_CL = 0.213, // Slope of PNA
IBU_CL = 0.838, // Effect of Ibuprofen Co-administration
CW_V = 0.919 // Exponent for Current Weight
// Default Covariate Values
BW = 1750, // g
CW = 1760, // g
PNA = 2, // days
NSAID = 0 // no co-administration of ibuprofen
$OMEGA // Population Parameter Variability
@labels ETA1
0.0899
// Added variance terms for second assay
$SIGMA // Residual Unexplained Variability
@labels ERRPROP1 ERRADD1 ERRPROP2 ERRADD2
0.0614 0.267 0.02 0.1
$MAIN // Individual Parameter Values
double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) * exp(ETA1);
double V1 = (POPV1 * pow(CW/1760,CW_V));
double V2 = V1;
double Q = (POPQFRAC * CL);
if(NSAID == 1) {
CL = CL * IBU_CL;
Q = Q * IBU_CL;
}
$ODE // Differential Equations
double k10 = CL/V1;
double k12 = Q/V1;
double k21 = Q/V2;
dxdt_CMT1 = - k12*CMT1 + k21*CMT2 - k10*CMT1;
dxdt_CMT2 = k12*CMT1 - k21*CMT2;
$TABLE //Determines Values and Includes in Output
double IPRED = CMT1/V1;
// Capturing DV for both assays below
double DV_ASSAY1 = IPRED*(1 + ERRPROP1) + ERRADD1;
while(DV_ASSAY1 < 0) {
simeps();
DV_ASSAY1 = IPRED*(1 + ERRPROP1) + ERRADD1; // keep concentrations positive while using ERRADD
}
double DV_ASSAY2 = IPRED*(1 + ERRPROP2) + ERRADD2;
while(DV_ASSAY2 < 0) {
simeps();
DV_ASSAY2 = IPRED*(1 + ERRPROP2) + ERRADD2; // keep concentrations positive while using ERRADD
}
$CAPTURE
CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV_ASSAY1 DV_ASSAY2
'
Compile and run model
mod_decock <- mcode("AmikacinPopPKModel",code)
## Building AmikacinPopPKModel ... done.
# Run Model
res <- mod_decock %>%
data_set(ev(ID = 1:100, amt = 500)) %>%
obsonly()%>%
mrgsim(end = 24) %>%
as.data.frame() %>%
gather(key= "Assay", value = "DV", DV_ASSAY1:DV_ASSAY2)
Plot concentration distribution to show effect of Assay
ggplot(data = filter(res)) +
geom_density(aes(y = DV, fill = Assay), alpha = 0.6, color = NA) +
facet_wrap("time", scales ="free_x", nrow=1, switch = "x") +
labs(x = "Time point (hrs.)", y = "Concentration",
title = "Distribution of Amikacin Concentration") +
scale_fill_manual(values = c( "grey40", "purple2")) +
theme(axis.text.x = element_blank(), axis.ticks.x = element_blank()) +
theme(strip.background =element_rect(fill="white", color = "white"))
[image: unnamed-chunk-3-1]
<https://user-images.githubusercontent.com/27981719/169620898-a9a3b9f7-f19c-4f18-9210-67f9142baaf7.png>
Of course there’s other ways to accomplish the same thing. If you have
many different sets of RUV values of interest, the above method may not be
the most convenient. In that case you can check out the smat() function
in mrgsolve which lets you update the sigma terms within your rscript, and
set.seed() to make sure the same IIV terms are sampled between runs.
—
Reply to this email directly, view it on GitHub
<#983 (comment)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AVWTFRZVYUDPIBOE72SYXJTVLAINFANCNFSM5WLBGZVQ>
.
You are receiving this because you were mentioned.Message ID:
***@***.***
.com>
|
Beta Was this translation helpful? Give feedback.
-
Hi @thinguyen94. I am still not sure I completely understand what you’re trying to accomplish, but I’ll take another stab at it. My interpretation is that you would like to run multiple simulations of the same design, and for each simulation keep the individuals exactly the same (covariates and random ETA terms), but randomly sample RUV between runs. Essentially you would like to keep the seed for randomly sampled ETAs the same between runs, but use a different seed for EPS for each run. See an example of this below. The trick that I’m using here is to resample the sigma terms n times depending on the simulation number with simeps(), guarenteeing the random samples will be unique for each iteration (you can see this at the start of the TABLE record). This feels incredibly hacky, and @kylebaron might know a better way, but it works for now 😏. library(mrgsolve)
library(ggplot2)
library(purrr)
library(dplyr)
theme_set(theme_bw())
code <- '
$INIT // Initial Conditions for Compartments
CMT1 = 0, // Central Compartment
CMT2 = 0, // Peripheral Compartment
$SET // Set Differential Equation Solver Options
atol = 1e-8, rtol = 1e-8
maxsteps = 100000
$PARAM // Population parameters
POPCL = 0.0493, // Clearance
POPV1 = 0.833, // Central Volume
POPV2 = 0.833, // Peripheral Volume
POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
// Covariate effects
BW_CL = 1.340, // Exponent for Birthweight
PNA_CL = 0.213, // Slope of PNA
IBU_CL = 0.838, // Effect of Ibuprofen Co-administration
CW_V = 0.919 // Exponent for Current Weight
// Default Covariate Values
BW = 1750, // g
CW = 1760, // g
PNA = 2, // days
NSAID = 0 // no co-administration of ibuprofen
// Simulation number which we will use to make epsilons unique
SIM_I = NA_real_
$OMEGA // Population Parameter Variability
@labels ETA1
0.0899
$SIGMA // Residual Unexplained Variability
@labels ERRPROP ERRADD
0.0614 0.267
$MAIN // Individual Parameter Values
double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) * exp(ETA1);
double V1 = (POPV1 * pow(CW/1760,CW_V));
double V2 = V1;
double Q = (POPQFRAC * CL);
if(NSAID == 1) {
CL = CL * IBU_CL;
Q = Q * IBU_CL;
}
$ODE // Differential Equations
double k10 = CL/V1;
double k12 = Q/V1;
double k21 = Q/V2;
dxdt_CMT1 = - k12*CMT1 + k21*CMT2 - k10*CMT1;
dxdt_CMT2 = k12*CMT1 - k21*CMT2;
$TABLE //Determines Values and Includes in Output
// resample epsilons once for first simulation, twice for second, etc.
double simn = 0;
while(simn < SIM_I){
simeps();
simn++;
}
double IPRED = CMT1/V1;
double DV = IPRED*(1 + ERRPROP) + ERRADD;
while(DV < 0) {
simeps();
DV = IPRED*(1 + ERRPROP) + ERRADD; // keep concentrations positive while using ERRADD
}
$CAPTURE
CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV SIM_I
' Compile and define simulation functionmod_decock <- mcode("AmikacinPopPKModel",code)
#simi = simulation number
simulate <- function(simi){
# setting seed so each simulation iteration samples same ETAS
set.seed(1)
dosing <- ev(ID = 1:10, amt = 500, SIM_I = simi)
mrgsim_df(mod_decock, dosing)
} Run simulations## Run ten simulations and bind rows together
res <- map_df(seq(10), simulate) Plot IPRED (black line) to show consistency between runs, and DVs (colored points) to show changing between runsggplot(res, aes(x = time, color = as.factor(SIM_I))) +
geom_point(aes(y = DV)) +
geom_line(aes(y = IPRED), color = "black", size = 1) +
facet_wrap("ID") +
labs(color = "Simulation Number") Another way you could accomplish the same thing is to use mrgsolve just to get the IPRED (run a single simulation), then tack on RUV directly in R using rnorm, where |
Beta Was this translation helpful? Give feedback.
-
Hi Mohamed,
This is exactly what I want. Thanks so so much for your help - much more
than I expected 😊
Kind regards,
Thi.
On Sun, 22 May 2022 at 04:18, Mohamed Ismail ***@***.***> wrote:
Hi @thinguyen94 <https://github.com/thinguyen94>. I am still not sure I
completely understand what you’re trying to accomplish, but I’ll take
another stab at it. My interpretation is that you would like to run
multiple simulations of the same design, and for each simulation keep the
individuals exactly the same (covariates and random ETA terms), but
randomly sample RUV between runs. Essentially you would like to keep the
seed for randomly sampled ETAs the same between runs, but use a different
seed for EPS for each run.
See an example of this below. The trick that I’m using here is to resample
the sigma terms n times depending on the simulation number with simeps(),
guarenteeing the random samples will be unique for each iteration (you can
see this at the start of the TABLE record). This feels incredibly hacky,
and @kylebaron <https://github.com/kylebaron> might know a better way,
but it works for now 😏.
library(mrgsolve)
library(ggplot2)
library(purrr)
library(dplyr)
theme_set(theme_bw())
code <- '
$INIT // Initial Conditions for Compartments
CMT1 = 0, // Central Compartment
CMT2 = 0, // Peripheral Compartment
$SET // Set Differential Equation Solver Options
atol = 1e-8, rtol = 1e-8
maxsteps = 100000
$PARAM // Population parameters
POPCL = 0.0493, // Clearance
POPV1 = 0.833, // Central Volume
POPV2 = 0.833, // Peripheral Volume
POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
// Covariate effects
BW_CL = 1.340, // Exponent for Birthweight
PNA_CL = 0.213, // Slope of PNA
IBU_CL = 0.838, // Effect of Ibuprofen Co-administration
CW_V = 0.919 // Exponent for Current Weight
// Default Covariate Values
BW = 1750, // g
CW = 1760, // g
PNA = 2, // days
NSAID = 0 // no co-administration of ibuprofen
// Simulation number which we will use to make epsilons unique
SIM_I = NA_real_
$OMEGA // Population Parameter Variability
@labels ETA1
0.0899
$SIGMA // Residual Unexplained Variability
@labels ERRPROP ERRADD
0.0614 0.267
$MAIN // Individual Parameter Values
double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) * exp(ETA1);
double V1 = (POPV1 * pow(CW/1760,CW_V));
double V2 = V1;
double Q = (POPQFRAC * CL);
if(NSAID == 1) {
CL = CL * IBU_CL;
Q = Q * IBU_CL;
}
$ODE // Differential Equations
double k10 = CL/V1;
double k12 = Q/V1;
double k21 = Q/V2;
dxdt_CMT1 = - k12*CMT1 + k21*CMT2 - k10*CMT1;
dxdt_CMT2 = k12*CMT1 - k21*CMT2;
$TABLE //Determines Values and Includes in Output
// resample epsilons once for first simulation, twice for second, etc.
double simn = 0;
while(simn < SIM_I){
simeps();
simn++;
}
double IPRED = CMT1/V1;
double DV = IPRED*(1 + ERRPROP) + ERRADD;
while(DV < 0) {
simeps();
DV = IPRED*(1 + ERRPROP) + ERRADD; // keep concentrations positive while using ERRADD
}
$CAPTURE
CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV SIM_I
'
Compile and define simulation function
mod_decock <- mcode("AmikacinPopPKModel",code)
## Building AmikacinPopPKModel ... done.
#simi = simulation number
simulate <- function(simi){
# setting seed so each simulation iteration samples same ETAS
set.seed(1)
dosing <- ev(ID = 1:10, amt = 500, SIM_I = simi)
mrgsim_df(mod_decock, dosing)
}
Run simulations
## Run ten simulations and bind rows together
res <- map_df(seq(10), simulate)
Plot IPRED (black line) to show consistency between runs, and DVs (colored
points) to show changing between runs
ggplot(res, aes(x = time, color = as.factor(SIM_I))) +
geom_point(aes(y = DV)) +
geom_line(aes(y = IPRED), color = "black", size = 1) +
facet_wrap("ID") +
labs(color = "Simulation Number")
[image: unnamed-chunk-4-1]
<https://user-images.githubusercontent.com/27981719/169669193-47242910-65f0-42a3-80d7-0ad47a265dd2.png>
Another way you could accomplish the same thing is to use mrgsolve just to
get the IPRED (run a single simulation), then tack on RUV directly in R
using rnorm, where DV = IPRED*(1 + rnorm(1, mean = 0, sd =
sqrt(ERRPROP))) + rnorm(1, mean = 0, sd = sqrt(ERRADD).
—
Reply to this email directly, view it on GitHub
<#983 (comment)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AVWTFR2HDJ7FV3USSH376F3VLFHKBANCNFSM5WLBGZVQ>
.
You are receiving this because you were mentioned.Message ID:
***@***.***
.com>
--
Kind regards,
Thi.
|
Beta Was this translation helpful? Give feedback.
-
Hi Mohamed,
Just one more question regarding your second option. I edited my codes as
below (in red), but DVs came out with NaNs. I cannot figure out where it
falls through. Could you please give me another hand?🙏🙏 🙏THANK YOU
code <- '
$INIT // Initial Conditions for Compartments
CMT1 = 0, // Central Compartment
CMT2 = 0, // Peripheral Compartment
$SET // Set Differential Equation Solver Options
atol = 1e-8, rtol = 1e-8
maxsteps = 100000
$PARAM // Population parameters
POPCL = 0.0493, // Clearance
POPV1 = 0.833, // Central Volume
POPV2 = 0.833, // Peripheral Volume
POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
// Covariate effects
BW_CL = 1.340, // Exponent for Birthweight
PNA_CL = 0.213, // Slope of PNA
IBU_CL = 0.838, // Effect of Ibuprofen Co-administration
CW_V = 0.919 // Exponent for Current Weight
// Default Covariate Values
BW = 1750, // g
CW = 1760, // g
PNA = 2, // days
NSAID = 0 // no co-administration of ibuprofen
$OMEGA // Population Parameter Variability
@labels ETA1
0.0899
$SIGMA // Residual Unexplained Variability
@labels ERRPROP ERRADD
0.0614 0.267
$MAIN // Individual Parameter Values
double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) *
exp(ETA1);
double V1 = (POPV1 * pow(CW/1760,CW_V));
double V2 = V1;
double Q = (POPQFRAC * CL);
if(NSAID == 1) {
CL = CL * IBU_CL;
Q = Q * IBU_CL;
}
$ODE // Differential Equations
double k10 = CL/V1;
double k12 = Q/V1;
double k21 = Q/V2;
dxdt_CMT1 = - k12*CMT1 + k21*CMT2 - k10*CMT1;
dxdt_CMT2 = k12*CMT1 - k21*CMT2;
$PLUGIN Rcpp
$TABLE //Determines Values and Includes in Output
double IPRED = CMT1/V1;
double DV = IPRED*(1 + R::rnorm(0, sqrt(ERRPROP))) + R::rnorm(0,
sqrt(ERRADD));
while(DV < 0) {
simeps();
DV = IPRED*(1 + R::rnorm(0, sqrt(ERRPROP))) + R::rnorm(0, sqrt(ERRADD)); //
keep concentrations positive while using ERRADD
}
$CAPTURE
CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV
'
mod_decock <- mcode("AmikacinPopPKModel",code)
## Building AmikacinPopPKModel ... done.
# Run Model
res <- mod_decock %>%
data_set(ev(ID = 1:100, amt = 500)) %>%
obsonly()%>%
mrgsim(end = 24)%>%
as.data.frame()
head(res)
ID time CMT1 CMT2 CL V1 V2 Q ETA1 BW
CW PNA NSAID IPRED DV
1 1 0 0.0000 0.00000 0.03069394 0.833 0.833 0.01273798 -0.6669555
1750 1760 2 0 0.0000 NaN
2 1 1 474.6546 7.39307 0.03069394 0.833 0.833 0.01273798 -0.6669555
1750 1760 2 0 569.8134 582.6754
3 1 2 450.7032 14.30003 0.03069394 0.833 0.833 0.01273798 -0.6669555
1750 1760 2 0 541.0603 NaN
4 1 3 428.0681 20.74881 0.03069394 0.833 0.833 0.01273798 -0.6669555
1750 1760 2 0 513.8873 NaN
Kind regards,
Thi.
…On Sun, 22 May 2022 at 09:24, Thi Nguyen ***@***.***> wrote:
Hi Mohamed,
This is exactly what I want. Thanks so so much for your help - much more
than I expected 😊
Kind regards,
Thi.
On Sun, 22 May 2022 at 04:18, Mohamed Ismail ***@***.***>
wrote:
> Hi @thinguyen94 <https://github.com/thinguyen94>. I am still not sure I
> completely understand what you’re trying to accomplish, but I’ll take
> another stab at it. My interpretation is that you would like to run
> multiple simulations of the same design, and for each simulation keep the
> individuals exactly the same (covariates and random ETA terms), but
> randomly sample RUV between runs. Essentially you would like to keep the
> seed for randomly sampled ETAs the same between runs, but use a different
> seed for EPS for each run.
>
> See an example of this below. The trick that I’m using here is to
> resample the sigma terms n times depending on the simulation number with
> simeps(), guarenteeing the random samples will be unique for each iteration
> (you can see this at the start of the TABLE record). This feels incredibly
> hacky, and @kylebaron <https://github.com/kylebaron> might know a better
> way, but it works for now 😏.
>
> library(mrgsolve)
>
> library(ggplot2)
>
> library(purrr)
>
> library(dplyr)
>
> theme_set(theme_bw())
>
>
>
>
> code <- '
> $INIT // Initial Conditions for Compartments
> CMT1 = 0, // Central Compartment
> CMT2 = 0, // Peripheral Compartment
>
> $SET // Set Differential Equation Solver Options
> atol = 1e-8, rtol = 1e-8
> maxsteps = 100000
>
> $PARAM // Population parameters
> POPCL = 0.0493, // Clearance
> POPV1 = 0.833, // Central Volume
> POPV2 = 0.833, // Peripheral Volume
> POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
>
> // Covariate effects
> BW_CL = 1.340, // Exponent for Birthweight
> PNA_CL = 0.213, // Slope of PNA
> IBU_CL = 0.838, // Effect of Ibuprofen Co-administration
> CW_V = 0.919 // Exponent for Current Weight
>
> // Default Covariate Values
> BW = 1750, // g
> CW = 1760, // g
> PNA = 2, // days
> NSAID = 0 // no co-administration of ibuprofen
>
> // Simulation number which we will use to make epsilons unique
> SIM_I = NA_real_
>
>
> $OMEGA // Population Parameter Variability
> @labels ETA1
> 0.0899
>
> $SIGMA // Residual Unexplained Variability
> @labels ERRPROP ERRADD
> 0.0614 0.267
>
> $MAIN // Individual Parameter Values
> double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) * exp(ETA1);
> double V1 = (POPV1 * pow(CW/1760,CW_V));
> double V2 = V1;
> double Q = (POPQFRAC * CL);
> if(NSAID == 1) {
> CL = CL * IBU_CL;
> Q = Q * IBU_CL;
> }
>
> $ODE // Differential Equations
> double k10 = CL/V1;
> double k12 = Q/V1;
> double k21 = Q/V2;
> dxdt_CMT1 = - k12*CMT1 + k21*CMT2 - k10*CMT1;
> dxdt_CMT2 = k12*CMT1 - k21*CMT2;
>
> $TABLE //Determines Values and Includes in Output
>
> // resample epsilons once for first simulation, twice for second, etc.
> double simn = 0;
> while(simn < SIM_I){
> simeps();
> simn++;
> }
>
>
> double IPRED = CMT1/V1;
> double DV = IPRED*(1 + ERRPROP) + ERRADD;
> while(DV < 0) {
> simeps();
> DV = IPRED*(1 + ERRPROP) + ERRADD; // keep concentrations positive while using ERRADD
> }
> $CAPTURE
> CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV SIM_I
> '
>
> Compile and define simulation function
>
> mod_decock <- mcode("AmikacinPopPKModel",code)
>
> ## Building AmikacinPopPKModel ... done.
>
>
> #simi = simulation number
> simulate <- function(simi){
>
> # setting seed so each simulation iteration samples same ETAS
>
> set.seed(1)
>
>
>
> dosing <- ev(ID = 1:10, amt = 500, SIM_I = simi)
>
> mrgsim_df(mod_decock, dosing)
>
> }
>
> Run simulations
>
> ## Run ten simulations and bind rows together
> res <- map_df(seq(10), simulate)
>
> Plot IPRED (black line) to show consistency between runs, and DVs
> (colored points) to show changing between runs
>
> ggplot(res, aes(x = time, color = as.factor(SIM_I))) +
>
> geom_point(aes(y = DV)) +
>
> geom_line(aes(y = IPRED), color = "black", size = 1) +
>
> facet_wrap("ID") +
>
> labs(color = "Simulation Number")
>
> [image: unnamed-chunk-4-1]
> <https://user-images.githubusercontent.com/27981719/169669193-47242910-65f0-42a3-80d7-0ad47a265dd2.png>
>
> Another way you could accomplish the same thing is to use mrgsolve just
> to get the IPRED (run a single simulation), then tack on RUV directly in R
> using rnorm, where DV = IPRED*(1 + rnorm(1, mean = 0, sd =
> sqrt(ERRPROP))) + rnorm(1, mean = 0, sd = sqrt(ERRADD).
>
> —
> Reply to this email directly, view it on GitHub
> <#983 (comment)>,
> or unsubscribe
> <https://github.com/notifications/unsubscribe-auth/AVWTFR2HDJ7FV3USSH376F3VLFHKBANCNFSM5WLBGZVQ>
> .
> You are receiving this because you were mentioned.Message ID:
> <metrumresearchgroup/mrgsolve/repo-discussions/983/comments/2797230@
> github.com>
>
--
Kind regards,
Thi.
|
Beta Was this translation helpful? Give feedback.
-
Hi @thinguyen94. For that option, you should add RUV in R after the simulation is done, not within the mrgsolve code itself. You can run the mrgsolve simulation per your original code, then use the IPREDs and rnorm to get DV values using R. ## Set up function to add RUV to the IPRED values
addRUV <- function(iteration, sim_results){
ERRPROP <- 0.0614
ERRADD <- 0.267
sim_results %>%
mutate(DV = IPRED*(1 + rnorm(n(), mean = 0, sd = sqrt(ERRPROP))) +
rnorm(n(), mean = 0, sd = sqrt(ERRADD)),
SIM_I = iteration)
}
## Run simulation
dosing <- ev(ID = 1:10, amt = 500)
## Now we're only running the simulation once with mrgsolve
sim_results <- mrgsim_df(mod_decock, dosing)
## Adding RUV in R post simulation, 10 iterations
res <- map_df(seq(10), ~addRUV(.x, sim_results))
## Plotting
ggplot(res, aes(x = time, color = as.factor(SIM_I))) +
geom_point(aes(y = DV)) +
geom_line(aes(y = IPRED), color = "black", size = 1) +
facet_wrap("ID") +
labs(color = "Simulation Number")
|
Beta Was this translation helpful? Give feedback.
-
Hi Mohamed,
Fantastic! Thank you 😉
Wishing you a nice week ahead.
Kind regards,
Thi.
…On Sun, 22 May 2022 at 22:19, Mohamed Ismail ***@***.***> wrote:
Hi @thinguyen94 <https://github.com/thinguyen94>. For that option, you
should add RUV in R after the simulation is done, not within the mrgsolve
code itself. You can run the mrgsolve simulation per your original code,
then use the IPREDs and rnorm to get DV values using R.
## Set up function to add RUV to the IPRED valuesaddRUV <- function(iteration, sim_results){
ERRPROP <- 0.0614
ERRADD <- 0.267
sim_results %>%
mutate(DV = IPRED*(1 + rnorm(n(), mean = 0, sd = sqrt(ERRPROP))) +
rnorm(n(), mean = 0, sd = sqrt(ERRADD)),
SIM_I = iteration)
}
## Run simulationdosing <- ev(ID = 1:10, amt = 500)
## Now we're only running the simulation once with mrgsolvesim_results <- mrgsim_df(mod_decock, dosing)
## Adding RUV in R post simulation, 10 iterationsres <- map_df(seq(10), ~addRUV(.x, sim_results))
## Plotting
ggplot(res, aes(x = time, color = as.factor(SIM_I))) +
geom_point(aes(y = DV)) +
geom_line(aes(y = IPRED), color = "black", size = 1) +
facet_wrap("ID") +
labs(color = "Simulation Number")
[image: image]
<https://user-images.githubusercontent.com/27981719/169702429-3e3cb4b9-39b1-4b44-bfd2-5c3fbaaf0c91.png>
—
Reply to this email directly, view it on GitHub
<#983 (comment)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AVWTFR7D4IOYRUU6354ACR3VLJF7LANCNFSM5WLBGZVQ>
.
You are receiving this because you were mentioned.Message ID:
***@***.***
.com>
|
Beta Was this translation helpful? Give feedback.
-
Hi @kylebaron,
Hope that everything is going well at your end.
I am working on a project aimed to assess the impact of bioanalytical variation on dose decisions in neonates using simulation. I am trying to run a sensitivity analysis of RUV representing intra-individual variability on dose decisions. So I would like to randomly generate RUV while maintaining IIV and all the remaining variables (clinical characteristics). But I couldn't find a way yet by reading through your tutorials.
What I have done so far is randomly generate both RUV and IIV while maintaining the rest simply by rerunning mrgsim() function. I am very much appreciative if you can guide me on how to do it.
Thanks so much and I look forward to your response.
Kind regards,
Thi
Below is my PK model code:
code <- '
$INIT // Initial Conditions for Compartments
CMT1 = 0, // Central Compartment
CMT2 = 0, // Peripheral Compartment
$SET // Set Differential Equation Solver Options
atol = 1e-8, rtol = 1e-8
maxsteps = 100000
$PARAM // Population parameters
POPCL = 0.0493, // Clearance
POPV1 = 0.833, // Central Volume
POPV2 = 0.833, // Peripheral Volume
POPQFRAC = 0.415, // Fraction of Inter-compartmental Clearance
$OMEGA // Population Parameter Variability
@labels ETA1
0.0899
$SIGMA // Residual Unexplained Variability
@labels ERRPROP ERRADD
0.0614 0.267
$MAIN // Individual Parameter Values
double CL = (POPCL * pow(BW/1750, BW_CL) * (1 + PNA_CL*(PNA/2))) * exp(ETA1);
double V1 = (POPV1 * pow(CW/1760,CW_V));
double V2 = V1;
double Q = (POPQFRAC * CL);
if(NSAID == 1) {
CL = CL * IBU_CL;
Q = Q * IBU_CL;
}
$ODE // Differential Equations
double k10 = CL/V1;
double k12 = Q/V1;
double k21 = Q/V2;
dxdt_CMT1 = - k12CMT1 + k21CMT2 - k10CMT1;
dxdt_CMT2 = k12CMT1 - k21*CMT2;
$TABLE //Determines Values and Includes in Output
double IPRED = CMT1/V1;
double DV = IPRED*(1 + ERRPROP) + ERRADD;
while(DV < 0) {
simeps();
DV = IPRED*(1 + ERRPROP) + ERRADD; // keep concentrations positive while using ERRADD
}
$CAPTURE
CL V1 V2 Q ETA1 BW CW PNA NSAID IPRED DV
'
mod_decock <- mcode("AmikacinPopPKModel",code)
Beta Was this translation helpful? Give feedback.
All reactions