library(haven)
data <- read_dta("Data1_R.dta")
View(data)Clase 11 de mayo
Encuesta Nacional de Salud y Nutrición 2018
Cargar las librerías necesarias
library(dplyr) # Para manipular los datos
library(ggplot2) # Para crear gráficos
library(“foreign”) #Cargas datos .dta
library(“msm”)
library(“car”)
require(plotrix)
library(“plotrix”)
Cambiar dirección de trabajo: setwd(“C:/Users/yomara.ruiz/Desktop/Clase demost. 11”)
head(data)# A tibble: 6 × 50
area empleo region edad t_hijos nac_vivo_murieron mortinato_2
<dbl+lbl> <dbl+lbl> <dbl+l> <dbl> <dbl> <dbl+lbl> <dbl+lbl>
1 1 [Urbano] 1 [Trabajó al… 1 [Sie… 19 1 0 [No] 0 [No]
2 1 [Urbano] 0 [No trabajó] 1 [Sie… 23 1 0 [No] 0 [No]
3 1 [Urbano] 1 [Trabajó al… 1 [Sie… 38 5 0 [No] 0 [No]
4 1 [Urbano] 0 [No trabajó] 1 [Sie… 18 1 0 [No] 0 [No]
5 1 [Urbano] 0 [No trabajó] 1 [Sie… 21 1 0 [No] 0 [No]
6 1 [Urbano] 1 [Trabajó al… 1 [Sie… 22 1 0 [No] 0 [No]
# ℹ 43 more variables: depresion_pp <dbl+lbl>, intensidad_dpp <dbl+lbl>,
# etnia <dbl+lbl>, f2_s2_216_1 <dbl+lbl>, f2_s2_216_2 <dbl>,
# f2_s2_218_1_a <dbl+lbl>, tiempo_dpp <dbl+lbl>, f2_s5_504a_1 <dbl+lbl>,
# f2_s5_504b_1 <dbl+lbl>, f2_s5_504c_1 <dbl+lbl>, f2_s5_504d_1 <dbl+lbl>,
# f2_s5_504e_1 <dbl+lbl>, f2_s5_504f_1 <dbl+lbl>, f2_s5_504g_1 <dbl+lbl>,
# f2_s5_504h_1 <dbl+lbl>, f2_s5_504i_1 <dbl+lbl>, f2_s5_504j_1 <dbl+lbl>,
# f2_s5_504k_1 <dbl+lbl>, est_civil <dbl+lbl>, q_usted <dbl+lbl>, …
media <- mean(data$ingrl, na.rm = TRUE) # Calcula la media
print(media)[1] 162.6633
desviacion <- sd(data$ingrl, na.rm = TRUE) # Calcula la desviación estándar
print(desviacion)[1] 329.8832
n <- length(data$ingrl) # Número de observaciones
print(n)[1] 16451
summary(data$ingrl) Min. 1st Qu. Median Mean 3rd Qu. Max.
0.0 0.0 0.0 162.7 200.0 3000.0
error_estandar <- desviacion / sqrt(n)
print(error_estandar)[1] 2.571959
z <- qnorm(0.975) # Z-score para un intervalo de confianza del 95%
print(z)[1] 1.959964
z2 <- qnorm(0.95) # Z-score para un intervalo de confianza del 90%
print(z2)[1] 1.644854
z3 <- qnorm(0.995) # Z-score para un intervalo de confianza del 99%
print(z3)[1] 2.575829
margen_error <- z * error_estandar # Margen de error
print(margen_error)[1] 5.040948
Crear intervalo de confianza
IC_inferior <- media - margen_error
print(IC_inferior)[1] 157.6223
IC_superior <- media + margen_error
print(IC_superior)[1] 167.7042
cat("El intervalo de confianza para la media de 'ingrl' es: [", IC_inferior, ",", IC_superior, "]\n")El intervalo de confianza para la media de 'ingrl' es: [ 157.6223 , 167.7042 ]
intervalo95<-cbind(IC_inferior,media,IC_superior);intervalo95 #Para transformar a matriz IC_inferior media IC_superior
[1,] 157.6223 162.6633 167.7042
colnames(intervalo95)<-c("IC_low","media","IC_high");intervalo95 #Para cambiar nombres más cortos de columnas IC_low media IC_high
[1,] 157.6223 162.6633 167.7042
#Análisis
Otra forma de calcular intervalos de confianza
media_test <- t.test(data$ingrl, conf.level = 0.95)
print(media_test$conf.int)[1] 157.6220 167.7046
attr(,"conf.level")
[1] 0.95
labels<-c("CI_95")
names(labels)<-labelslabelsre1<-round(c(media),2)
labelsre2<-round(c(IC_inferior),2)
labelsre3<-round(c(IC_superior),2) names(labelsre1)<-labelsre1xre1<-c(media)
y1<-(4)
lre1<-c(IC_inferior)
ure1<-c(IC_superior)options(repos = c(CRAN = "https://cloud.r-project.org"))
install.packages("plotrix")Installing package into 'C:/Users/yomara.ruiz/AppData/Local/R/win-library/4.5'
(as 'lib' is unspecified)
package 'plotrix' successfully unpacked and MD5 sums checked
The downloaded binary packages are in
C:\Users\yomara.ruiz\AppData\Local\Temp\Rtmpe6n9qG\downloaded_packages
library(plotrix)
plotCI(xre1,y1, ui=ure1,li=lre1,col="#898989",scol="#898989",
err="x",
axes=FALSE, ## disable axes (including tick labels)
pch=21,
pt.bg=16,
cex = 0.30,
slty = 1,
lwd=1,
xlab="",
ylab="",
ylim=c(1,7),
xlim=c(154,169), ## suppress x-axis label
main="Intervalo de Confianza para el Ingreso Promedio",
cex.main=0.85,
font.main = 1)
abline(v = xre1, col = "red",lty = 2)
#text(-200,"Media",col="blue")
axis(side=1,cex.axis=0.65) ## add default y-axis (ticks+labels)
axis(side=2,at=4, ## add custom x-axis
labels= labels,cex.axis=0.30)
box(bty="l") ## add box
text(xre1, y1+0.03, labelsre1, cex = 0.65, pos= 3)
text(lre1, y1+0.03, labelsre2, cex = 0.65, pos= 3)
text(ure1, y1+0.03, labelsre3, cex = 0.65, pos= 3)Ejemplo de 2
install.packages("haven")Warning: package 'haven' is in use and will not be installed
library(haven)
p <- sum(as_factor(data$depresion_pp) == "Si", na.rm = TRUE) / sum(!is.na(data$depresion_pp))
print(p)[1] 0.2202298
n <- length(data$depresion_pp)# Nivel de confianza (95%) -> Z-value 1.96
z <- qnorm(0.975) # Z-score para un intervalo de confianza del 95%
z2<- qnorm(0.95) # Z-score para un intervalo de confianza del 90%
z3<- qnorm(0.995) # Z-score para un intervalo de confianza del 99%error_estandarp <- sqrt((p * (1 - p)) / n)
print(error_estandarp)[1] 0.003230912
margen_errorp <- z * error_estandarp
print(margen_errorp)[1] 0.006332472
IC_inferiorp <- p - margen_errorp
print(IC_inferiorp)[1] 0.2138973
IC_superiorp <- p + margen_errorp
print(IC_superiorp)[1] 0.2265622
cat("El intervalo de confianza para la media de 'mujeres con depresión post partol' es: [", IC_inferiorp, ",", IC_superiorp, "]")El intervalo de confianza para la media de 'mujeres con depresión post partol' es: [ 0.2138973 , 0.2265622 ]
intervalo95p<-cbind(IC_inferiorp,p,IC_superiorp);intervalo95p #Para transformar a matriz IC_inferiorp p IC_superiorp
[1,] 0.2138973 0.2202298 0.2265622
colnames(intervalo95p)<-c("IC_low","p","IC_high");intervalo95p #Para cambiar nombres más cortos de columnas IC_low p IC_high
[1,] 0.2138973 0.2202298 0.2265622
Una forma más simplificada de calcular un intervalo de confianza para proporciones
library(haven)
depresion_factor <- as_factor(data$depresion_pp)
prop_test <- prop.test(
sum(depresion_factor == "Si", na.rm = TRUE), # éxitos
sum(!is.na(depresion_factor)), # total válido
conf.level = 0.95
)
print(prop_test$conf.int)[1] 0.2139329 0.2266579
attr(,"conf.level")
[1] 0.95
PRUEBAS DE HIPÓTESIS
Supongamos que queremos realizar una prueba de hipótesis para la media de ‘ingrl’
H0: El ingreso promedio poblacional es igual a $450 HA: El ingreso promedio poblacional es diferente a $450 Recordar que la HA nunca contiene los signos “=” , “≤” o “≥”.
t_prueba <- t.test(data$ingrl, mu = 450, conf.level = 0.95)
print(t_prueba)
One Sample t-test
data: data$ingrl
t = -111.72, df = 16450, p-value < 2.2e-16
alternative hypothesis: true mean is not equal to 450
95 percent confidence interval:
157.6220 167.7046
sample estimates:
mean of x
162.6633
if (t_prueba$p.value < 0.05) {
cat("Rechazamos la hipótesis nula: Entonces tenemos evidencia de que la media de 'ingrl' es diferente de 450.\n")
} else {
cat("No rechazamos la hipótesis nula: Entonces tenemos evidencia de que la media de 'ingrl' es igual a 450.\n")
}Rechazamos la hipótesis nula: Entonces tenemos evidencia de que la media de 'ingrl' es diferente de 450.
Supongamos que queremos realizar una prueba de hipótesis para la media de ‘ingrl’
H0: El ingreso promedio poblacional es mayor e igual a $450
HA: El ingreso promedio poblacional es menor a $450
Recordar que la HA nunca contiene los signos “=” , “≤” o “≥”.
t_prueba2 <- t.test(data$ingrl, mu = 450, alternative = "less", conf.level = 0.95)
print(t_prueba2)
One Sample t-test
data: data$ingrl
t = -111.72, df = 16450, p-value < 2.2e-16
alternative hypothesis: true mean is less than 450
95 percent confidence interval:
-Inf 166.894
sample estimates:
mean of x
162.6633
# Si el valor p es menor a 0.05, rechazamos la hipótesis nula
if (t_prueba$p.value < 0.05) {
cat("Rechazamos la hipótesis nula: Entonces tenemos evidencia de que la media poblacional de 'ingrl' es menor a 450.\n")
} else {
cat("No rechazamos la hipótesis nula: Entonces tenemos evidencia de que la media poblacional de 'ingrl' es mayor e igual a 450.\n")
}Rechazamos la hipótesis nula: Entonces tenemos evidencia de que la media poblacional de 'ingrl' es menor a 450.
Prueba de hipótesis de proporciones
Supongamos que queremos realizar una prueba de hipótesis para proporción de ‘depresion_pp’
H0: La proporción poblacional es igual a 0.10
HA: La proporción poblacional es diferente a 0.10
Recordar que la HA nunca contiene los signos “=” , “≤” o “≥”.
library(haven)
depresion <- as_factor(data$depresion_pp)
prop_test <- prop.test(
sum(depresion == "Si", na.rm = TRUE), # cantidad de "Si"
sum(!is.na(depresion)), # total válido
p = 0.10, # proporción esperada
conf.level = 0.95
)
print(prop_test)
1-sample proportions test with continuity correction
data: sum(depresion == "Si", na.rm = TRUE) out of sum(!is.na(depresion)), null probability 0.1
X-squared = 2640.9, df = 1, p-value < 2.2e-16
alternative hypothesis: true p is not equal to 0.1
95 percent confidence interval:
0.2139329 0.2266579
sample estimates:
p
0.2202298
if (prop_test$p.value < 0.05) {
cat("Rechazamos la hipótesis nula: Entonces tenemos evidencia de que la proporción de 'depresion_pp' es diferente a 0.10.\n")
} else {
cat("No rechazamos la hipótesis nula: Entonces tenemos evidencia de que la proporción de 'depresion_pp' es igual a 0.10.\n")
}Rechazamos la hipótesis nula: Entonces tenemos evidencia de que la proporción de 'depresion_pp' es diferente a 0.10.