Korelasi di R: cor(), cor.test(), dan matriks korelasi

Korelasi antara dua variabel: cor()

Korelasi bertanya: saat satu variabel naik, apakah variabel lainnya cenderung ikut naik (positif), turun (negatif), atau bergerak sendiri (mendekati nol)? Di R, cukup satu pemanggilan:

print(cor(mtcars$wt, mtcars$mpg))

Jawabannya sekitar −0,87: mobil yang lebih berat lebih boros bahan bakar, dan hubungannya kuat. Angka tunggal itu adalah koefisien korelasi Pearson, yang di mana-mana ditulis r.

Membaca r

Koefisiennya selalu berada di antara −1 dan +1. Tandanya menunjukkan arah; besarnya menunjukkan kekuatan:

| |r| | Cara membaca yang umum | | --- | --- | | 0,0 sampai 0,2 | dapat diabaikan | | 0,2 sampai 0,4 | lemah | | 0,4 sampai 0,6 | sedang | | 0,6 sampai 0,8 | kuat | | 0,8 sampai 1,0 | sangat kuat |

Anggap rentang ini sebagai bahan diskusi, bukan hukum: di fisika, r sebesar 0,6 mengecewakan, sedangkan di psikologi itu bisa jadi puncak karier. Dua sifat yang perlu dipahami betul: r tidak punya satuan (mengorelasikan berat dalam ton dengan mpg menghasilkan r yang sama dengan berat dalam kilogram dengan mpg, karena r dihitung dari nilai yang sudah distandarkan), dan r hanya mengukur hubungan linear. Hubungan berbentuk U yang sempurna bisa punya r ≈ 0.

Dan kalimat yang wajib diucapkan: korelasi bukan kausalitas. Penjualan es krim dan kematian karena tenggelam berkorelasi kuat sepanjang bulan-bulan dalam setahun, bukan karena es krim menenggelamkan orang, tapi karena musim panas mendorong keduanya. Korelasi memberi tahu bahwa dua variabel bergerak bersama; korelasi diam soal mengapa. Mungkin x mendorong y, mungkin y mendorong x, mungkin ada hal ketiga (dalam kasus ini, musim) yang mendorong keduanya. Memutuskan di antara kemungkinan itu butuh eksperimen atau penalaran sebab-akibat yang cermat, bukan r yang lebih besar.

Spearman dan Kendall: saat Pearson bukan alat yang tepat

Pearson bekerja pada nilai mentah, sehingga peka terhadap pencilan dan tidak bisa melihat hubungan yang melengkung. Argumen method beralih ke alternatif berbasis peringkat:

x <- c(1, 2, 3, 4, 5, 6, 7, 8)
y <- c(1, 4, 9, 16, 25, 36, 49, 300)  # increasing, with a wild last value

print(cor(x, y))                        # Pearson - dragged by the outlier
print(cor(x, y, method = "spearman"))   # ranks: perfectly monotonic
print(cor(x, y, method = "kendall"))

Spearman mengganti setiap nilai dengan peringkatnya lalu menghitung Pearson pada peringkat itu. Karena y selalu naik saat x naik, setiap peringkat sejajar dan Spearman melaporkan tepat 1: besarnya pencilan tidak lagi berpengaruh, hanya posisinya. Pilih Spearman saat datanya ordinal (skala survei), sangat miring, atau saat hubungannya monoton tapi tidak lurus. Kendall menjawab pertanyaan serupa dari pasangan yang searah dan berlawanan arah (concordant/discordant); hasilnya lebih tahan pada sampel kecil tapi lebih lambat, dan Spearman adalah default yang umum.

Matriks korelasi

Untuk memeriksa hubungan di antara banyak variabel sekaligus, berikan beberapa kolom numerik ke cor():

m <- cor(mtcars[, c("mpg", "wt", "hp", "disp")])
print(round(m, 2))

Setiap variabel dibandingkan dengan setiap variabel lain, dengan angka 1 di diagonal (setiap variabel berkorelasi sempurna dengan dirinya sendiri) dan bayangan cermin di kedua sisinya. Membulatkan ke dua desimal lebih penting daripada kedengarannya: matriks yang tidak dibulatkan adalah tembok angka, padahal gunanya matriks adalah untuk dipindai. Di sini pemindaian menunjukkan mpg berkorelasi negatif dengan ketiganya (mobil yang lebih berat, lebih bertenaga, dan bermesin lebih besar memakai lebih banyak bahan bakar), sedangkan wt, hp, dan disp semuanya berkorelasi positif kuat satu sama lain: sekelompok variabel "mobil besar" yang akan berpengaruh saat kamu sampai di regresi linear dan sakit kepala multikolinearitasnya.

