Modul 2 Analisis Multivariat
Jeika Antama (011), Ayda’ Nur Salzabillah (017), Grace Pinkkan Ladyna(137)
2026-04-02
1. Persiapan Data
1.1 Load data & Library
Mengimpor data dan library yang dibutuhkan, lalu melihat tipe setiap variabel.
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.0 ✔ readr 2.1.6
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.2 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(dplyr)
library(tidyr)
library(ggplot2)
library(psych)##
## Attaching package: 'psych'
##
## The following objects are masked from 'package:ggplot2':
##
## %+%, alphalibrary(kableExtra)##
## Attaching package: 'kableExtra'
##
## The following object is masked from 'package:dplyr':
##
## group_rowslibrary(knitr)
library(car)## Loading required package: carData
##
## Attaching package: 'car'
##
## The following object is masked from 'package:psych':
##
## logit
##
## The following object is masked from 'package:dplyr':
##
## recode
##
## The following object is masked from 'package:purrr':
##
## somelibrary(biotools)## Warning: package 'biotools' was built under R version 4.5.3## Loading required package: MASS
##
## Attaching package: 'MASS'
##
## The following object is masked from 'package:dplyr':
##
## select
##
## ---
## biotools version 4.3library(MVN)## Warning: package 'MVN' was built under R version 4.5.3##
## Attaching package: 'MVN'
##
## The following object is masked from 'package:psych':
##
## mardiadata <- read.csv("hsb2.csv")
head(data)## id gender race ses schtyp prog read write math science socst
## 1 70 male white low public general 57 52 41 47 57
## 2 121 female white middle public vocational 68 59 53 63 61
## 3 86 male white high public general 44 33 54 58 31
## 4 141 male white high public vocational 63 44 47 53 56
## 5 172 male white middle public academic 47 52 57 53 61
## 6 113 male white middle public academic 44 52 51 63 61str(data)## 'data.frame': 200 obs. of 11 variables:
## $ id : int 70 121 86 141 172 113 50 11 84 48 ...
## $ gender : chr "male" "female" "male" "male" ...
## $ race : chr "white" "white" "white" "white" ...
## $ ses : chr "low" "middle" "high" "high" ...
## $ schtyp : chr "public" "public" "public" "public" ...
## $ prog : chr "general" "vocational" "general" "vocational" ...
## $ read : int 57 68 44 63 47 44 50 34 63 57 ...
## $ write : int 52 59 33 44 52 52 59 46 57 55 ...
## $ math : int 41 53 54 47 57 51 42 45 54 52 ...
## $ science: int 47 63 58 53 53 63 53 39 58 50 ...
## $ socst : int 57 61 31 56 61 61 61 36 51 51 ...1.2 Deskripsi Data
No. <- c(1:11)
Variabel <- c("id",
"gender",
"race",
"ses",
"schtyp",
"prog",
"read",
"write",
"math",
"science",
"socst")
Deskripsi <- c("ID Pelajar.",
"Jenis kelamin siswa (L/P)",
"Ras / etnis siswa",
"Status sosial ekonomi (low, middle, high)",
"Tipe sekolah (negeri/swasta)",
"Program pendidikan (general, academic, vocational)",
"Nilai membaca siswa",
"Nilai menulis siswa",
"Nilai matematika siswa",
"Nilai sains siswa",
"Nilai social studies (mata pelajaran sosial)")
data.frame(No., Variabel, Deskripsi) %>%
kbl() %>%
kable_styling(bootstrap_options = c("bordered", "stripped", "hover"),
full_width = T)| No. | Variabel | Deskripsi |
|---|---|---|
| 1 | id | ID Pelajar. |
| 2 | gender | Jenis kelamin siswa (L/P) |
| 3 | race | Ras / etnis siswa |
| 4 | ses | Status sosial ekonomi (low, middle, high) |
| 5 | schtyp | Tipe sekolah (negeri/swasta) |
| 6 | prog | Program pendidikan (general, academic, vocational) |
| 7 | read | Nilai membaca siswa |
| 8 | write | Nilai menulis siswa |
| 9 | math | Nilai matematika siswa |
| 10 | science | Nilai sains siswa |
| 11 | socst | Nilai social studies (mata pelajaran sosial) |
1.3 Menggunakan Variabel id
Variabel id digunakan sebagai identitas unik tiap observasi dan tetap dianalisis secara deskriptif.
head(data$id)## [1] 70 121 86 141 172 113summary(data$id)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.00 50.75 100.50 100.50 150.25 200.002. Karakteristik Data
2.1 Ubah Variabel Kategori
Mengubah semua variabel kategori menjadi factor agar bisa dianalisis sebagai kelompok.
data$gender <- as.factor(data$gender)
data$race <- as.factor(data$race)
data$ses <- as.factor(data$ses)
data$prog <- as.factor(data$prog)
data$schtyp <- as.factor(data$schtyp)2.2 Statistik Deskriptif Variabel Numerik
numeric_vars <- data[, c("id", "read", "write", "math", "science", "socst")]
describe(numeric_vars)## vars n mean sd median trimmed mad min max range skew kurtosis
## id 1 200 100.50 57.88 100.5 100.50 74.13 1 200 199 0.00 -1.22
## read 2 200 52.23 10.25 50.0 52.03 10.38 28 76 48 0.19 -0.66
## write 3 200 52.77 9.48 54.0 53.36 11.86 31 67 36 -0.47 -0.78
## math 4 200 52.65 9.37 52.0 52.23 10.38 33 75 42 0.28 -0.69
## science 5 200 51.85 9.90 53.0 52.02 11.86 26 74 48 -0.19 -0.60
## socst 6 200 52.41 10.74 52.0 52.99 13.34 26 71 45 -0.38 -0.57
## se
## id 4.09
## read 0.72
## write 0.67
## math 0.66
## science 0.70
## socst 0.76Berdasarkan statistik deskriptif yang diperoleh, variabel id hanya berfungsi sebagai identitas tiap siswa, sehingga tidak memberikan informasi analisis kecuali jumlah siswa yang tercatat yaitu 200. Nilai rata-rata untuk variabel akademik menunjukkan bahwa nilai read, write, math, science, dan socst cenderung berada di kisaran 50–53, dengan penyebaran nilai (standar deviasi) sekitar 9–11 poin, yang menandakan variasi kemampuan siswa yang relatif moderat. Median dan nilai trimmed hampir serupa dengan mean, menunjukkan tidak ada nilai ekstrem yang sangat memengaruhi distribusi data. Nilai skewness sebagian besar mendekati nol, yang menandakan distribusi data cukup simetris, meskipun variabel write menunjukkan sedikit skew negatif. Begitu pula, nilai kurtosis negatif pada sebagian besar variabel menandakan distribusi yang agak lebih datar dibandingkan distribusi normal. Secara keseluruhan, data akademik siswa relatif merata, dengan tidak ada outlier ekstrim yang signifikan, sehingga dapat dijadikan dasar analisis lanjutan
2.3 Statistik Deskriptif Variabel Kategori
library(dplyr)
# Fungsi untuk membuat tabel frekuensi + persen
freq_table <- function(x, var_name) {
data.frame(
Kategori = names(table(x)),
Frekuensi = as.vector(table(x))
) %>%
mutate(Persen = round(Frekuensi / sum(Frekuensi) * 100, 2)) %>%
rename(!!var_name := Kategori)
}
# Gender
freq_table(data$gender, "Gender")## Gender Frekuensi Persen
## 1 female 109 54.5
## 2 male 91 45.5# Race
freq_table(data$race, "Race")## Race Frekuensi Persen
## 1 african american 20 10.0
## 2 asian 11 5.5
## 3 hispanic 24 12.0
## 4 white 145 72.5# SES
freq_table(data$ses, "SES")## SES Frekuensi Persen
## 1 high 58 29.0
## 2 low 47 23.5
## 3 middle 95 47.5# Program
freq_table(data$prog, "Program")## Program Frekuensi Persen
## 1 academic 105 52.5
## 2 general 45 22.5
## 3 vocational 50 25.0# School Type
freq_table(data$schtyp, "School_Type")## School_Type Frekuensi Persen
## 1 private 32 16
## 2 public 168 84Berdasarkan distribusi kategori, terlihat bahwa gender siswa relatif seimbang, dengan 54,5% perempuan dan 45,5% laki-laki. Untuk ras/etnis, mayoritas siswa adalah white (72,5%), diikuti oleh hispanic (12%), african american (10%), dan asian (5,5%), menunjukkan dominasi satu kelompok etnis dalam sampel. Dari sisi status sosial ekonomi (SES), hampir setengah siswa berada di kategori middle (47,5%), dengan proporsi high (29%) dan low (23,5%) lebih kecil, menandakan variasi moderat dalam latar belakang ekonomi. Pada program belajar, lebih dari setengah siswa mengikuti program academic (52,5%), sementara vocational (25%) dan general (22,5%) menempati porsi lebih kecil. Terakhir, untuk jenis sekolah, mayoritas siswa berasal dari sekolah public (84%), sedangkan private (16%) hanya sebagian kecil, menggambarkan representasi sekolah publik yang dominan dalam sampel.
2.4 Cek Multikolinearitas
cor(data[, c("read","write","math","science","socst")])## read write math science socst
## read 1.0000000 0.5967765 0.6622801 0.6301579 0.6214843
## write 0.5967765 1.0000000 0.6174493 0.5704416 0.6047932
## math 0.6622801 0.6174493 1.0000000 0.6307332 0.5444803
## science 0.6301579 0.5704416 0.6307332 1.0000000 0.4651060
## socst 0.6214843 0.6047932 0.5444803 0.4651060 1.0000000Berdasarkan hasil analisis, korelasi antar variabel dependen yaitu read, write, math, dan science berada pada kisaran 0,57 hingga 0,66, sedangkan korelasi antara variabel kovariat (socst) dengan variabel dependen berkisar antara 0,54 hingga 0,62. Nilai korelasi tersebut tergolong sedang dan tidak terlalu tinggi karena masih berada di bawah 0,9. Hal ini menunjukkan bahwa tidak terdapat masalah multikolinearitas yang serius antar variabel. Dengan demikian, seluruh variabel tersebut masih layak dan dapat digunakan secara bersamaan dalam analisis MANOVA maupun MANCOVA.
2.5 Visualisasi
2.5.1 Histogram untuk Distribusi Variabel Numerik
numeric_vars <- c("read", "write", "math", "science", "socst")
# Ubah data ke format long agar bisa dipakai facet_wrap
data_long <- data %>%
dplyr::select(all_of(numeric_vars)) %>% # pakai dplyr::select untuk menghindari konflik
pivot_longer(cols = everything(), names_to = "Variable", values_to = "Value")
# Plot histogram semua variabel dalam 1 figure
ggplot(data_long, aes(x = Value, fill = Variable)) +
geom_histogram(color = "black", bins = 20, alpha = 0.6) +
facet_wrap(~ Variable, scales = "free") +
labs(title = "Distribusi Semua Variabel Numerik",
x = "Nilai", y = "Frekuensi") +
theme_minimal() +
theme(legend.position = "none")
Berdasarkan grafik distribusi histogram, nilai pada semua variabel numerik (math, read, science, socst, dan write) umumnya tersebar di kisaran 40–60 dan cenderung membentuk pola yang cukup simetris, meskipun ada sedikit perbedaan. Variabel read dan write tampak memiliki beberapa nilai tinggi, sedangkan math dan science relatif lebih merata. Sementara itu, socst terlihat sedikit lebih menyebar ke nilai yang lebih tinggi dibanding yang lain.
2.5.2 Boxplot untuk Melihat Sebaran & Outlier
ggplot(data_long, aes(x = Variable, y = Value, fill = Variable)) +
geom_boxplot(alpha = 0.6) +
labs(title = "Boxplot Semua Variabel Numerik",
x = "Variabel", y = "Nilai") +
theme_minimal() +
theme(legend.position = "none")
Dari boxplot, median semua variabel berada di sekitar angka 50–55, sehingga dapat dikatakan bahwa nilai tengah antar variabel tidak jauh berbeda. Selain itu, tidak terlihat adanya titik di luar whisker, yang berarti tidak terdapat outlier pada data. Secara keseluruhan, data menunjukkan pola yang cukup konsisten antar variabel dengan sebaran yang masih dalam batas wajar dan tidak ada penyimpangan ekstrem.
3. Uji Asumsi
3.1 Normalitas Multivariat
Melihat apakah gabungan semua variabel dependen mengikuti distribusi normal.
numeric_vars <- c("read", "write", "math", "science")
# Ekstrak data numerik
data_numeric <- data[, numeric_vars]
# Normalitas multivariat menggunakan Mardia
mardia_result <- psych::mardia(data_numeric)
# Lihat hasil
mardia_result## Call: psych::mardia(x = data_numeric)
##
## Mardia tests of multivariate skew and kurtosis
## Use describe(x) the to get univariate tests
## n.obs = 200 num.vars = 4
## b1p = 1.13 skew = 37.52 with probability <= 0.01
## small sample skew = 38.31 with probability <= 0.0081
## b2p = 23.16 kurtosis = -0.86 with probability <= 0.39Berdasarkan hasil Mardia, skewness multivariat (b1p) sebesar 1,13 dengan p < 0,01 menunjukkan ada sedikit kemiringan pada distribusi gabungan semua variabel dependen (read, write, math, science). Sementara kurtosis (b2p) sebesar 23,16 dengan p = 0,39 menunjukkan distribusi relatif normal, tidak terlalu “peaked” atau datar. Artinya, secara keseluruhan, distribusi gabungan data cukup mendekati normal meski ada skew minor. Dengan sampel 200, ini masih aman untuk analisis MANOVA/MANCOVA.
3.1.1 Normalitas Univariat
# Normalitas per variabel numerik
numeric_vars <- c("read","write","math","science")
for (var in numeric_vars) {
shapiro <- shapiro.test(data[[var]])
print(paste("Shapiro-Wilk test for", var))
print(shapiro)
}## [1] "Shapiro-Wilk test for read"
##
## Shapiro-Wilk normality test
##
## data: data[[var]]
## W = 0.97979, p-value = 0.005553
##
## [1] "Shapiro-Wilk test for write"
##
## Shapiro-Wilk normality test
##
## data: data[[var]]
## W = 0.94703, p-value = 9.867e-07
##
## [1] "Shapiro-Wilk test for math"
##
## Shapiro-Wilk normality test
##
## data: data[[var]]
## W = 0.97681, p-value = 0.002145
##
## [1] "Shapiro-Wilk test for science"
##
## Shapiro-Wilk normality test
##
## data: data[[var]]
## W = 0.98525, p-value = 0.03478Uji Shapiro-Wilk untuk tiap variabel dependen menunjukkan p-value < 0,05. Ini berarti secara statistik tiap variabel tidak normal sempurna. Namun, ini normal terjadi pada data nyata, terutama dengan ukuran sampel besar, dan tidak terlalu mengganggu analisis. Jadi, setiap variabel tetap bisa digunakan untuk MANOVA/MANCOVA.
3.1.2 Cek Skewness
library(moments)
numeric_vars <- c("read","write","math","science","socst")
sapply(data[, numeric_vars], skewness)## read write math science socst
## 0.1948373 -0.4784158 0.2844115 -0.1872277 -0.3786624Hasil skewness tiap variabel: read = 0,19 → hampir simetris write = -0,48 → agak miring ke kiri tapi masih moderat math = 0,28 → cukup simetris science = -0,19 → hampir simetris socst = -0,38 → agak miring ke kiri tapi tidak ekstrem Artinya, tidak ada variabel yang miring secara ekstrem, distribusi data relatif seimbang.
Kesimpulan: Meskipun hasil uji formal menunjukkan data tidak normal sempurna, distribusi tiap variabel cukup seimbang dan kemiringannya (skewness) kecil. Distribusi gabungan semua variabel dependen juga cukup mendekati normal, terutama dengan jumlah sampel 200. Jadi, data aman digunakan langsung untuk MANOVA atau MANCOVA tanpa perlu diubah, dan disarankan memakai Pillai’s Trace karena tetap kuat meski ada pelanggaran normalitas minor.
3.2 Homogenitas Varians-Kovarians (Box’s M Test)
Untuk memastikan variansi dan kovarians antar kelompok faktor sama.
library(heplots)## Warning: package 'heplots' was built under R version 4.5.3## Loading required package: broom##
## Attaching package: 'heplots'## The following object is masked from 'package:biotools':
##
## boxMnumeric_vars <- c("read", "write", "math", "science")
boxm_result <- boxM(data[, numeric_vars], group = data$prog)
boxm_result##
## Box's M-test for Homogeneity of Covariance Matrices
##
## data: data[, numeric_vars] by data$prog
## Chi-Sq (approx.) = 22.9113, df = 20, p-value = 0.2932Berdasarkan hasil Box’s M Test, nilai Chi-square sebesar 22,91 dengan p-value = 0,2932. Karena p-value lebih besar dari 0,05, maka homogenitas varians-kovarians antar kelompok prog terpenuhi. Artinya, variasi dan hubungan antar variabel dependen (read, write, math, science) relatif sama di setiap kelompok program. Dengan kata lain, asumsi homogenitas varians-kovarians untuk MANOVA/MANCOVA sudah aman, sehingga analisis dapat dilanjutkan.
3.3 Linearitas dengan Scatterplot
Covariate harus linear dengan variabel dependen.
numeric_vars <- c("read", "write", "math", "science")
covariate <- data$socst
pairs(data[, c(numeric_vars, "socst")],
main = "Linearitas Covariate (socst) vs Variabel Dependen")
Berdasarkan scatterplot matrix, hubungan antara variabel kovariat (socst) dengan masing-masing variabel dependen (read, write, math, dan science) menunjukkan pola yang cenderung membentuk garis lurus ke atas (tren positif). Titik-titik data tersebar mengikuti arah diagonal dari kiri bawah ke kanan atas, tanpa pola melengkung atau menyebar secara acak yang ekstrem. Hal ini menandakan bahwa hubungan antara socst dan setiap variabel dependen bersifat linear. Dengan demikian, asumsi linearitas antara kovariat (socst) dan variabel dependen telah terpenuhi, sehingga analisis lanjutan seperti MANCOVA dapat dilakukan dengan lebih valid.
3.4 Homogenitas Slope
Pastikan interaksi antara faktor dan covariate tidak signifikan (p > 0,05)
lm_read <- lm(read ~ socst * prog, data = data)
anova(lm_read)## Analysis of Variance Table
##
## Response: read
## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 8080.0 8080.0 129.5542 < 2e-16 ***
## prog 2 525.7 262.9 4.2145 0.01615 *
## socst:prog 2 214.4 107.2 1.7192 0.18193
## Residuals 194 12099.3 62.4
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1lm_write <- lm(write ~ socst * prog, data = data)
anova(lm_write)## Analysis of Variance Table
##
## Response: write
## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 6539.6 6539.6 118.8515 < 2e-16 ***
## prog 2 473.0 236.5 4.2979 0.01491 *
## socst:prog 2 191.7 95.8 1.7417 0.17795
## Residuals 194 10674.6 55.0
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1lm_math <- lm(math ~ socst * prog, data = data)
anova(lm_math)## Analysis of Variance Table
##
## Response: math
## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 5177.9 5177.9 91.2979 < 2.2e-16 ***
## prog 2 1186.6 593.3 10.4611 4.847e-05 ***
## socst:prog 2 98.8 49.4 0.8707 0.4203
## Residuals 194 11002.6 56.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1lm_science <- lm(science ~ socst * prog, data = data)
anova(lm_science)## Analysis of Variance Table
##
## Response: science
## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 4219.9 4219.9 54.4553 4.542e-12 ***
## prog 2 205.1 102.5 1.3233 0.2686
## socst:prog 2 48.7 24.4 0.3144 0.7306
## Residuals 194 15033.7 77.5
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Karena semua p-value interaksi lebih besar dari 0,05, dapat disimpulkan bahwa slope covariate (socst) sama di semua kelompok program (prog) untuk tiap variabel dependen. Dengan kata lain, asumsi homogenitas slope terpenuhi, sehingga MANCOVA dapat dijalankan tanpa masalah terkait interaksi covariate × grup.
Selain itu, terlihat bahwa covariate socst berpengaruh signifikan terhadap semua variabel dependen (p < 0,001 untuk semua), yang berarti covariate memang relevan untuk dimasukkan dalam model MANCOVA. Sementara itu, efek prog sendiri signifikan untuk beberapa variabel (read, write, math) tapi tidak untuk science.
4. MANOVA
4.1 Variabel dependen
Variabel numerik disiapkan sebagai dependen
numeric_vars <- c("read", "write", "math", "science")4.2 Model MANOVA dengan semua faktor menggunakan Pillai’s Trace
manova_all <- manova(cbind(read, write, math, science) ~
prog + gender + race + ses + schtyp,
data = data)
summary(manova_all, test = "Pillai")## Df Pillai approx F num Df den Df Pr(>F)
## prog 2 0.34188 9.6907 8 376 3.052e-12 ***
## gender 1 0.20458 12.0243 4 187 1.022e-08 ***
## race 3 0.23282 3.9754 12 567 6.398e-06 ***
## ses 2 0.03368 0.8051 8 376 0.5984
## schtyp 1 0.00543 0.2551 4 187 0.9063
## Residuals 190
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Berdasarkan hasil uji MANOVA menggunakan Pillai’s Trace, diperoleh bahwa variabel prog, gender, dan race memiliki nilai p-value < 0,05. Hal ini menunjukkan bahwa ketiga faktor tersebut berpengaruh signifikan terhadap kombinasi variabel dependen (read, write, math, dan science). Sementara itu, variabel ses dan schtyp memiliki p-value > 0,05, sehingga tidak berpengaruh signifikan secara simultan terhadap variabel dependen.
5. MANCOVA
Analisis MANCOVA dilakukan untuk menguji pengaruh variabel program, gender, ras, tipe sekolah, dan status sosial ekonomi terhadap kemampuan membaca, menulis, matematika, dan sains secara simultan, dengan mengontrol nilai sosial studi sebagai kovariat. Pengujian dilakukan menggunakan statistik uji Pillai.
model_mancova <- manova(
cbind(read, write, math, science) ~
prog + gender + race + ses + schtyp + socst,
data = data
)
summary(model_mancova, test = "Pillai")## Df Pillai approx F num Df den Df Pr(>F)
## prog 2 0.42073 12.4546 8 374 7.329e-16 ***
## gender 1 0.20987 12.3510 4 186 6.282e-09 ***
## race 3 0.27579 4.7581 12 564 1.913e-07 ***
## ses 2 0.04639 1.1101 8 374 0.3554
## schtyp 1 0.00558 0.2609 4 186 0.9027
## socst 1 0.33566 23.4945 4 186 9.778e-16 ***
## Residuals 189
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Hasil analisis MANCOVA menggunakan uji Pillai’s Trace menunjukkan bahwa variabel program (prog), gender, dan ras (race) memiliki pengaruh signifikan terhadap variabel dependen secara simultan yaitu kemampuan membaca, menulis, matematika, dan sains, dengan nilai signifikansi masing-masing kurang dari 0,05. Selain itu, kovariat nilai sosial studi (socst) juga menunjukkan pengaruh yang signifikan (p < 0,05).
Sementara itu, variabel status sosial ekonomi (ses) dan tipe sekolah (schtyp) tidak menunjukkan pengaruh yang signifikan (p > 0,05). Hal ini menunjukkan bahwa perbedaan kemampuan akademik siswa lebih dipengaruhi oleh program, gender, ras, dan nilai sosial studi, sedangkan status sosial ekonomi dan tipe sekolah tidak memberikan pengaruh yang berarti dalam model ini.
6. ANOVA
Analisis ANOVA dilakukan untuk menguji pengaruh variabel program (prog), gender, ras (race), tipe sekolah (schtyp), dan status sosial ekonomi (ses) terhadap masing-masing variabel dependen, yaitu kemampuan membaca (read), menulis (write), matematika (math), dan sains (science) secara terpisah, dengan menggunakan uji F.
# ANOVA untuk kemampuan membaca
anova_read <- aov(read ~ prog + gender + race + schtyp + ses, data = data)
summary(anova_read)## Df Sum Sq Mean Sq F value Pr(>F)
## prog 2 3717 1858.4 22.996 1.14e-09 ***
## gender 1 72 72.4 0.896 0.34502
## race 3 1327 442.4 5.474 0.00125 **
## schtyp 1 15 14.8 0.183 0.66954
## ses 2 433 216.6 2.680 0.07119 .
## Residuals 190 15355 80.8
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# ANOVA untuk menulis
anova_write <- aov(write ~ prog + gender + race + schtyp + ses, data = data)
summary(anova_write)## Df Sum Sq Mean Sq F value Pr(>F)
## prog 2 3176 1587.8 25.069 2.18e-10 ***
## gender 1 1129 1128.7 17.820 3.76e-05 ***
## race 3 1364 454.7 7.179 0.000137 ***
## schtyp 1 5 4.9 0.077 0.781287
## ses 2 171 85.5 1.350 0.261730
## Residuals 190 12034 63.3
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# ANOVA untuk matematika
anova_math <- aov(math ~ prog + gender + race + schtyp + ses, data = data)
summary(anova_math)## Df Sum Sq Mean Sq F value Pr(>F)
## prog 2 4002 2001.1 32.066 1e-12 ***
## gender 1 23 22.5 0.361 0.548843
## race 3 1377 458.9 7.353 0.000109 ***
## schtyp 1 7 7.0 0.113 0.737400
## ses 2 201 100.3 1.607 0.203140
## Residuals 190 11857 62.4
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# ANOVA untuk sains
anova_science <- aov(science ~ prog + gender + race + schtyp + ses, data = data)
summary(anova_science)## Df Sum Sq Mean Sq F value Pr(>F)
## prog 2 1487 743.5 9.869 8.36e-05 ***
## gender 1 330 330.5 4.387 0.0375 *
## race 3 3007 1002.3 13.305 6.40e-08 ***
## schtyp 1 13 12.7 0.169 0.6818
## ses 2 356 178.2 2.366 0.0966 .
## Residuals 190 14314 75.3
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Hasil analisis ANOVA menunjukkan bahwa variabel program (prog) berpengaruh signifikan terhadap seluruh kemampuan akademik, yaitu membaca, menulis, matematika, dan sains (p < 0,05). Variabel ras (race) juga menunjukkan pengaruh signifikan terhadap semua variabel dependen (p < 0,05). Sementara itu, variabel gender hanya berpengaruh signifikan terhadap kemampuan menulis dan sains (p < 0,05), tetapi tidak berpengaruh terhadap kemampuan membaca dan matematika (p > 0,05).
Variabel status sosial ekonomi (ses) tidak menunjukkan pengaruh yang signifikan terhadap seluruh kemampuan akademik (p > 0,05), meskipun pada kemampuan membaca menunjukkan kecenderungan mendekati signifikan. Selain itu, tipe sekolah (schtyp) juga tidak berpengaruh signifikan terhadap semua variabel dependen (p > 0,05). Dengan demikian, dapat disimpulkan bahwa perbedaan kemampuan akademik siswa lebih dipengaruhi oleh program, ras, dan sebagian oleh gender, sedangkan status sosial ekonomi dan tipe sekolah tidak memberikan pengaruh yang berarti.
7. UJI POST HOC ANOVA
TukeyHSD(anova_read, "prog")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = read ~ prog + gender + race + schtyp + ses, data = data)
##
## $prog
## diff lwr upr p adj
## general-academic -6.406349 -10.19005 -2.6226504 0.0002668
## vocational-academic -9.961905 -13.61077 -6.3130378 0.0000000
## vocational-general -3.555556 -7.91913 0.8080187 0.1344998TukeyHSD(anova_read, "race")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = read ~ prog + gender + race + schtyp + ses, data = data)
##
## $race
## diff lwr upr p adj
## asian-african american 3.847019 -4.8995937 12.593632 0.6652624
## hispanic-african american -0.152494 -7.2071544 6.902166 0.9999367
## white-african american 6.138254 0.5803206 11.696188 0.0239716
## hispanic-asian -3.999513 -12.4835482 4.484522 0.6137467
## white-asian 2.291235 -4.9958248 9.578295 0.8473497
## white-hispanic 6.290748 1.1559413 11.425555 0.0093828TukeyHSD(anova_write, "prog")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = write ~ prog + gender + race + schtyp + ses, data = data)
##
## $prog
## diff lwr upr p adj
## general-academic -4.923810 -8.273481 -1.5741378 0.0018409
## vocational-academic -9.497143 -12.727449 -6.2668365 0.0000000
## vocational-general -4.573333 -8.436363 -0.7103035 0.0156484TukeyHSD(anova_write, "race")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = write ~ prog + gender + race + schtyp + ses, data = data)
##
## $race
## diff lwr upr p adj
## asian-african american 8.0070973 0.2638061 15.750388 0.0396526
## hispanic-african american -0.5220734 -6.7674956 5.723349 0.9963973
## white-african american 5.6627708 0.7423861 10.583155 0.0168897
## hispanic-asian -8.5291707 -16.0400043 -1.018337 0.0189824
## white-asian -2.3443265 -8.7954898 4.106837 0.7823036
## white-hispanic 6.1848442 1.6390495 10.730639 0.0029448TukeyHSD(anova_math, "prog")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = math ~ prog + gender + race + schtyp + ses, data = data)
##
## $prog
## diff lwr upr p adj
## general-academic -6.711111 -10.035985 -3.3862376 0.0000110
## vocational-academic -10.313333 -13.519725 -7.1069415 0.0000000
## vocational-general -3.602222 -7.436653 0.2322089 0.0705274TukeyHSD(anova_math, "race")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = math ~ prog + gender + race + schtyp + ses, data = data)
##
## $race
## diff lwr upr p adj
## asian-african american 9.1807073 1.494741 16.8666735 0.0120092
## hispanic-african american 0.7513248 -5.447861 6.9505111 0.9892329
## white-african american 6.2689299 1.384972 11.1528881 0.0057669
## hispanic-asian -8.4293824 -15.884612 -0.9741528 0.0197022
## white-asian -2.9117774 -9.315182 3.4916269 0.6410953
## white-hispanic 5.5176050 1.005464 10.0297464 0.0095589TukeyHSD(anova_science, "prog")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = science ~ prog + gender + race + schtyp + ses, data = data)
##
## $prog
## diff lwr upr p adj
## general-academic -1.355556 -5.008709 2.297598 0.6557515
## vocational-academic -6.580000 -10.102974 -3.057026 0.0000508
## vocational-general -5.224444 -9.437467 -1.011422 0.0106030TukeyHSD(anova_science, "race")## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = science ~ prog + gender + race + schtyp + ses, data = data)
##
## $race
## diff lwr upr p adj
## asian-african american 7.631169 -0.8136684 16.076006 0.0922458
## hispanic-african american 2.459006 -4.3522543 9.270267 0.7856738
## white-african american 10.618493 5.2523192 15.984667 0.0000043
## hispanic-asian -5.172163 -13.3634817 3.019157 0.3606336
## white-asian 2.987324 -4.0483179 10.022966 0.6896481
## white-hispanic 8.159487 3.2018409 13.117133 0.0001829Berdasarkan output Tukey HSD, keputusan didasarkan pada nilai p adj (< 0,05 menunjukkan signifikan). Pada variabel read, write, dan math, program academic berbeda signifikan dengan general dan vocational, menunjukkan bahwa program academic memiliki rata-rata lebih tinggi. Perbedaan antara general dan vocational umumnya tidak signifikan. Pada variabel science, perbedaan signifikan terutama terjadi antara vocational dengan dua program lainnya. Untuk variabel race, terdapat beberapa perbedaan signifikan, terutama melibatkan kelompok white dan asian dibandingkan kelompok lain pada beberapa variabel.
8. ANCOVA
Analisis kovarians (ANCOVA) dilakukan untuk menguji pengaruh variabel program (prog) terhadap kemampuan membaca, menulis, matematika, dan sains dengan mengontrol variabel kovariat, yaitu status sosial (socst), serta untuk menguji interaksi antara program dan status sosial menggunakan uji F.
ancova_read <- aov(read ~ socst + prog + gender + race + schtyp + ses, data = data)
summary(ancova_read)## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 8080 8080 132.737 <2e-16 ***
## prog 2 526 263 4.318 0.0147 *
## gender 1 149 149 2.443 0.1197
## race 3 602 201 3.298 0.0216 *
## schtyp 1 6 6 0.092 0.7615
## ses 2 52 26 0.429 0.6516
## Residuals 189 11505 61
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ancova_write <- aov(write ~ socst + prog + gender + race + schtyp + ses, data = data)
summary(ancova_write)## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 6540 6540 137.058 < 2e-16 ***
## prog 2 473 236 4.956 0.007981 **
## gender 1 923 923 19.354 1.81e-05 ***
## race 3 865 288 6.042 0.000598 ***
## schtyp 1 12 12 0.254 0.615120
## ses 2 48 24 0.501 0.606559
## Residuals 189 9018 48
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ancova_math <- aov(math ~ socst + prog + gender + race + schtyp + ses, data = data)
summary(ancova_math)## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 5178 5178 97.393 < 2e-16 ***
## prog 2 1187 593 11.159 2.62e-05 ***
## gender 1 53 53 1.003 0.317756
## race 3 985 328 6.178 0.000501 ***
## schtyp 1 3 3 0.051 0.822103
## ses 2 12 6 0.111 0.895379
## Residuals 189 10048 53
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ancova_science <- aov(science ~ socst + prog + gender + race + schtyp + ses, data = data)
summary(ancova_science)## Df Sum Sq Mean Sq F value Pr(>F)
## socst 1 4220 4220 64.161 1.15e-13 ***
## prog 2 205 103 1.559 0.2130
## gender 1 443 443 6.738 0.0102 *
## race 3 2136 712 10.828 1.34e-06 ***
## schtyp 1 6 6 0.096 0.7574
## ses 2 66 33 0.501 0.6068
## Residuals 189 12431 66
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analisis kovarians (ANCOVA) menunjukkan bahwa variabel socst berpengaruh signifikan terhadap seluruh kemampuan akademik, yaitu membaca, menulis, matematika, dan sains (p < 0,001). Variabel prog juga berpengaruh signifikan terhadap kemampuan membaca, menulis, dan matematika (p < 0,05), namun tidak menunjukkan pengaruh yang signifikan terhadap kemampuan sains (p > 0,05). Sementara itu, hasil uji interaksi antara socst dan prog tidak signifikan pada semua variabel dependen, yang mengindikasikan bahwa tidak terdapat perbedaan pengaruh socst pada masing-masing kategori program. Dengan demikian, dapat disimpulkan bahwa socst merupakan kovariat yang konsisten memengaruhi seluruh kemampuan akademik, sedangkan pengaruh prog terbatas pada beberapa aspek saja.
9. UJI POST HOC ANCOVA
#install.packages("emmeans")
library(emmeans)## Warning: package 'emmeans' was built under R version 4.5.3## Welcome to emmeans.
## Caution: You lose important information if you filter this package's results.
## See '? untidy'emmeans(ancova_read, pairwise ~ prog, adjust = "tukey")## $emmeans
## prog emmean SE df lower.CL upper.CL
## academic 52.6 1.15 189 50.4 54.9
## general 49.3 1.49 189 46.4 52.3
## vocational 49.0 1.66 189 45.7 52.3
##
## Results are averaged over the levels of: gender, race, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## academic - general 3.274 1.46 189 2.242 0.0669
## academic - vocational 3.620 1.57 189 2.309 0.0570
## general - vocational 0.347 1.69 189 0.206 0.9770
##
## Results are averaged over the levels of: gender, race, schtyp, ses
## P value adjustment: tukey method for comparing a family of 3 estimatesemmeans(ancova_read, pairwise ~ race, adjust = "tukey")## $emmeans
## race emmean SE df lower.CL upper.CL
## african american 48.1 1.90 189 44.3 51.8
## asian 52.1 2.53 189 47.2 57.1
## hispanic 48.5 1.78 189 45.0 52.0
## white 52.6 0.94 189 50.7 54.4
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## african american - asian -4.073 2.98 189 -1.368 0.5212
## african american - hispanic -0.415 2.39 189 -0.174 0.9981
## african american - white -4.524 1.94 189 -2.333 0.0943
## asian - hispanic 3.658 2.89 189 1.264 0.5872
## asian - white -0.451 2.48 189 -0.182 0.9979
## hispanic - white -4.109 1.77 189 -2.316 0.0980
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## P value adjustment: tukey method for comparing a family of 4 estimatesemmeans(ancova_write, pairwise ~ prog, adjust = "tukey")## $emmeans
## prog emmean SE df lower.CL upper.CL
## academic 54.0 1.02 189 52.0 56.0
## general 51.7 1.32 189 49.1 54.3
## vocational 50.5 1.47 189 47.6 53.4
##
## Results are averaged over the levels of: gender, race, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## academic - general 2.30 1.29 189 1.776 0.1804
## academic - vocational 3.49 1.39 189 2.511 0.0343
## general - vocational 1.19 1.49 189 0.796 0.7058
##
## Results are averaged over the levels of: gender, race, schtyp, ses
## P value adjustment: tukey method for comparing a family of 3 estimatesemmeans(ancova_write, pairwise ~ race, adjust = "tukey")## $emmeans
## race emmean SE df lower.CL upper.CL
## african american 48.8 1.690 189 45.5 52.1
## asian 57.4 2.240 189 53.0 61.8
## hispanic 48.8 1.580 189 45.7 51.9
## white 53.3 0.832 189 51.6 54.9
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## african american - asian -8.6232 2.64 189 -3.271 0.0069
## african american - hispanic -0.0436 2.11 189 -0.021 1.0000
## african american - white -4.4934 1.72 189 -2.617 0.0468
## asian - hispanic 8.5797 2.56 189 3.348 0.0054
## asian - white 4.1299 2.19 189 1.883 0.2389
## hispanic - white -4.4498 1.57 189 -2.833 0.0261
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## P value adjustment: tukey method for comparing a family of 4 estimatesemmeans(ancova_math, pairwise ~ prog, adjust = "tukey")## $emmeans
## prog emmean SE df lower.CL upper.CL
## academic 54.6 1.07 189 52.5 56.7
## general 49.9 1.39 189 47.2 52.7
## vocational 48.9 1.55 189 45.9 52.0
##
## Results are averaged over the levels of: gender, race, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## academic - general 4.69 1.36 189 3.433 0.0021
## academic - vocational 5.70 1.47 189 3.891 0.0004
## general - vocational 1.02 1.58 189 0.645 0.7954
##
## Results are averaged over the levels of: gender, race, schtyp, ses
## P value adjustment: tukey method for comparing a family of 3 estimatesemmeans(ancova_math, pairwise ~ race, adjust = "tukey")## $emmeans
## race emmean SE df lower.CL upper.CL
## african american 47.3 1.780 189 43.8 50.8
## asian 56.7 2.360 189 52.0 61.3
## hispanic 48.3 1.670 189 45.0 51.6
## white 52.4 0.878 189 50.7 54.1
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## african american - asian -9.36 2.78 189 -3.363 0.0051
## african american - hispanic -1.03 2.23 189 -0.463 0.9669
## african american - white -5.10 1.81 189 -2.812 0.0277
## asian - hispanic 8.33 2.71 189 3.078 0.0127
## asian - white 4.26 2.32 189 1.842 0.2570
## hispanic - white -4.06 1.66 189 -2.450 0.0714
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## P value adjustment: tukey method for comparing a family of 4 estimatesemmeans(ancova_science, pairwise ~ race, adjust = "tukey")## $emmeans
## race emmean SE df lower.CL upper.CL
## african american 44.3 1.980 189 40.4 48.2
## asian 52.0 2.630 189 46.8 57.1
## hispanic 47.1 1.850 189 43.4 50.8
## white 53.6 0.977 189 51.7 55.5
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## african american - asian -7.65 3.10 189 -2.472 0.0676
## african american - hispanic -2.79 2.48 189 -1.124 0.6750
## african american - white -9.31 2.02 189 -4.621 <0.0001
## asian - hispanic 4.86 3.01 189 1.616 0.3722
## asian - white -1.66 2.58 189 -0.645 0.9170
## hispanic - white -6.52 1.84 189 -3.538 0.0028
##
## Results are averaged over the levels of: prog, gender, schtyp, ses
## P value adjustment: tukey method for comparing a family of 4 estimatesPada hasil post hoc ANCOVA menggunakan emmeans, diperoleh bahwa setelah mengontrol kovariat (socst), beberapa perbedaan antar kelompok menjadi tidak signifikan. Hal ini terlihat pada variabel read dan write, di mana perbedaan antar program tidak lagi signifikan (p-value > 0,05). Namun, pada variabel math dan science, beberapa perbedaan masih tetap signifikan, khususnya antara program academic dan vocational. Secara keseluruhan, output R menunjukkan bahwa perbedaan utama terdapat pada program pendidikan dan ras, namun pengaruh tersebut dapat berubah setelah mempertimbangkan variabel kovariat.