Kontras polinomial untuk regresi

17

Saya tidak dapat memahami penggunaan kontras polinomial dalam penyesuaian regresi. Secara khusus, saya mengacu pada pengkodean yang digunakan Runtuk mengekspresikan variabel interval (variabel ordinal dengan tingkat spasi yang sama), dijelaskan di halaman ini .

Dalam contoh halaman itu , jika saya mengerti dengan benar, R cocok dengan model untuk variabel interval, mengembalikan beberapa koefisien yang menimbang tren linier, kuadrat, atau kubik. Oleh karena itu, model yang dipasang harus:

wrsayate=52.7870+14.2587X-0,9680X2-0,1554X3,

di mana harus mengambil nilai , , , atau sesuai dengan level variabel interval yang berbeda.X1234

Apakah ini benar? Dan, jika demikian, apa tujuan dari kontras polinomial?

Pippo
sumber
7
Tidak, koefisien tersebut adalah untuk istilah polinom ortogonal : Anda telah menulis model untuk istilah polinomial mentah . Ganti X , X2 , & X3 , dengan nilai L. , Q , & C masing-masing (dari tabel pencarian).
Scortchi
1
Yang terhormat @Scortchi, terima kasih atas balasan Anda. Saya kira mengerti maksud Anda, tetapi kemudian saya tidak benar-benar mengerti bagaimana istilah polinomial ortogonal ini bekerja. : P
Pippo
1
Sebagai masalah notasi, apa yang Anda miliki bukan model yang pas. Anda juga perlu 'topi' lebih dari menulis (atau E [tulis]), yang berarti nilai prediksi dari penulisan atau nilai yang diharapkan dari penulisan; atau Anda memerlukan '+ e' di akhir untuk mengindikasikan residu.
gung - Reinstate Monica
@ Scortchi Apa itu, atau bagaimana Anda dapat menemukan, "tabel pencarian"?
Antoni Parellada
2
@AntoniParellada: Ini adalah tabel di halaman yang ditautkan oleh OP: ats.ucla.edu/stat/r/library/contrast_coding.htm#ORTHOGONAL . & didapat contr.polydi R.
Scortchi - Reinstate Monica

Jawaban:

29

Sekadar rekap (dan jika hyperlink OP gagal di masa mendatang), kami melihat dataset hsb2seperti:

   id     female race ses schtyp prog read write math science socst
1  70        0    4   1      1    1   57    52   41      47    57
2 121        1    4   2      1    3   68    59   53      63    61
...
199 118      1    4   2      1    1   55    62   58      58    61
200 137      1    4   3      1    2   63    65   65      53    61

yang bisa diimpor di sini .

Kami mengubah variabel readmenjadi dan memerintahkan / variabel ordinal:

hsb2$readcat<-cut(hsb2$read, 4, ordered = TRUE)
(means = tapply(hsb2$write, hsb2$readcat, mean))
 (28,40]  (40,52]  (52,64]  (64,76] 
42.77273 49.97849 56.56364 61.83333 

Sekarang kita siap untuk menjalankan ANOVA biasa - ya, itu adalah R, dan kita pada dasarnya memiliki variabel dependen kontinu write,, dan variabel penjelas dengan beberapa level,readcat ,. Di R bisa kita gunakanlm(write ~ readcat, hsb2)


1. Menghasilkan matriks kontras:

Ada empat level berbeda untuk variabel yang diurutkan readcat, jadi kita akan memiliki kontras.n-1=3

table(hsb2$readcat)

(28,40] (40,52] (52,64] (64,76] 
     22      93      55      30 

Pertama, mari kita cari uang, dan lihat fungsi R bawaan:

contr.poly(4)
             .L   .Q         .C
[1,] -0.6708204  0.5 -0.2236068
[2,] -0.2236068 -0.5  0.6708204
[3,]  0.2236068 -0.5 -0.6708204
[4,]  0.6708204  0.5  0.2236068

Sekarang mari kita membedah apa yang terjadi di bawah tenda:

scores = 1:4  # 1 2 3 4 These are the four levels of the explanatory variable.
y = scores - mean(scores) # scores - 2.5

y=[-1.5,-0,5,0,5,1.5]

seq_len (n) - 1=[0,1,2,3]

n = 4; X <- outer(y, seq_len(n) - 1, "^") # n = 4 in this case

[1-1.52.25-3.3751-0,50,25-0,12510,50,250,12511.52.253.375]

Apa yang terjadi disana? yang outer(a, b, "^")menaikkan elemen ake elemen b, sehingga kolom pertama dihasilkan dari operasi, , ( - 0.5 ) 0 , 0.5 0 dan 1.5 0 ; kolom kedua dari ( - 1.5 ) 1 , ( - 0.5 ) 1 , 0.5 1 dan 1.5 1 ; yang ketiga dari ( - 1.5 ) 2 = 2.25(-1.5)0(-0,5)00,501.50(-1.5)1(-0,5)10,511.51(-1.5)2=2.25 , , 0,5 2 = 0,25 dan 1,5 2 = 2,25 ; dan yang keempat, ( - 1.5 ) 3 = - 3.375 , ( - 0.5 ) 3 = - 0.125 , 0.5 3 = 0.125 dan 1.5 3 = 3.375 .(-0,5)2=0,250,52=0,251.52=2.25(-1.5)3=-3.375(-0,5)3=-0,1250,53=0,1251.53=3.375

