1 Introducción

La presión arterial es la fuerza que ejerce contra la pared arterial la sangre que circula por las arterias. La presión o tensión arterial incluye dos mediciones: la presión sistólica, que se mide durante el latido del corazón (momento de presión máxima), y la presión diastólica, que se mide durante el descanso entre dos latidos (momento de presión mínima). En el resultado, primero se presenta la presión sistólica y luego la diastólica, por ejemplo: 120/80.

Consideremos el caso en donde se estudia la reducción de la presión arterial diastólica (DBP, por sus siglas en inglés) en 40 pacientes. Fue realizada una medición al inicio del estudio (DBP_Inicial) en donde el paciente se colocaba en posición supina (acostado boca arriba) para tomar su presión diastólica. Después de la medición inicial, se realizó la asignación aleatoria de dos tratamientos objeto de estudio (A y B), a 20 participantes se les asignó el tratamiento A y los restantes 20 el tratamiento B (ensayo clínico con brazos de igual tamaño). A los cuatro meses después de haber recibido los tratamientos, en donde hubo supervisión continua de los participantes por parte de los investigadores, se realizaron nuevamente mediciones de DBP (DBP_Final).

2 Objetivo

Investigar si el tratamiento A (nuevo medicamento) tiene mayor eficiencia en reducir la DPB comparado con el tratamiento B (medicamento estándar o control).

3 Materiales

3.1 Instalación y activación de librerías

# Instalación
list.of.packages <- c("tidyverse", "dplyr","forcats", "ggplot2","car", "ggpubr","datarium", "rstatix", "nortest", "lmtest", "knitr")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)) install.packages(new.packages)

# Activación
library(tidyverse)
library(dplyr)
library(forcats)
library(ggplot2)
library(car)
library(ggpubr)
library(datarium)
library(rstatix)
library(nortest)
library(lmtest)
library(knitr )

3.2 Datos DBP

Fueron incluidos 40 pacientes hipertensos. La DBP (en mmHg) fue medida en los 40 sujetos al comienzo del estudio (\(y_1\)). Posteriormente, se realizó la asignación aleatoria de los tratamientos, los brazos fueron del mismo tamaño, es decir, a 20 pacientes se les asigno el tratamiento A y a los restantes 20 el tratamiento B. A los cuatro meses de administrados los tratamientos fue realizada una medición final (\(y_2\)).

A continuación son ingresados y presentados los datos del ensayo clínico aleatorizado (ECA) para DBP.

Sujeto <- c(1:40)
Tratamiento <- c("A","A","A","A","A","A","A","A","A","A","A","A","A","A","A","A","A","A","A","A","B","B","B","B","B","B","B","B","B","B","B","B","B","B","B","B","B","B","B","B")
DBP_Inicial <- c(114,116,119,115,116,117,118,120,114,115,117,116,119,118,115,114,117,120,114,117,114,116,114,114,116,114,119,118,114,120,117,118,121,116,118,119,116,116,117,118)
DBP_Final <- c(105,101,98,101,105,102,99,102,103,97,101,102,104,99,102,100,102,103,100,101,113,110,109,115,109,110,115,112,108,113,115,110,115,111,112,111,109,112,115,115)
di <- DBP_Inicial-DBP_Final

data <- tibble(Sujeto,Tratamiento,DBP_Inicial,DBP_Final,di)
data
## # A tibble: 40 × 5
##    Sujeto Tratamiento DBP_Inicial DBP_Final    di
##     <int> <chr>             <dbl>     <dbl> <dbl>
##  1      1 A                   114       105     9
##  2      2 A                   116       101    15
##  3      3 A                   119        98    21
##  4      4 A                   115       101    14
##  5      5 A                   116       105    11
##  6      6 A                   117       102    15
##  7      7 A                   118        99    19
##  8      8 A                   120       102    18
##  9      9 A                   114       103    11
## 10     10 A                   115        97    18
## # ℹ 30 more rows

4 Métodos

La prueba t de student es usada para comparar la diferencia media entre los dos grupos generados por los tratamientos (A y B), se requiere que la variable respuesta sea continua y con distribución aproximadamente normal. Utilizamos la prueba clásica t de student para datos pareados.

Estrategia de análisis:

Calcular la diferencia en el tratamiento A:

\[d_{Ai} = y_{1Ai}-y_{2Ai}, i = 1,2,...,20\] Dónde:

  • \(y_{1Ai}\) denota la medición DBP en la i-ésima persona al inicio del estudio y antes de consumir el medicamento A.

  • \(y_{2Ai}\) denota la medición DBP en la i-ésima persona después de los cuatro meses y tras haber consumido el tratamiento A.

  • \(d_{Ai}\) denota la diferencia o cambio (reducción) entre \(y_{1Ai}\) y \(y_{2Ai}\).

y la , diferencia en el tratamiento B:

\[d_{Bi} = y_{1Bi}-y_{2Bi}, i = 21,22,...,40\] Dónde:

  • \(y_{1Bi}\) denota la medición DBP en la i-ésima persona al inicio del estudio y antes de consumir el medicamento B.

  • \(y_{2Bi}\) denota la medición DBP en la i-ésima persona después de los cuatro meses y tras haber consumido el tratamiento B.

  • \(d_{Bi}\) denota la diferencia o cambio (reducción) entre \(y_{1Bi}\) y \(y_{2Bi}\).

Se utilizan \(d_{Ai}\) y \(d_{Bi}\) para comparar las medias de los brazos. En otras palabras, con la prueba t se compara el promedio de las diferencias o reducción de DBP para los dos tratamientos en cuestión.

5 Análisis estadístico

5.1 Hipótesis estadística

En primer, construimos el contraste de hipótesis como sigue:

\(H_0: \text{δ}=\mu_A - \mu_B = 0\) Vs \(H_1: \text{δ}=\mu_A - \mu_B > 0\)

\(\text{Equivalente a:}\)

\(H_0: \mu_A = \mu_B\) Vs \(H_1: \mu_A > \mu_B\)

Dónde:

  • \(\mu_A\): es la media de la diferencia en DBP bajo el tratamiento A
  • \(\mu_B\): es la media de la diferencia en DBP bajo el tratamiento B
De modo que, al NO RECHAZAR la hipótesis nula (\(H_0\)), la prueba indicaría que no hay evidencia suficiente para determinar que el tratamiento A genera una mayor reducción de la DBP comparado con el tratamiento B.

5.2 Estadísticos de centralidad y dispersión

A continuación presentamos el análisis exploratorio usando la media, varianza y desviación estándar.

ta <- data %>%
  group_by(Tratamiento) %>%
  summarise(prom_DBP_Inicial = mean(DBP_Inicial),
            prom_DBP_Final   = mean(DBP_Final),
            std_DBP_Inicial  = sd(DBP_Inicial),
            std_DBP_Final    = sd(DBP_Final),
            var_DBP_Inicial  = var(DBP_Inicial),
            var_DBP_Final    = var(DBP_Final),
            prom_Diferencia  = mean(DBP_Inicial-DBP_Final),
            std_Diferencia   = sd(DBP_Final-DBP_Inicial),
            var_Diferencia   = var(DBP_Final-DBP_Inicial),) 

ta
## # A tibble: 2 × 10
##   Tratamiento prom_DBP_Inicial prom_DBP_Final std_DBP_Inicial std_DBP_Final
##   <chr>                  <dbl>          <dbl>           <dbl>         <dbl>
## 1 A                       117.           101.            1.99          2.13
## 2 B                       117.           112.            2.12          2.44
## # ℹ 5 more variables: var_DBP_Inicial <dbl>, var_DBP_Final <dbl>,
## #   prom_Diferencia <dbl>, std_Diferencia <dbl>, var_Diferencia <dbl>

Se aprecia que el tratamiento A reduce en alrededor de 10 unidades la DBP en comparación con el tratamiento B; es decir, se tienen indicios de reducción DBP con el nuevo tratamiento. Se observa que las varianzas en ambos tratamientos al inicio del ensayo son similares y a los 4 meses el tratamiento B presenta una varianza mayor.

5.3 Visualización

Se exploran las distribuciones de DBP_Inicial y DBP_Final, para cada uno de los dos grupos de tratamiento A y B. 
Datos1 <- data %>%
  tidyr::gather(key = "Tiempo", value = "Mediciones", DBP_Inicial, DBP_Final)
# Distribución de mediciones 
ggplot(
  data = Datos1,
  aes(x = factor(Tiempo, level = c("DBP_Inicial", "DBP_Final")), y = Mediciones, fill=Tratamiento)) +
  geom_boxplot(stat = "boxplot",
  position = "dodge2") +
  scale_x_discrete(name = "") +
  facet_wrap(~Tratamiento)

