Какова мощность регрессионного теста F?

11

Классический F-тест для подмножеств переменных в полилинейной регрессии имеет вид где - сумма квадратов ошибок в «уменьшенной» модели, которая вложена в «большую» модель , а - степени свободы две модели. При нулевой гипотезе, что дополнительные переменные в «большой» модели не имеют линейной объяснительной силы, статистика распределяется как F с степенями свободы и .

F=(SSE(R)SSE(B))/(dfRdfB)SSE(B)/dfB,
SSE(R)BdfdfRdfBdfB

Какое распределение, однако, под альтернативой? Я предполагаю, что это нецентральный F (надеюсь, не вдвойне нецентральный), но я не могу найти ссылку на то, что именно является параметром нецентральности. Я собираюсь догадаться, что это зависит от истинных коэффициентов регрессии , и, вероятно, от проектной матрицы , но в этом я не уверен.βX

shabbychef
источник

Ответы:

9

Параметр нецентральности равен , проекция для ограниченной модели - , - вектор истинных параметров, - матрица проектирования для неограниченной (истинной) модели,это норма:δ2PrβX||x||

δ2=||XβPrXβ||2σ2

Вы можете прочитать формулу как это: вектор ожидаемых значений обусловливающие конструкции матрицы . Если вы рассматриваете как вектор эмпирических данных , то его проекция на подпространство ограниченной модели будет , что дает вам прогноз из ограниченной модели для этих «данных». Следовательно, аналогичен и дает вам ошибку этого прогноза. Отсюда дает сумму квадратов этой ошибки. Если ограниченная модель верна, тоE(y|X)=XβXXβyPrXβy^XβPrXβyy^||XβPrXβ||2Xβуже находится в подпространстве, определенном и , так что параметр нецентральности равен .XrPrXβ=Xβ0

Вы должны найти это в Мардии, Кенте и Бибби. (1980). Многомерный анализ.

каракал
источник
большой! должна ли норма быть в квадрате? В противном случае кажется, что единицы имеют значение? Вы утверждаете, что это «сумма квадратов», поэтому я думаю, что это норма в квадрате ..
Шаббычеф,
@shabbychef Конечно, вы правы, спасибо, что поймали это!
Каракал
7

Я подтвердил ответ @ Каракала экспериментом Монте-Карло. Я генерировал случайные экземпляры из линейной модели (со случайным размером), вычислял F-статистику и вычислял p-значение, используя нецентральный параметр затем я построил эмпирический cdf этих p-значений. Если параметр нецентральности (и код!) Верен, я должен получить почти одинаковый cdf, что имеет место:

δ2=||Xβ1Xβ2||2σ2,

эмпирический CDF о том, что должно быть нормальным

Вот код R (простите за стиль, я все еще учусь):

#sum of squares
sum2 <- function(x) { return(sum(x * x)) }
#random integer between n and 2n
rint <- function(n) { return(ceiling(runif(1,min=n,max=2*n))) }
#generate random instance from linear model plus noise.
#n observations of p2 vector
#regress against all variables and against a subset of p1 of them
#compute the F-statistic for the test of the p2-p1 marginal variables
#compute the p-value under the putative non-centrality parameter
gend <- function(n,p1,p2,sig = 1) {
 beta2 <- matrix(rnorm(p2,sd=0.1),nrow=p2)
 beta1 <- matrix(beta2[1:p1],nrow=p1)
 X <- matrix(rnorm(n*p2),nrow=n,ncol=p2)
 yt1 <- X[,1:p1] %*% beta1
 yt2 <- X %*% beta2
 y <- yt2 + matrix(rnorm(n,mean=0,sd=sig),nrow=n)
 ncp <- (sum2(yt2 - yt1)) / (sig ** 2)
 bhat2 <- lm(y ~ X - 1)
 bhat1 <- lm(y ~ X[,1:p1] - 1)
 SSE1 <- sum2(bhat1$residual)
 SSE2 <- sum2(bhat2$residual)
 df1 <- bhat1$df.residual
 df2 <- bhat2$df.residual
 Fstat <- ((SSE1 - SSE2) / (df1 - df2)) / (SSE2 / bhat2$df.residual)
 pval <- pf(Fstat,df=df1-df2,df2=df2,ncp=ncp)
 return(pval)
}
#call the above function, but randomize the problem size (within reason)
genr <- function(n,p1,p2,sig=1) {
 use.p1 <- rint(p1)
 use.p2 <- use.p1 + rint(p2 - p1)
 return(gend(n=rint(n),p1=use.p1,p2=use.p2,sig=sig+runif(1)))
}
ntrial <- 4096
ssize <- 256
z <- replicate(ntrial,genr(ssize,p1=4,p2=10))
plot(ecdf(z))
shabbychef
источник
2
+1 за продолжение с кодом. Всегда приятно это видеть.
mpiktas