En esta parte del curso examinaremos las siguientes herramientas de inferencia estadística

• Cálculo de intervalo de confianza
• Representación gráfica de intervalos de confianza

Base de datos para estos ejercicios: Familia y roles de género 2012, a descargar de:

http://iop-data.pucp.edu.pe/busqueda/encuesta/71?

Se sugiere descargar también el cuestionario para utilizarlo como referencia de libro de códigos. Descomprimir y grabar el archivo SPSS en el directorio de trabajo de R

# Importar la base de datos del SPSS a un data frame de R
library(foreign)
genero <- as.data.frame(read.spss("IOP_1212_01_B.sav"))

## re-encoding from UTF-8

En estos ejemplos trabajaremos con las siguientes variables:

• Edad ideal para que una mujer se case (P1)
• Edad ideal para que un hombre se case (P2)
• Edad en la que el entrevistado se casó o empezó a convivir (P23)
• Número de hijos (P44A)
• Horas dedicadas al trabajo doméstico (P19A)
• Horas que la pareja le dedica al trabajo doméstico (P28A)

Antes de utilizar estas variables es necesario evaluar si necesitan algún tipo de acondicionamiento (identificar valores perdidos o "raros", recodificar, etc)

Las siguientes transformaciones nos permitirán acondicionar las variables para utilizarlas en el análisis:

genero$p1r <- genero$P1
genero$p2r <- genero$P2
genero$p1r[genero$P1 > 40] <- NA
genero$p2r[genero$P2 == 99] <- NA

genero$p19ar <- genero$P19A
genero$p19ar[genero$P19A >= 140] <- NA

genero$p28ar <- genero$P28A
genero$p28ar[genero$P28A > 120] <- NA

Para calcular el intervalo de confianza de la media de una variable que proviene de una muestra simple al azar usamos la siguiente fórmula: \[ \begin{aligned} IC_{1-\alpha}=\bar{x} \pm Z\sigma_{\bar{x}} \end{aligned} \]

Donde:

• \(1-\alpha\) es el nivel de confianza deseado
• Z es el valor de la distribución normal estándar que representa el nivel de confianza
• \(\sigma_{\bar{x}}\) es el error estándar de la media que se obtiene mediante:

\[ \begin{aligned} \sigma_{\bar{x}} = \frac{S_{x}}{\sqrt{n}} \end{aligned} \]

Entonces para calcular el intervalo de confianza para una media podemos proceder de la siguiente manera:

media <- mean(na.omit(genero$p19ar))
desv <- sd(na.omit(genero$p19ar))
N <- length(na.omit(genero$p19ar))
error.est <- desv/sqrt(N)
error <- 2*error.est # Fijamos Z=2 para indicar un nivel de confianza de 95%
lim.inf <- media-error
lim.sup <- media+error
resultado1 <- data.frame(media, desv, N, error.est, error, lim.inf, lim.sup)
resultado1

##   media  desv    N error.est error lim.inf lim.sup
## 1 22.61 19.95 1196    0.5769 1.154   21.46   23.76

media <- mean(na.omit(genero$p19ar))
desv <- sd(na.omit(genero$p19ar))
N <- length(na.omit(genero$p19ar))
error.est <- desv/sqrt(N)
error <- qt(0.975, df= N-1) * error.est # Usar el cuantil 0.975 de t
lim.inf <- media-error
lim.sup <- media+error
resultado2 <- data.frame(media, desv, N, error.est, error, lim.inf, lim.sup)
resultado2

##   media  desv    N error.est error lim.inf lim.sup
## 1 22.61 19.95 1196    0.5769 1.132   21.48   23.74

resultado1

##   media  desv    N error.est error lim.inf lim.sup
## 1 22.61 19.95 1196    0.5769 1.154   21.46   23.76

Generamos la función:

int.conf <- function(x, ic = 95) {
  y <- na.omit(x)
  Media <- mean(y)
  Desv.Est <- sd(y)
  N <- length(y)
  Error.Est <- Desv.Est/sqrt(N)
  ci.y <- 1-(((100-ic)/100)/2)
  Error <- Error.Est * qt(ci.y, df = N-1)
  Lim.Inf <- Media - Error
  Lim.Sup <- Media + Error
  result <- data.frame(Media, Desv.Est, N, Error.Est, Error, Lim.Inf, Lim.Sup)
  return(result)
  }

