Como se faz uma ANOVA SS tipo III em R com códigos de contraste?

26

Forneça um código R que permita realizar uma ANOVA entre sujeitos com -3, -1, 1, 3 contrastes. Entendo que há um debate sobre o tipo apropriado de soma de quadrados (SS) para essa análise. No entanto, como o tipo padrão de SS usado no SAS e SPSS (Tipo III) é considerado o padrão na minha área. Assim, eu gostaria que os resultados dessa análise correspondessem perfeitamente ao que é gerado por esses programas estatísticos. Para ser aceita, uma resposta deve chamar diretamente aov (), mas outras respostas podem ser votadas (especialmente se forem fáceis de entender / usar).

sample.data <- data.frame(IV=rep(1:4,each=20),DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))

Edit: Por favor, note que o contraste que estou solicitando não é um simples contraste linear ou polinomial, mas é um contraste derivado de uma previsão teórica, ou seja, o tipo de contraste discutido por Rosenthal e Rosnow.

russellpierce
fonte
5
Entendo que você precisa da soma do Tipo III, mas este artigo ( stats.ox.ac.uk/pub/MASS3/Exegeses.pdf ) é uma boa leitura. Ilustra alguns pontos interessantes.
suncoolsu
Em relação à sua pergunta, você pode estar interessado na seguinte discussão: stats.stackexchange.com/questions/60362/… A escolha entre ANOVA tipo I, II e III não é tão fácil quanto parece.
phx
Promoveu sua pergunta como útil, pois provocou várias respostas aprendidas, mas também notei que você concordou com o entrevistado, que basicamente disse que a premissa da pergunta estava incorreta. Espero estar resumindo a posição de StaGuy como dizendo que os contrastes definidos foram por definição "tipo I" e a discussão de outros tipos só se tornou relevante ao avaliar estatísticas de regressão parcial, presumivelmente mais importantes ao deixar "a máquina dirigir" usando métodos automatizados.
Dwin
@ DWin: Não tenho certeza de segui-lo inteiramente. Pode-se legitimamente usar outros tipos de SS sem deixar a "máquina dirigir" (pelo menos como eu entendo essa frase). Eu posso estar um pouco enferrujado aqui, mas se a memória servir, outros tipos podem ser relevantes quando não estiver usando regressão parcial. Por exemplo, o SS do Tipo III não parcial dos principais efeitos da interação. A distinção entre os tipos é importante justamente porque o Tipo III não é parcial, enquanto o Tipo I é. O problema, como afirmado, incluía apenas um único contraste e, portanto, a distinção entre os tipos de SS era / é discutível.
117915 russellpierce
Meu entendimento era de que a justificativa dada pelo SAS para a escolha do SSS do tipo III (e esse parece ser o motivo pelo qual as pessoas pensam que o tipo III é o preferido) é que ele suporta melhor o processo de seleção para trás e para frente.
Dwin

Respostas:

22

A soma dos quadrados do tipo III para ANOVA está prontamente disponível através da Anova()função do pacote veicular .

A codificação de contraste pode ser feita de várias maneiras, usando C()a contr.*família (como indicado por @nico) ou diretamente a contrasts()função / argumento. Isso está detalhado no §6.2 (pp. 144-151) da Estatística Moderna Aplicada com S (Springer, 2002, 4ª ed.). Observe que aov()é apenas uma função de invólucro para a lm()função. É interessante quando se deseja controlar o termo de erro do modelo (como em um projeto dentro do assunto), mas, caso contrário, ambos produzem os mesmos resultados (e seja qual for a forma como você se encaixa no seu modelo, você ainda pode gerar ANOVA ou LM- como resumos com summary.aovou summary.lm).

Não tenho o SPSS para comparar as duas saídas, mas algo como

> library(car)
> sample.data <- data.frame(IV=factor(rep(1:4,each=20)),
                            DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))
