---
title: "Case study: zebrafish thermal tolerance under hypoxia, normoxia, and hyperoxia"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Case study: zebrafish thermal tolerance under hypoxia, normoxia, and hyperoxia}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 6,
  fig.height = 4,
  dpi = 96
)
```

This article fits a lethal thermal-death-time (TDT) model to zebrafish
(*Danio rerio*) larvae assayed under three oxygen treatments — **hypoxia,
normoxia, and hyperoxia** — and asks a question motivated by the **oxygen- and
capacity-limited thermal tolerance** (OCLTT) hypothesis (Pörtner 2001): does
oxygen availability shift the upper thermal limit? OCLTT is actively debated
rather than settled (see, e.g., Jutfelt et al. 2018), so the framing here is
hypothesis-testing, not mechanism confirmation: if oxygen matters, the critical
temperature `CTmax` and the thermal-sensitivity slope `z` should shift with the
oxygen treatment. It mirrors the zebrafish case
study in [`bayesTLS`](https://github.com/daniel1noble/bayesTLS), fit here by
maximum likelihood instead of MCMC.

The fit runs live (TMB, **no Stan**) and the treatment effect is the direct
CTmax/z parameterisation — one `CTmax` and one `z` per oxygen level. Uncertainty
is summarised by **confidence intervals**; for this well-identified three-group
fit we use the fast Wald intervals, while the profile-likelihood intervals that
are the `freqTLS` default are showcased on the single-fit studies (e.g.
`vignette("case-study-shrimp")`).

```{r setup}
library(freqTLS)
```

## The dataset and the question

`zebrafish_o2` is vendored from `bayesTLS` (Saruhashi et al. 2026) under
CC BY 4.0; cite the source and `citation("freqTLS")` when you use it
(see `?zebrafish_o2`). It records survival counts in a temperature × duration
grid for diploid and triploid larvae across the three oxygen treatments. We take
the **diploid** larvae — the focal group; triploids are a separate extension —
and let `CTmax` and `z` depend on the oxygen treatment.

```{r data}
data(zebrafish_o2)
zf <- subset(zebrafish_o2, ploidy == "diploid")
table(zf$oxygen)
```

```{r standardize}
std <- standardize_data(
  zf,
  temp = "temp", duration = "duration_min",
  n_total = "n_total", n_surv = "n_surv",
  duration_unit = "minutes"
)
```

The design spans assay temperatures from `r min(zf$temp)` to `r max(zf$temp)` °C
and exposures up to `r max(zf$duration_min)` minutes, replicated across the three
oxygen treatments — `r nrow(zf)` rows in total. `oxygen` is the covariate of
interest.

## The grouped fit (live)

A single `fit_4pl()` call fits all three oxygen treatments at once, each with its
own `CTmax` and `z` but a shared constant 4PL shape (the temperature effect runs
through the midpoint). `CTmax` is reported at a one-hour reference exposure
(`t_ref = 60` minutes).

```{r fit}
fit <- suppressWarnings(fit_4pl(
  std,
  ctmax = ~ 0 + oxygen,
  z     = ~ 0 + oxygen,
  t_ref = 60,
  family = "beta_binomial"
))
fit
```

Per-treatment `CTmax` and `z` with confidence intervals:

```{r tls}
ox <- tls(fit, by = "oxygen", method = "wald")
ox$summary
```

```{r ranking, echo = FALSE}
ct <- ox$summary[ox$summary$quantity == "CTmax", ]
ct <- ct[order(ct$median), ]
```

Evaluate the OCLTT prediction with `CTmax` (the one-hour critical temperature): the
treatments order from lowest to highest as
**`r paste(ct$oxygen, sprintf("(%.1f °C)", ct$median), collapse = " < ")`**.
Under the OCLTT hypothesis, restricting oxygen should *lower* the thermal limit,
so hypoxia is expected at the cold end and hyperoxia at the warm end; the fitted
ordering and the overlap (or separation) of the intervals are what the data say
about that prediction.

## A Confidence Eye per treatment

The default `freqTLS` uncertainty visual is the **Confidence Eye** — a confidence
lens with a hollow point estimate, carrying no prior and making no probability
statement about the parameter. One lens per oxygen treatment makes the `CTmax`
comparison legible; the parameter names are read off the fit so they always match
the fitted labels.

```{r eye-ctmax, fig.height = 4.5, fig.alt = "Confidence Eyes for zebrafish CTmax by oxygen treatment: three pale confidence lenses with hollow point estimates, one per treatment, the freqTLS uncertainty display rather than a posterior density."}
plot_confidence_eye(fit, parm = get_ctmax(fit)$parameter, method = "wald")
```

The same display for the thermal-sensitivity slope `z`:

```{r eye-z, fig.height = 4.5, fig.alt = "Confidence Eyes for zebrafish z (thermal sensitivity) by oxygen treatment: three pale confidence lenses with hollow point estimates."}
plot_confidence_eye(fit, parm = get_z(fit)$parameter, method = "wald")
```

## Do the treatments differ? Reading the intervals

The OCLTT prediction is a statement about *differences* between treatments. The
honest frequentist reading is the overlap of the per-treatment confidence
intervals above: where two treatments' `CTmax` intervals are disjoint, the data
reject equality at the 95% level; where they overlap, the data are consistent
with no difference. `freqTLS` also supports an explicit contrast (the difference
in `CTmax` or `z` between two treatments, with its own confidence interval and a
likelihood-ratio test) — see `vignette("frequentist-and-bayesian")` and
`?confint.profile_tls` for the contrast interface. Because the relative threshold
is the curve midpoint, these comparisons are prior-free and respect the
asymmetry of the likelihood.

## Seeing the fit

The fitted survival surface, one curve per oxygen treatment and assay
temperature, with the observed proportions overlaid:

```{r survival, fig.height = 4.5, fig.alt = "Fitted zebrafish survival curves under three oxygen treatments: survival probability declining with exposure duration, observed proportions overlaid."}
plot_survival_curves(fit)
```

## Boundary: what this case study does not cover

We fit only the **diploid** larvae; the triploid larvae are a separate ploidy
contrast left as an extension. The fit holds the 4PL **shape** constant across
treatments (the temperature effect runs through the midpoint), matching the
`bayesTLS` configuration; relaxing that to a per-group shape (letting `low`, `up`,
`k` vary) is the shape-covariate feature demonstrated in
`vignette("comparing-to-bayesTLS")`. As a lethal endpoint, the survival-count 4PL
applies directly; the sublethal loss-of-equilibrium endpoint is out of scope for
`freqTLS`.

**Where to next:** the aphid case study adds formal species contrasts;
`vignette("random-effects")` covers hierarchical (many-group) designs.