Las dos muestras son semejantes al inicio, ya que las distribuciones, alrededor de la media, de las DBP inicial son casi iguales, con esto observamos que tenemos dos muestras comparables.

Después de la aplicación de los tratamientos, notamos que el grupo A tiene un valor menor de centralidad y dispersión en comparación con el grupo B, lo cual indicaría que el tratamiento otorgado al grupo A tiene mayor efecto en disminuir la DBP comparado con el tratamiento A.

6 Revisión de supuestos

6.1 Normalidad

Se utiliza la prueba estadística de Shapiro-Wilk (\(H_0\): Hay normalidad) para la diferencia promedio de DBP en cada uno de los grupos. 


Normalidad <- data %>%
  group_by(Tratamiento) %>%
  shapiro_test(di)                                    # Prueba de Shapiro-Wilk
# Tabla  de bondad de ajuste
cumple <- NULL          # Variable  para cumplimiento o no cumplimiento
for(i in 1:2){
if(Normalidad[i,4] > 0.05){ cumple1 = "Cumple normalidad"      # Se cumple si valor p > 0.05
}else
{cumple1 = "No cumple normalidad"}
 cumple[i] = cumple1 }
Normalidad$Interpretación <- cumple
Normalidad[]        #Visualizamos la tabla
## # A tibble: 2 × 5
##   Tratamiento variable statistic     p Interpretación   
##   <chr>       <chr>        <dbl> <dbl> <chr>            
## 1 A           di           0.974 0.838 Cumple normalidad
## 2 B           di           0.931 0.163 Cumple normalidad
Observamos que para ambos grupos el valor p es mayor a 0.05 por lo cual no se rechaza \(H_0\) y, se concluye que existe compatibilidad con una distribución normal. En caso de no normalidad explorar una prueba no paramétrica.


6.2 Homocedasticidad

Existen diferentes métodos para evaluar la homocedasticidad entre grupos: Prueba F, Bartlett o Levene, donde todas estudian si existe una igualdad de varianzas entre dos o más grupos de estudio.

Se utiliza la prueba estadística de Bartlett (\(H_0: \sigma^2_A = \sigma^2_B\)) para estudiar la igualdad de varianzas en los grupos de tratamientos (A y B), siempre y cuando los datos provengan de una distribución normal, y dado que en el apartado anterior se validó el supuesto, utilizaremos la función bartlett.test para realizar la prueba con R. 


# Homocedasticidad con prueba de Bartlett

bartlett.test((data$DBP_Final-data$DBP_Inicial) ~ data$Tratamiento)    # Prueba de Barlett
## 
##  Bartlett test of homogeneity of variances
## 
## data:  (data$DBP_Final - data$DBP_Inicial) by data$Tratamiento
## Bartlett's K-squared = 0.76462, df = 1, p-value = 0.3819

Dado que el valor p es mayor a 0.05, no se rechaza \(H_0\), es decir, que existe igualdad entre las varianzas para la diferencia de DBP del grupo A y B.

Cuando los datos NO provienen de una distribución normal, una alternativa para evualuar dicho supuesto es la prueba de Levene (\(H_0: \sigma^2_A = \sigma^2_B\)) que se considera más robusta.

Utilizando la función leveneTest, tenemos:

# Homocedasticidad con prueba de Levene

leveneTest(di ~ Tratamiento, data = Datos1)              # Prueba de Levene
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  1   0.585 0.4467
##       78
El resultado de la prueba arroja un valor p (Pr) mayor a 0.05, lo cual indica que homogeneidad de varianzas.

7 Error tipo I, error tipo II y potencia de la prueba

7.1 Error tipo I

El error tipo I es:

\[ \text{P(rechazar } H_0 | H_0 \text{ es verdadera})= P(\overline{X}_1-\overline{X}_2 > x^*| \mu_A-\mu_B \leq 0)\] Dónde:

  • \(x^*\) denota el cuantil que define la región de rechazo considerando como cierta \(H_0: \mu_A-\mu_B = 0\).
# Definición de variables que se utilizarán para calcular el cuantil de la zona de rechazo

