---
title: "Tutorial"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Tutorial}
  %\VignetteEncoding{UTF-8}
  %\VignetteEngine{knitr::rmarkdown}
editor_options: 
  chunk_output_type: console
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

```{r loadPackage}
library(rbioacc)
# library(ggplot2)
```

## A Simple Example: *Male Gammarus Single*

### Data

Load the data set with the function `data()`, define the duration of the exposure `time_accumulation`, and check if the data set is correctly imported with the function `modelData()`. Here the data set is called *Male Gammarus Single*

```{r dataMGS}
data("Male_Gammarus_Single")
```

### Inference

The function `fitTK()` performs the inference process.

```{r fitMGS, cache=TRUE, results="hide"}
modelData_MGS <- modelData(Male_Gammarus_Single, time_accumulation = 4)
fit_MGS <- fitTK(modelData_MGS, iter = 10000)
```

### Results

#### TK parameters

The 4 MCMC are stored in the object `fitMCMC`. The quantiles for each TK parameter can be obtained with the `quantile()` function.

```{r statsMGS}
quantile_table(fit_MGS)
```

```{r plotMGS, fig.height=4, fig.width=5}
plot(fit_MGS)
```

```{r ppcMGS, fig.height=4, fig.width=5}
ppc(fit_MGS)
```

## Male Gammarus Merged

```{r fitMGM708, cache=TRUE, results="hide", eval=FALSE}
data("Male_Gammarus_Merged")
data_MGM708 <- Male_Gammarus_Merged[Male_Gammarus_Merged$expw == 7.08021e-05, ]
modelData_MGM708 <- modelData(data_MGM708, time_accumulation = 4)
fit_MGM708 <- fitTK(modelData_MGM708, iter = 10000)
```

```{r statMGM708, eval=FALSE}
quantile_table(fit_MGM708)
```

```{r plotMGM708, fig.height=4, fig.width=5, eval=FALSE}
plot(fit_MGM708)
```

```{r ppcMGM708, fig.height=4, fig.width=5, eval=FALSE}
ppc(fit_MGM708)
```

```{r fitMGM141, cache=TRUE, results="hide", eval=FALSE}
data_MGM141 <- Male_Gammarus_Merged[Male_Gammarus_Merged$expw == 1.41604e-04, ]
modelData_MGM141 <- modelData(data_MGM141, time_accumulation = 7)
fit_MGM141 <- fitTK(modelData_MGM141, iter = 20000)
```

```{r statMGM141, eval=FALSE}
quantile_table(fit_MGM141)
```

```{r plotMGM141, fig.height=4, fig.width=5, eval=FALSE}
plot(fit_MGM141)
```

```{r ppcMGM141, fig.height=4, fig.width=5, eval=FALSE}
ppc(fit_MGM141)
```

```{r fitMGM283, cache=TRUE, results="hide", eval=FALSE}
data_MGM283 <- Male_Gammarus_Merged[Male_Gammarus_Merged$expw == 2.83208e-04, ]
modelData_MGM283 <- modelData(data_MGM283, time_accumulation = 4)
fit_MGM283 <- fitTK(modelData_MGM283, iter = 10000)
```

```{r statMGM283, eval=FALSE}
quantile_table(fit_MGM283)
```

```{r plotMGM283, fig.height=4, fig.width=5, eval=FALSE}
plot(fit_MGM283)
```

```{r ppcMGM283, fig.height=4, fig.width=5, eval=FALSE}
ppc(fit_MGM283)
```

## Male Gammarus seanine with growth

```{r fitMGSG, eval=FALSE}
data("Male_Gammarus_seanine_growth")
modelData_MGSG <- modelData(Male_Gammarus_seanine_growth, time_accumulation = 1.417)
fit_MGSG <- fitTK(modelData_MGSG, iter = 10000)
```


```{r statsMGSG, eval=FALSE}
quantile_table(fit_MGSG)
```