Berikutnya kita melakukan ortonormal dekomposisi matriks ini dan mengambil representasi kompak dari Q ( ). Beberapa cara kerja fungsi yang digunakan dalam faktorisasi QR dalam R yang digunakan dalam posting ini dijelaskan lebih lanjut di sini .QRc_Q = qr(X)$qr

[-20-2.500,5-2.2360-4.5840,50,447200,50,894-0,9296-1.342]

... yang kami simpan hanya diagonal ( z = c_Q * (row(c_Q) == col(c_Q))). Apa yang terletak pada diagonal: Hanya "bawah" entri dari bagian dari Q R dekomposisi. Hanya? baik, tidak ... Ternyata diagonal dari matriks segitiga atas berisi nilai eigen dari matriks!RQR

Selanjutnya kita memanggil fungsi berikut:, raw = qr.qy(qr(X), z)hasil yang dapat direplikasi "secara manual" oleh dua operasi: 1. Mengubah bentuk kompak dari , yaitu , menjadi Q , transformasi yang dapat dicapai dengan , dan 2. Melakukan matriks perkalian Q z , seperti dalam .Qqr(X)$qrQQ = qr.Q(qr(X))QzQ %*% z

Krusial, mengalikan dengan nilai eigen dari R tidak mengubah orthogonality dari vektor-vektor kolom konstituen, tetapi mengingat bahwa nilai absolut dari nilai eigen muncul dalam urutan menurun dari kiri atas ke kanan bawah, perbanyakan Q z akan cenderung menurunkan nilai dalam kolom polinomial orde tinggi:QRQz

Matrix of Eigenvalues of R
     [,1]      [,2] [,3]      [,4]
[1,]   -2  0.000000    0  0.000000
[2,]    0 -2.236068    0  0.000000
[3,]    0  0.000000    2  0.000000
[4,]    0  0.000000    0 -1.341641

Membandingkan nilai-nilai dalam vektor kolom kemudian (kuadrat dan kubik) sebelum dan sesudah operasi faktorisasi, dan dengan tidak terpengaruh pertama dua kolom.QR

Before QR factorization operations (orthogonal col. vec.)
     [,1] [,2] [,3]   [,4]
[1,]    1 -1.5 2.25 -3.375
[2,]    1 -0.5 0.25 -0.125
[3,]    1  0.5 0.25  0.125
[4,]    1  1.5 2.25  3.375


After QR operations (equally orthogonal col. vec.)
     [,1] [,2] [,3]   [,4]
[1,]    1 -1.5    1 -0.295
[2,]    1 -0.5   -1  0.885
[3,]    1  0.5   -1 -0.885
[4,]    1  1.5    1  0.295

Akhirnya kita sebut (Z <- sweep(raw, 2L, apply(raw, 2L, function(x) sqrt(sum(x^2))), "/", check.margin = FALSE))mengubah matriks rawmenjadi vektor ortonormal :

Orthonormal vectors (orthonormal basis of R^4)
     [,1]       [,2] [,3]       [,4]
[1,]  0.5 -0.6708204  0.5 -0.2236068
[2,]  0.5 -0.2236068 -0.5  0.6708204
[3,]  0.5  0.2236068 -0.5 -0.6708204
[4,]  0.5  0.6708204  0.5  0.2236068

"/"col.xsaya2(saya) apply(raw, 2, function(x)sqrt(sum(x^2)))2 2.236 2 1.341(ii)(saya) .

R4contr.poly(4)

[-0,67082040,5-0,2236068-0,2236068-0,50,67082040,2236068-0,5-0,67082040,67082040,50,2236068]

(sum(Z[,3]^2))^(1/4) = 1z[,3]%*%z[,4] = 0skor - rata-rata123


2. Kontras (kolom) mana yang berkontribusi signifikan untuk menjelaskan perbedaan antar level dalam variabel penjelas?

Kami hanya dapat menjalankan ANOVA dan melihat ringkasan ...

summary(lm(write ~ readcat, hsb2))

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  52.7870     0.6339  83.268   <2e-16 ***
readcat.L    14.2587     1.4841   9.607   <2e-16 ***
readcat.Q    -0.9680     1.2679  -0.764    0.446    
readcat.C    -0.1554     1.0062  -0.154    0.877 

... untuk melihat bahwa ada efek linear readcataktif write, sehingga nilai asli (dalam potongan kode ketiga di awal posting) dapat direproduksi sebagai:

coeff = coefficients(lm(write ~ readcat, hsb2))
C = contr.poly(4)
(recovered = c(coeff %*% c(1, C[1,]),
               coeff %*% c(1, C[2,]),
               coeff %*% c(1, C[3,]),
               coeff %*% c(1, C[4,])))
