Hayır, veriler heteroskedastik değildir (onları nasıl simüle ettiğinize göre). Testin 0 serbestlik derecesini fark ettiniz mi? Bu, burada bir şeylerin ters gittiğine dair bir ipucu. B-P testi, modeldeki kalıntıların karesini alır ve modeldeki yordayıcıların (veya belirlediğiniz diğer tahmin edicilerin) bu değerlerde önemli miktarda değişkenliği hesaba katıp katamayacağını test eder. Modelde yalnızca kesişme noktası bulunduğundan, tanım gereği herhangi bir değişkenliği hesaba katamaz.
Bir göz atın: http://en.wikipedia.org/wiki/Breusch-Pagan_test
Ayrıca, yardım (bptest)
bölümünü okuduğunuzdan emin olun. Bu, bazı şeyleri açıklığa kavuşturmaya yardımcı olmalı.
Burada yanlış giden bir şey, görünüşe göre bptest ()
işlevinin bu hatalı durumu test etmemesi ve küçük bir p değeri. Aslında, bptest ()
işlevinin altında yatan koda dikkatlice bakarsanız, esasen şu oluyor:
format.pval (pchisq (0,0), digits = 4)
"< 2.2e-16"
değerini verir. Dolayısıyla, pchisq (0,0)
0
döndürür ve bu format.pval tarafından "< 2.2e-16"
haline getirilir. ()
. Bir bakıma, hepsi doğrudur, ancak bu tür bir karışıklığı önlemek için bptest ()
içinde sıfır dfs'yi test etmek muhtemelen yardımcı olacaktır.
EDIT
Bu soruyla ilgili hala çok fazla kafa karışıklığı var. Belki de B-P testinin gerçekte ne yaptığını göstermenin gerçekten faydası vardır. İşte bir örnek. İlk olarak, homoskedastik olan bazı verileri simüle edelim. Sonra iki yordayıcıya sahip bir regresyon modeli uydururuz. Daha sonra bptest ()
işlevi ile BP testini gerçekleştiriyoruz.
library (lmtest) n <- 100 x1i <- rnorm (n) x2i <- rnorm (n) yi <- rnorm (n) mod <- lm (yi ~ x1i + x2i) bptest (mod)
Peki gerçekten ne oluyor? İlk olarak, regresyon modeline göre kalan kareleri alın. Daha sonra, orijinal modele dahil edilen tahmin edicilerdeki bu kare artıkları geri getirirken $ n \ times R ^ 2 $ alın ( bptest ()
işlevinin orijinal modeldeki aynı öngörücüleri kullandığını unutmayın, ancak, heteroskedisitenin diğer değişkenlerin bir fonksiyonu olduğundan şüpheleniliyorsa, burada başka yordayıcılar da kullanılabilir). B-P testi için test istatistiği budur. Eşit varyansın boş hipotezi altında, bu test istatistiği, testte kullanılan yordayıcıların sayısına eşit serbestlik derecelerine sahip bir ki-kare dağılımını takip eder (kesişme hariç). Öyleyse, aynı sonuçları alıp alamayacağımızı görelim:
e2 <- resid (mod) ^ 2bp <- summary (lm (e2 ~ x1i + x2i)) $ r.squared * nbppchisq (bp, df = 2, lower.tail = FALSE)
Evet, işe yarıyor. Şans eseri, yukarıdaki test önemli olabilir (simüle edilen veriler homoskedastik olduğundan bu bir Tip I hatasıdır), ancak çoğu durumda önemli olmayacaktır.