Transcript Lecture 3: Resemblance Between Relatives

Matrix Algebra
and
Linear Models
Matrices: An array of elements
Vectors
0
1written in bold lowercase
Usually
12
a = @13 A
(3 x 1)
47
b = ( 2 0 5 21)
Column vector Row vector (1 x 4)
0
1
0 bold1 uppercase
General Matrices Usually
written
in
3 1 2
0 1
C = @2 5 4 A
1 1 2
(3 x 3)
Square matrix
D = @3 4 A
2 9
(3 x 2)
Dimensionality of a matrix: r x c (rows x columns)
think of Railroad Car
A matrix is defined by a list of its elements.
B has ij-th element Bij -- the element in row i
and column j
Addition and Subtraction of Matrices
If two matrices have the same dimension (same number
of rows and columns), then matrix addition and subtraction
simply follows by adding (or subtracting) on an element by
element basis
Matrix addition: (A+B)ij = A ij + B
ij
Matrix subtraction: (A-B)ij = A ij - B ij
Examples:
µ
A =
3 0
1 2
µ
C= A+ B =
∂
µ
and B =
4 2
3 3
∂
1 2
2 1
∂
and D = A °- B =
µ
2 °- 2
°- 1
1
∂
Partitioned Matrices
It will often prove useful to divide (or partition) the
elements of a matrix into a matrix whose elements are
itself matrices.
0
..
0
1
3
.
B…
3 1 2
B ¢¢¢ ¢¢¢
C = @2 5 4 A = B
..
B 2
.
@
1 1 2
..
1
.
1
1
2
C
…
…
¢¢¢ ¢¢¢C
5
1
µ
a = ( 3) ;
b = ( 1 2) ;
d=
µ
C=
4 C
A
2
∂
2
;
1
a b
d B
µ
B =
∂
5 4
1 2
∂
• One
A row
vector
whose
elements
are column
vectors
partition
a whose
matrix
in many
ways
A can
column
vector
elements
are row
vectors
0 1
0 1
0 1
1 3
1
2
b
1
@
A
@
A
@
where
c1 = 2 ; c2 =
5 ; c3 =
4A
C = @b 2 A 1
1
2
b3
0
C = ( c1
c2
c3 )
b1 = ( 3 1
2);
b2 = ( 2 5
4);
b3 = ( 1 1
2)
Towards Matrix Multiplication: dot products
The dot (or inner) product of two vectors (both of
length n) is defined as follows:
.
a ¢b =
Xn
ai bi
i=1
Example:
.
0 11
2C
a= B
@3 A
4
and b = ( 4 5 7 9)
a b = 1*4 + 2*5 + 3*7 + 4*9 = 60
Towards Matrix Multiplication: Systems of Equations
The
least-squares
solution
for the
Matrix
multiplication
is defined
in linear
such amodel
way so as
to compactly represent and…solve linear systems of
y = π + Ø1x 1 + ¢¢¢Øn x n
equations.
yields the following system of equations for the bi
æ(y; z1 ) = Ø1 æ2 (z1 )
+ Ø2 æ(z1; z2) + ¢¢¢+ Øn æ(z1 ; zn )
æ(y; z2 ) = Ø1 æ(z1 ; z2 ) + Ø2 æ2(z2 ) +
..
..
..
.
.
.
æ(y; zn )= Ø1æ(z1 ; zn ) + Ø2 æ(z2; zn ) +
¢¢¢+ Øn æ(z2 ; zn )
..
..
.
.
¢¢¢ + Øn æ2 (zn )
This can be more compactly written in matrix form as
0
æ2 (z1 )
B æ(z1 ; z2 )
B
..
@
.
æ(z1 ; z2 )
æ2 (z2 )
..
.
æ(z1 ; zn ) æ(z2 ; zn )
: : : æ(z1 ; zn ) 1 0 Ø1
: : : æ(z2 ; zn ) C B Ø2
CB .
..
..
A @ ..
.
.
: : : æ2 (zn )
Øn
1
0 æ(y; z ) 1
1
C
B æ(y; z2 ) C
C = B
C
..
A @
A
.
æ(y; zn )
Matrix Multiplication: Quick overview
The order in which matrices are multiplied effects
the matrix product, e.g. AB = BA
For the product of two matrices to exist, the matrices
must conform. For AB, the number of columns of A must
equal the number of rows of B.
The matrix C = AB has the same number of rows as A
and the same number of columns as B.
C(rxc) = A(rxk) B(kxc)
Elements
in the
Elements
in the
ith
ij-th element of C is given by
column A
of B
rowjth
of matrix
indices
dimensions of
Inner
indices
mustgiven
match
XkOuter
Ci j columns
= resulting
A i lof
B l jA =matrix,
with
rows of
B r rows (A)
l = and
c columns (B)
1
More formally, consider the product L = MN
Express the matrix M as a column vector of row vectors
0 m 1
1
B m2 C
C
M = B
@ ... A
where
…
mi = ( M i1 M i2
¢¢¢ M i c )
mr
Likewise express N as a row vector of
column vectors
…
N = (n
n
¢¢¢ n )
where
1
2
b
0 N 1j
B N 2j
nj = B
@ ...
The ij-th element of L is the inner product
of M’s row i with N’s column j
…
.
0 m .¢n
m 1 ¢n 2 ¢¢¢ m 1 .¢n b 1
1
1
…
.
.
B m 2 .¢n1 m 2 ¢n 2 ¢¢¢ m 2 ¢n b C
C
L= B
.
.
.
.
.
.
.
.
@ .
A
.
.
.
… m ¢n
.
m .¢n
m .¢n
¢¢¢
r
1
r
2
r
b
N cj
1
C
C
A
Examples
µ
AB =
(
a b
c d
∂ µ
)(
e f
g h
∂
)=(
µ
ae + bg af + bh
ce+ dg cf + dh
Likewise
µ
BA =
(
ae + cf
ga + ch
eb + df
gd + dh
∂
)
∂
)
The Transpose of a Matrix
The transpose of a matrix exchanges the
rows and columns, ATij = Aji
Useful identities
(AB)T = BT AT
(ABC)T = CT BT AT
0 a 1
1
.
a = @ .. A
an
0 b 1
1
.
b = @ .. A
bn
Inner product = aTb = aT(1 X n) b (n X 1)
0 b 1
1
Xn
.
T
¢¢¢ an ) @ .. A = a b =
ai bi
…
Note that bTa = (aTb)T = aTb
i= 1
bn
Indices
Dimension
match,
of resulting
matricesproduct
conform
T
T
Outer product = ab =isa 1(nXX 1)
1b
(i.e.(1 Xa n)scalar)
( a1
0 a1 1
@ ... A ( b1
an
0 a1 b1 a1 b2
a2 b1 a2 b2
… bn ) = ab T = B
B .
¢¢¢
Resulting product
@ .. is ...an
an b1 an b2
¢¢¢ a1 bn 1
¢¢¢ a2 bn C
.. C
..
nxn
matrix
.
. A
¢¢¢ an bn
The Identity Matrix, I
The identity matrix serves the role of the
number 1 in matrix multiplication: AI =A, IA = A
I is a diagonal matrix, with all diagonal elements
being one, all off-diagonal elements zero.
Iij =
1 for i = j
0 otherwise
0
I3x 3
1 0
= @0 1
0 0
1
0
0A
1
The Inverse Matrix,
-1
A
For a square matrix A, define the Inverse of A, A-1, as
the matrix satisfying
-1
A A = AA
For
µ
A =
(
a b
c d
∂
A
° 1
-1
= I
1
=
ad ° bc
µ
d °b
°c a
If det(A) is not zero, A-1 exists and A is said to be
non-singular. If det(A) = 0, A is singular, and no
If this quantity (the determinant)
unique inverse exists (generalized inverses do)
is zero, the inverse does not exist.
∂
Useful identities
(AB)-1 = B-1 A-1
(AT)-1 = (A-1)T
Key use of inverses --- solving systems of equations
Ax = c
Vector
of known
constants
Matrix
of known
constants
wish
to=solve
for
-1
A Ax
Ix = x
Hence, x = A-1 c
unknowns
we
-1 c.
A-1Vector
Ax = Aof
Note that
0 æ2 (z1 )
B æ(z1 ; z2 )
B
..
@
.
æ(z1; z2 )
æ2 (z2 )
..
.
æ(z1 ; zn ) æ(z2 ; zn )
0 æ(y; z1 ) 1
: : : æ(z1 ; zn ) 1 0 Ø1 1
: : : æ(z2 ; zn ) C B Ø2 C B æ(y; z2 ) C
CB . C = B
C
..
.
..
A @ .. A @ ..
A
.
.
: : : æ2 (zn )
Øn
æ(y; zn )
Vb = c ---> b = V-1 c
Example: solve the OLS for b in y = a + b1z1 + b2z2 + e
b = V-1 c
0
1
æ(y; z1)
A
c= @
æ(y; z2)
0
æ2 (z1)
V = @
æ(z1 ; z2 )
æ(z1; z2)
1
A
æ2(z2 )
It is more compact to use æ(z1 ; z2 ) = Ω12 æ(z1 )æ(z2)
V
0
°- 1
0
2
æ (z2 )
1
@
= 2
2
2
æ (z1)æ (z2) (1 ° Ω12) °- æ(z ; z )
1 2
-° æ(z1 ; z2) 1
A
æ2(z1 )
1
0
1 0
1
Ø1
æ2 (z2 )
° æ(z1; z2 )
æ(y; z1 )
1
@ A = 2
@
A@
A
2
2
æ (z1 )æ (z2 ) (1° Ω12 )
Ø2
° æ(z1; z2 )
æ2 (z1 )
æ(y; z2 )
∑
∏
1
æ(y; z1 )
æ(y; z2 )
Ø1 = - 2
°- Ω12
2
1 ° Ω12 æ (z1 )
æ(z1)æ(z2 )
∑
∏
1
æ(y; z2 ) æ(y; z1 )
Ø2 = - 2
° Ω12
1 ° Ω12 æ2(z2 )
æ(z1)æ(z2 )
If r12 = 0, these reduce to the two univariate slopes,
æ(y; z1 )
Ø1 = 2
æ (z1 )
æ(y; z2 )
and Ø2 = 2
æ (z2 )
Quadratic and Bilinear
Forms
Quadratic product: for An x n and xn x 1
T
x Ax =
Xn Xn
ai j x i x j
Scalar (1 x 1)
i= 1 j = 1
Bilinear Form (generalization of quadratic product)
for Am x n, an x 1, bm x1 their bilinear form is bT1 x m Am x n an x 1
T
b Aa =
Xm Xn
A i j bi aj
i= 1 j = 1
Note that bTA a = aTAT b
Covariance Matrices for
Transformed Variables
What is the variance of the linear combination,
c1x1 + c2x2 + … + cnxn ? (note this is a scalar)
2° T
√
(
¢
( )
2
æ c x = æ
=
!
n
X
ci x i
0
n
i= 1
j=1
1
X
X
@
=æ
ci xi ;
cj xj A
i= 1
Xn Xn
n
æ(ci x i ; cj x j ) =
i= 1 j = 1
T
Xn Xn
ci cj æ(x i ; x j )
i= 1 j = 1
= c Vc
Likewise, the covariance between two linear combinations
can be expressed as a bilinear form,
T
T
T
æ(a x ; b x) = a V b
Now suppose we transform one vector of random
variables into another vector of random variables
Transform x into
(i) yk x 1 = Ak x n xn x 1
(ii) zm x 1 = Bm x n xn x 1
The covariance between the elements of these
two transformed vectors of the original is a
k x m covariance matrix = AVBT
For example, the covariance between yi and yj
is given by the ij-th element of AVAT
Likewise, the covariance between yi and zj
is given by the ij-th element of AVBT
The Multivariate Normal
Distribution (MVN)
Consider the pdf for n independent normal
random variables, the ith of which has mean
mi and variance s2i
p(x) =
Yn
2∂
(xi ° πi )
-° =
(2º ) 1 2 æi exp -°
2
°- 1
i= 1
-
= (2º ) ° n =2
√
! °- 1
n
Y
i=1
µ
æi
2æi
√
n
X (xi °- πi ) 2
exp -°
2æi2
i= 1
!
This can be expressed more compactly in matrix form
Define the covariance matrix V for the vector x of
the n normal random variable by
0 æ2
1
B 0
V = B .
@ ..
0
0
æ22
..
.
…
¢¢¢ 0 1
Fun fact! For
a
n
Y 2
diagonal
matrix
jV j =
æi D, the
Det = | D i|= 1= product
of the diagonal elements
¢¢¢ 0 C
.. C
..
. . A
…
¢¢¢ ¢¢¢ æn2
Define the mean vector m by
0 π 1
1
B π2 C
π= B . C
@ .. A
Xn (xi °- πi )2
- π) T V °- 1 (x °- π)
=
(x
°
2
æ
i
i=1
Hence in matrix from the MVN pdf becomes
πn
p(x) = (2º )
-° n =2
jV j
°- 1=2
∑
∏
1
exp -° (x °- π) T V ° 1 (x °- π)
2
Now notice that this holds for any vector m and symmetric
(positive-definite) matrix V (pd -> aTVa = Var(aTx) > 0 for all nonzero a)
Properties of the MVN - I
1) If x is MVN, any subset of the variables in x is also MVN
2) If x is MVN, any linear combination of the
elements of x is also MVN. If x is MVN(m,V)
for y = x + a;
y is MVNn (π + a; V )
n
for y = aT x =
X
y is N(a T π; aT V a)
ai x i ;
k= 1
≥
for
y = A x;
T
y is MVNm A π; A V A
¥
Properties of the MVN - II
3) Conditional distributions are also MVN. Partition x
into two components, x1 (m dimensional column vector)
and x2 ( n-m dimensional column vector)
µ
x=
∂
( )
x1
x2
µ
π=
( )
π1
π2
0
∂
V x1 x1
V = @
T
V x1 x2
and
V x1 x2
1
A
V x2 x2
x1 | x2 is MVN with m-dimensional mean vector
°- 1
π x 1 j x 2 = π1 + V x 1 x 2 V x 2 x 2 (x 2 -° π2 )
and m x m covariance matrix
°- 1
V x 1 j x 2 = V x 1 x 1 ° V x 1 x 2 V x 2 x 2 V xT1 x 2
Properties of the MVN - III
4) If x is MVN, the regression of any subset of
X on another subset is linear
linear and homoscedastic
x 1 = πx 1 j x 2 + e
°- 1
= π 1 + V x 1 x 2 V x 2 x 2 (x 2 °- π2 ) + e
Where e is MVN with mean vector 0 and
Variance-covariance matrix V
x1 jx2
The regression is homoscedastic
because
the function
variance
linear because it
is a linear
is
that does not change with the particular
ofaxconstant
2
value of the random vector x2
Example: Regression of Offspring value on Parental values
0
1
20
1
0
1 3
2 and the
2
Assume
the
vector
of
offspring
value
zo
πo
1
h =2 h =2 values of
A is
4@
zs parents
ª MVNis
πs A ; Then
æ2z @h2from
=2
1the correlations
0 A5
both @
its
MVN.
2 =2
z
π
h
0
1
d
d
among (outbred) relatives,
µ
Let x 1 = ( zo ) ;
V x1 ;x1 =
Hence,
æ2z ;
x2 =
zs
zd
∂
( )
h2æz2
V x 1 ;x 2 =
( 1 1) ;
2
V x 2 ; x 2 = æz2
( )
1 0
0 1
µ
∂µ
∂
h2 æz2
-° 2 1 0
zs -° πs
zo = πo +
( 1 1) æz
+e
0
1
z
°
π
d
d
2
2
h2
h
= πo +
(zs °- πs ) +
(zd °- πd ) + e
2
2
(
(
Where e is normal with mean zero and variance
°- 1
- hπ
2
2
2 æ2) + e
1 x 2 V - 2 21(x0
=
π
+
V
°
h
æ
x
2
2z 1
z
°x2 x
2
2 - 1
æe = æz °
( 1 1 ) æz
µ
= æz2
2
h4
1°
2
∂
(0
1
2
(1
Hence, the regression of offspring trait value given
the trait values of its parents is
zo = mo + h2/2(zs- ms) + h2/2(zd- md) + e
where the residual e is normal with mean zero and
Var(e) = sz2(1-h4/2)
Similar logic gives the regression of offspring breeding
value on parental breeding value as
Ao = mo + (As- ms)/2 + (Ad- md)/2 + e
= As/2 + Ad/2 + e
where the residual e is normal with mean zero and
Var(e) = sA2/2
Linear Models
Quick Review of the Major Points
The general linear model can be written as
y = Xb + e
Vector
Vector
of
of
observed
residuals.
dependent
Usually
assumed
values
MVN with
mean
Design
Vector
matrix
of
--unknown
observations
parents
of
to
the
estimate
variables
in
the
If the
residuals
are
homoscedastic
and
uncorrelated,
zero linear model
assumed
so that we can write the cov matrix of e as Cov(e) = s2I
then use the OLS estimate, OLS(b) = (XTX)-1 XTy
If the residuals are heteroscedastic and/or dependent,
Cov(e) = V and we must use GLS solution for b,
GLS(b) = (XT V-1 X)-1 V-1 XTy
Linear Models
One tries to explain a dependent variable y as a linear
function of a number of independent (or predictor)
variables.
A multiple regression is a typical linear model,
…
y = π + Ø1x 1 + Ø2 x2 + ¢¢¢+
Øn x x + e
Here
e is the regression
residual, orcoefficients
deviation between
the true
The (partial)
are interpreted
value
observed
andchange
the value
by the linear
as follows.
A unit
in xpredicted
i while holding all
model.
other variables constant results in a change of bi in y
As with a univariate regression (y =a + bx + e), the model
parameters are typically chosen by least squares,
wherein they are chosen to minimize the sum of
squared residuals, S ei2
Predictor and Indicator Variables
In a regression, the predictor variables are
typically continuous, although they need not be.
Suppose we measuring the offspring of p sires. One
linear model would be
yij = m + si + eij
Trait
valuefor
of offspring
jmean
fromof
sire
Overall
The
effect
mean.
sire
term
(the
invariables
included
to
itsigive
offspring).
the
sisindicator
terms
Recall
Note
We
can
The
that
write
deviation
theThis
this
predictor
in
ofilinear
the
jth
model
offspring
form
here
by
from
areusing
the
the
, something
ifamily
a that
mean
thatmean
variance
value
zero,
in the
i.e.,
sestimate
they
are of
expressed
Cov(half
as deviations
= VA/2the
i estimates
variables,
we
are
ofof
trying
sire
i. to
The
variance
the e’ssibs)
estimates
Ω
from the
mean
within-family
variance.
1 if sire k = i
xi k =
0 ot herwise
Giving the above model as
yi j = š +
Xp
s k x i k + ei j
k= 1
Models consisting entirely of indicator variables
are typically called ANOVA, or analysis of variance
models
Models that contain no indicator variables (other than
for the mean), but rather consist of observed value of
continuous or discrete values are typically called
regression models
Both are special cases of the General Linear Model
(or GLM)
yijk = m + si + dij + bxijk + eijk
Regression
model
ANOVA
model
Example: Nested
half sib/full
sib design with an
age correction b on the trait
Linear Models in Matrix Form
Suppose we have 3 variables in a multiple regression,
with four (y,x) vectors of observations.
y i = π + Ø1 x i 1 + Ø2x i 2 + Ø2 x i 2 + ei
In matrix form,
y = X Ø+ e
0 y 1
0 1 x
0 π1
1
11
y2 C
B Ø1 C X = B 1 x 21
Ø
=
y= B
@Ø A
@y A
@1 x 31
2
3
Ø3
y4
1 x 41
x 12
x 22
x 32
x 42
0 e 1
1
x 13 1
x 23 C
B e2 C
e
=
@e3 A
x 33 A
x 43
e4
Vector of
parameters
The
design matrix X. Details of both
to estimate
the experimental design and the
observed values of the predictor
variables all reside solely in X
The Sire Model in Matrix Form
The model here is yij = m + si + eij
Consider three sires. Sire 1 has 2 offspring, sire 2 has
one and sire 3 has three. The GLM form is
0 y 1
11
B y 12 C
By C
B
C
y = B 21 C ;
B y 31 C
@y A
32
y 33
0 1
B1
B1
B
X = B
B1
@1
1
1
1
0
0
0
0
0 01
0 0C
1 0C
C
C
0 1C
0 1A
0 1
0 π1
Bs C
Ø= @ 1 A ;
s2
s3
0 e 1
11
B e12 C
B
C
e
B
C
and e = B 21 C
B e31 C
@e A
32
e33
Ordinary Least Squares (OLS)
When the covariance structure of the residuals has a
certain form, we solve for the vector b using OLS
IfIfresiduals
followare
a MVN
distribution,
OLS
= ML solution
the residuals
homoscedastic
and
uncorrelated,
s2(ei) = se2, s(ei,ej) = 0. Hence, each residual is equally
weighted,
Sum of squared
residuals can
be written as
n
X
- X b) T (y °- X b)
2
e
=
e
e
=
(y
°
bi
b b
T
i= 1
Predicted value of the y’s
T
° 1 T
Taking (matrix) derivatives shows this is minimized by
Thisvariance-covariance
is the OLS estimateestimate
of the vector
b sample
Xn
1
The
for
the
estimates
2
c2 =
æ
e
b
e
i
is
n °- ran k(X ) i
T
°- 1 2
=1
V b = (X X ) æe
b = (X X )
X y
Example: Regression Through the Origin
yi = bxi + ei
0 y 1
1
B y2 C
y= B . C
@. A
.
yn
0 x 1
1
B x2 C
C
X = B
@ ... A
Here
xn
T
X X =
Xn
2
xi
T
X y=
b=
¥-° 1
(X X )
T
1
æ2 (b) = n° 1
xi yi
i= 1
i= 1
≥
Xn
Ø = ( Ø)
P
x i yi
X y= P 2
xi
T
P
(yi °- bx i ) 2
P 2
xi
2
æ (b) =
æe2
=
≥
T
¥-° 1
1
X
(X X )
n- 1
æ2e
æe2
= P 2
xi
(yi - bx i ) 2
Polynomial Regressions and
Interaction Effects
GLM
can easily
handle any function of the observed
Quadratic
regression:
predictor variables, provided the parameters to estimate
2
are still linear,
= 1a x+b
f(x)
+
b
yi =e.g.
Æ+Y Ø
+
Ø
x
+ e+i … + e
i 1
2 i1g(x)
0
1
Æ
Ø = @Ø1 A
Ø2
0 1 x
1
B 1 x2
X = B
..
@ ...
.
1 xn
x 21 1
x 22 C
.. C
. A
x 2n
Interaction terms (e.g. sex x age) are handled similarly
yi = Æ+ Ø1 x i 1 + Ø2 x i 2 + Ø3 x i 1 xi 2 + ei
0 1 x
x12 x11 x 12 1
0
1
11
Æ
With
x held
constant,
change in
x bchanges
Likewise,
a unit
change ainunit
x changes
y by
+b x y
1
by b2 + b3x1
1
Ø1 C
Ø= B
@Ø A
2
Ø3
B 1 x 212 x122 x321 x222 C
X = B
..
..
.. C
@ ...
.
.
. A
1 xn 1 x n 2 x n 1 x n 2
Fixed vs. Random Effects
In linear models are are trying to accomplish two goals:
estimation the values of model parameters and estimate
any appropriate variances.
For example in the simplest regression model,
y = a + bx + e, we estimate the values for a and b and
also the variance of e. We, of course, can also
estimate the ei = yi - (a + bxi )
Note that a/b are fixed constants are we trying to
estimate (fixed factors or fixed effects), while the
ei values are drawn from some probability distribution
(typically Normal with mean 0, variance s2e). The
ei are random effects.
This distinction between fixed and random effects is
extremely important in terms of how we analyzed a model.
If a parameter is a fixed constant we wish to estimate,
it is a fixed effect. If a parameter are drawn from
some probability distribution and we are trying to make
inference on either the distribution and/or specific
realizations from this distribution, it is a random effect.
We generally speak of estimating fixed factors and
predicting random effects.
“Mixed” models contain both fixed and random factors
Example: Sire model
yij = m + si + eij
It
Conversely,
depends.
Is the sire
if
Ifwe
effect
we
Fixed
are
have
not
fixed
effect
(say)
Random
only
or10
interested
random?
sires,
effect
if we
in these
are ONLY
interested
10 particular
in the
sires
values
but also
of these
wish to
particular
make some
10 sires and
don’t
inference
care to
about
make
the
any
population
other inferences
from which
about
they
the
population
were drawn
from
(such
which
as the
theadditive
sires are
variance,
drawn, then
sincewe can
treat
fixedthe
effects.
In the effects.
case, the In
model
s2A =them
4s2s, as
), then
si are random
thisis
fully
covariance
for the residuals.
case specified
we wish tothe
estimate
m andstructure
the variances
2 . Since
Thus,
s2s and
wesneed
to estimate
2si also m,
estimates
s1 to s10 (or
andpredicts)
s2e, and we
w
2(e
2 .
the breeding
value
wish
to
write
the model
as for
yij =sire
m + si,i we
+ ealso
,
s
)
=
s
ij
ij
e
estimate (predict) these as well. Under a
Random-effects interpretation, we write the model as
yij = m + si + eij, s2(eij) = s2e, s2(si) = s2s
Generalized Least Squares (GLS)
Suppose the residuals no longer have the same
variance (i.e., display heteroscedasticity). Clearly
we do not wish to minimize the unweighted sum
of squared residuals, because those residuals with
smaller variance should receive more weight.
Likewise in the event the residuals are correlated,
we also wish to take this into account (i.e., perform
a suitable transformation to remove the correlations)
before minimizing the sum of squares.
Either of the above settings leads to a GLS solution
in place of an OLS solution.
In the GLS setting, the covariance matrix for the
vector e of residuals is written as R where
Rij = s(ei,ej)
The linear model becomes y = Xb + e, cov(e) = R
The GLS solution for b is
b=
(X
T
R
°- 1
X)
°- 1
T
X R
°-1
y
The variance-covariance of the estimated model
parameters is given by
≥
¥°-1
T -° 1
Vb = X R X
æe2
(
)
Example
One common setting is where the residuals are uncorrelated
but heteroscedastic, with s2(ei) = s2e/wi.
For example, sample i is the mean value of ni individuals,
with s2(ei) = se2/ni. Here wi = ni
R=
° -1
°-1
°- 1
Diag(w1 ; w2 ; : : : ; wn )
Consider the model yi = a + bxi
µ
Ø=
∂
Æ
Ø
0 1 x 1
1
.. A
X = @ ...
.
1 xn
Here
R
°- 1
= Diag(w1; w2; : : : ; wn )
giving
T
X R
-
° 1
0
-
yw
1
XTR
A
y = w@
xy w
0
°- 1
1
X = w@
A
xw
where
n
w=
X
i=1
n
wi ;
xw =
X
i=1
Xn wi yi
yw =
;
w
i= 1
wi xi
;
w
xw
n
x2 w
=
X
i= 1
Xn wi xi yi
xyw =
w
i= 1
1
wi x2i
w
x2 w
This gives the GLS estimators of a and b as
a = yw °- bx w
b=
xyw -° x w y w
x 2w °- x w2
Likewise, the resulting variances and covariance for
these estimators is
2
æ (a) =
æ2e ¢x 2 w
w ( x 2 w °- x w2 )
æe2
æ (b) =
w ( x 2 w °- x w2 )
2
æ(a; b) =
-° æ2e x w
w ( x 2 w -° x w2 )
Chi-square and F distributions
Let Ui be a unit normal
The sum U12 + U22 + … + Un2 is a chi-square random
variable with k degrees of freedom
Under appropriate normality assumptions, the
sums of squares that appear in linear models
are also chi-square distributed. In particular,
Xn
2
(x i °- x ) ~
ª ¬ 2n °- 1
i =1
The ratio of two chi-squares is an F distribution
In particular, an F distribution with k numerator
degrees of freedom, and n denominator degrees
of freedom is given by
¬ 2k =k ~
ª Fk ;n
2
¬ n =n
Thus, F distributions frequently arrive in tests
of linear models, as these usually involve ratios
of sums of squares.
Sums of Squares in linear models
The total sums of squares SST of a linear model
can be written as the sum of the error (or residual)
sum of squares and the model (or regression) sum
of squares
SST = SSM + SSE
X
(yi °- y)
2
X
( ybi °°- y)
2
X
(yi °- ybi ) 2
r2, the coefficient
SSof
is the
SSE
M determination,
2
=
=1fraction of rvariation
accounted
SS
SSfor by the model
T
T
Sums of Squares are quadratic products
n
SST =
X
( yi °- y ) =
2
i =1
n
X
n
X
yi -° y =
2
2
i =1
√
yi °2
i =1
1
n2
(
!
n
X
2
yi
i =1
We can write this as a quadratic product as
1
SST = y T y °- y T Jy = y T
n
µ
∂
1
I °- J y
n
(
)
Where J is a matrix all of whose elements are 1’s
SSE =
Xn
2
( yi °- ybi ) =
i =1
T
≥
¥ °-1
( I °- X (X X )
T
SSM = SST °- SSE = y
2
e
bi
i =1
µ
SSE = y
Xn
µ
T
≥
X
T
)y
¥ °- 1
( X (X X )
T
∂
X
T
∂
1
°- J y
n
)
Hypothesis testing
Provided the residual errors in the model are
MVN, then for a model with n observations
and p estimated parameters,
SSE
2
ª
¬
~
n °- p
æe2
Consider the comparison of a full (p parameters)
and reduced (q < p) models
µ difference
∂ ¡ µ
∂
µsum
∂ µsquares ∂
The
of
SSE r -° SSE f in the
SSEerror
n
°
p
SSE r f
=
° 1
for thepfull
°- q and reduced
n °- p model
p °-provided
q
SSE f a test
for whether the model fit is the same
(
This ratio follows an Fp-q,n-p distribution
)
Does our model account for a significant fraction of
the variation?
Here the reduced model is just yi = u + ei
In this case, the error sum of squares for the
reduced model is just the total sum of squares
and the F test ratio becomes
µ
n °- p
p °- 1
∂µ
(
∂
µ
∂ µ
∂
SST
n °- p
r2
°- 1 =
* 1 °- r 2
SSE f
p °- 1
)
This ratio follows an Fp-1,n-p distribution