[1] 42.77273 49.97849 56.56364 61.83333

... atau...

masukkan deskripsi gambar di sini

... atau jauh lebih baik ...

masukkan deskripsi gambar di sini

saya=1tSebuahsaya=0Sebuah1,,Sebuaht

masukkan deskripsi gambar di sini

Gagasan di balik kontras ortogonal adalah bahwa kesimpulan yang dapat kita tarik (dalam hal ini menghasilkan koefisien melalui regresi linier) akan menjadi hasil dari aspek independen dari data. Ini tidak akan terjadi jika kita hanya menggunakanX0,X1,.Xn sebagai kontras.

Secara grafis, ini jauh lebih mudah dipahami. Bandingkan sarana aktual oleh kelompok-kelompok dalam blok hitam persegi besar dengan nilai-nilai yang ditentukan sebelumnya, dan lihat mengapa perkiraan garis lurus dengan kontribusi minimal polinomial kuadratik dan kubik (dengan kurva yang hanya didekati dengan loess) adalah optimal:

masukkan deskripsi gambar di sini

Jika, hanya untuk efek, koefisien ANOVA telah sebesar untuk kontras linier untuk perkiraan lainnya (kuadrat dan kubik), plot nonsensik yang mengikuti akan menggambarkan lebih jelas plot polinomial dari setiap "kontribusi":

masukkan deskripsi gambar di sini

Kodenya ada di sini .

Antoni Parellada
sumber
+1 Wow. Bisakah jawaban ini (saya belum membacanya sampai akhir ini) dilihat sebagai jawaban untuk pertanyaan saya yang lama, lupa juga stats.stackexchange.com/q/63639/3277 ?
ttnphns
(+1) @ttnphns: Bisa dibilang lebih cocok di sana.
Scortchi
Hanya sebuah tip: Anda mungkin ingin berkomentar di sana dengan tautan ke sini; atau berikan jawaban di sana - yang kemungkinan akan saya terima.
ttnphns
1
@ttnphns dan @Scortchi Terima kasih! Saya menghabiskan cukup banyak waktu untuk mencoba memahami konsep-konsep ini, dan tidak berharap banyak reaksi. Jadi ini kejutan yang sangat positif. Saya pikir ada beberapa kerutan yang harus diselesaikan sehubungan dengan menjelaskan qr.qy()fungsi, tetapi saya pasti akan mencoba untuk melihat apakah saya dapat mengatakan sesuatu yang minimal koheren tentang pertanyaan Anda segera setelah saya punya waktu.
Antoni Parellada
1
@ Elvis saya memang mencoba memilih kalimat ringkasan yang bagus dan menempatkannya di suatu tempat di pos. Saya pikir ini adalah poin yang bagus, dan membutuhkan penjelasan matematis yang bagus, tetapi mungkin terlalu banyak untuk membahas lebih lanjut.
Antoni Parellada
5

Saya akan menggunakan contoh Anda untuk menjelaskan cara kerjanya. Menggunakan kontras polinomial dengan empat kelompok menghasilkan berikut.

Ewrsayate1=μ-0,67L.+0,5Q-0,22CEwrsayate2=μ-0,22L.-0,5Q+0,67CEwrsayate3=μ+0,22L.-0,5Q-0,67CEwrsayate4=μ+0,67L.+0,5Q+0,22C

Di mana persamaan pertama bekerja untuk kelompok skor bacaan terendah dan yang keempat untuk kelompok skor bacaan terbaik. kita dapat membandingkan persamaan ini dengan yang diberikan menggunakan regresi linear normal (seandainyareSebuahdsaya terus menerus)

Ewrsayatesaya=μ+reSebuahdsayaL.+reSebuahdsaya2Q+reSebuahdsaya3C

Biasanya bukan L.,Q,C kamu akan memiliki β1,β2,β3dan ditulis pada posisi pertama. Tetapi tulisan ini menyerupai yang dengan kontras polinomial. Jadi angka di depanL.,Q,C sebenarnya bukan reSebuahdsaya,reSebuahdsaya2,reSebuahdsaya3. Anda dapat melihat koefisien itu sebelumnyaL. memiliki tren linier, sebelumnya Q kuadrat dan sebelumnya C kubik.

Kemudian R memperkirakan parameter μ,L.,Q,C dan memberimu

μ^=52.79,L.^=14.26,Q^=-0,97,C^=-0,16
Dimana μ^=14saya=14Ewrsayatesaya dan estimasi koefisien μ^,L.^,Q^,C^adalah sesuatu seperti perkiraan pada regresi linier normal. Jadi dari output Anda dapat melihat apakah koefisien yang diperkirakan secara signifikan berbeda dari nol, sehingga Anda dapat mengantisipasi beberapa jenis tren linier, kuadratik atau kubik.

Dalam contoh itu secara signifikan bukan nol saja L.^. Jadi kesimpulan Anda bisa jadi: Kami melihat bahwa skor yang lebih baik dalam penulisan tergantung secara linear pada skor membaca, tetapi tidak ada efek kuadratik atau kubik yang signifikan.

Fimba
sumber