Usamos la función para calcular un intervalo de confianza al 95% (ic por defecto de la función)

int.conf(genero$p19ar)

##   Media Desv.Est    N Error.Est Error Lim.Inf Lim.Sup
## 1 22.61    19.95 1196    0.5769 1.132   21.48   23.74

En vez de un intervalo de 95% podemos usar otros niveles de confianza:

int.conf(genero$p19ar, ic=99) # Intervalo de 99% de confianza

##   Media Desv.Est    N Error.Est Error Lim.Inf Lim.Sup
## 1 22.61    19.95 1196    0.5769 1.488   21.12    24.1

int.conf(genero$p19ar, ic=90) # Intervalo de 90% de confianza

##   Media Desv.Est    N Error.Est  Error Lim.Inf Lim.Sup
## 1 22.61    19.95 1196    0.5769 0.9496   21.66   23.56

También podemos calcular un intervalo de confianza para una proporción a partir de una muestra grande.

\[ \begin{aligned} IC_{1-\alpha}= p \pm Z\sigma_{p} \end{aligned} \]

Donde:

• \(1-\alpha\) es el nivel de confianza deseado
• Z es el valor de la distribución normal estándar que representa el nivel de confianza
• \(\sigma_{p}\) es el error estándar de la propoción que se obtiene mediante:

\[ \begin{aligned} \sigma_{p} = \sqrt{\frac{pq}{n}} \end{aligned} \]

Donde:

• p es la proporción de la categoría de análisis en nuestra muestra
• q es igual a 1 - p

Vamos a usar la variable P32

library(descr)
freq(genero$P32, plot = FALSE)

## genero$P32 
##                                           Frequency Percent Valid Percent
## Hago mucho más de lo que me corresponde         147  12.219        21.304
## Hago algo más de lo que me corresponde          141  11.721        20.435
## Hago más o menos lo que me corresponde          256  21.280        37.101
## Hago algo menos de lo que me corresponde         65   5.403         9.420
## Hago mucho menos de lo que me corresponde        53   4.406         7.681
## No contesta                                      28   2.328         4.058
## NA's                                            513  42.643              
## Total                                          1203 100.000       100.000

Preparamos los datos para el cálculo del intervalo de confianza:

p32r <- na.omit(genero$P32)
cat <- ifelse(p32r=="Hago más o menos lo que me corresponde", 1, 0)
prop.table(table(cat))

## cat
##     0     1 
## 0.629 0.371

p <- mean(cat)
p

## [1] 0.371

Preparamos el resto de la información

n <- length(cat)
error.est.p <- sqrt((p*(1-p))/n)
error.p <- 2 * error.est.p # Usamos Z = 2 para indicar un nivel de confianza del 95%
lim.inf.p <- p - error.p
lim.sup.p <- p + error.p
result.p <- data.frame(p, n, error.est.p, error.p, lim.inf.p, lim.sup.p)
result.p

##       p   n error.est.p error.p lim.inf.p lim.sup.p
## 1 0.371 690     0.01839 0.03678    0.3342    0.4078

Generamos la función

int.conf.p <- function(x, cat, ic = 95){
  vcat <- na.omit(x)
  p <- mean(ifelse(vcat==cat, 1, 0))
  n <- length(vcat)
  error.est.p <- sqrt((p*(1-p))/n)
  beta <- 1-(((100-ic)/100)/2)
  error <- error.est.p * qt(beta, df = n-1)
  lim.inf.p <- p - error
  lim.sup.p <- p + error
  result.p <- data.frame(p, n, error.est.p, error, lim.inf.p, lim.sup.p)
  return(result.p)
}

Usamos la función

int.conf.p(genero$P32, "Hago más o menos lo que me corresponde")# IC95%

##       p   n error.est.p   error lim.inf.p lim.sup.p
## 1 0.371 690     0.01839 0.03611    0.3349    0.4071

int.conf.p(genero$P32, "Hago más o menos lo que me corresponde", ic=99)# IC99%

##       p   n error.est.p  error lim.inf.p lim.sup.p
## 1 0.371 690     0.01839 0.0475    0.3235    0.4185

int.conf.p(genero$P32, "Hago más o menos lo que me corresponde", ic=90)# IC90%

##       p   n error.est.p   error lim.inf.p lim.sup.p
## 1 0.371 690     0.01839 0.03029    0.3407    0.4013

