Prueba de hipótesis para la diferencia de medias

Suponga que se tienen dos muestras aleatorias que provienen de poblaciones normales así:

$n_1$ observaciones $X_{11}, X_{12}, \ldots, X_{1,n1}$ de una población I con media $\mu_1$ y varianza $\sigma^2_1$,
$n_2$ observaciones $X_{21}, X_{22}, \ldots, X_{2,n2}$ de una población II con media $\mu_2$ y varianza $\sigma^2_2$,
ambas muestras son independientes entre sí.

En este problema se quiere estudiar la hipótesis nula $H_0: \mu_1 - \mu_2 = \delta_0$ y se sospecha que la diferencia de medias $\mu_1 - \mu_2$ podría estar en alguna de las siguientes situaciones:

$H_a: \mu_1 - \mu_2 < \delta_0$
$H_a: \mu_1 - \mu_2 \neq \delta_0$
$H_a: \mu_1 - \mu_2 > \delta_0$

El valor $\delta_0$ corresponde a la diferencia entre las medias poblacionales; cuando $\delta_0=0$, la hipótesis nula se convierte en $H_0: \mu_1 = \mu_2$ lo que significa que se está probando igualdad de medias.

Para realizar esta prueba de hipótesis se deben diferenciar dos casos, uno en el que las varianzas poblacionales son iguales y otro caso en el que las varianzas poblaciones son diferentes, esto se puede chequear utilizando la prueba de comparación de varianzas. Para cada uno de los casos hay un estadístico de prueba, una distribución asociada y un valor-$P$, a continuación se presentan los dos casos en detalle.

Caso 1: varianzas poblacionales iguales $\sigma_1^2 = \sigma_2^2$

En este caso el estadístico para realizar la prueba es: \[ T_0=\frac{\bar{X}_1 - \bar{X}_2 - \delta_0}{S_p \sqrt{\frac{1}{n_1} + \frac{1}{n_2}}}, \] donde $\bar{X}_1$ y $\bar{X}_2$ son las medias de las muestras I y II respectivamente; la cantidad $S_p^2$ es una varianza combinada y se calcula como: \[S_p^2=\frac{(n_1-1)S_1^2+(n_2-1)S_2^2}{n_1+n_2-2},\] donde $S_1^2$ y $S_2^2$ son las varianzas de las muestras I y II respectivamente. En este caso el estadístico $T_0$, bajo la suposición de que $H_0$ es verdadera, tiene distribución $t$-student con $n_1+n_2-2$ grados de libertad .

Si el valor calculado para el estadístico dado en la ecuación se denota por $t_0$, entonces el valor-$P$ de la prueba se calcula de acuerdo a la hipótesis alterna $H_a$ así:

Si $H_a: \mu_1 - \mu_2 < \delta_0$ entonces valor-$P$=$P(t_{n_1+n_2-2} < t_0)$.
Si $H_a: \mu_1 - \mu_2 \neq \delta_0$ entonces valor-$P$=$2 \times P(t_{n_1+n_2-2} > \lvert t_0 \rvert)$.
Si $H_a: \mu_1 - \mu_2 > \delta_0$ entonces valor-$P$=$P(t_{n_1+n_2-2} > t_0)$.

Caso 2: varianzas poblacionales diferentes $\sigma_1^2 \neq \sigma_2^2$

En este caso el estadístico para realizar la prueba es: \[ T_0=\frac{\bar{X}_1 - \bar{X}_2 - \delta_0}{\sqrt{\frac{S_1^2}{n_1} + \frac{S_2^2}{n_2}}}. \] Bajo la suposición de que $H_0$ es verdadera, $T_0$ tiene distribución $t$-student con $v$ grados de libertad , en donde $v$ se calcula como: \[v=\frac{ \left( \frac{S_1^2}{n_1} + \frac{S_2^2}{n_2} \right)^2 }{ \frac{(S_1^2/n_1)^2}{n_1-1} + \frac{(S_2^2/n_2)^2}{n_2-1}}\]

Si el valor calculado para el estadístico dado en la ecuación se denota por $t_0$, entonces el valor-$P$ de la prueba se calcula de acuerdo a la hipótesis alterna $H_a$ así:

Si $H_a: \mu_1 - \mu_2 < \delta_0$ entonces valor-$P$=$P(t_{v} < t_0)$.
Si $H_a: \mu_1 - \mu_2 \neq \delta_0$ entonces valor-$P$=$2 \times P(t_{v} > \lvert t_0 \rvert)$.
Si $H_a: \mu_1 - \mu_2 > \delta_0$ entonces valor-$P$=$P(t_{v} > t_0)$.