> Anova(lm1 <- lm(DV ~ IV, data=sample.data, 
                  contrasts=list(IV=contr.poly)), type="III")
Anova Table (Type III tests)

Response: DV
            Sum Sq Df F value    Pr(>F)    
(Intercept)  18.08  1  21.815  1.27e-05 ***
IV          567.05  3 228.046 < 2.2e-16 ***
Residuals    62.99 76                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1 

vale a pena tentar em primeira instância.

Sobre a codificação de fatores no R vs. SAS: R considera o nível de linha de base ou referência como o primeiro nível em ordem lexicográfica, enquanto o SAS considera o último. Portanto, para obter resultados comparáveis, você precisa usar contr.SAS()ou ao relevel()seu fator R.

chl
fonte
11
Eu não acho que esta resposta use o contraste -3, -1,1,3 que eu especifiquei, nem parece fornecer um teste de 1 df do contraste.
russellpierce
@drknexus Sim, você está certo. Escreveu muito rapidamente. Algo como Anova(lm(DV ~ C(IV, c(-3,-1,1,3),1), data=sample.data), type="III")deveria ser melhor. Por favor, deixe-me saber se isso está bem com você.
chl
Obrigado! Parece bom. Eu validarei o SPSS amanhã e volto para você.
russellpierce
11
BTW, dê uma olhada no pacote ez ( cran.r-project.org/web/packages/ez/index.html ) para agrupar o código Anova ...
Tal Galili
2
@drknexus: Se ao menos houvesse uma página de apresentação pedido de recurso e questões para ez ... github.com/mike-lawrence/ez/issues :)
Mike Lawrence
11

Isso pode parecer um pouco de autopromoção (e acho que sim). Mas desenvolvi um pacote lsmeans para R (disponível no CRAN), projetado para lidar exatamente com esse tipo de situação. Aqui está como isso funciona para o seu exemplo:

> sample.data <- data.frame(IV=rep(1:4,each=20),DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))
> sample.aov <- aov(DV ~ factor(IV), data = sample.data)

> library("lsmeans")
> (sample.lsm <- lsmeans(sample.aov, "IV"))
 IV    lsmean        SE df   lower.CL  upper.CL
  1 -3.009669 0.2237448 76 -3.4552957 -2.564043
  2 -3.046072 0.2237448 76 -3.4916980 -2.600445
  3  1.147080 0.2237448 76  0.7014539  1.592707
  4  3.049153 0.2237448 76  2.6035264  3.494779

> contrast(sample.lsm, list(mycon = c(-3,-1,1,3)))
 contrast estimate       SE df t.ratio p.value
 mycon    22.36962 1.000617 76  22.356  <.0001

Você pode especificar contrastes adicionais na lista, se quiser. Neste exemplo, você obterá os mesmos resultados com o contraste polinomial linear interno:

> con <- contrast(sample.lsm, "poly")
> con
 contrast   estimate        SE df t.ratio p.value
 linear    22.369618 1.0006172 76  22.356  <.0001
 quadratic  1.938475 0.4474896 76   4.332  <.0001
 cubic     -6.520633 1.0006172 76  -6.517  <.0001

Para confirmar isso, observe que a "poly"especificação direciona a chamada poly.lsmc, o que produz estes resultados:

> poly.lsmc(1:4)
  linear quadratic cubic
1     -3         1    -1
2     -1        -1     3
3      1        -1    -3
4      3         1     1

Se você deseja fazer um teste conjunto de vários contrastes, use a testfunção com joint = TRUE. Por exemplo,

> test(con, joint = TRUE)

Isso produzirá um teste "tipo III". Ao contrário car::Anova(), ele será feito corretamente, independentemente da codificação de contraste usada no estágio de ajuste do modelo. Isso ocorre porque as funções lineares que estão sendo testadas são especificadas diretamente, e não implicitamente, através da redução do modelo. Uma característica adicional é que um caso em que os contrastes que estão sendo testados são linearmente dependentes é detectado e a estatística de teste correta e os graus de liberdade são produzidos.

