#LIBRARIES library(ggplot2) library(dplyr) library(tidyr) library(car) # For Levene’s Test #install.packages(“ggpubr”) #install.packages(“pwr”) #to check the power analysis library(ggpubr) # For ggqqplot and anova visualization library(DescTools) # For post hoc test if needed library(multcompView) # For letters on boxplots library(AgroR) library(agricolae) library(tidyverse) library(reshape2) library(emmeans) install.packages(“emmeans”) #install.packages(“multcompView”) library(multcomp)

#load data data<- read.csv(“DATA.csv”) view(data) # Convert Treatment to factor for ANOVA data\(Treatment <- as.factor(data\)Treatment) data\(location<- as.factor(data\)location) view(data) summary(data) boxplot(data\(Grain.Yield~data\)Treatment) hist(data$Grain.Yield)

#anova grain yield anova <- aov(Grain.Yield ~ Treatment * location, data = data) summary(anova)

#mean comparison emm<- emmeans(anova, pairwise ~ Treatment | location)

Summarize means and SE for plotting

summary_data <- data %>% group_by(location, Treatment) %>% summarise(mean_yield = mean(Grain.Yield, na.rm = TRUE), se_yield = sd(Grain.Yield, na.rm = TRUE)/sqrt(n()))

#Visualizatiojn

ggplot(summary_data, aes(x = Treatment, y = mean_yield, fill = location)) + geom_bar(stat = “identity”, position = position_dodge(0.9)) + geom_errorbar(aes(ymin = mean_yield - se_yield, ymax = mean_yield + se_yield), position = position_dodge(0.9), width = 0.25) + labs(x = “Treatment”, y = “Grain Yield (t/ha)”) + theme_minimal(base_size = 14) + theme(axis.text.x = element_text(angle = 45, hjust = 1))

#anova plant heigth 2 weeks aph <- aov(data$PH2WAS ~ Treatment * location, data = data) summary(aph)

#mean comparison emmeans(aph, pairwise ~ Treatment | location)

Summarize means and SE for plotting

summary_data <- data %>% group_by(location, Treatment) %>% summarise(mean_ph = mean(PH2WAS, na.rm = TRUE), se_ph = sd(PH2WAS, na.rm = TRUE)/sqrt(n()))

#Visualizatiojn

ggplot(summary_data, aes(x = Treatment, y = mean_ph, fill = location)) + geom_bar(stat = “identity”, position = position_dodge(0.9)) + geom_errorbar(aes(ymin = mean_ph - se_ph, ymax = mean_ph + se_ph), position = position_dodge(0.9), width = 0.25) + labs(x = “Treatment”, y = “Plant heigth 2WAS”) + theme_minimal(base_size = 14) + theme(axis.text.x = element_text(angle = 45, hjust = 1))

#anova plant heigth 4 weeks aph4 <- aov(data$PH4WAS ~ Treatment * location, data = data) summary(aph4)

#mean comparison emmeans(aph4, pairwise ~ Treatment | location)

Summarize means and SE for plotting

summary_data4 <- data %>% group_by(location, Treatment) %>% summarise(mean_ph4 = mean(PH4WAS, na.rm = TRUE), se_ph4 = sd(PH4WAS, na.rm = TRUE)/sqrt(n()))

#Visualizatiojn

ggplot(summary_data4, aes(x = Treatment, y = mean_ph4, fill = location)) + geom_bar(stat = “identity”, position = position_dodge(0.9)) + geom_errorbar(aes(ymin = mean_ph4 - se_ph4, ymax = mean_ph4 + se_ph4), position = position_dodge(0.9), width = 0.25) + labs(x = “Treatment”, y = “Plant heigth 4WAS”) + theme_minimal(base_size = 14) + theme(axis.text.x = element_text(angle = 45, hjust = 1))

#anova plant heigth 6 weeks aph6 <- aov(data$PH6WAS ~ Treatment * location, data = data) summary(aph6)

#mean comparison emmeans(aph6, pairwise ~ Treatment | location)

Summarize means and SE for plotting

summary_data <- data %>% group_by(location, Treatment) %>% summarise(mean_ph6 = mean(PH6WAS, na.rm = TRUE), se_ph6 = sd(PH6WAS, na.rm = TRUE)/sqrt(n()))

#Visualizatiojn

ggplot(summary_data, aes(x = Treatment, y = mean_ph6, fill = location)) + geom_bar(stat = “identity”, position = position_dodge(0.9)) + geom_errorbar(aes(ymin = mean_ph6 - se_ph6, ymax = mean_ph6 + se_ph6), position = position_dodge(0.9), width = 0.25) + labs(x = “Treatment”, y = “Plant heigth 6WAS”) + theme_minimal(base_size = 14) + theme(axis.text.x = element_text(angle = 45, hjust = 1))