AGRON Info-Tech
November 02, 2021
Introduction
Analysis of variance (ANOVA) is commonly used to evaluate treatment effects across biological and agricultural experiments. Preparing a clear table requires more than formatting: the model must reflect the design, and the reported degrees of freedom, mean squares, test statistics, and p-values must come from the correct analysis.
This tutorial builds a reusable R workflow that fits the specified ANOVA model for multiple response variables and assembles the results into a publication-ready table. The workflow reduces repetitive manual steps, but it does not choose the experimental model for you: verify the design terms and error structure before interpreting or reporting any result.
So let’s get started!
Load the packages used to import the workbook, manipulate data, and format the final table. Install any package that is not already available in your R library before calling library().
library(readxl)
library(dplyr)
library(tibble)
library(flextable)The example workbook contains a replication variable, two treatment factors, and several numeric response variables. The code below assumes that the first three columns identify replication and treatments, while columns 4 onward contain outcomes. Check this structure against your own file before running the analysis. The example fits a two-factor randomized complete block design; adapt the model terms if your design differs.
data = readxl::read_excel(
path = "Data.xlsx",
col_names = TRUE
)
head(data)# # A tibble: 6 x 10
# Rep Water Priming `Plant height` Spikelets `Spike length` `Grain per spike`
# <dbl> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
# 1 1 CONV NP 100 204 12 35
# 2 1 CONV HP 103 293 15.3 59
# 3 1 CONV OP 95 322 17.5 64
# 4 1 FLOOD NP 99 284 12 34
# 5 1 FLOOD HP 86 304 12.2 44
# 6 1 FLOOD OP 88 326 16.8 69
# # ... with 3 more variables: grain weight <dbl>, Biological yield <dbl>,
# # Grain yield <dbl>Spreadsheet imports can store treatment labels as character strings. Convert the design variables to factors so that R treats them as categorical predictors in the model. Keep response columns numeric and inspect the structure after conversion.
# Convert categorical variables to factor variables
data$Rep = as.factor(data$Rep)
data$Water = as.factor(data$Water)
data$Priming = as.factor(data$Priming)
str(data)# tibble [27 x 10] (S3: tbl_df/tbl/data.frame)
# $ Rep : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 2 ...
# $ Water : Factor w/ 3 levels "CONV","FLOOD",..: 1 1 1 2 2 2 3 3 3 1 ...
# $ Priming : Factor w/ 3 levels "HP","NP","OP": 2 1 3 2 1 3 2 1 3 2 ...
# $ Plant height : num [1:27] 100 103 95 99 86 88 91 89 95 97 ...
# $ Spikelets : num [1:27] 204 293 322 284 304 326 274 291 313 252 ...
# $ Spike length : num [1:27] 12 15.3 17.5 12 12.2 ...
# $ Grain per spike : num [1:27] 35 59 64 34 44 69 32 37 48 32 ...
# $ grain weight : num [1:27] 37 46 56.5 37.5 49.8 59 39.5 44.8 50.6 40.5 ...
# $ Biological yield: num [1:27] 13.3 16.5 17.4 12.3 15 ...
# $ Grain yield : num [1:27] 4.2 5.3 5.85 4.5 5.3 6.1 4.8 5.2 6.3 4.03 ...The example uses attach() so column names can be referenced directly in later commands. For reusable analysis scripts, explicit data arguments (for example, data = data) are safer because they avoid name conflicts and make each variable's source clear.
attach(data)Creating for loop function
The next section wraps the repeated analysis steps in a loop. For each response variable, the workflow fits the specified model, selects the required ANOVA columns, applies consistent labels and rounding, and stores the result. Review each step before applying the loop to a new design.
Assign column names
The response variables begin in column 4 in this workbook. The index vector therefore selects columns 4 through the final column. If your workbook places outcomes elsewhere, update this range; do not rely on column position if the file structure may change.
cols <- names(data)[4:ncol(data)]
cols# [1] "Plant height" "Spikelets" "Spike length" "Grain per spike"
# [5] "grain weight" "Biological yield" "Grain yield"Analysis of variance model
For each outcome, the lapply() call applies the model-fitting operation and returns a list of results. The model formula is supplied through the term labels shown in the code. Replace those terms with a formula justified by your experimental design, including the correct blocking or error terms where required.
aov.model <- lapply(cols, FUN = function(x)
aov(reformulate(termlabels = "Rep + Water*Priming",
response = x),
data = data))
anova(aov.model[[2]])# Analysis of Variance Table
#
# Response: Spikelets
# Df Sum Sq Mean Sq F value Pr(>F)
# Rep 2 15894.9 7947.5 11.2003 0.0009083 ***
# Water 2 533.1 266.5 0.3756 0.6927466
# Priming 2 28332.7 14166.3 19.9646 4.486e-05 ***
# Water:Priming 4 918.1 229.5 0.3235 0.8581023
# Residuals 16 11353.2 709.6
# ---
# Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1The ANOVA table contains degrees of freedom, sums of squares, mean squares, F statistics, and p-values. This workflow retains the degrees of freedom, mean squares, and p-values for display. If you choose to report F statistics instead, select the corresponding column and label it accurately.
# print df, MS and Pvalue
final = anova(aov.model[[2]])[,c(1,3,5)]
final# Df Mean Sq Pr(>F)
# Rep 2 7947.5 0.0009083 ***
# Water 2 266.5 0.6927466
# Priming 2 14166.3 4.486e-05 ***
# Water:Priming 4 229.5 0.8581023
# Residuals 16 709.6
# ---
# Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Set the row labels to the model terms so that the final table identifies each source of variation. Confirm that the labels correspond to the terms actually fitted, especially when the model is changed for a different design.
# Getting rownames
rnames = rownames(final)
rnames# [1] "Rep" "Water" "Priming" "Water:Priming"
# [5] "Residuals"Setting column names
Assign concise column headings that identify the response variable and statistic. Use labels that remain unambiguous when the tables are combined across outcomes.
# Setting column names
colnames(final) = c("DF", "MS", "P-value")
colnames(final)[2] = cols[2]
final# DF Spikelets P-value
# Rep 2 7947.5 0.00091
# Water 2 266.5 0.69275
# Priming 2 14166.3 0.00004
# Water:Priming 4 229.5 0.85810
# Residuals 16 709.6Rounding values
Round displayed values to three decimal places for readability. Keep full-precision values in the analysis object; rounding is a presentation choice and should not affect the fitted tests or conclusions.
# Rounding values to 2 decimal place
final = as.data.frame(round(final, digits = 2))
final# DF Spikelets P-value
# Rep 2 7947.47 0.00
# Water 2 266.55 0.69
# Priming 2 14166.35 0.00
# Water:Priming 4 229.53 0.86
# Residuals 16 709.58 NAAssign astericks
Rather than adding a separate p-value column for every response, this example maps p-values to significance labels: * for p < 0.05, ** for p < 0.01, and ns for p ≥ 0.05. These thresholds are conventions, not measures of effect size. Retain exact p-values where they are needed for transparent reporting, and state the thresholds in the table note.
# Assign astericks according to p values
final$sign[final$`P-value` < 0.05] <- "*"
final$sign[final$`P-value` < 0.01] <- "**"
final$sign[final$`P-value` > 0.05] <- "ns"
# Merge MS and significance column together
final[[2]] = paste(final[[2]], ifelse(is.na(final[[4]]), "", final[[4]]))
final# DF Spikelets P-value sign
# Rep 2 7947.47 ** 0.00 **
# Water 2 266.55 ns 0.69 ns
# Priming 2 14166.35 ** 0.00 **
# Water:Priming 4 229.53 ns 0.86 ns
# Residuals 16 709.58 NA <NA>After creating the significance labels, combine them with the mean-square values for compact display. Check the resulting column order and labels so readers can distinguish the statistic from its significance indicator.
final = final[-c(3,4)]
final# DF Spikelets
# Rep 2 7947.47 **
# Water 2 266.55 ns
# Priming 2 14166.35 **
# Water:Priming 4 229.53 ns
# Residuals 16 709.58Create list of tables in excel sheets
The script writes each response-specific table to a separate workbook sheet. This intermediate export is used to collect the tables in the next step; it is not a substitute for checking the final table for completeness, consistent rounding, and correct labels.
anova = writexl::write_xlsx(final, path = paste(cols[i], '-ANOVA.xlsx'))Print final ANOVA object
The following code reads the generated sheets back into R. Confirm that the file path, sheet names, and imported column types are correct before combining the results.
# Print final ANOVA object
file.list <- list.files(pattern='*-ANOVA.xlsx')
df.list <- lapply(X = file.list, FUN = read_excel)
df.list# [[1]]
# # A tibble: 5 x 2
# DF `Biological yield`
# <dbl> <chr>
# 1 2 1.93 ns
# 2 2 0.2 ns
# 3 2 46.74 **
# 4 4 0.18 ns
# 5 16 0.7
#
# [[2]]
# # A tibble: 5 x 2
# DF `Grain per spike`
# <dbl> <chr>
# 1 2 24.41 ns
# 2 2 21.12 ns
# 3 2 1374.04 **
# 4 4 30.41 ns
# 5 16 42.05
#
# [[3]]
# # A tibble: 5 x 2
# DF `grain weight`
# <dbl> <chr>
# 1 2 79.42 **
# 2 2 3.37 ns
# 3 2 836.57 **
# 4 4 7 ns
# 5 16 5.23
#
# [[4]]
# # A tibble: 5 x 2
# DF `Grain yield`
# <dbl> <chr>
# 1 2 0.27 **
# 2 2 0.2 *
# 3 2 8.58 **
# 4 4 0.02 ns
# 5 16 0.03
#
# [[5]]
# # A tibble: 5 x 2
# DF `Plant height`
# <dbl> <chr>
# 1 2 0.45 ns
# 2 2 29.93 ns
# 3 2 0.5 ns
# 4 4 15.85 ns
# 5 16 15.5
#
# [[6]]
# # A tibble: 5 x 2
# DF `Spike length`
# <dbl> <chr>
# 1 2 8.71 *
# 2 2 1.1 ns
# 3 2 40.3 **
# 4 4 1.44 ns
# 5 16 1.27
#
# [[7]]
# # A tibble: 5 x 2
# DF Spikelets
# <dbl> <chr>
# 1 2 7947.47 **
# 2 2 266.55 ns
# 3 2 14166.35 **
# 4 4 229.53 ns
# 5 16 709.58Combine all variables
The list.cbind() step places response-specific results side by side. Because degrees of freedom are common to the model terms, repeated copies may appear; the code removes duplicate columns before producing the combined table. Verify that this operation removes only redundant columns and preserves all outcome-specific statistics.
# Combined ANOVA table for all variables
aov.table = rlist::list.cbind(df.list)
# Remove duplicate columns for DF
dup.cols = which(duplicated(names(aov.table)))
aov.table = aov.table[,-dup.cols]
aov.table# DF Biological yield Grain per spike grain weight Grain yield Plant height
# 1 2 1.93 ns 24.41 ns 79.42 ** 0.27 ** 0.45 ns
# 2 2 0.2 ns 21.12 ns 3.37 ns 0.2 * 29.93 ns
# 3 2 46.74 ** 1374.04 ** 836.57 ** 8.58 ** 0.5 ns
# 4 4 0.18 ns 30.41 ns 7 ns 0.02 ns 15.85 ns
# 5 16 0.7 42.05 5.23 0.03 15.5
# Spike length Spikelets
# 1 8.71 * 7947.47 **
# 2 1.1 ns 266.55 ns
# 3 40.3 ** 14166.35 **
# 4 1.44 ns 229.53 ns
# 5 1.27 709.58Restore clear row labels for the sources of variation in the combined table. These labels should match the fitted model and the order of rows in the ANOVA output.
# Names for sources of variation in ANOVA
rownames(aov.table) = rnames
aov.table# DF Biological yield Grain per spike grain weight Grain yield
# Rep 2 1.93 ns 24.41 ns 79.42 ** 0.27 **
# Water 2 0.2 ns 21.12 ns 3.37 ns 0.2 *
# Priming 2 46.74 ** 1374.04 ** 836.57 ** 8.58 **
# Water:Priming 4 0.18 ns 30.41 ns 7 ns 0.02 ns
# Residuals 16 0.7 42.05 5.23 0.03
# Plant height Spike length Spikelets
# Rep 0.45 ns 8.71 * 7947.47 **
# Water 29.93 ns 1.1 ns 266.55 ns
# Priming 0.5 ns 40.3 ** 14166.35 **
# Water:Priming 15.85 ns 1.44 ns 229.53 ns
# Residuals 15.5 1.27 709.58The loop repeats the workflow over the selected response columns. The example excludes the first three design columns and processes the remaining outcomes. Before running it on another data set, confirm the selected range, response types, missing-data handling, and model formula.
Printing publication ready ANOVA table
Use flextable() to format the assembled results for a Word report. The code places the row labels in a dedicated column and emphasizes the header. Review table width, decimal precision, footnotes, and page fit before including the table in a manuscript or report.
Render the report in the required output format. Inspect the generated document to ensure that the table is not clipped and that the model terms, labels, and significance notes are legible.
The final table combines results across the response variables. Treat it as a formatted presentation of the fitted models: statistical interpretation still depends on the design, model assumptions, diagnostics, and planned comparisons relevant to the experiment.
For another experiment, update the factor conversions, response-column range, and model formula to match the data and design. Then verify the ANOVA output for each response before exporting. A loop can make repetitive analysis consistent, but it cannot replace design-aware model specification or diagnostic checks.
Please comment below if you have any questions.
Download data file — Click_here
Download Rscript — Click_here
Download R program —
Click_hereDownload R studio —
Click_here
good
ReplyDeleteThanks for your support
DeleteThank you very much for the vital information. May you demystify on the ANOVA using alpha lattice design.
ReplyDeleteGood web site! I really love how it is simple on my eyes and the data are well written. I am wondering how I might be notified whenever a new post has been made. I've subscribed to your RSS feed which must do the trick! Have a great day!quite gaming mouse
ReplyDeleteHello Dear
ReplyDeleteIt is Good article Blog Thanks for sharing.
affordable video production
brand marketing video production washington dc
video production companies in washington dc
video production service washington dc
video production washington dc
affordable video marketing
affordable video production services
corporate video production dc
corporate video production washington dc
video production services dc
Thanks for sharing such an informative Article. It will be beneficial to those who seek information. Continue to share your knowledge through articles like these.
ReplyDeleteData Engineering Solutions
Artificial Intelligence Services
Data Analytics Services
Data Modernization Services
I appreciate you taking the time and effort to share your knowledge. This material proved to be really efficient and beneficial to me. Thank you very much for providing this information. Continue to write your blog.
ReplyDeleteData Engineering Services
Machine Learning Solutions
Data Analytics Solutions
Data Modernization Services