library(readxl)
CRD <- read_excel("D:/MARV BS MATH/Marv 4th year, 1st sem/Regression Analysis/CRDEqual.xlsx", 
    col_types = c("numeric", "numeric"))
View(CRD)

Fitting of linear model

model <-lm(CRD$Yield ~ CRD$Treatment)

Obtains R Square and other statistics of fitted model

summary <-summary(model)

Null Hypothesis : All the five treatment means are the same.

Alternative Hypothesis: At least one of the treatment is different.

Carryout ANOVA

anova <-anova(model)
anova
Analysis of Variance Table

Response: CRD$Yield
              Df Sum Sq Mean Sq F value    Pr(>F)    
CRD$Treatment  1 220.90 220.900  19.335 0.0003477 ***
Residuals     18 205.65  11.425                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Interpretation

The treatment source has p-value less than 0.01, so it is significant at 1% level of significance. This means that one or more treatments means are unequal so we reject our null hypothesis. Now it’s a green signal for multiple mean comparison test like LSD (least significant difference).

Below codes are used to obtain plots of fitted vs Residuals and Normal QQ plots

par(mfrow=c(1,2))
plot(model, which=1)
plot(model, which=2)

Load the package

library(agricolae)
Warning: package 'agricolae' was built under R version 4.2.2

Carry out LSD test

LSD <-LSD.test(CRD$Yield,CRD$Treatment,anova$`Df`[2],anova$`Mean Sq`[2])
LSD
$statistics
  MSerror Df  Mean       CV  t.value      LSD
   11.425 18 19.85 17.02815 2.100922 5.021379

$parameters
        test p.ajusted        name.t ntr alpha
  Fisher-LSD      none CRD$Treatment   5  0.05

$means
  CRD$Yield      std r       LCL      UCL Min Max   Q25  Q50   Q75
1     21.25 2.872281 4 17.699349 24.80065  18  25 20.25 21.0 22.00
2     25.50 1.732051 4 21.949349 29.05065  24  28 24.75 25.0 25.75
3     21.25 3.774917 4 17.699349 24.80065  16  24 19.75 22.5 24.00
4     18.00 1.825742 4 14.449349 21.55065  16  20 16.75 18.0 19.25
5     13.25 1.707825 4  9.699349 16.80065  11  15 12.50 13.5 14.25

$comparison
NULL

$groups
  CRD$Yield groups
2     25.50      a
1     21.25     ab
3     21.25     ab
4     18.00     bc
5     13.25      c

attr(,"class")
[1] "group"

Treatment 2 has the highest mean which was at par with Treatment 1 and 3.

Generate the txt file of analysis

sink("crdanalysis.txt")
print("ANOVA of CRD")
## [1] "ANOVA of CRD"
print(anova)
## Analysis of Variance Table
## 
## Response: CRD$Yield
##               Df Sum Sq Mean Sq F value    Pr(>F)    
## CRD$Treatment  1 220.90 220.900  19.335 0.0003477 ***
## Residuals     18 205.65  11.425                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
print("LSD ANALYSIS")
[1] "LSD ANALYSIS"
print(LSD$statistics)
  MSerror Df  Mean       CV  t.value      LSD
   11.425 18 19.85 17.02815 2.100922 5.021379
print(LSD$groups)
  CRD$Yield groups
2     25.50      a
1     21.25     ab
3     21.25     ab
4     18.00     bc
5     13.25      c
sink()