Jednym ze sposobów podejścia do tego pytania jest spojrzenie na to w odwrotny sposób: jak moglibyśmy zacząć od normalnie rozmieszczonych reszt i ustawić je tak, aby były heteroscedastyczne? Z tego punktu widzenia odpowiedź staje się oczywista: powiązać mniejsze reszty z mniejszymi przewidywanymi wartościami.
Aby to zilustrować, oto wyraźna konstrukcja.
Dane po lewej stronie są wyraźnie heteroscedastyczne względem dopasowania liniowego (pokazane na czerwono). Jest to zależne od reszty w porównaniu do przewidywanego wykresu po prawej stronie. Ale - z założenia - nieuporządkowany zbiór reszt jest prawie normalnie rozłożony, jak pokazuje ich histogram na środku. (Wartość p w teście normalności Shapiro-Wilka wynosi 0,60, uzyskana za pomocą R
polecenia shapiro.test(residuals(fit))
wydanego po uruchomieniu poniższego kodu.)
Tak też mogą wyglądać prawdziwe dane. Morał polega na tym, że heteroscedastyczność charakteryzuje związek między wielkością resztkową a przewidywaniami, podczas gdy normalność nie mówi nam nic o tym, jak resztki odnoszą się do czegokolwiek innego.
Oto R
kod tej konstrukcji.
set.seed(17)
n <- 256
x <- (1:n)/n # The set of x values
e <- rnorm(n, sd=1) # A set of *normally distributed* values
i <- order(runif(n, max=dnorm(e))) # Put the larger ones towards the end on average
y <- 1 + 5 * x + e[rev(i)] # Generate some y values plus "error" `e`.
fit <- lm(y ~ x) # Regress `y` against `x`.
par(mfrow=c(1,3)) # Set up the plots ...
plot(x,y, main="Data", cex=0.8)
abline(coef(fit), col="Red")
hist(residuals(fit), main="Residuals")
plot(predict(fit), residuals(fit), cex=0.8, main="Residuals vs. Predicted")
ncvTest
funkcji pakietu samochodowego doR
przeprowadzenia formalnego testu dla heteroskedastyczności. W przykładzie Whubera poleceniencvTest(fit)
daje wartość która jest prawie równa zero i dostarcza mocnych dowodów przeciwko stałej wariancji błędu (czego oczekiwano oczywiście).