---
title: "Compare EC50 models"
author: "Kaique S Alves"
date: "`r Sys.Date()`"
output:
  rmarkdown::html_vignette:
    toc: true
    toc_depth: 2
vignette: >
  %\VignetteIndexEntry{Compare EC50 models}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

## When to Use This Workflow

Use `ec50_multimodel()` when you need to compare candidate dose-response models
before reporting EC50 estimates. The function fits each model within each
isolate and stratum, stores the fitted `drc` objects, and returns a data frame
with EC50 estimates plus model-selection statistics.

## Packages and Data

```{r message=FALSE, warning=FALSE}
library(ec50estimator)
library(drc)
library(ggplot2)
```

```{r}
data(multi_isolate)

example_data <- subset(
  multi_isolate,
  isolate %in% 1:5 & fungicida == "Fungicide A"
)

head(example_data)
```

## Inspect the Raw Pattern

Repeated observations at the same dose should be shown as points, not connected
with lines. Lines are reserved for fitted model curves. The exploratory plot
below uses positive doses because the x-axis is log-scaled.

```{r fig.alt="Jittered raw dose-response observations for five isolates, faceted by field system."}
raw_plot_data <- subset(example_data, dose > 0)

ggplot(raw_plot_data, aes(dose, growth, color = factor(isolate))) +
  geom_jitter(width = 0.08, height = 0, alpha = 0.85, size = 2) +
  facet_wrap(~field) +
  scale_x_log10() +
  labs(x = "Dose", y = "Growth", color = "Isolate") +
  theme_light()
```

## Check the Data

Use `check_ec50_data()` before model fitting. The returned flags help separate
data problems from model-fitting problems.

```{r}
check_ec50_data(
  example_data,
  response = "growth",
  dose = "dose",
  isolate = "isolate",
  strata = "field"
)
```

## Fit Candidate Models

Pass a named or unnamed list of `drc` model functions to `fct`. The fitted object
keeps the EC50 estimates as rows and stores the fitted models for downstream
helpers.

```{r}
fit <- ec50_multimodel(
  growth ~ dose,
  data = example_data,
  isolate_col = "isolate",
  strata_col = "field",
  fct = list(drc::LL.3(), drc::LL.4(), drc::W2.3()),
  interval = "delta"
)

head(ec50_estimates(fit))
```

## Diagnose Fits

`fit_quality()` summarizes successful model fits. `fit_failures()` returns
failed group/model combinations as data.

```{r}
fit_quality(fit)

fit_failures(fit)
```

## Rank Candidate Models

`model_selection()` ranks models inside each isolate and field. The default
criterion is `IC`, which is returned by `ec50_multimodel()`. Smaller values are
better.

```{r}
selection <- model_selection(fit)

selection[, c("ID", "field", "model", "IC", "delta", "weight", "rank")]
```

`best_model()` keeps the top-ranked model in each group.

```{r}
best <- best_model(fit)

best[, c("ID", "field", "model", "Estimate", "Lower", "Upper", "IC")]
```

Model weights are relative to the candidate set supplied by the user. They are
useful for ranking fitted models, but they do not prove that a model is
biologically correct.

## Plot All Models or Best Models

Use `plot_EC50_curves()` directly on the fitted object. With multiple models,
line type identifies the model.

```{r fig.alt="Faceted dose-response plot with raw observations and all candidate fitted curves for each isolate and field."}
plot_EC50_curves(fit)
```

For a cleaner reporting figure, plot only the best-supported model in each
isolate and stratum.

```{r fig.alt="Faceted dose-response plot with raw observations and the best-supported fitted curve for each isolate and field."}
plot_EC50_curves(fit, models = "best")
```

## Customize with Curve Data

Use `curve_data()` when you want to build a custom `ggplot2` figure. This helper
uses the stored models and does not refit them.

```{r}
curves <- curve_data(fit)

head(curves)
```

```{r fig.alt="Custom ggplot2 curve plot built from fitted EC50 curve coordinates."}
ggplot(curves, aes(dose, growth, color = factor(isolate), linetype = model)) +
  geom_line(linewidth = 1) +
  facet_wrap(~field) +
  scale_x_log10() +
  labs(x = "Dose", y = "Growth", color = "Isolate") +
  theme_light()
```

## Predict and Report Selected Models

Use `models = "best"` when you want predictions or tables from only the
best-supported model in each group.

```{r}
predict_ec50(
  fit,
  dose = c(0.001, 0.01, 0.1),
  models = "best"
)
```

```{r}
report_ec50(fit, models = "best")
```

## Residual Diagnostics

Residuals can be returned as data or plotted with `ggplot2`.

```{r}
head(residual_data(fit, models = "best"))
```

```{r fig.alt="Residual diagnostic plot for best-supported EC50 models."}
plot_residuals(fit, models = "best")
```