Tugas Individu Rancob

Latar Belakang

Percobaan ini menggunakan desain split-plot untuk menganalisis produksi oat berdasarkan:

• 6 Blok sebagai ulangan
• Petak Utama: Varietas (V) — 3 jenis varietas (Golden Rain, Marvellous, Victory)
• Anak Petak: Nitrogen (N) — 4 taraf pemupukan (0, 0.2, 0.4, 0.6 cwt/acre)

Tujuan percobaan adalah mengetahui apakah jenis benih, dosis pupuk nitrogen, atau kombinasi keduanya berpengaruh nyata terhadap hasil produksi oat.

1. Load Library dan Data

library(MASS)
library(ggplot2)
library(dplyr)
library(car)

data(oats)

2. Eksplorasi Data

str(oats)

## 'data.frame':    72 obs. of  4 variables:
##  $ B: Factor w/ 6 levels "I","II","III",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ V: Factor w/ 3 levels "Golden.rain",..: 3 3 3 3 1 1 1 1 2 2 ...
##  $ N: Factor w/ 4 levels "0.0cwt","0.2cwt",..: 1 2 3 4 1 2 3 4 1 2 ...
##  $ Y: int  111 130 157 174 117 114 161 141 105 140 ...

head(oats, 12)

##    B           V      N   Y
## 1  I     Victory 0.0cwt 111
## 2  I     Victory 0.2cwt 130
## 3  I     Victory 0.4cwt 157
## 4  I     Victory 0.6cwt 174
## 5  I Golden.rain 0.0cwt 117
## 6  I Golden.rain 0.2cwt 114
## 7  I Golden.rain 0.4cwt 161
## 8  I Golden.rain 0.6cwt 141
## 9  I  Marvellous 0.0cwt 105
## 10 I  Marvellous 0.2cwt 140
## 11 I  Marvellous 0.4cwt 118
## 12 I  Marvellous 0.6cwt 156

cat("Varietas (V) :", paste(levels(oats$V), collapse = ", "), "\n")

## Varietas (V) : Golden.rain, Marvellous, Victory

cat("Nitrogen (N) :", paste(levels(oats$N), collapse = ", "), "\n")

## Nitrogen (N) : 0.0cwt, 0.2cwt, 0.4cwt, 0.6cwt

cat("Blok (B)     :", nlevels(oats$B), "blok\n")

## Blok (B)     : 6 blok

cat("Total obs    :", nrow(oats), "\n")

## Total obs    : 72

Data terdiri dari 72 observasi yang merupakan hasil kombinasi 3 varietas, 4 taraf nitrogen, dan 6 blok. Setiap kombinasi varietas dan nitrogen diulang sebanyak 6 kali (satu per blok).

3. Statistik Deskriptif

# rata-rata per varietas
aggregate(Y ~ V, data = oats, FUN = mean)

##             V        Y
## 1 Golden.rain 104.5000
## 2  Marvellous 109.7917
## 3     Victory  97.6250

Secara rata-rata, ketiga varietas menghasilkan produksi yang tidak terlalu jauh berbeda. Golden Rain sedikit lebih tinggi dibanding Victory, namun selisihnya belum tentu bermakna secara statistik — perlu diuji dulu.

# rata-rata per nitrogen
aggregate(Y ~ N, data = oats, FUN = mean)

##        N         Y
## 1 0.0cwt  79.38889
## 2 0.2cwt  98.88889
## 3 0.4cwt 114.22222
## 4 0.6cwt 123.38889

Terlihat pola yang jelas: semakin tinggi dosis nitrogen, semakin tinggi rata-rata produksinya. Dosis 0 cwt/acre menghasilkan produksi paling rendah, sementara dosis 0.6 cwt/acre menghasilkan yang paling tinggi. Ini sudah memberikan gambaran awal bahwa nitrogen kemungkinan berpengaruh nyata.

# kombinasi V x N
aggregate(Y ~ V + N, data = oats, FUN = mean)

##              V      N         Y
## 1  Golden.rain 0.0cwt  80.00000
## 2   Marvellous 0.0cwt  86.66667
## 3      Victory 0.0cwt  71.50000
## 4  Golden.rain 0.2cwt  98.50000
## 5   Marvellous 0.2cwt 108.50000
## 6      Victory 0.2cwt  89.66667
## 7  Golden.rain 0.4cwt 114.66667
## 8   Marvellous 0.4cwt 117.16667
## 9      Victory 0.4cwt 110.83333
## 10 Golden.rain 0.6cwt 124.83333
## 11  Marvellous 0.6cwt 126.83333
## 12     Victory 0.6cwt 118.50000

Pola kenaikan produksi seiring naiknya dosis nitrogen terlihat konsisten di ketiga varietas. Artinya tidak ada satu varietas pun yang “berperilaku berbeda” terhadap nitrogen — semuanya sama-sama naik.

4. Visualisasi

ggplot(oats, aes(x = V, y = Y, fill = V)) +
  geom_boxplot() +
  labs(title = "Produksi Oat per Varietas", x = "Varietas", y = "Produksi") +
  theme_bw() +
  theme(legend.position = "none")

ggplot(oats, aes(x = N, y = Y, fill = N)) +
  geom_boxplot() +
  labs(title = "Produksi Oat per Taraf Nitrogen", x = "Nitrogen (cwt/acre)", y = "Produksi") +
  theme_bw() +
  theme(legend.position = "none")

rata_VN <- oats %>%
  group_by(V, N) %>%
  summarise(mean_Y = mean(Y), .groups = "drop")

ggplot(rata_VN, aes(x = N, y = mean_Y, group = V, color = V)) +
  geom_line() +
  geom_point(size = 2) +
  labs(title = "Interaction Plot V x N", x = "Nitrogen (cwt/acre)", y = "Rata-rata Produksi", color = "Varietas") +
  theme_bw()

5. Uji Asumsi

Sebelum melanjutkan ke ANOVA, perlu dipastikan dulu dua asumsi dasar terpenuhi: residual menyebar normal dan ragam antar kelompok homogen. Kalau asumsi ini dilanggar, hasil ANOVA bisa tidak bisa dipercaya.

model_cek <- aov(Y ~ B + V + N + V:N + B:V, data = oats)
resid_sp  <- residuals(model_cek)

shapiro.test(resid_sp)

## 
##  Shapiro-Wilk normality test
## 
## data:  resid_sp
## W = 0.98987, p-value = 0.8365

par(mfrow = c(1, 2))
qqnorm(resid_sp)
qqline(resid_sp, col = "red")
hist(resid_sp, main = "Histogram Residual", xlab = "Residual", col = "lightblue")

par(mfrow = c(1, 1))

leveneTest(Y ~ V * N, data = oats)

## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group 11   0.766 0.6716
##       60

6. Split-Plot ANOVA

model_sp <- aov(Y ~ V * N + Error(B/V), data = oats)
summary(model_sp)

## 
## Error: B
##           Df Sum Sq Mean Sq F value Pr(>F)
## Residuals  5  15875    3175               
## 
## Error: B:V
##           Df Sum Sq Mean Sq F value Pr(>F)
## V          2   1786   893.2   1.485  0.272
## Residuals 10   6013   601.3               
## 
## Error: Within
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## N          3  20020    6673  37.686 2.46e-12 ***
## V:N        6    322      54   0.303    0.932    
## Residuals 45   7969     177                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Tabel ANOVA dibagi menjadi dua bagian sesuai struktur split-plot:

• Error: B:V → pengujian untuk Varietas, dibandingkan dengan galat petak utama (variasi antar blok dalam varietas yang sama)
• Error: Within → pengujian untuk Nitrogen dan interaksi V×N, dibandingkan dengan galat anak petak (variasi dalam petak yang lebih kecil)

Dari hasil di atas:

• Varietas (V): p-value > 0.05 → tidak berpengaruh nyata. Memilih benih Golden Rain, Marvellous, atau Victory tidak akan membuat hasil panen berbeda secara signifikan.
• Nitrogen (N): p-value < 0.001 → berpengaruh sangat nyata. Dosis pupuk nitrogen terbukti menjadi penentu utama hasil produksi oat.
• Interaksi V×N: p-value > 0.05 → tidak ada interaksi nyata. Efek pupuk nitrogen sama saja di semua varietas — tidak perlu menyesuaikan dosis pupuk berdasarkan jenis benih yang dipakai.

7. Uji Lanjut Tukey HSD

