1. Montar Tabela de hepatectomia parcial (benigna)

hepatectomia_benigna <- data.frame(
  regiao = c("Norte", "Nordeste", "Sudeste", "Sul", "Centro-Oeste"),
  AIH_aprovadas = c(405, 1191, 2604, 1403, 409),
  obitos = c(41, 98, 212, 88, 45)
)

# Calcular percentual de óbitos
hepatectomia_benigna$percentual_obitos <- round(100 * hepatectomia_benigna$obitos / hepatectomia_benigna$AIH_aprovadas, 2)

hepatectomia_benigna
##         regiao AIH_aprovadas obitos percentual_obitos
## 1        Norte           405     41             10.12
## 2     Nordeste          1191     98              8.23
## 3      Sudeste          2604    212              8.14
## 4          Sul          1403     88              6.27
## 5 Centro-Oeste           409     45             11.00

2. Montar Tabela de hepatectomia parcial em oncologia (maligna)

hepatectomia_maligna <- data.frame(
  regiao = c("Norte", "Nordeste", "Sudeste", "Sul", "Centro-Oeste"),
  AIH_aprovadas = c(420, 1777, 6415, 3066, 690),
  obitos = c(50, 137, 417, 183, 53)
)

# Calcular percentual de óbitos
hepatectomia_maligna$percentual_obitos <- round(100 * hepatectomia_maligna$obitos / hepatectomia_maligna$AIH_aprovadas, 2)

hepatectomia_maligna
##         regiao AIH_aprovadas obitos percentual_obitos
## 1        Norte           420     50             11.90
## 2     Nordeste          1777    137              7.71
## 3      Sudeste          6415    417              6.50
## 4          Sul          3066    183              5.97
## 5 Centro-Oeste           690     53              7.68

3. Juntar os dados para comparação

dados_comparados <- hepatectomia_benigna %>%
  rename(AIH_benigna = AIH_aprovadas, obitos_benigna = obitos, perc_benigna = percentual_obitos) %>%
  inner_join(
    hepatectomia_maligna %>%
      rename(AIH_maligna = AIH_aprovadas, obitos_maligna = obitos, perc_maligna = percentual_obitos),
    by = "regiao"
  )

dados_comparados
##         regiao AIH_benigna obitos_benigna perc_benigna AIH_maligna
## 1        Norte         405             41        10.12         420
## 2     Nordeste        1191             98         8.23        1777
## 3      Sudeste        2604            212         8.14        6415
## 4          Sul        1403             88         6.27        3066
## 5 Centro-Oeste         409             45        11.00         690
##   obitos_maligna perc_maligna
## 1             50        11.90
## 2            137         7.71
## 3            417         6.50
## 4            183         5.97
## 5             53         7.68

4. Testes de proporções por região

# Norte - p-valor = 0.4806. Não rejeita H0. Evidência sugere que as proporções são iguais
prop.test(x = c(41, 50), n = c(405, 420))
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(41, 50) out of c(405, 420)
## X-squared = 0.49749, df = 1, p-value = 0.4806
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.06292574  0.02729964
## sample estimates:
##    prop 1    prop 2 
## 0.1012346 0.1190476
# Nordeste - p-valor = 0.6575. Não rejeita H0. Evidência sugere que as proporções são iguais
prop.test(x = c(98, 137), n = c(1191, 1777))
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(98, 137) out of c(1191, 1777)
## X-squared = 0.19686, df = 1, p-value = 0.6573
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.01544790  0.02582303
## sample estimates:
##     prop 1     prop 2 
## 0.08228380 0.07709623
# Sudeste - p-valor = 0.006392. Rejeita H0. Evidência sugere que as proporções são diferentes.
prop.test(x = c(212, 417), n = c(2604, 6415))
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(212, 417) out of c(2604, 6415)
## X-squared = 7.4363, df = 1, p-value = 0.006392
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.00402657 0.02879206
## sample estimates:
##     prop 1     prop 2 
## 0.08141321 0.06500390
# Sul - p-value = 0.7436. Não rejeita H0. Evidência sugere que as proporções são iguais
prop.test(x = c(88, 183), n = c(1403, 3066))
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(88, 183) out of c(1403, 3066)
## X-squared = 0.107, df = 1, p-value = 0.7436
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.01269165  0.01876335
## sample estimates:
##     prop 1     prop 2 
## 0.06272274 0.05968689
# Centro-Oeste - p-value = 0.07874. Rejeita H0 somente a 10%. Mas 0 está contido no intervalor de confiança. Resultado é analisar com cautela.
prop.test(x = c(45, 53), n = c(409, 690))
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(45, 53) out of c(409, 690)
## X-squared = 3.0906, df = 1, p-value = 0.07874
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.004989972  0.071415683
## sample estimates:
##     prop 1     prop 2 
## 0.11002445 0.07681159

5. Ajustar as bases para estimar o modelo GLM

library(dplyr)


benigna <- hepatectomia_benigna %>%
  mutate(tipo = "Benigna",
         vivos = AIH_aprovadas - obitos)

maligna <- hepatectomia_maligna %>%
  mutate(tipo = "Maligna",
         vivos = AIH_aprovadas - obitos)

# Unir as duas bases
hepatectomia_total <- bind_rows(benigna, maligna)

# Verificar a base final
print(hepatectomia_total)
##          regiao AIH_aprovadas obitos percentual_obitos    tipo vivos
## 1         Norte           405     41             10.12 Benigna   364
## 2      Nordeste          1191     98              8.23 Benigna  1093
## 3       Sudeste          2604    212              8.14 Benigna  2392
## 4           Sul          1403     88              6.27 Benigna  1315
## 5  Centro-Oeste           409     45             11.00 Benigna   364
## 6         Norte           420     50             11.90 Maligna   370
## 7      Nordeste          1777    137              7.71 Maligna  1640
## 8       Sudeste          6415    417              6.50 Maligna  5998
## 9           Sul          3066    183              5.97 Maligna  2883
## 10 Centro-Oeste           690     53              7.68 Maligna   637

