Model linier: membandingkan daya prediksi dua metode pengukuran yang berbeda

9

Saya tertarik untuk memprediksi Ydan saya mempelajari dua teknik pengukuran yang berbeda X1dan X2. Bisa jadi misalnya saya ingin memprediksi rasa pisang, baik dengan mengukur berapa lama sudah tergeletak di atas meja, atau dengan mengukur jumlah bintik-bintik coklat pada pisang.

Saya ingin tahu mana salah satu teknik pengukuran yang lebih baik, haruskah saya memilih untuk melakukan hanya satu.

Saya dapat membuat model linier di R:

m1 = lm(Y ~ X1)
m2 = lm(Y ~ X2)

Sekarang katakanlah X1adalah prediktor unggul dari rasa pisang daripada X2. Saat menghitung dari dua model, model jelas lebih tinggi dari model . Sebelum menulis makalah tentang bagaimana metode lebih baik daripada , saya ingin memiliki semacam indikasi bahwa perbedaannya bukan karena kebetulan, mungkin dalam bentuk nilai-p.R2R2m1m2X1X2

Bagaimana orang akan melakukan ini? Bagaimana cara melakukannya ketika saya menggunakan berbagai merek pisang dan beralih ke model Linear Mixed Effect yang memasukkan merek pisang sebagai efek acak?

Rodin
sumber
Bisakah Anda menjelaskan mengapa Anda tidak bisa memasukkan kedua prediktor dalam model? Dalam kasus Anda, X1dan X2mungkin akan dikorelasikan, karena bintik-bintik coklat mungkin meningkat dengan meningkatnya waktu berbaring di atas meja.
COOLSerdash
Saya menarik dalam menguji apakah X1 atau X2 adalah metode pengukuran yang lebih baik. Jika memasukkan keduanya dalam satu model dapat menjawab pertanyaan itu, tidak ada masalah untuk melakukannya. Jelas keduanya berkorelasi karena mereka mengukur hal yang sama.
Rodin
Saya ingin mengatakan: ketika mencoba mengukur rasa pisang, mengukur berapa lama itu tergeletak di atas meja adalah cara yang lebih baik untuk menentukan ini daripada menghitung jumlah bintik-bintik cokelat (p <0,05).
Rodin

Jawaban:

7

Kemudian

Satu hal yang ingin saya tambahkan setelah mendengar bahwa Anda memiliki model efek campuran linier: dan masih dapat digunakan untuk membandingkan model. Lihat tulisan ini , misalnya. Dari pertanyaan serupa lainnya di situs, tampaknya makalah ini sangat penting.AIC,AICcBIC


Jawaban asli

Pada dasarnya yang Anda inginkan adalah membandingkan dua model yang tidak bersarang. Seleksi Model Burnham dan Anderson dan inferensi multimodel membahas hal ini dan merekomendasikan penggunaan , atau dll. Karena uji rasio kemungkinan tradisional hanya berlaku pada model bersarang. Mereka secara eksplisit menyatakan bahwa kriteria informasi-teoretis seperti dll. Bukan tes dan bahwa kata "signifikan" harus dihindari ketika melaporkan hasil.AICAICcBICAIC,AICc,BIC

Berdasarkan ini dan ini jawaban, saya sarankan ini pendekatan:

  1. Membuat matriks sebar (SPLOM) dari dataset Anda termasuk smoothers: pairs(Y~X1+X2, panel = panel.smooth, lwd = 2, cex = 1.5, col = "steelblue", pch=16). Periksa apakah garis (smoothers) kompatibel dengan hubungan linier. Sempurnakan model jika perlu.
  2. Hitung model m1dan m2. Lakukan beberapa pengecekan model (residu, dll.): plot(m1)Dan plot(m2).
  3. Hitung ( dikoreksi untuk ukuran sampel kecil) untuk kedua model dan hitung perbedaan absolut antara kedua s. The paket menyediakan fungsi untuk ini: . Jika perbedaan absolut ini lebih kecil dari 2, kedua model pada dasarnya tidak bisa dibedakan. Kalau tidak, pilih model dengan lebih rendah .AICcAICAICcR psclAICcabs(AICc(m1)-AICc(m2))AICc
  4. Hitung tes rasio kemungkinan untuk model non-bersarang. The R paketlmtest memiliki fungsi coxtest(uji Cox), jtest(uji Davidson-MacKinnon J) dan encomptest(uji meliputi dari Davidson & MacKinnon).

