Użyję małych liter dla wektorów i wielkich liter dla matryc.
W przypadku modelu liniowego formy:
y=Xβ+ε
gdzie jest macierzą rangi , i przyjmujemy . n × ( k + 1 ) k + 1 ≤ n ε ∼ N ( 0 , σ 2 )Xn×(k+1)k+1≤nε∼N(0,σ2)
Możemy oszacować przez , ponieważ istnieje odwrotność .β^(X⊤X)−1X⊤yX⊤X
Teraz w przypadku ANOVA mamy do czynienia z tym, że nie ma już pełnej rangi. Oznacza to, że nie mamy i musimy zadowolić się uogólnionym odwrotnością .X(X⊤X)−1(X⊤X)−
Jednym z problemów korzystania z tej uogólnionej odwrotności jest to, że nie jest ona unikalna. Innym problemem jest to, że nie możemy znaleźć obiektywnego estymatora dla , ponieważ
β
β^=(X⊤X)−X⊤y⟹E(β^)=(X⊤X)−X⊤Xβ.
Dlatego nie możemy oszacować . Ale czy możemy oszacować liniową kombinację ?ββ
Mamy liniową kombinację , powiedzmy , można oszacować, jeśli istnieje wektor taki, że .βg⊤βaE(a⊤y)=g⊤β
W kontrasty są szczególnym przypadku funkcji oszacowania, w których suma współczynników jest równa zeru.g
I pojawiają się kontrasty w kontekście predyktorów jakościowych w modelu liniowym. (jeśli sprawdzisz instrukcję połączoną przez @amoeba, zobaczysz, że wszystkie ich kodowanie kontrastu jest powiązane ze zmiennymi kategorialnymi). Następnie, odpowiadając na @Curious i @amoeba, widzimy, że powstają one w ANOVA, ale nie w „czystym” modelu regresji z tylko ciągłymi predyktorami (możemy również mówić o kontrastach w ANCOVA, ponieważ mamy w sobie pewne zmienne kategoryczne).
Teraz w modelu gdzie nie ma pełnej rangi, a , funkcja liniowa jest możliwa do oszacowania, jeśli istnieje wektor taki, że . Oznacza to, że jest liniową kombinacją wierszy . Ponadto istnieje wiele opcji wektora , na przykład , jak widać w poniższym przykładzie.
y=Xβ+ε
XE(y)=X⊤βg⊤βaa⊤X=g⊤g⊤Xaa⊤X=g⊤
Przykład 1
Rozważmy model jednokierunkowy:
yij=μ+αi+εij,i=1,2,j=1,2,3.
X=⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢111111111000000111⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥,β=⎡⎣⎢μτ1τ2⎤⎦⎥
Załóżmy, że , więc chcemy oszacować .g⊤=[0,1,−1][0,1,−1]β=τ1−τ2
Widzimy, że istnieją różne opcje wektora które dają : take ; lub ; lub .aa⊤X=g⊤a⊤=[0,0,1,−1,0,0]a⊤=[1,0,0,0,0,−1]a⊤=[2,−1,0,0,1,−2]
Przykład 2
Weźmy dwukierunkowy model:
.
yij=μ+αi+βj+εij,i=1,2,j=1,2
X=⎡⎣⎢⎢⎢11111100001110100101⎤⎦⎥⎥⎥,β=⎡⎣⎢⎢⎢⎢⎢⎢μα1α2β1β2⎤⎦⎥⎥⎥⎥⎥⎥
Możemy zdefiniować przewidywalne funkcje, biorąc liniowe kombinacje wierszy .X
Odejmowanie wiersza 1 od wierszy 2, 3 i 4 (z ):
X
⎡⎣⎢⎢⎢1000−10−1−10011−1−1−0−10101⎤⎦⎥⎥⎥
I biorąc rzędy 2 i 3 z czwartego rzędu:
⎡⎣⎢⎢⎢1000−10−1−00010−1−1−0−00100⎤⎦⎥⎥⎥
Pomnożenie tego przez daje:
β
g⊤1βg⊤2βg⊤3β=μ+α1+β1=β2−β1=α2−α1
Mamy więc trzy liniowo niezależne, estymowalne funkcje. Teraz tylko i można uznać za kontrasty, ponieważ suma jego współczynników (lub wiersza suma odpowiedniego wektora ) jest równa zero.g⊤2βg⊤3βg
Wracając do modelu zrównoważonego w jedną stronę
yij=μ+αi+εij,i=1,2,…,k,j=1,2,…,n.
Załóżmy, że chcemy przetestować hipotezę .H0:α1=…=αk
W tym ustawieniu macierz nie ma pełnej rangi, więc nie jest unikalny i nie do oszacowania. Aby to oszacować, możemy pomnożyć przez , o ile . Innymi słowy, można iff .Xβ=(μ,α1,…,αk)⊤βg⊤∑igi=0∑igiαi∑igi=0
Dlaczego to prawda?
Wiemy, że można jeśli istnieje wektor taki, że . Biorąc różne wiersze i , a następnie:
g⊤β=(0,g1,…,gk)β=∑igiαiag⊤=a⊤XXa⊤=[a1,…,ak]
[0,g1,…,gk]=g⊤=a⊤X=(∑iai,a1,…,ak)
Wynik jest następujący.
Jeśli chcielibyśmy przetestować konkretny kontrast, nasza hipoteza to . Na przykład: , którą można zapisać jako , więc porównujemy ze średnią i .H0:∑giαi=0H0:2α1=α2+α3H0:α1=α2+α32α1α2α3
Hipotezę tę można wyrazić jako , gdzie . W tym przypadku i testujemy tę hipotezę za pomocą następującej statystyki:
H0:g⊤β=0g⊤=(0,g1,g2,…,gk)q=1
F=[g⊤β^]⊤[g⊤(X⊤X)−g]−1g⊤β^SSE/k(n−1).
Jeśli jest wyrażone jako gdzie wiersze macierzy
są wzajemnie ortogonalnymi kontrastami ( ), a następnie możemy przetestować przy użyciu statystyki , gdzieH0:α1=α2=…=αkGβ=0
G=⎡⎣⎢⎢⎢⎢⎢g⊤1g⊤2⋮g⊤k⎤⎦⎥⎥⎥⎥⎥
g⊤igj=0H0:Gβ=0F=SSHrank(G)SSEk(n−1)SSH=[Gβ^]⊤[G(X⊤X)−1G⊤]−1Gβ^.
Przykład 3
Aby to lepiej zrozumieć, i załóżmy, że chcemy przetestować które można wyrazić jako
k=4H0:α1=α2=α3=α4,
H0:⎡⎣⎢α1−α2α1−α3α1−α4⎤⎦⎥=⎡⎣⎢000⎤⎦⎥
Lub, jako :
H0:Gβ=0
H0:⎡⎣⎢000111−1−0−0−0−1−1−0−0−1⎤⎦⎥G,our contrast matrix⎡⎣⎢⎢⎢⎢⎢⎢μα1α2α3α4⎤⎦⎥⎥⎥⎥⎥⎥=⎡⎣⎢000⎤⎦⎥
Widzimy więc, że trzy rzędy naszej matrycy kontrastu są określone przez współczynniki kontrastów będących przedmiotem zainteresowania. Każda kolumna podaje poziom współczynnika, którego używamy w naszym porównaniu.
Prawie wszystko, co napisałem, zostało wzięte \ skopiowane (bezwstydnie) z Rencher & Schaalje, „Modele liniowe w statystyce”, rozdziały 8 i 13 (przykłady, sformułowanie twierdzeń, niektóre interpretacje), ale inne rzeczy, takie jak termin „matryca kontrastowa” „(które tak naprawdę nie pojawia się w tej książce) i podana tutaj definicja była moja.
Powiązanie macierzy kontrastu OP z moją odpowiedzią
Jedna z macierzy OP (które można również znaleźć w tym podręczniku ) to:
> contr.treatment(4)
2 3 4
1 0 0 0
2 1 0 0
3 0 1 0
4 0 0 1
W tym przypadku nasz współczynnik ma 4 poziomy i możemy napisać model w następujący sposób: Można to zapisać w postaci macierzy jako:
⎡⎣⎢⎢⎢y11y21y31y41⎤⎦⎥⎥⎥=⎡⎣⎢⎢⎢⎢μμμμ⎤⎦⎥⎥⎥⎥+⎡⎣⎢⎢⎢a1a2a3a4⎤⎦⎥⎥⎥+⎡⎣⎢⎢⎢ε11ε21ε31ε41⎤⎦⎥⎥⎥
Lub
⎡⎣⎢⎢⎢y11y21y31y41⎤⎦⎥⎥⎥=⎡⎣⎢⎢⎢11111000010000100001⎤⎦⎥⎥⎥X⎡⎣⎢⎢⎢⎢⎢⎢μa1a2a3a4⎤⎦⎥⎥⎥⎥⎥⎥β+⎡⎣⎢⎢⎢ε11ε21ε31ε41⎤⎦⎥⎥⎥
Teraz, dla przykładu kodowania fikcyjnego w tej samej instrukcji, używają jako grupy referencyjnej. W ten sposób odejmujemy wiersz 1 od każdego innego wiersza w macierzy , co daje :a1XX˜
⎡⎣⎢⎢⎢1000−1−1−1−1010000100001⎤⎦⎥⎥⎥
Jeśli zaobserwujesz numerację wierszy i kolumn w macierzy contr.treatment (4), zobaczysz, że uwzględniają one wszystkie wiersze i tylko kolumny związane z czynnikami 2, 3 i 4. Jeśli zrobimy to samo w powyższa macierz daje:
⎡⎣⎢⎢⎢010000100001⎤⎦⎥⎥⎥
W ten sposób macierz contr.treatment (4) mówi nam, że porównują czynniki 2, 3 i 4 ze współczynnikiem 1 i porównują czynnik 1 ze stałą (oto moje rozumienie powyższego).
I, definiując (tj. Biorąc tylko wiersze, które sumują się do 0 w powyższej macierzy):
G
⎡⎣⎢000−1−1−1100010001⎤⎦⎥
Możemy przetestować i znaleźć szacunki kontrastów.H0:Gβ=0
hsb2 = read.table('http://www.ats.ucla.edu/stat/data/hsb2.csv', header=T, sep=",")
y<-hsb2$write
dummies <- model.matrix(~factor(hsb2$race)+0)
X<-cbind(1,dummies)
# Defining G, what I call contrast matrix
G<-matrix(0,3,5)
G[1,]<-c(0,-1,1,0,0)
G[2,]<-c(0,-1,0,1,0)
G[3,]<-c(0,-1,0,0,1)
G
[,1] [,2] [,3] [,4] [,5]
[1,] 0 -1 1 0 0
[2,] 0 -1 0 1 0
[3,] 0 -1 0 0 1
# Estimating Beta
X.X<-t(X)%*%X
X.y<-t(X)%*%y
library(MASS)
Betas<-ginv(X.X)%*%X.y
# Final estimators:
G%*%Betas
[,1]
[1,] 11.541667
[2,] 1.741667
[3,] 7.596839
Szacunki są takie same.
Odnosząc odpowiedź @ttnphns do mojej.
W pierwszym przykładzie konfiguracja ma kategoryczny współczynnik A mający trzy poziomy. Możemy to zapisać jako model (załóżmy dla uproszczenia, że ):
j=1
yij=μ+ai+εij,for i=1,2,3
Załóżmy, że chcemy przetestować lub , przy czym jest naszą grupą odniesienia / czynnikiem.H0:a1=a2=a3H0:a1−a3=a2−a3=0a3
Można to zapisać w formie macierzy jako:
⎡⎣⎢y11y21y31⎤⎦⎥=⎡⎣⎢μμμ⎤⎦⎥+⎡⎣⎢a1a2a3⎤⎦⎥+⎡⎣⎢ε11ε21ε31⎤⎦⎥
Lub
⎡⎣⎢y11y21y31⎤⎦⎥=⎡⎣⎢111100010001⎤⎦⎥X⎡⎣⎢⎢⎢μa1a2a3⎤⎦⎥⎥⎥β+⎡⎣⎢ε11ε21ε31⎤⎦⎥
Teraz, jeśli odejmiemy wiersz 3 od rzędu 1 i rzędu 2, otrzymamy, że stanie się (nazywam to :XX˜
X˜=⎡⎣⎢001100010−1−1−1⎤⎦⎥
Porównaj 3 ostatnie kolumny powyższej macierzy z macierzą @ttnphns . Pomimo kolejności są dość podobne. Rzeczywiście, jeśli pomnożymy , otrzymamy:LX˜β
⎡⎣⎢001100010−1−1−1⎤⎦⎥⎡⎣⎢⎢⎢μa1a2a3⎤⎦⎥⎥⎥=⎡⎣⎢a1−a3a2−a3μ+a3⎤⎦⎥
Mamy więc przewidywalne funkcje: ; ; .c⊤1β=a1−a3c⊤2β=a2−a3c⊤3β=μ+a3
Ponieważ , z powyższego wynika, że porównujemy naszą stałą ze współczynnikiem dla grupy odniesienia (a_3); współczynnik grupy 1 do współczynnika grupy 3; oraz współczynnik grupy 2 do grupy3. Lub, jak powiedział @ttnphns: „Od razu widzimy, zgodnie ze współczynnikami, że oszacowana stała będzie równa średniej Y w grupie odniesienia; ten parametr b1 (tj. Zmiennej obojętnej A1) będzie równy różnicy: średnia Y w grupie 1 minus Średnia Y w grupie 3; parametr b2 jest różnicą: średnia w grupie 2 minus średnia w grupie 3. ”H0:c⊤iβ=0
Ponadto zauważ, że (zgodnie z definicją kontrastu: funkcja przewidywalna + suma wiersza = 0), że wektory i są kontrastami. A jeśli utworzymy macierz kontrastów, otrzymamy:c1c2G
G=[001001−1−1]
Nasza matryca kontrastu do testowaniaH0:Gβ=0
Przykład
Użyjemy tych samych danych, co „Przykład kontrastu zdefiniowany przez użytkownika” @ttnphns (chciałbym wspomnieć, że teoria, którą tu napisałem, wymaga kilku modyfikacji w celu uwzględnienia modeli z interakcjami, dlatego wybrałem ten przykład. Jednak , definicje kontrastów i - jak to nazywam - matryca kontrastu pozostają takie same).
Y<-c(0.226,0.6836,-1.772,-0.5085,1.1836,0.5633,0.8709,0.2858,0.4057,-1.156,1.5199,
-0.1388,0.4865,-0.7653,0.3418,-1.273,1.4042,-0.1622,0.3347,-0.4576,0.7585,0.4084,
1.4165,-0.5138,0.9725,0.2373,-1.562,1.3985,0.0397,-0.4689,-1.499,-0.7654,0.1442,
-1.404,-0.2201,-1.166,0.7282,0.9524,-1.462,-0.3478,0.5679,0.5608,1.0338,-1.161,
-0.1037,2.047,2.3613,0.1222)
F_<-c(1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,4,4,
5,5,5,5,5,5,5,5,5,5,5)
dummies.F<-model.matrix(~as.factor(F_)+0)
X_F<-cbind(1,dummies.F)
G_F<-matrix(0,4,6)
G_F[1,]<-c(0,3,3,-2,-2,-2)
G_F[2,]<-c(0,1,-1,0,0,0)
G_F[3,]<-c(0,0,0,2,-1,-1)
G_F[4,]<-c(0,0,0,0,1,-1)
G
[,1] [,2] [,3] [,4] [,5] [,6]
[1,] 0 3 3 -2 -2 -2
[2,] 0 1 -1 0 0 0
[3,] 0 0 0 2 -1 -1
[4,] 0 0 0 0 1 -1
# Estimating Beta
X_F.X_F<-t(X_F)%*%X_F
X_F.Y<-t(X_F)%*%Y
Betas_F<-ginv(X_F.X_F)%*%X_F.Y
# Final estimators:
G_F%*%Betas_F
[,1]
[1,] 0.5888183
[2,] -0.1468029
[3,] 0.6115212
[4,] -0.9279030
Mamy więc te same wyniki.
Wniosek
Wydaje mi się, że nie ma jednego definiującego pojęcia, czym jest matryca kontrastu.
Jeśli weźmiesz definicję kontrastu podaną przez Scheffe („Analiza wariancji”, strona 66), zobaczysz, że jest to przewidywalna funkcja, której współczynniki sumują się do zera. Jeśli więc chcemy przetestować różne liniowe kombinacje współczynników naszych zmiennych kategorialnych, używamy macierzy . Jest to macierz, w której rzędy sumują się do zera, której używamy do pomnożenia naszej macierzy współczynników, aby umożliwić oszacowanie tych współczynników. Jego rzędy wskazują różne liniowe kombinacje kontrastów, które testujemy, a kolumny wskazują, które czynniki (współczynniki) są porównywane.G
Ponieważ powyższa macierz jest skonstruowana w taki sposób, że każdy jej wiersz składa się z wektora kontrastu (który sumuje się do 0), dla mnie sensowne jest nazywanie „macierzą kontrastu” ( Monahan - „Elementarz w modelach liniowych” - również używa tej terminologii).GG
Jednak, jak to pięknie wyjaśnia @ttnphns, oprogramowanie nazywa coś innego „matrycą kontrastu” i nie mogłem znaleźć bezpośredniego związku między macierzą a wbudowanymi poleceniami / macierzami SPSS (@ttnphns ) lub R (pytanie OP), tylko podobieństwa. Ale wierzę, że miła dyskusja / współpraca przedstawione tutaj pomogą wyjaśnić takie pojęcia i definicje.G
