Multivariate linear mixed model using lmer

Question

I want to detect differences between 4 treatments of a food product based on the results from a sensory evaluation, where 9 panelists assessed the 4 products (i.e. treatments) following a number of sensory attributes. The output variables are scores (on a scale from 0 to 100) reflecting the intensity of each attribute. For each attribute I want to know how the 4 products differ (or not) from each other. To achieve this, I use a multivariate linear mixed model using the lmer() function. Here’s how the data looks like:

panelist treatment repetition Odor1 Odor2 Odor3 Flavor1 Flavor2 Flavor3 Flavor4 Texture1
    1     CTRL           1    35    25     0      64      83      21      20       84
    1     T1             1    29    19    48      29      38      53      53       58
    1     T2             1    28    22    38      36      44      51      51       39
    1     T3             1    28    28    67      27      66      46      45       72
    1     CTRL           2    19    32     0      42      43      22      21       42
    1     T1             2     6    32    11      65      57      79      79       71
    1     T2             2    46    32    78      38      63      63      63       63
    1     T3             2    31    31    78      85      31      61      60       48
    1     CTRL           3    64    51    50      50      67      67      67       76
    1     T1             3    42    37    59      82      58      35      34       67
    1     T2             3     0    26    58      57      70      21      21       54
    1     T3             3     0    41    69      49      70      24      24       57

In order to perform the multivariate analysis with lmer(), I reshaped the dataset to a long-format table as explained here (and there) so that the scores are grouped in one column and indexed to the "attribute" variable (listing all sensory attributes):

panelist treatment repetition attribute value
    1    CTRL          1        Odor1     35
    1    T1            1        Odor1     29
    1    T2            1        Odor1     28
    1    T3            1        Odor1     28
    1    CTRL          2        Odor1     19
    1    T1            2        Odor1     6
    1    T2            2        Odor1     46 
    1    T3            2        Odor1     31
    1    CTRL          3        Odor1     64
    ...

"Panelist" is treated as a random factor to adjust for different uses of the scale among the assessors. I then used the following model as a starting point:

mod <- lmer(score ~ treatment + (1|panelist) + (1|attribute), data)

I feel that treating “attribute” as a random effect is the only way to go and the results obtained make sense with my data as opposed to treating it as fixed effect. But it only returns the significance between the 4 groups based on the 12 attributes all together. So my questions are:

• How can I figure out the significance (or not) among treatment groups for each attribute?

• Is this model the right choice in this case or should I use a generalized linear mixed model instead?

Answer 1

With regard to your first question, you should include attribute also as a fixed effect, and potentially also consider its interaction with treatment.

The choice of linear mixed model vs a generalized linear mixed model depends on the outcome variable. Namely, linear mixed models assume that the distribution of your outcome conditional on the random effects is normal. A generalized mixed model would be used if you would have, for example, a binary outcome (Binomial distribution) or a count (Poisson distribution).

Answer 2

I actually found the solution using the SensMixed package (cran.r-project.org/web/packages/SensMixed/SensMixed.pdf) with no need to reshape the dataset.

Pretty good explainations are also available from the package author's PhD thesis