Sedikit inkonsistensi antara fungsi R bawaan Kruskal-Wallis dan perhitungan manual

Saya bingung dengan yang berikut ini, dan saya belum bisa menggali jawabannya di tempat lain.

Saya mencoba mempelajari R sambil melakukan beberapa statistik, dan, sebagai latihan, saya mencoba memeriksa ulang hasil fungsi R bawaan dengan juga melakukan ini 'dengan tangan', seperti di R. , untuk tes Kruskal-Wallis saya terus mendapatkan hasil yang berbeda, dan saya tidak tahu mengapa.

Sebagai contoh, saya melihat data berikut yang diberikan dalam latihan

activity <- c(2, 4, 3, 2, 3, 3, 4, 0, 4, 3, 4, 0, 0, 1, 3, 1, 2, 0, 3, 1, 0, 3, 4, 0, 1, 2, 2, 2, 3, 2) 
group <- c(rep("A", 11), rep("B", 10), rep("C", 9))
group <- factor(group)
data.raw <- data.frame(activity, group)

Dan saya ingin menganalisis aktivitas berdasarkan kelompok. Pertama saya menjalankan tes Kruskal-Wallis menggunakan fungsi R built-in

kruskal.test(activity ~ group, data = data.raw)

Yang mengembalikan . $H = 8.9056$

Untuk mengecek, saya mencoba melakukan hal yang sama 'dengan tangan' di R, dengan kode berikut (tidak diragukan lagi tak berdaya)

rank <- rank(activity)
data.rank <- data.frame(rank, group)
rank.sum <- aggregate(rank ~ group, data = data.rank, sum)

x <- rank.sum[1,2]^2 / 11 + rank.sum[2,2]^2 / 10 + rank.sum[3,2]^2 / 9
H <- (12 / (length(activity) * (length(activity) + 1))) * x - 3 * (length(activity) + 1)
H

Yang dimaksudkan untuk mencerminkan rumus berikut:

H = \frac{12}{N (N + 1)} \sum_{i = 1}^{g} (\frac{R_{i}^{2}}{n_{i}}) - 3 (N + 1)

$H =\frac{12}{N(N+1)}\sum_{i = 1}^g \left(\frac{R^2_i}{n_i} \right) - 3(N + 1)$

Di mana adalah jumlah total pengamatan, adalah jumlah kelompok, adalah jumlah pengamatan dalam kelompok ke- , dan adalah jumlah peringkat dari kelompok ke- . $N$ $g$ $n_i$ $i$ $R_i$ $i$

Dan sekarang saya mendapatkan , yang, menambah kebingungan saya, juga merupakan jawaban yang diberikan untuk latihan tersebut. Saya sudah mencoba ini untuk beberapa set data yang berbeda, dan saya cenderung mendapatkan nilai yang sedikit lebih tinggi untuk menggunakan fungsi built-in. $H = 8.499$ $H$

Saya sudah mencoba mencari tahu apa yang saya lakukan salah atau gagal untuk mengerti, tetapi tidak berhasil. Adakah yang bisa membantu saya memahami mengapa kruskal.testfungsi inbuilt mengembalikan nilai yang berbeda dari yang saya dapatkan dengan mengeja semuanya?

r kruskal-wallis ties MSR
sumber

kruskal.testmenerapkan koreksi untuk ikatan seperti yang dijelaskan dalam artikel Wikipedia ini (poin 4):

Koreksi untuk ikatan jika menggunakan rumus pintasan yang dijelaskan pada poin sebelumnya dapat dilakukan dengan membagi H dengan , ... $1 - \frac{\sum_{i=1}^G (t_i^3 - t_i)}{N^3-N}$

Melanjutkan dari kode Anda:

TIES <- table(activity)
H / (1 - sum(TIES^3 - TIES)/(length(activity)^3 - length(activity)))
#[1] 8.9056

Anda dapat mengetahui apa fungsi R dengan mempelajari kode secara hati-hati, yang dapat Anda lihat menggunakan getAnywhere(kruskal.test.default).

Roland
sumber

@MichaelChernick Tidak, tidak. Intinya adalah bahwa OP telah diajarkan penyederhanaan tes yang harus digunakan hanya jika tidak ada ikatan.

Roland

@MichaelChernick Saya tidak mengatakan itu tidak muat di Stack Overflow. Tapi saya berpendapat bahwa itu cocok juga di CV. Jelas, akan sangat membantu jika OP tidak hanya membagikan kode mereka tetapi juga formula yang mereka gunakan.

Roland

@Michael Status utas ini adalah panggilan yang mudah: tepat di lingkup kami karena berusaha memahami uji statistik.

whuber

Diedit untuk memasukkan rumus yang tercermin dalam kode. Seharusnya berpikir untuk melakukannya pertama kali. Permintaan maaf.

MSR

Lihat juga fungsi Hmiscpaket R spearman2yang menggunakan midranks untuk pengikatan dan Ftes untuk mendapatkan Kruskal-Wallis. Saya pikir ini lebih akurat daripada beberapa metode.

Frank Harrell

Sedikit inkonsistensi antara fungsi R bawaan Kruskal-Wallis dan perhitungan manual

Jawaban: