library(readr)
data1 <- read_csv("skulls.csv")
## Rows: 150 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): epoch
## dbl (5): rownames, mb, bh, bl, nh
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
head(data1)
## # A tibble: 6 × 6
## rownames epoch mb bh bl nh
## <dbl> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 1 c4000BC 131 138 89 49
## 2 2 c4000BC 125 131 92 48
## 3 3 c4000BC 131 132 99 50
## 4 4 c4000BC 119 132 96 44
## 5 5 c4000BC 136 143 100 54
## 6 6 c4000BC 138 137 89 56
Pada tes one way ANOVA digunakan data mengenai pengukuran tengkorak mesir dari lima zaman yang berbeda, yaitu tahun 4000, 3300, 1850, 200 sebelum masehi dan 150 tahun masehi. Terdapat 150 observasi dan 5 variabel, yaitu epoch (tengkorang mesir dari berbagai zaman), mb (napas maksimum tengkorak), bh (ketinggian basibergmatic tengkorak), bl (panjang basliveolar tengkorak), dan nh (ketinggian nasal tengkorak). Dalam tes one way ANOVA ini ingin dilihat apakah tengkorang mesir dari lima zaman yang berbeda memiliki ukuran-ukuran yang berbeda.
library(MVN)
## Warning: package 'MVN' was built under R version 4.4.1
data1_fix <- data1[,3:6]
head(data1_fix)
## # A tibble: 6 × 4
## mb bh bl nh
## <dbl> <dbl> <dbl> <dbl>
## 1 131 138 89 49
## 2 125 131 92 48
## 3 131 132 99 50
## 4 119 132 96 44
## 5 136 143 100 54
## 6 138 137 89 56
test = mvn(data1_fix, mvnTest = "mardia", univariateTest = "SW", multivariatePlot = "qq")
test
## $multivariateNormality
## Test Statistic p value Result
## 1 Mardia Skewness 25.2679738909465 0.191370740880957 YES
## 2 Mardia Kurtosis -0.0433467219770312 0.965425147447826 YES
## 3 MVN <NA> <NA> YES
##
## $univariateNormality
## Test Variable Statistic p value Normality
## 1 Shapiro-Wilk mb 0.9913 0.4925 YES
## 2 Shapiro-Wilk bh 0.9864 0.1503 YES
## 3 Shapiro-Wilk bl 0.9940 0.7887 YES
## 4 Shapiro-Wilk nh 0.9840 0.0804 YES
##
## $Descriptives
## n Mean Std.Dev Median Min Max 25th 75th Skew Kurtosis
## mb 150 133.97333 4.890680 134 119 148 131 137 -0.02808971 0.2162337
## bh 150 132.54667 4.939346 133 120 145 129 136 -0.17427582 -0.1485636
## bl 150 96.46000 5.377844 96 81 114 93 100 0.13808882 0.1432677
## nh 150 50.93333 3.207932 51 44 60 49 53 0.08172942 -0.2068898
Untuk melakukan uji One Way ANOVA, perlu dilakukan pengecekan terlebih dahulu apakah data yang digunakan telah berdistribusi normal multivariat. Maka, dilakukan uji normalitas multivariat dengan mardia test. Dengan Hipotesis,
H0 : Data berdistribusi normal multivariat
H1 : Data tidak berdistribusi normal multivariat
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| Mardia | P-Value | Alpha | Keputusan |
|---|---|---|---|
| Skewness | 0,19137 | 0,05 | Terima H0 |
| Kurtosis | 0,96542 | 0,05 | Terima H0 |
Data berdistribusi normal multivariat apabila P-Value > alpha (0,05). Dari hasil p-value di atas, dapat disimpulkan bahwa kedua nilai p-value > alpha (0,05), maka Terima H0. Sehingga, dapat disimpulkan bahwa data berdistribusi normal multivariat.
library(biotools)
## Warning: package 'biotools' was built under R version 4.4.1
## Loading required package: MASS
## ---
## biotools version 4.2
grup <- data1$epoch
head(grup)
## [1] "c4000BC" "c4000BC" "c4000BC" "c4000BC" "c4000BC" "c4000BC"
boxM(data = data1_fix, grouping = grup)
##
## Box's M-test for Homogeneity of Covariance Matrices
##
## data: data1_fix
## Chi-Sq (approx.) = 45.667, df = 40, p-value = 0.2483
Asumsi selanjutnya yang harus terpenuhi pada data adalah homogenitas. Maka, perlu dilakukan pengecekan lagi terhadap kehomegenitas data. Dengan hipotesis,
H0 : Matriks kovarians grup adalah sama
H1 : Minimal ada satu matriks kovarians grup yang berbeda
Kriteria Uji:
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| P-Value BoxM Test | 0,2483 | 0,05 | Terima H0 |
Data memenuhi kehomogenitas apabila P-Value > alpha (0,05). Dari hasil p-value di atas, dapat disimpulkan bahwa kedua nilai p-value > alpha (0,05), maka Terima H0. Sehingga, dapat disimpulkan bahwa data memenuhi kehomgenitasan.
owm = manova(cbind(data1$mb, data1$bh, data1$bl, data1$nh)~data1$epoch)
summary(owm)
## Df Pillai approx F num Df den Df Pr(>F)
## data1$epoch 4 0.35331 3.512 16 580 4.675e-06 ***
## Residuals 145
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Karena sudah memenuhi asumsi normal dan homogenitas, maka bisa dilakukan tes one way ANOVA. Dengan hipotesis sebagai berikut:
H0 : Faktor epoch tidak memberikan hasil pengukuran tengkorak mesir dari lima zaman berbeda
H1 : Setidaknya ada satu vektor yang dipengaruhi signifikan oleh faktor epoch
Kriteria Uji:
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| P-Value One Way ANOVA | 4,675e-06 | 0,05 | Tolak H0 |
summary.aov(owm)
## Response 1 :
## Df Sum Sq Mean Sq F value Pr(>F)
## data1$epoch 4 502.83 125.707 5.9546 0.0001826 ***
## Residuals 145 3061.07 21.111
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Response 2 :
## Df Sum Sq Mean Sq F value Pr(>F)
## data1$epoch 4 229.9 57.477 2.4474 0.04897 *
## Residuals 145 3405.3 23.485
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Response 3 :
## Df Sum Sq Mean Sq F value Pr(>F)
## data1$epoch 4 803.3 200.823 8.3057 4.636e-06 ***
## Residuals 145 3506.0 24.179
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Response 4 :
## Df Sum Sq Mean Sq F value Pr(>F)
## data1$epoch 4 61.2 15.300 1.507 0.2032
## Residuals 145 1472.1 10.153
Karena hasil tes one way MANOVA menunjukan bahwa faktor epoch memberikan hasil pengukuran tengkorak mesir dari lima zaman yang berbeda. Maka, ingin diketahui lebih lanjut faktor epoch memberikan hasil yang berbeda terhadap variabel vektor yang mana. Dilakukan uji lanjut Post-Hoc Test dengan hipotesis sebagai berikut:
H0 : Faktor epoch tidak memberikan hasil pengukuran tengkorak mesir dari lima zaman berbeda pada vektor yang diamati
H1 : Faktor epoch memberikan hasil pengukuran tengkorak mesir dari lima zaman yang berbeda pada vektor yang diamati
Kriteria Uji:
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| Response | P-Value | alpha | Keputusan |
|---|---|---|---|
| 1. Nafas Maksimum | 0,0001826 | 0,05 | Tolak H0 |
| 2. Ketinggian Basibregmatic | 0,04897 | 0,05 | Tolak H0 |
| 3. Panjang Basialiveolar | 4,636e-06 | 0,05 | Tolak H0 |
| 4. Ketinggian Nasal | 0,2032 | 0,05 | Terima H0 |
library(readr)
data2 <- read_csv("Salaries.csv")
## Rows: 397 Columns: 7
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (3): rank, discipline, sex
## dbl (4): rownames, yrs.since.phd, yrs.service, salary
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
head(data2)
## # A tibble: 6 × 7
## rownames rank discipline yrs.since.phd yrs.service sex salary
## <dbl> <chr> <chr> <dbl> <dbl> <chr> <dbl>
## 1 1 Prof B 19 18 Male 139750
## 2 2 Prof B 20 16 Male 173200
## 3 3 AsstProf B 4 3 Male 79750
## 4 4 Prof B 45 39 Male 115000
## 5 5 Prof B 40 41 Male 141500
## 6 6 AssocProf B 6 6 Male 97000
Pada tes one way ANOVA digunakan data mengenai professor, professor madya, dan asisten professor dari sebuah kampus di Amerika Serikat pada tahun 2008-2009. Terdapat 397 observasi dan 7 variabel, yaitu rank (tingkatan profesor), discipline (bagian mengajar), yrs.since.phd (tahun semenjak Phd), yrs.of.service (masa kerja), sex (gender), dan salary (gaji). Dalam tes two way ANOVA ini ingin dilihat apakah tingkatan professor dan bagian mengajar dari sebuah kampus di Amerika Serikat pada tahun 2008-2009 berpengaruh secara signifikan terhadap tahun semenjak Phd, masa kerja, dan gaji.
## Uji Normalitas Multivariate
x1 <- data2[,4]
x2 <- data2[,5]
x3 <- data2[,7]
data2_fix <- data.frame(x1=x1, x2=x2, x3=x3)
library(MVN)
test = mvn(data2_fix, mvnTest = "mardia", univariateTest = "SW", multivariatePlot = "qq")
test
## $multivariateNormality
## Test Statistic p value Result
## 1 Mardia Skewness 292.53660804796 5.87259177177525e-57 NO
## 2 Mardia Kurtosis 5.93488545071384 2.94050717108973e-09 NO
## 3 MVN <NA> <NA> NO
##
## $univariateNormality
## Test Variable Statistic p value Normality
## 1 Shapiro-Wilk yrs.since.phd 0.9696 <0.001 NO
## 2 Shapiro-Wilk yrs.service 0.9418 <0.001 NO
## 3 Shapiro-Wilk salary 0.9599 <0.001 NO
##
## $Descriptives
## n Mean Std.Dev Median Min Max 25th 75th
## yrs.since.phd 397 22.31486 12.88700 21 1 56 12 32
## yrs.service 397 17.61461 13.00602 16 0 60 7 27
## salary 397 113706.45844 30289.03869 107300 57800 231545 91000 134185
## Skew Kurtosis
## yrs.since.phd 0.2986097 -0.8105942
## yrs.service 0.6456616 -0.3360614
## salary 0.7091778 0.1814838
Untuk melakukan uji Two Way MANOVA, perlu dilakukan pengecekan terlebih dahulu apakah data yang digunakan telah berdistribusi normal multivariat. Maka, dilakukan uji normalitas multivariat dengan mardia test. Dengan Hipotesis,
H0 : Data berdistribusi normal multivariat
H1 : Data tidak berdistribusi normal multivariat
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| Mardia | P-Value | Alpha | Keputusan |
|---|---|---|---|
| Skewness | 5,87259 | 0,05 | Tolak H0 |
| Kurtosis | 2,94051 | 0,05 | Tolak H0 |
Data berdistribusi normal multivariat apabila Terima H0,P-Value > alpha (0,05).
Dari hasil p-value di atas, dapat disimpulkan bahwa kedua nilai p-value < alpha (0,05), maka Tolak H0. Sehingga, dapat disimpulkan bahwa data tidak berdistribusi normal multivariat. Dengan tujuan pembelajaran, akan dilakukan pemotongan terhadap data yang awalnya 397 menjadi 30 (85:114) agar data menjadi berdistribusi normal multivariat.
## Potong Data
data2_p <- data2_fix[85:114,]
dim(data2_p)
## [1] 30 3
x1_p <- data2[85:114,4]
x2_p <- data2[85:114,5]
x3_p <- data2[85:114,7]
## Uji Normalitas Lagi
library(MVN)
test = mvn(data2_p, mvnTest = "mardia", univariateTest = "SW", multivariatePlot = "qq")
test
## $multivariateNormality
## Test Statistic p value Result
## 1 Mardia Skewness 18.1303365648699 0.052804447146691 YES
## 2 Mardia Kurtosis 0.453356519733072 0.650292040267642 YES
## 3 MVN <NA> <NA> YES
##
## $univariateNormality
## Test Variable Statistic p value Normality
## 1 Shapiro-Wilk yrs.since.phd 0.9302 0.0496 NO
## 2 Shapiro-Wilk yrs.service 0.9258 0.0381 NO
## 3 Shapiro-Wilk salary 0.9461 0.1332 YES
##
## $Descriptives
## n Mean Std.Dev Median Min Max 25th
## yrs.since.phd 30 20.30000 11.95437 17.5 2 42 10.25
## yrs.service 30 15.63333 11.08115 14.0 0 38 7.00
## salary 30 112631.73333 26428.28317 106848.5 72500 172272 90202.75
## 75th Skew Kurtosis
## yrs.since.phd 29.5 0.4133458 -1.0711037
## yrs.service 24.5 0.5740170 -0.8426250
## salary 129135.0 0.5074847 -0.6278946
Setelah dilakukan pemotongan data, dilakukan kembali uji normlaitas pada data. Dengan Hipotesis,
H0 : Data berdistribusi normal multivariat
H1 : Data tidak berdistribusi normal multivariat
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| Mardia | P-Value | Alpha | Keputusan |
|---|---|---|---|
| Skewness | 0,0528 | 0,05 | Terima H0 |
| Kurtosis | 0,65029 | 0,05 | Terima H0 |
Data berdistribusi normal multivariat apabila Terima H0,P-Value > alpha (0,05).
Dari hasil p-value di atas, dapat disimpulkan bahwa kedua nilai p-value > alpha (0,05), maka Tolak H0. Sehingga, dapat disimpulkan bahwa data tidak berdistribusi normal multivariat.
library(biotools)
grup1 <- ifelse(data2[85:114,2] == "Prof",1, ifelse(data2[85:114,2] == "AsstProf",2,3))
head(data2[85:114,2])
## # A tibble: 6 × 1
## rank
## <chr>
## 1 Prof
## 2 Prof
## 3 Prof
## 4 AsstProf
## 5 Prof
## 6 AssocProf
boxM(data = data2_p, grouping = grup1)
##
## Box's M-test for Homogeneity of Covariance Matrices
##
## data: data2_p
## Chi-Sq (approx.) = 26.216, df = 12, p-value = 0.01
grup2 <- ifelse(data2[85:114,3] == "Prof",1,2)
head(data2[85:114,3])
## # A tibble: 6 × 1
## discipline
## <chr>
## 1 B
## 2 B
## 3 B
## 4 B
## 5 B
## 6 B
boxM(data = data2_p, grouping = grup2)
##
## Box's M-test for Homogeneity of Covariance Matrices
##
## data: data2_p
## Chi-Sq (approx.) = NaN, df = 0, p-value = NA
Asumsi selanjutnya yang harus terpenuhi pada data adalah homogenitas. Maka, perlu dilakukan pengecekan lagi terhadap kehomegenitas data. Dengan hipotesis,
H0 : Matriks kovarians grup adalah sama
H1 : Minimal ada satu matriks kovarians grup yang berbeda
Kriteria Uji:
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| P-Value BoxM Test (Rank) | 0,01 | 0,05 | Tolak H0 |
| P-Value BoxM Test (Discipline) | 0,168 | 0,05 | Terima H0 |
Data memenuhi kehomogenitas apabila P-Value > alpha (0,05). Dari hasil p-value di atas, dapat disimpulkan bahwa grup rank nilai p-value < alpha (0,05), maka Tolak H0. Sehingga, dapat disimpulkan bahwa data grup rank tidak memenuhi kehomgenitasan. Dari hasil p-value di atas, dapat disimpulkan bahwa grup discipline nilai p-value > alpha (0,05), maka Terima H0. Sehingga, dapat disimpulkan bahwa data grup discipline memenuhi kehomgenitasan. Namun, dengan tujuan pembelajaran untuk grup rank diasumsikan telah memenuhi asumsi.
## Two Way MANOVA
rank <- as.factor(data2$rank[85:114])
discipline <- as.factor(data2$discipline[85:114])
x1_p<-unlist(x1_p)
x2_p<-unlist(x2_p)
x3_p<-unlist(x3_p)
manova <- manova(cbind(x1_p, x2_p, x3_p) ~ discipline * rank, data = data2_p)
# Menampilkan hasil
summary(manova)
## Df Pillai approx F num Df den Df Pr(>F)
## discipline 1 0.53385 8.3985 3 22 0.0006607 ***
## rank 2 0.75378 4.6372 6 46 0.0009205 ***
## discipline:rank 2 0.06562 0.2601 6 46 0.9525695
## Residuals 24
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Karena sudah memenuhi asumsi normal dan homogenitas, maka bisa dilakukan tes Two Way MANOVA. Dengan hipotesis sebagai berikut:
Rank:
H0 : Faktor rank tidak berpengaruh signifikan terhadap vektor
H1 : Setidaknya ada satu vektor yang dipengaruhi signifikan oleh faktor rank
Discipline:
H0 : Faktor discipline tidak berpengaruh signifikan terhadap vektor
H1 : Setidaknya ada satu vektor yang dipengaruhi signifikan oleh faktor discipline
Rank & Discipline:
H0 : Interaksi antara faktor rank & discipline tidak berpengaruh signifikan terhadap vektor
H1 : Setidaknya ada satu vektor yang dipengaruhi signifikan oleh interaksi antara faktor rank & discipline
Kriteria Uji:
Tolak H0 apabila P-value < alpha (0,05).
Melalui perhitungan syntax R di atas, didapatkan p-value sebagai berikut:
| P-Value Rank | 0,00067 | 0,05 | Tolak H0 |
| P-Value Discipline | 0,0009205 | 0,05 | Tolak H0 |
| P-Value Rank & Discipline | 0,95257 | 0,05 | Terima H0 |
summary.aov(manova)
## Response x1_p :
## Df Sum Sq Mean Sq F value Pr(>F)
## discipline 1 0.15 0.15 0.0022 0.9634
## rank 2 2432.88 1216.44 17.4615 2.085e-05 ***
## discipline:rank 2 39.32 19.66 0.2822 0.7566
## Residuals 24 1671.94 69.66
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Response x2_p :
## Df Sum Sq Mean Sq F value Pr(>F)
## discipline 1 6.67 6.67 0.1049 0.7488
## rank 2 1959.32 979.66 15.4212 4.933e-05 ***
## discipline:rank 2 70.33 35.17 0.5535 0.5821
## Residuals 24 1524.65 63.53
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Response x3_p :
## Df Sum Sq Mean Sq F value Pr(>F)
## discipline 1 4871472763 4871472763 19.9437 0.0001616 ***
## rank 2 9513287357 4756643678 19.4736 9.437e-06 ***
## discipline:rank 2 8133878 4066939 0.0166 0.9834993
## Residuals 24 5862276395 244261516
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Kriteria Uji:
Tolak H0 apabila P-value < alpha (0,05).
| Response | Group | P-Value | alpha | Keputusan |
|---|---|---|---|---|
| Rank | Years since Phd. | 2,085e-05 | 0,05 | Tolak H0 |
| Rank | Years of service | 4,933e-05 | 0,05 | Tolak H0 |
| Rank | Salary | 9,437e-06 | 0,05 | Tolak H0 |
| Discipline | Years since Phd. | 0,9643 | 0,05 | Terima H0 |
| Discipline | Years of service | 0,7488 | 0,05 | Terima H0 |
| Discipline | Salary | 0,0001616 | 0,05 | Tolak H0 |
| Rank & Discipline | Years since Phd. | 0,7566 | 0,05 | Terima H0 |
| Rank & Discipline | Years of service | 0,5821 | 0,05 | Terima H0 |
| Rank & Discipline | Salary | 0,9835 | 0,05 | Terima H0 |
Kesimpulan:
Faktor rank berpengaruh signifikan terhadap tahun semenjak Phd
Faktor rank berpengaruh signifikan terhadap masa kerja
Faktor rank berpengaruh signifikan terhadap gaji
Faktor discipline tidak berpengaruh signifikan terhadap tahun semenjak Phd
Faktor discipline tidak berpengaruh signifikan terhadap masa kerja
Faktor discipline berpengaruh signifikan terhadap gaji
Faktor interaksi antara rank dan discipline tidak berpengaruh signifikan terhadap tahun semenjak Phd
Faktor interaksi rank dan discipline tidak berpengaruh signifikan terhadap masa kerja
Faktor interaksi antara rank dan discipline tidak berpengaruh signifikan terhadap gaji