Beberapa pemikiran: Jika dua ukuran pisang benar - benar mengukur hal yang sama, keduanya mungkin sama - sama cocok untuk prediksi dan mungkin tidak ada model "terbaik".

Makalah ini mungkin juga bermanfaat.

Berikut adalah contoh dalam R:

#==============================================================================
# Generate correlated variables
#==============================================================================

set.seed(123)

R <- matrix(cbind(
  1   , 0.8 , 0.2,
  0.8 , 1   , 0.4,
  0.2 , 0.4 , 1),nrow=3) # correlation matrix
U <- t(chol(R))
nvars <- dim(U)[1]
numobs <- 500
set.seed(1)
random.normal <- matrix(rnorm(nvars*numobs,0,1), nrow=nvars, ncol=numobs);
X <- U %*% random.normal
newX <- t(X)
raw <- as.data.frame(newX)
names(raw) <- c("response","predictor1","predictor2")

#==============================================================================
# Check the graphic
#==============================================================================

par(bg="white", cex=1.2)
pairs(response~predictor1+predictor2, data=raw, panel = panel.smooth,
      lwd = 2, cex = 1.5, col = "steelblue", pch=16, las=1)

SPLOM

Para smoothers mengkonfirmasi hubungan linier. Ini dimaksudkan, tentu saja.

#==============================================================================
# Calculate the regression models and AICcs
#==============================================================================

library(pscl)

m1 <- lm(response~predictor1, data=raw)
m2 <- lm(response~predictor2, data=raw)

summary(m1)

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.004332   0.027292  -0.159    0.874    
predictor1   0.820150   0.026677  30.743   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.6102 on 498 degrees of freedom
Multiple R-squared:  0.6549,    Adjusted R-squared:  0.6542 
F-statistic: 945.2 on 1 and 498 DF,  p-value: < 2.2e-16

summary(m2)

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.01650    0.04567  -0.361    0.718    
predictor2   0.18282    0.04406   4.150 3.91e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.021 on 498 degrees of freedom
Multiple R-squared:  0.03342,   Adjusted R-squared:  0.03148 
F-statistic: 17.22 on 1 and 498 DF,  p-value: 3.913e-05

AICc(m1)
[1] 928.9961

AICc(m2)
[1] 1443.994

abs(AICc(m1)-AICc(m2))
[1] 514.9977

#==============================================================================
# Calculate the Cox test and Davidson-MacKinnon J test
#==============================================================================

library(lmtest)

coxtest(m1, m2)

Cox test

Model 1: response ~ predictor1
Model 2: response ~ predictor2
                Estimate Std. Error   z value  Pr(>|z|)    
fitted(M1) ~ M2   17.102     4.1890    4.0826 4.454e-05 ***
fitted(M2) ~ M1 -264.753     1.4368 -184.2652 < 2.2e-16 ***

jtest(m1, m2)

J test

Model 1: response ~ predictor1
Model 2: response ~ predictor2
                Estimate Std. Error t value  Pr(>|t|)    
M1 + fitted(M2)  -0.8298   0.151702  -5.470 7.143e-08 ***
M2 + fitted(M1)   1.0723   0.034271  31.288 < 2.2e-16 ***

The dari model pertama jelas lebih rendah dan jauh lebih tinggi.AICcm1R2

Penting: Dalam model linier dengan kompleksitas yang sama dan distribusi kesalahan Gaussian , dan harus memberikan jawaban yang sama (lihat posting ini ). Dalam model nonlinear , penggunaan untuk kinerja model (goodness of fit) dan pemilihan model harus dihindari: lihat posting ini dan makalah ini , misalnya.R2,AICBICR2

COOLSerdash
sumber
Ini adalah jawaban praktis yang menjelaskan dengan baik apa yang bisa dilakukan. Apa yang lebih mengejutkan adalah contoh di mana plot pencar, dan menghasilkan jawaban yang berbeda . R2AICc
Nick Cox
@NickCox Itu memang menarik! Jika modelnya linier dan kompleksitasnya sama, dan kesalahan memiliki distribusi Gaussian, dan harus memberikan jawaban yang sama. Yang benar adalah bahwa saya tidak tahu bagaimana menghasilkan contoh ad-hoc di mana scatterplots, dan mengarah pada jawaban yang berbeda (mungkin seseorang dapat membantu saya?). Maaf. R2,AICBICR2AICc
COOLSerdash
Anda tidak perlu meminta maaf; Saya melihatnya sebagai kabar baik setiap kali metode yang masuk akal berbeda menyiratkan kesimpulan yang sama!
Nick Cox
Luar biasa. Tes Cox persis seperti yang saya inginkan. Sayangnya, model saya adalah model linear-efek campuran yang dilengkapi dengan paket lme4, yang tidak didukung secara langsung oleh lmtestpaket tersebut. Sebelum saya terjun ke literatur yang menulis tentang tes mirip cox dengan LME, apakah ada yang tahu paket R yang tersedia untuk melakukannya?
Rodin
@Rodin Tidak, saya tidak tahu Rpaket apa pun yang bisa melakukan itu. Mungkin posting ini bisa memberi Anda panduan lebih lanjut.
COOLSerdash
3

