tdmore-dev / tdmore

In silico evaluation of precision dosing
https://tdmore-dev.github.io/tdmore/dev
GNU Affero General Public License v3.0
3 stars 1 forks source link

Standardize `predict` and `stat_predict` functions #64

Closed rfaelens closed 5 years ago

rfaelens commented 5 years ago

We should make sure that the tdmore prediction output is standard. This is not easy, because we have both parameter uncertainty as well as residual error.

Looking at the other modeling libraries, there are some that use the se.fit argument to add a SE column. This does not work here, as our problem is non-linear.

There are some that define the interval to show: interval = c("none", "confidence", "prediction"). This does not really match with what we are doing... Furthermore, we would like to combine this!

So instead, we have tried to define exhaustively all the prediction scenario's.

rfaelens commented 5 years ago

Key thought

We use predict for monte carlo simulations, and summarizing many samples of monte carlo simulations.

We use model.frame to add on residual error.

How should predict work?

predict(tdmore, regimen, newdata=seq(0, 24), samples=0, level=NA) ID EKA TIME CONC PD 0 0 0 10 3 0 0 1 14 2

predict(tdmore, regimen, newdata=seq(0, 24), samples=500, level=NA) ID EKA TIME CONC PD 1 0.32 0 10 3 1 0.32 1 14 2 1 ... 2 -0.28 0 8 1 2 -0.28 1 8.3 1.4 2 ...

predict(tdmore, regimen, newdata=seq(0, 24), samples=500, level=0.95) TIME CONC CONC.median CONC.lower CONC.upper PD PD.median ... 0 10 3 1 14 2 ...

predict(tdmore, regimen, newdata=seq(0, 24), samples=500, level=0.95) %>% meltPredictions() TIME variable value value.median value.lower value.upper 0 CONC 10 10.2 8.2 13.2 1 CONC ...

How should model.frame work?

model.frame(tdmorefit) TIME CONC 23 80

model.frame(tdmorefit, level=0.80) TIME CONC CONC.lower CONC.upper 23 80 79 81

Idea: let 'predict' return a 'tdmore.prediction' that you can directly pass to model.frame without specifying tdmore again

predict(tdmore, regimen, newdata=seq(0, 24)) %>% model.frame(tdmore, level=0.80, newdata=., id.vars="TIME") TIME CONC PD CONC.upper CONC.lower 0 10 3
1 14 2 ...

predict(tdmore, regimen, newdata=seq(0, 24), samples=500) %>% model.frame(tdmore, level=0.80, newdata=., id.vars=c("ID","TIME","EKA")) ID EKA TIME CONC PD CONC.upper CONC.lower 1 0.32 0 10 3
1 0.32 1 14 2 ... 1 -0.28 0 10 3
1 -0.28 1 14 2 ...

predict(tdmore, regimen, newdata=seq(0, 24), samples=500, level=0.95) %>% moel.frame(tdmore, level=0.80, newdata=., id.vars=c("TIME")) TIME CONC.lower CONC.median CONC.upper CONC.median.lower CONC.median.upper CONC.median.median CONC.upper.lower CONC.upper.median CONC.upper.upper 0 10 3 1 14 2 ...

predict(tdmore, regimen, newdata=seq(0, 24), samples=500, level=0.95) %>% melt() %>% model.frame(tdmore, level=0.80, newdata=., id.vars=c("TIME")) TIME variable value.median value.lower value.upper 0 CONC 10 9.9 10.1 0 CONC.lower 8.2 8.0 8.4 0 CONC.median 10.2 10.1 10.3 0 CONC.upper 12.3 12.1 12.5 ...

rfaelens commented 5 years ago

So how should stat_ipred work

stat_ipred( geom=NULL, tdmore, regimen, newdata, mapping, samples=0, ci=NA, re=NA )

method samples ci re usage
ipred 0 NA NA the typical value, as CONC
ipred 0 0.3 NA does not make sense; error
ipred 0 NA 0.3 add RE to typical value, as CONC.upper / CONC.lower
ipred 10 NA NA spaghetti plot, with extra ID column
ipred 10 0.3 NA ribbon, with CONC.upper / CONC.lower
ipred 10 NA 0.3 spaghetti plot, with CONC.upper, CONC.lower
ipred 10 0.3 0.3 ##huge amount of columns, with CONC.ci.upper.re.upper, CONC.ci.upper.re.lower, etc.
## DECISION: we fail to see how this would be routinely used
## We decided to consider this as an ERROR
rfaelens commented 5 years ago

Using stat_ipredmelt generates the same, but with a different variable for every model output.

We should still think about the name of this function. I like stat_predict a lot, as it is a direct analogue of predict. This also allows stat_fitted and stat_model.frame. However, it is a bit difficult to find a proper name for the a_prioriprediction in that case...

rfaelens commented 5 years ago

predict and model.frame are now quite stable, and the interface works well. stat_predict can be done later; most plotting is sufficient with either plot() or by doing a manual predict() and then including it in a plot.