Bastão de Asclépio & Distribuição Normal

invisible(Sys.setlocale("LC_CTYPE", "pt_BR.UTF-8"))
invisible(Sys.setlocale("LC_ALL", "pt_BR.UTF-8"))

options(warn=-1)
suppressMessages(library(bootES, 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(ggplot2, warn.conflicts = FALSE))
suppressMessages(library(gplots, warn.conflicts = FALSE))
suppressMessages(library(Hmisc, warn.conflicts=FALSE))
suppressMessages(library(jmv, warn.conflicts=FALSE))
suppressMessages(library(knitr, warn.conflicts=FALSE))
suppressMessages(library(lmerTest, warn.conflicts=FALSE))
suppressMessages(library(MBESS, warn.conflicts=FALSE))
suppressMessages(library(multcomp, warn.conflicts=FALSE))
suppressMessages(library(psych, warn.conflicts=FALSE))
suppressMessages(library(RcmdrMisc, warn.conflicts=FALSE))
suppressMessages(library(readxl, warn.conflicts=FALSE))
suppressMessages(library(Rmisc, warn.conflicts=FALSE))
suppressMessages(library(rstatix, warn.conflicts=FALSE))
suppressMessages(library(sjPlot, warn.conflicts=FALSE))
suppressMessages(library(WRS2, warn.conflicts=FALSE))
source("summarySEwithin2.R")
# source("eiras.OneWayANOVA.Welch.WORD.R")
options(warn=0)

Material

HTML de Rmarkdown em RPubs

Introdução

Esta prática é relativa ao capítulo 10 do livro-texto.

Os seguintes arquivos de dados serão usados nesta aula:

ANOVA unifatorial independente de Welch

Supondo que IMC é uma medida melhor que a MCT isoladamente, pois leva a estatura também em consideração, verifique com a ANOVA de Welch se entre os estudantes de medicina de 2015-2018 há efeito de atividade física sobre IMC. Sugerimos que analise cada sexo separadamente.

Solução em ANOVAWelch.R

ANOVA unifatorial independente de Welch sem dados brutos

Podemos realizar os testes omnibus e post hoc da ANOVA unifatorial de Welch quando temos apenas as médias, os desvios padrão e os tamanhos de cada grupo amostrado. Desenvolvemos uma função que lida com esta situação e demonstra-se seu funcionamento colocando-se as informações que obtivemos do arquivo de Biometria em ExemploWelchSDB_Biometria.xlsx. Verificar que são obtidos os mesmos resultados que conseguimos com os dados brutos.

Solução em ANOVAWelch_SemDadosBrutos.R

Análise de poder retrospectivo (plug in) de ANOVA unifatorial independente balanceada de Fisher (Gerard et al., 1998)

Gerard, PD et al. (1998) Limits of retrospective power analysis. Journal of Wildlife Management 62(2): 801-7. https://doi.org/10.2307/3802357

A análise de poder retrospectivo, de maneira geral, para o modelo linear misto geral (GLMM), que inclui ANOVA, é inútil e conceitualmente equivocada.

Exemplo em ANOVApoderretrospectivo.R

ANOVA unifatorial relacionada

pulse

Utilize os dados de NewDrug.rds para testar o efeito do fator droga considerando as três observações de pulso (variáveis pulse1, pulse2 e pulse3).

Solução em NewDrug_ANOVAintra_pulse.R

resp

Repita a análise para resp.

Solução em NewDrug_ANOVAintra_resp.R

com uma terceira condição de tratamento

Suponha que um terceiro grupo de pacientes recebeu um tratamento tradicional. Utilize os dados de NewDrug3long.rds para testar o efeito do fator droga considerando as três medidas de pulse e resp. A transformação de NewDrug3.xlsx para NewDrug3long.rds está em NewDrug3xlsxToNewDrug3longrds.R.

Solução em NewDrug_ANOVAintra_3cond.R

quatro drogas sem medidas repetidas

Utilize os dados de Drug4.xlsx para testar o efeito do fator droga na resposta (variável response).

Solução em NewDrug_ANOVAintra_drug4.R

Pacientes diferentes?

Suponha que foram observadas em seis dias consecutivos frequência cardíaca (FC) em 3 pacientes.

Paciente	Dia1	Dia2	Dia3	Dia4	Dia5	Dia6
1	70	72	74	69	71	73
2	75	76	78	74	77	79
3	65	66	68	64	67	69

Série temporal é uma sequência em ordem cronológica equiespaçada.

Neste caso, há 3 séries temporais com 6 valores diários de FC.

Questão

Os três pacientes têm FC médias iguais?

Solução em PacienteDiaFC_ANOVAintra.R

Cinco condições independentes: tamanho de amostra muito grande

ANOVA de Fisher: idade de óbito e raça/cor

ANOVAuni_Fisher_DOBR2020.R

Ler o arquivo DOBR2020.rds.

Dados <- readRDS("DOBR2020.rds")
# sjPlot::view_df(Dados)

Questões

Usar ANOVA de Fisher testar o efeito de raça/cor sobre a idade de óbito de adultos de adultos (18+).
Interpretar o valor p e d Cohen. Qual faz mais sentido: valor p ou d de Cohen?

boxplot(data=Dados, 
        IdadeObito~RACACOR)

print(psych::describeBy(IdadeObito~RACACOR,
                        mat=2,
                        digits=6,
                        data=Dados[Dados$IdadeObito >= 18,]))

            item   group1 vars      n     mean       sd median  trimmed     mad
IdadeObito1    1   Branca    1 746441 71.41314 16.67186     74 72.90715 16.3086
IdadeObito2    2    Preta    1 129774 65.61251 17.89297     67 66.62632 17.7912
IdadeObito3    3  Amarela    1   9323 74.39730 15.73094     77 75.89516 14.8260
IdadeObito4    4    Parda    1 576586 64.46338 19.26617     67 65.65579 19.2738
IdadeObito5    5 Indígena    1   4360 64.41193 21.56470     68 65.46359 22.2390
            min max range      skew  kurtosis       se
IdadeObito1  18 125   107 -0.794643  0.390440 0.019297
IdadeObito2  18 117    99 -0.464506 -0.231533 0.049669
IdadeObito3  18 111    93 -0.933965  0.916803 0.162921
IdadeObito4  18 120   102 -0.482287 -0.441473 0.025372
IdadeObito5  18 117    99 -0.380515 -0.727123 0.326588

alfa <- 0.05
gplots::plotmeans(data=Dados,
                  subset=IdadeObito >= 18,
                  IdadeObito~RACACOR,
                  p=1-alfa/5,
                  connect=FALSE,
                  barcol="black",
                  main="IC95% Bonferroni")

Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, :
zero-length arrow is of indeterminate angle and so skipped
Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, :
zero-length arrow is of indeterminate angle and so skipped
Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, :
zero-length arrow is of indeterminate angle and so skipped

Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, :
zero-length arrow is of indeterminate angle and so skipped
Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, :
zero-length arrow is of indeterminate angle and so skipped
Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, :
zero-length arrow is of indeterminate angle and so skipped

cat("\nFisher's One-way ANOVA ")


Fisher's One-way ANOVA

cat("\nStatistical analysis: omnibus test")


Statistical analysis: omnibus test

modelo <- lm(IdadeObito~RACACOR, 
             data=Dados[Dados$IdadeObito >= 18,])
# ANOVA da one-way ANOVA
cat("\nANOVA")


ANOVA

print(anv <- car::Anova(modelo))

Anova Table (Type II tests)

Response: IdadeObito
             Sum Sq      Df F value    Pr(>F)    
RACACOR    17044665       4   13370 < 2.2e-16 ***
Residuals 467375277 1466479                      
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

print(reg <- summary(modelo))


Call:
lm(formula = IdadeObito ~ RACACOR, data = Dados[Dados$IdadeObito >= 
    18, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-56.397 -10.463   2.587  13.537  55.537 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)     71.41314    0.02066 3456.06   <2e-16 ***
RACACORPreta    -5.80063    0.05369 -108.03   <2e-16 ***
RACACORAmarela   2.98416    0.18604   16.04   <2e-16 ***
RACACORParda    -6.94976    0.03130 -222.03   <2e-16 ***
RACACORIndígena -7.00121    0.27115  -25.82   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 17.85 on 1466479 degrees of freedom
  (40494 observations deleted due to missingness)
Multiple R-squared:  0.03519,   Adjusted R-squared:  0.03518 
F-statistic: 1.337e+04 on 4 and 1466479 DF,  p-value: < 2.2e-16

cat("R^2 = eta^2 =", reg$r.squared)

R^2 = eta^2 = 0.03518572

eta2 <- effectsize::eta_squared(modelo,
                                partial=FALSE,
                                generalized=FALSE,
                                ci=1-alfa,
                                alternative="two.sided",
                                verbose=TRUE)
eta2$interpret <- effectsize::interpret_eta_squared(eta2$Eta2)
print(eta2, digits=4)

# Effect Size for ANOVA (Type I)

Parameter |   Eta2 |           95% CI | interpret
-------------------------------------------------
RACACOR   | 0.0352 | [0.0346, 0.0358] |     small

# gráfico do test F
reg <- summary(modelo)
fobs <- 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="One-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: igualdade das médias populacionais",
         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, 
       bty="n")

cat("\nFisher's One-way ANOVA")


Fisher's One-way ANOVA

cat("\nStatistical analysis: Post hoc tests")


Statistical analysis: Post hoc tests

emm <- emmeans::emmeans(modelo, 
                        pairwise~"RACACOR", 
                        adjust="holm",
                        level=1-alfa)
print(emm)

$emmeans
 RACACOR  emmean     SE      df lower.CL upper.CL
 Branca    71.41 0.0207 1466479    71.37    71.45
 Preta     65.61 0.0496 1466479    65.52    65.71
 Amarela   74.40 0.1850 1466479    74.03    74.76
 Parda     64.46 0.0235 1466479    64.42    64.51
 Indígena  64.41 0.2700 1466479    63.88    64.94

Confidence level used: 0.95 

$contrasts
 contrast           estimate     SE      df t.ratio p.value
 Branca - Preta       5.8006 0.0537 1466479 108.035  <.0001
 Branca - Amarela    -2.9842 0.1860 1466479 -16.040  <.0001
 Branca - Parda       6.9498 0.0313 1466479 222.034  <.0001
 Branca - Indígena    7.0012 0.2710 1466479  25.820  <.0001
 Preta - Amarela     -8.7848 0.1910 1466479 -45.893  <.0001
 Preta - Parda        1.1491 0.0549 1466479  20.950  <.0001
 Preta - Indígena     1.2006 0.2750 1466479   4.368  <.0001
 Amarela - Parda      9.9339 0.1860 1466479  53.299  <.0001
 Amarela - Indígena   9.9854 0.3280 1466479  30.486  <.0001
 Parda - Indígena     0.0515 0.2710 1466479   0.190  0.8496

P value adjustment: holm method for 10 tests

print(plot(emm$emmeans, 
           colors="black"))

print(plot(emm$contrasts, 
           colors="black"))

print(multcomp::cld(object=emm$emmeans,
                    level=1-alfa,
                    adjust="holm",
                    Letters=letters,
                    alpha=alfa))

 RACACOR  emmean     SE      df lower.CL upper.CL .group
 Indígena  64.41 0.2700 1466479    63.72    65.11  a    
 Parda     64.46 0.0235 1466479    64.40    64.52  a    
 Preta     65.61 0.0496 1466479    65.48    65.74   b   
 Branca    71.41 0.0207 1466479    71.36    71.47    c  
 Amarela   74.40 0.1850 1466479    73.92    74.87     d 

Confidence level used: 0.95 
Conf-level adjustment: bonferroni method for 5 estimates 
P value adjustment: holm method for 10 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.

Quando o arquivo de dados é enorme, as estatisticas de significancia estatística (valor p) e de significancia pratica (eta^2) tornam-se desprovidas de significado para tomada de decisão sobre a hipotese nula. Algoritmos de aprendizagem de maquina podem ser uma alternativa?

Sim, é um ponto válido. Quando você tem um conjunto de dados muito grande, até pequenas diferenças podem ser estatisticamente significativas, mesmo que não sejam praticamente significativas. Por exemplo, uma diferença muito pequena entre dois grupos pode ter um valor-p muito baixo (indicando significância estatística) em um conjunto de dados muito grande, mas essa diferença pode não ser relevante ou útil na prática.

A significância prática, frequentemente representada por medidas como o tamanho do efeito (por exemplo, \(\eta^2\) para ANOVA), é crucial para entender se uma descoberta é apenas estatisticamente significativa ou também relevante na prática.

Prática R: ANOVA Unifatorial

MSP4112: Estatística Aplicada à Medicina

José Siqueira (silveira@usp.br)
Paulo Silveira (silveira@usp.br)

19 agosto 2025 18:44h

Material

Introdução

ANOVA unifatorial independente de Welch

ANOVA unifatorial independente de Welch sem dados brutos

Análise de poder retrospectivo (plug in) de ANOVA unifatorial independente balanceada de Fisher (Gerard et al., 1998)

ANOVA unifatorial relacionada

pulse

resp

com uma terceira condição de tratamento

quatro drogas sem medidas repetidas

Pacientes diferentes?

Questão

Cinco condições independentes: tamanho de amostra muito grande

ANOVA de Fisher: idade de óbito e raça/cor

Questões

Prática R: ANOVA Unifatorial

MSP4112: Estatística Aplicada à Medicina

José Siqueira (silveira@usp.br) Paulo Silveira (silveira@usp.br)

19 agosto 2025 18:44h

Material

Introdução

ANOVA unifatorial independente de Welch

ANOVA unifatorial independente de Welch sem dados brutos

Análise de poder retrospectivo (plug in) de ANOVA unifatorial independente balanceada de Fisher (Gerard et al., 1998)

ANOVA unifatorial relacionada

pulse

resp

com uma terceira condição de tratamento

quatro drogas sem medidas repetidas

Pacientes diferentes?

Questão

Cinco condições independentes: tamanho de amostra muito grande

ANOVA de Fisher: idade de óbito e raça/cor

Questões

José Siqueira (silveira@usp.br)
Paulo Silveira (silveira@usp.br)