freq(genero$P10, plot=FALSE)

## genero$P10 
##             Frequency   Percent
## Si               1070  88.94431
## No                132  10.97257
## No contesta         1   0.08313
## Total            1203 100.00000

int.conf.p(genero$P10, "Si")

##        p    n error.est.p   error lim.inf.p lim.sup.p
## 1 0.8894 1203    0.009041 0.01774    0.8717    0.9072

Podemos crear un gráfico del intervalo de confianza de la media de una variable. En este ejemplo usaremos el gráfico de puntos del paquete ggplot2. Vamos a representar gráficamente el promedio de horas semanales dedicadas a labores domésticas que le dedican los hombres y las mujeres.

Primero vamos a generar un data frame que contenga los estadísticos que necesitamos graficar (la media, los límites inferior y superior del intervalo de confianza). Para ello usaremos la función "summarySE" del paquete "Rmisc" (que debe instalarse en su sistema). Esta función summaySE hace algo parecido a la función que hemos creado en una de las diapositivas anteriores.

library(Rmisc)

## Loading required package: lattice
## Loading required package: plyr

df <- summarySE(genero, measurevar="p19ar", groupvars="SEXO", na.rm=T)
df

##        SEXO   N p19ar    sd     se    ci
## 1 Masculino 587  13.7 13.50 0.5574 1.095
## 2  Femenino 609  31.2 21.35 0.8653 1.699

Un gráfico de puntos que muestra las medias de p19ar según sexo

library(ggplot2)
graf.punto1 <- ggplot(df, aes(x=SEXO, y=p19ar)) + 
  geom_point(size = 3) + ylim(0, 40)

graf.punto1

Ahora añadimos las barras que indican el intervalo de confianza, con títulos

graf.punto1 +  geom_errorbar(aes(ymin=p19ar-ci, ymax=p19ar+ci), width = 0.2) + 
  ylab("Horas semanales") + 
  ggtitle("Intervalo de confianza al 95% para la media de las horas semanales\n dedicadas al trabajo doméstico, según sexo")

Podemos graficar el mismo concepto con barras:

graf.barra1 <- ggplot(df, aes(x=SEXO, y=p19ar)) + 
  geom_bar(width=0.5, fill="white", colour="black", stat="identity") + 
  geom_errorbar(aes(ymin=p19ar-ci, ymax=p19ar+ci), width = 0.2) +
  ylab("Horas semanales") + ggtitle("Intervalo de confianza al 95% para la media de las horas semanales\n dedicadas al trabajo doméstico, según sexo")

graf.barra1

Según ámbito y sexo del entrevistado

df2 <- summarySE(genero, measurevar="p19ar", groupvars=c("SEXO", "Ambito"), 
                 na.rm=T)
df2

##        SEXO          Ambito   N p19ar    sd     se    ci
## 1 Masculino     Lima-Callao 218 14.49 14.27 0.9666 1.905
## 2 Masculino Interior Urbano 246 12.68 12.26 0.7819 1.540
## 3 Masculino  Interior Rural 123 14.34 14.43 1.3007 2.575
## 4  Femenino     Lima-Callao 229 30.16 19.83 1.3107 2.583
## 5  Femenino Interior Urbano 264 29.62 19.59 1.2057 2.374
## 6  Femenino  Interior Rural 116 36.83 26.71 2.4795 4.912

graf.punto2 <- ggplot(df2, aes(x=Ambito, y=p19ar, colour=SEXO)) +
  geom_point(size = 3) + ylim(0, 50) + 
  geom_errorbar(aes(ymin=p19ar-ci, ymax=p19ar+ci), width = 0.2) + 
  ylab("Horas semanales") + 
  ggtitle("Intervalo de confianza al 95% para la media de las horas semanales\n dedicadas al trabajo doméstico, según ámbito, por sexo")

graf.punto2

graf.bar2 <- ggplot(df2, aes(x=Ambito, y=p19ar, fill=SEXO)) +
  geom_bar(position="dodge", stat = "identity") +
  geom_errorbar(aes(ymin=p19ar-ci, ymax=p19ar+ci), width = 0.2, 
  position=position_dodge(.9)) + ylab("Horas semanales") + 
  ggtitle("Intervalo de confianza al 95% para la media de las horas semanales\n dedicadas al trabajo doméstico, según ámbito, por sexo")

graf.bar2