Summary and Analysis of Extension Program Evaluation in R

Salvatore S. Mangiafico

Least Square Means for Multiple Comparisons

Comparisons of values across groups in cumulative link models can be conducted easily with the lsmeans package.  The package can be used with a variety of other model types as well, and with fairly complex models.  Importantly, it can make comparisons among interactions of factors.


This will often be the easiest way to compare among treatment means for CLM.


L.S. means stands for least square means.  Least square means are means for treatment levels that are adjusted for means of other factors in the model.  For a more complete explanation, see the What are least square means? chapter.


p-value adjustments for multiple comparisons

Adjustment of p-values for multiple comparisons is indicated with the adjust= option in the lsmeans function.  Options are "tukey", "scheffe", "sidak", "bonferroni", "dunnettx", "mvt", and "none".  For further information on these options, see the summary function section in:


Lenth, R.V. and M. Hervé. 2015. Package ‘lsmeans’. cran.r-project.org/web/packages/lsmeans/lsmeans.pdf.


Note that the adjust= option should also be applied to the cld function if a compact letter display is requested.


Ignore values of lsmeans for clm and clmm models

Typically you should ignore the values of the LS means themselves (lsmeans) when using them with clm and clmm models.  With default settings, the values of the LS means and the values of differences among the LS means are not easy to interpret.  See “clm” in help(models, package="lsmeans") for details and for other options.


Packages used in this chapter


The packages used in this chapter include:

•  psych

•  ordinal

•  car

•  RVAideMemoire

•  lsmeans

•  multcompView


The following commands will install these packages if they are not already installed:


Multiple comparisons with lsmeans


This example uses data set and model from the One-way Ordinal ANOVA with CLM chapter.

Input =("
 Speaker  Likert
 Pooh      3
 Pooh      5
 Pooh      4
 Pooh      4
 Pooh      4
 Pooh      4
 Pooh      4
 Pooh      4
 Pooh      5
 Pooh      5
 Piglet    2
 Piglet    4
 Piglet    2
 Piglet    2
 Piglet    1
 Piglet    2
 Piglet    3
 Piglet    2
 Piglet    2
 Piglet    3
 Tigger    4
 Tigger    4
 Tigger    4
 Tigger    4
 Tigger    5
 Tigger    3
 Tigger    5
 Tigger    4
 Tigger    4
 Tigger    3

Data = read.table(textConnection(Input),header=TRUE)

### Order levels of the factor; otherwise R will alphabetize them

Data$Speaker = factor(Data$Speaker,

### Create a new variable which is the likert scores as an ordered factor

Data$Likert.f = factor(Data$Likert,
                       ordered = TRUE)

###  Check the data frame





### Remove unnecessary objects


Define model and conduct analysis of deviance


model = clm(Likert.f ~ Speaker,
            data = Data)



      type = "II")

Analysis of Deviance Table (Type II tests)

        LR Chisq Df Pr(>Chisq)   
Speaker   23.395  2  8.315e-06 ***

Group separations by lsmeans in table format

In the lsmeans function, model specifies the model object that was previously fitted.  Note the specialized formula where pairwise indicates that all pairwise comparisons should be conducted, and Speaker indicates the variable whose levels will be compared.


Remember to ignore the estimate, SE, and lsmean values when using lsmeans with a clm or clmm model object, unless specific options in lsmeans are selected.


        pairwise ~ Speaker,
        adjust="tukey")       ### Tukey-adjusted comparisons

 contrast         estimate        SE df    z.ratio p.value
 Pooh - Piglet    4.943822 1.3764706 NA  3.5916658  0.0010
 Pooh - Tigger    0.633731 0.9055691 NA  0.6998152  0.7636
 Piglet - Tigger -4.310091 1.3173294 NA -3.2718403  0.0031

P value adjustment: tukey method for comparing a family of 3 estimates

   ### Remember to ignore “estimate” and “SE” of differences with CLM,
   ###   as well as “lsmeans” in the output not shown here

Compact letter display

Here we’ll create an object of the lsmeans output called leastsquare.  Then we’ll pass this object to the cld function to create a compact letter display.



leastsquare = lsmeans(model,
                      pairwise ~ Speaker,
                      adjust="tukey")      ###  Tukey-adjusted comparisons

    Letters=letters,      ### Use lower-case letters for .group
    adjust="tukey")       ### Tukey-adjusted comparisons 

Speaker    lsmean        SE df  asymp.LCL  asymp.UCL .group
 Piglet  -1.783599 0.7755958 NA -3.6355182 0.06832068  a   
 Tigger   2.526492 0.8413200 NA  0.5176407 4.53534421   b  
 Pooh     3.160223 0.8912968 NA  1.0320404 5.28840661   b  

Confidence level used: 0.95
Conf-level adjustment: sidak method for 3 estimates
P value adjustment: tukey method for comparing a family of 3 estimates
significance level used: alpha = 0.05

   ### Groups sharing a letter in .group are not significantly different

   ### Remember to ignore “lsmean”, “SE”, “LCL”, and “UCL” for CLM
   ###   unless certain options are used.