Nilai yang hilang: argumen use

Saat ada data yang hilang, default cor() adalah mengembalikan NA alih-alih menebak:

d <- mtcars[, c("mpg", "wt", "hp")]
d$wt[c(3, 7)] <- NA  # knock two values out

print(cor(d)[1, 2])
print(cor(d, use = "complete.obs")[1, 2])
print(cor(d, use = "pairwise.complete.obs")[1, 2])

• use = "complete.obs" membuang setiap baris yang berisi NA apa pun, lalu menghitung seluruh matriks dari baris yang tersisa. Konsisten, tapi boros: wt yang hilang juga menghapus baris itu dari pasangan mpg dan hp.
• use = "pairwise.complete.obs" menghitung setiap sel dari semua baris tempat pasangan itu ada. Lebih banyak data yang dipertahankan, tapi sel yang berbeda bertumpu pada subset yang berbeda, yang sesekali bisa menghasilkan matriks yang tidak konsisten secara internal.

Untuk satu atau dua NA yang nyasar, keduanya tidak masalah; cukup sebutkan mana yang kamu pakai.

Apakah signifikan? cor.test()

cor() memberikan angka, tapi tidak memberi gambaran apakah angka itu bisa jadi hanya derau. cor.test() menjalankan uji hipotesisnya:

print(cor.test(mtcars$wt, mtcars$mpg))

Telusuri output-nya blok demi blok:

• t = −9.56, df = 30: statistik ujinya. Hipotesis nolnya adalah korelasi sebenarnya nol; r yang teramati diubah menjadi statistik t dengan n − 2 derajat kebebasan (32 mobil − 2).
• p-value = 1.29e-10: kalau korelasi sebenarnya nol, peluang melihat r sejauh ini dari nol pada sampel 32 mobil adalah sekitar 0,0000000001. Itu bukti yang sangat kuat bahwa hubungannya nyata, tapi ingat, p-value berbicara tentang apakah r berbeda dari nol, bukan apakah hubungannya besar atau kausal.
• 95 percent confidence interval: −0.93 to −0.74: rentang yang masuk akal untuk korelasi sebenarnya. Sering kali lebih berguna daripada p-value: bahkan ujung yang paling optimis dari interval ini pun masih korelasi negatif yang kuat.
• sample estimates: cor = −0.87: angka yang sama dengan yang diberikan cor().

Sampel kecil butuh kewaspadaan ekstra di sini: dengan n = 10, korelasi ±0,5 muncul secara kebetulan dengan frekuensi yang mengkhawatirkan, dan interval kepercayaan yang lebar akan memberitahumu hal itu. Laporkan intervalnya, bukan hanya p-value-nya.

Melihatnya: selalu buat plot

Koefisien korelasi memadatkan seluruh hubungan menjadi satu angka, dan pemadatan itu bisa menyembunyikan lengkungan, gerombolan, atau satu titik yang melakukan semua pekerjaannya. Sebelum memercayai r mana pun, lihat scatter plot-nya:

plot(mtcars$wt, mtcars$mpg)              # one pair
pairs(mtcars[, c("mpg", "wt", "hp", "disp")])  # every pair in the matrix

pairs() menggambar kisi scatter plot yang sesuai dengan matriks korelasimu, cara tercepat untuk memastikan angka-angkanya berarti seperti yang kamu kira. Untuk grafik matriks bergaya heatmap yang rapi, package corrplot (install.packages("corrplot"), lalu corrplot(cor(m))) adalah alat standarnya.

Yang kamu bawa pulang

• cor(x, y) memberikan r Pearson: tanda adalah arah, besar adalah kekuatan, selalu di [−1, 1], tanpa satuan.
• Korelasi mengukur gerak bersama yang linear dan tidak mengatakan apa pun tentang sebab-akibat: variabel ketiga yang tersembunyi selalu jadi kandidat.
• method = "spearman" untuk peringkat: data ordinal, pencilan, hubungan yang monoton tapi melengkung.
• cor(df) pada kolom numerik memberikan matriksnya; bulatkan dengan round(, 2), dan perhatikan argumen use = saat ada data yang hilang.
• cor.test(x, y) menambahkan p-value dan interval kepercayaan; laporkan intervalnya.
• Selalu lihat scatter plot sebelum memercayai angkanya.

Berikutnya: saat pertanyaannya menajam dari "apakah keduanya bergerak bersama?" menjadi "apakah rata-rata kelompok ini berbeda dari kelompok itu?", yaitu uji t.