Ada jawaban abad ke-19 yang bagus yang berisiko diabaikan di sini. Untuk membandingkan dua garis lurus yang berbeda, plot data dan garis yang cocok dan pikirkan apa yang Anda lihat. Sangat mungkin bahwa satu model akan jelas lebih baik, dan itu tidak berarti lebih tinggi . Sebagai contoh, ada kemungkinan bahwa model garis lurus secara kualitatif salah dalam satu atau lain kasus. Bahkan lebih baik, data dan kecocokan mungkin menyarankan model yang lebih baik. Jika kedua model muncul tentang sama-sama baik atau buruk, itu jawaban lain.R2

Contoh pisang mungkin bercanda di sini, tapi saya tidak berharap garis lurus cocok untuk bekerja sama sekali ....

Mesin inferensial yang dibuat oleh orang lain sebagai jawaban adalah sesuatu yang indah secara intelektual, tetapi kadang-kadang Anda tidak perlu palu godam canggih untuk memecahkan masalah. Kadang-kadang tampak bahwa siapa pun yang menerbitkan malam itu lebih gelap daripada siang akan selalu memiliki seseorang yang bertanya, "Apakah Anda mengujinya secara formal? Apa nilai P Anda?".

Nick Cox
sumber
+1, poin bagus. Saya telah secara eksplisit menyatakan bahwa seseorang tidak boleh menggunakan pengujian signifikan ketika berhadapan dengan dan sama. AIC
COOLSerdash
1
Selalu baik untuk mengambil langkah mundur dan mengingat ini jadi +1 tetapi saya bertanya-tanya apakah ini mungkin hanya kasus aneh di mana ini tidak terlalu berguna. Apakah benar-benar sangat mungkin bahwa satu model akan jelas lebih baik daripada yang lain ketika dua prediktor tentatif adalah ukuran yang berbeda dari hal yang sama sebagai lawan variabel yang berbeda secara substantif? Mengesampingkan pisang selama satu menit, pikirkan kuesioner yang sedikit berbeda atau penggaris vs. pengintai laser. Satu ukuran harus sangat kurang untuk non-linearitas muncul dalam satu kasus tetapi tidak yang lain.
Gala
2
Saya setuju, tetapi saya merasa senang untuk memperluas jawaban saya ke plot residu, dll. Lebih umum, cara apa pun untuk menilai ini dapat berjalan ke keberatan yang sama, misalnya tes signifikansi mungkin gagal untuk menolak kinerja yang sama, mungkin bisa dibilang identik. Setelah mempublikasikan perbandingan metode pengukuran yang berbeda dengan data nyata, saya menegaskan bahwa pilihan metode biasanya akan menemukan metode terbaik, tetapi jika dua metode tampak sama bagusnya, maka itulah jawabannya, bukan masalah yang tidak terpecahkan yang menyebabkan kesedihan metodologis. R2
Nick Cox
Saya setuju bahwa cara terbaik adalah dengan jelas menyajikan data dan kecocokan oleh model. Dengan cara ini pembaca dapat memutuskan sendiri apakah akan menerima pernyataan tentang seseorang yang lebih baik atau tidak. Namun, saya khawatir pengulas akan menuntut tes signifikan, murni karena reaksi spontan. Komentar Gaël tentang siang dan malam tidak begitu jauh.
Rodin
1
Itu saya pada malam dan siang ....
Nick Cox
2

Lakukan tes Cox untuk model yang tidak bersarang.

y <- rnorm( 10 )
x1 <- y + rnorm( 10 ) / 2
x2 <- y + rnorm( 10 )

lm1 <- lm( y ~ x1 )
lm2 <- lm( y ~ x2 )

library( lmtest )

coxtest( lm1, lm2 )
?coxtest

(Anda akan menemukan referensi untuk tes lain).

Lihat juga komentar ini dan pertanyaan ini . Khususnya, pertimbangkan untuk menggunakan AIC / BIC.

Januari
sumber