How to conduct Least Significant difference (LSD) test using R STUDIO?
For the mean comparison among variables, Least Significant difference (LSD) test is the most common method. Today I’ll introduce LSD test using R Studio.
Here is one data.
Cultivar=rep(c("CV1"), each=12)
Nitrogen=rep(rep(c("N0","N1","N2","N3"), each=3))
Block=rep(c("I","II","III"),4)
Yield=c (99, 109, 89, 115, 142, 133, 121, 157, 142, 125, 150, 139)
dataA=data.frame(Cultivar, Nitrogen, Block, Yield)
dataA
Cultivar Nitrogen Block Yield
1 CV1 N0 I 99
2 CV1 N0 II 109
3 CV1 N0 III 89
4 CV1 N1 I 115
5 CV1 N1 II 142
6 CV1 N1 III 133
7 CV1 N2 I 121
8 CV1 N2 II 157
9 CV1 N2 III 142
10 CV1 N3 I 125
11 CV1 N3 II 150
12 CV1 N3 III 139This data is about the yield difference of CV1 in response to 4 different nitrogen fertilizer (N0 ,N1, N2, N3). First of all, let’s check the mean per each nitrogen fertilizer.
mean(subset(dataA, Nitrogen=="N0")$ Yield) #99
mean(subset(dataA, Nitrogen=="N1")$ Yield) #130
mean(subset(dataA, Nitrogen=="N2")$ Yield) #140
mean(subset(dataA, Nitrogen=="N3")$ Yield) #138It seems that yield is different from nitrogen fertilizers, but we need to confirm it statistically. First, I’ll run One-Way ANOVA to check there is a difference in final yield among nitrogen fertilizer through the below codes.
Fert_ANOVA= aov(Yield ~ Nitrogen + factor(Block), data=dataA)
summary(Fert_ANOVA)
Df Sum Sq Mean Sq F value Pr(>F)
Nitrogen 3 3248 1082.7 19.14 0.00179 **
factor(Block) 2 1207 603.3 10.66 0.01059 *
Residuals 6 339 56.6
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1The result shows that different nitrogen fertilizer affected final yield (p-value: 0.00179), indicating that yield is different from N0 to N3. Now I’d like to know which nitrogen fertilizer brings about the greatest yield.
To compare the mean of yield per each nitrogen fertilizer, let’s install agricolae package and run LSD test through the below codes.
install.packages("agricolae")
library(agricolae)
LSD_Test= LSD.test(Fert_ANOVA,"Nitrogen")
LSD_Test
$statistics
MSerror Df Mean CV t.value LSD
56.58333 6 126.75 5.934666 2.446912 15.02855
$parameters
test p.ajusted name.t ntr alpha
Fisher-LSD none Nitrogen 4 0.05
$means
Yield std r LCL UCL Min Max Q25 Q50 Q75
N0 99 10.00000 3 88.37321 109.6268 89 109 94.0 99 104.0
N1 130 13.74773 3 119.37321 140.6268 115 142 124.0 133 137.5
N2 140 18.08314 3 129.37321 150.6268 121 157 131.5 142 149.5
N3 138 12.52996 3 127.37321 148.6268 125 150 132.0 139 144.5
$comparison
NULL
$groups
Yield groups
N2 140 a
N3 138 a
N1 130 a
N0 99 bThe LSD Test indicates that yield was not different from N1, N2, N3 fertilizer applications, while N0 application resulted in the lowest yield (If the letter is the same, it means the mean is the same). If N0 is 0kg nitrogen application, it makes sense that yield is lower than under nitrogen applications.
What if Two-Way ANOVA? If we want to compare the mean about interaction between two variables, how can we do that?
Here is one data.
Cultivar=rep(c("CV1","CV2"),each=12)
Nitrogen=rep(rep(c("N0","N1","N2","N3"), each=3),2)
Block=rep(c("I","II","III"),8)
Yield=c (99, 109, 89, 115, 142, 133, 121, 157, 142, 125, 150, 139, 82, 104, 99, 117, 125, 127, 145, 154, 154, 151, 166, 175)
dataA=data.frame(Cultivar,Nitrogen,Block,Yield)
dataA
Cultivar Nitrogen Block Yield
1 CV1 N0 I 99
2 CV1 N0 II 109
3 CV1 N0 III 89
4 CV1 N1 I 115
5 CV1 N1 II 142
6 CV1 N1 III 133
7 CV1 N2 I 121
8 CV1 N2 II 157
9 CV1 N2 III 142
10 CV1 N3 I 125
11 CV1 N3 II 150
12 CV1 N3 III 139
13 CV2 N0 I 82
14 CV2 N0 II 104
15 CV2 N0 III 99
16 CV2 N1 I 117
17 CV2 N1 II 125
18 CV2 N1 III 127
19 CV2 N2 I 145
20 CV2 N2 II 154
21 CV2 N2 III 154
22 CV2 N3 I 151
23 CV2 N3 II 166
24 CV2 N3 III 175This data is a factorial experiment between two genotype (CV1, CV2) and 4 different nitrogen fertilizer (N0, N1, N2, N3), and final yield was measured. First, I’ll run Two-Way ANOVA to check there is an interaction between genotype and nitrogen fertilizer about yield through the below codes.
Fert_ANOVA= aov(Yield~Cultivar+Nitrogen+Cultivar:Nitrogen+factor(Block), data=dataA)
summary(Fert_ANOVA)
Df Sum Sq Mean Sq F value Pr(>F)
Cultivar 1 253 253 4.99 0.042327 *
Nitrogen 3 10695 3565 70.17 1.12e-08 ***
factor(Block) 2 1505 752 14.81 0.000351 ***
Cultivar:Nitrogen 3 1040 347 6.82 0.004604 **
Residuals 14 711 51
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1There is an interaction between genotype and nitrogen fertilizer about final yield (p-value: 0.004604). Now, I’ll compare the means among the interactions through the below codes.
LSD_Test= LSD.test(Fert_ANOVA, c("Cultivar", "Nitrogen"))
LSD_Test
test p.ajusted name. t ntr alpha
Fisher-LSD none Cultivar:Nitrogen 8 0.05
$means
Yield std r LCL UCL Min Max Q25 Q50 Q75
CV1:N0 99 10.000000 3 90.17386 107.8261 89 109 94.0 99 104.0
CV1:N1 130 13.747727 3 121.17386 138.8261 115 142 124.0 133 137.5
CV1:N2 140 18.083141 3 131.17386 148.8261 121 157 131.5 142 149.5
CV1:N3 138 12.529964 3 129.17386 146.8261 125 150 132.0 139 144.5
CV2:N0 95 11.532563 3 86.17386 103.8261 82 104 90.5 99 101.5
CV2:N1 123 5.291503 3 114.17386 131.8261 117 127 121.0 125 126.0
CV2:N2 151 5.196152 3 142.17386 159.8261 145 154 149.5 154 154.0
CV2:N3 164 12.124356 3 155.17386 172.8261 151 175 158.5 166 170.5
$comparison
NULL
$groups
Yield groups
CV2:N3 164 a
CV2:N2 151 b
CV1:N2 140 bc
CV1:N3 138 c
CV1:N1 130 cd
CV2:N1 123 d
CV1:N0 99 e
CV2:N0 95 eUnder the N0, yield was the lowest in both CV1 and CV2, and yield was the greatest in CV2 when N3 was applied.
full code: https://github.com/agronomy4future/r_code/blob/main/How_to_conduct_Least_Significant_difference_(LSD)_test_using_R_STUDIO.ipynb
