Why does Type II/III ANOVA on cumulative link model (R-package ordinal) give different results depending on factor order?

Question

From what I understand, type II and type III ANOVA should give the same result irrespective of the order of the factors in the formula because they calculate:

Type II: 
    SS(A | B) for factor A. 
    SS(B | A) for factor B. 
    SS(A:B|A,B) for the interaction. 
Type III:
    SS(A | B, A:B) for factor A. 
    SS(B | A, A:B) for factor B. 
    SS(A:B|A,B) for the interaction.

But when I run clm from the ordinal R-package on my data and then perform type II (or type III) ANOVA using the Anova function from the car R-package, I obtain completely different results depending on the order of my factors.

Edit 2: After a little bit of trial and error, I managed to make better example data for my issue:

dat<-data.frame(Treatment=c("Control","TreatmentA","TreatmentB","Control","TreatmentA","TreatmentB"),
Judge=c("C","C","C","E","E","E"),
VeryLow=c(2,3,2,2,3,2),
Low=c(1,3,6,1,3,6),
High=c(3,3,5,3,3,5),
VeryHigh=c(68,76,72,68,76,72))

# Some code to make the count data compatible with the clm package
orderedLevels<-c("VeryLow","Low","High","VeryHigh")
longDat<-reshape(dat, varying=orderedLevels, v.names="rating_count",timevar="rating",times=orderedLevels,direction="long")
indivDat<-longDat[rep(seq(1, nrow(longDat)), longDat$rating_count),]
indivDat$rating<-ordered(indivDat$rating, levels=orderedLevels)
summary(indivDat)

# The clm models without interactions
library(ordinal);library(car)
fm1 <- clm(rating ~ Treatment + Judge, data=indivDat)
fm2 <- clm(rating ~ Judge + Treatment, data=indivDat)
Anova(fm1, type="II"); Anova(fm2, type="II")

# The lm equivalent models without interactions
indivDat$rating_num<- as.numeric(indivDat$rating)
fm3 <- lm(rating_num ~ Treatment + Judge, data=indivDat)
fm4 <- lm(rating_num ~ Judge + Treatment, data=indivDat)
Anova(fm3, type="II"); Anova(fm4, type="II")

The results for this new example are:

> Anova(fm1, type="II")
Analysis of Deviance Table (Type II tests)

Response: rating
          Df  Chisq Pr(>Chisq)    
Treatment  2 92.725    < 2e-16 ***
Judge      1 53.275    2.9e-13 ***
---

> Anova(fm2, type="II")
Analysis of Deviance Table (Type II tests)

Response: rating
          Df  Chisq Pr(>Chisq)    
Judge      1 92.426  < 2.2e-16 ***
Treatment  2 76.882  < 2.2e-16 ***
---

> Anova(fm3, type="II")
Anova Table (Type II tests)

Response: rating_num
           Sum Sq  Df F value Pr(>F)
Treatment   1.177   2  1.3917 0.2497
Judge       0.000   1  0.0000 1.0000
Residuals 204.659 484   

> Anova(fm4, type="II")
Anova Table (Type II tests)

Response: rating_num
           Sum Sq  Df F value Pr(>F)
Judge       0.000   1  0.0000 1.0000
Treatment   1.177   2  1.3917 0.2497
Residuals 204.659 484  

How can the Judge effect be significant with the clm package when the data for each judge is just a copy of each other? And why does the order of the factor in the cum link model matters?

What is going on here? Did I misunderstood the type II/III Anova? And what is the best course of action when analyzing such a data set?

Original example (without results): Here is a reproducible example for Type II using the wine dataset:

library(ordinal)
fm1 <- clm(rating ~ contact * bottle, data=wine)
fm2 <- clm(rating ~ bottle * contact, data=wine)

library(car)
Anova(fm1, type="II"); Anova(fm2, type="II")

Edit 1: Same model without interactions.

> fm1 <- clm(rating ~ contact + bottle, data=wine)
> Anova(fm1, type="II")
Analysis of Deviance Table (Type II tests)

Response: rating
        Df   Chisq Pr(>Chisq)    
contact  1  1.8182     0.1775    
bottle   7 53.1100  3.526e-09 ***

> fm2 <- clm(rating ~ bottle + contact, data=wine)
> Anova(fm2, type="II")
Analysis of Deviance Table (Type II tests)

Response: rating
        Df  Chisq Pr(>Chisq)    
bottle   7 71.534  7.231e-13 ***
contact  1  6.847   0.008879 ** 

Answer 1

The documentation for car::Anova does not indicate that it works with clm or clmm objects. From your examples, it appears it can have some unusual behavior.

My solution is to use RVAideMemoire::Anova.clm.

So, for your EDIT 2 example:

if(!require(RVAideMemoire)){install.packages("RVAideMemoire")}
if(!require(car)){install.packages("car")}

library(car)

library(RVAideMemoire)

Anova(fm1, type="II")

### Analysis of Deviance Table (Type II tests)

### Response: rating
###           LR Chisq Df Pr(>Chisq)
### Treatment   4.0507  2     0.1319
### Judge       0.0000  1     1.0000

Anova(fm2, type="II")

### Analysis of Deviance Table (Type II tests)

### Response: rating
###           LR Chisq Df Pr(>Chisq)
### Judge       0.0000  1     1.0000
### Treatment   4.0507  2     0.1319

Note that this produces the same p-values as the anova function.

fm1 <- clm(rating ~ Treatment + Judge, data=indivDat)
fmx <- clm(rating ~ Treatment, data=indivDat)
fmy <- clm(rating ~ Judge, data=indivDat)

anova(fm1, fmx)

### Likelihood ratio tests of cumulative link models:

###     formula:                   link: threshold:
### fmx rating ~ Treatment         logit flexible  
### fm1 rating ~ Treatment + Judge logit flexible  

###     no.par    AIC  logLik LR.stat df Pr(>Chisq)
### fmx      5 474.85 -232.42                      
### fm1      6 476.85 -232.42       0  1  

anova(fm1, fmy)

### Likelihood ratio tests of cumulative link models:

###     formula:                   link: threshold:
### fmy rating ~ Judge             logit flexible  
### fm1 rating ~ Treatment + Judge logit flexible  

###     no.par    AIC  logLik LR.stat df Pr(>Chisq)
### fmy      4 476.90 -234.45                      
### fm1      6 476.85 -232.42  4.0507  2     0.1319