i am trying to use matinv from the Design package
to compute the generalized inverse of the normal equations
of a 3x3 design via the sweep operator.
That is, for the linear model
y = ยต + x1 + x2 + x1*x2
where x1, x2 are 3-level factors and dummy coding is being used
the matrix to be inverted is
X'X =
9 3 3 3 3 3 3 1 1 1 1 1 1 1 1 1
3 3 0 0 1 1 1 1 0 0 1 0 0 1 0 0
3 0 3 0 1 1 1 0 1 0 0 1 0 0 1 0
3 0 0 3 1 1 1 0 0 1 0 0 1 0 0 1
3 1 1 1 3 0 0 1 1 1 0 0 0 0 0 0
3 1 1 1 0 3 0 0 0 0 1 1 1 0 0 0
3 1 1 1 0 0 3 0 0 0 0 0 0 1 1 1
1 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0
1 0 1 0 1 0 0 0 1 0 0 0 0 0 0 0
1 0 0 1 1 0 0 0 0 1 0 0 0 0 0 0
1 1 0 0 0 1 0 0 0 0 1 0 0 0 0 0
1 0 1 0 0 1 0 0 0 0 0 1 0 0 0 0
1 0 0 1 0 1 0 0 0 0 0 0 1 0 0 0
1 1 0 0 0 0 1 0 0 0 0 0 0 1 0 0
1 0 1 0 0 0 1 0 0 0 0 0 0 0 1 0
1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 1
this matrix has rank=9
however, matinv(X'X) falsely returns rank=4,
no matter what the tolerance threshold eps is set to.
also the defining property of the
generalized inverse
_
X'X %*% (X'X) %*% X'X = X'X
is not satisfied.
if i use qr (from the base package) the rank is
correctly determined as 9.
any ideas?
thank you
--
View this message in context:
http://r.789695.n4.nabble.com/generalized-inverse-using-matinv-Design-tp3747337p3747337.html
Sent from the R help mailing list archive at Nabble.com.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
