Nilai P sama dengan 0 dalam uji permutasi

15

Saya memiliki dua set data dan saya ingin tahu apakah mereka berbeda secara signifikan atau tidak (ini berasal dari " Dua kelompok berbeda secara signifikan? Tes untuk digunakan ").

Saya memutuskan untuk menggunakan tes permutasi, melakukan hal berikut di R:

permutation.test <- function(coding, lncrna) {
    coding <- coding[,1] # dataset1
    lncrna <- lncrna[,1] # dataset2

    ### Under null hyphotesis, both datasets would be the same. So:
    d <- c(coding, lncrna)

    # Observed difference
    diff.observed = mean(coding) - mean(lncrna)
    number_of_permutations = 5000
    diff.random = NULL

    for (i in 1:number_of_permutations) {
        # Sample from the combined dataset
        a.random = sample (d, length(coding), TRUE)
        b.random = sample (d, length(lncrna), TRUE)
        # Null (permuated) difference
        diff.random[i] = mean(b.random) - mean(a.random)
    }

    # P-value is the fraction of how many times the permuted difference is equal or more extreme than the observed difference
    pvalue = sum(abs(diff.random) >= abs(diff.observed)) / number_of_permutations
    pvalue
}

Namun demikian, nilai-p tidak boleh 0 menurut makalah ini: http://www.statsci.org/smyth/pubs/permp.pdf

Apa yang Anda rekomendasikan untuk saya lakukan? Apakah ini cara untuk menghitung nilai p:

pvalue = sum(abs(diff.random) >= abs(diff.observed)) / number_of_permutations

cara yang baik? Atau lebih baik melakukan hal berikut?

pvalue = sum(abs(diff.random) >= abs(diff.observed)) + 1 / number_of_permutations + 1
pengguna2886545
sumber
(1) Baris terakhir dalam pertanyaan salah karena tidak termasuk tanda kurung yang diperlukan untuk menjalankan perhitungan yang dimaksud. (Dijamin menghasilkan hasil lebih besar dari , yang tidak mungkin untuk nilai p apa pun.) (2) Anda tidak benar-benar melakukan tes permutasi: dua sampel dan jarang akan terdiri dari partisi acak data tetapi biasanya akan tumpang tindih secara substansial. Sebaliknya, hitung sebagai pelengkap dari dalam persatuan dan . 1a.randomb.randomb.randoma.randomcodinglncrna
whuber
Karena nilai-p adalah himpunan nilai setidaknya yang ekstrim seperti yang diamati, jika seseorang mengevaluasi distribusi permutasi, statistik yang diamati adalah dalam "permutasi" yang dihitung. Ketika melakukan pengacakan, itu biasa untuk menghitung statistik yang diamati di antara statistik permutasi yang dipertimbangkan (untuk alasan yang sama).
Glen_b -Reinstate Monica

Jawaban:

15

Diskusi

Tes permutasi menghasilkan semua permutasi yang relevan dari suatu dataset, menghitung statistik uji yang ditunjuk untuk setiap permutasi tersebut, dan menilai statistik pengujian aktual dalam konteks distribusi permutasi yang dihasilkan dari statistik. Cara yang umum untuk menilai itu adalah melaporkan proporsi statistik yang (dalam beberapa hal) "sebagai atau lebih ekstrem" daripada statistik aktual. Ini sering disebut "nilai-p".

Karena dataset aktual adalah salah satu permutasi itu, statistiknya tentu akan berada di antara yang ditemukan dalam distribusi permutasi. Oleh karena itu, nilai-p tidak pernah nol.

Kecuali jika dataset sangat kecil (kurang dari sekitar 20-30 jumlah total, biasanya) atau statistik uji memiliki bentuk matematika yang sangat bagus, tidak praktis untuk menghasilkan semua permutasi. (Contoh di mana semua permutasi dihasilkan muncul di Uji Permutasi di R. ) Oleh karena itu implementasi komputer dari tes permutasi biasanya sampel dari distribusi permutasi. Mereka melakukannya dengan menghasilkan beberapa permutasi acak independen dan berharap bahwa hasilnya adalah sampel representatif dari semua permutasi.

Oleh karena itu, angka apa pun (seperti "nilai-p") yang berasal dari sampel semacam itu hanyalah penaksir properti dari distribusi permutasi. Sangat mungkin - dan sering terjadi ketika efeknya besar - bahwa nilai p yang diperkirakan adalah nol. Tidak ada yang salah dengan hal itu, tetapi hal itu segera menimbulkan masalah yang sebelumnya diabaikan tentang seberapa besar estimasi nilai p berbeda dari yang benar? Karena distribusi sampling proporsi (seperti estimasi nilai-p) adalah Binomial, ketidakpastian ini dapat diatasi dengan interval kepercayaan Binomial .


Arsitektur

Implementasi yang dibangun dengan baik akan mengikuti diskusi dengan cermat dalam segala hal. Ini akan dimulai dengan rutin untuk menghitung statistik tes, karena ini untuk membandingkan cara dua kelompok:

diff.means <- function(control, treatment) mean(treatment) - mean(control)

Tulis rutin lain untuk menghasilkan permutasi acak dataset dan menerapkan statistik uji. Antarmuka yang satu ini memungkinkan penelepon untuk menyediakan statistik uji sebagai argumen. Ini akan membandingkan melemen pertama dari sebuah array (dianggap sebagai grup referensi) dengan elemen lainnya (grup "perawatan").