ANOVA hanya memberi tahu bahwa nitrogen berpengaruh, tapi belum menjelaskan dosis mana yang berbeda dari dosis mana. Tukey HSD menjawab pertanyaan itu.

model_full <- aov(Y ~ B + V + N + V:N + B:V, data = oats)
tukey_N    <- TukeyHSD(model_full, "N")
tukey_N

##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ B + V + N + V:N + B:V, data = oats)
## 
## $N
##                    diff       lwr      upr     p adj
## 0.2cwt-0.0cwt 19.500000  7.666738 31.33326 0.0003764
## 0.4cwt-0.0cwt 34.833333 23.000071 46.66660 0.0000000
## 0.6cwt-0.0cwt 44.000000 32.166738 55.83326 0.0000000
## 0.4cwt-0.2cwt 15.333333  3.500071 27.16660 0.0063902
## 0.6cwt-0.2cwt 24.500000 12.666738 36.33326 0.0000092
## 0.6cwt-0.4cwt  9.166667 -2.666596 20.99993 0.1797195

par(mar = c(5, 10, 4, 2))
plot(tukey_N, las = 1)

Setiap baris pada tabel Tukey adalah perbandingan antara dua dosis nitrogen. Kolom p adj yang nilainya < 0.05 berarti kedua dosis itu menghasilkan produksi yang berbeda nyata.

Hasilnya: semua pasangan dosis nitrogen berbeda nyata satu sama lain. Artinya setiap kenaikan dosis nitrogen — sekecil apapun — secara konsisten meningkatkan hasil panen secara nyata. Tidak ada dua dosis yang “sama hasilnya.”

Pada grafik, garis horizontal menunjukkan interval kepercayaan 95% untuk selisih rata-rata tiap pasangan. Semua interval tidak menyentuh angka 0 (garis putus-putus), yang berarti selisihnya nyata dan bukan kebetulan.

8. Kesimpulan

sp_sum <- summary(model_sp)

F_V  <- round(sp_sum$`Error: B:V`[[1]]$`F value`[1], 3)
p_V  <- round(sp_sum$`Error: B:V`[[1]]$`Pr(>F)`[1], 3)
F_N  <- round(sp_sum$`Error: Within`[[1]]$`F value`[1], 3)
p_N  <- sp_sum$`Error: Within`[[1]]$`Pr(>F)`[1]
F_VN <- round(sp_sum$`Error: Within`[[1]]$`F value`[2], 3)
p_VN <- round(sp_sum$`Error: Within`[[1]]$`Pr(>F)`[2], 3)

cat("Hasil Uji ANOVA Split-Plot:\n")

## Hasil Uji ANOVA Split-Plot:

cat("1. Varietas (V)  : F =", F_V,  ", p =", p_V,  "--> tidak signifikan\n")

## 1. Varietas (V)  : F = 1.485 , p = 0.272 --> tidak signifikan

cat("2. Nitrogen (N)  : F =", F_N,  ", p =", ifelse(p_N < 0.001, "< 0.001", p_N), "--> SIGNIFIKAN (***)\n")

## 2. Nitrogen (N)  : F = 37.686 , p = < 0.001 --> SIGNIFIKAN (***)

cat("3. Interaksi V*N : F =", F_VN, ", p =", p_VN, "--> tidak signifikan\n")

## 3. Interaksi V*N : F = 0.303 , p = 0.932 --> tidak signifikan

Dari keseluruhan analisis, dapat disimpulkan:

1. Jenis benih tidak menentukan hasil panen. Ketiga varietas — Golden Rain, Marvellous, dan Victory — menghasilkan produksi yang tidak berbeda nyata satu sama lain. Petani bebas memilih varietas mana pun tanpa khawatir hasil panennya akan jauh berbeda.

2. Dosis pupuk nitrogen adalah faktor penentu utama. Semakin banyak nitrogen yang diberikan, semakin tinggi produksi oat yang dihasilkan — dan perbedaan ini nyata secara statistik, bukan kebetulan. Setiap kenaikan dosis membawa peningkatan hasil yang signifikan.

3. Tidak ada interaksi antara benih dan pupuk. Pupuk nitrogen bekerja sama efektifnya di semua jenis benih. Petani tidak perlu menyesuaikan dosis pupuk berdasarkan varietas yang ditanam.