6. Rodar modelo GLM com regiao como preditor (sem separar por tipo)

#O problema desse modelo é que ele não compara as regiões entre si. Todas as regiões são comparadas com Centro Oeste.
modelo_geral <- glm(cbind(obitos, vivos) ~ regiao,
                    data = hepatectomia_total,
                    family = binomial)

summary(modelo_geral)
## 
## Call:
## glm(formula = cbind(obitos, vivos) ~ regiao, family = binomial, 
##     data = hepatectomia_total)
## 
## Coefficients:
##                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)     -2.3238     0.1058 -21.955  < 2e-16 ***
## regiaoNordeste  -0.1298     0.1258  -1.032  0.30221    
## regiaoNorte      0.2361     0.1535   1.539  0.12390    
## regiaoSudeste   -0.2669     0.1136  -2.349  0.01884 *  
## regiaoSul       -0.4165     0.1230  -3.386  0.00071 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 44.261  on 9  degrees of freedom
## Residual deviance: 11.967  on 5  degrees of freedom
## AIC: 85.477
## 
## Number of Fisher Scoring iterations: 4
exp(coef(modelo_geral))  # Odds ratios
##    (Intercept) regiaoNordeste    regiaoNorte  regiaoSudeste      regiaoSul 
##      0.0979021      0.8782866      1.2663488      0.7657671      0.6593786

7. Fazer a comparação entre todas as regiões

# Transformar a variável 'regiao' em fator
hepatectomia_total$regiao <- as.factor(hepatectomia_total$regiao)

# Rodar o modelo novamente
modelo_regiao <- glm(cbind(obitos, vivos) ~ regiao,
                     data = hepatectomia_total,
                     family = binomial)

# Comparações múltiplas entre todos os pares de regiões
library(multcomp)
## Carregando pacotes exigidos: mvtnorm
## Carregando pacotes exigidos: survival
## Carregando pacotes exigidos: TH.data
## Carregando pacotes exigidos: MASS
## 
## Anexando pacote: 'MASS'
## O seguinte objeto é mascarado por 'package:dplyr':
## 
##     select
## 
## Anexando pacote: 'TH.data'
## O seguinte objeto é mascarado por 'package:MASS':
## 
##     geyser
comparacoes <- glht(modelo_regiao, linfct = mcp(regiao = "Tukey"))

# Resultados
summary(comparacoes) # Resumo do modelo
## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## 
## Fit: glm(formula = cbind(obitos, vivos) ~ regiao, family = binomial, 
##     data = hepatectomia_total)
## 
## Linear Hypotheses:
##                              Estimate Std. Error z value Pr(>|z|)    
## Nordeste - Centro-Oeste == 0 -0.12978    0.12579  -1.032  0.83161    
## Norte - Centro-Oeste == 0     0.23614    0.15347   1.539  0.52087    
## Sudeste - Centro-Oeste == 0  -0.26688    0.11363  -2.349  0.12147    
## Sul - Centro-Oeste == 0      -0.41646    0.12301  -3.386  0.00578 ** 
## Norte - Nordeste == 0         0.36592    0.13028   2.809  0.03678 *  
## Sudeste - Nordeste == 0      -0.13709    0.07956  -1.723  0.40330    
## Sul - Nordeste == 0          -0.28668    0.09246  -3.100  0.01527 *  
## Sudeste - Norte == 0         -0.50302    0.11858  -4.242  < 0.001 ***
## Sul - Norte == 0             -0.65260    0.12759  -5.115  < 0.001 ***
## Sul - Sudeste == 0           -0.14958    0.07508  -1.992  0.25607    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
confint(comparacoes)  # Intervalos de confiança
## 
##   Simultaneous Confidence Intervals
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## 
## Fit: glm(formula = cbind(obitos, vivos) ~ regiao, family = binomial, 
##     data = hepatectomia_total)
## 
## Quantile = 2.701
## 95% family-wise confidence level
##  
## 
## Linear Hypotheses:
##                              Estimate lwr      upr     
## Nordeste - Centro-Oeste == 0 -0.12978 -0.46956  0.20999
## Norte - Centro-Oeste == 0     0.23614 -0.17840  0.65068
## Sudeste - Centro-Oeste == 0  -0.26688 -0.57380  0.04004
## Sul - Centro-Oeste == 0      -0.41646 -0.74871 -0.08420
## Norte - Nordeste == 0         0.36592  0.01403  0.71781
## Sudeste - Nordeste == 0      -0.13709 -0.35200  0.07781
## Sul - Nordeste == 0          -0.28668 -0.53642 -0.03693
## Sudeste - Norte == 0         -0.50302 -0.82329 -0.18274
## Sul - Norte == 0             -0.65260 -0.99722 -0.30797
## Sul - Sudeste == 0           -0.14958 -0.35238  0.05322
exp(coef(comparacoes)) # OR
## Nordeste - Centro-Oeste    Norte - Centro-Oeste  Sudeste - Centro-Oeste 
##               0.8782866               1.2663488               0.7657671 
##      Sul - Centro-Oeste        Norte - Nordeste      Sudeste - Nordeste 
##               0.6593786               1.4418401               0.8718875 
##          Sul - Nordeste         Sudeste - Norte             Sul - Norte 
##               0.7507557               0.6047047               0.5206927 
##           Sul - Sudeste 
##               0.8610694