X_A <- ta[1,8]  # Diferencia promedio en el tratamiento A
X_B <- ta[2,8]  # Diferencia promedio en el tratamiento B
S_A <- ta[1,10] # varianza de diferencia promedio en el tratamiento A
S_B <- ta[2,10] # varianza de diferencia promedio en el tratamiento B
sum1 <- tapply(data$di, data$Tratamiento,length) # Tamaño de ambos grupos
n_A <- sum1[1] # Tamaño de grupo A
n_B <- sum1[2] # Tamaño de grupo B

Una vez definidas las variables, se calcula la varianza conjunta estimada (\(S_p^2\)), dada por:

\[S_p^2 = \frac{(n_A-1)s_A^2 + (n_B-1)s_B^2}{(n_A + n_B - 2)}\] donde:

  • \(n_A\): denota el tamaño de grupo del tratamiento A
  • \(n_B\): denota el tamaño de grupo del tratamiento B
  • \(s_A^2\): denota la varianza de la diferencia en el tratamiento A
  • \(s_B^2\): denota la varianza de la diferencia en el tratamiento B
# Varianza conjunta estimada

Sp <- (((n_A-1)*S_A)+((n_B-1)*S_B))/(n_A+n_B-2)
Sp <- as.numeric((Sp))       # Resultado como valor numérico
(Sp <- round(Sp,2))          # Resultado redondeado a dos decimales
## [1] 7.33
Ahora, calculamos la zona de rechazo empleando la distribución t:

\[x^* = t_{0.05,20+20-2}*S_p*\sqrt{\frac{1}{n_A}+\frac{1}{n_B}}\]

# Calculamos el cuantil para una distribución t con una probabilidad de 0.05 y 38 grados de libertad

t <- qt(0.05, 38, lower.tail = F)          # Cuantil de distribución t
(t<-round(t,digits=2))                     # Redondeamos cuantil t a dos decimales
## [1] 1.69
# Calculamos la zona de rechazo

X_ <- t*sqrt(Sp)*sqrt((1/n_A)+(1/n_B))
X_<-round(X_,digits=2)                    # Resultado redondeado a dos decimales
(X_<- as.numeric(X_))                     # Resultado como valor numérico
## [1] 1.45
#Creamos gráfico para la prueba unilateral de distribución t de student

plot(x_1,y, type = "l", lwd = 2, axes = FALSE, xlab = "", ylab = "",main="Prueba unilateral distribución t_student")
axis(1, at = -3:3, labels = c("-3s", "-2s", "-1s", "Media", "1s", "2s", "3s"))
legend(x = 1.5, y = 0.3,   
       legend = c("Zona de rechazo H0","Zona de NO rechazo H0"),
       fill = c("red","white"),       
       cex = 0.8,
       border = "black")

arrows(x0 = X_,
       x1 = X_+0.5,
       y0 = dt(x = X_, df = 20),
       y1 = dt(x = X_, df = 20),
       length = 0.1) 

text(x = X_+1, y = dt(x = X_, df = 38),      # Coordenadas del texto
     label = "t = 1.45")

text(x = X_+1, y = dt(x = X_, df = 38)-0.08, # Coordenadas del texto
     label = "5%")

text(x = X_-1.5, y = dt(x = X_, df = 38)-0.08, # Coordenadas del texto
     label = "95%")
       
valor <- X_
valor2 <- seq(X_, 4, length= length(x_1[x_1 >= valor]))

polygon(c( x_1[x_1 >= valor],valor),
        c( dt(x = valor, df = 38), dt(x = valor2, df = 38)),
        col = "red",
        border = 1)

\[\text{Error tipo I = P}(\overline{X}_1-\overline{X}_2 > 1.45| \mu_A-\mu_B = 0)\]

# Error tipo 1
# La probabilidad que la variable t de Student de 38 grados de libertad deja a la derecha de 1.69

# P[t38 > 1.69]
# lower.tail = TRUE (por default)

t_3 <- pt(1.69, 38, lower.tail = FALSE)         # Probabilidad a la derecha de t=1.69
                                                
round(t_3,digits=2)                             # Resultado redondeado a dos decimales
## [1] 0.05
La probabilidad de cometer el error tipo I, es decir, la probabilidad de rechazar la hipótesis nula siendo ésta verdadera, es 5%.

7.2 Error tipo II

El error tipo II para el ejemplo que se ha venido desarrollando es:

