I'm conducting a multivariate multiple regression. As response variable, I have a function [ y(t) ] that I have discretized on a grid of 27 points, and 3 scalar regressors (x1,x2,x3). I have replaced the response function with an nxq (q=27) matrix Y and what I need to solve is:
Y = XB + E
where X [nxp (p=3)] is x1,x2,x3 column-stacked, B is a p×q matrix of regression coefficients and E is an n×q matrix of errors.
What I have done up to now is calling lm: mylm<-lm(Y ~ X)
which regresses each dependent variable separately on the predictors.
Now I want to determine whether a predictor jointly contributes to all the 27 models I get, but I don't know how to overcome the errors I get.
When I call Anova
this is what I get
Anova(mylm)
Error in eigen(qr.coef(if (repeated) qr(x$SSPE[[term]]) else SSPE.qr, :
infinite or missing values in 'x'
and, if for example, I want to test if x2 is statistically different from 0, I get the following
hyp =c(0,0,1,0)
rhs =rep(0,27)
lh.out <- linearHypothesis(mylm, hyp,rhs)
Error in linearHypothesis.mlm(mylm, hyp, rhs) :
The error SSP matrix is apparently of deficient rank = 25 < 27
If the problem is related to the non-singularity of a matrix, how can I ask to use the pseudoinverse?
EDIT--------------
I put your data in this gist, so it can be easily sourced in. I actually get the Anova()
function from car
to work. I think the difference is putting all of the data in a data frame (or tibble in this case) and then using cbind()
to create the multivariate DV.
source("https://gist.githubusercontent.com/davidaarmstrong/469e2159d4802bae6fa09bad34527df0/raw/a25a02280841f4df3c6d8c10720f404586498b13/github_data1.r")
mydat <- cbind(Y, X)
mylm <- lm(cbind(V1, V2, V3, V4, V5, V6, V7,
V8, V9, V10, V11, V12, V13,
V14, V15, V16, V17, V18, V19,
V20, V21, V22, V23, V24, V25) ~
PC1D + PC2D + PC1H, data=mydat)
Anova(mylm)
# Type II MANOVA Tests: Pillai test statistic
# Df test stat approx F num Df den Df Pr(>F)
# PC1D 1 0.20085 4.424 25 440 3.981e-11 ***
# PC2D 1 0.18632 4.030 25 440 9.065e-10 ***
# PC1H 1 0.71205 43.522 25 440 < 2.2e-16 ***
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1