La base de datos puede descargarse desde:
https://www.dropbox.com/s/rfyc4guf1hha130/wvs.rdata?dl=0
Luego de descargar el archivo, éste debe copiarse en el directorio de trabajo
load("wvs.rdata")
library(car)
# Sexo del entrevistado
wvs$sex <- wvs$X001
wvs$sex[wvs$sex < 1] <- NA
wvs$sex <- factor(wvs$sex)
levels(wvs$sex) <- c("Masculino", "Femenino")
table(wvs$sex)
##
## Masculino Femenino
## 160465 172298
# Grupos de edad
table(wvs$X003R2)
##
## -5 -4 -3 -2 -1 1 2 3 13 14
## 124 3437 23 440 183 101279 135649 96342 2 12
wvs$gedad <- recode(wvs$X003R2, "lo:0= NA; 4:hi=NA")
wvs$gedad <- factor(wvs$gedad)
levels(wvs$gedad) <- c("15-29", "30-49", "50+")
table(wvs$gedad)
##
## 15-29 30-49 50+
## 101279 135649 96342
# Clase social subjetiva
table(wvs$X045)
##
## -5 -4 -3 -2 -1 1 2 3 4 5
## 239 45770 3 2584 8101 5204 54140 105332 78649 37469
wvs$classoc <- recode(wvs$X045, "1:2=1; 3=2; 4=3; 5=4; else=NA")
wvs$classoc <- factor(wvs$classoc)
levels(wvs$classoc) <- c("Alta/Media alta", "Media baja", "Trabajadora", "Baja")
table(wvs$classoc)
##
## Alta/Media alta Media baja Trabajadora Baja
## 59344 105332 78649 37469
# Nivel Educativo
table(wvs$X025)
##
## -5 -4 -3 -2 -1 1 2 3 4 5 6 7
## 461 25415 17125 1442 656 23722 42073 22637 55033 26742 51786 22730
## 8
## 47668
wvs$nivedu <- factor(recode(wvs$X025, "1:3=1; 4:5=2; 6=3; 7=4; 8=5; else=NA"))
levels(wvs$nivedu) <- c("Secundaria incompleta o menos", "Secundaria completa", "Técnica completa", "Univ. incompleta", "Univ. completa")
# Satisfacción con la vida
wvs$satvid <- wvs$A170
wvs$satvid[wvs$satvid < 0] <- NA
table(wvs$satvid)
##
## 1 2 3 4 5 6 7 8 9 10
## 14566 10067 16817 19263 47067 36679 49091 60399 33873 44275
# Felicidad
table(wvs$oleada, wvs$A008)
##
## -5 -4 -2 -1 1 2 3 4
## 1981-1984 0 970 230 0 2155 5497 1274 181
## 1990-1994 0 0 60 529 5304 12333 5527 805
## 1995-1998 0 3029 30 877 17531 37193 14570 2379
## 1999-2004 0 0 54 517 17024 28599 8788 1788
## 2005-2009 21 0 269 588 23326 44742 12593 2436
## 2010-2014 4 0 237 491 27906 44051 11041 2542
wvs$feliz1 <- wvs$A008
wvs$feliz1[wvs$feliz1 < 0] <- NA
wvs$feliz2 <- recode(wvs$feliz1, "1:2=1; 3:4=2")
wvs$feliz1 <- factor(wvs$feliz1)
wvs$feliz2 <- factor(wvs$feliz2)
levels(wvs$feliz1) <- c("Muy feliz", "Bastante feliz", "Poco feliz", "Nada Feliz")
levels(wvs$feliz2) <- c("Felices", "Poco/Nada felices")
table(wvs$feliz1)
##
## Muy feliz Bastante feliz Poco feliz Nada Feliz
## 93246 172415 53793 10131
table(wvs$feliz2)
##
## Felices Poco/Nada felices
## 265661 63924
# Asistencia a Servicios religiosos
wvs$serv.rel <- factor(recode(wvs$F028, "1:2 = 1; 3=2; 4:5=3; 6:8=4; else=NA"))
levels(wvs$serv.rel) <- c("1 vez x semana o más", "1 vez x mes", "En ocasiones especiales", "1 vez x año o menos")
table(wvs$serv.rel)
##
## 1 vez x semana o más 1 vez x mes En ocasiones especiales
## 105808 35223 50672
## 1 vez x año o menos
## 122364
# Se considera una persona religiosa
wvs$perso.rel <- factor(recode(wvs$F034, "1=1; 2:3=2; else=NA"))
levels(wvs$perso.rel) <- c("Es una persona religiosa", "No es una persona religiosa o es ateo")
table(wvs$perso.rel)
##
## Es una persona religiosa
## 216543
## No es una persona religiosa o es ateo
## 90948
# Importancia de Dios en la vida
wvs$imp.dios <- wvs$F063
wvs$imp.dios[wvs$imp.dios < 1] <- NA
table(wvs$imp.dios)
##
## 1 2 3 4 5 6 7 8 9 10
## 24590 10014 10230 7506 20513 15488 20278 24956 21720 160075
Base de datos de Perú
wvs.peru <- subset(wvs, S003 == 604)
# Se convierte nuevamente en factor la variable "oleada" para eliminar las categorías vacías
wvs.peru$oleada <- factor(wvs.peru$oleada)
# Subconjunto de datos para el año 2012
peru.12 <- subset(wvs.peru, oleada=="2010-2014")
En esta sección se trabaja con la función “tabla.cont” que puede descargarse desde el siguiente enlace:
https://sites.google.com/a/pucp.pe/data_est/archivos/Funcion_tabla_cont.R?attredirects=0&d=1
Deben seguirse las instruccione.s para cargar la función tabla.cont antes de proceder con los pasos siguientes
Un tutorial que explica el manejo de esta función puede verse en este video de Youtube:
https://youtu.be/7-sDCocwjII?list=PL5BNuzISw7jxg9mx6t_Xt9L76Sb87mp17
Funcion para generar tablas de contingencia. Requiere los paquetes “vcd”, “vcdExtra” y “data.table”
library(vcd)
## Loading required package: grid
library(vcdExtra)
## Loading required package: gnm
##
## Attaching package: 'vcdExtra'
## The following object is masked from 'package:car':
##
## Burt
library(data.table)
tabla.cont <- function(df, x, y, pc="s", asoc = FALSE){
library(vcd)
library(vcdExtra)
library(data.table)
as <- c()
gam <- c()
misvars <- c(y, x)
data <- na.omit(df[misvars])
tab.0 <- table(data)
asoc.0 <- assocstats(tab.0)
gam.0 <- GKgamma(tab.0)
if(pc=="s"){
tab.1 <- addmargins(tab.0, 1, FUN = list(list(TOTALc = sum)))
tab.1 <- addmargins(tab.1, 2, FUN = list(list(TOTALf = sum)))
if(asoc==TRUE){
as <- asoc.0
gam <- gam.0}
return(list(tab.1, as, gam))}
if (pc=="col"){
tab.2 <- addmargins(tab.0, 2, FUN = list(list(TOTAL = sum)))
tab.2 <- round(prop.table(tab.2, 2)*100,2)
tab.2 <- round(addmargins(tab.2, 1, FUN = list(list(TOTAL = sum))),1)
tab.0 <- tab.0
tab.0 <- addmargins(tab.0, 2, FUN = list(list(TOTAL = sum)))
tab.0 <- addmargins(tab.0, 1, FUN = list(list(Nvalid = sum)))
tab.2 <- rbind(tab.2, tail(tab.0, 1))
if(asoc==TRUE){
as <- asoc.0
gam <- gam.0
}
return(list(tab.2, as, gam))}
if (pc=="fila"){
tab.2 <- addmargins(tab.0, 1, FUN = list(list(TOTAL = sum)))
tab.2 <- round(prop.table(tab.2, 1)*100,2)
tab.2 <- round(addmargins(tab.2, 2, FUN = list(list(TOTAL = sum))),1)
tab.0 <- tab.0
tab.0.1 <- addmargins(tab.0, 1, FUN = list(list(TOTAL = sum)))
tab.0.1 <- addmargins(tab.0.1, 2, FUN = list(list(Nvalid = sum)))
tab.2 <- cbind(tab.2, Nvalid=tab.0.1[, "Nvalid"])
if(asoc==TRUE){
as <- asoc.0
gam <- gam.0
}
return(list(tab.2, as, gam))}
}
Las variables sociodemográficas que estamos considerando como variables independientes son: sexo, grupo de edad, clase social subjetiva y nivel educativo del entrevistado.
Paso 1: Se generan las tablas de contingencia correspondientes usando la función “tabla.cont”. En este caso estamos colocando las categorías de la variable dependiente en las columnas y las categorías de las variables independientes en las filas, por ello se pide que se calculen los % para cada fila (% horizontales). Se pide además los estadísticos de asociación
t1.1 <- tabla.cont(peru.12, "perso.rel","sex", pc="fila", asoc=T)
t1.2 <- tabla.cont(peru.12, "perso.rel","gedad", pc="fila", asoc=T)
t1.3 <- tabla.cont(peru.12, "perso.rel","classoc", pc="fila", asoc=T)
t1.4 <- tabla.cont(peru.12, "perso.rel", "nivedu", pc="fila", asoc=T)
Paso 2: Observamos los resultados de cada tabla y sus respectivos estadísticos de asociación
t1.1
## [[1]]
## Es una persona religiosa No es una persona religiosa o es ateo
## Masculino 79.6 20.4
## Femenino 89.2 10.8
## TOTAL 84.3 15.7
## TOTAL Nvalid
## Masculino 100 587
## Femenino 100 581
## TOTAL 100 1168
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 20.664 1 5.4738e-06
## Pearson 20.364 1 6.4009e-06
##
## Phi-Coefficient : 0.132
## Contingency Coeff.: 0.131
## Cramer's V : 0.132
##
## [[3]]
## gamma : -0.357
## std. error : 0.073
## CI : -0.501 -0.214
t1.2
## [[1]]
## Es una persona religiosa No es una persona religiosa o es ateo TOTAL
## 15-29 78.5 21.5 100
## 30-49 84.9 15.1 100
## 50+ 91.0 9.0 100
## TOTAL 84.3 15.7 100
## Nvalid
## 15-29 391
## 30-49 477
## 50+ 300
## TOTAL 1168
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 20.892 2 2.9064e-05
## Pearson 20.222 2 4.0640e-05
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.13
## Cramer's V : 0.132
##
## [[3]]
## gamma : -0.299
## std. error : 0.062
## CI : -0.421 -0.177
t1.3
## [[1]]
## Es una persona religiosa
## Alta/Media alta 79.1
## Media baja 85.2
## Trabajadora 85.6
## Baja 86.9
## TOTAL 84.2
## No es una persona religiosa o es ateo TOTAL Nvalid
## Alta/Media alta 20.9 100 239
## Media baja 14.8 100 404
## Trabajadora 14.4 100 368
## Baja 13.1 100 122
## TOTAL 15.8 100 1133
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 5.8823 3 0.11748
## Pearson 6.1857 3 0.10292
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.074
## Cramer's V : 0.074
##
## [[3]]
## gamma : -0.131
## std. error : 0.064
## CI : -0.255 -0.006
t1.4
## [[1]]
## Es una persona religiosa
## Secundaria incompleta o menos 89.2
## Secundaria completa 82.1
## Técnica completa 83.2
## Univ. incompleta 81.4
## Univ. completa 82.0
## TOTAL 84.1
## No es una persona religiosa o es ateo TOTAL
## Secundaria incompleta o menos 10.8 100
## Secundaria completa 17.9 100
## Técnica completa 16.8 100
## Univ. incompleta 18.6 100
## Univ. completa 18.0 100
## TOTAL 15.9 100
## Nvalid
## Secundaria incompleta o menos 306
## Secundaria completa 392
## Técnica completa 161
## Univ. incompleta 129
## Univ. completa 150
## TOTAL 1138
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 8.9502 4 0.062357
## Pearson 8.4010 4 0.077946
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.086
## Cramer's V : 0.086
##
## [[3]]
## gamma : 0.137
## std. error : 0.056
## CI : 0.028 0.247
Paso 3: Utilizamos la función xtable para exportar las tablas archivos en formato html que luego podemos abrir con el Excel para diseñar las tablas en un formato publicable. Los archivos generados se guardarán en el directorio de trabajo.
library(xtable)
print(xtable(t1.1[[1]]), type = "html", file="tab1_1.html")
print(xtable(t1.2[[1]]), type = "html", file="tab1_2.html")
print(xtable(t1.3[[1]]), type = "html", file="tab1_3.html")
print(xtable(t1.4[[1]]), type = "html", file="tab1_4.html")
Los estadísticos de asociación (que pueden verse en el paso 2) se copian manualmente al Excel para mostrarlos en tablas publicables.
Paso 1:
t2.1 <- tabla.cont(peru.12, "serv.rel","sex", pc="fila", asoc=T)
t2.2 <- tabla.cont(peru.12, "serv.rel","gedad", pc="fila", asoc=T)
t2.3 <- tabla.cont(peru.12, "serv.rel","classoc", pc="fila", asoc=T)
t2.4 <- tabla.cont(peru.12, "serv.rel", "nivedu", pc="fila", asoc=T)
Paso 2:
t2.1
## [[1]]
## 1 vez x semana o más 1 vez x mes En ocasiones especiales
## Masculino 33.7 21.6 2.0
## Femenino 48.7 17.6 1.9
## TOTAL 41.2 19.6 2.0
## 1 vez x año o menos TOTAL Nvalid
## Masculino 42.6 100 584
## Femenino 31.8 100 585
## TOTAL 37.2 100 1169
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 27.670 3 4.2605e-06
## Pearson 27.543 3 4.5290e-06
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.152
## Cramer's V : 0.153
##
## [[3]]
## gamma : -0.238
## std. error : 0.046
## CI : -0.328 -0.148
t2.2
## [[1]]
## 1 vez x semana o más 1 vez x mes En ocasiones especiales
## 15-29 32.3 19.1 2.5
## 30-49 42.0 20.7 2.1
## 50+ 51.9 18.5 1.0
## TOTAL 41.2 19.6 2.0
## 1 vez x año o menos TOTAL Nvalid
## 15-29 46.1 100 393
## 30-49 35.3 100 479
## 50+ 28.6 100 297
## TOTAL 37.2 100 1169
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 33.165 6 9.7463e-06
## Pearson 32.990 6 1.0532e-05
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.166
## Cramer's V : 0.119
##
## [[3]]
## gamma : -0.221
## std. error : 0.039
## CI : -0.297 -0.146
t2.3
## [[1]]
## 1 vez x semana o más 1 vez x mes En ocasiones especiales
## Alta/Media alta 41.6 22.0 2.5
## Media baja 42.9 17.0 2.0
## Trabajadora 36.0 21.7 1.9
## Baja 50.4 15.4 1.6
## TOTAL 41.2 19.4 2.0
## 1 vez x año o menos TOTAL Nvalid
## Alta/Media alta 33.9 100 245
## Media baja 38.1 100 399
## Trabajadora 40.4 100 369
## Baja 32.5 100 123
## TOTAL 37.3 100 1136
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 12.075 9 0.20914
## Pearson 12.051 9 0.21043
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.102
## Cramer's V : 0.059
##
## [[3]]
## gamma : 0.014
## std. error : 0.038
## CI : -0.06 0.088
t2.4
## [[1]]
## 1 vez x semana o más 1 vez x mes
## Secundaria incompleta o menos 46.0 19.3
## Secundaria completa 39.1 19.9
## Técnica completa 40.9 20.1
## Univ. incompleta 32.8 21.1
## Univ. completa 40.4 19.9
## TOTAL 40.6 19.9
## En ocasiones especiales 1 vez x año o menos
## Secundaria incompleta o menos 0.7 34.0
## Secundaria completa 2.0 38.9
## Técnica completa 1.8 37.2
## Univ. incompleta 4.7 41.4
## Univ. completa 2.6 37.1
## TOTAL 2.0 37.4
## TOTAL Nvalid
## Secundaria incompleta o menos 100 300
## Secundaria completa 100 396
## Técnica completa 100 164
## Univ. incompleta 100 128
## Univ. completa 100 151
## TOTAL 100 1139
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 13.452 12 0.33707
## Pearson 13.690 12 0.32093
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.109
## Cramer's V : 0.063
##
## [[3]]
## gamma : 0.062
## std. error : 0.036
## CI : -0.008 0.132
Paso 3:
print(xtable(t2.1[[1]]), type = "html", file="tab2_1.html")
print(xtable(t2.2[[1]]), type = "html", file="tab2_2.html")
print(xtable(t2.3[[1]]), type = "html", file="tab2_3.html")
print(xtable(t2.4[[1]]), type = "html", file="tab2_4.html")
Además de las variables sociodemográficas, en este grupo de tablas se incluyen dos variables sobre religiosidad
Paso 1:
t3.1 <- tabla.cont(peru.12, "feliz2","sex", pc="fila", asoc=T)
t3.2 <- tabla.cont(peru.12, "feliz2","gedad", pc="fila", asoc=T)
t3.3 <- tabla.cont(peru.12, "feliz2","classoc", pc="fila", asoc=T)
t3.4 <- tabla.cont(peru.12, "feliz2", "nivedu", pc="fila", asoc=T)
t3.5 <- tabla.cont(peru.12, "feliz2", "perso.rel", pc="fila", asoc=T)
t3.6 <- tabla.cont(peru.12, "feliz2", "serv.rel", pc="fila", asoc=T)
Paso 2:
t3.1
## [[1]]
## Felices Poco/Nada felices TOTAL Nvalid
## Masculino 77.1 22.9 100 603
## Femenino 75.6 24.4 100 599
## TOTAL 76.4 23.6 100 1202
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 0.36894 1 0.54358
## Pearson 0.36891 1 0.54360
##
## Phi-Coefficient : 0.018
## Contingency Coeff.: 0.018
## Cramer's V : 0.018
##
## [[3]]
## gamma : 0.041
## std. error : 0.068
## CI : -0.092 0.174
t3.2
## [[1]]
## Felices Poco/Nada felices TOTAL Nvalid
## 15-29 82.4 17.6 100 410
## 30-49 76.0 24.0 100 492
## 50+ 68.7 31.3 100 300
## TOTAL 76.4 23.6 100 1202
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 18.223 2 0.00011038
## Pearson 18.269 2 0.00010789
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.122
## Cramer's V : 0.123
##
## [[3]]
## gamma : 0.236
## std. error : 0.053
## CI : 0.131 0.34
t3.3
## [[1]]
## Felices Poco/Nada felices TOTAL Nvalid
## Alta/Media alta 88.6 11.4 100 246
## Media baja 78.0 22.0 100 414
## Trabajadora 70.0 30.0 100 383
## Baja 68.5 31.4 100 124
## TOTAL 76.6 23.4 100 1167
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 36.692 3 5.3461e-08
## Pearson 34.160 3 1.8327e-07
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.169
## Cramer's V : 0.171
##
## [[3]]
## gamma : 0.303
## std. error : 0.048
## CI : 0.208 0.397
t3.4
## [[1]]
## Felices Poco/Nada felices TOTAL Nvalid
## Secundaria incompleta o menos 64.4 35.6 100 309
## Secundaria completa 75.7 24.3 100 408
## Técnica completa 83.1 16.9 100 166
## Univ. incompleta 88.8 11.2 100 134
## Univ. completa 85.8 14.2 100 155
## TOTAL 76.6 23.4 100 1172
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 48.954 4 5.9694e-10
## Pearson 48.273 4 8.2790e-10
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.199
## Cramer's V : 0.203
##
## [[3]]
## gamma : -0.347
## std. error : 0.048
## CI : -0.441 -0.252
t3.5
## [[1]]
## Felices Poco/Nada felices TOTAL
## Es una persona religiosa 76.6 23.4 100
## No es una persona religiosa o es ateo 76.4 23.6 100
## TOTAL 76.6 23.4 100
## Nvalid
## Es una persona religiosa 979
## No es una persona religiosa o es ateo 182
## TOTAL 1161
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 0.0047233 1 0.94521
## Pearson 0.0047308 1 0.94516
##
## Phi-Coefficient : 0.002
## Contingency Coeff.: 0.002
## Cramer's V : 0.002
##
## [[3]]
## gamma : 0.007
## std. error : 0.095
## CI : -0.18 0.193
t3.6
## [[1]]
## Felices Poco/Nada felices TOTAL Nvalid
## 1 vez x semana o más 79.1 20.9 100 478
## 1 vez x mes 73.6 26.4 100 227
## En ocasiones especiales 69.6 30.4 100 23
## 1 vez x año o menos 75.3 24.6 100 434
## TOTAL 76.4 23.6 100 1162
##
## [[2]]
## X^2 df P(> X^2)
## Likelihood Ratio 3.7665 3 0.28781
## Pearson 3.7784 3 0.28641
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.057
## Cramer's V : 0.057
##
## [[3]]
## gamma : 0.079
## std. error : 0.056
## CI : -0.03 0.188
Paso 3:
print(xtable(t3.1[[1]]), type = "html", file="tab3_1.html")
print(xtable(t3.2[[1]]), type = "html", file="tab3_2.html")
print(xtable(t3.3[[1]]), type = "html", file="tab3_3.html")
print(xtable(t3.4[[1]]), type = "html", file="tab3_4.html")
print(xtable(t3.5[[1]]), type = "html", file="tab3_5.html")
print(xtable(t3.6[[1]]), type = "html", file="tab3_6.html")
En este caso el objetivo es evaluar si las medias de una variable dependiente de intervalo o razón varían significativamente entre los niveles de una variable independiente cualitativa que tenga tres o más categorías.
En esta sección estamos generando tablas de estadísticos descriptivos con la función summarySE del paquete Rmisc, así como gráficos descriptivos con el paquete ggplot.
Los resultados serán exportados o copiados al Excel para diseñar las tablas publicables. Para ello utilizamos la función xtable para exportar las tablas de estadísticos descriptivos y de ANOVA en archivos de formato hmlt que luego podemos abrir y editar en Excel.
library(Rmisc)
## Loading required package: lattice
##
## Attaching package: 'lattice'
## The following object is masked from 'package:gnm':
##
## barley
## Loading required package: plyr
##
## Attaching package: 'plyr'
## The following object is masked from 'package:vcdExtra':
##
## summarise
library(ggplot2)
¿El nivel de importancia de Dios en la vida del entrevistado ha variado significativamente entre las oleadas de la Encuesta Mundial de Valores?
Paso 1: Generar los estadísticos descriptivos usando summarySE
desc1 <- summarySE(wvs.peru, measurevar = "imp.dios", groupvars = "oleada",
na.rm = T)
print(xtable(desc1), type = "html", file="descr1.html")
desc1
## oleada N imp.dios sd se ci
## 1 1995-1998 1198 9.017529 1.793107 0.05180572 0.10164012
## 2 1999-2004 1490 9.175839 1.732231 0.04487585 0.08802661
## 3 2005-2009 1484 8.978437 1.865199 0.04841816 0.09497536
## 4 2010-2014 1207 8.913007 1.855891 0.05341939 0.10480527
Paso 2: Generar un gráfico que muestre las medias de cada grupo con sus respectivos intervalos de confianza al 95%
graf3.1 <- ggplot(desc1, aes(x=oleada, y=imp.dios)) + geom_point() +
geom_errorbar(aes(ymin=imp.dios-ci, ymax=imp.dios+ci), width = 0.2) +
ylim(7,10) + xlab("Oleada") + ylab("Importancia de Dios") +
ggtitle("Importancia de Dios según oleada: Media e \n Intervalo de Confianza al 95%") + theme_bw()
graf3.1
jpeg("graf3_1.jpeg")
graf3.1
dev.off()
## quartz_off_screen
## 2
Paso 3: Pedir la tabla de anova de importancia de Dios según oleada
anova1 <- aov(imp.dios~oleada, data=wvs.peru)
summary(anova1)
## Df Sum Sq Mean Sq F value Pr(>F)
## oleada 3 52 17.43 5.315 0.00118 **
## Residuals 5375 17630 3.28
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 43 observations deleted due to missingness
print(xtable(summary(anova1)), type="html", file="anova1.html")
Paso 4: Si se rechaza la hipótesis cero en la prueba de ANOVA se realizan las pruebas de comparaciones múltiples de Tukey:
TukeyHSD(anova1)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = imp.dios ~ oleada, data = wvs.peru)
##
## $oleada
## diff lwr upr p adj
## 1999-2004-1995-1998 0.15830971 -0.02229626 0.33891568 0.1095130
## 2005-2009-1995-1998 -0.03909256 -0.21986118 0.14167606 0.9450086
## 2010-2014-1995-1998 -0.10452176 -0.29432976 0.08528624 0.4898522
## 2005-2009-1999-2004 -0.19740227 -0.36808886 -0.02671567 0.0157259
## 2010-2014-1999-2004 -0.26283147 -0.44306381 -0.08259913 0.0010355
## 2010-2014-2005-2009 -0.06542920 -0.24582453 0.11496613 0.7876176
Paso 1:
desc2 <- summarySE(peru.12, measurevar = "imp.dios", groupvars = "gedad",
na.rm = T)
print(xtable(desc2), type = "html", file="descr2.html")
desc2
## gedad N imp.dios sd se ci
## 1 15-29 410 8.680488 2.122843 0.10483977 0.2060920
## 2 30-49 491 8.997963 1.757200 0.07930131 0.1558126
## 3 50+ 306 9.088235 1.581261 0.09039466 0.1778761
Paso 2:
graf3.2 <- ggplot(desc2, aes(x=gedad, y=imp.dios)) + geom_point() +
geom_errorbar(aes(ymin=imp.dios-ci, ymax=imp.dios+ci), width = 0.2) +
ylim(7,10) + xlab("Grupos de edad") + ylab("Importancia de Dios") +
ggtitle("Perú 2012: Importancia de Dios según grupos de edad: \nMedia e Intervalo de Confianza al 95%") + theme_bw()
graf3.2
jpeg("graf3_2.jpeg")
graf3.2
dev.off()
## quartz_off_screen
## 2
Paso 3:
anova2 <- aov(imp.dios~gedad, data=peru.12)
summary(anova2)
## Df Sum Sq Mean Sq F value Pr(>F)
## gedad 2 35 17.553 5.131 0.00604 **
## Residuals 1204 4119 3.421
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 3 observations deleted due to missingness
print(xtable(summary(anova2)), type="html", file="anova2.html")
Paso 4:
TukeyHSD(anova2)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = imp.dios ~ gedad, data = peru.12)
##
## $gedad
## diff lwr upr p adj
## 30-49-15-29 0.31747554 0.02711321 0.6078379 0.0280794
## 50+-15-29 0.40774749 0.07986788 0.7356271 0.0100251
## 50+-30-49 0.09027195 -0.22583805 0.4063820 0.7809076
Paso 1
desc3 <- summarySE(peru.12, measurevar = "imp.dios", groupvars = "nivedu",
na.rm = T)
print(xtable(desc3), type = "html", file="descr3.html")
desc3
## nivedu N imp.dios sd se ci
## 1 Secundaria incompleta o menos 310 9.170968 1.617000 0.09183942 0.1807097
## 2 Secundaria completa 408 8.928922 1.832641 0.09072923 0.1783564
## 3 Técnica completa 167 9.035928 1.443274 0.11168387 0.2205039
## 4 Univ. incompleta 134 8.559701 2.252855 0.19461704 0.3849450
## 5 Univ. completa 157 8.464968 2.254706 0.17994516 0.3554434
## 6 <NA> 31 9.258065 1.731430 0.31097399 0.6350936
Paso 2:
graf3.3 <- ggplot(na.omit(desc3), aes(x=nivedu, y=imp.dios)) + geom_point() +
geom_errorbar(aes(ymin=imp.dios-ci, ymax=imp.dios+ci), width = 0.2) +
ylim(7,10) + xlab("Nivel Educativo") + ylab("Importancia de Dios") +
ggtitle("Perú 2012: Importancia de Dios según nivel educativo: \nMedia e Intervalo de Confianza al 95%") + coord_flip() + theme_bw()
graf3.3
jpeg("graf3_3.jpeg", width=600, height=400)
graf3.3
dev.off()
Paso 3:
anova3 <- aov(imp.dios~nivedu, data=peru.12)
summary(anova3)
## Df Sum Sq Mean Sq F value Pr(>F)
## nivedu 4 71 17.850 5.24 0.000348 ***
## Residuals 1171 3989 3.406
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 34 observations deleted due to missingness
print(xtable(summary(anova3)), type="html", file="anova3.html")
TukeyHSD(anova3)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = imp.dios ~ nivedu, data = peru.12)
##
## $nivedu
## diff lwr
## Secundaria completa-Secundaria incompleta o menos -0.24204617 -0.6219525
## Técnica completa-Secundaria incompleta o menos -0.13503960 -0.6190396
## Univ. incompleta-Secundaria incompleta o menos -0.61126625 -1.1325611
## Univ. completa-Secundaria incompleta o menos -0.70599959 -1.1999155
## Técnica completa-Secundaria completa 0.10700658 -0.3561959
## Univ. incompleta-Secundaria completa -0.36922008 -0.8712648
## Univ. completa-Secundaria completa -0.46395342 -0.9375075
## Univ. incompleta-Técnica completa -0.47622665 -1.0610138
## Univ. completa-Técnica completa -0.57095999 -1.1314781
## Univ. completa-Univ. incompleta -0.09473334 -0.6877535
## upr p adj
## Secundaria completa-Secundaria incompleta o menos 0.137860147 0.4093549
## Técnica completa-Secundaria incompleta o menos 0.348960369 0.9412880
## Univ. incompleta-Secundaria incompleta o menos -0.089971364 0.0121358
## Univ. completa-Secundaria incompleta o menos -0.212083706 0.0009434
## Técnica completa-Secundaria completa 0.570209006 0.9700344
## Univ. incompleta-Secundaria completa 0.132824625 0.2620950
## Univ. completa-Secundaria completa 0.009600619 0.0580639
## Univ. incompleta-Técnica completa 0.108560522 0.1711845
## Univ. completa-Técnica completa -0.010441856 0.0434432
## Univ. completa-Univ. incompleta 0.498286870 0.9924640
Paso 1:
desc4 <- summarySE(peru.12, measurevar = "satvid", groupvars = "serv.rel",
na.rm = T)
print(xtable(desc4), type = "html", file="descr4.html")
desc4
## serv.rel N satvid sd se ci
## 1 1 vez x semana o más 481 7.203742 2.232354 0.1017866 0.2000023
## 2 1 vez x mes 229 7.222707 2.102256 0.1389210 0.2737332
## 3 En ocasiones especiales 23 6.739130 1.888178 0.3937124 0.8165095
## 4 1 vez x año o menos 432 7.115741 2.165010 0.1041641 0.2047328
## 5 <NA> 41 6.243902 2.299788 0.3591665 0.7259026
Paso 2:
graf3.4 <- ggplot(na.omit(desc4), aes(x=serv.rel, y=satvid)) + geom_point() +
geom_errorbar(aes(ymin=satvid-ci, ymax=satvid+ci), width = 0.2) +
ylim(2,10) + xlab("Asist. Serv. Rel.") + ylab("Satisfacción con la vida") +
ggtitle("Perú 2012: Satisfacción con la vida según asistencia a\n servicios religiosos: Media e Intervalo de Confianza al 95%") + coord_flip() + theme_bw()
graf3.4
jpeg("graf3_4.jpeg", width=600, height=400)
graf3.4
dev.off()
## quartz_off_screen
## 2
Paso 3:
anova4 <- aov(satvid~serv.rel, data=peru.12)
summary(anova4)
## Df Sum Sq Mean Sq F value Pr(>F)
## serv.rel 3 7 2.235 0.472 0.702
## Residuals 1161 5498 4.736
## 45 observations deleted due to missingness
print(xtable(summary(anova4)), type="html", file="anova4.html")