Alumno: Daniel Lorenzo Medina Salcedo
Correo: daniel.medina@cuantico.com
Asignatura: Estadística para investigadores
Profesor: Dra. Maria Luisa Díez Platas
Doctorado en Ciencias de la Computación
UNIR
Abril 2025
Este documento se puede consulta en línea en el siguiente link:
Una administración de lotería emite boletos rasca y gana, de los cuales el 2 % deben tener premio. Una asociación de consumidores estudia los resultados de una muestra de 500 boletos para asegurarse de que no se estén emitiendo menos premios de los debidos obteniendo un porcentaje de premios del \(\hat{p} = 1,6 \%\). Para saber si debemos seguir confiando en que se están emitiendo un 2 % de boletos premiados, se plantea el siguiente contraste de hipótesis:
Calcular la desviación estándar de la proporción de la hipótesis nula \(\sigma_0 = \sqrt{\frac{p_0(1 - p_0)}{n}}\), la región de aceptación de \(H_0\) y valor crítico \(c_\alpha\) con nivel de confianza del 90 %. Para ello sabemos que:
\[ RA_z(H_0; p_0, \sigma_0) = [p_0 - |z_\alpha| \cdot SD(p_0), 1] = [c_l, 1] \]
# Parámetros dados
n <- 500
p_0 <- 0.02
p_hat <- 0.016
conf_level <- 0.90
alpha <- 1 - conf_level
# Cálculo de la desviación estándar bajo la hipótesis nula
sd_0 <- sqrt(p_0 * (1 - p_0) / n)
# Z crítico para el nivel de confianza
z_alpha <- qnorm(alpha)
# Cálculo de la región de aceptación
cl <- p_0 - abs(z_alpha) * sd_0
sd_0
## [1] 0.00626099
z_alpha
## [1] -1.281552
cl
## [1] 0.01197622
# Crear valores para la distribución
x <- seq(0, 0.05, length.out = 1000)
y <- dnorm(x, mean = p_0, sd = sd_0)
# Graficar la distribución normal con color morado
plot(x, y, type = "l", col = "purple",
xlab = "Proporción de boletos premiados",
ylab = "Densidad de probabilidad",
main = "Desviación Estándar, Región de Aceptación y Valor Crítico",
lwd = 2)
# Agregar la región crítica (zona de rechazo) con color rojo
polygon(c(x[x < p_0 - abs(z_alpha) * sd_0], p_0 - abs(z_alpha) * sd_0),
c(y[x < p_0 - abs(z_alpha) * sd_0], 0), col = rgb(1, 0, 0, alpha = 0.5)) # Rojo
# Agregar la zona de aceptación con color gris claro
x_acc <- x[x >= cl]
y_acc <- y[x >= cl]
polygon(c(cl, x_acc, max(x_acc), cl),
c(0, y_acc, 0, 0),
col = rgb(0.8, 0.8, 0.8, 0.5), border = NA)
# Marcar p_hat
abline(v = p_hat, col = "green", lty = 2, lwd = 2)
# Marcar p_0 con una línea vertical
abline(v = p_0, col = "blue", lty = 2, lwd = 2)
# Marcar z_alpha con una línea vertical
abline(v = p_0 - abs(z_alpha) * sd_0, col = "red", lty = 2, lwd = 2)
# Agregar leyenda
legend("topright", legend = c(paste("p^ = ", p_hat), paste("Desviación estándar =", round(sd_0, 4)),
paste("p0 = 0.02"), paste("z_alpha =", round(z_alpha, 4)), "Región Crítica", "Región de Aceptación"),
col = c("green", "purple", "blue", "red", "red", "gray"), lty = c(2, 1, 1, 2, 1, 1), lwd = 2, pch = c(NA, NA, NA, NA, 15, 15), pt.cex = c(1, 1, 1, 1, 2, 2))
pnorm() para
obtener el p-valor de \(H_0\).
Otra forma de obtener el p-valor de \(H_0\) es con el comando:
prop.test(x, n, p, alternative="less", conf.level, correct=FALSE)
Donde x es el número de boletos premiados, n el número de boletos de la muestra, p la proporción de boletos premiados en la hipótesis nula, alternative=“less” indica que la hipótesis alternativa es de cola izquierda, conf.level es el nivel de confianza, y correct=FALSE evita la corrección de continuidad que se hace en las distribuciones discretas.
Ejecute este comando con los parámetros adecuados y compruebe que obtiene el mismo p-valor ¿Debemos seguir confiando en que se están emitiendo un 2% de boletos premiados?
Nota: El comando prop.test() nos devuelve el \(IC(\hat{p})\) utilizando un método distinto
al visto en clase llamado el método Wilson.
# Calcular el p-valor usando pnorm
z_score <- (p_hat - p_0) / sd_0
p_value <- pnorm(z_score)
p_value
## [1] 0.2614516
# 1. Dibujar la curva
curve(dnorm(x, mean = p_0, sd = sd_0),
from = 0, to = 0.05,
col = "blue", lwd = 2,
xlab = "Proporción",
ylab = "Densidad",
main = "P - Valor")
# 2. Rellenar la región crítica (x < c_l) en rojo
x_crit <- seq(0, cl, length.out = 200)
y_crit <- dnorm(x_crit, mean = p_0, sd = sd_0)
polygon(c(x_crit, rev(x_crit)),
c(y_crit, rep(0, length(y_crit))),
col = rgb(1, 0, 0, 0.5), border = NA)
# 3. Marcar p-hat con una línea verde
abline(v = p_hat, col = "green", lty = 2, lwd = 2) # p̂
abline(v = cl, col = "red", lty = 2, lwd = 2) # c_l (valor crítico)
# 4. Anotación del p‑valor en centro vertical
y_centro <- max(y) / 2
text(x = p_hat, y = y_centro,
labels = bquote(p-valor == .(round(p_value, 4))),
col = "green4", cex = 1)
# 5. Leyenda
legend("topright",
legend = c(
bquote(p-valor == .(round(p_value,4))),
bquote(z[alpha] == .(round(z_alpha, 4))),
"Región Crítica"
),
col = c("green", "red", rgb(1,0,0,0.5)),
lty = c(2, 2, NA),
lwd = c(2, 2, NA),
pch = c(NA, NA, 15),
pt.cex = c(NA, NA, 2),
pt.bg = c(NA, NA, rgb(1,0,0,0.5))
)
# Número de boletos premiados
x <- p_hat * n
# Realizar el test de proporciones
test_result <- prop.test(x, n, p_0, alternative = "less", conf.level = conf_level, correct = FALSE)
test_result$p.value
## [1] 0.2614516
## **No se rechaza H₀**
## – Estadístico de prueba: z = -0.6389
## – p‑valor = 0.2615 > 0.1
## – \(\hat{p} = 0.016\) cae en la zona de aceptación (\(c_l \approx 0.012\)).
CONCLUSIONES:
Los estudios poblacionales dicen que el vello corporal femenino crece de media 13,1 mm/mes. Un laboratorio anuncia que su medicamento «reduce significativamente la velocidad de crecimiento del vello en mujeres». Nuestra asociación de consumidores ha realizado un estudio sobre una muestra de n=30 mujeres a las que ha suministrado el medicamento obteniendo los siguientes crecimientos en mm/mes:
12.9, 2.92, 14.38, 9.03, 13.17, 19.98, 15.4, 7.41, 6.19, 10.32, 16.9, 11.33, 11.31, 12.33, 15.28, -1.53, 20.87, 11.99, 14.23, 4.23, 10.38, 11.39, 6.65, 9.96, 8.58, 14.39, 10.73, 16.61, 21.6, 12.74
En este ejercicio desconocemos la varianza poblacional \(\sigma\), por lo que debemos realizar un t‑test.
Recordar que, tal como muestra la figura, el intervalo de confianza se calcula sobre la distribución muestral centrada en \(\hat\mu\).
Con significancia \(\alpha = 0{,}05\) calcule la media muestral estimada \(\hat\mu\), la desviación estándar muestral estimada \(\hat\sigma\), el error estándar de la muestra \(\hat\sigma_{\bar X} = SE(\hat\mu)\), y el intervalo de confianza para \(\mu\). Para hacer este cálculo sabemos que:
\[ IC_t(\mu) \;=\; \bigl[\hat\mu \;-\; \lvert t_{\alpha/2,\,df}\rvert \,SE(\hat\mu),\;\hat\mu \;+\; \lvert t_{\alpha/2,\,df}\rvert \,SE(\hat\mu)\bigr] \]
También podemos calcular el \(IC_t(\mu)\) usando el siguiente comando,
donde x es la muestra, y conf.level indica el
nivel de confianza.
t.test(x, conf.level = 1 - alpha)
Ejecute este comando con los parámetros adecuados e indique si los resultados obtenidos son consistentes con los cálculos anteriores. De no ser así, explique a qué se debe la diferencia.
# Crecimientos en mm/mes
x <- c(12.9, 2.92, 14.38, 9.03, 13.17, 19.98, 15.4, 7.41, 6.19, 10.32, 16.9, 11.33, 11.31, 12.33, 15.28, -1.53, 20.87, 11.99, 14.23, 4.23, 10.38, 11.39, 6.65, 9.96, 8.58, 14.39, 10.73, 16.61, 21.6, 12.74)
n <- length(x)
mu_0 <- 13.1
alpha <- 0.05
# Media muestral
mu_hat <- mean(x)
# Desviación estándar muestral
sd_hat <- sd(x)
# Error estándar
se_x_bar <- sd_hat / sqrt(n)
mu_hat
## [1] 11.72233
sd_hat
## [1] 5.123752
se_x_bar
## [1] 0.9354649
# Graficar la distribución muestral
hist(x, col = "lightblue", main = "Distribución de Crecimientos de Vello",
xlab = "Crecimiento (mm/mes)", ylab = "Frecuencia", breaks = 10)
abline(v = mu_hat, col = "red", lwd = 2)
legend("topright", legend = paste("Media muestral:", round(mu_hat, 2)), col = "red", lwd = 2)
# Cálculo del intervalo de confianza usando t.test
t_test_result <- t.test(x, conf.level = 1 - alpha)
t_test_result$conf.int
## [1] 9.809093 13.635574
## attr(,"conf.level")
## [1] 0.95
# Graficar el intervalo de confianza
hist(x, col = "lightblue", main = "Intervalo de Confianza para la Media",
xlab = "Crecimiento (mm/mes)", ylab = "Frecuencia", breaks = 10)
abline(v = t_test_result$conf.int, col = "green", lwd = 2)
legend("topright", legend = paste("IC:", round(t_test_result$conf.int[1], 2),
"a", round(t_test_result$conf.int[2], 2)), col = "green", lwd = 2)
Con significancia \(\alpha = 0{,}05\) calcule el valor crítico \(c_l\) para \(H_0\). Para hacer este cálculo sabemos que:
\[ RA_t\bigl(H_0; \mu_0, \hat\sigma_{\bar X}\bigr) \;=\; \bigl[\mu_0 \;-\; \lvert t_{\alpha,\,df}\rvert \,SE(\hat\mu),\;\infty\bigr] \;=\; [c_l,\;\infty] \]
# Estadísticos muestrales
mu_hat <- mean(x)
sd_hat <- sd(x)
se_x_bar <- sd_hat / sqrt(n)
# t crítico y valor crítico c_l (una cola izquierda)
t_alpha <- qt(alpha, df = n - 1) # t_(α, df)
c_l <- mu_0 + t_alpha * se_x_bar # límite inferior de aceptación
# Crear secuencia para la distribución de la media muestral
x_seq <- seq(mu_0 - 4*se_x_bar, mu_0 + 4*se_x_bar, length.out = 1000)
# Densidad de la t-Student escalada para la media
y_seq <- dt((x_seq - mu_0)/se_x_bar, df = n - 1) / se_x_bar
# Graficar la curva de la distribución muestral
plot(x_seq, y_seq, type = "l", lwd = 2, col = "purple",
xlab = expression(bar(X)),
ylab = "Densidad",
main = "Distribución de la media muestral y región crítica")
# Rellenar región crítica (cola izquierda)
crit_x <- x_seq[x_seq < c_l]
crit_y <- y_seq[x_seq < c_l]
polygon(c(crit_x, rev(crit_x)),
c(crit_y, rep(0, length(crit_y))),
col = rgb(1, 0, 0, 0.5), border = NA)
# Líneas de referencia
abline(v = mu_0, col = "blue", lty = 3, lwd = 2) # mu_0
abline(v = mu_hat, col = "green", lty = 2, lwd = 2) # media muestral
abline(v = c_l, col = "red", lty = 2, lwd = 2) # valor crítico
# Leyenda
legend("topright",
legend = c(
bquote(mu[0] == .(mu_0)),
bquote(bar(X) == .(round(mu_hat, 2))),
bquote(c[l] == .(round(c_l, 2))),
"Región crítica (cola izquierda)"
),
col = c("blue", "green", "red", rgb(1,0,0,0.5)),
lty = c(3, 2, 2, NA),
lwd = c(2, 2, 2, NA),
pch = c(NA, NA, NA, 15),
pt.cex = c(NA, NA, NA, 2),
pt.bg = c(NA, NA, NA, rgb(1,0,0,0.5))
)
A continuación calcule el p‑valor con el siguiente comando:
pt(t_score, df = n - 1, lower.tail = TRUE)
También podemos calcular el t‑score y p‑valor usando el siguiente
comando, donde x es la muestra, mu indica la
media poblacional de \(H_0\),
alternative indica cuál es la hipótesis alternativa, y
conf.level indica el nivel de confianza:
t.test(x, mu = mu_0,alternative = "less",conf.level = 1 - alpha)
Ejecute este comando con los parámetros adecuados e indique si el
t-score para \(\hat\mu\) y el
p‑valor obtenidos por t.test() son consistentes
con los cálculos anteriores. De no ser así, explique a qué se debe.
Indique también si refutaría,
aceptaría o dejaría en duda \(H_0\).
# ---- Cálculo t-score y p-valor manual ----
# ---- Cálculo de estadísticas ----
t_score <- (mu_hat - mu_0) / se_x_bar
p_value <- pt(t_score, df = n - 1, lower.tail = TRUE)
# Valor t crítico para cola izquierda
t_alpha <- qt(alpha, df = n - 1)
# ---- Comparativa con t.test() ----
res_t <- t.test(x, mu = mu_0, alternative = "less", conf.level = 1 - alpha)
t_ttest <- as.numeric(res_t$statistic)
p_ttest <- res_t$p.value
# ---- Preparar distribución t ----
x_t <- seq(-4, 4, length.out = 500)
y_t <- dt(x_t, df = n - 1)
# ---- Gráfico ----
plot(x_t, y_t, type = "l", lwd = 2, col = "darkblue",
xlab = "t", ylab = "Densidad",
main = "Cálculo t-score, p-valor y región crítica")
# ---- Región crítica (t <= t_alpha) ----
crit_x <- x_t[x_t <= t_alpha]
crit_y <- y_t[x_t <= t_alpha]
polygon(c(crit_x, rev(crit_x)),
c(crit_y, rep(0, length(crit_y))),
col = rgb(1, 0, 0, 0.5), border = NA)
# ---- Líneas verticales ----
abline(v = t_score, col = "red", lty = 2, lwd = 2) # t-score observado
abline(v = t_alpha, col = "purple", lty = 3, lwd = 2) # t crítico
# ---- Anotación p‑valor ----
y_centro <- max(y_t) / 2
text(x = t_score, y = y_centro,
bquote(p-valor == .(round(p_value,4))),
col = "red4", cex = 1)
# ---- Leyenda ----
legend("topright",
legend = c(
bquote("t-score (manual) ="~.(round(t_score, 4))),
bquote("p-valor (manual) ="~.(round(p_value, 4))),
bquote("t-score (t.test)"~.(round(t_ttest, 4))),
bquote("p-valor (t.test)"~.(round(p_ttest, 4))),
bquote(alpha == .(alpha)),
"Región crítica"
),
col = c("red", "red4", "purple", "black", "black", rgb(1,0,0,0.5)),
lty = c(2, NA, 3, NA, NA, NA),
lwd = c(2, NA, 2, NA, NA, NA),
pch = c(NA, NA, NA, NA, NA, 15),
pt.cex = c(NA, NA, NA, NA, NA, 2),
pt.bg = c(NA, NA, NA, NA, NA, rgb(1,0,0,0.5)),
bty = "n"
)
CONCLUSIONES:
t.test().
Imaginemos que ahora descubrimos un estudio realizado con 5000 mujeres donde aparece que la desviación estándar de crecimiento de vello en mujeres es 3,2 mm/mes. Dado que el estudio se ha hecho sobre una muestra grande, podemos asumir que la desviación estándar poblacional es \(\sigma = 3{,}2\).
Usando el método del valor crítico y \(\alpha = 0{,}05\), calcule para nuestra muestra de \(n = 30\) mujeres la desviación estándar de la hipótesis nula \(\sigma_0 = SD(\mu_0)\), el z‑score \(z_\alpha\) y el valor crítico \(c_l\). Para ello nos dicen que la región de aceptación para la hipótesis nula es:
\[ RA_z\bigl(H_0; \mu_0, \sigma_0\bigr) =\bigl[\mu_0 - z_\alpha\,SD(\mu_0), \;\infty\bigr] =[c_l, \infty] \]
Calcule también el z‑score de \(\hat\mu\) y el p‑valor. Para ello sabemos que:
\[ z \;=\;\frac{\hat\mu - \mu_0}{\sigma_0} \]
Para validar los resultados anteriores, podemos utilizar el siguiente
comando para calcular la \(RA_z(H_0)\)
y el p‑valor, donde mu es la media poblacional y
sd la desviación estándar poblacional:
Ejemplo de comando z.test()
install.packages("TeachingDemos")
library(TeachingDemos)
# Datos y parámetros
# --- Datos y parámetros ---
x <- c(12.9, 2.92, 14.38, 9.03, 13.17, 19.98, 15.4, 7.41, 6.19,
10.32, 16.9, 11.33, 11.31, 12.33, 15.28, -1.53, 20.87,
11.99, 14.23, 4.23, 10.38, 11.39, 6.65, 9.96, 8.58,
14.39, 10.73, 16.61, 21.6, 12.74)
n <- length(x)
mu_0 <- 13.1
sigma_p <- 3.2 # σ poblacional conocido
alpha <- 0.05
# --- Estadísticos muestrales ---
mu_hat <- mean(x)
# --- Desviación estándar de la media bajo H0 ---
sigma_0 <- sigma_p / sqrt(n)
# --- Estadístico Z correcto ---
z_score <- (mu_hat - mu_0) / sigma_0
# --- P‑valor para prueba unilateral izquierda ---
p_value <- pnorm(z_score, lower.tail = TRUE)
# --- Comparativa con z.test() de TeachingDemos ---
library(TeachingDemos)
res_z <- z.test(x, mu = mu_0, sd = sigma_p, alternative = "less", conf.level = 1 - alpha)
z_ttest <- res_z$statistic
p_ttest <- res_z$p.value
# --- Valores críticos ---
z_crit <- qnorm(1 - alpha) # cuantil positivo
c_l <- mu_0 - z_crit * sigma_0 # límite de la región crítica
# --- Gráfico ---
x_norm <- seq(mu_0 - 4*sigma_0, mu_0 + 4*sigma_0, length.out = 500)
y_norm <- dnorm(x_norm, mean = mu_0, sd = sigma_0)
plot(x_norm, y_norm, type = "l", lwd = 2, col = "darkblue",
xlab = expression(bar(X)),
ylab = "Densidad",
main = "Contraste Z test con σ poblacional conocida")
# Región crítica (cola izquierda)
crit <- x_norm < c_l
polygon(c(x_norm[crit], rev(x_norm[crit])),
c(y_norm[crit], rep(0, sum(crit))),
col = rgb(1, 0, 0, 0.5), border = NA)
# Líneas de referencia
abline(v = mu_0, col = "blue", lty = 3, lwd = 2) # μ₀
abline(v = mu_hat, col = "green", lty = 2, lwd = 2) # μ̂
abline(v = c_l, col = "purple", lty = 2, lwd = 2) # c_l
# Anotación del p-valor
text(x = mu_hat, y = max(y_norm)*0.6,
labels = bquote(p-valor == .(round(p_value,4))),
col = "darkgreen", cex = 1)
# Leyenda arriba derecha en varias líneas
legend("topright",
legend = c(
bquote(mu[0] == .(mu_0)),
bquote(hat(mu) == .(round(mu_hat,2))),
bquote(c[l] == .(round(c_l,2))),
bquote("Z-score (manual) ="~.(round(z_score,4))),
bquote("p-valor (manual) <"~alpha),
bquote("Z-score (z.test)"~.(round(z_ttest,4))),
bquote("p-valor (z.test) <"~alpha)
),
col = c("blue", "green", "purple", "black", "darkgreen", "black", "darkgreen"),
lty = c(3, 2, 2, 1, NA, 1, NA),
lwd = c(2, 2, 2, 2, NA, 2, NA),
pch = c(NA, NA, NA, NA, 15, NA, 15),
pt.cex = c(NA, NA, NA, NA, 2, NA, 2),
pt.bg = c(NA, NA, NA, NA, rgb(0,1,0,0.5), NA, rgb(0,1,0,0.5)),
bty = "n"
)
CONCLUSIONES:
Z.test().
Un laboratorio desarrolla un nuevo test de embarazo y se realizan 50 experimentos que revelan los siguientes resultados, donde 0 significa «no embarazo» y 1 significa «embarazo».
decision <- c(0, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 0, 0)
realidad <- c(0, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 1, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 0, 0, 1, 1, 0, 0)
Calcula la significancia \(\alpha\), tolerancia \(\beta\), especificidad \(1 - \alpha\) y sensibilidad \(1 - \beta\) de este test de embarazo.
# Datos
decision <- c(0,0,1,1,1,1,0,1,1,0,1,0,1,0,0,0,0,0,0,0,1,0,1,0,1,0,0,1,1,0,1,1,0,0,1,0,0,1,1,1,1,1,1,0,0,0,1,1,0,0)
realidad <- c(0,1,1,1,1,1,0,1,1,0,1,0,0,0,0,0,0,0,0,1,1,0,1,0,1,0,0,1,1,0,1,1,0,0,1,0,1,1,1,0,1,1,1,0,0,0,1,1,0,0)
# Construir matriz de confusión
conf_mat <- table(Predicción = decision, Realidad = realidad)
# Calcular métricas
TP <- conf_mat["1","1"]
TN <- conf_mat["0","0"]
FP <- conf_mat["1","0"]
FN <- conf_mat["0","1"]
sensibilidad <- TP / (TP + FN)
especificidad <- TN / (TN + FP)
tolerancia <- 1 - sensibilidad
ylim_max <- max(conf_mat) * 1.5
# Gráfico de barras
bp <- barplot(conf_mat, beside = TRUE,
col = c("lightblue","lightgreen"),
ylim = c(0, ylim_max),
main = "Matriz de Confusión",
ylab = "Número de casos",
xlab = "Clase Real",
names.arg = c("No Embarazo","Embarazo"))
# Anotar conteos sobre cada barra
for(i in seq_len(nrow(conf_mat))) {
for(j in seq_len(ncol(conf_mat))) {
text(x = bp[j,i], y = conf_mat[i,j] + 2,
labels = conf_mat[i,j], cex = 0.9)
}
}
# Anotaciones de métricas un poco más abajo
usr <- par("usr")
x_pos <- mean(usr[1:2])
y_pos <- usr[4] * 0.83
text(x = x_pos, y = y_pos,
labels = paste0("Sensibilidad = ", round(sensibilidad, 2),
" | Especificidad = ", round(especificidad, 2),
" | Tolerancia = ", round(tolerancia, 2)),
cex = 0.9, font = 2)
# Leyenda centrada arriba, ligeramente dentro del gráfico
par(xpd = TRUE)
legend(x = 1.3, y = ylim_max ,
legend = c("Predicción = No Embarazo","Predicción = Embarazo"),
fill = c("lightblue","lightgreen"),
horiz = TRUE,
bty = "n")
par(xpd = FALSE)