# Para generar base de datos DIFERENCIA DE MEDIAS ------------------------- dt <- data.frame(y=c(rnorm(n=100, mean=170, sd=5), rnorm(n=150, mean=173, sd=5)), group=rep(c('Grupo 1', 'Grupo 2'), times=c(100, 150))) write(x=t(dt), file='equal_var_data.txt', ncolumns=2, sep='\t') dt <- data.frame(y=c(rnorm(n=100, mean=65, sd=10), rnorm(n=150, mean=70, sd=5)), group=rep(c('Grupo 1', 'Grupo 2'), times=c(100, 150))) write(x=t(dt), file='unequal_var_data.txt', ncolumns=2, sep='\t') # Para generar base de datos MEDIAS --------------------------------------- dt <- data.frame(y=rnorm(n=100, mean=60, sd=5), x=rnorm(n=100, mean=170, sd=15)) write(x=t(dt), file='means_data.txt', ncolumns=2, sep='\t')

library(shiny) shinyServer(function(input,output,session){ observe({ inFile <- input$file1 if(is.null(inFile)) dt <- read.table('unequal_var_data.txt', header=T, sep='\t') else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) updateSelectInput(session, "variable1", choices=names(dt[!sapply(dt, is.factor)])) # Asegurar cuanti updateSelectInput(session, "variable2", choices=names(dt[sapply(dt, is.factor)])) # Asegurar cuali }) output$summary <- renderTable({ inFile <- input$file1 if(is.null(inFile)) dt <- read.table('unequal_var_data.txt', header=T, sep='\t') else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) dt <- na.omit(dt) # Para eliminar obs con NA dt }) output$statistic <- renderTable({ inFile <- input$file1 if(is.null(inFile)) dt <- read.table('unequal_var_data.txt', header=T, sep='\t') else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) dt <- na.omit(dt) # Para eliminar obs con NA x <- dt[, input$variable1] # Variable de interes group <- dt[, input$variable2] # Variable de clasificacion if (nlevels(group) != 2) group <- dt[, sapply(dt, is.factor)] xx <- split(x, group) # Lista con variable interes resumen <- function(x) c(mean(x), var(x), length(x)) res <- sapply(xx, resumen) rownames(res) <- c('Media', 'Varianza', 'Número obs') t(res) }, rownames = TRUE, align='c', bordered = TRUE) # Para obtener tabla con rownames output$appPlot <- renderPlot({ inFile <- input$file1 if(is.null(inFile)) dt <- read.table('unequal_var_data.txt', header=T, sep='\t') else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) # Aqui inicia la figura par(mfrow=c(1, 2), bg='gray98') dt <- na.omit(dt) # Para eliminar obs con NA x <- dt[, input$variable1] group <- dt[, input$variable2] if (nlevels(group) != 2) group <- dt[, sapply(dt, is.factor)] # Para dibujar las densidades xx <- split(x, group) den <- lapply(xx, density) plot(den[[1]], lwd=4, col='deepskyblue3', main='Densidad', las=1, xlab=as.character(input$variable1), ylab='Densidad', xlim=range(range(den[[1]]$x), range(den[[2]]$x)), ylim=c(0, max(c(den[[1]]$y, den[[2]]$y)))) lines(den[[2]], lwd=4, col='firebrick3') # Leyenda para distiguir las densidades legend('topright', bty='n', lwd=4, col=c('deepskyblue3', 'firebrick3'), legend=unique(group)) # Para dibujar los qqplot qq1 <- qqnorm(xx[[1]], plot.it=FALSE) qq2 <- qqnorm(xx[[2]], plot.it=FALSE) plot(qq1, las=1, main='QQplot', pch=19, col='deepskyblue3', xlim=range(c(qq1$x, qq2$x)), ylim=range(c(qq1$y, qq2$y)), xlab='Cuantiles teóricos N(0, 1)', ylab=as.character(input$variable1)) points(qq2, pch=19, col='firebrick3') # Para construir los qqplot qqline(xx[[1]], col='deepskyblue3') qqline(xx[[2]], col='firebrick3') # Para incluir el valor P de Shapiro shapi <- lapply(xx, shapiro.test) leyenda <- c(paste('Valor P=', round(shapi[[1]]$p.value, 2)), paste('Valor P=', round(shapi[[2]]$p.value, 2))) legend('topleft', bty='n', text.col=c('deepskyblue3', 'firebrick3'), legend=leyenda) }) output$resul1 <- renderText({ inFile <- input$file1 if(is.null(inFile)) dt <- read.table('unequal_var_data.txt', header=T, sep='\t') else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) dt <- na.omit(dt) # Para eliminar obs con NA x <- dt[, input$variable1] group <- dt[, input$variable2] if (nlevels(group) != 2) group <- dt[, sapply(dt, is.factor)] xx <- split(x, group) ph <- t.test(x=xx[[1]], y=xx[[2]], alternative=input$h0, mu=input$delta0, conf.level=input$alfa, var.equal=input$var.equal) paste0('El estadístico de prueba es t0=', round(ph$statistic, 4), ' con un valor-P de ', round(ph$p.value, 4), '.') }) output$resul2 <- renderText({ inFile <- input$file1 if(is.null(inFile)) dt <- read.table('unequal_var_data.txt', header=T, sep='\t') else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) dt <- na.omit(dt) # Para eliminar obs con NA x <- dt[, input$variable1] group <- dt[, input$variable2] if (nlevels(group) != 2) group <- dt[, sapply(dt, is.factor)] xx <- split(x, group) ph <- t.test(x=xx[[1]], y=xx[[2]], alternative=input$h0, mu=input$delta0, conf.level=input$alfa, var.equal=input$var.equal) intervalo <- paste("(", round(ph$conf.int[1], digits=4), ", ", round(ph$conf.int[2], digits=4), ").", sep='') paste0('El intervalo de confianza del ', 100*input$alfa, '% para la diferencia de medias poblacionales es ', intervalo) }) output$miteoria <- renderUI({ HTML(markdown::markdownToHTML(knit(input='teoria.md', quiet = TRUE))) }) })

library(shiny) library(markdown) shinyUI(pageWithSidebar( headerPanel(title=HTML("Prueba de hipótesis para la diferencia de medias μ1 - μ2"), windowTitle="PH dif medias"), sidebarPanel( h5('Esta aplicación realiza la prueba de hipótesis para la diferencia de medias de variables cuantitativas que tengan distribución normal.'), h6('La aplicación usa una base de datos de ejemplo pero el usuario puede cargar su propia base de datos.'), fileInput(inputId='file1', label='Use el botón siguiente para cargar su base de datos.', accept = c( 'text/csv', 'text/comma-separated-values', 'text/tab-separated-values', 'text/plain', '.csv', '.tsv' )), checkboxInput(inputId='header', label='¿Tiene encabezado la base de datos?', value=TRUE), selectInput(inputId="sep", label = "¿Cuál es la separación de los datos?", choices = list(Tab='\t', Comma=',', Semicolon=';', 'Space'=' '), selected = ';'), selectInput(inputId="variable1", label=p("Elija la variable", span("cualitativa", style = "color:red"), "para realizar la prueba de hipótesis."), choices=""), selectInput(inputId="variable2", label=p("Elija la variable", span("cualitativa", style = "color:blue"), "de agrupacion, DEBE tener 2 niveles y ser un factor."), choices=""), numericInput(inputId='delta0', label=HTML("Ingrese el valor de referencia δ0 para la probar H0: μ1 - μ2 = δ0"), value=0), selectInput(inputId="h0", label=HTML("Elija el tipo de hipotesis alterna < , ≠ o >"), choices=list("Menor" = "less", "Diferente" = "two.sided", "Mayor" = "greater"), selected = "two.sided"), checkboxInput(inputId='var.equal', label='Marque la casilla si las varianzas poblacionales son iguales', value=FALSE), sliderInput(inputId='alfa', label=HTML("Opcional: elija un nivel de confianza para construir el intervalo de confianza para la diferencia μ1 - μ2"), min=0.90, max=0.99, value=0.95, step=0.01), img(src="logo.png", height = 60, width = 120), img(src="udea.png", height = 25, width = 70), img(src="cua.png", height = 40, width = 110), br(), tags$a(href="https://srunal.github.io", "https://srunal.github.io") ), mainPanel( tabsetPanel(type = "pills", tabPanel(title="Resultados", h5('A continuación el histograma, la densidad, el QQplot y valor-P para la prueba de normalidad Shapiro-Wilk para cada una de las dos muestras.'), plotOutput("appPlot", width='500px', height='300px'), h4("- Tabla de resumen con los estadísticos muestrales:"), tableOutput('statistic'), h4("- Resultados de la prueba de hipótesis:"), textOutput("resul1"), h4(HTML("- Intervalo de confianza para la diferencia de medias μ1 - μ2:")), textOutput("resul2")), tabPanel("Datos", "A continuación los datos que está usando la aplicación.", uiOutput('summary')), tabPanel("Teoría", includeHTML("include.html")) ) ) ))

Prueba de hipótesis para la diferencia de medias μ₁ - μ₂

Esta aplicación realiza la prueba de hipótesis para la diferencia de medias de variables cuantitativas que tengan distribución normal.

La aplicación usa una base de datos de ejemplo pero el usuario puede cargar su propia base de datos.

A continuación el histograma, la densidad, el QQplot y valor-P para la prueba de normalidad Shapiro-Wilk para cada una de las dos muestras.

- Tabla de resumen con los estadísticos muestrales:

- Resultados de la prueba de hipótesis:

- Intervalo de confianza para la diferencia de medias μ₁ - μ₂:

Prueba de hipótesis para la diferencia de medias

Caso 1: varianzas poblacionales iguales \(\sigma_1^2 = \sigma_2^2\)

Caso 2: varianzas poblacionales diferentes \(\sigma_1^2 \neq \sigma_2^2\)