\[ \text{P(No rechazar } H_0 | H_0 \text{ es falsa}) = P(\overline{X}_1-\overline{X}_2 ≤ x^*| \mu_A-\mu_B > 0)\]
Si empleamos la región de rechazo anterior, cuando la diferencia de media poblacional “suponemos” tiene un valor de 2 (δ = \(\mu_A-\mu_B=2)\) en el valor de DBP, entonces la probabilidad cometer error tipo II es:

\[\text{Error tipo II = P}(\overline{X}_1-\overline{X}_2 < 1.45| \mu_A-\mu_B = 2)\]

# Estadístico de prueba t student

t_1 <- (X_- 2)/(sqrt(Sp)*sqrt((1/n_A)+(1/n_B)))
t_1 <- as.numeric(t_1)                                # Resultado de estadístico como valor numérico
(t_1<-round(t_1,digits=2))                            # Resultado de estadístico redondeado a dos decimales
## [1] -0.64
# Error tipo II
# La probabilidad de que una variable t de Student de 38 grados de libertad deja a la izquierda de -0.65:

# P[t38 < -0.65]

t_3 <- pt(t_1, 38)                              # lower.tail = TRUE (por default)
round(t_3,digits=2)
## [1] 0.26
La probabilidad de cometer el error tipo II, es decir, la probabilidad de NO rechazar la hipótesis nula siendo ésta falsa ya que en realidad hay una diferencia de 2 unidades de DBP, es de 26%.

7.3 Potencia estadística

\[Poder_{prueba} = 1- ErrorTipoII\]
Cuando la diferencia de las medias de los tratamientos disminuye la probabilidad de cometer error tipo II aumenta y la potencia estadística disminuye. 
# Beta (Potencia estadística) para una diferencia de 2 unidades en DBP

beta_1 <- 1-t_3
(round(beta_1,digits=2))
## [1] 0.74
A continuación se explora la potencia estadística (poder) para distintos valores de diferencia (δ) en la hipótesis alternativa: 
# Beta para distintos valores de diferencia en la hipótesis alternativa

f <- -0.5
Potencia <- NULL
t_6 <- NULL

for (i in 1:6){

f <- f+1

t_4 <- (X_- f)/(sqrt(Sp)*sqrt((1/n_A)+(1/n_B)))
t_4 <- as.numeric(t_4)
(t_4<-round(t_4,digits=2))

t_5 <- pt(-t_4, 38,lower.tail = FALSE) 
round(t_5,digits=2)

t_6 [i] <- t_5
Potencia[i] <- 1-t_5
}
#Creamos una tabla para exponer los resultados

datosTabla <- 
data.frame(
  PoderPrueba= round(Potencia,digits=3),
  Error_Tipo_II = round(t_6,digits=3),
  Diferencia_en_Medias = c(0.5:5.5)) 

datosTabla
##   PoderPrueba Error_Tipo_II Diferencia_en_Medias
## 1       0.137         0.863                  0.5
## 2       0.524         0.476                  1.5
## 3       0.887         0.113                  2.5
## 4       0.989         0.011                  3.5
## 5       0.999         0.001                  4.5
## 6       1.000         0.000                  5.5
# Gráfico para visualizar el aumento de la potencia estadística

ggplot(
  data = datosTabla, aes(x = Diferencia_en_Medias, y = PoderPrueba)) +
  geom_line(aes(color='red')) +
  geom_point(aes(color='red')) +
scale_y_continuous(
  limits = c(0, 1),
  breaks = seq(0, 1, 0.1),
  name = "Potencia (beta)") +
  scale_x_continuous(
  limits = c(0.5, 5.5),
  breaks = seq(0.5, 5.5, 1),
  name = "Diferencia de medias en hipótesis alternativa") +
  ggtitle("Potencia estadística") + theme(plot.title = element_text(hjust = 0.5))

8 Tamaño de muestra

Para el cálculo del tamaño de muestra, se debe de contar con una diferencia mínima entre la media del valor de DBP del grupo de sujetos con el tratamiento A y la del grupo de sujetos del tratamiento B, esta diferencia será denotada con la letra δ (\(\mu_A-\mu_B = δ =1\)). La varianza se debe de conocer o estimar, por lo general se obtiene a partir de información de algún estudio previo desarrollado en una población objetivo similar a la que se tiene en el actual estudio. Supongamos que “cuidadosamente” de un estudio piloto se obtuvo una estimación de la desviación estándar igual a \(\sigma\)= 8 mmHg. Considerando un valor del error tipo I igual a \(\alpha\)= 0.05 y una potencia de prueba 1 − \(\beta\) = 0.80, los cuales son lo más comunes en investigaciones médicas, se tiene:

Para una prueba de hipótesis unilateral (una cola) del tipo:

\[H_0: \mu_A - \mu_B = 0 \text{ } \text{ Vs } \text{ } H_1: \mu_A - \mu_B > 0\]

que corresponde a una prueba de superioridad, el tamaño de muestra está dado por:

\[ \text{n} \ge \frac{2*[Z_{1-\alpha}+Z_{1-\beta}]^2}{(\frac{\mu_1 - \mu_2}{\sigma})^2} \]

donde:

  • \(\sigma\): es una estimación de la desviación estándar (variabilidad)
  • \(δ = (\mu_1 - \mu_2)\): es la diferencia a detectar entre grupos (clínicamente relevante)
  • \(Z_{1-\alpha}\): es el punto crítico de la distribución normal estándar para el error tipo I cuando la prueba es unilateral
  • \(Z_{1-\beta}\): es el punto crítico de la distribución normal estándar para el error tipo II cuando la prueba es unilateral

Para el tipo de hipótesis unilateral o de una cola planteada en el problema tenemos:

\(Z_{1-\alpha} = Z_{1-0.05} = Z_{0.95}=1.64\)

# Cuantil de una distribución normal con área de 0.95 a la izquierda

Z_1_alfa <- qnorm(0.95)
(round(Z_1_alfa,digits=2))
## [1] 1.64

\(Z_{1-\beta} = Z_{1-0.20} = Z_{0.80}\)

# Cuantil de una distribución normal con área de 0.80 a la izquierda

Z_1_beta <- qnorm(0.80)
(round(Z_1_beta,digits=2))
## [1] 0.84

sustituyendo valores

\[ \text{n} \ge \frac{2*[Z_{0.95}+Z_{0.80}]^2}{(\frac{1}{8})^2} \]

# Determinamos los valores dados por el enunciado

desvEsta <- 8
mA_mB <- 1
# Tamaño de muestra

n <- (2*(Z_1_alfa+Z_1_beta)^2)/((mA_mB)/desvEsta)^2
(ceiling(n))
## [1] 792
# Tamaño por grupo

n_2 <- n/2
ceiling(n_2)
## [1] 396
En cada grupo se requerirán 396 participantes para lograr detectar una diferencia de δ = 1, con un poder de la prueba de 0.8 (80%) y una probabilidad de cometer error tipo I igual a 0.05 (5%).

9 Prueba de comparación: t de student

Utilizamos la prueba t de student (muestras pareadas) para comparar las medias de los dos grupos de tratamientos (A y B). Otro ejemplo ver: https://rpubs.com/Humberto_Mtz_Bta/906503.

\(H_0: \text{δ}=\mu_A - \mu_B = 0\) Vs \(H_1: \text{δ}=\mu_A - \mu_B > 0\) 


# Prueba t para la diferencia promedio en DBP entre los tratamientos

Prueba <- t.test(data$di ~ data$Tratamiento, alternative = "greater", paired= TRUE)
Prueba
## 
##  Paired t-test
## 
## data:  data$di by data$Tratamiento
## t = 12.701, df = 19, p-value = 4.939e-11
## alternative hypothesis: true mean difference is greater than 0
## 95 percent confidence interval:
##  8.984088      Inf
## sample estimates:
## mean difference 
##            10.4
Dado que el valor p es menor que 0.05, se encuentra que existe una diferencia estadísticamente significativa entre la diferencia promedio de los tratamientos A y B. 


# Distribución de mediciones 

ggplot(
  data = Datos1,
  aes(x = , y = di, fill=Tratamiento)) +
  geom_boxplot(stat = "boxplot", position = "dodge2") +
  scale_x_discrete(name = "") +
  labs ( y='Diferencia promedio de DBP')+
  ggtitle("Distribución de la diferencia en DBP") + 
  theme(plot.title = element_text(hjust = 0.5))+
  facet_wrap(~Tratamiento)

10 Conclusión

El tratamiento A presenta una mayor disminución en la presión arterial diastólica en comparación con el tratamiento B (p < 0.001), con base en la prueba t de comparación.