rvl
fonte
7

Quando você está fazendo contrastes, está fazendo uma combinação linear específica declarada de médias de células no contexto do termo de erro apropriado. Como tal, o conceito de "Tipo de SS" não tem significado para contrastes. Cada contraste é essencialmente o primeiro efeito usando um SS Tipo I. "Tipo de SS" tem a ver com o que é parcial ou explicado pelos outros termos. Para contrastes, nada é dividido ou contabilizado. O contraste permanece por si só.

StatGuy
fonte
Você está absolutamente certo.
russellpierce
3

O fato de os testes do tipo III serem usados ​​no seu local de trabalho é o mais fraco dos motivos para continuar usando-os. O SAS causou grandes danos às estatísticas nesse sentido. A exegese de Bill Venables, mencionada acima, é um ótimo recurso. Apenas diga não ao tipo III; baseia-se em uma noção defeituosa de equilíbrio e tem menor poder devido ao peso tolo das células no caso de desequilíbrio.

Uma maneira mais natural e menos propensa a erros de obter contrastes gerais e de poder descrever o que você fez é fornecida pela função do rmspacote R. contrast.rmsOs contrastes podem ser muito complexos, mas para o usuário são muito simples porque são declarados em termos de diferenças nos valores preditivos. Testes e contrastes simultâneos são suportados. Ele lida com efeitos de regressão não lineares, efeitos de interação não lineares, contrastes parciais, todo tipo de coisa.

Frank Harrell
fonte
Isso é bom e bom para você, como uma pessoa estabelecida de renome a dizer. Outros não têm influência para discordar dos revisores. Como as interpretações das estatísticas diferem, você solicita que as pessoas novas sigam os princípios e incorram em um custo irracional. Eu digo isso como alguém que morreu minha parte de vezes no topo desta (e similar) colinas. A mudança da OMI nesta frente é de responsabilidade dos porteiros, ou seja, editores e revisores.
russellpierce
As pessoas que são realmente boas com dados têm uma grande variedade de empregos e podem ter a opção de trabalhar em áreas onde suas habilidades e opiniões são respeitadas.
Frank # # # # Harrell Harrell
11
... e é isso que eu faço agora. Mas as pessoas que estão chegando a essa pergunta freqüentemente não pertencem a essa classe. Assim como eu não era há 7 anos. Só advogo um pouco de empatia pelo novato e suas circunstâncias.
Russellpierce 5/05
2

Experimente o comando Anova na biblioteca de carros. Use o argumento type = "III", pois o padrão é o tipo II. Por exemplo:

library(car)
mod <- lm(conformity ~ fcategory*partner.status, data=Moore, contrasts=list(fcategory=contr.sum, partner.status=contr.sum))
Anova(mod, type="III")
jebyrnes
fonte
3
Sei que Moore está na biblioteca de carros, mas Quando os dados de amostra são fornecidos, é mais fácil para o solicitante de entender sua resposta se você usar os dados de amostra.
russellpierce
0

Também autopromotora, escrevi uma função para exatamente isso: https://github.com/samuelfranssens/type3anova

Instale da seguinte maneira:

library(devtools)
install_github(samuelfranssens/type3anova)
library(type3anova)

sample.data <- data.frame(IV=rep(1:4,each=20),DV=rep(c(-3,-3,1,3),each=20)+rnorm(80))

type3anova(lm(DV ~ IV, data = sample.data))

Você também precisará ter o carpacote instalado.

sam_f
fonte
Como você aplicaria isso à parte contrastante da pergunta?
russellpierce
11
Desculpas, eu não li a pergunta corretamente. Minha função simplificará apenas a realização de Anova tipo III. Como o StatGuy acima, não vejo onde o SS entra em jogo ao testar contrastes específicos.
sam_f