```{r plotMGSG, fig.height=6, fig.width=7, eval=FALSE}
plot(fit_MGSG)
```

```{r ppcMGSG, fig.height=6, fig.width=7, eval=FALSE}
ppc(fit_MGSG)
```

## Oncorhynchus

```{r fitOT440, cache=TRUE, results="hide", eval=FALSE}
data("Oncorhynchus_two") #  Pimephales_two
data_OT440 = Oncorhynchus_two[Oncorhynchus_two$expw == 0.00440,]
modelData_OT440 <- modelData(data_OT440, time_accumulation = 49)
fit_OT440 <- fitTK(modelData_OT440, iter = 10000)
```

```{r statOT440, eval=FALSE}
quantile_table(fit_OT440)
```

```{r plotOT440, fig.height=4, fig.width=5, eval=FALSE}
plot(fit_OT440)
```

```{r ppcOT440, fig.height=4, fig.width=5, eval=FALSE}
ppc(fit_OT440)
```

```{r fitOT041, cache=TRUE, results="hide", eval=FALSE}
data_OT041 <- Oncorhynchus_two[Oncorhynchus_two$expw == 0.00041,]
modelData_OT041 <- modelData(data_OT041, time_accumulation = 49)
fit_OT041 <- fitTK(modelData_OT041, iter = 10000)
```

```{r statOT041, eval=FALSE}
quantile_table(fit_OT041)
```

```{r plotOT041, fig.height=4, fig.width=5, eval=FALSE}
plot(fit_OT041)
```

```{r ppcOT041, fig.height=4, fig.width=5, eval=FALSE}
ppc(fit_OT041)
```

## Chironomus benzo-a-pyrene

```{r fitCB, cache=TRUE, results="hide", eval=FALSE}
data("Chironomus_benzoapyrene")
modelData_CB <- modelData(Chironomus_benzoapyrene, time_accumulation = 3)
modelData_CB$unifMax = modelData_CB$unifMax  * 100
fit_CB <- fitTK(modelData_CB, iter = 10000)
```
 
```{r statCB, eval=FALSE}
quantile_table(fit_CB)
```

```{r plotCB, fig.height=4, fig.width=5, eval=FALSE}
plot(fit_CB)
```

```{r ppcCB, fig.height=4, fig.width=5, eval=FALSE}
ppc(fit_CB)
```

# Prediction

```{r predictMGS, eval=FALSE}
data("Male_Gammarus_Single")
modelData_MGS <- modelData(Male_Gammarus_Single, time_accumulation = 4)
fit_MGS <- fitTK(modelData_MGS, iter = 5000, chains = 3)

# Data 4 prediction should respect the exposure routes
data_4pred <- data.frame( time = 1:25, expw = 4e-5 )
predict_MGS <- predict(fit_MGS, data_4pred)
plot(predict_MGS)
```

```{r predictMGSG, eval=FALSE}
# data("Male_Gammarus_seanine_growth")
# modelData_MGSG <- modelData(Male_Gammarus_seanine_growth, time_accumulation = 4)
# fit_MGSG <- fitTK(modelData_MGSG, iter = 5000, chains = 3)
# 
# # Data 4 prediction should respect the exposure routes
# data_4pred <- data.frame( time = 1:25, expw = 18 )
# predict_MGSG <- predict(fit_MGSG, data_4pred)
# plot(predict_MGSG)
```

```{r predictCC, eval=FALSE}
data("Chiro_Creuzot")
Chiro_Creuzot <- Chiro_Creuzot[Chiro_Creuzot$replicate == 1,]
modelData_CC <- modelData(Chiro_Creuzot, time_accumulation = 1.0)
fit_CC <- fitTK(modelData_CC, iter = 5000, chains = 3)
# --------
quantile_table(fit_CC)

# Data 4 prediction should respect the exposure routes
data_4pred <- data.frame( time = 1:25, expw = 18, exps = 1200, exppw = 15 )
predict_CC <- predict(fit_CC, data_4pred)
plot(predict_CC)
```
