1. Draw a plot to show the experiment design in the field

Varieties: A, B, C, D

replicates: 5

arrangement: random

x=rep(1:5,4)
y=rep(1:4,5)
z=sample(rep(c("A", "B", "C","D"), 5))
z
##  [1] "A" "D" "C" "C" "A" "B" "B" "B" "C" "D" "B" "B" "D" "A" "C" "C" "A" "D" "A"
## [20] "D"
design = data.frame(x=factor(x), y=factor(y), z = z)
design
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.0.5
ggplot(data=design,
       aes(x=x,y=y))+
  geom_tile(fill = y,
            alpha = 0.2,
            col = "blue")+
  geom_text(aes(label = z))+
  theme_void()
The design of experiment

The design of experiment

2. Randomly generate the yield data following normal distribution

for A (mean = 50, sd = 3), B (mean = 45, sd =2),

C (mean = 40, sd = 5), and D (mean = 52, sd = 2).

yield = seq(30, 60,0.5)
y_A = dnorm(yield, mean = 50, sd = 3)
y_B = dnorm (yield,mean = 45, sd =2)
y_C = dnorm (yield,mean = 40, sd =5)
y_D = dnorm (yield,mean = 52, sd =2)

mydata = data.frame(yield,y_A,y_B,y_C, y_D)

ggplot(data = mydata,
       aes(x=yield, y = y_A))+
  geom_line(col = "green")+
  geom_vline(xintercept = 50,col = "green")+
  geom_line(aes(y=y_B),
            col="red")+
  geom_vline(xintercept=45, col = "red")+
  geom_line(aes(y=y_C),
            col="blue")+
  geom_vline(xintercept=40, col = "blue")+
  geom_line(aes(y=y_D),
            col="purple")+
  geom_vline(xintercept=52, col = "purple")

Yield_A = rnorm(5, mean = 50, sd =3)
Yield_B = rnorm(5, mean = 45, sd = 2)
Yield_C = rnorm(5, mean = 40, sd =5)
Yield_D = rnorm(5, mean = 52, sd = 2)

Yield = c(Yield_A, Yield_B,Yield_C, Yield_D)
Yield
##  [1] 49.64649 51.53185 47.35195 50.90054 54.21135 48.78155 44.52975 48.04881
##  [9] 47.18506 46.85433 50.75428 42.68644 37.59236 44.91129 41.71280 52.19397
## [17] 53.41457 53.01384 49.84794 54.22150
Variety = rep(c("A", "B", "C", "D"), each = 5)
Variety
##  [1] "A" "A" "A" "A" "A" "B" "B" "B" "B" "B" "C" "C" "C" "C" "C" "D" "D" "D" "D"
## [20] "D"
Yield_data = data.frame(Yield,Variety)
Yield_data
summary(Yield_data)
##      Yield         Variety         
##  Min.   :37.59   Length:20         
##  1st Qu.:46.37   Class :character  
##  Median :49.21   Mode  :character  
##  Mean   :48.47                     
##  3rd Qu.:51.70                     
##  Max.   :54.22
mean = with(Yield_data, tapply(Yield, Variety,mean))
mean
##        A        B        C        D 
## 50.72843 47.07990 43.53144 52.53837
sd = with(Yield_data, tapply(Yield, Variety,sd))
sd
##        A        B        C        D 
## 2.518784 1.612422 4.831346 1.672105

3. Test whether there are significant differences

on yield among the 4 varieties and draw a plot to

show the results as cool as you can

res_aov = aov(Yield~Variety, data = Yield_data)
summary(res_aov)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## Variety      3  239.9   79.96   9.117 0.000942 ***
## Residuals   16  140.3    8.77                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(res_aov)
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Yield ~ Variety, data = Yield_data)
## 
## $Variety
##          diff          lwr       upr     p adj
## B-A -3.648535  -9.00728079  1.710210 0.2480842
## C-A -7.196999 -12.55574453 -1.838253 0.0070606
## D-A  1.809932  -3.54881375  7.168677 0.7700310
## C-B -3.548464  -8.90720930  1.810282 0.2690422
## D-B  5.458467   0.09972148 10.817213 0.0451330
## D-C  9.006931   3.64818521 14.365676 0.0009975
plot(TukeyHSD(res_aov))

library(ggplot2)
ggplot(Yield_data,
       aes(x=Variety,y=Yield), show.legend = T)+
  geom_boxplot(aes(fill = Variety))+
  scale_fill_manual (values = cm.colors(4))+
  ylim(25, 60)+
  labs (title = "The differences among varieties") +
  theme(plot.title = element_text(face = "italic", color = "red", size = 14))+
  xlab ("Varieties")+
  ylab ("Yield (Ton)")+
  theme(axis.title = element_text(face = "bold", color = "blue", size =14))+
  scale_x_discrete(labels = c("Type A", "Type B", "Type C", "Type D")) +
  theme(axis.text = element_text(face = "bold",
                                 color = "Black", size = 12, angle = 0,
                                 hjust = 1))+
  geom_segment(x=3,
               y=55,
               xend = 3,
               yend = 58)+
  geom_segment(x=4,
               y=58,
               xend = 4,
               yend = 59)+
  geom_segment(x=3,
               y=59,
               xend = 4,
               yend = 59)+
  annotate(geom = "text",
           x=3.5,
           y=59.5,
           label = "*")+
  annotate(geom = "text",
           x=2.5,
           y=57,
           label = "P=2.13e-05")