---
title: "Methodology"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Methodology}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

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

```{=html}
<!--
  STATUS: STUB. Each section below is a placeholder that lays out the
  intended scope of the technical manual. Expand the prose with material
  from the manuscript (Córdoba, Paccioretti, Balzarini, in review) and
  the IBC 2024 abstract / presentation.
-->
```

This vignette describes the statistical protocol implemented by the
`ofemeantest` package, intended for readers who want to understand
*what* the function computes — beyond the user-facing recipe in
`vignette("getting-started")`.

## 1. Problem setting

On-farm experiments (OFE) on commercial fields are typically laid out
as side-by-side strips with **no replication** of treatments and no
randomisation. Yield-monitor data over these strips are dense, irregular
and **strongly spatially autocorrelated**, so classical ANOVA on the
raw observations grossly inflates the effective sample size and breaks
its independence assumption.

The protocol implemented here adapts ANOVA to this setting by (i)
aggregating the raw observations into a regular grid, (ii) estimating
the autocorrelation-corrected effective sample size, and (iii)
inferring treatment effects from a distribution of permutational
*p*-values computed on random subsamples of that ESS.

## 2. Grid construction

`make_ofe_grid()` builds a regular grid (square cells, optional rotation
and origin shift) over the bounding box of the trial data and assigns
each yield-monitor point to a cell. A cell is **selected for analysis**
if and only if:

- it contains at least `min_per_cell` observations, and
- all those observations come from a single treatment.

Mixed-treatment cells (those that straddle the boundary between strips)
are dropped. The per-cell response used downstream is the **median**
of the raw observations falling in that cell, which is robust to the
typical outliers of yield-monitor data.

*Implementation details:* the grid origin and cell size are passed
through deterministically; the same arguments always produce the same
cells, so `ofemt(grid = NULL, ...)` and
`ofemt(grid = make_ofe_grid(...))` are interchangeable.

## 3. Effective sample size (ESS)

A one-way ANOVA on cell medians yields residuals \\(e_i\\). Their
spatial autocorrelation is estimated with the approximate profile
likelihood estimator (APLE) of \\(\\rho\\) (Li, Calder & Cressie, 2007),
and the **effective sample size** is computed from Griffith (2005):

\\[
n^{*} = n \\left(1 - a \\cdot \\frac{n-1}{n} \\cdot \\left(1 - \\exp(-b\\rho + c\\sqrt{\\rho})\\right)\\right)
\\]

with the coefficients given in the published formula. We also report
**Moran's *I*** of the residuals as a sanity check on the magnitude of
the autocorrelation. \\(\\rho\\) is constrained to \\([0, 1]\\) before
entering the ESS formula.

## 4. Permutational ANOVA on subsamples

Given the ESS, `ofemt()` runs `n_s` independent permutational ANOVAs.
Each run:

1. samples \\(\\lceil n^{*}/k \\rceil\\) cells **without replacement**
   from each of the \\(k\\) treatments;
2. fits a one-way ANOVA on the subsample;
3. permutes the treatment labels `n_p` times and records the proportion
   of permuted *F*-statistics that exceed the observed one (this is the
   per-run *p*-value).

The reported test statistic for the comparison is the **median** of the
`n_s` per-run *p*-values. The full distribution is retained in
`res$perm_runs` and can be inspected with `plot_pvalue_hist()`.

## 5. Multiple comparisons

When the trial contains more than two treatments, every pairwise
comparison is tested as above and the results are adjusted for
multiplicity via `stats::p.adjust()` (Bonferroni, Holm or
Benjamini–Hochberg, selected by the `p_adjust_method` argument).

The adjustment is applied **within each run**, across the
\\(\\binom{k}{2}\\) comparisons produced by that run, and the reported
`p_adj` is the *median of the adjusted values*. Adjusting inside the run
keeps the two summaries on the same footing: `p_value` and `p_adj` are
both medians of a genuine empirical distribution over the `n_s` runs,
and the median line drawn by `plot_pvalue_hist()` is exactly the value
reported in the table. (Adjusting the median instead — the reverse
order — would mix a summary statistic with a correction computed from a
single vector of medians, and the histogram could then disagree with the
table.) Both per-run quantities are kept in `res$perm_runs` as
`p_value` and `p_adj`.

A compact letter display is then computed from the adjusted *p*-values
with `multcompView::multcompLetters()`. Treatments are ordered by
decreasing median response before the letters are assigned, so `"a"`
marks the highest-yielding group and the letters read monotonically down
the `Means comparison` table. Treatment labels are swapped for internal
placeholders before the letters are computed and restored afterwards, so
labels containing spaces, `+`, `-`, parentheses or accents are reported
verbatim.

## 6. Assumptions and limitations

```{=html}
<!-- TODO: expand with the assumptions listed in the manuscript. -->
```

- The protocol assumes the trial layout produces enough cells per
  treatment for the subsampling step to be meaningful.
- The ESS formula is an approximation valid for approximately
  normally distributed data with non-negative autocorrelation.
- The current implementation supports a **single categorical
  treatment**. Multi-factor designs, factor interactions and
  continuous covariates are not yet supported.

## References

- Córdoba, M., Paccioretti, P., & Balzarini, M. (2025). A new method to compare
  treatments in unreplicated on-farm experimentation. Precision Agriculture,
  26(1), 4. https://doi.org/10.1007/s11119-024-10206-0
- Griffith, D.A. (2005). Effective geographic sample size in the
  presence of spatial autocorrelation. *Annals of the Association of
  American Geographers* 95(4): 740–760.
- Li, H., Calder, C.A. & Cressie, N. (2007). Beyond Moran's *I*:
  Testing for spatial dependence based on the spatial autoregressive
  model. *Geographical Analysis* 39(4): 357–375.
- Vega A., Córdoba M., Balzarini M. (2019). Protocol for automating
  error removal from yield maps. *Precision Agriculture* 20:
  1030–1044.
