-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
194 lines (140 loc) · 8.01 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
---
output: github_document
---
# campsismisc <img src='man/figures/r_package_campsismisc.png' align="right" width="120"/>
<!-- badges: start -->
[](https://github.com/Calvagone/campsismisc/actions)
[](https://app.codecov.io/gh/Calvagone/campsismisc)
<!-- badges: end -->
## Forest plots
Assume the following estimated 2-compartment model:
```{r, message=FALSE, warning=FALSE}
library(campsismisc)
library(progressr)
model <- read.campsis("tests/testthat/campsis_models/metaboliser_effect_on_cl/")
model
```
This model contains:
- a metabolism effect on the clearance: `THETA_METAB_CL`
- allometric scaling on central and peripheral volumes with a fixed exponent of 1
- a variance-covariance matrix
### Configure your simulation settings
```{r, message=F}
# Configure hardware
settings <- Settings(Hardware(cpu=8, replicate_parallel=TRUE))
dest <- "mrgsolve" # Much faster than RxODE
```
### Effect of metabolism effect and weight on model parameters
Let's show how the metabolism effect and the weight influence the clearance or the volume. This can be achieved as follows:
```{r, message=F}
outputs <- OATOutputs() %>%
add(ModelParameterOutput("CL")) %>%
add(ModelParameterOutput("VC"))
fp <- ForestPlot(model=model, outputs=outputs, replicates=100, dest=dest, settings=settings) %>%
add(CategoricalLabeledCovariate(name="METAB", default_value=0, label="Metaboliser", categories=c(Slow=0, Fast=1))) %>%
add(LabeledCovariate(name="WT", default_value=70, label="Weight", unit="kg")) %>%
add(ForestPlotItem(Covariate("METAB", 0))) %>%
add(ForestPlotItem(Covariate("METAB", 1))) %>%
add(ForestPlotItem(Covariate("WT", 60))) %>%
add(ForestPlotItem(Covariate("WT", 80)))
fp <- with_progress({fp %>% prepare()})
```
The forest plot for clearance can be otained as follows (use `index=1`):
```{r fp_metabolism_effect_on_cl, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=1) +
ggplot2::scale_y_continuous(breaks=c(0.7,0.8,1,1.25,1.4), limits=c(0.5, 1.5))
```
Similarly, the forest plot for volume can be obtained as follows (use `index=2`):
```{r fp_metabolism_effect_on_vc, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=2) +
ggplot2::scale_y_continuous(breaks=c(0.7,0.8,1,1.25,1.4), limits=c(0.5, 1.5))
```
### Effect of metabolism and weight on PK metrics
Let's show how the metabolism effect and the weight influence PK metrics. For that, we need to create first a dataset. Assume we are interested to compute AUC/Cmax on day 1 and day 7. We give the drug every day and we observe on day 1 and day 7.
```{r}
dataset <- Dataset(1) %>%
add(Infusion(time=0, amount=1000, compartment=1, ii=24, addl=6)) %>%
add(Observations(times=c(0:24, 144:168)))
```
Instead of creating a `ModelParameterOutput` object, we create a `NcaMetricOutput` for each metric we want to compute (by referring to the proper `campsisnca` metric). Each output is added to the list of outputs as follows:
```{r}
outputs <- OATOutputs() %>%
add(NcaMetricOutput(AUC(variable="CONC"), filter=~campsisnca::timerange(.x, min=0, max=24))) %>%
add(NcaMetricOutput(AUC(variable="CONC"), filter=~campsisnca::timerange(.x, min=144, max=168), suffix="Day 7")) %>%
add(NcaMetricOutput(Cmax(variable="CONC"), filter=~campsisnca::timerange(.x, min=0, max=24))) %>%
add(NcaMetricOutput(Cmax(variable="CONC"), filter=~campsisnca::timerange(.x, min=144, max=168), suffix="Day 7"))
```
A forest plot object can be instantiated and run as follows:
```{r, message=F}
fp <- ForestPlot(model=model, dataset=dataset, outputs=outputs,
replicates=100, dest=dest, settings=settings) %>%
add(CategoricalLabeledCovariate(name="METAB", default_value=0, label="Metaboliser", categories=c(Slow=0, Fast=1))) %>%
add(LabeledCovariate(name="WT", default_value=70, label="Weight", unit="kg")) %>%
add(ForestPlotItem(Covariate("METAB", 0))) %>%
add(ForestPlotItem(Covariate("METAB", 1))) %>%
add(ForestPlotItem(Covariate("WT", 60))) %>%
add(ForestPlotItem(Covariate("WT", 80)))
fp <- with_progress({fp %>% prepare()})
```
We can look at AUC day 1 by doing:
```{r fp_metabolism_effect_on_auc_d1, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=1) +
ggplot2::scale_y_continuous(breaks=c(0.7,0.8,1,1.25,1.4), limits=c(0.5, 1.5))
```
If we are interested to see the absolute change in AUC, argument relative can be set to FALSE as follows.
```{r fp_metabolism_effect_on_auc_d1_absolute, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=1, relative=FALSE)
```
To look at AUC on day 7, we proceed as follows:
```{r fp_metabolism_effect_on_auc_d7, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=2) +
ggplot2::scale_y_continuous(breaks=c(0.7,0.8,1,1.25,1.4), limits=c(0.5, 1.5))
```
Again, we can look at the absolute values.
```{r fp_metabolism_effect_on_auc_d7_absolute, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=2, relative=FALSE)
```
Let's produce one more plot with Cmax on day 1.
```{r fp_metabolism_effect_on_cmax_d1, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=3) +
ggplot2::scale_y_continuous(breaks=c(0.7,0.8,1,1.25,1.4), limits=c(0.5, 1.5))
```
And a final one with Cmax on day 7.
```{r fp_metabolism_effect_on_cmax_d7, fig.align='center', fig.height=4, fig.width=8, message=F}
fp %>% getForestPlot(index=4) +
ggplot2::scale_y_continuous(breaks=c(0.7,0.8,1,1.25,1.4), limits=c(0.5, 1.5))
```
Finally, you need to know that you can combine model parameter ouputs together with NCA metric outputs, if you want! This will make your simulations even faster. In that specific case, Campsis will check that your model parameters do not vary over time (since several observations are provided in the dataset).
## OAT-method based sensitivity analysis
### Effect of model parameters on PK metrics
Assume the same 2-compartment model. Let's see how the model parameters influence AUC or Cmax if they are multiplied by two. This can be achieved with a one-at-a-time (OAT)-method-based sensitivity analysis. Most of the time, these analyses are done without parameter uncertainty (replicates=1, by default) and represented with tornado plots, as shown below. However if parameter uncertainty is used, you can still call `getForestPlot()` on this new type of object.
```{r sensibility_analysis_auc_example, fig.align='center', fig.height=3, fig.width=8, message=F}
outputs <- OATOutputs() %>%
add(NcaMetricOutput(AUC(variable="CONC"))) %>%
add(NcaMetricOutput(Cmax(variable="CONC")))
dataset <- Dataset(1) %>%
add(Infusion(time=0, amount=1000, compartment=1)) %>%
add(Observations(times=seq(0, 24, by=0.1))) %>%
add(Covariate("METAB", 0)) %>%
add(Covariate("WT", 70))
object <- SensitivityAnalysis(model=model, dataset=dataset, outputs=outputs) %>%
add(LabeledParameters(c(DUR="Duration", VC="Central volume", VP="Peripheral volume", Q="Inter-compartmental\nclearance", CL="Clearance"))) %>%
add(SensitivityAnalysisItem(Change("DUR", up=2, down=2, log=TRUE))) %>%
add(SensitivityAnalysisItem(Change("VC", up=2, down=2, log=TRUE))) %>%
add(SensitivityAnalysisItem(Change("VP", up=2, down=2, log=TRUE))) %>%
add(SensitivityAnalysisItem(Change("Q", up=2, down=2, log=TRUE))) %>%
add(SensitivityAnalysisItem(Change("CL", up=2, down=2, log=TRUE)))
object <- object %>% prepare()
object %>% getTornadoPlot(index=1)
```
Absolute values may also be shown:
```{r sensibility_analysis_auc_example_absolute, fig.align='center', fig.height=3, fig.width=8, message=F}
object %>% getTornadoPlot(index=1, relative=FALSE)
```
Cmax can also be shown using `index=2`.
```{r sensibility_analysis_cmax_example, fig.align='center', fig.height=3, fig.width=8, message=F}
object %>% getTornadoPlot(index=2)
```
```{r, echo=FALSE, message=F, results='hide'}
setupPlanSequential()
```