f <- function(..., sample, m, statistic) {
  s <- sample(sample)
  statistic(s[1:m], s[-(1:m)])
}

Tes permutasi dilakukan pertama-tama dengan menemukan statistik untuk data aktual (diasumsikan di sini untuk disimpan dalam dua array controldan treatment) dan kemudian menemukan statistik untuk banyak permutasi acak independen daripadanya:

z <- stat(control, treatment) # Test statistic for the observed data
sim<- sapply(1:1e4, f, sample=c(control,treatment), m=length(control), statistic=diff.means)

Sekarang hitung estimasi binomial dari nilai-p dan interval kepercayaan untuknya. Satu metode menggunakan binconfprosedur bawaan dalam HMiscpaket:

require(Hmisc)                                    # Exports `binconf`
k <- sum(abs(sim) >= abs(z))                      # Two-tailed test
zapsmall(binconf(k, length(sim), method='exact')) # 95% CI by default

Ini bukan ide yang buruk untuk membandingkan hasilnya dengan tes lain, bahkan jika itu diketahui tidak cukup berlaku: setidaknya Anda mungkin mendapatkan urutan besarnya di mana hasilnya seharusnya terletak. Dalam contoh ini (alat pembanding), Student t-test biasanya memberikan hasil yang baik:

t.test(treatment, control)

Arsitektur ini diilustrasikan dalam situasi yang lebih kompleks, dengan Rkode kerja , di Test Apakah Variabel Ikuti Distribusi yang Sama .


Contoh

100201.5

set.seed(17)
control <- rnorm(10)
treatment <- rnorm(20, 1.5)

Setelah menggunakan kode sebelumnya untuk menjalankan tes permutasi, saya merencanakan sampel distribusi permutasi bersama dengan garis merah vertikal untuk menandai statistik aktual:

h <- hist(c(z, sim), plot=FALSE)
hist(sim, breaks=h$breaks)
abline(v = stat(control, treatment), col="Red")

Angka

Perhitungan batas kepercayaan binomial menghasilkan

 PointEst Lower        Upper
        0     0 0.0003688199

00,000373.16e-050,000370,000370,050,010,001


Komentar

kN k/N(k+1)/(N+1)N

10102=1000,0000051.611.7bagian per juta: sedikit lebih kecil dari uji-t Student yang dilaporkan. Meskipun data dihasilkan dengan generator angka acak normal, yang akan membenarkan menggunakan uji-t Student, hasil tes permutasi berbeda dari hasil uji-t Student karena distribusi dalam setiap kelompok pengamatan tidak normal normal.

whuber
sumber
Makalah oleh Smyth & Phipson yang dikutip di atas dengan jelas menunjukkan mengapa k / N adalah pilihan yang buruk untuk penduga nilai-p. Singkatnya, untuk tingkat signifikansi yang relevan seperti alpha = 0,05, P ((k / N) <alpha | H0) dapat secara mengejutkan lebih besar dari alpha. Ini berarti bahwa uji permutasi acak menggunakan k / N sebagai penaksir nilai-p dan 0,05 sebagai ambang penolakannya akan menolak hipotesis nol lebih dari 5% kali! Nilai nol p adalah kasus ekstrim dari masalah ini - dengan kriteria alpha = 0 kami berharap untuk tidak pernah menolak nol, namun b / m dapat sama dengan nol di bawah nol, yang mengarah ke penolakan palsu.
Trisoloriansunscreen
1
@Tal "Pilihan yang buruk" untuk tujuan tertentu. Yang membedakan kami sebagai ahli statistik dari yang lain adalah pemahaman kami tentang peran variabilitas dalam analisis data dan pengambilan keputusan, bersama dengan kemampuan kami untuk mengukur variabilitas itu dengan tepat. Itulah pendekatan yang dicontohkan (dan secara tersirat dianjurkan) dalam jawaban saya di sini. Ketika itu dilakukan tidak ada masalah seperti yang Anda gambarkan, karena pengguna dari prosedur permutasi diarahkan untuk memahami keterbatasan dan kekuatannya dan akan memiliki kebebasan untuk bertindak sesuai dengan tujuannya.
whuber
13

Karena estimasi nilai-p digunakan untuk memutuskan apakah akan menolak hipotesis nol, penting untuk mempertimbangkan bagaimana pilihan estimator memengaruhi probabilitas penolakan palsu. Makalah yang dikutip oleh Smyth & Phipson menunjukkan bahwa penduga tidak bias (BM.) gagal mengendalikan tingkat kesalahan tipe-I dengan benar. Sebaliknya, (B+1M.+1) adalah penaksir nilai-nilai yang valid (tapi konservatif) - tidak menyebabkan penolakan berlebih dari nol.

(B adalah jumlah permutasi acak di mana statistik lebih besar atau sama dengan yang diamati diperoleh dan M adalah jumlah total permutasi acak sampel).

Smyth & Phipson juga menunjukkan bahwa ketidakabsahan (BM.) menjadi kritis dalam beberapa pengaturan perbandingan, di mana estimasi nilai p yang sangat kecil diturunkan dan kemudian dikoreksi dengan perkalian dengan suatu faktor. Perkiraan nilai p nol di bawah nol sangat berbahaya dalam pengaturan ini, karena tetap nol terlepas dari koreksi yang diterapkan.

Trisoloriansunscreen
sumber
1
+1 Ini adalah ringkasan bagus dari poin utama makalah ini. Saya terutama menghargai perhatian Anda pada perbedaan antara nilai p yang diperkirakan dan nilai p permutasi yang sebenarnya.
whuber