#PACKAGES
library(car)
## Warning: package 'car' was built under R version 4.5.3
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.5.3
library(agricolae)
## Warning: package 'agricolae' was built under R version 4.5.3
#Memanggil Data
data = read.table("dataTMKOMSTAT.txt", header = TRUE)
data
##   MerekTablet JenisPelarut  Data
## 1           1            1 46.18
## 2           1            2 54.54
## 3           1            3 56.63
## 4           2            1 57.69
## 5           2            2 71.57
## 6           2            3 81.51
## 7           3            1 48.61
## 8           3            2 51.01
## 9           3            3 80.02
#Pendefinisian Variabel
data$MerekTablet=factor(data$MerekTablet,levels=1:3,labels=c("CDR","Redoxon","Enervon-C"))
data$JenisPelarut=factor(data$JenisPelarut,levels=1:3,labels=c("AirBiasa","AirSoda","AirKelapa"))
data
##   MerekTablet JenisPelarut  Data
## 1         CDR     AirBiasa 46.18
## 2         CDR      AirSoda 54.54
## 3         CDR    AirKelapa 56.63
## 4     Redoxon     AirBiasa 57.69
## 5     Redoxon      AirSoda 71.57
## 6     Redoxon    AirKelapa 81.51
## 7   Enervon-C     AirBiasa 48.61
## 8   Enervon-C      AirSoda 51.01
## 9   Enervon-C    AirKelapa 80.02
#Normalitas per Merek Tablet
for(i in levels(data$MerekTablet)) {
  cat("\n========================================\n")
  cat(paste("Uji Normalitas Merek Tablet:", i), "\n")
  cat("========================================\n")
  print(shapiro.test(data$Data[data$MerekTablet == i]))
}
## 
## ========================================
## Uji Normalitas Merek Tablet: CDR 
## ========================================
## 
##  Shapiro-Wilk normality test
## 
## data:  data$Data[data$MerekTablet == i]
## W = 0.89286, p-value = 0.3631
## 
## 
## ========================================
## Uji Normalitas Merek Tablet: Redoxon 
## ========================================
## 
##  Shapiro-Wilk normality test
## 
## data:  data$Data[data$MerekTablet == i]
## W = 0.99096, p-value = 0.8182
## 
## 
## ========================================
## Uji Normalitas Merek Tablet: Enervon-C 
## ========================================
## 
##  Shapiro-Wilk normality test
## 
## data:  data$Data[data$MerekTablet == i]
## W = 0.80695, p-value = 0.1312
#Normalitas per Jumlah Pelarut
for(j in levels(data$JenisPelarut)) {
  cat("\n========================================\n")
  cat(paste("Uji Normalitas Jenis Pelarut:", j), "\n")
  cat("========================================\n")
  print(shapiro.test(data$Data[data$JenisPelarut == j]))
}
## 
## ========================================
## Uji Normalitas Jenis Pelarut: AirBiasa 
## ========================================
## 
##  Shapiro-Wilk normality test
## 
## data:  data$Data[data$JenisPelarut == j]
## W = 0.89987, p-value = 0.3851
## 
## 
## ========================================
## Uji Normalitas Jenis Pelarut: AirSoda 
## ========================================
## 
##  Shapiro-Wilk normality test
## 
## data:  data$Data[data$JenisPelarut == j]
## W = 0.87434, p-value = 0.3079
## 
## 
## ========================================
## Uji Normalitas Jenis Pelarut: AirKelapa 
## ========================================
## 
##  Shapiro-Wilk normality test
## 
## data:  data$Data[data$JenisPelarut == j]
## W = 0.79474, p-value = 0.102
#Homogenitas variansi#
leveneTest(Data ~ MerekTablet, center = mean, data = data)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  2  2.4902 0.1632
##        6
leveneTest(Data ~ JenisPelarut, center = mean, data = data)
## Levene's Test for Homogeneity of Variance (center = mean)
##       Df F value Pr(>F)
## group  2  1.9806 0.2185
##        6
#Uji ANOVA#
Y <- matrix(data$Data, nrow = 3, byrow = TRUE,
            dimnames = list(Merek = c("CDR", "Redoxon", "Enervon-C"),
                            Pelarut = c("AirBiasa", "AirSoda", "AirKelapa")))

anava2arah_tanpa <- function(Y, alpha = 0.05){
  a <- nrow(Y)
  b <- ncol(Y)
  N <- a * b
  FK <- sum(Y)^2 / N
  JKT <- sum(Y^2) - FK
  JK_A <- sum(rowSums(Y)^2) / b - FK
  JK_B <- sum(colSums(Y)^2) / a - FK
  JKG <- JKT - JK_A - JK_B
  
  db <- c(a - 1, b - 1, (a - 1) * (b - 1), N - 1)
  JK <- c(JK_A, JK_B, JKG, JKT)
  KT <- JK / db
  F_hit <- c(KT[1:2] / KT[3], NA, NA)
  F_tabel <- c(qf(1 - alpha, db[1:2], db[3]), NA, NA)
  p_value <- c(pf(F_hit[1:2], db[1:2], db[3], lower.tail = FALSE), NA, NA)
  
  hasil <- data.frame(
    Sumber = c("Faktor Baris (MerekTablet)", "Faktor Kolom (JenisPelarut)", "Galat", "Total"),
    db, JK, KT, F_hit, F_tabel, p_value
  )
  print(hasil, digits = 5, na.print = "", row.names = FALSE)
}

anava2arah_tanpa(Y, alpha = 0.05)
##                       Sumber db      JK      KT  F_hit F_tabel  p_value
##   Faktor Baris (MerekTablet)  2  479.96 239.979 4.2696  6.9443 0.101762
##  Faktor Kolom (JenisPelarut)  2  733.92 366.960 6.5287  6.9443 0.054991
##                        Galat  4  224.83  56.207     NA      NA       NA
##                        Total  8 1438.70 179.838     NA      NA       NA
uji_anava = aov(Data ~ MerekTablet + JenisPelarut, data = data)
anova(uji_anava)
## Analysis of Variance Table
## 
## Response: Data
##              Df Sum Sq Mean Sq F value  Pr(>F)  
## MerekTablet   2 479.96  239.98  4.2696 0.10176  
## JenisPelarut  2 733.92  366.96  6.5287 0.05499 .
## Residuals     4 224.83   56.21                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1