Quadratsummen vom Typ III

Ich habe ein lineares Regressionsmodell mit einer kategorischen Variable (männlich und weiblich) und einer kontinuierlichen Variable .$A$$B$

Ich habe Kontrastcodes in R mit eingerichtet options(contrasts=c("contr.sum","contr.poly")). Und jetzt habe ich Quadratsummen vom Typ III für , und deren Interaktion (A: B) mit .$A$$B$drop1(model, .~., test="F")

Ich bleibe dabei, wie die Quadratsummen für berechnet werden . Ich$B$ denke es ist sum((predicted y of the full model - predicted y of the reduced model)^2). Das reduzierte Modell würde so aussehen y~A+A:B. Aber wenn ich benutze predict(y~A+A:B), gibt R vorhergesagte Werte zurück, die mit den vorhergesagten Werten des vollständigen Modells übereinstimmen. Daher wären die Quadratsummen 0.

(Für die Quadratsummen von ich ein reduziertes Modell von verwendet , das dasselbe ist wie .)$A$y~B+A:By~A:B

Hier ist ein Beispielcode für zufällig generierte Daten:

A<-as.factor(rep(c("male","female"), each=5))
set.seed(1)
B<-runif(10)
set.seed(5)
y<-runif(10)

model<-lm(y~A+B+A:B)

options(contrasts = c("contr.sum","contr.poly"))

#type3 sums of squares
drop1(model, .~., test="F")
#or same result:
library(car)
Anova(lm(y~A+B+A:B),type="III")

#full model
predFull<-predict(model)

#Calculate sum of squares
#SS(A|B,AB)
predA<-predict(lm(y~B+A:B))
sum((predFull-predA)^2) 

#SS(B|A,AB) (???)
predB<-predict(lm(y~A+A:B))
sum((predFull-predB)^2) 
#Sums of squares should be 0.15075 (according to anova table)
#but calculated to be 2.5e-31

#SS(AB|A,B)
predAB<-predict(lm(y~A+B))
sum((predFull-predAB)^2)


#Anova Table (Type III tests)
#Response: y
#             Sum Sq Df F value Pr(>F)
#(Intercept) 0.16074  1  1.3598 0.2878
#A           0.00148  1  0.0125 0.9145
#B           0.15075  1  1.2753 0.3019
#A:B         0.01628  1  0.1377 0.7233
#Residuals   0.70926  6    

Ich habe Unterschiede in der Schätzung der Regressoren zwischen R 2.15.1 und SAS 9.2 festgestellt, aber nach der Aktualisierung von R auf die Version 3.0.1 waren die Ergebnisse dieselben. Zuerst empfehle ich Ihnen, R auf die neueste Version zu aktualisieren.

Sie verwenden den falschen Ansatz, weil Sie die Quadratsumme für zwei verschiedene Modelle berechnen, was zwei verschiedene Entwurfsmatrizen impliziert. Dies führt zu einer völlig anderen Schätzung der Regressoren, die von lm () zur Berechnung der vorhergesagten Werte verwendet werden (Sie verwenden Regressoren mit unterschiedlichen Werten zwischen den beiden Modellen). SS3 wird basierend auf einem Hypotesetest berechnet, wobei angenommen wird, dass alle Konditionierungsregressoren gleich Null sind, während der konditionierte Regressor gleich 1 ist. Für die Berechnungen verwenden Sie dieselbe Entwurfsmatrix, die zur Schätzung des vollständigen Modells verwendet wird, wie für den vollständig geschätzten Regressor Modell. Denken Sie daran, dass die SS3s nicht vollständig additiv sind. Dies bedeutet, dass Sie, wenn Sie die geschätzte SS3 summieren, das Modell SS (SSM) nicht erhalten.

Hier schlage ich eine R-Implementierung der Mathematik vor, die den GLS-Algorithmus implementiert, der zur Schätzung von SS3 und Regressoren verwendet wird.

Die von diesem Code generierten Werte sind genau die gleichen, die mit SAS 9.2 generiert wurden, wie für die Ergebnisse, die Sie in Ihrem Code angegeben haben, während der SS3 (B | A, AB) 0,167486 anstelle von 0,15075 beträgt. Aus diesem Grund empfehle ich erneut, Ihre R-Version auf die neueste verfügbare Version zu aktualisieren.

Hoffe das hilft :)

A<-as.factor(rep(c("male","female"), each=5))
set.seed(1)
B<-runif(10)
set.seed(5)
y<-runif(10)


# Create a dummy vector of 0s and 1s
dummy <- as.numeric(A=="male")

# Create the design matrix
R <- cbind(rep(1, length(y)), dummy, B, dummy*B)

# Estimate the regressors
bhat <- solve(t(R) %*% R) %*% t(R) %*% y
yhat <- R %*% bhat
ehat <- y - yhat

# Sum of Squares Total
# SST <- t(y)%*%y - length(y)*mean(y)**2
# Sum of Squares Error
# SSE <- t(ehat) %*% ehat
# Sum of Squares Model
# SSM <- SST - SSE

# used for ginv()
library(MASS)

# Returns the Sum of Squares of the hypotesis test contained in the C matrix
SSH_estimate <- function(C)
{
    teta <- C%*%bhat
    M <- C %*% ginv(t(R)%*%R) %*% t(C)
    SSH <- t(teta) %*% ginv(M) %*% teta
    SSH
}

# SS(A|B,AB)
# 0.001481682
SSH_estimate(matrix(c(0, 1, 0, 0), nrow=1, ncol=4))
# SS(B|A,AB)
# 0.167486
SSH_estimate(matrix(c(0, 0, 1, 0), nrow=1, ncol=4))
# SS(AB|A,B)
# 0.01627824
SSH_estimate(matrix(c(0, 0, 0, 1), nrow=1, ncol=4))