[R] a mixed effects model with collinear parameters.

[ivan.borozan at utoronto.ca]

Hi all,

I have the following unbalanced data set :

F1<-rep("F1",29)
F2<-rep("F2",83)
F3<-rep("F3",18)
F4<-rep("F4",11)
F5<-rep("F5",25)
F6<-rep("F6",10)

all.fac<-c(F1,F2,F3,F4,F5,F6)

batch1<-rep("B1",29)
batch2<-rep("B2",83)
batch3<-rep("B2",18)
batch4<-rep("B1",11)
batch5<-rep("B2",25)
batch6<-rep("B2",10)


batch<-c(batch1,batch2,batch3,batch4,batch5,batch6)


sample0[batch == "B1"]<-rnorm(length(sample0[batch == "B1"]),5,1)
sample0[batch == "B2"]<-rnorm(length(sample0[batch == "B2"]),2,1)


dataA <- data.frame(samples = sample0, batches=factor(batch),  
soil=factor(all.fac))
dataA$soil <- relevel(dataA$soil,ref="F2")

X<-matrix
XX<-t(X)%*%X
solve(XX)
Error in solve.default(XX) : system is computationally singular:  
reciprocal condition number = 3.74708e-18


Meaning there is collinearity between predictors. (F1 and B1) and (F4 and B1)


I have tried a mixed effects model.


> amod <- lme(fixed = samples ~ soil, random = ~ 1|batches, dataA, method="ML")
> summary(amod)
Linear mixed-effects model fit by maximum likelihood
  Data: dataA
        AIC      BIC    logLik
   507.0955 532.4594 -245.5478

Random effects:
  Formula: ~1 | batches
         (Intercept)  Residual
StdDev: 2.02532e-05 0.9764997

Fixed effects: samples ~ soil
                  Value Std.Error  DF   t-value p-value
(Intercept)  2.1540927 0.1090599 169 19.751470  0.0000
soilF1       3.1223796 0.2143260 169 14.568363  0.0000
soilF3      -0.2070161 0.2583386 169 -0.801336  0.4241
soilF4       2.7412508 0.3188104 169  8.598372  0.0000
soilF5      -0.1668292 0.2266767 169 -0.735979  0.4628
soilF6      -0.2522510 0.3325879 169 -0.758449  0.4492
  Correlation:
        (Intr) soilF1 soilF3 soilF4 soilF5
soilF1 -0.509
soilF3 -0.422  0.215
soilF4 -0.342  0.174  0.144
soilF5 -0.481  0.245  0.203  0.165
soilF6 -0.328  0.167  0.138  0.112  0.158

Standardized Within-Group Residuals:
         Min          Q1         Med          Q3         Max
-3.02335670 -0.59885469  0.03137624  0.65601671  3.71106690

Number of Observations: 176
Number of Groups: 2

> ranef(amod)
     (Intercept)
B1 5.944455e-23
B2 1.795242e-23


The fixed effects part seems to be right, should I just ignore the  
random effects part ?