1
$\begingroup$

I performed mixed model ANOVAs using the lmer function from lme4.

My model is as follows.

mixedhrt<-lmer(HR_COR~Optimel+Order+Gender+SipSize +Hunger + GHI+(1|Participant),data=Taste2)

I want to see the effect of package condition, the order of my blocks, gender, the sip size, hunger ratings and general health interest on heart rate.

The structure of the data:

'data.frame':   498 obs. of  8 variables:
 $ Participant     : Factor w/ 98 levels "2","3","4","5",..: 1 1 1 1 1 2 2 2 2 2 ...
 $ Optimel         : Factor w/ 5 levels "BHL","BLH","Practice",..: 1 2 3 4 5 1 2 3 4 5 ...
 $ HR_COR          : num  6.091 5.855 -0.773 NA 2.676 ...
 $ Hunger          : num  1.5 1.5 1.5 1.5 1.5 13.7 13.7 13.7 13.7 13.7 ...
 $ Order           : Factor w/ 5 levels "0","1","2","4",..: 2 3 1 4 5 3 5 1 2 4 ...
 $ SipSize         : num  12.6 13.8 11.3 13.2 13.0 11.6  9.9 10.7 12.3 13.1 ...
 $ Gender          : Factor w/ 2 levels "male","female": 2 2 2 2 2 2 2 2 2 2 ...
 $ GHI             : num  6 6 6 6 6 5 5 5 5 5 ...

When running anova(mixedhrt) I get the following output:

fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
Analysis of Variance Table of type II  with  Satterthwaite approximation for degrees of freedom
     Sum Sq Mean Sq NumDF  DenDF F.value    Pr(>F)    
Optimel  35.840   8.960     4 328.90  1.1218   0.34608    
Order   183.865  61.288     3 328.22  7.6733 5.698e-05 ***
Gender   10.714  10.714     1  85.53  1.3414   0.25001    
Hunger    2.288   2.288     1  85.72  0.2865   0.59388    
GHI       6.168   6.168     1  84.32  0.7723   0.38201    
SipSize  33.142  33.142     1 392.39  4.1494   0.04232 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

I am interested in which order factors are significantly different from each other (pairwise comparison). When I perform the lsmeans function lsmeans(mixedhrt, pairwise ~ Order) on the lmer model, I only receive the contrasts for order factors 1,2,4,5 and no contrasts with factor 0.

$lsmeans
Order lsmean SE df asymp.LCL asymp.UCL
 0         NA NA NA        NA        NA
 1         NA NA NA        NA        NA
 2         NA NA NA        NA        NA
 4         NA NA NA        NA        NA
 5         NA NA NA        NA        NA

Results are averaged over the levels of: Optimel, Gender 
Degrees-of-freedom method: satterthwaite 
Confidence level used: 0.95 

$contrasts
 contrast    estimate        SE     df t.ratio p.value
 0 - 1             NA        NA     NA      NA      NA
 0 - 2             NA        NA     NA      NA      NA
 0 - 4             NA        NA     NA      NA      NA
 0 - 5             NA        NA     NA      NA      NA
 1 - 2     1.32707259 0.4381103 327.41   3.029  0.0221
 1 - 4     1.24663731 0.4434918 328.08   2.811  0.0415
 1 - 5     2.06891027 0.4379091 328.43   4.725  <.0001
 2 - 4    -0.08043529 0.4380253 328.92  -0.184  0.9997
 2 - 5     0.74183767 0.4390263 328.85   1.690  0.4417
 4 - 5     0.82227296 0.4433427 327.58   1.855  0.3441

Results are averaged over the levels of: Optimel, Gender 
P value adjustment: tukey method for comparing a family of 5 estimates 

My guess is, that it has to do with the rank deficiency, which in turn has to do with the nature of my data. I have 5 order time points (0,1,2,4,5). During time points 1,2,4,5 I tested one of 4 different conditions, whereas time point 0 was a practice run and was therefore always the same condition. The conditions are classified through the variable Optimel. Therefore, there seem to be issues with the calculations, as the anova also only returns 3 NumDF, where there should be 4. I'm afraid I'm lacking necessary statistical and R knowledge to understand how to deal with this. Therefore, I would be very grateful for any help on possibilities to calculate all contrasts with my current data structure.

$\endgroup$
  • $\begingroup$ I think your diagnosis is correct. Why not just ignore all the comparisons with 0 if they make no sense? $\endgroup$ – mdewey Jan 4 '17 at 13:20
  • $\begingroup$ Because I am also interested in how the averages of the conditions at the different time points compare to each other, including time point o. $\endgroup$ – Aurelia Jan 5 '17 at 8:32
2
$\begingroup$

The issue is definitely related to the rank deficiency.

The OP seems most concerned about the comparisons, but first and foremost, one should be concerned about the fact that none of the lsmeans for Order are estimable. That means that we really don't know what we are discussing when we compare them.

Often, a rank deficiency comes about when a model includes too many interactions and there are empty cells in the design. The fact that the fixed part of the model is additive, and you still have rank deficiencies, means that the deficiency is particularly egregious.

To work around this, either a simpler model (omitting one or more factors) or subsetting the dataset (e.g., omitting time point 0) will be necessary. Another possibility may be to analyze gain scores (difference from time point 0), which of course also excludes that time point as a factor. If time 0 is really a baseline measurement, that could make a lot of sense, and also incorporates a certain kind of comparison with time 0, which the OP states is of interest. If you continue to get non-estimable lsmeans, you're not done.

After fitting a more suitable model, use appropriate diagnostics to check that the model fits the data reasonably. This is an important step and it may be wise to seek consulting help to make sure this is done correctly. Statistics is more than just calculations; it is fitting models that actually do a reasonable job of explaining things and follow the distributional assumptions that underlie it.

Remember that the GIGO (garbage in, garbage out) principle always applies. If you can't fit a sensible model and do the comparisons you need to do, then there is no getting around the fact that the data are unsuitable for the inferences you need -- as much as one may wish otherwise.

$\endgroup$
  • $\begingroup$ Thank you for your answer. It was starting to dawn on me that the comparisons were just a symptom and not the real problem. I will talk to my supervisor and see which one of your suggested options works best in this case. I thought I had checked the diagnostics correctly, but I will see that I find some help on that as well. $\endgroup$ – Aurelia Jan 6 '17 at 9:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.