ANOVA Bifatorial em R
MSP4112: Estatística Aplicada à Medicina
28 outubro 2025 17:25h
Bastão de Asclépio & Distribuição Normal
suppressMessages(library(bruceR, warn.conflicts=FALSE))
suppressMessages(library(car, warn.conflicts=FALSE))
suppressMessages(library(effectsize, warn.conflicts=FALSE))
suppressMessages(library(emmeans, warn.conflicts=FALSE))
suppressMessages(library(estimatr, warn.conflicts=FALSE))
suppressMessages(library(FSA, warn.conflicts=FALSE))
suppressMessages(library(ggplot2, warn.conflicts=FALSE))
suppressMessages(library(gplots, warn.conflicts=FALSE))
suppressMessages(library(jmv, warn.conflicts=FALSE))
suppressMessages(library(lmboot, warn.conflicts=FALSE))
suppressMessages(library(lmerTest, warn.conflicts=FALSE))
suppressMessages(library(multcomp, warn.conflicts=FALSE))
suppressMessages(library(MuMIn, warn.conflicts=FALSE))
suppressMessages(library(patchwork, warn.conflicts=FALSE))
suppressMessages(library(phia, warn.conflicts=FALSE))
suppressMessages(library(psych, warn.conflicts=FALSE))
suppressMessages(library(pwr2, warn.conflicts=FALSE))
suppressMessages(library(RcmdrMisc, warn.conflicts=FALSE))
suppressMessages(library(readxl, warn.conflicts=FALSE))
suppressMessages(library(twowaytests, warn.conflicts=FALSE))
suppressMessages(library(WebPower, warn.conflicts=FALSE))Adicionar
bruceRbruceR::cor_diffbruceR::CorrbruceR::DescribebruceR::RECODEbruceR::model_summarybruceR::regressbruceR::EMMEANSbruceR::GLM_summarybruceR::HLM_summarybruceR::TTESTbruceR::MANOVA: https://psychbruce.github.io/bruceR/reference/MANOVA.html
apaTables::apa.aov.table
apaTables::apa.aov.table(fit,
conf.level=1-alfa,
type=3)
Material
- HTML de Rmarkdown em
RPubs
Objetivos
- Discorrer sobre a extensão da ANOVA unifatorial que inclui dois ou mais fatores.
- Descrever e indicar três delineamentos diferentes da ANOVA
bifatorial:
- independente: entre participantes,
- relacionada: intraparticipantes,
- mista: um fator entre e um fator intraparticipantes.
- Proceder com a análise descritiva numérica e gráfica dos dados.
- Formular e implementar os modelos de ANOVA bifatorial adequados a cada situação.
- Determinar as significâncias estatística e prática dos efeitos de interação e principais da ANOVA bifatorial.
- Executar e interpretar a análise de efeito principal simples.
- Executar e interpretar o teste post hoc, quando adequados.
Introdução
ANOVA bifatorial testa efeitos principais e de interação entre fatores numa única análise.
As variáveis envolvidas no teste estatístico são:
VD (variável dependente ou de desfecho) intervalar com distribuição normal por condição: univariada
duas ou mais VI (variável independente) com dois ou mais níveis (fatores nominais dicotômicos ou politômicos): multifatorial
ANOVA multifatorial pode ser:
- independente (fator entre participantes)
- relacionada (fator intraparticipantes)
- mista (fatores entre e intraparticipantes)
ANOVA multifatorial é uma extensão do ANOVA unifatorial.
As suposições de normalidade homocedasticidade são condições suficientes e a suposição de independência é necessária para ANOVA multifatorial.
ANOVA multifatorial em R
- Planejamento de ANOVA
- independente bifatorial de Fisher
pwr2::pwr.2way
- Multifatorial Independente de Fisher
WebPower::wp.kanova
- ANOVA multifatorial independente
- Fisher (homocedástico)
- teste omnibus:
lm,car::Anova,effectsize::eta_squared- teste de efeito principal simples:
phia::testInteractions- teste post hoc:
emmeans::emmeans,multcomp::cld,emmeans::emmip
- teste post hoc:
- teste de efeito principal simples:
- testes omnibus e post hoc:
jmv::ANOVA
- Fisher-White (heterocedástico)
- teste omnibus:
lm,car::Anova(white.adjust="hc2", ...)),estimatr::lm_robust(se_type="HC2", ...)- teste de efeito principal simples:
phia::testInteractions- teste post hoc:
sandwich::vcovHC,emmeans::emmeans,multcomp::cld,emmeans::emmip
- teste post hoc:
- teste de efeito principal simples:
- Parametric Bootstrap based Generalized Test (robusto à heterocedasticidade e à não normalidade)
- teste omnibus:
twowaytests::gpTwoWay- teste post hoc:
twowaytests::paircompTwoWay
- teste post hoc:
- ANOVA multifatorial relacionada ou para medidas repetidas
- teste omnibus:
lmerTest::lmerecar::Anova(test.statistic="F", ...),effectsize::eta_squared- teste de efeito principal simples:
phia::testInteractions- teste post hoc:
emmeans::emmeans,multcomp::cld,emmeans::emmip
- teste post hoc:
- teste de efeito principal simples:
- ANOVA multifatorial mista
- teste omnibus:
lmerTest::lmerecar::Anova(test.statistic="F", ...),effectsize::eta_squared- teste de efeito principal simples:
phia::testInteractions- teste post hoc:
emmeans::emmeans,multcomp::cld,emmeans::emmip
- teste post hoc:
- teste de efeito principal simples:
Modelo linear geral (GLM) univariado multifatorial
O modelo linear geral (GLM) univariado multifatorial nada mais é do que uma forma unificada de tratar regressão e ANOVA. Ele descreve a relação entre uma variável dependente (VD) intervalar e uma ou mais variáveis independentes (VI) (que podem ser intervalares e/ou nominais). A VD também pode ser chamada de variável de desfecho, medida (measure) ou variável de resposta. A VI pode ser chamada de fator (VI nominal) ou covariável (VI intervalar).
Em palavras simples: você observa uma VD e verificar se ela está relacionada com uma VI interesse do pesquisador.
GLM permite o controle estatístico das variáveis de confusão na análise da relação entre as variáveis dependentes (VD) e a variável independente (VI) de interesse. Esse controle é feito pela inclusão das potenciais variáveis de confusão como covariáveis e fatores no modelo, de modo que o efeito estimado da VI sobre a VD seja ajustado e condicional às demais variáveis incluídas. Assim, o coeficiente associado à VI reflete sua associação com a VD mantendo fixos as covariáveis e os fatores de controle.
A equação básica é:
\[ Y = \beta_0 + \beta_1X_1 + \beta_2X_2 + \dots + \varepsilon \]
Aqui, cada \(\beta\) mede o efeito de uma VI, mantendo as demais fixas (ceteris paribus).
Exemplo: se a VD é pressão arterial e a VI de interesse é o fator uso de medicamento (sim/não), com as VI de controle idade (covariável) e sexo (fator), o GLM testa estatisticamente o efeito do fator medicamento sobre a pressão arterial controlando por idade e sexo.
Ou seja: GLM é um modelo que responde à pergunta “qual é o efeito da VI de interesse sobre a VD, quando controlo estatisticamente por meio de outras VI (variáveis de confusão)?”.
Controlar estatisticamente significa excluir os efeitos de sexo e idade da relação entre pressão arterial e medicamento. Este controle estatistico das variáveis de confusão sexo e idade consiste numa maneira de purificar a relação entre pressão arterial e medicamento.
Portanto, GLM é um modelo explicativo em contraposição ao modelo preditivo.
Além do controle estatístico, em um estudo empírico também se pode recorrer ao controle experimental, que visa aumentar a plausibilidade da relação causal entre VD e VI. Esse controle se subdivide em estratégias voltadas à validade externa e à validade interna. A validade externa está relacionada à possibilidade de generalizar os resultados, dependendo de fatores como o método de amostragem, o processo de recrutamento e a definição clara da população-alvo por meio de critérios de inclusão e exclusão. Já a validade interna diz respeito à capacidade de assegurar que a relação observada entre VD e VI não seja espúria, o que pode ser alcançado por técnicas como a randomização, o matching (pareamento) e o controle direto de variáveis de confusão.
Assim, a combinação entre o controle estatístico proporcionado pelo modelo linear geral e os procedimentos de controle experimental fortalece tanto a validade interna quanto a externa do estudo, oferecendo maior veracidade às conclusões sobre a relação entre VD e VI de interesse.
ANOVA independente bifatorial
Suponha que conduzimos um estudo com delineamento totalmente entre participantes, balanceado e não-randomizado com 48 participantes para investigar o efeito da cafeína na habilidade de dirigir sob o efeito do álcool.
O estudo analisa os efeitos de presença e ausência de álcool e de cafeína na quantidade de infrações (erros) cometidas ao dirigir um simulador de direção veicular durante 20 minutos.
Dada a antiga crença de que o café auxilia a mantermo-nos alertas, podemos fazer algumas previsões para esse experimento:
- Altos níveis de álcool diminuem a capacidade de dirigir.
- Altos níveis de cafeína podem melhorar a habilidade de dirigir devido ao efeito estimulante.
- Um aumento de cafeína reduz a influência do álcool na habilidade de dirigir.
Se beber, tomar café para dirigir?
Conforme Howland et al. (2010),
“RESUMO: Objetivos A comercialização que promove a mistura de bebidas alcoólicas com bebidas energéticas cafeinadas (por exemplo, Red Bull com vodka) tem como alvo jovens consumidores e transmite a expectativa de que a cafeína compensará os efeitos sedativos do álcool e aumentará o estado de alerta. Tais crenças podem resultar em comportamentos de risco injustificados (por exemplo, dirigir sob efeito de álcool). O objetivo deste estudo foi avaliar os efeitos agudos de bebidas alcoólicas cafeinadas versus não cafeinadas em uma tarefa simulada de direção e no tempo de atenção/reação. Delineamento Conduzimos um ensaio randomizado entre grupos em esquema 2 × 2, no qual os participantes foram alocados em uma de quatro condições: cerveja e cerveja não alcoólica, com e sem adição de cafeína. A cafeína foi adicionada na mesma proporção encontrada em uma cerveja cafeinada comercialmente disponível (69 mg/355 ml de cerveja a 4,8% de álcool por volume). Participantes Os participantes foram 127 jovens adultos bebedores episódicos pesados, não dependentes (idade entre 21 e 30 anos), estudantes universitários ou recém-formados. O nível-alvo de álcool no ar expirado foi de 0.12% g%. Medidas O desempenho na direção foi avaliado com um simulador de direção; a atenção sustentada/reação foi avaliada com o Psychomotor Vigilance Task (PVT). Resultados Em toda a tarefa de direção e no tempo de atenção/reação encontramos efeitos principais do álcool, que prejudicaram significativamente a direção e a atenção sustentada/tempo de reação, com efeitos estatísticos em geral grandes; no entanto, a adição de cafeína não apresentou efeitos principais nem de interação sobre o desempenho. Conclusão A adição de cafeína ao álcool não parece melhorar o desempenho na direção ou no tempo de atenção/reação sustentada em comparação com o álcool isoladamente.”
Conforme Gulick & Gould (2009),
“Resumo:
Introdução A cafeína é frequentemente consumida simultaneamente ou imediatamente após o consumo de etanol. Identificar como a cafeína e o etanol interagem para modular o comportamento é essencial para compreender o uso concomitante dessas drogas. A tarefa de esquiva discriminativa em labirinto em cruz elevado (PMDAT) permite a medição intra-sujeito de aprendizagem, ansiedade e locomoção.
Métodos Para o treinamento, cada camundongo foi colocado no centro do labirinto em cruz elevado por cinco minutos e, a cada vez que o animal entrava no braço fechado aversivo, uma luz e um ruído branco eram ligados. No teste, cada camundongo era retornado ao centro do labirinto por três minutos. Nenhum estímulo foi apresentado durante o teste.
Resultados O etanol (1.0–1.4 g/kg) diminuiu de forma dose-dependente a ansiedade e a aprendizagem, e aumentou a locomoção. A cafeína (5.0–40.0 mg/kg) aumentou de forma dose-dependente a ansiedade e diminuiu a locomoção e a aprendizagem. A cafeína não conseguiu reverter os déficits de aprendizagem induzidos pelo etanol. No entanto, 1.4 g/kg de etanol bloqueou o efeito ansiogênico da cafeína.
Discussão Embora a cafeína e o etanol interajam para modular o comportamento na PMDAT, a cafeína não reverte os déficits de aprendizagem induzidos pelo etanol. A ansiólise induzida pelo etanol pode contribuir para o consumo de álcool, enquanto o bloqueio pelo etanol da ansiogênese induzida pela cafeína pode contribuir para o uso concomitante.”
Os dados estão no arquivo CafEntre_AlcEntre.rds. Foram obtidos e
adaptados de Dancey & Reidy (2006).
Dados <- data.frame(readxl::read_excel("CafEntre_AlcEntre.xlsx"))
Dados$UE <- factor(Dados$UE)
Dados$Cafeina <- factor(Dados$Cafeina,
levels=c("SemC","ComC"))
Dados$Alcool <- factor(Dados$Alcool,
levels=c("SemA","ComA"))
saveRDS(Dados, "CafEntre_AlcEntre.rds") UE Cafeina Alcool NumErros
1 1 SemC SemA 4
2 2 SemC SemA 9
3 3 SemC SemA 10
4 4 SemC SemA 8
5 5 SemC SemA 6
6 6 SemC SemA 11
7 7 SemC SemA 2
8 8 SemC SemA 11
9 9 SemC SemA 11
10 10 SemC SemA 10
11 11 SemC SemA 3
12 12 SemC SemA 10
13 13 SemC ComA 28
14 14 SemC ComA 22
15 15 SemC ComA 21
16 16 SemC ComA 27
17 17 SemC ComA 21
18 18 SemC ComA 20
19 19 SemC ComA 19
20 20 SemC ComA 16
21 21 SemC ComA 25
22 22 SemC ComA 17
23 23 SemC ComA 19
24 24 SemC ComA 20
25 25 ComC SemA 15
26 26 ComC SemA 8
27 27 ComC SemA 17
28 28 ComC SemA 0
29 29 ComC SemA 15
30 30 ComC SemA 11
31 31 ComC SemA 11
32 32 ComC SemA 6
33 33 ComC SemA 0
34 34 ComC SemA 15
35 35 ComC SemA 17
36 36 ComC SemA 15
37 37 ComC ComA 10
38 38 ComC ComA 11
39 39 ComC ComA 27
40 40 ComC ComA 15
41 41 ComC ComA 27
42 42 ComC ComA 10
43 43 ComC ComA 21
44 44 ComC ComA 15
45 45 ComC ComA 19
46 46 ComC ComA 21
47 47 ComC ComA 15
48 48 ComC ComA 15 Alcool
Cafeina SemA ComA
SemC 12 12
ComC 12 12 Alcool SemA ComA
UE Cafeina
1 SemC 1 0
ComC 0 0
2 SemC 1 0
ComC 0 0
3 SemC 1 0
ComC 0 0
4 SemC 1 0
ComC 0 0
5 SemC 1 0
ComC 0 0
6 SemC 1 0
ComC 0 0
7 SemC 1 0
ComC 0 0
8 SemC 1 0
ComC 0 0
9 SemC 1 0
ComC 0 0
10 SemC 1 0
ComC 0 0
11 SemC 1 0
ComC 0 0
12 SemC 1 0
ComC 0 0
13 SemC 0 1
ComC 0 0
14 SemC 0 1
ComC 0 0
15 SemC 0 1
ComC 0 0
16 SemC 0 1
ComC 0 0
17 SemC 0 1
ComC 0 0
18 SemC 0 1
ComC 0 0
19 SemC 0 1
ComC 0 0
20 SemC 0 1
ComC 0 0
21 SemC 0 1
ComC 0 0
22 SemC 0 1
ComC 0 0
23 SemC 0 1
ComC 0 0
24 SemC 0 1
ComC 0 0
25 SemC 0 0
ComC 1 0
26 SemC 0 0
ComC 1 0
27 SemC 0 0
ComC 1 0
28 SemC 0 0
ComC 1 0
29 SemC 0 0
ComC 1 0
30 SemC 0 0
ComC 1 0
31 SemC 0 0
ComC 1 0
32 SemC 0 0
ComC 1 0
33 SemC 0 0
ComC 1 0
34 SemC 0 0
ComC 1 0
35 SemC 0 0
ComC 1 0
36 SemC 0 0
ComC 1 0
37 SemC 0 0
ComC 0 1
38 SemC 0 0
ComC 0 1
39 SemC 0 0
ComC 0 1
40 SemC 0 0
ComC 0 1
41 SemC 0 0
ComC 0 1
42 SemC 0 0
ComC 0 1
43 SemC 0 0
ComC 0 1
44 SemC 0 0
ComC 0 1
45 SemC 0 0
ComC 0 1
46 SemC 0 0
ComC 0 1
47 SemC 0 0
ComC 0 1
48 SemC 0 0
ComC 0 1 Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 12 7.9 3.3 2 5.5 9.5 10 11 0
2 ComA SemC 12 21.2 3.7 16 19.0 20.5 23 28 0
3 SemA ComC 12 10.8 6.1 0 7.5 13.0 15 17 17
4 ComA ComC 12 17.2 5.9 10 14.0 15.0 21 27 0Este é um formato de arquivo wide, no qual cada linha contêm todos os dados de cada unidade experimental (UE). A identificação de cada unidade experimental, que no caso é um participante, é o elemento-chave que distingue os delineamentos. Neste caso, há 48 participantes e 48 linhas no arquivo.
Modelagem da ANOVA
Para abordar este tipo de situação, precisamos testar todos os efeitos principais e de interação numa única análise. Utilizaremos ANOVA independente bifatorial (two-way ANOVA), uma extensão da ANOVA independente unifatorial na qual:
- a variável dependente (VD) é intervalar.
- as variáveis independentes (VI) são fatores nominais com pelo menos dois níveis.
Dependendo do delineamento do estudo, três formas de executar a ANOVA são necessárias (independente, relacionada e mista).
Neste exemplo aplica-se a ANOVA independente bifatorial, na qual cada unidade observacional/experimental é alocada em apenas um tratamento. Cada tratamento é uma combinação dos níveis dos fatores. No exemplo acima, há quatro tratamentos:
- com cafeína e com álcool,
- com cafeína e sem álcool,
- sem cafeína e com álcool,
- sem cafeína e sem álcool.
Quando existem dois fatores há efeitos de interação e efeitos principais (de cada um dos fatores em separado). Portanto, o propósito da ANOVA é testar todos os efeitos dos fatores, que são fontes de variação para identificar quanto da variação total nos escores da VD pode ser atribuída a cada um destes efeitos.
Neste exemplo de fatores entre participantes, há quatro efeitos fixos:
- interação entre Álcool e Cafeína,
- principal de Álcool,
- principal de Cafeína,
- outros fatores não considerados no modelo, chamado de termo de erro.
Note que a forma mais correta é dizer que testamos efeitos dos fatores e não, como coloquialmente costuma-se dizer, que testam-se os fatores.
O caso deste exemplo é um modelo completo (full model) porque os quatro efeitos são considerados.
| Em capítulos subsequentes, abordaremos os modelos parciais (custom models) nos quais podemos testar apenas o efeito de interação, os efeitos principais, parte dos efeitos de interação (quando há três ou mais fatores) ou parte dos efeitos principais de interesse do pesquisador. No modelo parcial, o efeito do erro sempre aparece, contendo tudo o que não é incorporado ao modelo. |
Hipótese nula e teste F
A hipótese nula omnibus consiste na ausência de efeitos principais e de interação populacionais.
O teste do modelo como um todo (omnibus) implica em testar implicitamente as três hipóteses nulas envolvidas na ANOVA bifatorial:
\[ \begin{cases} H_{0}^\text{Alcool}: &\mu_{\text{SemA}} = \mu_{\text{ComA}}\\ H_{0}^\text{Cafeina}: &\mu_{\text{SemC}} = \mu_{\text{ComC}}\\ H_{0}^\text{Alcool:Cafeina}: &\text{ausência de efeito de interação populacional} \end{cases}\\ \alpha=5\% \]
Mais formalmente:
\[ \begin{cases} H_{0}^\text{Alcool}: &\mu_{\text{SemA}} - \mu_{\text{ComA}}=0\\ H_{0}^\text{Cafeina}: &\mu_{\text{SemC}} - \mu_{\text{ComC}}=0\\ H_{0}^\text{Alcool:Cafeina}: &\text{ausência de efeito de interação populacional} \end{cases}\\ \alpha=5\% \] Se o valor p omnibus é menor do que \(\alpha\), o modelo proposto é significante, i.e., é plausível estatisticamente.
Se o modelo é significante, a análise de significância estatística
dos efeitos do modelo é hierárquica. Ele tem início pelo nível mais alto
de efeito de interação. O próximo passo é testar o de efeito de
interação dos fatores (Alcool:Cafeina). Se o seu valor
p é menor do que \(\alpha\), a
hipótese nula de ausência de efeito populacional de interação é
rejeitada. Se o efeito de interação é significante, os valores
p dos efeitos principais não devem ser considerados, pois os
efeitos dos fatores Alcool e Cafeina não são
independentes (aditivos).
Se o seu valor p é maior ou igual do que \(\alpha\), a hipótese nula de ausência de
efeito populacional de interação não é rejeitada. Se o efeito de
interação não é significante, os valores p dos efeitos
principais podem ser considerados, pois os efeitos dos fatores são
independentes (aditivos). Nesta situação, os valores p dos
efeitos principais dos fatores Alcool e
Cafeina podem ser comparados com o nível de significância
adotado. Se o fator tem três ou mais níveis, podem ser realizados os
testes post hoc.
Podemos, apenas, considerar seus tamanhos de efeito, os \(\eta^2\) parciais, que são as porcentagens da variância da VD que são explicadas pelos respectivos efeitos. Cada \(\eta\) parcial é uma correlação de Pearson parcial absoluta implícita entre o respectivo efeito e a VD. Como há interação de efeitos, existem redundâncias que fazem com que a soma dos \(\eta^2\) parciais não seja necessariamente igual ao \(\eta^2\) omnibus.
Efeito principal simples
Se o efeito de interação Alcool:Cafeina é significante,
precisamos estudar os efeitos principais simples, que é a verificação da
significância do efeito de um fator em todos os níveis do outro
fator.
Se houve efeito de interação, a verificação dos efeitos principais
simples de Alcool e Cafeina é feita com função
R que fornece os valores p ajustados.
Portanto, podemos compará-los diretamente com \(\alpha\).
Neste exemplo, verificamos que o efeito do Alcool é
significante nos dois níveis de Cafeina. Portanto, podemos
dizer que o efeito principal do Alcool é genuíno, i.e., tem
existência per se. Neste caso, o teste post hoc pode
ser realizado sobre o efeito principal simples do fator
Alcool.
Ao contrário, o efeito da Cafeina foi não significante
nas condições com ou sem Alcool (bastaria ser não
significante em apenas uma condição do fator Alcool).
Portanto, não se trata de um efeito principal genuíno. Isto equivale a
dizer Cafeina não tem existência per se. Neste
caso, o teste post hoc não pode ser realizado sobre o efeito
principal do fator Cafeina. Pode ser realizado o teste
post hoc sobre Cafeina|Alcool:Cafeina.
Os gráficos dos intervalos de confiança de \(1-\alpha\) com as devidas correções são
pós-modelo (observe que as amplitudes são iguais). Se há efeito de
interação de Alcool e Cafeina, o gráfico de
perfis de médias não são paralelos populacionalmente.
O teste de efeito principal simples é realizado com a função
phia::testInteractions.
Tabela ANOVA
Como existe modelo, podemos verificar os efeitos na saída ANOVA, para
explicitamente podermos definir a rejeição ou não rejeição de cada uma
das hipóteses. Utilizamos o resultado de lm no formato de
ANOVA exibido com a função car::Anova, que denominaremos,
aqui, de “ANOVA da ANOVA”.
Há parâmetros para a função car::Anova. O parâmetro
type=2 é o default. Estes tipos são maneiras
diferentes de particionar a variância total da VD entre os efeitos
elencados no modelo.
Existe também o parâmetro white.adjust=TRUE serve para
executar método mais robusto, com correção para heterocedasticidade (a
mesma que utilizamos anteriormente com ANOVA unifatorial).
Existe a função anova nativa no R. Ela é mais limitada
que a car::Anova (por não oferecer correção para
heterocedasticidade) e que também utiliza tipo de teste 2 por
default (portanto, neste exemplo, mostraria resultados
idênticos). Cuidado com a função anova, pois os valores
p podem variar em função da ordem dos fatores na fórmula do
modelo. A função car::Anova faz os valores p serem
invariantes com relação à ordem dos fatores na fórmula do modelo.
Modelo Linear Geral (GLM) em R
A função lm implementa o GLM Neste caso é um General
Linear Model (GLM) porque testa efeitos fixos de fatores, tratando
NumErros em função (~) de Alcool
interagindo (*) com Cafeina.
NumErros é a VD, intervalar. Alcool e
Cafeina são os fatores nominais.
Esta forma de escrever a interação é equivalente à sintaxe na qual os efeitos são explicitamente exibidos (exceto o termo de erro, que é inerente ao modelo):
modelo <- lm(NumErros ~ Alcool+Cafeina+Alcool:Cafeina,
data=Dados)
modelo <- lm(NumErros ~ Alcool*Cafeina,
data=Dados)
modelo <- lm(NumErros ~ (Alcool+Cafeina)^2,
data=Dados)Portanto, NumErros é tratado em função dos efeitos
principais do Alcool e da Cafeina, e também em
função do efeito de interação Alcool:Cafeina.
Este é o teste omnibus, e utilizamos apenas o valor
p associado à estatística \(F\) da saída da função lm
(neste caso, rejeita-se a hipótese nula de ausência de modelo).
O valor de \(R^2\) que aparece nesta saída corresponde à estimativa pontual do tamanho de efeito omnibus a posteriori.
Observe que obtivemos esta saída estabelecendo o modelo com
lm e exibindo no formato de regressão com
summary. Vamos denominar este procedimento de “Regressão da
ANOVA”.
A regressão funciona como um teste omnibus, cuja hipótese nula é a não existência de um modelo.
Para o nível de significância adotado (\(\alpha=0.05\)), o valor p rejeita a hipótese nula (a hipótese alternativa é a violação de qualquer igualdade), portanto existe modelo, significando que a estrutura proposta detecta a existência de associação entre os efeitos fixos dos fatores e a VD.
\(R^2\) (Multiple R-squared) é a porcentagem da variância da VD explicada pelos efeitos do modelo e, portanto, é um tamanho de efeito que equivale ao \(\eta^2\) global.
|
Em geral, os trabalhos publicados não têm mais do que três fatores. Imagine um estudo com Álcool, Cafeína e Glicose, cada um com dois níveis. O total de efeitos é \(2×2×2=8\):
Genericamente, uma ANOVA \(a_1×a_2×\cdots×a_k\) tem \(2^k\) efeitos, sendo que \(a_i\) é o número de níveis do fator \(i\). Assim, \(a_1×a_2×\cdots×a_k\) é número de tratamentos/condições experimentais, divididos em:
…
Como regra heurística em ANOVA com mais de três fatores analisam-se os efeitos de interação até ordem 3, sendo que os efeitos de interação podem ser escolhidos para análise. Como foi visto, os testes de efeitos de ordem inferior só podem ser considerados se não há efeitos de ordem superior estatisticamente significantes mas podem ser considerados diretamente se o efeito de interação não é estatisticamente significante.
|
ANOVA independente bifatorial de Fisher
Testes:
- Omnibus:
lm,car::Anova,effectsize::eta_squared- Efeito principal simples:
phia::testInteractions,emmeans::emmip- Post hoc:
emmeans::emmeans,multcomp::cld,emmeans::emmip
- Post hoc:
- Efeito principal simples:
- testes omnibus e post hoc:
jmv::ANOVA
Modelo incompleto: efeitos principais
O modelo incompleto de ANOVA independente bifatorial de Fisher supõe que o efeito de interação é nulo. Dessa forma, os efeitos principais são considerados independentes e aditivos.
Portanto, o modelo incompleto é:
lm(NumErros ~ Alcool + Cafeina, data=Dados)
Descriptive Statistics
By Alcool Alcool n mean sd min Q1 median Q3 max percZero
1 SemA 24 9.38 5.04 0 6 10.0 12.00 17 8.33
2 ComA 24 19.21 5.27 10 15 19.5 21.25 28 0.00
By Cafeina Cafeina n mean sd min Q1 median Q3 max percZero
1 SemC 24 14.58 7.63 2 9.75 13.5 20.25 28 0.00
2 ComC 24 14.00 6.72 0 10.75 15.0 17.00 27 8.33
alfa <- 0.05
g1 <- nlevels(Dados$Alcool)
g2 <- nlevels(Dados$Cafeina)
alfaBonf1 <- alfa/g1
alfaBonf2 <- alfa/g2
gplots::plotmeans(NumErros~Alcool,
data=Dados,
p=1-alfaBonf1,
main="CI95% Bonferroni",
barcol="black",
connect=FALSE)
gplots::plotmeans(NumErros~Cafeina,
data=Dados,
p=1-alfaBonf2,
main="CI95% Bonferroni",
barcol="black",
connect=FALSE)
Assumptions
Normalitynorm.test <- by(data=Dados$NumErros,
INDICES=list("Alcool" = Dados$Alcool,
"Cafeina" = Dados$Cafeina),
FUN=shapiro.test)
print(norm.test, digits=5)Alcool: SemA
Cafeina: SemC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.84, p-value = 0.028
------------------------------------------------------------
Alcool: ComA
Cafeina: SemC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.93, p-value = 0.38
------------------------------------------------------------
Alcool: SemA
Cafeina: ComC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.847, p-value = 0.034
------------------------------------------------------------
Alcool: ComA
Cafeina: ComC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.9, p-value = 0.16
HomoscedasticityRegistered S3 methods overwritten by 'car':
method from
hist.boot FSA
confint.boot FSA Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 3 1.36 0.27
44 # ANOVA bifatorial independente de Fisher
cat("\nmain effects Fisher two-factor between-subjects ANOVA")
main effects Fisher two-factor between-subjects ANOVA
ANOVA tableAnova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
Alcool 1160.33 1 42.8871 4.769e-08 ***
Cafeina 4.08 1 0.1509 0.6995
Residuals 1217.50 45
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = NumErros ~ Alcool + Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-9.0833 -3.7292 0.3333 2.1875 8.5000
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.6667 1.3004 7.434 2.33e-09 ***
AlcoolComA 9.8333 1.5015 6.549 4.77e-08 ***
CafeinaComC -0.5833 1.5015 -0.388 0.699
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.201 on 45 degrees of freedom
Multiple R-squared: 0.4889, Adjusted R-squared: 0.4661
F-statistic: 21.52 on 2 and 45 DF, p-value: 2.768e-07
Call:
lm(formula = NumErros ~ Alcool + Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-9.0833 -3.7292 0.3333 2.1875 8.5000
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.6667 1.3004 7.434 2.33e-09 ***
AlcoolComA 9.8333 1.5015 6.549 4.77e-08 ***
CafeinaComC -0.5833 1.5015 -0.388 0.699
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.201 on 45 degrees of freedom
Multiple R-squared: 0.4889, Adjusted R-squared: 0.4661
F-statistic: 21.52 on 2 and 45 DF, p-value: 2.768e-07fobs <- reg$fstatistic[1]
dfn <- reg$fstatistic[2]
dfd <- reg$fstatistic[3]
fc <- qf(1-alfa, dfn, dfd)
p <- 1-pf(fobs, dfn, dfd)
if (p < 1e-4)
{
p <- sprintf("%.2e",p)
} else
{
p <- sprintf("%.4f",p)
}
f <- seq(0,1.4*max(fc,fobs),length.out=300)
densf <- df(f, dfn, dfd)
plot(f, densf,
main="Fisher Two-way ANOVA: Omnibus test",
xlab="F", ylab="Density",
lwd=2, type="l")
abline(v=fc, lty=3)
abline(v=fobs, lty=4)
legend("topright",
c("H0: there is not model",
paste("F(",dfn,",",dfd,",",1-alfa,") = ",round(fc,3),sep=""),
paste("F(",dfn,",",dfd,") = ",round(fobs,3),"\n",
"p = ",p,sep="")
),
lwd=c(2,1,1), lty=c(1,3,4),
cex=0.8)
Effect size analysisR^2 = eta^2 omnibus = 0.4889[1] "large"
(Rules: cohen1992)nef <- length(sort(anv$`F value`))
eta2 <- effectsize::eta_squared(anv,
partial=TRUE,
generalized=FALSE,
ci=1-alfa/nef,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 (partial) | 97.5% CI | interpret
----------------------------------------------------------
Alcool | 0.4880 | [0.2437, 0.6553] | large
Cafeina | 0.0033 | [0.0000, 0.1294] | very small
Gráfico de perfil de médiasg3 <- emmeans::emmip(modelo,
~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE) +
ggplot2::theme_bw()
g4 <- emmeans::emmip(modelo,
~Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE) +
ggplot2::theme_bw()
p <- (g1 | g2)
print(p)[1] TRUE
Post hoc test se efeito de interação Alcool:Cafeina não é significante
AlcoolEMM.A <- emmeans::emmeans(modelo,
specs=pairwise~"Alcool",
level=1-alfa,
adjust="holm")
print(summary(EMM.A$emmeans)) Alcool emmean SE df lower.CL upper.CL
SemA 9.38 1.06 45 7.24 11.5
ComA 19.21 1.06 45 17.07 21.3
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.5 45 -12.9 -6.81 -6.549 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95 
print(multcomp::cld(object=EMM.A$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.06 45 6.91 11.8 a
ComA 19.21 1.06 45 16.75 21.7 b
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
CafeinaEMM.C <- emmeans::emmeans(modelo,
specs=pairwise~"Cafeina",
level=1-alfa,
adjust="holm")
print(summary(EMM.C$emmeans)) Cafeina emmean SE df lower.CL upper.CL
SemC 14.6 1.06 45 12.4 16.7
ComC 14.0 1.06 45 11.9 16.1
Results are averaged over the levels of: Alcool
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 0.583 1.5 45 -2.44 3.61 0.388 0.6995
Results are averaged over the levels of: Alcool
Confidence level used: 0.95 
print(multcomp::cld(object=EMM.C$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Cafeina emmean SE df lower.CL upper.CL .group
ComC 14.0 1.06 45 11.5 16.5 a
SemC 14.6 1.06 45 12.1 17.0 a
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. Para \(\alpha=5\%\), a conclusão é que o efeito de álcool é significante, mas o efeito de cafeína não é significante.
Note que o efeito de interação entre álcool e cafeína foi desprezado indevidamente nessa análise, sendo que ele pode ser testado estatisticamente.
Modelo completo (full model): efeitos principais e de interação
O modelo completo de ANOVA independente bifatorial de Fisher supõe que o efeito de interação não é nulo, em princípio. Dessa forma, os efeitos principais podem ser dependentes.
Portanto, o modelo completo é:
lm(NumErros ~ Alcool + Cafeina + Alcool:Cafeina, data=Dados)
ou, mais consisamente:
lm(NumErros ~ Alcool*Cafeina, data=Dados)
Descriptive Statistics
By Alcool Alcool n mean sd min Q1 median Q3 max percZero
1 SemA 24 9.38 5.04 0 6 10.0 12.00 17 8.33
2 ComA 24 19.21 5.27 10 15 19.5 21.25 28 0.00
By Cafeina Cafeina n mean sd min Q1 median Q3 max percZero
1 SemC 24 14.58 7.63 2 9.75 13.5 20.25 28 0.00
2 ComC 24 14.00 6.72 0 10.75 15.0 17.00 27 8.33
By Alcool and Cafeina Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 12 7.92 3.32 2 5.5 9.5 10.25 11 0.00
2 ComA SemC 12 21.25 3.72 16 19.0 20.5 22.75 28 0.00
3 SemA ComC 12 10.83 6.12 0 7.5 13.0 15.00 17 16.67
4 ComA ComC 12 17.17 5.92 10 14.0 15.0 21.00 27 0.00
alfa <- 0.05
g1 <- nlevels(Dados$Alcool)
g2 <- nlevels(Dados$Cafeina)
alfaBonf <- alfa/(g1*g2)
gplots::plotmeans(NumErros~interaction(Cafeina, Alcool),
data=Dados,
p=1-alfaBonf,
main="CI95% Bonferroni",
barcol="black",
connect=FALSE)
Assumptions
Normalitynorm.test <- by(data=Dados$NumErros,
INDICES=list("Alcool" = Dados$Alcool,
"Cafeina" = Dados$Cafeina),
FUN=shapiro.test)
print(norm.test, digits=5)Alcool: SemA
Cafeina: SemC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.84, p-value = 0.028
------------------------------------------------------------
Alcool: ComA
Cafeina: SemC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.93, p-value = 0.38
------------------------------------------------------------
Alcool: SemA
Cafeina: ComC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.847, p-value = 0.034
------------------------------------------------------------
Alcool: ComA
Cafeina: ComC
Shapiro-Wilk normality test
data: dd[x, ]
W = 0.9, p-value = 0.16
HomoscedasticityLevene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 3 1.36 0.27
44
Fisher two-factor between-subjects ANOVA
ANOVA tableAnova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
Alcool 1160.33 1 47.6924 1.57e-08 ***
Cafeina 4.08 1 0.1678 0.68403
Alcool:Cafeina 147.00 1 6.0420 0.01798 *
Residuals 1070.50 44
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1modelo_aditivo <- lm(NumErros ~ Alcool+Cafeina,
data=Dados)
print(car::Anova(modelo_aditivo), digits=5)Anova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
Alcool 1160.33 1 42.8871 4.769e-08 ***
Cafeina 4.08 1 0.1509 0.6995
Residuals 1217.50 45
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analysis of Variance Table
Model 1: NumErros ~ Alcool + Cafeina
Model 2: NumErros ~ Alcool * Cafeina
Res.Df RSS Df Sum of Sq F Pr(>F)
1 45 1217.5
2 44 1070.5 1 147 6.042 0.01798 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = NumErros ~ Alcool * Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-10.833 -2.396 0.125 3.771 9.833
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.917 1.424 5.560 1.49e-06 ***
AlcoolComA 13.333 2.014 6.621 4.11e-08 ***
CafeinaComC 2.917 2.014 1.448 0.155
AlcoolComA:CafeinaComC -7.000 2.848 -2.458 0.018 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 4.932 on 44 degrees of freedom
Multiple R-squared: 0.5506, Adjusted R-squared: 0.5199
F-statistic: 17.97 on 3 and 44 DF, p-value: 9.277e-08
Call:
lm(formula = NumErros ~ Alcool * Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-10.833 -2.396 0.125 3.771 9.833
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.917 1.424 5.560 1.49e-06 ***
AlcoolComA 13.333 2.014 6.621 4.11e-08 ***
CafeinaComC 2.917 2.014 1.448 0.155
AlcoolComA:CafeinaComC -7.000 2.848 -2.458 0.018 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 4.932 on 44 degrees of freedom
Multiple R-squared: 0.5506, Adjusted R-squared: 0.5199
F-statistic: 17.97 on 3 and 44 DF, p-value: 9.277e-08fobs <- reg$fstatistic[1]
dfn <- reg$fstatistic[2]
dfd <- reg$fstatistic[3]
fc <- qf(1-alfa, dfn, dfd)
p <- 1-pf(fobs, dfn, dfd)
if (p < 1e-4)
{
p <- sprintf("%.2e",p)
} else
{
p <- sprintf("%.4f",p)
}
f <- seq(0,1.4*max(fc,fobs),length.out=300)
densf <- df(f, dfn, dfd)
plot(f, densf,
main="Fisher Two-way ANOVA: Omnibus test",
xlab="F", ylab="Density",
lwd=2, type="l")
abline(v=fc, lty=3)
abline(v=fobs, lty=4)
legend("topright",
c("H0: there is not model",
paste("F(",dfn,",",dfd,",",1-alfa,") = ",round(fc,3),sep=""),
paste("F(",dfn,",",dfd,") = ",round(fobs,3),"\n",
"p = ",p,sep="")
),
lwd=c(2,1,1), lty=c(1,3,4),
cex=0.8)
Effect size analysisR^2 = eta^2 omnibus = 0.5506[1] "large"
(Rules: cohen1992)nef <- length(sort(anv$`F value`))
eta2 <- effectsize::eta_squared(anv,
partial=TRUE,
generalized=FALSE,
ci=1-alfa/nef,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 (partial) | 98.3333% CI | interpret
---------------------------------------------------------------
Alcool | 0.5201 | [0.2588, 0.6884] | large
Cafeina | 0.0038 | [0.0000, 0.0077] | very small
Alcool:Cafeina | 0.1207 | [0.0000, 0.3549] | medium
Gráfico de perfil de médiasg1 <- emmeans::emmip(object=modelo,
formula=~Cafeina ~Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
dodge=.3) +
ggplot2::theme_bw()
g2 <- emmeans::emmip(object=modelo,
formula=~Alcool ~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
dodge=.3) +
ggplot2::theme_bw()
g3 <- emmeans::emmip(modelo,
~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE) +
ggplot2::theme_bw()NOTE: Results may be misleading due to involvement in interactionsNOTE: Results may be misleading due to involvement in interactions
Post hoc test se efeito de interação Alcool:Cafeina não é significante
AlcoolNOTE: Results may be misleading due to involvement in interactions Alcool emmean SE df lower.CL upper.CL
SemA 9.38 1.01 44 7.35 11.4
ComA 19.21 1.01 44 17.18 21.2
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.42 44 -12.7 -6.96 -6.906 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95 
print(multcomp::cld(object=EMM.A$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.01 44 7.04 11.7 a
ComA 19.21 1.01 44 16.87 21.5 b
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
CafeinaNOTE: Results may be misleading due to involvement in interactions Cafeina emmean SE df lower.CL upper.CL
SemC 14.6 1.01 44 12.6 16.6
ComC 14.0 1.01 44 12.0 16.0
Results are averaged over the levels of: Alcool
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 0.583 1.42 44 -2.29 3.45 0.410 0.6840
Results are averaged over the levels of: Alcool
Confidence level used: 0.95 
print(multcomp::cld(object=EMM.C$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Cafeina emmean SE df lower.CL upper.CL .group
ComC 14.0 1.01 44 11.7 16.3 a
SemC 14.6 1.01 44 12.2 16.9 a
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Post hoc teste se efeito de interação Alcool:Cafeina é significante
Post hoc test: interaction Cafeina:Alcoolplot(emmeans::emmip(modelo,
~Cafeina:Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE) +
ggplot2::theme_bw())
EMM.AC <- emmeans::emmeans(modelo,
specs=pairwise~"Alcool:Cafeina",
level=1-alfa,
adjust="holm")
print(summary(EMM.AC$emmeans)) Alcool Cafeina emmean SE df lower.CL upper.CL
SemA SemC 7.92 1.42 44 5.05 10.8
ComA SemC 21.25 1.42 44 18.38 24.1
SemA ComC 10.83 1.42 44 7.96 13.7
ComA ComC 17.17 1.42 44 14.30 20.0
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA SemC - ComA SemC -13.33 2.01 44 -18.90 -7.77 -6.621 <.0001
SemA SemC - SemA ComC -2.92 2.01 44 -8.48 2.65 -1.448 0.1546
SemA SemC - ComA ComC -9.25 2.01 44 -14.81 -3.69 -4.594 0.0001
ComA SemC - SemA ComC 10.42 2.01 44 4.85 15.98 5.173 <.0001
ComA SemC - ComA ComC 4.08 2.01 44 -1.48 9.65 2.028 0.0973
SemA ComC - ComA ComC -6.33 2.01 44 -11.90 -0.77 -3.145 0.0089
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 6 estimates
P value adjustment: holm method for 6 tests 
print(multcomp::cld(object=EMM.AC$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool Cafeina emmean SE df lower.CL upper.CL .group
SemA SemC 7.92 1.42 44 4.21 11.6 a
SemA ComC 10.83 1.42 44 7.12 14.5 a
ComA ComC 17.17 1.42 44 13.46 20.9 b
ComA SemC 21.25 1.42 44 17.54 25.0 b
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 4 estimates
P value adjustment: holm method for 6 tests
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Alcool effect by Cafeina levelsF Test:
P-value adjustment method: holm
Value SE Df Sum of Sq F Pr(>F)
SemC -13.3333 2.0137 1 1066.67 43.842 8.225e-08 ***
ComC -6.3333 2.0137 1 240.67 9.892 0.002973 **
Residuals 44 1070.50
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(emmeans::emmip(modelo,
~Alcool|Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE) +
ggplot2::theme_bw()) 
Post hoc test of the simple main effect of AlcoolNOTE: Results may be misleading due to involvement in interactions Alcool emmean SE df lower.CL upper.CL t.ratio p.value
SemA 9.38 1.01 44 7.35 11.4 9.311 <.0001
ComA 19.21 1.01 44 17.18 21.2 19.078 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.42 44 -12.7 -6.96 -6.906 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
print(multcomp::cld(object=model_means$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.01 44 7.04 11.7 a
ComA 19.21 1.01 44 16.87 21.5 b
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Cafeina effect by Alcool levelsF Test:
P-value adjustment method: holm
Value SE Df Sum of Sq F Pr(>F)
SemA -2.9167 2.0137 1 51.04 2.0979 0.15459
ComA 4.0833 2.0137 1 100.04 4.1119 0.09733 .
Residuals 44 1070.50
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(emmeans::emmip(modelo,
~Cafeina|Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE) +
ggplot2::theme_bw())
Post hoc test of the simple main effect of Cafeinamodel_means <- emmeans::emmeans(object=modelo,
specs=pairwise~Cafeina|Alcool:Cafeina,
level=1-alfa,
adjust="holm")
print(summary(model_means$emmeans, infer=TRUE))Alcool = SemA:
Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 7.92 1.42 44 5.05 10.8 5.560 <.0001
ComC 10.83 1.42 44 7.96 13.7 7.608 <.0001
Alcool = ComA:
Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 21.25 1.42 44 18.38 24.1 14.924 <.0001
ComC 17.17 1.42 44 14.30 20.0 12.056 <.0001
Confidence level used: 0.95 Alcool = SemA:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC -2.92 2.01 44 -6.975 1.14 -1.448 0.1546
Alcool = ComA:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 4.08 2.01 44 0.025 8.14 2.028 0.0487
Confidence level used: 0.95 
print(multcomp::cld(object=model_means$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters))Alcool = SemA:
Cafeina emmean SE df lower.CL upper.CL .group
SemC 7.92 1.42 44 4.61 11.2 a
ComC 10.83 1.42 44 7.53 14.1 a
Alcool = ComA:
Cafeina emmean SE df lower.CL upper.CL .group
ComC 17.17 1.42 44 13.86 20.5 a
SemC 21.25 1.42 44 17.95 24.6 b
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. # jmv::ANOVA: Fisher
print(jmv::ANOVA(dep = "NumErros",
factors = c("Alcool", "Cafeina"),
modelTes = TRUE,
postHoc= NumErros ~ Alcool+Cafeina+Alcool:Cafeina,
emMeans = NumErros ~ Alcool+Cafeina+Alcool:Cafeina,
emmPlots = TRUE,
emmTables = TRUE,
data = Dados))NOTE: Results may be misleading due to involvement in interactions
NOTE: Results may be misleading due to involvement in interactions
ANOVA
ANOVA - NumErros
─────────────────────────────────────────────────────────────────────────────────────
Sum of Squares df Mean Square F p
─────────────────────────────────────────────────────────────────────────────────────
Overall model 1311.416667 3 437.138889 17.9674088 < .0000001
Alcool 1160.333333 1 1160.333333 47.6923556 < .0000001
Cafeina 4.083333 1 4.083333 0.1678343 0.6840314
Alcool:Cafeina 147.000000 1 147.000000 6.0420364 0.0179772
Residuals 1070.500000 44 24.329545
─────────────────────────────────────────────────────────────────────────────────────
POST HOC TESTS
Post Hoc Comparisons - Alcool
───────────────────────────────────────────────────────────────────────────────────────────────
Alcool Alcool Mean Difference SE df t p-tukey
───────────────────────────────────────────────────────────────────────────────────────────────
SemA - ComA -9.833333 1.423890 44.00000 -6.905965 < .0000001
───────────────────────────────────────────────────────────────────────────────────────────────
Note. Comparisons are based on estimated marginal means
Post Hoc Comparisons - Cafeina
────────────────────────────────────────────────────────────────────────────────────────────────
Cafeina Cafeina Mean Difference SE df t p-tukey
────────────────────────────────────────────────────────────────────────────────────────────────
SemC - ComC 0.5833333 1.423890 44.00000 0.4096759 0.6840314
────────────────────────────────────────────────────────────────────────────────────────────────
Note. Comparisons are based on estimated marginal means
Post Hoc Comparisons - Alcool:Cafeina
────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Alcool Cafeina Alcool Cafeina Mean Difference SE df t p-tukey
────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
SemA SemC - SemA ComC -2.916667 2.013684 44.00000 -1.448423 0.4767730
- ComA SemC -13.333333 2.013684 44.00000 -6.621362 0.0000002
- ComA ComC -9.250000 2.013684 44.00000 -4.593570 0.0002075
ComC - ComA SemC -10.416667 2.013684 44.00000 -5.172939 0.0000315
- ComA ComC -6.333333 2.013684 44.00000 -3.145147 0.0151866
ComA SemC - ComA ComC 4.083333 2.013684 44.00000 2.027792 0.1934857
────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Note. Comparisons are based on estimated marginal means
ESTIMATED MARGINAL MEANS
ALCOOL
Estimated Marginal Means - Alcool
────────────────────────────────────────────────────────────
Alcool Mean SE Lower Upper
────────────────────────────────────────────────────────────
SemA 9.375000 1.006842 7.345843 11.40416
ComA 19.208333 1.006842 17.179176 21.23749
────────────────────────────────────────────────────────────
CAFEINA
Estimated Marginal Means - Cafeina
───────────────────────────────────────────────────────────
Cafeina Mean SE Lower Upper
───────────────────────────────────────────────────────────
SemC 14.58333 1.006842 12.55418 16.61249
ComC 14.00000 1.006842 11.97084 16.02916
───────────────────────────────────────────────────────────
ALCOOL:CAFEINA
Estimated Marginal Means - Alcool:Cafeina
───────────────────────────────────────────────────────────────────────
Cafeina Alcool Mean SE Lower Upper
───────────────────────────────────────────────────────────────────────
SemC SemA 7.916667 1.423890 5.047005 10.78633
ComA 21.250000 1.423890 18.380339 24.11966
ComC SemA 10.833333 1.423890 7.963672 13.70299
ComA 17.166667 1.423890 14.297005 20.03633
─────────────────────────────────────────────────────────────────────── 
Conclusão:
“Se beber, tomar café para dirigir?”
Adotando \(\alpha=5\%\), temos os seguintes resultados.
O efeito de interação é significante.
O efeito principal simples de Alcool é significante.
Portanto, o resultado do teste post hoc de Alcool
é:
Alcool emmean SE df lower.CL upper.CL .group SemA 9.38 1.01 44 7.04 11.7 a ComA 19.21 1.01 44 16.87 21.5 b
O efeito principal simples de Cafeina não é
significante. Portanto, o resultado do teste post hoc de
Cafeina|Alcool:Cafeina é:
Alcool = SemA: Cafeina emmean SE df lower.CL upper.CL .group SemC 7.92 1.42 44 4.61 11.2 a ComC 10.83 1.42 44 7.53 14.1 a Alcool = ComA: Cafeina emmean SE df lower.CL upper.CL .group ComC 17.17 1.42 44 13.86 20.5 a SemC 21.25 1.42 44 17.95 24.6 b
Em função dos resultados acima, pode-se concluir que:
Sem álcool:
- Médias: 7.92 (SemC) vs 10.83 (ComC).
- Mesma letra (“a”): logo, não há diferença significante entre tomar ou não café: portanto, café sozinho não altera o desempenho.
Com álcool:
- Médias: 21.25 (SemC) vs 17.17 (ComC).
- Letras diferentes (“b” vs “a”): logo, há diferença significante: portanto, com café, o desempenho foi melhor (menor prejuízo).
Comparando situações globais:
- Todos os valores com álcool são piores do que sem álcool.
- O café atenua o prejuízo do álcool, mas não restaura o desempenho ao nível sem álcool.
Portanto, tomar café depois de beber melhora um pouco o desempenho, mas não resolve: dirigir após álcool continua significantemente pior do que dirigir sem álcool. O café não é antídoto, apenas reduz parcialmente os efeitos negativos.
Teste do efeito de interação por bootstrapping
Dados <- readRDS("CafEntre_AlcEntre.rds")
alfa <- 0.05
B <- 1e5
modelo_boot <- lmboot::ANOVA.boot(NumErros ~ Alcool*Cafeina,
B=B,
type="residual",
wild.dist="normal",
seed=123,
data=Dados,
keep.boot.resp=FALSE)
Fintp <- ncol(modelo_boot$bootFStats)
Fc <- as.numeric(quantile(modelo_boot$bootFStats[Fintp],probs = 1-alfa))
Fobs <- qf(1-modelo_boot$`p-values`[Fintp], modelo_boot$df[Fintp],
modelo_boot$df[Fintp+1])
d <- density(modelo_boot$bootFStats)
plot(d,
main=paste("Independent two-way ANOVA: interaction effect\n",B,
" replicates", sep=""),
xlim=c(0, max(Fobs, Fc)*1.1),
xlab="F", ylab="Density", lwd=2)
abline(v=c(Fobs,Fc), lty=c(4,3))
legend("topright",
c(expression(H[0]),
paste("F(",modelo_boot$df[Fintp],", ",modelo_boot$df[Fintp+1],
"; ",1-alfa,") = ",round(Fc,2),sep=""),
paste("\nF(",modelo_boot$df[Fintp],", ",modelo_boot$df[Fintp+1],") = ",
round(Fobs,2),"\n",
"p = ",round(modelo_boot$`p-values`[Fintp],5),sep="")
),
lwd=c(2,1,1), lty=c(1,3,4))
cat(paste("Interaction effect test:\nF(",modelo_boot$df[Fintp],", ",modelo_boot$df[Fintp+1],") = ",
round(Fobs,2),
", p = ",round(modelo_boot$`p-values`[Fintp],5),"\n", sep=""))Interaction effect test:
F(1, 44) = 5.96, p = 0.01868(1e+05 bootstrap samples)Effect size analysiseta2 <- modelo_boot$df[Fintp]*Fobs/(modelo_boot$df[Fintp]*Fobs+modelo_boot$df[Fintp+1])
cat("\neta^2 = ", round(eta2,2), sep="")
eta^2 = 0.12es <- effectsize::interpret_eta_squared(eta2)
names(es) <- c("\nTamanho de efeito: estimativa pontual")
print(es)\nTamanho de efeito: estimativa pontual
"medium"
(Rules: field2013)ANOVA independente bifatorial de Fisher desbalanceada
UE Cafeina Alcool NumErros
2 2 SemC SemA 9
3 3 SemC SemA 10
4 4 SemC SemA 8
5 5 SemC SemA 6
6 6 SemC SemA 11
7 7 SemC SemA 2
8 8 SemC SemA 11
9 9 SemC SemA 11
10 10 SemC SemA 10
11 11 SemC SemA 3
12 12 SemC SemA 10
13 13 SemC ComA 28
14 14 SemC ComA 22
15 15 SemC ComA 21
16 16 SemC ComA 27
17 17 SemC ComA 21
18 18 SemC ComA 20
19 19 SemC ComA 19
20 20 SemC ComA 16
21 21 SemC ComA 25
22 22 SemC ComA 17
23 23 SemC ComA 19
24 24 SemC ComA 20
25 25 ComC SemA 15
26 26 ComC SemA 8
27 27 ComC SemA 17
28 28 ComC SemA 0
29 29 ComC SemA 15
30 30 ComC SemA 11
31 31 ComC SemA 11
32 32 ComC SemA 6
33 33 ComC SemA 0
34 34 ComC SemA 15
35 35 ComC SemA 17
36 36 ComC SemA 15
37 37 ComC ComA 10
38 38 ComC ComA 11
39 39 ComC ComA 27
40 40 ComC ComA 15
41 41 ComC ComA 27
42 42 ComC ComA 10
43 43 ComC ComA 21
44 44 ComC ComA 15
45 45 ComC ComA 19
46 46 ComC ComA 21
47 47 ComC ComA 15
48 48 ComC ComA 15 Alcool
Cafeina SemA ComA
SemC 11 12
ComC 12 12 Alcool SemA ComA
UE Cafeina
1 SemC 0 0
ComC 0 0
2 SemC 1 0
ComC 0 0
3 SemC 1 0
ComC 0 0
4 SemC 1 0
ComC 0 0
5 SemC 1 0
ComC 0 0
6 SemC 1 0
ComC 0 0
7 SemC 1 0
ComC 0 0
8 SemC 1 0
ComC 0 0
9 SemC 1 0
ComC 0 0
10 SemC 1 0
ComC 0 0
11 SemC 1 0
ComC 0 0
12 SemC 1 0
ComC 0 0
13 SemC 0 1
ComC 0 0
14 SemC 0 1
ComC 0 0
15 SemC 0 1
ComC 0 0
16 SemC 0 1
ComC 0 0
17 SemC 0 1
ComC 0 0
18 SemC 0 1
ComC 0 0
19 SemC 0 1
ComC 0 0
20 SemC 0 1
ComC 0 0
21 SemC 0 1
ComC 0 0
22 SemC 0 1
ComC 0 0
23 SemC 0 1
ComC 0 0
24 SemC 0 1
ComC 0 0
25 SemC 0 0
ComC 1 0
26 SemC 0 0
ComC 1 0
27 SemC 0 0
ComC 1 0
28 SemC 0 0
ComC 1 0
29 SemC 0 0
ComC 1 0
30 SemC 0 0
ComC 1 0
31 SemC 0 0
ComC 1 0
32 SemC 0 0
ComC 1 0
33 SemC 0 0
ComC 1 0
34 SemC 0 0
ComC 1 0
35 SemC 0 0
ComC 1 0
36 SemC 0 0
ComC 1 0
37 SemC 0 0
ComC 0 1
38 SemC 0 0
ComC 0 1
39 SemC 0 0
ComC 0 1
40 SemC 0 0
ComC 0 1
41 SemC 0 0
ComC 0 1
42 SemC 0 0
ComC 0 1
43 SemC 0 0
ComC 0 1
44 SemC 0 0
ComC 0 1
45 SemC 0 0
ComC 0 1
46 SemC 0 0
ComC 0 1
47 SemC 0 0
ComC 0 1
48 SemC 0 0
ComC 0 1 Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 11 8.3 3.2 2 7.0 10 10 11 0
2 ComA SemC 12 21.2 3.7 16 19.0 20 23 28 0
3 SemA ComC 12 10.8 6.1 0 7.5 13 15 17 17
4 ComA ComC 12 17.2 5.9 10 14.0 15 21 27 0
Fisher two-factor between-subjects ANOVA# lm: Fisher
modelo <- lm(NumErros ~ Alcool*Cafeina,
data=Dados)
print(anv <- car::Anova(modelo), digits=3)Anova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
Alcool 1078 1 43.98 4.4e-08 ***
Cafeina 8 1 0.33 0.566
Alcool:Cafeina 129 1 5.28 0.026 *
Residuals 1054 43
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = NumErros ~ Alcool * Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-10.8333 -2.2614 0.1667 3.7917 9.8333
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 8.273 1.493 5.543 1.69e-06 ***
AlcoolComA 12.977 2.066 6.280 1.43e-07 ***
CafeinaComC 2.561 2.066 1.239 0.2220
AlcoolComA:CafeinaComC -6.644 2.890 -2.299 0.0264 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 4.95 on 43 degrees of freedom
Multiple R-squared: 0.5366, Adjusted R-squared: 0.5042
F-statistic: 16.59 on 3 and 43 DF, p-value: 2.619e-07Planejamento de ANOVA independente bifatorial de Fisher
Nas condições deste exemplo (serve apenas para ANOVA independente
bifatorial), podemos saber o poder a priori usando a função
pwr2::pwr.2way:
Balanced two-way analysis of variance power calculation
a = 2
b = 2
n.A = 24
n.B = 24
sig.level = 0.05
power.A = 0.6787416
power.B = 0.6787416
power = 0.6787416
NOTE: power is the minimum power among two factorsMatemática da ANOVA independente bifatorial de Fisher
O objetivo é produzir o teste F omnibus na ANOVA independente bifatorial com interação de Fisher.
Para testar todos os efeitos principais e de interações populacionais produzidos pelos dois fatores entre participantes, com base no teste F, os graus de liberdade do numerador dependem do número de condições independentes e os do denominador dependem do tamanho da amostra do número de condições independentes.
Modelo bifatorial balanceado com interação (efeitos fixos):
\[ Y_{ijk} = \mu + \tau^{A}_{i} + \tau^{B}_{j} + \tau^{AB}_{ij} + \varepsilon_{ijk} \\\\ i=1,2,\ldots,a\\ j=1,2,\ldots,b\\ k=1,2,\ldots,n\\ n=\sum_{i=1}^{a}{\sum_{j=1}^{b}{n_{ij}}} \]
onde \(\mu\) é a média geral, \(\tau^{A}_{i}\) é o efeito do nível \(i\) do fator \(A\), \(\tau^{B}_{j}\) é o efeito do nível \(j\) do fator \(B\), \(\tau^{AB}_{ij}\) é o efeito de interação de \(A \times B\), e \(\varepsilon_{ijk} \sim \mathcal{N}\text{IID}(0,\sigma^2)\).
Restrições usuais:
\[ \sum_{i=1}^a\tau^{A}_{i}=0,\qquad \sum_{j=1}^b\tau^{B}_{j}=0,\qquad \sum_{i=1}^a\tau^{AB}_{ij}=0=\sum_{j=1}^b\tau^{AB}_{ij},\quad \forall i,j \]
Hipótese nula omnibus:
\[ H_0^{\text{omnibus}}: \tau^{A}_{\forall i} = \tau^{B}_{\forall j} = \tau^{AB}_{\forall ij} = 0 \]
Ou seja, todos os efeitos populacionais do GLM univariado (principais e interação) são nulos.
Estatística de teste:
A soma de quadrados associada ao modelo completo (sem intercepto) é
\[ \text{SQ}_{\text{modelo}} = \text{SQ}_A + \text{SQ}_B + \text{SQ}_{AB} \]
Logo, a estatística de teste F omnibus é
\[ F^{\text{omnibus}} = \dfrac{\dfrac{\text{SQ}_A + \text{SQ}_B + \text{SQ}_{AB}}{df_A+df_B+df_{AB}}} {\dfrac{\text{SQ}_E}{df_E}} \]
Sendo que os números de graus de liberdade são:
\[ \begin{align} df_A &= a-1\\ df_B &= b-1\\ df_{AB} &= (a-1)(b-1)\\ df_T &= n - 1\\ df_E &= df_T -(df_A+df_B+df_{AB})=n - ab \end{align} \]
Portanto,
\[ \begin{align} F^{\text{omnibus}} & \underset{H_0}{\sim} F_{\,df_A+df_B+df_{AB},\,df_E}\\ F^{\text{omnibus}} &\sim F_{\,ab-1,\ n-ab} \end{align} \]
Somas de quadrados (balanceado com \(n\) por célula; use médias amostrais \(\bar Y_{ij\cdot}\), \(\bar Y_{i\cdot\cdot}\), \(\bar Y_{\cdot j\cdot}\), \(\bar Y_{\cdot\cdot\cdot}\)):
\[ \text{SQ}_A = nb\sum_{i=1}^a\left(\bar Y_{i\cdot\cdot}-\bar Y_{\cdot\cdot\cdot}\right)^2\qquad \text{SQ}_B = na\sum_{j=1}^b\left(\bar Y_{\cdot j\cdot}-\bar Y_{\cdot\cdot\cdot}\right)^2 \]
\[ \text{SQ}_{AB} = n\sum_{i=1}^a\sum_{j=1}^b\left(\bar Y_{ij\cdot}-\bar Y_{i\cdot\cdot}-\bar Y_{\cdot j\cdot}+\bar Y_{\cdot\cdot\cdot}\right)^2\qquad \text{SQ}_E = \sum_{i,j}\sum_{k=1}^{n}\left(Y_{ijk}-\bar Y_{ij\cdot}\right)^2 \]
Valores esperados dos quadrados médios em termos de efeitos populacionais:
\[ \begin{align} \mathbb{E}[\text{MS}_E]&=\sigma^2_{Y}\\ \mathbb{E}[\text{MS}_A]&=\sigma^2_Y+\frac{nb}{a-1}\sum_{i=1}^{a}{\left(\tau^{A}_{i}\right)^2}\\ \mathbb{E}[\text{MS}_B]&=\sigma^2_Y+\frac{na}{b-1}\sum_{j=1}^{b}{\left({\tau^{B}_{j}}\right)^2}\\ \mathbb{E}[\text{MS}_{AB}]&=\sigma^2_Y+\frac{n}{(a-1)(b-1)}\sum_{i=1}^a\sum_{j=1}^{b}{\left(\tau^{AB}_{ij}\right)^2} \end{align} \]
Se \(df_E \ge 30\), pode-se recorrer ao Teorema Central do Limite, dispensando-se a suposição de multinormalidade para a validade do teste de hipótese nula omnibus.
Se \(df_E < 30\), pode-se recorrer ao argumento de que a distribuição multinormal em cada condição independente é uma condição suficiente para a validade do teste de hipótese nula omnibus. Nenhum teste estatístico de normalidade multivariada ou univariada é prova de normalidade, sendo possível apenas desprovar a normalidade por meio de teste estatístico. Dessa forma, mesmo em situação de amostra pequena, o teste omnibus pode ser válido.
Dados <- readRDS("CafEntre_AlcEntre.rds")
modelo <- lm(NumErros ~ Alcool*Cafeina,
data=Dados)
print(anv <- car::Anova(modelo), digits=4)Anova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
Alcool 1160.3 1 47.692 1.57e-08 ***
Cafeina 4.1 1 0.168 0.684
Alcool:Cafeina 147.0 1 6.042 0.018 *
Residuals 1070.5 44
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1nef <- nrow(anv)
SS_M <- sum(anv$`Sum Sq`[1:(nef-1)])
df_M <- sum(anv$Df[1:(nef-1)])
MS_M <- SS_M/df_M
SS_E <- sum(anv$`Sum Sq`[nef])
df_E <- sum(anv$Df[nef])
MS_E <- SS_E/df_E
F_omnibus <- MS_M/MS_E
cat("\nomnibus F(", df_M, ", ", df_E, ") = ", round(F_omnibus,2), sep="")
omnibus F(3, 44) = 17.97
Call:
lm(formula = NumErros ~ Alcool * Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-10.833 -2.396 0.125 3.771 9.833
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.917 1.424 5.560 1.49e-06 ***
AlcoolComA 13.333 2.014 6.621 4.11e-08 ***
CafeinaComC 2.917 2.014 1.448 0.155
AlcoolComA:CafeinaComC -7.000 2.848 -2.458 0.018 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 4.932 on 44 degrees of freedom
Multiple R-squared: 0.5506, Adjusted R-squared: 0.5199
F-statistic: 17.97 on 3 and 44 DF, p-value: 9.277e-08ANOVA independente bifatorial de Fisher-White
Testes:
- teste omnibus:
lm,car::Anova(white.adjust="hc2", ...)),estimatr::lm_robust(se_type="HC2", ...),effectsize::eta_squared- testes post hoc:
sandwich::vcovHC,emmeans::emmeans,emmeans::contrast,multcomp::cld
- testes post hoc:
ANOVA independente bifatorial de Fisher-White é um teste robusto à heterocedasticidade. Não há testes de efeito principal simples e post hoc em R.
O código completo está implementado em demo_ANOVA2f_Indep_Fisher-White_CafAlc.R.
Fisher-White two-factor between-subjects ANOVA
Statistical analysis: omnibus test
ANOVACoefficient covariances computed by hccm()Analysis of Deviance Table (Type II tests)
Response: NumErros
Df F Pr(>F)
Alcool 1 86.504 5.951e-12 ***
Cafeina 1 0.158 0.69297
Alcool:Cafeina 1 6.042 0.01798 *
Residuals 44
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
estimatr::lm_robust(formula = NumErros ~ Alcool * Cafeina, data = Dados,
se_type = "HC2")
Standard error type: HC2
Coefficients:
Estimate Std. Error t value Pr(>|t|) CI Lower CI Upper
(Intercept) 7.92 0.957 8.27 1.65e-10 5.99 9.85
AlcoolComA 13.33 1.439 9.27 6.59e-12 10.43 16.23
CafeinaComC 2.92 2.009 1.45 1.54e-01 -1.13 6.96
AlcoolComA:CafeinaComC -7.00 2.848 -2.46 1.80e-02 -12.74 -1.26
DF
(Intercept) 44
AlcoolComA 44
CafeinaComC 44
AlcoolComA:CafeinaComC 44
Multiple R-squared: 0.551 , Adjusted R-squared: 0.52
F-statistic: 30.9 on 3 and 44 DF, p-value: 6.73e-11
R^2 = eta^2 omnibus = 0.55[1] "large"
(Rules: field2013)
# Effect Size for ANOVA (Type I)
Parameter | Eta2 (partial) | 98.3333% CI
--------------------------------------------------
Alcool | 0.5201 | [0.2588, 0.6884]
Cafeina | 0.0038 | [0.0000, 0.0077]
Alcool:Cafeina | 0.1207 | [0.0000, 0.3549]
Gráfico de perfil de médiasNOTE: Results may be misleading due to involvement in interactionsNOTE: Results may be misleading due to involvement in interactions
Post hoc test se efeito de interação Alcool:Cafeina não é significante
AlcoolNOTE: Results may be misleading due to involvement in interactions Alcool emmean SE df lower.CL upper.CL
SemA 9.38 1.00 44 7.35 11.4
ComA 19.21 1.01 44 17.17 21.2
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.42 44 -12.7 -6.96 -6.906 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95 
Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.00 44 7.04 11.7 a
ComA 19.21 1.01 44 16.87 21.6 b
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
CafeinaNOTE: Results may be misleading due to involvement in interactions Cafeina emmean SE df lower.CL upper.CL
SemC 14.6 0.719 44 13.1 16.0
ComC 14.0 1.230 44 11.5 16.5
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 0.583 1.42 44 -2.29 3.45 0.410 0.6840
Results are averaged over the levels of: Alcool
Confidence level used: 0.95 
Cafeina emmean SE df lower.CL upper.CL .group
ComC 14.0 1.230 44 11.1 16.9 a
SemC 14.6 0.719 44 12.9 16.3 a
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Post hoc teste se efeito de interação Alcool:Cafeina é significante
Post hoc test: interaction Cafeina:Alcool
Alcool Cafeina emmean SE df lower.CL upper.CL
SemA SemC 7.92 0.957 44 5.99 9.85
ComA SemC 21.25 1.070 44 19.09 23.41
SemA ComC 10.83 1.770 44 7.27 14.39
ComA ComC 17.17 1.710 44 13.72 20.61
Confidence level used: 0.95
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA,SemC - ComA,SemC -13.33 1.44 44 -17.31 -9.359 -9.269 <.0001
SemA,SemC - SemA,ComC -2.92 2.01 44 -8.47 2.633 -1.452 0.1536
SemA,SemC - ComA,ComC -9.25 1.96 44 -14.66 -3.838 -4.722 0.0001
ComA,SemC - SemA,ComC 10.42 2.07 44 4.71 16.127 5.040 <.0001
ComA,SemC - ComA,ComC 4.08 2.02 44 -1.49 9.661 2.023 0.0984
SemA,ComC - ComA,ComC -6.33 2.46 44 -13.12 0.457 -2.577 0.0402
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 6 estimates
P value adjustment: holm method for 6 tests 
Alcool Cafeina emmean SE df lower.CL upper.CL .group
SemA SemC 7.92 0.957 44 5.42 10.4 a
SemA ComC 10.83 1.770 44 6.23 15.4 a
ComA ComC 17.17 1.710 44 12.71 21.6 b
ComA SemC 21.25 1.070 44 18.45 24.0 b
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 4 estimates
P value adjustment: holm method for 6 tests
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Alcool effect by Cafeina levels
F Test:
P-value adjustment method: holm
Value SE Df Sum of Sq F Pr(>F)
SemC -13.3333 2.0137 1 1066.67 43.842 8.225e-08 ***
ComC -6.3333 2.0137 1 240.67 9.892 0.002973 **
Residuals 44 1070.50
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Post hoc test of the simple main effect of AlcoolNOTE: Results may be misleading due to involvement in interactions Alcool emmean SE df lower.CL upper.CL t.ratio p.value
SemA 9.38 1.00 44 7.35 11.4 9.335 <.0001
ComA 19.21 1.01 44 17.17 21.2 19.030 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.42 44 -12.7 -6.96 -6.906 <.0001
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.00 44 7.04 11.7 a
ComA 19.21 1.01 44 16.87 21.6 b
Results are averaged over the levels of: Cafeina
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Cafeina effect by Alcool levels
F Test:
P-value adjustment method: holm
Value SE Df Sum of Sq F Pr(>F)
SemA -2.9167 2.0137 1 51.04 2.0979 0.15459
ComA 4.0833 2.0137 1 100.04 4.1119 0.09733 .
Residuals 44 1070.50
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Post hoc test of the simple main effect of CafeinaNOTE: Results may be misleading due to involvement in interactions Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 14.6 0.719 44 13.1 16.0 20.275 <.0001
ComC 14.0 1.230 44 11.5 16.5 11.393 <.0001
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 0.583 1.42 44 -2.29 3.45 0.410 0.6840
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
Cafeina emmean SE df lower.CL upper.CL .group
ComC 14.0 1.230 44 11.1 16.9 a
SemC 14.6 0.719 44 12.9 16.3 a
Results are averaged over the levels of: Alcool
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. ANOVA independente bifatorial por Parametric Bootstrap based Generalized Test
Testes:
- teste omnibus:
twowaytests::gpTwoWay- testes post hoc:
twowaytests::paircompTwoWay
- testes post hoc:
ANOVA independente bifatorial por Parametric Bootstrap based Generalized Test é um teste robusto à normalidade e heterocedasticadade. Não há teste post hoc em R.
O código completo está implementado em demo_ANOVA2f_Indep_PBoot_CafAlc.R.
alfa <- 0.05
Dados <- readRDS("CafEntre_AlcEntre.rds")
# ANOVA multifatorial independente por Parametric Bootstrap based Generalized Test:
omni_Caf <- twowaytests::gpTwoWay(NumErros ~ Cafeina*Alcool,
alpha=alfa,
data=Dados)
Generalized p-value by PB method (alpha = 0.05)
--------------------------------------------------
Factor P.value Result
Cafeina 0.0681 Not reject
Alcool 0.0000 Reject
Cafeina:Alcool 0.0213 Reject
--------------------------------------------------
Holm correction (alpha = 0.05) for subgroups of each Cafeina level
--------------------------------------------------------------------------------
Cafeina Alcool(a) Alcool(b) P.value No difference
1 SemC SemA ComA 5.316791e-09 Reject
2 ComC SemA ComA 1.721035e-02 Reject
-------------------------------------------------------------------------------- # Teste de efeito simples de Cafeina
omni_Alc <- twowaytests::gpTwoWay(NumErros ~ Alcool*Cafeina,
alpha=alfa,
data=Dados)
Generalized p-value by PB method (alpha = 0.05)
--------------------------------------------------
Factor P.value Result
Alcool 0.0000 Reject
Cafeina 0.0684 Not reject
Alcool:Cafeina 0.0205 Reject
--------------------------------------------------# Teste de efeito simples de Alcool
twowaytests::paircompTwoWay(omni_Alc,
alpha=alfa,
adjust.method="holm")
Holm correction (alpha = 0.05) for subgroups of each Alcool level
--------------------------------------------------------------------------------
Alcool Cafeina(a) Cafeina(b) P.value No difference
1 SemA SemC ComC 0.16475000 Not reject
2 ComA SemC ComC 0.05778413 Not reject
-------------------------------------------------------------------------------- ANOVA relacionada bifatorial
Testes:
- Omnibus:
lmerTest::lmerecar::Anova(test.statistic="F", ...),effectsize::eta_squared- Efeito principal simples:
phia::testInteractions- Post hoc:
emmeans::emmeans,multcomp::cld,emmeans::emmip
- Post hoc:
- Efeito principal simples:
Usando os mesmos valores do exemplo anterior (para comparação dos resultados), imagine que tenhamos apenas 12 participantes, cada um deles sendo submetido às quatro condições experimentais.
Neste caso, é usado General Linear Mixed Model (GLMM) por
meio da função lmerTest::lmer.
“Vale a pena observar que a MANOVA é ainda um teste válido, mesmo com modestas violações na suposição de normalidade, particularmente quando os tamanhos dos grupos são iguais e existe um número razoável de participantes em cada grupo; por “razoável” entendemos que [para um delineamento] completamente intraparticipantes, [deve haver] pelo menos 22 participantes ao todo.”
Dancey & Reidy, 2019, p. 472
Neste exemplo aplica-se a ANOVA relacionada bifatorial. Os tratamentos são os mesmos, com a diferença importante de que todos os participantes são submetidos às quatro condições:
- com cafeína e com álcool,
- com cafeína e sem álcool,
- sem cafeína e com álcool,
- sem cafeína e sem álcool.
Dizemos, portanto, que ambos os fatores são intra participantes. Nesta modelagem da ANOVA é preciso controlar pelo participante, incluindo a variável de identificação como um efeito aleatório.
|
Adiante veremos ANOVA relacionada (com dois fatores intraparticipantes) e mista (com um fator intra e outro entre participantes):
|
Os dados estão em CafIntra_AlcIntra.rds.
Dados <- data.frame(readxl::read_excel("CafIntra_AlcIntra.xlsx"))
Dados$UE <- factor(Dados$UE)
Dados$Cafeina <- factor(Dados$Cafeina,
levels=c("SemC","ComC"))
Dados$Alcool <- factor(Dados$Alcool,
levels=c("SemA","ComA"))
saveRDS(Dados, "CafIntra_AlcIntra.rds") UE Alcool Cafeina NumErros
1 SemA SemC 4
1 SemA ComC 15
1 ComA SemC 28
1 ComA ComC 10
2 SemA SemC 9
2 SemA ComC 8
2 ComA SemC 22
2 ComA ComC 11
3 SemA SemC 10
3 SemA ComC 17
3 ComA SemC 21
3 ComA ComC 27
4 SemA SemC 8
4 SemA ComC 0
4 ComA SemC 27
4 ComA ComC 15
5 SemA SemC 6
5 SemA ComC 15
5 ComA SemC 21
5 ComA ComC 27
6 SemA SemC 11
6 SemA ComC 11
6 ComA SemC 20
6 ComA ComC 10
7 SemA SemC 2
7 SemA ComC 11
7 ComA SemC 19
7 ComA ComC 21
8 SemA SemC 11
8 SemA ComC 6
8 ComA SemC 16
8 ComA ComC 15
9 SemA SemC 11
9 SemA ComC 0
9 ComA SemC 25
9 ComA ComC 19
10 SemA SemC 10
10 SemA ComC 15
10 ComA SemC 17
10 ComA ComC 21
11 SemA SemC 3
11 SemA ComC 17
11 ComA SemC 19
11 ComA ComC 15
12 SemA SemC 10
12 SemA ComC 15
12 ComA SemC 20
12 ComA ComC 15 Alcool Cafeina NumErros
SemA:24 SemC:24 Min. : 0.00
ComA:24 ComC:24 1st Qu.:10.00
Median :15.00
Mean :14.29
3rd Qu.:19.25
Max. :28.00 Alcool
Cafeina SemA ComA
SemC 12 12
ComC 12 12 Cafeina SemC ComC
UE Alcool
1 SemA 1 1
ComA 1 1
2 SemA 1 1
ComA 1 1
3 SemA 1 1
ComA 1 1
4 SemA 1 1
ComA 1 1
5 SemA 1 1
ComA 1 1
6 SemA 1 1
ComA 1 1
7 SemA 1 1
ComA 1 1
8 SemA 1 1
ComA 1 1
9 SemA 1 1
ComA 1 1
10 SemA 1 1
ComA 1 1
11 SemA 1 1
ComA 1 1
12 SemA 1 1
ComA 1 1Estatística descritiva
Os dados estão no formato long. Observe a coluna das unidades experimentais (UE). Cada observação aparece em uma linha, mas a identificação das unidades experimentais é repetida, de forma que há apenas 12 participantes, cada um deles submetido aos quatro tratamentos.
Os gráficos que foram utilizados para a condição com dois fatores independentes não devem ser usados neste contexto, pois os intervalos de confiança dependem dos desvios-padrão e estes, por sua vez, pressupõem observações independentes. O gráfico de intervalo de confiança apresentado é o mais adequado para medidas repetidas.
ANOVA relacionada bifatorial em R
A linha que define este modelo é
Note o termo (1|UE). Em delineamento com medidas
repetidas, esta é forma de incluir a identificação do participante como
um controle por meio de efeito aleatório.
A ANOVA é computada por
O tamanho de efeito global (\(\eta^2\)) é dado por
Estas funções são adequadas para o participante modelado como efeito aleatório.
ANOVA relacionada bifatorial mostra que há efeito de interação entre Álcool e Cafeína para \(\alpha=5\%\).
Falta determinar o valor p omnibus por meio de estatística de teste F do GLMM e do GLM.
GLM:
- \(a=2\): níveis de Álcool (SemA, ComA)
- \(b=2\): níveis de Cafeína (SemC, ComC)
- \(c=12\): unidades experimentais (UE)
- \(r=1\): observação por célula
- \(n = a \cdot b \cdot c \cdot r = 2 \cdot 2 \cdot 12 \cdot 1 = 48\)
Graus de liberdade
Total: \(df_T = N - 1 = 48 - 1 = 47\)
Fator A (Álcool): \(df_A = a - 1 = 1\)
Fator B (Cafeína): \(df_B = b - 1 = 1\)
Interação A × B: \(df_{AB} = (a-1)(b-1) = 1\)
Bloco (UE): \(df_{UE} = c - 1 = 11\)
Resíduo: \(df_E = df_T - (df_A+df_B+df_{AB}+df_{UE})= 47 - (1+1+1+11)= 33\)
A estatística F omnibus do GLM sob \(H_0\):
\[ \begin{align} F^{\text{omnibus}}_{\text{GLM}}&\sim F_{ab-1,\,(n-1)-ab-c}\\ F^{\text{omnibus}}_{\text{GLM}}&\sim F_{3,\,33}\\ \end{align} \]
GLMM:
A estatística F omnibus do GLM é igual à F omnibus do GLMM sob \(H_0\) se os dados são perfeitamente balanceados:
\[ \begin{align} F^{\text{omnibus}}_{\text{GLMM}}&\sim F_{ab-1,\,(n-1)-ab-c}\\ F^{\text{omnibus}}_{\text{GLMM}}&\sim F_{3,\,33}\\ \end{align} \]
# valores da tabela ANOVA do GLM
SQ_A <- 1160.3
SQ_B <- 4.1
SQ_AB <- 147.0
SQ_E <- 884.1
df_A <- 1
df_B <- 1
df_AB <- 1
df_E <- 33
# graus de liberdade do numerador
df1 <- df_A + df_B + df_AB
df2 <- df_E
# estatística F omnibus
MS_modelo <- (SQ_A + SQ_B + SQ_AB) / df1
MS_E <- SQ_E / df2
F_omnibus <- MS_modelo / MS_E
# valor crítico de 95%
F_crit <- qf(0.95, df1, df2)
# valor-p
p_val <- pf(F_omnibus, df1, df2, lower.tail = FALSE)
cat("F omnibus = ", round(F_omnibus,2), "\n", sep="")F omnibus = 16.32F(3, 33; 0.95) = 2.9 p = 1.12e-06Conclusão: rejeitamos \(H_0\) omnibus para qualquer nível de significância razoável. Há efeito global significativo de Álcool, Cafeína ou de interação.
Análise dos efeitos principais simples
A análise, neste caso, utiliza as mesmas funções do caso
independente. Porém, o objeto veio de lmerTest::lmer.
Por ser um objeto diferente, a estatística de teste mudou, bem como os respectivos valores p. A conclusão ainda é a mesma: há efeito de álcool nos dois níveis de cafeína, mas não há efeito de cafeína nos dois níveis de álcool. O gráfico das médias marginais estimadas é praticamente o mesmo.
Código R completo
Da mesma forma que fizemos no exemplo anterior, o código completo
está em demo_ANOVA2f_CafIntra_AlcIntra.R.
source("summarySEwithin2.R")
Dados <- readRDS("CafIntra_AlcIntra.rds")
print(FSA::Summarize(NumErros ~ Alcool*Cafeina,
data=Dados),
digits=3) Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 12 7.92 3.32 2 5.5 9.5 10.2 11 0.0
2 ComA SemC 12 21.25 3.72 16 19.0 20.5 22.8 28 0.0
3 SemA ComC 12 10.83 6.12 0 7.5 13.0 15.0 17 16.7
4 ComA ComC 12 17.17 5.92 10 14.0 15.0 21.0 27 0.0alfa <- 0.05
g1 <- nlevels(Dados$Alcool)
g2 <- nlevels(Dados$Cafeina)
alfaBonf <- alfa/(g1*g2)
ic <- summarySEwithin2(Dados,
measurevar="NumErros",
withinvars=c("Alcool","Cafeina"),
idvar="UE",
na.rm=TRUE,
conf.interval=1-alfaBonf)
print(ic, digits=3) Alcool Cafeina NumErros n_obs.x NumErrosNormed sd n_obs.y se ci
1 ComA ComC 17.17 12 17.17 5.23 12 1.51 4.50
2 ComA SemC 21.25 12 21.25 5.03 12 1.45 4.33
3 SemA ComC 10.83 12 10.83 5.91 12 1.71 5.09
4 SemA SemC 7.92 12 7.92 4.41 12 1.27 3.79grf <- ggplot2::ggplot(ic,
ggplot2::aes(x=Alcool,
y=NumErros,
colour=Cafeina)) +
ggplot2::geom_errorbar(position=ggplot2::position_dodge(.9),
width=.1,
ggplot2::aes(ymin=NumErros-ci,
ymax=NumErros+ci)) +
ggplot2::geom_point(shape=21,
size=3,
fill="white",
position=ggplot2::position_dodge(.9)) +
ggplot2::ylab("NumErros") +
ggplot2::ggtitle("Álcool & Cafeína: Número de Erros\nWithin-subject CI95% Bonferroni") +
ggplot2::theme_bw()
print(grf)
Two-factor within-subjects ANOVA
GLMMboundary (singular) fit: see help('isSingular')
ANOVA tableAnalysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)
Response: NumErros
F Df Df.res Pr(>F)
Alcool 47.692 1 33 6.89e-08 ***
Cafeina 0.168 1 33 0.6847
Alcool:Cafeina 6.042 1 33 0.0194 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: NumErros ~ Alcool * Cafeina + (1 | UE)
Data: Dados
REML criterion at convergence: 275.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.19632 -0.48572 0.02534 0.76449 1.99358
Random effects:
Groups Name Variance Std.Dev.
UE (Intercept) 0.00 0.000
Residual 24.33 4.932
Number of obs: 48, groups: UE, 12
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 7.917 1.424 44.000 5.560 1.49e-06 ***
AlcoolComA 13.333 2.014 44.000 6.621 4.11e-08 ***
CafeinaComC 2.917 2.014 44.000 1.448 0.155
AlcoolComA:CafeinaComC -7.000 2.848 44.000 -2.458 0.018 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
GLM
ANOVA tableAnova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
UE 186.4 11 0.633 0.7880
Alcool 1160.3 1 43.312 1.76e-07 ***
Cafeina 4.1 1 0.152 0.6987
Alcool:Cafeina 147.0 1 5.487 0.0253 *
Residuals 884.1 33
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = NumErros ~ UE + Alcool * Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-10.292 -3.021 0.375 3.021 7.542
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.875 2.893 2.722 0.0103 *
UE2 -1.750 3.660 -0.478 0.6357
UE3 4.500 3.660 1.230 0.2276
UE4 -1.750 3.660 -0.478 0.6357
UE5 3.000 3.660 0.820 0.4183
UE6 -1.250 3.660 -0.342 0.7349
UE7 -1.000 3.660 -0.273 0.7864
UE8 -2.250 3.660 -0.615 0.5429
UE9 -0.500 3.660 -0.137 0.8922
UE10 1.500 3.660 0.410 0.6846
UE11 -0.750 3.660 -0.205 0.8389
UE12 0.750 3.660 0.205 0.8389
AlcoolComA 13.333 2.113 6.310 3.89e-07 ***
CafeinaComC 2.917 2.113 1.380 0.1768
AlcoolComA:CafeinaComC -7.000 2.988 -2.342 0.0253 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.176 on 33 degrees of freedom
Multiple R-squared: 0.6288, Adjusted R-squared: 0.4714
F-statistic: 3.994 on 14 and 33 DF, p-value: 0.0005265
Effect size analysisRegistered S3 methods overwritten by 'MuMIn':
method from
nobs.multinom broom
nobs.fitdistr broom
Tamanho de efeito: eta^2 omnibus = 0.5342es <- effectsize::interpret_eta_squared(eta2g,
rules="cohen1992")
names(es) <- c("Tamanho de efeito omnibus: estimativa pontual")
print(es)Tamanho de efeito omnibus: estimativa pontual
"large"
(Rules: cohen1992)eta2 <- effectsize::eta_squared(anv,
partial=FALSE,
generalized=FALSE,
ci=1-alfa/3,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 | 98.3333% CI | interpret
-------------------------------------------------------
UE | 0.0783 | [0.0000, 0.0135] | medium
Alcool | 0.4871 | [0.1782, 0.6855] | large
Cafeina | 0.0017 | [0.0000, 0.0035] | very small
Alcool:Cafeina | 0.0617 | [0.0000, 0.3150] | medium
Gráfico de perfil de médiasg1 <- emmeans::emmip(object=modelo,
formula=~Cafeina ~Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
dodge=.3,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados)) +
ggplot2::theme_bw()
g2 <- emmeans::emmip(object=modelo,
formula=~Alcool ~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
dodge=.3,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados)) +
ggplot2::theme_bw()
g3 <- emmeans::emmip(modelo,
~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactionsg4 <- emmeans::emmip(modelo,
~Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactions
Post hoc test se efeito de interação Alcool:Cafeina não é significante
AlcoolEMM.A <- emmeans::emmeans(modelo,
specs=pairwise~"Alcool",
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactions Alcool emmean SE df lower.CL upper.CL
SemA 9.38 1.01 44 7.35 11.4
ComA 19.21 1.01 44 17.18 21.2
Results are averaged over the levels of: Cafeina
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.42 44 -12.7 -6.96 -6.906 <.0001
Results are averaged over the levels of: Cafeina
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 
print(multcomp::cld(object=EMM.A$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.01 44 7.04 11.7 a
ComA 19.21 1.01 44 16.87 21.5 b
Results are averaged over the levels of: Cafeina
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
CafeinaEMM.C <- emmeans::emmeans(modelo,
specs=pairwise~"Cafeina",
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactions Cafeina emmean SE df lower.CL upper.CL
SemC 14.6 1.01 44 12.6 16.6
ComC 14.0 1.01 44 12.0 16.0
Results are averaged over the levels of: Alcool
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 0.583 1.42 44 -2.29 3.45 0.410 0.6840
Results are averaged over the levels of: Alcool
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 
print(multcomp::cld(object=EMM.C$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Cafeina emmean SE df lower.CL upper.CL .group
ComC 14.0 1.01 44 11.7 16.3 a
SemC 14.6 1.01 44 12.2 16.9 a
Results are averaged over the levels of: Alcool
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Post hoc teste se efeito de interação Alcool:Cafeina é significante
Post hoc test: interaction Cafeina:Alcool
EMM.AC <- emmeans::emmeans(modelo,
specs=pairwise~"Alcool:Cafeina",
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))
print(summary(EMM.AC$emmeans)) Alcool Cafeina emmean SE df lower.CL upper.CL
SemA SemC 7.92 1.42 44 5.05 10.8
ComA SemC 21.25 1.42 44 18.38 24.1
SemA ComC 10.83 1.42 44 7.96 13.7
ComA ComC 17.17 1.42 44 14.30 20.0
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA,SemC - ComA,SemC -13.33 2.01 44 -18.90 -7.77 -6.621 <.0001
SemA,SemC - SemA,ComC -2.92 2.01 44 -8.48 2.65 -1.448 0.1546
SemA,SemC - ComA,ComC -9.25 2.01 44 -14.81 -3.69 -4.594 0.0001
ComA,SemC - SemA,ComC 10.42 2.01 44 4.85 15.98 5.173 <.0001
ComA,SemC - ComA,ComC 4.08 2.01 44 -1.48 9.65 2.028 0.0973
SemA,ComC - ComA,ComC -6.33 2.01 44 -11.90 -0.77 -3.145 0.0089
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 6 estimates
P value adjustment: holm method for 6 tests 
print(multcomp::cld(object=EMM.AC$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool Cafeina emmean SE df lower.CL upper.CL .group
SemA SemC 7.92 1.42 44 4.21 11.6 a
SemA ComC 10.83 1.42 44 7.12 14.5 a
ComA ComC 17.17 1.42 44 13.46 20.9 b
ComA SemC 21.25 1.42 44 17.54 25.0 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 4 estimates
P value adjustment: holm method for 6 tests
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Alcool effect by Cafeina levelsChisq Test:
P-value adjustment method: holm
Value SE Df Chisq Pr(>Chisq)
SemC -13.3333 2.0137 1 43.842 7.118e-11 ***
ComC -6.3333 2.0137 1 9.892 0.00166 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(emmeans::emmip(modelo,
~Alcool|Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))) 
Post hoc test of the simple main effect of Alcoolmodel_means <- emmeans::emmeans(object=modelo,
specs=pairwise~Alcool|Cafeina:Alcool,
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))
print(summary(model_means$emmeans, infer=TRUE))Cafeina = SemC:
Alcool emmean SE df lower.CL upper.CL t.ratio p.value
SemA 7.92 1.42 44 5.05 10.8 5.560 <.0001
ComA 21.25 1.42 44 18.38 24.1 14.924 <.0001
Cafeina = ComC:
Alcool emmean SE df lower.CL upper.CL t.ratio p.value
SemA 10.83 1.42 44 7.96 13.7 7.608 <.0001
ComA 17.17 1.42 44 14.30 20.0 12.056 <.0001
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 Cafeina = SemC:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -13.33 2.01 44 -17.4 -9.28 -6.621 <.0001
Cafeina = ComC:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -6.33 2.01 44 -10.4 -2.28 -3.145 0.0030
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 
print(multcomp::cld(object=model_means$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters))Cafeina = SemC:
Alcool emmean SE df lower.CL upper.CL .group
SemA 7.92 1.42 44 4.61 11.2 a
ComA 21.25 1.42 44 17.95 24.6 b
Cafeina = ComC:
Alcool emmean SE df lower.CL upper.CL .group
SemA 10.83 1.42 44 7.53 14.1 a
ComA 17.17 1.42 44 13.86 20.5 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Cafeina effect by Alcool levelsChisq Test:
P-value adjustment method: holm
Value SE Df Chisq Pr(>Chisq)
SemA -2.9167 2.0137 1 2.0979 0.14750
ComA 4.0833 2.0137 1 4.1119 0.08516 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(emmeans::emmip(modelo,
~Cafeina|Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))) 
Post hoc test of the simple main effect of Cafeinamodel_means <- emmeans::emmeans(object=modelo,
specs=pairwise~Cafeina|Cafeina:Alcool,
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))
print(summary(model_means$emmeans, infer=TRUE))Alcool = SemA:
Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 7.92 1.42 44 5.05 10.8 5.560 <.0001
ComC 10.83 1.42 44 7.96 13.7 7.608 <.0001
Alcool = ComA:
Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 21.25 1.42 44 18.38 24.1 14.924 <.0001
ComC 17.17 1.42 44 14.30 20.0 12.056 <.0001
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 Alcool = SemA:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC -2.92 2.01 44 -6.975 1.14 -1.448 0.1546
Alcool = ComA:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 4.08 2.01 44 0.025 8.14 2.028 0.0487
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 
print(multcomp::cld(object=model_means$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters))Alcool = SemA:
Cafeina emmean SE df lower.CL upper.CL .group
SemC 7.92 1.42 44 4.61 11.2 a
ComC 10.83 1.42 44 7.53 14.1 a
Alcool = ComA:
Cafeina emmean SE df lower.CL upper.CL .group
ComC 17.17 1.42 44 13.86 20.5 a
SemC 21.25 1.42 44 17.95 24.6 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. ANOVA relacionada bifatorial desbalanceada
UE Alcool Cafeina NumErros
2 1 SemA ComC 15
3 1 ComA SemC 28
4 1 ComA ComC 10
5 2 SemA SemC 9
6 2 SemA ComC 8
7 2 ComA SemC 22
8 2 ComA ComC 11
9 3 SemA SemC 10
10 3 SemA ComC 17
11 3 ComA SemC 21
12 3 ComA ComC 27
13 4 SemA SemC 8
14 4 SemA ComC 0
15 4 ComA SemC 27
16 4 ComA ComC 15
17 5 SemA SemC 6
18 5 SemA ComC 15
19 5 ComA SemC 21
20 5 ComA ComC 27
21 6 SemA SemC 11
22 6 SemA ComC 11
23 6 ComA SemC 20
24 6 ComA ComC 10
25 7 SemA SemC 2
26 7 SemA ComC 11
27 7 ComA SemC 19
28 7 ComA ComC 21
29 8 SemA SemC 11
30 8 SemA ComC 6
31 8 ComA SemC 16
32 8 ComA ComC 15
33 9 SemA SemC 11
34 9 SemA ComC 0
35 9 ComA SemC 25
36 9 ComA ComC 19
37 10 SemA SemC 10
38 10 SemA ComC 15
39 10 ComA SemC 17
40 10 ComA ComC 21
41 11 SemA SemC 3
42 11 SemA ComC 17
43 11 ComA SemC 19
44 11 ComA ComC 15
45 12 SemA SemC 10
46 12 SemA ComC 15
47 12 ComA SemC 20
48 12 ComA ComC 15 Alcool
Cafeina SemA ComA
SemC 11 12
ComC 12 12 Cafeina SemC ComC
UE Alcool
1 SemA 0 1
ComA 1 1
2 SemA 1 1
ComA 1 1
3 SemA 1 1
ComA 1 1
4 SemA 1 1
ComA 1 1
5 SemA 1 1
ComA 1 1
6 SemA 1 1
ComA 1 1
7 SemA 1 1
ComA 1 1
8 SemA 1 1
ComA 1 1
9 SemA 1 1
ComA 1 1
10 SemA 1 1
ComA 1 1
11 SemA 1 1
ComA 1 1
12 SemA 1 1
ComA 1 1 Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 11 8.3 3.2 2 7.0 10 10 11 0
2 ComA SemC 12 21.2 3.7 16 19.0 20 23 28 0
3 SemA ComC 12 10.8 6.1 0 7.5 13 15 17 17
4 ComA ComC 12 17.2 5.9 10 14.0 15 21 27 0
Two-factor within-subjects ANOVA
GLMMboundary (singular) fit: see help('isSingular')Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)
Response: NumErros
F Df Df.res Pr(>F)
Alcool 43.826 1 32.40 1.72e-07 ***
Cafeina 0.339 1 32.40 0.5645
Alcool:Cafeina 5.273 1 32.42 0.0283 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
GLMAnova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
UE 191.5 11 0.646 0.7760
Alcool 1052.6 1 39.065 5.27e-07 ***
Cafeina 9.7 1 0.362 0.5519
Alcool:Cafeina 124.2 1 4.608 0.0395 *
Residuals 862.2 32
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analysis of Variance Table
Response: NumErros
Df Sum Sq Mean Sq F value Pr(>F)
UE 11 218.3 19.8 0.737 0.6963
Alcool 1 1059.3 1059.3 39.312 4.98e-07 ***
Cafeina 1 9.7 9.7 0.362 0.5519
Alcool:Cafeina 1 124.2 124.2 4.608 0.0395 *
Residuals 32 862.2 26.9
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ANOVA bifatorial mista
Testes:
- Omnibus:
lmerTest::lmerecar::Anova(test.statistic="F", ...),effectsize::eta_squared- Efeito principal simples:
phia::testInteractions- Post hoc:
emmeans::emmeans,multcomp::cld,emmeans::emmip
- Post hoc:
- Efeito principal simples:
Na ANOVA bifatorial mista, a identificação de cada um dos participantes (UE) aparece duas vezes, de tal forma que há apenas 24 participantes no total, alocados nos quatro tratamentos da seguinte maneira:
- 12 participantes sem Álcool, submetidos às duas condições:
- com Cafeína,
- sem Cafeína
- 12 participantes com Álcool, submetidos às duas condições:
- com Cafeína,
- sem Cafeína.
Por isso, o Álcool é o fator entre participantes, mas a Cafeína é intra participantes.
Sinonímia:
- Two-factor mixed-design ANOVA
- Two-way mixed-design ANOVA
- Two-factor repeated measures ANOVA with split-plot structure
- Two-way ANOVA with split-plot design
- Two-way split-plot ANOVA
Neste caso, é usado General Linear Mixed Model (GLMM) por
meio da função lmerTest::lmer.
Os dados estão em CafIntra_AlcEntre.rds.
Dados <- data.frame(readxl::read_excel("CafIntra_AlcEntre.xlsx"))
Dados$UE <- factor(Dados$UE)
Dados$Cafeina <- factor(Dados$Cafeina,
levels=c("SemC","ComC"))
Dados$Alcool <- factor(Dados$Alcool,
levels=c("SemA","ComA"))
saveRDS(Dados, "CafIntra_AlcEntre.rds") UE NumErros Alcool Cafeina
1 1 4 SemA SemC
2 1 15 SemA ComC
3 2 9 SemA SemC
4 2 8 SemA ComC
5 3 10 SemA SemC
6 3 17 SemA ComC
7 4 8 SemA SemC
8 4 0 SemA ComC
9 5 6 SemA SemC
10 5 15 SemA ComC
11 6 11 SemA SemC
12 6 11 SemA ComC
13 7 2 SemA SemC
14 7 11 SemA ComC
15 8 11 SemA SemC
16 8 6 SemA ComC
17 9 11 SemA SemC
18 9 0 SemA ComC
19 10 10 SemA SemC
20 10 15 SemA ComC
21 11 3 SemA SemC
22 11 17 SemA ComC
23 12 10 SemA SemC
24 12 15 SemA ComC
25 13 28 ComA SemC
26 13 10 ComA ComC
27 14 22 ComA SemC
28 14 11 ComA ComC
29 15 21 ComA SemC
30 15 27 ComA ComC
31 16 27 ComA SemC
32 16 15 ComA ComC
33 17 21 ComA SemC
34 17 27 ComA ComC
35 18 20 ComA SemC
36 18 10 ComA ComC
37 19 19 ComA SemC
38 19 21 ComA ComC
39 20 16 ComA SemC
40 20 15 ComA ComC
41 21 25 ComA SemC
42 21 19 ComA ComC
43 22 17 ComA SemC
44 22 21 ComA ComC
45 23 19 ComA SemC
46 23 15 ComA ComC
47 24 20 ComA SemC
48 24 15 ComA ComC NumErros Alcool Cafeina
Min. : 0.00 SemA:24 SemC:24
1st Qu.:10.00 ComA:24 ComC:24
Median :15.00
Mean :14.29
3rd Qu.:19.25
Max. :28.00 Alcool
Cafeina SemA ComA
SemC 12 12
ComC 12 12 Cafeina SemC ComC
UE Alcool
1 SemA 1 1
ComA 0 0
2 SemA 1 1
ComA 0 0
3 SemA 1 1
ComA 0 0
4 SemA 1 1
ComA 0 0
5 SemA 1 1
ComA 0 0
6 SemA 1 1
ComA 0 0
7 SemA 1 1
ComA 0 0
8 SemA 1 1
ComA 0 0
9 SemA 1 1
ComA 0 0
10 SemA 1 1
ComA 0 0
11 SemA 1 1
ComA 0 0
12 SemA 1 1
ComA 0 0
13 SemA 0 0
ComA 1 1
14 SemA 0 0
ComA 1 1
15 SemA 0 0
ComA 1 1
16 SemA 0 0
ComA 1 1
17 SemA 0 0
ComA 1 1
18 SemA 0 0
ComA 1 1
19 SemA 0 0
ComA 1 1
20 SemA 0 0
ComA 1 1
21 SemA 0 0
ComA 1 1
22 SemA 0 0
ComA 1 1
23 SemA 0 0
ComA 1 1
24 SemA 0 0
ComA 1 1Implementado em demo_ANOVA2f_CafIntra_AlcEntre.R. Este
é o mesmo código R usado para a ANOVA relacionada bifatorial. Somente a
estrutura da planilha (a coluna com a variável UE é a diferença) basta
para que o comportamento da análise seja adaptado por
lmerTest::lmer:
source("summarySEwithin2.R")
Dados <- readRDS("CafIntra_AlcEntre.rds")
print(FSA::Summarize(NumErros ~ Alcool*Cafeina,
data=Dados),
digits=3) Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 12 7.92 3.32 2 5.5 9.5 10.2 11 0.0
2 ComA SemC 12 21.25 3.72 16 19.0 20.5 22.8 28 0.0
3 SemA ComC 12 10.83 6.12 0 7.5 13.0 15.0 17 16.7
4 ComA ComC 12 17.17 5.92 10 14.0 15.0 21.0 27 0.0alfa <- 0.05
g1 <- nlevels(Dados$Alcool)
g2 <- nlevels(Dados$Cafeina)
alfaBonf <- alfa/(g1*g2)
ic <- summarySEwithin2(Dados,
measurevar="NumErros",
withinvars=c("Alcool","Cafeina"),
idvar="UE",
na.rm=TRUE,
conf.interval=1-alfaBonf)
print(ic, digits=3) Alcool Cafeina NumErros n_obs.x NumErrosNormed sd n_obs.y se ci
1 ComA ComC 17.17 12 12.2 4.47 12 1.29 3.85
2 ComA SemC 21.25 12 16.3 4.47 12 1.29 3.85
3 SemA ComC 10.83 12 15.7 4.56 12 1.32 3.93
4 SemA SemC 7.92 12 12.8 4.56 12 1.32 3.93grf <- ggplot2::ggplot(ic,
ggplot2::aes(x=Alcool,
y=NumErros,
colour=Cafeina)) +
ggplot2::geom_errorbar(position=ggplot2::position_dodge(.9),
width=.1,
ggplot2::aes(ymin=NumErros-ci,
ymax=NumErros+ci)) +
ggplot2::geom_point(shape=21,
size=3,
fill="white",
position=ggplot2::position_dodge(.9)) +
ggplot2::ylab("NumErros") +
ggplot2::ggtitle("Álcool & Cafeína: Número de Erros\nMixed CI95% Bonferroni") +
ggplot2::theme_bw()
print(grf)
Two-factor mixed ANOVA
GLMMboundary (singular) fit: see help('isSingular')
ANOVA tableAnalysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)
Response: NumErros
F Df Df.res Pr(>F)
Alcool 47.692 1 22 6.19e-07 ***
Cafeina 0.168 1 22 0.6860
Alcool:Cafeina 6.042 1 22 0.0223 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: NumErros ~ Alcool * Cafeina + (1 | UE)
Data: Dados
REML criterion at convergence: 275.2
Scaled residuals:
Min 1Q Median 3Q Max
-2.19632 -0.48572 0.02534 0.76449 1.99358
Random effects:
Groups Name Variance Std.Dev.
UE (Intercept) 0.00 0.000
Residual 24.33 4.932
Number of obs: 48, groups: UE, 24
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 7.917 1.424 44.000 5.560 1.49e-06 ***
AlcoolComA 13.333 2.014 44.000 6.621 4.11e-08 ***
CafeinaComC 2.917 2.014 44.000 1.448 0.155
AlcoolComA:CafeinaComC -7.000 2.848 44.000 -2.458 0.018 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
GLM
ANOVA tableAnalysis of Variance Table
Response: NumErros
Df Sum Sq Mean Sq F value Pr(>F)
UE 23 1557.9 67.74 2.215 0.0334 *
Cafeina 1 4.1 4.08 0.133 0.7183
Alcool:Cafeina 1 147.0 147.00 4.806 0.0392 *
Residuals 22 672.9 30.59
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Note: model has aliased coefficients
sums of squares computed by model comparisonAnova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
UE 397.6 22 0.591 0.8875
Alcool 0
Cafeina 4.1 1 0.133 0.7183
Alcool:Cafeina 147.0 1 4.806 0.0392 *
Residuals 672.9 22
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = NumErros ~ UE + Alcool * Cafeina, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-6.958 -3.042 0.000 3.042 6.958
Coefficients: (1 not defined because of singularities)
Estimate Std. Error t value Pr(>|t|)
(Intercept) 8.042 4.070 1.976 0.06087 .
UE2 -1.000 5.531 -0.181 0.85817
UE3 4.000 5.531 0.723 0.47714
UE4 -5.500 5.531 -0.994 0.33081
UE5 1.000 5.531 0.181 0.85817
UE6 1.500 5.531 0.271 0.78875
UE7 -3.000 5.531 -0.542 0.59297
UE8 -1.000 5.531 -0.181 0.85817
UE9 -4.000 5.531 -0.723 0.47714
UE10 3.000 5.531 0.542 0.59297
UE11 0.500 5.531 0.090 0.92878
UE12 3.000 5.531 0.542 0.59297
UE13 13.000 5.756 2.258 0.03418 *
UE14 10.500 5.756 1.824 0.08176 .
UE15 18.000 5.756 3.127 0.00491 **
UE16 15.000 5.756 2.606 0.01614 *
UE17 18.000 5.756 3.127 0.00491 **
UE18 9.000 5.756 1.563 0.13221
UE19 14.000 5.756 2.432 0.02361 *
UE20 9.500 5.756 1.650 0.11308
UE21 16.000 5.756 2.780 0.01093 *
UE22 13.000 5.756 2.258 0.03418 *
UE23 11.000 5.756 1.911 0.06913 .
UE24 11.500 5.756 1.998 0.05825 .
AlcoolComA NA NA NA NA
CafeinaComC 2.917 2.258 1.292 0.20984
AlcoolComA:CafeinaComC -7.000 3.193 -2.192 0.03923 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.531 on 22 degrees of freedom
Multiple R-squared: 0.7175, Adjusted R-squared: 0.3965
F-statistic: 2.235 on 25 and 22 DF, p-value: 0.03035
Effect size analysiseta2g <- as.numeric(MuMIn::r.squaredGLMM(modelo)[1])
cat("\nTamanho de efeito: eta^2 omnibus =", round(eta2g,4))
Tamanho de efeito: eta^2 omnibus = 0.5342es <- effectsize::interpret_eta_squared(eta2g,
rules="cohen1992")
names(es) <- c("Tamanho de efeito omnibus: estimativa pontual")
print(es)Tamanho de efeito omnibus: estimativa pontual
"large"
(Rules: cohen1992)eta2 <- effectsize::eta_squared(anv,
partial=FALSE,
generalized=FALSE,
ci=1-alfa/3,
alternative="two.sided",
verbose=TRUE)Warning: Some CIs could not be estimated due to non-finite F, df, or df_error
values.# Effect Size for ANOVA (Type II)
Parameter | Eta2 | 98.3333% CI | interpret
-----------------------------------------------
UE | | |
Alcool | | |
Cafeina | | |
Alcool:Cafeina | | |
Gráfico de perfil de médiasg1 <- emmeans::emmip(object=modelo,
formula=~Cafeina ~Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
dodge=.3,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados)) +
ggplot2::theme_bw()
g2 <- emmeans::emmip(object=modelo,
formula=~Alcool ~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
dodge=.3,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados)) +
ggplot2::theme_bw()
g3 <- emmeans::emmip(modelo,
~Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactionsg4 <- emmeans::emmip(modelo,
~Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactions
Post hoc test se efeito de interação Alcool:Cafeina não é significante
AlcoolEMM.A <- emmeans::emmeans(modelo,
specs=pairwise~"Alcool",
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactions Alcool emmean SE df lower.CL upper.CL
SemA 9.38 1.01 44 7.35 11.4
ComA 19.21 1.01 44 17.18 21.2
Results are averaged over the levels of: Cafeina
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -9.83 1.42 44 -12.7 -6.96 -6.906 <.0001
Results are averaged over the levels of: Cafeina
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
print(multcomp::cld(object=EMM.A$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool emmean SE df lower.CL upper.CL .group
SemA 9.38 1.01 44 7.04 11.7 a
ComA 19.21 1.01 44 16.87 21.5 b
Results are averaged over the levels of: Cafeina
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
CafeinaEMM.C <- emmeans::emmeans(modelo,
specs=pairwise~"Cafeina",
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))NOTE: Results may be misleading due to involvement in interactions Cafeina emmean SE df lower.CL upper.CL
SemC 14.6 1.01 44 12.6 16.6
ComC 14.0 1.01 44 12.0 16.0
Results are averaged over the levels of: Alcool
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 0.583 1.42 44 -2.29 3.45 0.410 0.6840
Results are averaged over the levels of: Alcool
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
print(multcomp::cld(object=EMM.C$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Cafeina emmean SE df lower.CL upper.CL .group
ComC 14.0 1.01 44 11.7 16.3 a
SemC 14.6 1.01 44 12.2 16.9 a
Results are averaged over the levels of: Alcool
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Post hoc teste se efeito de interação Alcool:Cafeina é significante
Post hoc test: interaction Cafeina:Alcool
EMM.AC <- emmeans::emmeans(modelo,
specs=pairwise~"Alcool:Cafeina",
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))
print(summary(EMM.AC$emmeans)) Alcool Cafeina emmean SE df lower.CL upper.CL
SemA SemC 7.92 1.42 44 5.05 10.8
ComA SemC 21.25 1.42 44 18.38 24.1
SemA ComC 10.83 1.42 44 7.96 13.7
ComA ComC 17.17 1.42 44 14.30 20.0
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA,SemC - ComA,SemC -13.33 2.01 44 -18.90 -7.77 -6.621 <.0001
SemA,SemC - SemA,ComC -2.92 2.01 44 -8.48 2.65 -1.448 0.1546
SemA,SemC - ComA,ComC -9.25 2.01 44 -14.81 -3.69 -4.594 0.0001
ComA,SemC - SemA,ComC 10.42 2.01 44 4.85 15.98 5.173 <.0001
ComA,SemC - ComA,ComC 4.08 2.01 44 -1.48 9.65 2.028 0.0973
SemA,ComC - ComA,ComC -6.33 2.01 44 -11.90 -0.77 -3.145 0.0089
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 6 estimates
P value adjustment: holm method for 6 tests
print(multcomp::cld(object=EMM.AC$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters)) Alcool Cafeina emmean SE df lower.CL upper.CL .group
SemA SemC 7.92 1.42 44 4.21 11.6 a
SemA ComC 10.83 1.42 44 7.12 14.5 a
ComA ComC 17.17 1.42 44 13.46 20.9 b
ComA SemC 21.25 1.42 44 17.54 25.0 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 4 estimates
P value adjustment: holm method for 6 tests
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Alcool effect by Cafeina levelsChisq Test:
P-value adjustment method: holm
Value SE Df Chisq Pr(>Chisq)
SemC -13.3333 2.0137 1 43.842 7.118e-11 ***
ComC -6.3333 2.0137 1 9.892 0.00166 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(emmeans::emmip(modelo,
~Alcool|Cafeina,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados)))
Post hoc test of the simple main effect of Alcoolmodel_means <- emmeans::emmeans(object=modelo,
specs=pairwise~Alcool|Cafeina:Alcool,
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))
print(summary(model_means$emmeans, infer=TRUE))Cafeina = SemC:
Alcool emmean SE df lower.CL upper.CL t.ratio p.value
SemA 7.92 1.42 44 5.05 10.8 5.560 <.0001
ComA 21.25 1.42 44 18.38 24.1 14.924 <.0001
Cafeina = ComC:
Alcool emmean SE df lower.CL upper.CL t.ratio p.value
SemA 10.83 1.42 44 7.96 13.7 7.608 <.0001
ComA 17.17 1.42 44 14.30 20.0 12.056 <.0001
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 Cafeina = SemC:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -13.33 2.01 44 -17.4 -9.28 -6.621 <.0001
Cafeina = ComC:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemA - ComA -6.33 2.01 44 -10.4 -2.28 -3.145 0.0030
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
print(multcomp::cld(object=model_means$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters))Cafeina = SemC:
Alcool emmean SE df lower.CL upper.CL .group
SemA 7.92 1.42 44 4.61 11.2 a
ComA 21.25 1.42 44 17.95 24.6 b
Cafeina = ComC:
Alcool emmean SE df lower.CL upper.CL .group
SemA 10.83 1.42 44 7.53 14.1 a
ComA 17.17 1.42 44 13.86 20.5 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Simple main effect test: Cafeina effect by Alcool levelsChisq Test:
P-value adjustment method: holm
Value SE Df Chisq Pr(>Chisq)
SemA -2.9167 2.0137 1 2.0979 0.14750
ComA 4.0833 2.0137 1 4.1119 0.08516 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(emmeans::emmip(modelo,
~Cafeina|Alcool,
level=1-alfa,
adjust="holm",
CIs=TRUE,
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados)))
Post hoc test of the simple main effect of Cafeinamodel_means <- emmeans::emmeans(object=modelo,
specs=pairwise~Cafeina|Cafeina:Alcool,
level=1-alfa,
adjust="holm",
lmer.df="satterthwaite",
lmerTest.limit=nrow(Dados))
print(summary(model_means$emmeans, infer=TRUE))Alcool = SemA:
Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 7.92 1.42 44 5.05 10.8 5.560 <.0001
ComC 10.83 1.42 44 7.96 13.7 7.608 <.0001
Alcool = ComA:
Cafeina emmean SE df lower.CL upper.CL t.ratio p.value
SemC 21.25 1.42 44 18.38 24.1 14.924 <.0001
ComC 17.17 1.42 44 14.30 20.0 12.056 <.0001
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95 Alcool = SemA:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC -2.92 2.01 44 -6.975 1.14 -1.448 0.1546
Alcool = ComA:
contrast estimate SE df lower.CL upper.CL t.ratio p.value
SemC - ComC 4.08 2.01 44 0.025 8.14 2.028 0.0487
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
print(multcomp::cld(object=model_means$emmeans,
level=1-alfa,
alpha=alfa,
adjust="holm",
Letters=letters))Alcool = SemA:
Cafeina emmean SE df lower.CL upper.CL .group
SemC 7.92 1.42 44 4.61 11.2 a
ComC 10.83 1.42 44 7.53 14.1 a
Alcool = ComA:
Cafeina emmean SE df lower.CL upper.CL .group
ComC 17.17 1.42 44 13.86 20.5 a
SemC 21.25 1.42 44 17.95 24.6 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. ANOVA mista bifatorial desbalanceada
UE NumErros Alcool Cafeina
2 1 15 SemA ComC
3 2 9 SemA SemC
4 2 8 SemA ComC
5 3 10 SemA SemC
6 3 17 SemA ComC
7 4 8 SemA SemC
8 4 0 SemA ComC
9 5 6 SemA SemC
10 5 15 SemA ComC
11 6 11 SemA SemC
12 6 11 SemA ComC
13 7 2 SemA SemC
14 7 11 SemA ComC
15 8 11 SemA SemC
16 8 6 SemA ComC
17 9 11 SemA SemC
18 9 0 SemA ComC
19 10 10 SemA SemC
20 10 15 SemA ComC
21 11 3 SemA SemC
22 11 17 SemA ComC
23 12 10 SemA SemC
24 12 15 SemA ComC
25 13 28 ComA SemC
26 13 10 ComA ComC
27 14 22 ComA SemC
28 14 11 ComA ComC
29 15 21 ComA SemC
30 15 27 ComA ComC
31 16 27 ComA SemC
32 16 15 ComA ComC
33 17 21 ComA SemC
34 17 27 ComA ComC
35 18 20 ComA SemC
36 18 10 ComA ComC
37 19 19 ComA SemC
38 19 21 ComA ComC
39 20 16 ComA SemC
40 20 15 ComA ComC
41 21 25 ComA SemC
42 21 19 ComA ComC
43 22 17 ComA SemC
44 22 21 ComA ComC
45 23 19 ComA SemC
46 23 15 ComA ComC
47 24 20 ComA SemC
48 24 15 ComA ComC Alcool
Cafeina SemA ComA
SemC 11 12
ComC 12 12 Alcool Cafeina n mean sd min Q1 median Q3 max percZero
1 SemA SemC 11 8.3 3.2 2 7.0 10 10 11 0
2 ComA SemC 12 21.2 3.7 16 19.0 20 23 28 0
3 SemA ComC 12 10.8 6.1 0 7.5 13 15 17 17
4 ComA ComC 12 17.2 5.9 10 14.0 15 21 27 0
Two-factor within-subjects ANOVA
GLMMboundary (singular) fit: see help('isSingular')Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)
Response: NumErros
F Df Df.res Pr(>F)
Alcool 43.826 1 21.73 1.25e-06 ***
Cafeina 0.339 1 21.71 0.5664
Alcool:Cafeina 5.273 1 21.73 0.0317 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
GLMNote: model has aliased coefficients
sums of squares computed by model comparisonAnova Table (Type II tests)
Response: NumErros
Sum Sq Df F value Pr(>F)
UE 416.5 22 0.624 0.8603
Alcool 0
Cafeina 13.6 1 0.448 0.5107
Alcool:Cafeina 112.6 1 3.712 0.0677 .
Residuals 637.3 21
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analysis of Variance Table
Response: NumErros
Df Sum Sq Mean Sq F value Pr(>F)
UE 23 1510.2 65.66 2.164 0.0400 *
Cafeina 1 13.6 13.59 0.448 0.5107
Alcool:Cafeina 1 112.6 112.64 3.712 0.0677 .
Residuals 21 637.3 30.35
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1|
Na ANOVA independente é requerido pelo menos duas observações em cada um dos tratamentos para testar efeito de interação. Para ANOVA relacionada ou mista é requerido pelo menos uma observação de cada tratamento para testar efeito de interação. Quando um modelo completo como não pode ser testado ou não se deseja testar, uma alternativa é usar um modelo de efeitos principais: |
Exemplo: Estudo dos Efeitos da Vitamina C no Desenvolvimento Dentário de Porquinhos-da-Índia: Comparação entre Fontes Biológicas e Químicas
porquinho-da-índia
Uma preocupação durante a Segunda Guerra Mundial era o fornecimento de vitamina C aos soldados e, nesse amplo contexto, os efeitos do ácido ascórbico e do suco de laranja foram estudados em animais. Um desses estudos foi realizado por Crampton (1947).
Sessenta porquinhos-da-índia (com 28 dias de idade) receberam um suplemento dietético de vitamina C em uma de três doses (0.5, 1 ou 2 mg/dia) administradas de uma de duas maneiras (como ácido ascórbico ou em suco de laranja). Havia dez porquinhos-da-índia, cinco machos e cinco fêmeas, em cada combinação de dose e método de administração. Após 42 dias na dieta, os porquinhos-da-índia foram sacrificados; os incisivos foram removidos e seccionados para obter medidas do comprimento dos odontoblastos — células que são importantes para o desenvolvimento dos dentes. Poderia haver múltiplas medições (odontoblastos) por animal, então os comprimentos foram médios para fornecer um único valor para cada animal. O resultado de interesse é o comprimento médio dos odontoblastos incisivos; isso foi medido em micrômetros. O interesse é como a dose e o método de administração influenciam o resultado.
Crampton (1947, p. 495)
Você pode estar se perguntando o que um estudo de células em porquinhos-da-índia tem a ver com o fornecimento de vitamina C aos humanos. Os pesquisadores haviam estabelecido que o crescimento dos odontoblastos em porquinhos-da-índia era sensível à ingestão de vitamina C e, portanto, medir os odontoblastos poderia fornecer uma maneira de medir de forma confiável a absorção de vitamina C. Havia um debate sobre se as fontes ‘biológicas’ (por exemplo, suco de laranja) de vitamina C eram melhores do que as formas ‘químicas’ (por exemplo, ácido ascórbico), e isso foi explorado em um modelo animal.
ANOVA unifatorial independente de Fisher com o fator tipo de suplemento resulta na não significância deste fator (\(p=0.060\)) com \(\alpha=5\%\). Os resultados a seguir mostram que há o efeito de interação entre tipo de suplemento e dose (\(p=0.022\)) com \(\alpha=5\%\). Há efeitos simples de tipo de suplemento e dose. Então, suco de laranja tem efeito maior do que ácido ascórbico. As três doses têm efeitos progressivamente maiores.
suppressMessages(suppressWarnings(
invisible(Sys.setlocale("LC_ALL", "pt_BR.UTF-8"))))
# The Effect of Vitamin C on Tooth Growth in Guinea Pigs
# len: numeric Tooth length
# supp: factor Supplement type (VC or OJ).
# dose: numeric Dose in milligrams/day
Dados <- datasets::ToothGrowth
Dados$supp <- factor(Dados$supp,
levels=c("VC", "OJ"))
print(head(Dados)) len supp dose
1 4.2 VC 0.5
2 11.5 VC 0.5
3 7.3 VC 0.5
4 5.8 VC 0.5
5 6.4 VC 0.5
6 10.0 VC 0.5 len supp dose
55 24.8 OJ 2
56 30.9 OJ 2
57 26.4 OJ 2
58 27.3 OJ 2
59 29.4 OJ 2
60 23.0 OJ 2# print(skimr::skim_without_charts(Dados))
Dados$dose <- factor(Dados$dose)
print(ftable(Dados$supp, Dados$dose)) 0.5 1 2
VC 10 10 10
OJ 10 10 10 len supp dose
Min. : 4.20 VC:30 0.5:20
1st Qu.:13.07 OJ:30 1 :20
Median :19.25 2 :20
Mean :18.81
3rd Qu.:25.27
Max. :33.90 
alfa <- 0.05
g1 <- nlevels(Dados$supp)
alfaBonf.supp <- alfa/g1
g2 <- nlevels(Dados$dose)
alfaBonf.dose <- alfa/g2
gplots::plotmeans(len~supp,
connect=FALSE,
col="black",
barcol="black",
p=1-alfaBonf.supp,
main="Porquinho-da-Índia",
data=Dados)
Call:
lm(formula = len ~ supp, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-12.7633 -5.7633 0.4367 5.5867 16.9367
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 16.963 1.366 12.418 <2e-16 ***
suppOJ 3.700 1.932 1.915 0.0604 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 7.482 on 58 degrees of freedom
Multiple R-squared: 0.05948, Adjusted R-squared: 0.04327
F-statistic: 3.668 on 1 and 58 DF, p-value: 0.06039Anova Table (Type II tests)
Response: len
Sum Sq Df F value Pr(>F)
supp 205.4 1 3.6683 0.06039 .
Residuals 3246.9 58
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1alfaBonf <- alfa/(g1*g2)
gplots::plotmeans(len~interaction(dose,supp),
connect=FALSE,
col="black",
barcol="black",
p=1-alfaBonf,
main="Porquinho-da-Índia",
data=Dados)
Call:
lm(formula = len ~ supp * dose, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-8.20 -2.72 -0.27 2.65 8.27
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.980 1.148 6.949 4.98e-09 ***
suppOJ 5.250 1.624 3.233 0.00209 **
dose1 8.790 1.624 5.413 1.46e-06 ***
dose2 18.160 1.624 11.182 1.13e-15 ***
suppOJ:dose1 0.680 2.297 0.296 0.76831
suppOJ:dose2 -5.330 2.297 -2.321 0.02411 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 3.631 on 54 degrees of freedom
Multiple R-squared: 0.7937, Adjusted R-squared: 0.7746
F-statistic: 41.56 on 5 and 54 DF, p-value: < 2.2e-16Anova Table (Type II tests)
Response: len
Sum Sq Df F value Pr(>F)
supp 205.35 1 15.572 0.0002312 ***
dose 2426.43 2 92.000 < 2.2e-16 ***
supp:dose 108.32 2 4.107 0.0218603 *
Residuals 712.11 54
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1[1] 18.8 supp len
1 VC 16.96333
2 OJ 20.66333 dose len
1 0.5 10.605
2 1 19.735
3 2 26.100 dose supp len
1 0.5 VC 7.98
2 1 VC 16.77
3 2 VC 26.14
4 0.5 OJ 13.23
5 1 OJ 22.70
6 2 OJ 26.06# Somas de Quadrados
g <- nlevels(Dados$supp)
b <- nlevels(Dados$dose)
n <- nrow(Dados) / (g*b)
p <- 1
SS_total <- sum((Dados$len - X_bar)^2)
SS_A <- b * n * sum((X_bar_l$len - X_bar)^2)
SS_B <- g * n * sum((X_bar_k$len - X_bar)^2)
SS_error <- sum((Dados$len - rep(X_bar_lk$len, times=rep(n,g*b)))^2)
SS_AB <- SS_total - SS_A - SS_B - SS_error
# Graus de Liberdade
df_A <- g - 1
df_B <- b - 1
df_AB <- (g - 1)*(b - 1)
df_error <- g*b*(n - 1)
# Quadrado Médio
MS_A <- SS_A / df_A
MS_B <- SS_B / df_B
MS_AB <- SS_AB / df_AB
MS_error <- SS_error / df_error
# Estatística F
F_A <- MS_A / MS_error
F_B <- MS_B / MS_error
F_AB <- MS_AB / MS_error
# Valor-p (usando a distribuição F)
p_value_A <- formatC(1-pf(F_A, df_A, df_error),
format="e", digits=2)
p_value_B <- formatC(1-pf(F_B, df_B, df_error),
format="e", digits=2)
p_value_AB <- formatC(1-pf(F_AB, df_AB, df_error),
format="e", digits=2)
# Tabela ANOVA
anova_table <- data.frame(
Source = c("Factor A", "Factor B", "Interaction AB", "Error", "Total"),
SS = c(SS_A, SS_B, SS_AB, SS_error, SS_total),
df = c(df_A, df_B, df_AB, df_error, g*b*n - 1),
MS = c(MS_A, MS_B, MS_AB, MS_error, NA),
"F value" = c(F_A, F_B, F_AB, NA, NA),
"p value" = c(p_value_A, p_value_B, p_value_AB, NA, NA),
stringsAsFactors = FALSE,
check.names = FALSE
)
print(anova_table, digits=3, row.names=FALSE) Source SS df MS F value p value
Factor A 205 1 205.3 15.57 2.31e-04
Factor B 2426 2 1213.2 92.00 0.00e+00
Interaction AB 108 2 54.2 4.11 2.19e-02
Error 712 54 13.2 NA <NA>
Total 3452 59 NA NA <NA># Solução model.matrix e OLS
# https://www.maths.usyd.edu.au/u/UG/SM/STAT3022/r/current/Lecture/lecture18_2020JC.html#17
fit <- lm(len~supp*dose,
data=Dados)
print(anova(fit))Analysis of Variance Table
Response: len
Df Sum Sq Mean Sq F value Pr(>F)
supp 1 205.35 205.35 15.572 0.0002312 ***
dose 2 2426.43 1213.22 92.000 < 2.2e-16 ***
supp:dose 2 108.32 54.16 4.107 0.0218603 *
Residuals 54 712.11 13.19
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova Table (Type II tests)
Response: len
Sum Sq Df F value Pr(>F)
supp 205.35 1 15.572 0.0002312 ***
dose 2426.43 2 92.000 < 2.2e-16 ***
supp:dose 108.32 2 4.107 0.0218603 *
Residuals 712.11 54
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1y <- Dados$len
X0 <- model.matrix(~ 1, data=Dados)
XA <- model.matrix(~ -1 + supp, data=Dados)
XB <- model.matrix(~ -1 + dose, data=Dados)
XAB <- model.matrix(~ -1 + supp:dose, data=Dados)
PV0 <- X0 %*% solve(t(X0) %*% X0) %*% t(X0)
PWA <- XA %*% solve(t(XA) %*% XA) %*% t(XA) - PV0
PWB <- XB %*% solve(t(XB) %*% XB) %*% t(XB) - PV0
PVAB <- XAB %*% solve(t(XAB) %*% XAB) %*% t(XAB)
PWAB <- PVAB - PWA - PWB - PV0
PRes <- diag(nrow(Dados)) - PVAB
SSA <- t(y) %*% PWA %*% y
SSB <- t(y) %*% PWB %*% y
SSAB <- t(y) %*% PWAB %*% y
SSE <- t(y) %*% PRes %*% y
dfA <- sum(diag(PWA))
dfB <- sum(diag(PWB))
dfAB <- sum(diag(PWAB))
dfE <- sum(diag(PRes))
dfT <- dfA+dfB+dfAB+dfE
SST <- SSA+SSB+SSAB+SSE
R2A <- SSA/SST
R2B <- SSB/SST
R2AB <- SSAB/SST
R2E <- SSE/SST
R2 <- 1 - SSE/SST
eta2partialA <- R2A/(R2A+1-R2)
eta2partialB <- R2B/(R2B+1-R2)
eta2partialAB <- R2AB/(R2AB+1-R2)
FA <- (SSA/dfA)/(SSE/dfE)
pA <- formatC(1-pf(FA, dfA, dfE),
format="e", digits=2)
FB <- (SSB/dfB)/(SSE/dfE)
pB <- formatC(1-pf(FB, dfB, dfE),
format="e", digits=2)
FAB <- (SSAB/dfAB)/(SSE/dfE)
pAB <- formatC(1-pf(FAB, dfAB, dfE),
format="e", digits=2)
R2 <- 1 - SSE/SST
df1 <- dfT-dfE
df2 <- dfE
F_omnibus <- (R2/df1)/((1-R2)/df2) # Pestana & Gageiro, 2005, p. 77
pv <- 1-pf(F_omnibus, df1, df2)
if(pv<2.2e-16) p_omnibus <- "< 2.2e-16"
if(pv>2.2e-16) p_omnibus <- paste0("< ", formatC(pv, format="e", digits=2))
Fcrit <- qf(1-alfa, df1, df2)
cat("Omnibus F(",df1, ",", df2, ") = ", round(Fcrit,2),
", F = ", round(F_omnibus,2),
", p = ", p_omnibus, "\nR^2 = eta^2 = ", round(R2,4), sep="")Omnibus F(5,54) = 2.39, F = 41.56, p = < 2.2e-16
R^2 = eta^2 = 0.7937ANOVA2wtable <- data.frame("Source"=c("supp", "dose", "supp:dose", "Error", "Total"),
# "df"=c(qr(PWA)$rank, qr(PWB)$rank, qr(PWAB)$rank, qr(PRes)$rank),
"df"=c(dfA, dfB, dfAB, dfE, dfT),
"SS"=c(SSA, SSB, SSAB, SSE, SST),
"F value"=c(FA, FB, FAB, NA, NA),
"p value"=c(pA, pB, pAB, NA, NA),
"R^2"=c(R2A, R2B, R2AB, R2E, NA),
"eta^2 partial"=c(eta2partialA, eta2partialB, eta2partialAB, NA, NA),
stringsAsFactors = FALSE,
check.names = FALSE)
print(ANOVA2wtable, digits=2, row.names=FALSE) Source df SS F value p value R^2 eta^2 partial
supp 1 205 15.6 2.31e-04 0.059 0.22
dose 2 2426 92.0 0.00e+00 0.703 0.77
supp:dose 2 108 4.1 2.19e-02 0.031 0.13
Error 54 712 NA <NA> 0.206 NA
Total 59 3452 NA <NA> NA NAeta2 <- effectsize::eta_squared(anv,
partial=TRUE,
generalized=FALSE,
ci=1-alfa/3,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 (partial) | 98.3333% CI | interpret
---------------------------------------------------------
supp | 0.2238 | [0.0337, 0.4383] | large
dose | 0.7731 | [0.6299, 0.8496] | large
supp:dose | 0.1320 | [0.0000, 0.3346] | mediumeta2 <- effectsize::eta_squared(anv,
partial=FALSE,
generalized=FALSE,
ci=1-alfa/3,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 | 98.3333% CI | interpret
-------------------------------------------------
supp | 0.0595 | [0.0000, 0.2525] | small
dose | 0.7029 | [0.5236, 0.8020] | large
supp:dose | 0.0314 | [0.0000, 0.1848] | small# Regressão Linear Múltipla de two-way ANOVA
fit <- lm(len~supp*dose,
data=Dados)
print(car::Anova(fit))Anova Table (Type II tests)
Response: len
Sum Sq Df F value Pr(>F)
supp 205.35 1 15.572 0.0002312 ***
dose 2426.43 2 92.000 < 2.2e-16 ***
supp:dose 108.32 2 4.107 0.0218603 *
Residuals 712.11 54
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = len ~ supp * dose, data = Dados)
Residuals:
Min 1Q Median 3Q Max
-8.20 -2.72 -0.27 2.65 8.27
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 7.980 1.148 6.949 4.98e-09 ***
suppOJ 5.250 1.624 3.233 0.00209 **
dose1 8.790 1.624 5.413 1.46e-06 ***
dose2 18.160 1.624 11.182 1.13e-15 ***
suppOJ:dose1 0.680 2.297 0.296 0.76831
suppOJ:dose2 -5.330 2.297 -2.321 0.02411 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 3.631 on 54 degrees of freedom
Multiple R-squared: 0.7937, Adjusted R-squared: 0.7746
F-statistic: 41.56 on 5 and 54 DF, p-value: < 2.2e-16y <- Dados$len
X <- model.matrix(len~supp*dose, data=Dados)
coeff <- solve(t(X) %*% X) %*% t(X) %*% y # OLS
cat("\nCoefficients:\n")
Coefficients: Estimate
(Intercept) 7.98
suppOJ 5.25
dose1 8.79
dose2 18.16
suppOJ:dose1 0.68
suppOJ:dose2 -5.33[1] 2740.103[1] 712.106[1] 3452.209R^2 = eta^2 omnibus = 0.79df1 <- dim(X)[2] - 1
N <- dim(X)[1]
df2 <- N-df1-1
F_omnibus <- (R2/df1)/((1-R2)/df2) # Pestana & Gageiro, 2005, p. 77
pv <- 1-pf(F_omnibus, df1, df2)
if(pv<2.2e-16) p_omnibus <- "< 2.2e-16"
if(pv>2.2e-16) p_omnibus <- paste0("< ", formatC(pv, format="e", digits=2))
Fcrit <- qf(1-alfa, df1, df2)
cat("Omnibus F(",df1, ",", df2, ") = ", round(Fcrit,2),
", F = ", round(F_omnibus,2),
", p ", p_omnibus, "\nR^2 = eta^2 = ", round(R2,2), sep="")Omnibus F(5,54) = 2.39, F = 41.56, p < 2.2e-16
R^2 = eta^2 = 0.79# print(sjPlot::tab_model(fit,
# p.adjust="holm",
# encoding="Windows-1252",
# title="Effect of Vitamin C on Tooth Growth in Guinea Pigs: holm",
# show.reflvl=TRUE,
# show.ngroups=TRUE,
# show.stat=TRUE,
# show.df=TRUE,
# show.se=TRUE,
# digits=4,
# digits.p=4,
# digits.rsq=6,
# digits.re=6))
# Solução GLM
fit <- lm(len~supp,
data=Dados)
print(anova(fit))Analysis of Variance Table
Response: len
Df Sum Sq Mean Sq F value Pr(>F)
supp 1 205.4 205.35 3.6683 0.06039 .
Residuals 58 3246.9 55.98
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analysis of Variance Table
Response: len
Df Sum Sq Mean Sq F value Pr(>F)
supp 1 205.35 205.35 14.017 0.0004293 ***
dose 2 2426.43 1213.22 82.811 < 2.2e-16 ***
Residuals 56 820.43 14.65
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analysis of Variance Table
Response: len
Df Sum Sq Mean Sq F value Pr(>F)
supp 1 205.35 205.35 15.572 0.0002312 ***
dose 2 2426.43 1213.22 92.000 < 2.2e-16 ***
supp:dose 2 108.32 54.16 4.107 0.0218603 *
Residuals 54 712.11 13.19
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova Table (Type II tests)
Response: len
Sum Sq Df F value Pr(>F)
supp 205.35 1 15.572 0.0002312 ***
dose 2426.43 2 92.000 < 2.2e-16 ***
supp:dose 108.32 2 4.107 0.0218603 *
Residuals 712.11 54
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1eta2 <- effectsize::eta_squared(anv,
partial=TRUE,
generalized=FALSE,
ci=1-alfa/3,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 (partial) | 98.3333% CI | interpret
---------------------------------------------------------
supp | 0.2238 | [0.0337, 0.4383] | large
dose | 0.7731 | [0.6299, 0.8496] | large
supp:dose | 0.1320 | [0.0000, 0.3346] | mediumeta2 <- effectsize::eta_squared(anv,
partial=FALSE,
generalized=FALSE,
ci=1-alfa/3,
alternative="two.sided",
verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)# Effect Size for ANOVA (Type II)
Parameter | Eta2 | 98.3333% CI | interpret
-------------------------------------------------
supp | 0.0595 | [0.0000, 0.2525] | small
dose | 0.7029 | [0.5236, 0.8020] | large
supp:dose | 0.0314 | [0.0000, 0.1848] | smallo.par <- par()
alfaBonf_int <- alfa/(g1*g2)
fit.means <- phia::interactionMeans(fit)
plot(fit.means,
errorbar=paste0("ci",
round((1-alfaBonf)*100,4)),
abbrev.levels=FALSE)
par(o.par)Warning in par(o.par): parâmetro gráfico "cin" não pode ser especificadoWarning in par(o.par): parâmetro gráfico "cra" não pode ser especificadoWarning in par(o.par): parâmetro gráfico "csi" não pode ser especificadoWarning in par(o.par): parâmetro gráfico "cxy" não pode ser especificadoWarning in par(o.par): parâmetro gráfico "din" não pode ser especificado
Warning in par(o.par): parâmetro gráfico "page" não pode ser especificado# alfaBonf_int <- alfa/(length(unique(Dados$supp))*
# length((unique(Dados$dose))))
# plot(effects::effect(c("supp", "dose"), fit,
# confidence.level=1-alfaBonf_int),
# multiline=TRUE, ci.style="bars")
# alfaBonf_fat1 <- alfa/(length(unique(Dados$supp)))
# plot(effects::effect(c("supp"), fit, confidence.level=1-alfaBonf_fat1),
# ci.style="bars")
# alfaBonf_fat2 <- alfa/(length(unique(Dados$dose)))
# plot(effects::effect(c("dose"), fit, confidence.level=1-alfaBonf_fat2),
# ci.style="bars")
cat("supp simple main effect test in the presence of interaction with dose \n")supp simple main effect test in the presence of interaction with dose F Test:
P-value adjustment method: holm
Value SE Df Sum of Sq F Pr(>F)
0.5 -5.25 1.624 1 137.81 10.4505 0.004185 **
1 -5.93 1.624 1 175.82 13.3330 0.001769 **
2 0.08 1.624 1 0.03 0.0024 0.960893
Residuals 54 712.11
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Post hoc test: supp|supp:dose (there is no simple main effect)emm <- emmeans::emmeans(object=fit,
specs=pairwise~supp|supp:dose,
alpha=alfa,
adjust="holm")
print(emm)$emmeans
dose = 0.5:
supp emmean SE df lower.CL upper.CL
VC 7.98 1.15 54 5.68 10.3
OJ 13.23 1.15 54 10.93 15.5
dose = 1:
supp emmean SE df lower.CL upper.CL
VC 16.77 1.15 54 14.47 19.1
OJ 22.70 1.15 54 20.40 25.0
dose = 2:
supp emmean SE df lower.CL upper.CL
VC 26.14 1.15 54 23.84 28.4
OJ 26.06 1.15 54 23.76 28.4
Confidence level used: 0.95
$contrasts
dose = 0.5:
contrast estimate SE df t.ratio p.value
VC - OJ -5.25 1.62 54 -3.233 0.0021
dose = 1:
contrast estimate SE df t.ratio p.value
VC - OJ -5.93 1.62 54 -3.651 0.0006
dose = 2:
contrast estimate SE df t.ratio p.value
VC - OJ 0.08 1.62 54 0.049 0.9609
dose = 0.5:
supp emmean SE df lower.CL upper.CL .group
VC 7.98 1.15 54 5.33 10.6 a
OJ 13.23 1.15 54 10.58 15.9 b
dose = 1:
supp emmean SE df lower.CL upper.CL .group
VC 16.77 1.15 54 14.12 19.4 a
OJ 22.70 1.15 54 20.05 25.3 b
dose = 2:
supp emmean SE df lower.CL upper.CL .group
OJ 26.06 1.15 54 23.41 28.7 a
VC 26.14 1.15 54 23.49 28.8 a
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 2 estimates
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. dose simple main effect test in the presence of interaction with supp F Test:
P-value adjustment method: holm
dose1 dose2 SE1 SE2 Df Sum of Sq F Pr(>F)
VC -18.16 -9.37 1.624 1.624 2 1649.49 62.541 1.751e-14 ***
OJ -12.83 -3.36 1.624 1.624 2 885.26 33.565 3.363e-10 ***
Residuals 54 712.11
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Post hoc test: dose (there is simple main effect)NOTE: Results may be misleading due to involvement in interactions$emmeans
dose emmean SE df lower.CL upper.CL
0.5 10.6 0.812 54 8.98 12.2
1 19.7 0.812 54 18.11 21.4
2 26.1 0.812 54 24.47 27.7
Results are averaged over the levels of: supp
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
dose0.5 - dose1 -9.13 1.15 54 -7.951 <.0001
dose0.5 - dose2 -15.49 1.15 54 -13.493 <.0001
dose1 - dose2 -6.37 1.15 54 -5.543 <.0001
Results are averaged over the levels of: supp
P value adjustment: holm method for 3 tests 
dose emmean SE df lower.CL upper.CL .group
0.5 10.6 0.812 54 8.6 12.6 a
1 19.7 0.812 54 17.7 21.7 b
2 26.1 0.812 54 24.1 28.1 c
Results are averaged over the levels of: supp
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 3 estimates
P value adjustment: holm method for 3 tests
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same.
Post hoc test: supp:dose$emmeans
supp dose emmean SE df lower.CL upper.CL
VC 0.5 7.98 1.15 54 5.68 10.3
OJ 0.5 13.23 1.15 54 10.93 15.5
VC 1 16.77 1.15 54 14.47 19.1
OJ 1 22.70 1.15 54 20.40 25.0
VC 2 26.14 1.15 54 23.84 28.4
OJ 2 26.06 1.15 54 23.76 28.4
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
VC,dose0.5 - OJ,dose0.5 -5.25 1.62 54 -3.233 0.0105
VC,dose0.5 - VC,dose1 -8.79 1.62 54 -5.413 <.0001
VC,dose0.5 - OJ,dose1 -14.72 1.62 54 -9.064 <.0001
VC,dose0.5 - VC,dose2 -18.16 1.62 54 -11.182 <.0001
VC,dose0.5 - OJ,dose2 -18.08 1.62 54 -11.133 <.0001
OJ,dose0.5 - VC,dose1 -3.54 1.62 54 -2.180 0.1346
OJ,dose0.5 - OJ,dose1 -9.47 1.62 54 -5.831 <.0001
OJ,dose0.5 - VC,dose2 -12.91 1.62 54 -7.949 <.0001
OJ,dose0.5 - OJ,dose2 -12.83 1.62 54 -7.900 <.0001
VC,dose1 - OJ,dose1 -5.93 1.62 54 -3.651 0.0035
VC,dose1 - VC,dose2 -9.37 1.62 54 -5.770 <.0001
VC,dose1 - OJ,dose2 -9.29 1.62 54 -5.720 <.0001
OJ,dose1 - VC,dose2 -3.44 1.62 54 -2.118 0.1346
OJ,dose1 - OJ,dose2 -3.36 1.62 54 -2.069 0.1346
VC,dose2 - OJ,dose2 0.08 1.62 54 0.049 0.9609
P value adjustment: holm method for 15 tests 
supp dose emmean SE df lower.CL upper.CL .group
VC 0.5 7.98 1.15 54 4.83 11.1 a
OJ 0.5 13.23 1.15 54 10.08 16.4 b
VC 1 16.77 1.15 54 13.62 19.9 b
OJ 1 22.70 1.15 54 19.55 25.8 c
OJ 2 26.06 1.15 54 22.91 29.2 c
VC 2 26.14 1.15 54 22.99 29.3 c
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 6 estimates
P value adjustment: holm method for 15 tests
significance level used: alpha = 0.05
NOTE: If two or more means share the same grouping symbol,
then we cannot show them to be different.
But we also did not show them to be the same. Exemplo de GLM univariado multifatorial independente em Silveira, Tempski, Mayer, Enns, Peleias, Martins & Siqueira (2025)
Título: Comparação entre amostras aleatórias e de conveniência em levantamento multicêntrico para avaliar a qualidade de vida de estudantes de medicina.
Resumo: A avaliação da saúde mental e física de estudantes de medicina é dificultada pela obtenção de amostras aleatórias e pelas limitações de amostras de conveniência. Conduzimos um estudo multicêntrico para avaliar o ambiente educacional, a qualidade de vida e a competência emocional de estudantes de medicina. Entre 2011 e 2012, um total de 1.350 estudantes selecionados aleatoriamente de 22 escolas e 1.201 estudantes voluntários de 50 escolas em todo o Brasil completaram todos os questionários (WHOQOL-BREF, VERAS-Q, IRI, RS-14, BDI, PSQI, ESS, IDATE, MBI e DREEM). Monitoramento, apoio a pesquisadores locais e estratégias de feedback personalizado foram aplicados para assegurar a participação dos estudantes aleatorizados, alcançando taxa de resposta de 81,8%. A plataforma também ficou disponível aos voluntários. A análise estatística examinou o efeito dessas duas estratégias de recrutamento por meio de modelos lineares gerais, controlando por sexo, idade, massa corporal, ano do curso, atividade física e equivalentes metabólicos (METs), tipo de escola, população da cidade e localização. Adotou-se nível de significância de 5% e tamanhos de efeito estimados pelo eta-quadrado de Cohen, ambos com correção de Bonferroni nas variáveis de interesse. O grupo de voluntários apresentou mais mulheres, menos alunos dos anos finais do curso e maior proporção de estudantes de escolas privadas e de cidades maiores. Essas variáveis explicaram em grande parte as diferenças estatisticamente significantes e os tamanhos de efeito observados entre os grupos aleatorizado e voluntário. Em conclusão, embora tenham sido identificadas lições e estratégias motivacionais úteis, o esforço considerável necessário para obter alta adesão via busca ativa pode não se justificar, pois os resultados entre amostras aleatórias e de voluntários mostraram diferenças mínimas nas respostas aos questionários. Nossos achados sugerem que estudos futuros considerem estratégias motivacionais mais leves para reduzir a carga das equipes de pesquisa. Dados e scripts em R para replicar a análise estatística estão disponíveis no Harvard Dataverse em https://doi.org/10.7910/DVN/YECV8E.
Exemplo de GLM misto multifatorial com interação em Takayanagi, Siqueira, Silveira & Valentova (2024)
Título: O que Diferentes Pessoas Buscam em um Parceiro? Efeitos do Sexo, da Orientação Sexual e das Estratégias de Acasalamento nas Preferências de Parceiro
A artigo apresenta uma aplicação de GLMM (modelo linear misto geral) com fatores entre e intraparticipantes.
Os dados e código R estão disponíveis em OSF Home.
Resumo: As preferências de parceiro são um diferencial importante na formação de relacionamentos e na aptidão evolutiva, variando de acordo com fatores individuais, ecológicos e sociais. Neste estudo, avaliamos a variação na preferência por inteligência, bondade, atratividade física, saúde e nível socioeconômico entre indivíduos de diferentes sexos e orientações sexuais em uma amostra brasileira. Analisamos as pontuações de preferência de 778 homens e mulheres heterossexuais, bissexuais e homossexuais em três tarefas de orçamento de parceiro (baixo vs. médio vs. alto orçamento) e sua associação com sociosexualidade, estilos de apego, homogamia e disposição para engajar-se em relacionamentos de curto e longo prazo. Os resultados indicaram uma ordem global de preferência por traços, com a inteligência em primeiro lugar, seguida por bondade, atratividade física, saúde e, por último, status socioeconômico. Diferenças típicas de sexo foram observadas principalmente no grupo heterossexual, e combinações específicas de sexo e orientação sexual estiveram associadas a variações na preferência por atratividade física, bondade e status socioeconômico. Também encontramos associações únicas de outras variáveis com preferências de parceiro e com a disposição para engajar-se em relacionamentos de curto ou longo prazo. Ao explorar as preferências de parceiro de indivíduos não heterossexuais de um país latino-americano, um grupo sub-representado na pesquisa em psicologia evolucionista, nossos resultados ajudam a compreender os fatores universais e específicos que orientam as preferências de parceiro e o comportamento sexual humano.
Análises Estatísticas: Os dados foram analisados no software
R, versão 4.2.1 (Core Team, 2022). Computamos modelos lineares mistos
gerais (GLMM) utilizando a função lmer do pacote
lmerTest (Kuznetsova et al., 2017).
Na análise principal, utilizamos um GLMM para avaliar como o nível de preferência por cada traço variou em função das características individuais. Usamos a pontuação de preferência pelo traço como variável dependente; o identificador do participante como efeito aleatório; traço do parceiro, sexo e orientação sexual como variáveis fixas nominais; e tamanho do orçamento, sociosexualidade, apego ansioso, apego evitativo e pontuação de autoavaliação do traço como variáveis fixas intervalares. Idade, estado civil, renda mensal domiciliar per capita e nível educacional foram incluídos como variáveis de controle.
Nas análises adicionais, comparamos as pontuações de disposição para engajar-se em relacionamentos de curto e longo prazo entre condições de orçamento, sexo e orientações sexuais, utilizando uma análise de variância (ANOVA) de medidas repetidas com quatro fatores. Também usamos GLMM para explorar como características individuais e diferenças nas preferências de parceiro influenciaram a inclinação para formar relacionamentos de curto ou longo prazo com o parceiro hipotético. Usamos o nível de disposição para engajar-se em relacionamentos de curto ou longo prazo como variável dependente de cada modelo; o identificador do participante como efeito aleatório; sexo e orientação sexual como variáveis fixas nominais; e tamanho do orçamento, sociosexualidade, apego ansioso, apego evitativo e pontuações de preferência por traços como variáveis fixas intervalares. Mais uma vez, idade, estado civil, renda mensal domiciliar per capita e nível educacional foram incluídos como variáveis de controle.
Referências
- Crampton EW (1947) O crescimento do odontoblasto dos dentes incisivos como critério de ingestão de vitamina C do porquinho-da-índia. The Journal of Nutrition 33(5): 491–504.
- Dancey CP & Reidy J (2006) Estatística sem Matemática para Psicologia. 3a ed. Porto Alegre: Bookman.
- Dancey CP & Reidy J (2019) Estatística sem Matemática para Psicologia. 7a ed. Porto Alegre: PENSO.
- Gulick, D & Gould, TJ (2009) Effects of ethanol and caffeine on behavior in C57BL/6 mice in the plus-maze discriminative avoidance task. Behavioral Neuroscience 123(6): 1271–78. https://doi.org/10.1037/a0017610
- Howell GT & Lacroix GL (2012) Decomposing interactions using GLM in combination with the COMPARE, LMATRIX and MMATRIX subcommands in SPSS. Tutorials in Quantitative Methods for Psychology 8(1):1-22. DOI:10.20982/tqmp.08.1.p001
- Howland J et al. (2011) The acute effects of caffeinated versus non-caffeinated alcoholic beverage on driving performance and attention/reaction time. Addiction. 106(2):335-41. doi:10.1111/j.1360-0443.2010.03219.x
- Silveira PSP, Tempski PZ, Mayer FB, Enns SC, Peleias M, Martins MdA, Siqueira JO (2025) Comparison between random and convenience samples in a multicenter survey to evaluate medical students’ quality of life. PLoS One 20(10): e0332850. https://doi.org/10.1371/journal.pone.0332850
- Smith, CE & Cribbie, R (2014) Factorial ANOVA with unbalanced data: A fresh look at the types of sums of squares. Journal of Data Science 12(2): 385–404.
- Takayanagi, JFGB, Siqueira, JO, Silveira, PSP, Valentova, JV (2024) What Do Different People Look for in a Partner? Effects of Sex, Sexual Orientation, and Mating Strategies on Partner Preferences. Arch Sex Behav 53, 981–1000 (2024). https://doi.org/10.1007/s10508-023-02767-4. Dados e código R: https://osf.io/pc6rn/?view%2520only=38097a82dd274d34bbddc8dcc2258bb9