Práctica 5. Variables cuantitativas y variables nominales
AE
14/01/2020
Previo
¡Recuerda poner el directorio!
setwd("/Users/anaescoto/Dropbox/2020/R_invierno")Vamos a importar la base completa del cuestionario sociodemográfico de la Encuesta Nacional de Ocupación y Empleo, trimestre III de 2019
library(haven)
SDEMT319 <- read_dta("./datos/SDEMT319.dta")Gráficos de comparación cuanti y cuali
Los gráficos de comparaciones de cajas y distribuciones son muy útiles. Nos permiten ver cómo una variable cuantitativa se comporta entre la población dividida entre los grupos.
Vamos a volver a usar nuestros paquetes ggplot2, dplyr y sjlabelled.
library(ggplot2)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(sjlabelled)##
## Attaching package: 'sjlabelled'## The following object is masked from 'package:dplyr':
##
## as_label## The following objects are masked from 'package:haven':
##
## as_factor, read_sas, read_spss, read_stata, write_sas,
## zap_labelsBoxplot
Empecemos con el “boxplot” univariado y con la función que viene en “graphics”, sin necesidad de instalar paquetería
boxplot(SDEMT319[SDEMT319$clase2==1,]$ing_x_hrs)
La función boxplot, también sirve para obtener el resumen de cinco números.
a<-boxplot(SDEMT319[SDEMT319$clase2==1,]$ing_x_hrs)
a$stats## [,1]
## [1,] 0.00000
## [2,] 0.00000
## [3,] 21.53316
## [4,] 37.03704
## [5,] 92.59259Con formato dplyr
SDEMT319 %>%
filter(clase2==1) %>%
with(boxplot(ing_x_hrs))
Boxplot bivariado
Podemos comparar la distribución de una variable cuantitativa entre las categorías de una variable categórica
#otra forma
boxplot(ing_x_hrs ~ as_label(sex), data= SDEMT319[SDEMT319$clase2==1,])
Con dplyr:
SDEMT319 %>%
filter(clase2==1) %>%
with(boxplot(ing_x_hrs ~ as_label(sex)))
con ggplot
Para “producción”, quizás queremos un gráfico más acabado en ggplot. Como siempre, primero haremos nuestro lienzo, aprovecharemos para quitar el valor atípico para mejor visualización
lienzo <-SDEMT319 %>%
filter(clase2==1 & ing_x_hrs<30000) %>%
ggplot(aes(x=as_label(sex), y=ing_x_hrs))
lienzo
Hoy podemos agregar nuestra geometría
lienzo + geom_boxplot()
g<-lienzo + geom_boxplot(fill='#A4A4A4', color="black")+
theme_classic()
g
Para ponerlo de modeo
g + coord_flip()
Una tercera variable
g + facet_wrap(~ as_label(t_loc))
g2<-lienzo + geom_boxplot(aes(color=as_label(t_loc))) +
theme_classic()
g2
Gráfico de densidad bivariado
Para poder editar la transparencia de nuestro relleno podemos usar la opción alpha y ponerle un valor menor a 1.
lienzo_bi <-SDEMT319 %>%
filter(clase2==1) %>%
ggplot(aes(x=log(ing_x_hrs+1), fill=as_label(sex),
color=as_label(sex),
alpha=I(0.5)))
lienzo_bi + geom_density()
Vamos -sin usar el ggTheme Assistant a ponerle otros elementos a nuestro gráfico.
lienzo_bi + geom_density() +
ggtitle("Logaritmo de los ingresos por trabajo") +
xlab("Log(Dólares mensuales)") + ylab("Densidad") +
labs(fill = "sex", color="sex")
Análisis descriptivo
Un truco para hacer el análisis descriptivo es generar una sub-base con las variables que vamos analizar
mydata<- SDEMT319[SDEMT319$clase2==1, c("ing_x_hrs", "t_loc","eda", "sex", "anios_esc", "fac")]
tail(mydata)## # A tibble: 6 x 6
## ing_x_hrs t_loc eda sex anios_esc fac
## <dbl> <dbl+lbl> <dbl> <dbl+lbl> <dbl> <dbl>
## 1 37.5 4 [Localidades menores de 2 50… 33 1 [Hombr… 7 228
## 2 7.41 4 [Localidades menores de 2 50… 27 1 [Hombr… 9 308
## 3 25.7 4 [Localidades menores de 2 50… 57 1 [Hombr… 12 308
## 4 26.7 4 [Localidades menores de 2 50… 25 2 [Mujer] 10 306
## 5 0 4 [Localidades menores de 2 50… 70 2 [Mujer] 1 190
## 6 0 4 [Localidades menores de 2 50… 15 1 [Hombr… 9 252Si queremos obtener las medias de estas tres variables, usamos la función sapply. Esto nos permite repetir la función para todo el data frame Excluyendo missing values
sapply(as_label(mydata), mean, na.rm=TRUE) # checa los NA's## Warning in mean.default(X[[i]], ...): argument is not numeric or logical:
## returning NA
## Warning in mean.default(X[[i]], ...): argument is not numeric or logical:
## returning NA## ing_x_hrs t_loc eda sex anios_esc fac
## 28.18256 NA 39.48739 NA 10.60863 310.94767Además de media (mean) tenemos las siguientes opciones: mean, sd, var, min, max, median, range, and quantile.
Para obtener de una sola vez las funciones mean,median,25th and 75th quartiles,min,max, usamos la función “summary”
summary(as_label(mydata)) #total## ing_x_hrs t_loc
## Min. : 0.00 Localidades mayores de 100 000 habitantes:109291
## 1st Qu.: 0.00 Localidades de 15 000 a 99 999 habitantes: 23803
## Median : 21.53 Localidades de 2 500 a 14 999 habitantes : 21865
## Mean : 28.18 Localidades menores de 2 500 habitantes : 24333
## 3rd Qu.: 37.04
## Max. :4375.00
## eda sex anios_esc fac
## Min. :12.00 Hombre:106755 Min. : 0.00 Min. : 22.0
## 1st Qu.:28.00 Mujer : 72537 1st Qu.: 9.00 1st Qu.: 110.0
## Median :39.00 Median :10.00 Median : 166.0
## Mean :39.49 Mean :10.61 Mean : 310.9
## 3rd Qu.:50.00 3rd Qu.:12.00 3rd Qu.: 398.0
## Max. :98.00 Max. :99.00 Max. :22671.0summary(as_label(mydata[mydata$sex==1,])) #hombres## ing_x_hrs t_loc
## Min. : 0.00 Localidades mayores de 100 000 habitantes:62873
## 1st Qu.: 0.00 Localidades de 15 000 a 99 999 habitantes:13935
## Median : 22.50 Localidades de 2 500 a 14 999 habitantes :13356
## Mean : 28.57 Localidades menores de 2 500 habitantes :16591
## 3rd Qu.: 37.50
## Max. :4375.00
## eda sex anios_esc fac
## Min. :12.00 Hombre:106755 Min. : 0.00 Min. : 22.0
## 1st Qu.:27.00 Mujer : 0 1st Qu.: 8.00 1st Qu.: 113.0
## Median :38.00 Median : 9.00 Median : 174.0
## Mean :39.51 Mean :10.39 Mean : 318.9
## 3rd Qu.:50.00 3rd Qu.:12.00 3rd Qu.: 411.0
## Max. :98.00 Max. :99.00 Max. :12595.0summary(as_label(mydata[mydata$sex==2,])) #mujeres## ing_x_hrs t_loc
## Min. : 0.00 Localidades mayores de 100 000 habitantes:46418
## 1st Qu.: 0.00 Localidades de 15 000 a 99 999 habitantes: 9868
## Median : 20.41 Localidades de 2 500 a 14 999 habitantes : 8509
## Mean : 27.61 Localidades menores de 2 500 habitantes : 7742
## 3rd Qu.: 35.39
## Max. :2325.58
## eda sex anios_esc fac
## Min. :12.00 Hombre: 0 Min. : 0.00 Min. : 22.0
## 1st Qu.:28.00 Mujer :72537 1st Qu.: 9.00 1st Qu.: 106.0
## Median :39.00 Median :11.00 Median : 157.0
## Mean :39.45 Mean :10.93 Mean : 299.2
## 3rd Qu.:49.00 3rd Qu.:14.00 3rd Qu.: 375.0
## Max. :98.00 Max. :99.00 Max. :22671.0Tapply, se parece a sapply, nos permite hacer una función en un subconjunto que se deriva de un vector. Por ejemplo, a continuación podemos pedirle al programa que nos haga un summary del ingreso, tomando en cuenta un vector de factor
tapply(mydata$ing_x_hrs, as_label(mydata$sex), summary)## $Hombre
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 0.00 22.50 28.57 37.50 4375.00
##
## $Mujer
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 0.00 20.41 27.61 35.39 2325.58tapply(mydata$ing_x_hrs, as_label(mydata$sex), mean, na.rm=T)## Hombre Mujer
## 28.57272 27.60834Con dplyr
Con dplyr a veces es más sencillo por la función “as group” y summarise (o summarize)
SDEMT319 %>%
filter(SDEMT319$clase2==1 )%>%
summarise(avg_ing= mean(ing_x_hrs))## # A tibble: 1 x 1
## avg_ing
## <dbl>
## 1 28.2Con la función de grupos es más sencillo que tapply y sapply
SDEMT319 %>%
filter(SDEMT319$clase2==1 )%>%
group_by(as_label(sex)) %>%
summarise(avg_ing = mean(ing_x_hrs))## # A tibble: 2 x 2
## `as_label(sex)` avg_ing
## <fct> <dbl>
## 1 Hombre 28.6
## 2 Mujer 27.6SDEMT319 %>%
filter(SDEMT319$clase2==1 )%>%
group_by(as_label(sex)) %>%
tally() # Nos da el conteo de casos## # A tibble: 2 x 2
## `as_label(sex)` n
## <fct> <int>
## 1 Hombre 106755
## 2 Mujer 72537Y podemos “summarizar” los indicadores que querramos. Recuerda que estos resultados se pueden guardar en un objeto.
resultados1<-SDEMT319 %>%
filter(SDEMT319$clase2==1 )%>%
group_by(as_label(sex)) %>%
summarise(avg_ing = mean(ing_x_hrs),
avg_esc = mean(anios_esc))
resultados1## # A tibble: 2 x 3
## `as_label(sex)` avg_ing avg_esc
## <fct> <dbl> <dbl>
## 1 Hombre 28.6 10.4
## 2 Mujer 27.6 10.9Usando datos expandidos
Un elemento importante, sobre todo cuando estamos usando una encuesta con representatividad nacional como la ENOE es que podemos establecer los valores poblacionales. Para ello, los descriptivos se presentan con los datos de la población y no de la muestra. Esto se hace expandiendo los datos muestrales según el peso o factor de espansión.
Para tabulados podemos usar el comando tally
SDEMT319 %>%
filter(SDEMT319$clase2==1 )%>%
group_by(as_label(sex)) %>%
tally(fac) #SUMA TODOS LOS RENGLONES DE LA VARIABLE fac## # A tibble: 2 x 2
## `as_label(sex)` n
## <fct> <dbl>
## 1 Hombre 34046315
## 2 Mujer 21704114Para sacar medias, podemos usar la función de media expandida
SDEMT319 %>%
filter(SDEMT319$clase2==1) %>%
group_by(as_label(sex)) %>%
summarise(avg_ing =weighted.mean(ing_x_hrs, weights = fac,na.rm = T))## # A tibble: 2 x 2
## `as_label(sex)` avg_ing
## <fct> <dbl>
## 1 Hombre 28.6
## 2 Mujer 27.6También existen librerías para hacer tabulados de manera más tradicional. Es más fácil cuando queremos calcular porcentajes y proporciones. Usaremos la librería “questionr”, instalada en prácticas anteriores.
library(questionr)Para usar nuestro factor de expansión primero tenemos que verificar qué tipo de valor es en nuestro objeto.
class(mydata$fac)## [1] "numeric"Para variables cualitativas podemos hacer tablas con los factores de expansión
wtd.table(mydata$sex, weights = mydata$fac)## 1 2
## 34046315 21704114prop.table(wtd.table(mydata$sex, weights = mydata$fac))## 1 2
## 0.6106915 0.3893085addmargins(prop.table(wtd.table(mydata$sex, weights = mydata$fac)))*100## 1 2 Sum
## 61.06915 38.93085 100.00000Para una tabla cruzada
wtd.table(mydata$sex,mydata$t_loc, weights = mydata$fac)## 1 2 3 4
## 1 16152047 4896720 4876237 8121311
## 2 11600634 3453422 3003438 3646620Mi primera función
Unos de los elementos más poderosos de R es hacer nuestra propias funciones. Y cuando trabajamos bases de datos a veces hacemos un procedimiento que lo queremos repetir. Tal vez para una base de datos, tal vez para una variable, tal vez para grupos de casos.
Para ello haremos una función sencilla. Para sumarle un valor un 1
mi_funcion<-function(x) {
x<-x+1
x
}
mi_funcion(5)## [1] 6Qué tal hoy una una función donde podamos establecer qué valor sumamos
mi_funcion<-function(x,a) {
x<-x+a
x
}
mi_funcion(5,2)## [1] 7Mi función para expandir medias
expandir<- function(x){
weighted.mean(x, w = SDEMT319$fac)
}Entonces puedo ya sólo poner la variable
expandir(SDEMT319$eda)## [1] NAOtras funciones
percentage1<- function(x) {
addmargins(prop.table(table(x)))*100
}
percentage2<- function(x, y) {
addmargins(prop.table(table(x,y)))*100
}
percent_col<- function(x, y) {
addmargins(prop.table(table(x,y),2))*100}
percent_row<- function(x, y) {
addmargins(prop.table(table(x,y),1))*100
}Guardando todo nuestro Environment
save.image("ambiente_p5.RData")Ejercicio
- Tome la variable “t_loc” y estime la media y la mediana(funcion “median”) de la variable “ingocup” (Ingreso mensual), pero realícelo sólo para los jefes de hogar (ponga el filtro de par_c==101)
- Realice estos resultados también agrupando la variable por sexo, para revisar si los jefes ganan más que las jefas, pero también añada el promedio de miembros del hogar. ¿Qué opina de los ingresos en los hogares?
- Realice las operaciones correspondientes ¿Cuántos hogares (con estimaciones para la población total) tienen jefas?
- Compare la distribución de los ingresos por hora según sexo del trabajador
Envíe el código https://www.dropbox.com/request/eSRRn257s7g8R8R4OHoW