Anàlisi de dies d’espera segons especialitat, prioritat i sexe

L’objectiu del present exercisi comentat és basant-nos en un exemple d’exportacio procedent de BO-SAP, realitzar una análisi descriptiva, gràfica i inferencial sencilles per tal de respondre les següents preguntes:

Resumeix el numero de visites i dies d’espera per especialitat.
Existeixen diferencies significatives entre les diferents especialitats en termes de dies d’espera?
El nivell de prioritat assignat te repercusió significativa en els dies d’espera reals?
Existeix un bias de gènere en l’espera fins a la visita al especialista?
Hi ha relació entre l’edat i els dies d’espera?

# carreguem llibreries
library(readxl)
library(dplyr)
library(ggplot2)
library(janitor)
library(gt)
library(ggmosaic)

1 Lectura de dades i neteja

df <- read_excel("visites_prioritat.xlsx")
glimpse(df)

Rows: 300
Columns: 8
$ `ID Pacient`          <chr> "PAC0000", "PAC0001", "PAC0002", "PAC0003", "PAC…
$ `Especialitat mèdica` <chr> "Oftalmologia", "Traumatologia", "Pediatria", "O…
$ `Data Visita`         <dttm> 2024-02-29, 2024-02-15, 2024-03-05, 2024-02-25,…
$ `Data Sol·licitud`    <dttm> 2024-01-20, 2024-01-08, 2024-01-15, 2024-01-13,…
$ `Dies d'espera`       <dbl> 40, 38, 50, 43, 9, 20, 16, 46, 36, 39, 22, 46, 3…
$ Sexe                  <chr> "Dona", "Home", "Dona", "Dona", "Dona", "Home", …
$ Edat                  <dbl> 72, 37, 40, 5, 66, 55, 19, 18, 34, 64, 44, 35, 3…
$ `Nivell de prioritat` <chr> "Normal", "Baixa", "Normal", "Normal", "Alta", "…

# Netejem aquests noms horrorosos i recodifiquem com a factor si és necessari.

df<- df %>% janitor::clean_names()
df$nivell_de_prioritat <- factor(df$nivell_de_prioritat, levels = c("Alta", "Normal", "Baixa"))

glimpse(df)

Rows: 300
Columns: 8
$ id_pacient          <chr> "PAC0000", "PAC0001", "PAC0002", "PAC0003", "PAC00…
$ especialitat_medica <chr> "Oftalmologia", "Traumatologia", "Pediatria", "Oft…
$ data_visita         <dttm> 2024-02-29, 2024-02-15, 2024-03-05, 2024-02-25, 2…
$ data_sol_licitud    <dttm> 2024-01-20, 2024-01-08, 2024-01-15, 2024-01-13, 2…
$ dies_despera        <dbl> 40, 38, 50, 43, 9, 20, 16, 46, 36, 39, 22, 46, 35,…
$ sexe                <chr> "Dona", "Home", "Dona", "Dona", "Dona", "Home", "H…
$ edat                <dbl> 72, 37, 40, 5, 66, 55, 19, 18, 34, 64, 44, 35, 31,…
$ nivell_de_prioritat <fct> Normal, Baixa, Normal, Normal, Alta, Alta, Alta, B…

head(df)

# A tibble: 6 × 8
  id_pacient especialitat_medica data_visita         data_sol_licitud   
  <chr>      <chr>               <dttm>              <dttm>             
1 PAC0000    Oftalmologia        2024-02-29 00:00:00 2024-01-20 00:00:00
2 PAC0001    Traumatologia       2024-02-15 00:00:00 2024-01-08 00:00:00
3 PAC0002    Pediatria           2024-03-05 00:00:00 2024-01-15 00:00:00
4 PAC0003    Oftalmologia        2024-02-25 00:00:00 2024-01-13 00:00:00
5 PAC0004    Oftalmologia        2024-02-09 00:00:00 2024-01-31 00:00:00
6 PAC0005    Cardiologia         2024-01-30 00:00:00 2024-01-10 00:00:00
# ℹ 4 more variables: dies_despera <dbl>, sexe <chr>, edat <dbl>,
#   nivell_de_prioritat <fct>

