--- title: "Number of replications" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Number of replications} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, include = FALSE} knitr::opts_chunk$set(collapse = TRUE, comment = "#>", fig.width = 6.8, fig.height = 4.4, dpi = 120) ``` ```{r setup} library(trialSizing) ``` ## Theory Choosing the plot size settles one half of the design; the other half is how many times to replicate it. The two are linked: the plot size determines the experimental precision, expressed as $CV_{Xo}$, and that precision determines how many replications are needed to detect a given difference between treatment means. The starting point is the least significant difference of the Tukey test, expressed as a percentage of the experiment mean: $$ d = \frac{q_{\alpha(i;\,GLE)} \sqrt{QME / r}}{m} \times 100 $$ where $q_{\alpha(i;GLE)}$ is the Tukey critical value for $i$ treatments and $GLE$ error degrees of freedom, $QME$ is the error mean square, $r$ is the number of replications and $m$ the mean. Substituting the experimental coefficient of variation, $CV = 100\sqrt{QME}/m$, and solving for $r$: $$ r = \left(\frac{q_{\alpha(i;\,GLE)} \cdot CV}{d}\right)^{2} $$ with error degrees of freedom $$ GLE = i\,(r - 1) \ \text{ for CRD}, \qquad GLE = (i - 1)(r - 1) \ \text{ for RCBD} $$ Here $CV$ is $CV_{Xo}$: the CV expected once the experiment uses the optimal plot size. ### Why this needs iteration The expression is circular. The critical value $q_\alpha$ depends on the error degrees of freedom, which depend on $r$, which is what we are solving for. So $r$ is obtained iteratively: guess $r$, compute $GLE$ and $q_\alpha$, get a new $r$, repeat until it stops moving. `calc_replicates()` runs that fixed point and reports **two readings** of the result: - `r_continuous`: the converged value, e.g. 10.46. This is what the published tables print, and what to quote when comparing with an article. - `r_optimal`: the practical integer, `ceiling(r_continuous)`, floored at 2. You cannot run 10.46 replications, and no design runs fewer than 2. Smaller $d$ means finer resolution, and the cost is steep: because $r$ grows with the square of $1/d$, halving the difference you want to detect roughly quadruples the replications needed. ### References Method: Cargnelutti Filho, A., Alves, B. M., Toebe, M., Burin, C., Santos, G. O., Facco, G., Neu, I. M. M. & Stefanello, R. B. (2014). Tamanho de parcela e número de repetições em aveia preta. *Ciência Rural*, 44(10), 1732-1739. The examples below reproduce Tables 2 (CRD) and 3 (RCBD) of that article, which use black oat with $X_o = 4.14$ m² and $CV_{Xo} = 9.25\%$. ## Basic use ```{r basic} reps <- calc_replicates( treatments = c(3, 10, 50), cv_percent = 9.25, lsd_percent = c(10, 20, 30), design = "CRD" ) reps ``` Reading the first row: with 3 treatments and a CV of 9.25%, detecting a difference of 10% of the mean needs 10.46 replications in theory, so 11 in practice. The published table gives 10.46 for that cell. ```{r summary} summary(reps) ``` The full grid is in `$data`, one row per treatments × LSD combination: ```{r access} head(reps$data[, c("Treatments", "LSD_percent", "r_continuous", "r_optimal", "df_error", "q_tukey")]) ``` ## Design, significance level Randomized complete blocks need slightly more replications than a completely randomized design at the same precision, because blocking costs error degrees of freedom. The gap narrows as the number of treatments grows: ```{r design} rbind( CRD = calc_replicates(c(3, 50), 9.25, 10, design = "CRD")$data$r_continuous, RCBD = calc_replicates(c(3, 50), 9.25, 10, design = "RCBD")$data$r_continuous ) ``` With 3 treatments the difference is 10.46 against 10.95; with 50 it has practically vanished. Both match the published tables. A stricter significance level raises the critical value and therefore the replications: ```{r alpha} rbind( `alpha = 0.05` = calc_replicates(10, 9.25, 20, alpha = 0.05)$data$r_continuous, `alpha = 0.01` = calc_replicates(10, 9.25, 20, alpha = 0.01)$data$r_continuous ) ``` ## The full workflow In practice $CV_{Xo}$ comes from one of the plot-size methods rather than being typed in. Fitting the bundled simulated trial (`?uniformity_trial`) and feeding the result straight through: ```{r workflow, message = FALSE} grid1 <- as.matrix(uniformity_trial[uniformity_trial$trial == "T1", grep("^col", names(uniformity_trial))]) cv_tab <- calc_cv_shapes(grid1) lrp <- fit_lrp(cv_tab, x = "x", cv = "cv", step = 0.01) cvxo <- unname(lrp$parameters["Breakpoint_Response"]) cvxo calc_replicates(treatments = c(5, 10, 20), cv_percent = cvxo, lsd_percent = c(10, 20), design = "RCBD")$data[ , c("Treatments", "LSD_percent", "r_continuous", "r_optimal")] ``` So a trial designed with plots of about 9 m² would carry an expected CV near 7%, and the table reads off how many replications are needed to detect a given difference between treatment means. ## Plot Over a range of treatments the trade-off becomes visual: one line per LSD level, showing how the requirement climbs as the number of treatments grows and as the target difference tightens. ```{r plot} reps <- calc_replicates(treatments = 3:50, cv_percent = 9.25, lsd_percent = c(10, 20, 30), design = "CRD") plot(reps, title = "Replications needed (black oat)") ``` The default plots the integer `r_optimal`, which is why the lines are stepped. For the smooth theoretical curve, plot the continuous value: ```{r plot-continuous} plot(reps, y_var = "r_continuous", title = "Replications (continuous)") ``` ```{r save, eval = FALSE} plot(reps, title = "Replications needed", save = TRUE, file = "replications.tiff", format = "tiff", dpi = 300) ``` ## Reading the results The article's conclusion illustrates how these numbers are used in practice. With a CV of 9.25%, detecting differences of 10% of the mean would need around 27 replications for 50 treatments, which is not feasible in the field. Turning the question around and fixing $r = 4$, a common choice in the literature, the detectable difference becomes about 26.7% of the mean. That is the honest statement of what the experiment can and cannot resolve. The function does not pick a precision for you, and it should not: that decision depends on how large a difference matters agronomically and how much field area is available. What it does is make the trade-off explicit before the experiment is installed rather than after. ## Where to go next `vignette("lrp")`, `vignette("qrp")` and `vignette("mcm")` cover the CV-based plot-size methods that supply $CV_{Xo}$; `vignette("paranaiba")` derives it directly from the raw trial grid.