2 Resum per especialitat

df_resum <- df %>%
  filter(!is.na(especialitat_medica)) %>%
  group_by(especialitat_medica) %>%
  summarise(
    visites = n(),
    espera_mitja = round(mean(dies_despera, na.rm = TRUE),2)
  ) %>%
  arrange(desc(visites))

#formategem en una taula elegant
df_resum %>% gt()

especialitat_medica	visites	espera_mitja
Oftalmologia	56	35.73
Cardiologia	54	35.43
Reumatologia	53	32.25
Pediatria	51	31.69
Dermatologia	43	35.58
Traumatologia	43	29.81

3 Gràfic de barres (visites)

ggplot(df_resum, aes(x = reorder(especialitat_medica, visites), y = visites)) +
  geom_col() +
  coord_flip() +
  labs(
    title = "Nombre de visites per especialitat",
    x = "Especialitat",
    y = "Nombre de visites"
  )

4 Comparació de dies d’espera per especialitat

#graficament
ggplot(df, aes(x = especialitat_medica, y = dies_despera)) +
  geom_boxplot() +
  labs(
    title = "Dies d'espera per especialitat",
    x = "Especialitat",
    y = "Dies d'espera"
  ) +
  theme_minimal()

#numericament
anova_especialitat <- aov(dies_despera ~ especialitat_medica, data = df)
summary(anova_especialitat)

                     Df Sum Sq Mean Sq F value Pr(>F)  
especialitat_medica   5   1501   300.2   2.501 0.0308 *
Residuals           294  35288   120.0                 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

TukeyHSD(anova_especialitat)

  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = dies_despera ~ especialitat_medica, data = df)

$especialitat_medica
                                 diff        lwr       upr     p adj
Dermatologia-Cardiologia    0.1554694  -6.267955 6.5788940 0.9999998
Oftalmologia-Cardiologia    0.3062169  -5.687791 6.3002250 0.9999905
Pediatria-Cardiologia      -3.7396514  -9.876205 2.3969025 0.5011070
Reumatologia-Cardiologia   -3.1806429  -9.257359 2.8960732 0.6636105
Traumatologia-Cardiologia  -5.6119724 -12.035397 0.8114522 0.1255884
Oftalmologia-Dermatologia   0.1507475  -6.221626 6.5231211 0.9999998
Pediatria-Dermatologia     -3.8951208 -10.401756 2.6115147 0.5214518
Reumatologia-Dermatologia  -3.3361123  -9.786344 3.1141195 0.6750065
Traumatologia-Dermatologia -5.7674419 -12.545303 1.0104189 0.1457423
Pediatria-Oftalmologia     -4.0458683 -10.128964 2.0372275 0.3993331
Reumatologia-Oftalmologia  -3.4868598  -9.509587 2.5358671 0.5588881
Traumatologia-Oftalmologia -5.9181894 -12.290563 0.4541842 0.0857145
Reumatologia-Pediatria      0.5590085  -5.605600 6.7236173 0.9998374
Traumatologia-Pediatria    -1.8723210  -8.378957 4.6343145 0.9626544
Traumatologia-Reumatologia -2.4313295  -8.881561 4.0189023 0.8885813

Tot i que Trauma i Oft son bastant diferents a simple vista, no hi han diferències significatives.

5 Gràfic mosaic de sexe per especialitat

ggplot(data = df) +
  geom_mosaic(aes(x = product(sexe), fill = especialitat_medica)) +
  labs(title = "Distribució de sexe per especialitat")

6 Contrast d’hipòtesi: diferència de dies d’espera entre sexes

#graficament
ggplot(df, aes(x = sexe, y = dies_despera)) +
  geom_boxplot(fill = "gold") +
  labs(
    title = "Diferències de dies d'espera i gènere",
    x = "Genere",
    y = "Dies d'espera")

# analiticament t.test
t.test(dies_despera ~ sexe, data = df)


    Welch Two Sample t-test

data:  dies_despera by sexe
t = -1.3226, df = 294.71, p-value = 0.187
alternative hypothesis: true difference in means between group Dona and group Home is not equal to 0
95 percent confidence interval:
 -4.2172291  0.8272402
sample estimates:
mean in group Dona mean in group Home 
          32.62759           34.32258

7 Analisis dels dies d’espera per nivell de prioritat

df %>% group_by(nivell_de_prioritat) %>% 
  summarise(mitjana_dies = mean(dies_despera, na.rm=T))

# A tibble: 3 × 2
  nivell_de_prioritat mitjana_dies
  <fct>                      <dbl>
1 Alta                        19  
2 Normal                      39.6
3 Baixa                       41.7

ggplot(df, aes(x = nivell_de_prioritat, y = dies_despera, fill = nivell_de_prioritat)) +
  geom_violin(trim = FALSE, alpha = 0.6) +
  geom_boxplot(width = 0.1, fill = "white") +
  labs(
    title = "Dies d'espera segons el nivell de prioritat",
    x = "Prioritat",
    y = "Dies d'espera"
  ) +
  theme_minimal()

8 ANOVA: diferències en dies d’espera segons la prioritat i anàlisi post hoc.

anova_model <- aov(dies_despera ~ nivell_de_prioritat, data = df)
summary(anova_model)

                     Df Sum Sq Mean Sq F value Pr(>F)    
nivell_de_prioritat   2  29441   14720     595 <2e-16 ***
Residuals           297   7348      25                   
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

tukey_result <- TukeyHSD(anova_model)
tukey_result

  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = dies_despera ~ nivell_de_prioritat, data = df)

$nivell_de_prioritat
                 diff        lwr       upr     p adj
Normal-Alta  20.57042 19.0174486 22.123396 0.0000000
Baixa-Alta   22.69841 20.7947393 24.602086 0.0000000
Baixa-Normal  2.12799  0.3543788  3.901601 0.0138968

Aqui veiem diferencies significatives entre tots els nivells de prioritat. Tot i que Baixa i Normal tenen una diferència petita de 2 dies, el test diu que és significativa. Ho és també en termes pràctics? Reflexionar en entre significació estadística i importància clínica o operativa.

9 Relació entre edat i dies d’espera

ggplot(df, aes(x = edat, y = dies_despera)) +
  geom_point(alpha = 0.4, size = 2, color = "darkblue") +
  geom_smooth(method = "lm", se = TRUE, color = "red") +
  labs(
    title = "Relació entre edat i dies d'espera",
    x = "Edat (anys)",
    y = "Dies d'espera"
  ) +
  theme_minimal()

cor.test(df$edat, df$dies_despera, method = "pearson")


    Pearson's product-moment correlation

data:  df$edat and df$dies_despera
t = 0.066407, df = 298, p-value = 0.9471
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.1094417  0.1170367
sample estimates:
        cor 
0.003846854

model <- lm(dies_despera ~ edat, data = df)
summary(model)


Call:
lm(formula = dies_despera ~ edat, data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-25.473 -10.487   3.490   8.491  21.596 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 33.382174   1.933977  17.261   <2e-16 ***
edat         0.002451   0.036916   0.066    0.947    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 11.11 on 298 degrees of freedom
Multiple R-squared:  1.48e-05,  Adjusted R-squared:  -0.003341 
F-statistic: 0.00441 on 1 and 298 DF,  p-value: 0.9471

No hi ha correlació i el model de regresió lineal mostra un coeficient d’edat ≈ 0 → no hi ha pràcticament cap relació entre edat i dies d’espera. p-valor = 0.947 → no significativa. No podem rebutjar la hipòtesi nula: l’edat no té cap efecte sobre els dies d’espera.