Estimación puntual de estadísticos descriptivos usando R
Inferencia estadística univariada: de la estimación puntual al parámetro
¿Por qué R?
Directamente desde: The R Project for Statistical Computing
Ver también RStudio (RStudio 2022.12.0 ÚLTIMO)
https://www.filehorse.com/es/descargar-rstudio/
R install
Paquetes en R
Paquetes en R
> a = 1
> A = 2
> a == A
[1] FALSE
>help.start() # abre ventana con web para ayuda
>help("pp") # ayuda con la función "pp"
>?pp # ayuda con la función "pp"
>help.search("pp") # busca apariciones de la cadena "pp"
>??pp # busca apariciones de la cadena "pp"
>apropos("pp", mode="function") # lista funciones con "pp" en el nombre
>example(topic) # corre el código en R con ejemplos en un tópico
# determinado (por ejemplo “example(plot)”)
O usar directamente el buscador del sistema de ayuda de *RStudio*
>R.home() # vuelve al directorio raíz de R
[1] "/usr/lib64/R"
>getwd() # vuelve al directorio de trabajo
[1] "/home/user/R"
>setwd("/home/user/newRdir") # establece un nuevo directorio de trabajo
>dir() # muestra el contenido del directorio actual
>ls()... # muestra el contenido del espacio de trabajo
>history(n) # muestra los últimos n comandos
>source("filename.R") # ejecuta comandos en el script filename.R
>sink("register.txt") # envía la salida de R a un fichero externo
>sink() # termina el envío anterior (vuelve a terminal)
>save.image("filename”) # guarda espacio de trabajo
>load("filename”) # carga espacio de trabajo
>save(objectlist,file="filename”) # guarda objetos específicos
>options() # muestra las opciones
>options(opt=...) # fija una opción
R es un lenguaje orientado a objetos (un objeto es algo que puede ser asignado a una variable)
Vectores (unidimensionales)
Matrices (bidimensionales)
Arrays (multidimensionales)
Factores (vectores de variables categóricas, para agrupar los componentes de otro vector)
Listas (colección de objetos, cada uno puede ser de un tipo)
Data Frames (generalización de matrices; cada fila es un elemento, cada columna una variable de diferente tipo)
Funciones (objetos creados para hacer operaciones)
Vectores y Arrays
a<-1.7
a #Muestra el resultado
## [1] 1.7
a=1.7
a #Muestra el resultado
## [1] 1.7
assign("a",1.7)
a #Imprime el valor de "a"
## [1] 1.7
b<-c(10,11,12,15,19)
b #Imprime el valor de "b"
## [1] 10 11 12 15 19
b*b
## [1] 100 121 144 225 361
1/b
## [1] 0.10000000 0.09090909 0.08333333 0.06666667 0.05263158
c <- b-1
c
## [1] 9 10 11 14 18
x[n] # elemento en la posición n
x[c(n1,n2)] # elemento en las posiciones n1 y n2
x[n1:n2] # elemento en las posiciones de n1 a n2
2:10 # Imprime números del 2 al 10
## [1] 2 3 4 5 6 7 8 9 10
5:1 # Imprime números del 5 al 1
## [1] 5 4 3 2 1
seq(from=1, to=10, by=3) # Imprime números del 1 al 10 de 3 en 3
## [1] 1 4 7 10
seq(1, 10, 3) # Imprime números del 1 al 10 de 3 en 3
## [1] 1 4 7 10
seq(length=10, from=1, by=3)# Imprime 10 números empezando del 1 de 3 en 3
## [1] 1 4 7 10 13 16 19 22 25 28
a <- 1:3; b <- rep(a, times=3); c <- rep(a, each=3)
a #Imprime números del 1 al 3
## [1] 1 2 3
b #Imprime números del 1 al 3 dos veces
## [1] 1 2 3 1 2 3 1 2 3
c #Imprime números del 1 al 3 pero cada cifra 3 veces
## [1] 1 1 1 2 2 2 3 3 3
a <- seq(1:10)
a
## [1] 1 2 3 4 5 6 7 8 9 10
b <- (a>5) # asigna valores TRUE o FALSE a partir de una desigualdad
b
## [1] FALSE FALSE FALSE FALSE FALSE TRUE TRUE TRUE TRUE TRUE
a[b] # muestra los valores que cumplen una condición
## [1] 6 7 8 9 10
a[a>5] # igual, pero evitando variables intermedias
## [1] 6 7 8 9 10
a <- "Esto es un ejemplo"
a
## [1] "Esto es un ejemplo"
x <- 1.5
y <- -2.7
paste("los puntos son (",x,",",y,")", sep="") # concatenando x, y más un texto usando 'paste'
## [1] "los puntos son (1.5,-2.7)"
a <- "Esto es un ejemplo" # extrayendo una parte de la cadena
a <- matrix(1:12, nrow=3, ncol=4)
a
## [,1] [,2] [,3] [,4]
## [1,] 1 4 7 10
## [2,] 2 5 8 11
## [3,] 3 6 9 12
b <- matrix(1:8, nrow=4, ncol=4, byrow=TRUE)
b
## [,1] [,2] [,3] [,4]
## [1,] 1 2 3 4
## [2,] 5 6 7 8
## [3,] 1 2 3 4
## [4,] 5 6 7 8
Los elementos se reciclan si hay menos elementos que espacios a a rellenar.
Si no se dice lo contrario, se rellenan por columnas,para hacerlo por filas incluir “byrow=TRUE”
z <- array(1:24, dim=c(2,3,4))
z
## , , 1
##
## [,1] [,2] [,3]
## [1,] 1 3 5
## [2,] 2 4 6
##
## , , 2
##
## [,1] [,2] [,3]
## [1,] 7 9 11
## [2,] 8 10 12
##
## , , 3
##
## [,1] [,2] [,3]
## [1,] 13 15 17
## [2,] 14 16 18
##
## , , 4
##
## [,1] [,2] [,3]
## [1,] 19 21 23
## [2,] 20 22 24
dim(z)
## [1] 2 3 4
Vectores que contienen información categórica útil para agrupar los elementos de otros vectores del mismo tamaño
bv <- c(0.92,0.97,0.87,0.91,0.92,1.04,0.91,0.94,0.96,0.90,0.96,0.86,0.85)
# colores (B-V) para 13 galaxias
morfo <- c("Sab","E","Sab","S0","E","E","S0","S0","E", "Sab","E","Sab","S0") # información morfológica
length(morfo) # comprobamos que es el mismo tamaño
## [1] 13
fmorfo <- factor(morfo) # se crea el factor con 'factor()'
fmorfo
## [1] Sab E Sab S0 E E S0 S0 E Sab E Sab S0
## Levels: E S0 Sab
# muestra los diferentes factores
levels(fmorfo) # muestra los niveles
## [1] "E" "S0" "Sab"
bv[fmorfo=="E"]
## [1] 0.97 0.92 1.04 0.96 0.96
# separa un determinado nivel
mean(bv[fmorfo=="E"]) # hace una operación sobre un nivel
## [1] 0.97
Útil para identificar submuestras y realizar operaciones sólo sobre sus elementos
Colección ordenadas de objetos, donde se pueden agrupar objetos de diferentes tipos (por ejemplo una combinación de vectores, matrices, factores, otras listas, etc.)
gal <- list(name="NGC3379", morf="E", T.RC3=-5, colours=c(0.53,0.96))
gal
## $name
## [1] "NGC3379"
##
## $morf
## [1] "E"
##
## $T.RC3
## [1] -5
##
## $colours
## [1] 0.53 0.96
length(gal) # muestra cuántos elementos tiene
## [1] 4
names(gal) # muestra los nombres de los elementos
## [1] "name" "morf" "T.RC3" "colours"
gal$radio <- TRUE # pueden añadirse nuevos elementos
gal$redshift <- 0.002922 # add a numeric element
names(gal)
## [1] "name" "morf" "T.RC3" "colours" "radio" "redshift"
# se pueden concatenar varias listas para producir otra
Muy versátiles porque pueden almacenar cualquier tipo de información, pero pueden convertirse en estructuras muy complejas.
Tipo especial de lista de gran utilidad estadística. Sus componente son vectores, o factores, de igual longitud.
La información se organiza en una tabla. Típicamente, cada fila corresponde a un elemento de la muestra. Cada columna a una variable medida en toda (o parte de) la muestra.
Los elementos dentro de cada columna son del mismo tipo. Cada columna puede ser de un tipo diferente.
Cuidado: Al crearse, los caracteres se convierten automáticamente en factores. Para evitarlo usar: “options(stringsAsFactors = FALSE)”
Se crean con la función: “data.frame(…)”
Cada columna (variable) tiene un título o nombre
options(stringsAsFactors = FALSE)
df <- data.frame(numbers=c(10,20,30,40),text=c("a","b","c","a"))
df
## numbers text
## 1 10 a
## 2 20 b
## 3 30 c
## 4 40 a
options(stringsAsFactors = TRUE)
## Warning in options(stringsAsFactors = TRUE): 'options(stringsAsFactors = TRUE)'
## is deprecated and will be disabled
df <-data.frame(numbers=c(10,20,30,40),text=c("a","b","c","a"))
df$text
## [1] "a" "b" "c" "a"
x<-data.frame(vector1, vector2, vector3) #Vectores para cada columna
# Los nombres de cada variable (Columna) seran "vector1", etc.
Para referirse a las columnas,x[n],x[n1,n2], x[c(n1,n2)]
O por el nombre:x[c(“vector1”,“vector2”)], x$vector1
Se pueden generar factores a partir de una variable continua usando el comando cut. Si el parámetro break, es un número, indica el número de intervalos. Con table, se construye una tabla de frecuencias:
bv<-c(0.92,0.97,0.87,0.91,0.92,1.04,0.91,0.94,0.96,0.90,0.96,0.86,0.85) #Colores (B-V) para 13 galaxias
fbv<-cut(bv,breaks=3) #divide "bv" en 3 intervalos de igual longitud
fbv # Muestra el intervalo de cada galaxia
## [1] (0.913,0.977] (0.913,0.977] (0.85,0.913] (0.85,0.913] (0.913,0.977]
## [6] (0.977,1.04] (0.85,0.913] (0.913,0.977] (0.913,0.977] (0.85,0.913]
## [11] (0.913,0.977] (0.85,0.913] (0.85,0.913]
## Levels: (0.85,0.913] (0.913,0.977] (0.977,1.04]
table(fbv) #Genera una tabla contando el numero de galaxias
## fbv
## (0.85,0.913] (0.913,0.977] (0.977,1.04]
## 6 6 1
# (frecuencias) en cada intervalo
Si el parámetro break es un vector, indica los extremos de los intervalos:
ffbv<-cut(bv,breaks = c(0.80,0.90,1.00,1.10))
table(ffbv)
## ffbv
## (0.8,0.9] (0.9,1] (1,1.1]
## 4 8 1
Para dividir por cuantiles, usar la opción quantile:
ffffbv<-cut(bv,quantile(bv,(0:4)/4))
table(ffffbv)
## ffffbv
## (0.85,0.9] (0.9,0.92] (0.92,0.96] (0.96,1.04]
## 3 4 3 2
Se puede usar la opción pretty para dividir en n intervalos más o menos redondeados:
fffbv<-cut(bv,pretty(bv,3))
Construcción de tablas multidimensionales
heights<-c(1.64,1.76,1.79,1.65,1.68,1.65,1.86,1.82,1.73,1.75,1.59,1.87,1.73,1.57,1.63,1.71,1.68,1.73,1.53,1.82)
weights<-c(64,77,82,62,71,72,85,68,72,75,81,88,72,71,74,69,81,67,65,73)
ages<-c(12,34,23,53,23,12,53,38,83,28,28,58,38,63,72,44,33,27,32,38)
fheights<-cut(heights,c(1.50,1.60,1.70,1.80,1.90))
fweights<-cut(weights,c(60,70,80,90))
fages<-cut(ages,seq(10,90,10))
Tabla entre dos variables:
ta<-table(fheights,fweights);ta
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 1 1 1
## (1.6,1.7] 2 3 1
## (1.7,1.8] 2 4 1
## (1.8,1.9] 1 1 2
Añadiendo frecuencias marginales:
addmargins(ta)
## fweights
## fheights (60,70] (70,80] (80,90] Sum
## (1.5,1.6] 1 1 1 3
## (1.6,1.7] 2 3 1 6
## (1.7,1.8] 2 4 1 7
## (1.8,1.9] 1 1 2 4
## Sum 6 9 5 20
Frecuencias relativas (con el comando prop.table):
tta<-prop.table(ta)
addmargins(tta)
## fweights
## fheights (60,70] (70,80] (80,90] Sum
## (1.5,1.6] 0.05 0.05 0.05 0.15
## (1.6,1.7] 0.10 0.15 0.05 0.30
## (1.7,1.8] 0.10 0.20 0.05 0.35
## (1.8,1.9] 0.05 0.05 0.10 0.20
## Sum 0.30 0.45 0.25 1.00
De la misma forma se pueden crear tablas tridimensionales, etc
table(fheights,fweights,fages)
## , , fages = (10,20]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 0 0
## (1.6,1.7] 1 1 0
## (1.7,1.8] 0 0 0
## (1.8,1.9] 0 0 0
##
## , , fages = (20,30]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 0 1
## (1.6,1.7] 0 1 0
## (1.7,1.8] 1 1 1
## (1.8,1.9] 0 0 0
##
## , , fages = (30,40]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 1 0 0
## (1.6,1.7] 0 0 1
## (1.7,1.8] 0 2 0
## (1.8,1.9] 1 1 0
##
## , , fages = (40,50]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 0 0
## (1.6,1.7] 0 0 0
## (1.7,1.8] 1 0 0
## (1.8,1.9] 0 0 0
##
## , , fages = (50,60]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 0 0
## (1.6,1.7] 1 0 0
## (1.7,1.8] 0 0 0
## (1.8,1.9] 0 0 2
##
## , , fages = (60,70]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 1 0
## (1.6,1.7] 0 0 0
## (1.7,1.8] 0 0 0
## (1.8,1.9] 0 0 0
##
## , , fages = (70,80]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 0 0
## (1.6,1.7] 0 1 0
## (1.7,1.8] 0 0 0
## (1.8,1.9] 0 0 0
##
## , , fages = (80,90]
##
## fweights
## fheights (60,70] (70,80] (80,90]
## (1.5,1.6] 0 0 0
## (1.6,1.7] 0 0 0
## (1.7,1.8] 0 1 0
## (1.8,1.9] 0 0 0
Se pueden crear tablas bidimensionales a partir de matrices:
mtab<-matrix(c(30,12,47,58,25,32),ncol=2,byrow=TRUE)
colnames(mtab)<-c("ellipticals","spirals")
rownames(mtab)<-c("sample1","sample2","new sample")
mtab
## ellipticals spirals
## sample1 30 12
## sample2 47 58
## new sample 25 32
Pero no son “tablas reales”. Para transformarlas hacer:
rtab<-as.table(mtab)
class(mtab);class(rtab)
## [1] "matrix" "array"
## [1] "table"
El comando summary devuelve información diferente para una matriz y para una tabla:
summary(mtab)
## ellipticals spirals
## Min. :25.0 Min. :12
## 1st Qu.:27.5 1st Qu.:22
## Median :30.0 Median :32
## Mean :34.0 Mean :34
## 3rd Qu.:38.5 3rd Qu.:45
## Max. :47.0 Max. :58
summary(rtab)
## Number of cases in table: 204
## Number of factors: 2
## Test for independence of all factors:
## Chisq = 9.726, df = 2, p-value = 0.007726
Objetos que pueden ser creados por el usuario para hacer, y repetir, operaciones específicas:
Ejemplo: función para calcular la desviación típica de un vector:
stddev <- function(x){
res=sqrt(sum((x-mean(x))^2)/(length(x)-1))
return(res)
}
Se pueden usar y definir funciones dentro de funciones.
El valor devuelto por una función es el resultado de la última expresión evaluada o el especificado con el comando return
Los argumentos de las funciones pueden especificarse por su posición o por su nombre.
Puede haber argumentos con valores por defecto.
mynumbers<-c(1,2,3,4,5)
stddev(mynumbers)
## [1] 1.581139
stddev(x=mynumbers)
## [1] 1.581139
Ej. Función sd de R (calcula la desviación típica):
sd(x=mynumbers)
## [1] 1.581139
sd(x=mynumbers,na.rm=TRUE)
## [1] 1.581139
sd(mynumbers,na.rm=TRUE)
## [1] 1.581139
sd(na.rm=TRUE,x=mynumbers)
## [1] 1.581139
bv.list<-list(
colsSab=c(0.92,0.87,0.90,0.86),
colsE=c(0.97,0.92,1.04,0.96,0.96),
colsSO=c(0.91,0.91,0.94,0.85))
lapply(bv.list,mean)
## $colsSab
## [1] 0.8875
##
## $colsE
## [1] 0.97
##
## $colsSO
## [1] 0.9025
sapply(bv.list,mean)# devuelve un vector
## colsSab colsE colsSO
## 0.8875 0.9700 0.9025
a<-matrix(1:12,nrow=3,ncol=4)
a
## [,1] [,2] [,3] [,4]
## [1,] 1 4 7 10
## [2,] 2 5 8 11
## [3,] 3 6 9 12
apply(a,1,mean) #media por filas (“1”)
## [1] 5.5 6.5 7.5
rowMeans(a)
## [1] 5.5 6.5 7.5
apply(a,1,sum)
## [1] 22 26 30
rowSums(a)
## [1] 22 26 30
apply(a,2,mean) #medias por columnas (“2”)
## [1] 2 5 8 11
apply(a,2,sum)
## [1] 6 15 24 33
bv<-c(0.92,0.97,0.87,0.91,0.92,1.04,0.91,0.94,0.96,0.90,0.96,0.86,0.85)
morfo<-c("Sab","E","Sab","S0","E","E","S0","S0","E","Sab","E","Sab","S0")
fmorfo<-factor(morfo) # crea factor
tapply(bv,fmorfo,mean)
## E S0 Sab
## 0.9700 0.9025 0.8875
library(grid)
colors() # Lista de los 657 colores activos
## [1] "white" "aliceblue" "antiquewhite"
## [4] "antiquewhite1" "antiquewhite2" "antiquewhite3"
## [7] "antiquewhite4" "aquamarine" "aquamarine1"
## [10] "aquamarine2" "aquamarine3" "aquamarine4"
## [13] "azure" "azure1" "azure2"
## [16] "azure3" "azure4" "beige"
## [19] "bisque" "bisque1" "bisque2"
## [22] "bisque3" "bisque4" "black"
## [25] "blanchedalmond" "blue" "blue1"
## [28] "blue2" "blue3" "blue4"
## [31] "blueviolet" "brown" "brown1"
## [34] "brown2" "brown3" "brown4"
## [37] "burlywood" "burlywood1" "burlywood2"
## [40] "burlywood3" "burlywood4" "cadetblue"
## [43] "cadetblue1" "cadetblue2" "cadetblue3"
## [46] "cadetblue4" "chartreuse" "chartreuse1"
## [49] "chartreuse2" "chartreuse3" "chartreuse4"
## [52] "chocolate" "chocolate1" "chocolate2"
## [55] "chocolate3" "chocolate4" "coral"
## [58] "coral1" "coral2" "coral3"
## [61] "coral4" "cornflowerblue" "cornsilk"
## [64] "cornsilk1" "cornsilk2" "cornsilk3"
## [67] "cornsilk4" "cyan" "cyan1"
## [70] "cyan2" "cyan3" "cyan4"
## [73] "darkblue" "darkcyan" "darkgoldenrod"
## [76] "darkgoldenrod1" "darkgoldenrod2" "darkgoldenrod3"
## [79] "darkgoldenrod4" "darkgray" "darkgreen"
## [82] "darkgrey" "darkkhaki" "darkmagenta"
## [85] "darkolivegreen" "darkolivegreen1" "darkolivegreen2"
## [88] "darkolivegreen3" "darkolivegreen4" "darkorange"
## [91] "darkorange1" "darkorange2" "darkorange3"
## [94] "darkorange4" "darkorchid" "darkorchid1"
## [97] "darkorchid2" "darkorchid3" "darkorchid4"
## [100] "darkred" "darksalmon" "darkseagreen"
## [103] "darkseagreen1" "darkseagreen2" "darkseagreen3"
## [106] "darkseagreen4" "darkslateblue" "darkslategray"
## [109] "darkslategray1" "darkslategray2" "darkslategray3"
## [112] "darkslategray4" "darkslategrey" "darkturquoise"
## [115] "darkviolet" "deeppink" "deeppink1"
## [118] "deeppink2" "deeppink3" "deeppink4"
## [121] "deepskyblue" "deepskyblue1" "deepskyblue2"
## [124] "deepskyblue3" "deepskyblue4" "dimgray"
## [127] "dimgrey" "dodgerblue" "dodgerblue1"
## [130] "dodgerblue2" "dodgerblue3" "dodgerblue4"
## [133] "firebrick" "firebrick1" "firebrick2"
## [136] "firebrick3" "firebrick4" "floralwhite"
## [139] "forestgreen" "gainsboro" "ghostwhite"
## [142] "gold" "gold1" "gold2"
## [145] "gold3" "gold4" "goldenrod"
## [148] "goldenrod1" "goldenrod2" "goldenrod3"
## [151] "goldenrod4" "gray" "gray0"
## [154] "gray1" "gray2" "gray3"
## [157] "gray4" "gray5" "gray6"
## [160] "gray7" "gray8" "gray9"
## [163] "gray10" "gray11" "gray12"
## [166] "gray13" "gray14" "gray15"
## [169] "gray16" "gray17" "gray18"
## [172] "gray19" "gray20" "gray21"
## [175] "gray22" "gray23" "gray24"
## [178] "gray25" "gray26" "gray27"
## [181] "gray28" "gray29" "gray30"
## [184] "gray31" "gray32" "gray33"
## [187] "gray34" "gray35" "gray36"
## [190] "gray37" "gray38" "gray39"
## [193] "gray40" "gray41" "gray42"
## [196] "gray43" "gray44" "gray45"
## [199] "gray46" "gray47" "gray48"
## [202] "gray49" "gray50" "gray51"
## [205] "gray52" "gray53" "gray54"
## [208] "gray55" "gray56" "gray57"
## [211] "gray58" "gray59" "gray60"
## [214] "gray61" "gray62" "gray63"
## [217] "gray64" "gray65" "gray66"
## [220] "gray67" "gray68" "gray69"
## [223] "gray70" "gray71" "gray72"
## [226] "gray73" "gray74" "gray75"
## [229] "gray76" "gray77" "gray78"
## [232] "gray79" "gray80" "gray81"
## [235] "gray82" "gray83" "gray84"
## [238] "gray85" "gray86" "gray87"
## [241] "gray88" "gray89" "gray90"
## [244] "gray91" "gray92" "gray93"
## [247] "gray94" "gray95" "gray96"
## [250] "gray97" "gray98" "gray99"
## [253] "gray100" "green" "green1"
## [256] "green2" "green3" "green4"
## [259] "greenyellow" "grey" "grey0"
## [262] "grey1" "grey2" "grey3"
## [265] "grey4" "grey5" "grey6"
## [268] "grey7" "grey8" "grey9"
## [271] "grey10" "grey11" "grey12"
## [274] "grey13" "grey14" "grey15"
## [277] "grey16" "grey17" "grey18"
## [280] "grey19" "grey20" "grey21"
## [283] "grey22" "grey23" "grey24"
## [286] "grey25" "grey26" "grey27"
## [289] "grey28" "grey29" "grey30"
## [292] "grey31" "grey32" "grey33"
## [295] "grey34" "grey35" "grey36"
## [298] "grey37" "grey38" "grey39"
## [301] "grey40" "grey41" "grey42"
## [304] "grey43" "grey44" "grey45"
## [307] "grey46" "grey47" "grey48"
## [310] "grey49" "grey50" "grey51"
## [313] "grey52" "grey53" "grey54"
## [316] "grey55" "grey56" "grey57"
## [319] "grey58" "grey59" "grey60"
## [322] "grey61" "grey62" "grey63"
## [325] "grey64" "grey65" "grey66"
## [328] "grey67" "grey68" "grey69"
## [331] "grey70" "grey71" "grey72"
## [334] "grey73" "grey74" "grey75"
## [337] "grey76" "grey77" "grey78"
## [340] "grey79" "grey80" "grey81"
## [343] "grey82" "grey83" "grey84"
## [346] "grey85" "grey86" "grey87"
## [349] "grey88" "grey89" "grey90"
## [352] "grey91" "grey92" "grey93"
## [355] "grey94" "grey95" "grey96"
## [358] "grey97" "grey98" "grey99"
## [361] "grey100" "honeydew" "honeydew1"
## [364] "honeydew2" "honeydew3" "honeydew4"
## [367] "hotpink" "hotpink1" "hotpink2"
## [370] "hotpink3" "hotpink4" "indianred"
## [373] "indianred1" "indianred2" "indianred3"
## [376] "indianred4" "ivory" "ivory1"
## [379] "ivory2" "ivory3" "ivory4"
## [382] "khaki" "khaki1" "khaki2"
## [385] "khaki3" "khaki4" "lavender"
## [388] "lavenderblush" "lavenderblush1" "lavenderblush2"
## [391] "lavenderblush3" "lavenderblush4" "lawngreen"
## [394] "lemonchiffon" "lemonchiffon1" "lemonchiffon2"
## [397] "lemonchiffon3" "lemonchiffon4" "lightblue"
## [400] "lightblue1" "lightblue2" "lightblue3"
## [403] "lightblue4" "lightcoral" "lightcyan"
## [406] "lightcyan1" "lightcyan2" "lightcyan3"
## [409] "lightcyan4" "lightgoldenrod" "lightgoldenrod1"
## [412] "lightgoldenrod2" "lightgoldenrod3" "lightgoldenrod4"
## [415] "lightgoldenrodyellow" "lightgray" "lightgreen"
## [418] "lightgrey" "lightpink" "lightpink1"
## [421] "lightpink2" "lightpink3" "lightpink4"
## [424] "lightsalmon" "lightsalmon1" "lightsalmon2"
## [427] "lightsalmon3" "lightsalmon4" "lightseagreen"
## [430] "lightskyblue" "lightskyblue1" "lightskyblue2"
## [433] "lightskyblue3" "lightskyblue4" "lightslateblue"
## [436] "lightslategray" "lightslategrey" "lightsteelblue"
## [439] "lightsteelblue1" "lightsteelblue2" "lightsteelblue3"
## [442] "lightsteelblue4" "lightyellow" "lightyellow1"
## [445] "lightyellow2" "lightyellow3" "lightyellow4"
## [448] "limegreen" "linen" "magenta"
## [451] "magenta1" "magenta2" "magenta3"
## [454] "magenta4" "maroon" "maroon1"
## [457] "maroon2" "maroon3" "maroon4"
## [460] "mediumaquamarine" "mediumblue" "mediumorchid"
## [463] "mediumorchid1" "mediumorchid2" "mediumorchid3"
## [466] "mediumorchid4" "mediumpurple" "mediumpurple1"
## [469] "mediumpurple2" "mediumpurple3" "mediumpurple4"
## [472] "mediumseagreen" "mediumslateblue" "mediumspringgreen"
## [475] "mediumturquoise" "mediumvioletred" "midnightblue"
## [478] "mintcream" "mistyrose" "mistyrose1"
## [481] "mistyrose2" "mistyrose3" "mistyrose4"
## [484] "moccasin" "navajowhite" "navajowhite1"
## [487] "navajowhite2" "navajowhite3" "navajowhite4"
## [490] "navy" "navyblue" "oldlace"
## [493] "olivedrab" "olivedrab1" "olivedrab2"
## [496] "olivedrab3" "olivedrab4" "orange"
## [499] "orange1" "orange2" "orange3"
## [502] "orange4" "orangered" "orangered1"
## [505] "orangered2" "orangered3" "orangered4"
## [508] "orchid" "orchid1" "orchid2"
## [511] "orchid3" "orchid4" "palegoldenrod"
## [514] "palegreen" "palegreen1" "palegreen2"
## [517] "palegreen3" "palegreen4" "paleturquoise"
## [520] "paleturquoise1" "paleturquoise2" "paleturquoise3"
## [523] "paleturquoise4" "palevioletred" "palevioletred1"
## [526] "palevioletred2" "palevioletred3" "palevioletred4"
## [529] "papayawhip" "peachpuff" "peachpuff1"
## [532] "peachpuff2" "peachpuff3" "peachpuff4"
## [535] "peru" "pink" "pink1"
## [538] "pink2" "pink3" "pink4"
## [541] "plum" "plum1" "plum2"
## [544] "plum3" "plum4" "powderblue"
## [547] "purple" "purple1" "purple2"
## [550] "purple3" "purple4" "red"
## [553] "red1" "red2" "red3"
## [556] "red4" "rosybrown" "rosybrown1"
## [559] "rosybrown2" "rosybrown3" "rosybrown4"
## [562] "royalblue" "royalblue1" "royalblue2"
## [565] "royalblue3" "royalblue4" "saddlebrown"
## [568] "salmon" "salmon1" "salmon2"
## [571] "salmon3" "salmon4" "sandybrown"
## [574] "seagreen" "seagreen1" "seagreen2"
## [577] "seagreen3" "seagreen4" "seashell"
## [580] "seashell1" "seashell2" "seashell3"
## [583] "seashell4" "sienna" "sienna1"
## [586] "sienna2" "sienna3" "sienna4"
## [589] "skyblue" "skyblue1" "skyblue2"
## [592] "skyblue3" "skyblue4" "slateblue"
## [595] "slateblue1" "slateblue2" "slateblue3"
## [598] "slateblue4" "slategray" "slategray1"
## [601] "slategray2" "slategray3" "slategray4"
## [604] "slategrey" "snow" "snow1"
## [607] "snow2" "snow3" "snow4"
## [610] "springgreen" "springgreen1" "springgreen2"
## [613] "springgreen3" "springgreen4" "steelblue"
## [616] "steelblue1" "steelblue2" "steelblue3"
## [619] "steelblue4" "tan" "tan1"
## [622] "tan2" "tan3" "tan4"
## [625] "thistle" "thistle1" "thistle2"
## [628] "thistle3" "thistle4" "tomato"
## [631] "tomato1" "tomato2" "tomato3"
## [634] "tomato4" "turquoise" "turquoise1"
## [637] "turquoise2" "turquoise3" "turquoise4"
## [640] "violet" "violetred" "violetred1"
## [643] "violetred2" "violetred3" "violetred4"
## [646] "wheat" "wheat1" "wheat2"
## [649] "wheat3" "wheat4" "whitesmoke"
## [652] "yellow" "yellow1" "yellow2"
## [655] "yellow3" "yellow4" "yellowgreen"
demo(colors) # Demostración
##
##
## demo(colors)
## ---- ~~~~~~
##
## > ### ----------- Show (almost) all named colors ---------------------
## >
## > ## 1) with traditional 'graphics' package:
## > showCols1 <- function(bg = "gray", cex = 0.75, srt = 30) {
## + m <- ceiling(sqrt(n <- length(cl <- colors())))
## + length(cl) <- m*m; cm <- matrix(cl, m)
## + ##
## + require("graphics")
## + op <- par(mar=rep(0,4), ann=FALSE, bg = bg); on.exit(par(op))
## + plot(1:m,1:m, type="n", axes=FALSE)
## + text(col(cm), rev(row(cm)), cm, col = cl, cex=cex, srt=srt)
## + }
##
## > showCols1()
##
## > ## 2) with 'grid' package:
## > showCols2 <- function(bg = "grey", cex = 0.75, rot = 30) {
## + m <- ceiling(sqrt(n <- length(cl <- colors())))
## + length(cl) <- m*m; cm <- matrix(cl, m)
## + ##
## + require("grid")
## + grid.newpage(); vp <- viewport(width = .92, height = .92)
## + grid.rect(gp=gpar(fill=bg))
## + grid.text(cm, x = col(cm)/m, y = rev(row(cm))/m, rot = rot,
## + vp=vp, gp=gpar(cex = cex, col = cm))
## + }
##
## > showCols2()
##
## > showCols2(bg = "gray33")
##
## > ###
## >
## > ##' @title Comparing Colors
## > ##' @param col
## > ##' @param nrow
## > ##' @param ncol
## > ##' @param txt.col
## > ##' @return the grid layout, invisibly
## > ##' @author Marius Hofert, originally
## > plotCol <- function(col, nrow=1, ncol=ceiling(length(col) / nrow),
## + txt.col="black") {
## + stopifnot(nrow >= 1, ncol >= 1)
## + if(length(col) > nrow*ncol)
## + warning("some colors will not be shown")
## + require(grid)
## + grid.newpage()
## + gl <- grid.layout(nrow, ncol)
## + pushViewport(viewport(layout=gl))
## + ic <- 1
## + for(i in 1:nrow) {
## + for(j in 1:ncol) {
## + pushViewport(viewport(layout.pos.row=i, layout.pos.col=j))
## + grid.rect(gp= gpar(fill=col[ic]))
## + grid.text(col[ic], gp=gpar(col=txt.col))
## + upViewport()
## + ic <- ic+1
## + }
## + }
## + upViewport()
## + invisible(gl)
## + }
##
## > ## A Chocolate Bar of colors:
## > plotCol(c("#CC8C3C", paste0("chocolate", 2:4),
## + paste0("darkorange", c("",1:2)), paste0("darkgoldenrod", 1:2),
## + "orange", "orange1", "sandybrown", "tan1", "tan2"),
## + nrow=2)
##
## > ##' Find close R colors() to a given color {original by Marius Hofert)
## > ##' using Euclidean norm in (HSV / RGB / ...) color space
## > nearRcolor <- function(rgb, cSpace = c("hsv", "rgb255", "Luv", "Lab"),
## + dist = switch(cSpace, "hsv" = 0.10, "rgb255" = 30,
## + "Luv" = 15, "Lab" = 12))
## + {
## + if(is.character(rgb)) rgb <- col2rgb(rgb)
## + stopifnot(length(rgb <- as.vector(rgb)) == 3)
## + Rcol <- col2rgb(.cc <- colors())
## + uniqC <- !duplicated(t(Rcol)) # gray9 == grey9 (etc)
## + Rcol <- Rcol[, uniqC] ; .cc <- .cc[uniqC]
## + cSpace <- match.arg(cSpace)
## + convRGB2 <- function(Rgb, to)
## + t(convertColor(t(Rgb), from="sRGB", to=to, scale.in=255))
## + ## the transformation, rgb{0..255} --> cSpace :
## + TransF <- switch(cSpace,
## + "rgb255" = identity,
## + "hsv" = rgb2hsv,
## + "Luv" = function(RGB) convRGB2(RGB, "Luv"),
## + "Lab" = function(RGB) convRGB2(RGB, "Lab"))
## + d <- sqrt(colSums((TransF(Rcol) - as.vector(TransF(rgb)))^2))
## + iS <- sort.list(d[near <- d <= dist])# sorted: closest first
## + setNames(.cc[near][iS], format(zapsmall(d[near][iS]), digits=3))
## + }
##
## > nearRcolor(col2rgb("tan2"), "rgb")
## 0.0 21.1 25.8 29.5
## "tan2" "tan1" "sandybrown" "sienna1"
##
## > nearRcolor(col2rgb("tan2"), "hsv")
## 0.0000 0.0410 0.0618 0.0638 0.0667 0.0766
## "tan2" "sienna2" "coral2" "tomato2" "tan1" "coral"
## 0.0778 0.0900 0.0912 0.0918
## "sienna1" "sandybrown" "coral1" "tomato"
##
## > nearRcolor(col2rgb("tan2"), "Luv")
## 0.00 7.42 7.48 12.41 13.69
## "tan2" "tan1" "sandybrown" "orange3" "orange2"
##
## > nearRcolor(col2rgb("tan2"), "Lab")
## 0.00 5.56 8.08 11.31
## "tan2" "tan1" "sandybrown" "peru"
##
## > nearRcolor("#334455")
## 0.0867
## "darkslategray"
##
## > ## Now, consider choosing a color by looking in the
## > ## neighborhood of one you know :
## >
## > plotCol(nearRcolor("deepskyblue", "rgb", dist=50))
##
## > plotCol(nearRcolor("deepskyblue", dist=.1))
##
## > plotCol(nearRcolor("tomato", "rgb", dist= 50), nrow=3)
##
## > plotCol(nearRcolor("tomato", "hsv", dist=.12), nrow=3)
##
## > plotCol(nearRcolor("tomato", "Luv", dist= 25), nrow=3)
##
## > plotCol(nearRcolor("tomato", "Lab", dist= 18), nrow=3)
Especificar un color
col=92 #de 1 a 657
col="orange"
col=rgb(1,2,3)
col=hsv(299,92,552)
plot(x,y,col=)
Paletas de color con n colores contiguos:
n<-456
col=rainbow(n)
col=heat.colors(n)
col=terrain.colors(n)
col=topo.colors(n)
col=cm.colors(n)
col=gray(1:n/n)
demo(graphics)
##
##
## demo(graphics)
## ---- ~~~~~~~~
##
## > # Copyright (C) 1997-2009 The R Core Team
## >
## > require(datasets)
##
## > require(grDevices); require(graphics)
##
## > ## Here is some code which illustrates some of the differences between
## > ## R and S graphics capabilities. Note that colors are generally specified
## > ## by a character string name (taken from the X11 rgb.txt file) and that line
## > ## textures are given similarly. The parameter "bg" sets the background
## > ## parameter for the plot and there is also an "fg" parameter which sets
## > ## the foreground color.
## >
## >
## > x <- stats::rnorm(50)
##
## > opar <- par(bg = "white")
##
## > plot(x, ann = FALSE, type = "n")
##
## > abline(h = 0, col = gray(.90))
##
## > lines(x, col = "green4", lty = "dotted")
##
## > points(x, bg = "limegreen", pch = 21)
##
## > title(main = "Simple Use of Color In a Plot",
## + xlab = "Just a Whisper of a Label",
## + col.main = "blue", col.lab = gray(.8),
## + cex.main = 1.2, cex.lab = 1.0, font.main = 4, font.lab = 3)
##
## > ## A little color wheel. This code just plots equally spaced hues in
## > ## a pie chart. If you have a cheap SVGA monitor (like me) you will
## > ## probably find that numerically equispaced does not mean visually
## > ## equispaced. On my display at home, these colors tend to cluster at
## > ## the RGB primaries. On the other hand on the SGI Indy at work the
## > ## effect is near perfect.
## >
## > par(bg = "gray")
##
## > pie(rep(1,24), col = rainbow(24), radius = 0.9)
##
## > title(main = "A Sample Color Wheel", cex.main = 1.4, font.main = 3)
##
## > title(xlab = "(Use this as a test of monitor linearity)",
## + cex.lab = 0.8, font.lab = 3)
##
## > ## We have already confessed to having these. This is just showing off X11
## > ## color names (and the example (from the postscript manual) is pretty "cute".
## >
## > pie.sales <- c(0.12, 0.3, 0.26, 0.16, 0.04, 0.12)
##
## > names(pie.sales) <- c("Blueberry", "Cherry",
## + "Apple", "Boston Cream", "Other", "Vanilla Cream")
##
## > pie(pie.sales,
## + col = c("purple","violetred1","green3","cornsilk","cyan","white"))
##
## > title(main = "January Pie Sales", cex.main = 1.8, font.main = 1)
##
## > title(xlab = "(Don't try this at home kids)", cex.lab = 0.8, font.lab = 3)
##
## > ## Boxplots: I couldn't resist the capability for filling the "box".
## > ## The use of color seems like a useful addition, it focuses attention
## > ## on the central bulk of the data.
## >
## > par(bg="cornsilk")
##
## > n <- 10
##
## > g <- gl(n, 100, n*100)
##
## > x <- rnorm(n*100) + sqrt(as.numeric(g))
##
## > boxplot(split(x,g), col="lavender", notch=TRUE)
##
## > title(main="Notched Boxplots", xlab="Group", font.main=4, font.lab=1)
##
## > ## An example showing how to fill between curves.
## >
## > par(bg="white")
##
## > n <- 100
##
## > x <- c(0,cumsum(rnorm(n)))
##
## > y <- c(0,cumsum(rnorm(n)))
##
## > xx <- c(0:n, n:0)
##
## > yy <- c(x, rev(y))
##
## > plot(xx, yy, type="n", xlab="Time", ylab="Distance")
##
## > polygon(xx, yy, col="gray")
##
## > title("Distance Between Brownian Motions")
##
## > ## Colored plot margins, axis labels and titles. You do need to be
## > ## careful with these kinds of effects. It's easy to go completely
## > ## over the top and you can end up with your lunch all over the keyboard.
## > ## On the other hand, my market research clients love it.
## >
## > x <- c(0.00, 0.40, 0.86, 0.85, 0.69, 0.48, 0.54, 1.09, 1.11, 1.73, 2.05, 2.02)
##
## > par(bg="lightgray")
##
## > plot(x, type="n", axes=FALSE, ann=FALSE)
##
## > usr <- par("usr")
##
## > rect(usr[1], usr[3], usr[2], usr[4], col="cornsilk", border="black")
##
## > lines(x, col="blue")
##
## > points(x, pch=21, bg="lightcyan", cex=1.25)
##
## > axis(2, col.axis="blue", las=1)
##
## > axis(1, at=1:12, lab=month.abb, col.axis="blue")
##
## > box()
##
## > title(main= "The Level of Interest in R", font.main=4, col.main="red")
##
## > title(xlab= "1996", col.lab="red")
##
## > ## A filled histogram, showing how to change the font used for the
## > ## main title without changing the other annotation.
## >
## > par(bg="cornsilk")
##
## > x <- rnorm(1000)
##
## > hist(x, xlim=range(-4, 4, x), col="lavender", main="")
##
## > title(main="1000 Normal Random Variates", font.main=3)
##
## > ## A scatterplot matrix
## > ## The good old Iris data (yet again)
## >
## > pairs(iris[1:4], main="Edgar Anderson's Iris Data", font.main=4, pch=19)
##
## > pairs(iris[1:4], main="Edgar Anderson's Iris Data", pch=21,
## + bg = c("red", "green3", "blue")[unclass(iris$Species)])
##
## > ## Contour plotting
## > ## This produces a topographic map of one of Auckland's many volcanic "peaks".
## >
## > x <- 10*1:nrow(volcano)
##
## > y <- 10*1:ncol(volcano)
##
## > lev <- pretty(range(volcano), 10)
##
## > par(bg = "lightcyan")
##
## > pin <- par("pin")
##
## > xdelta <- diff(range(x))
##
## > ydelta <- diff(range(y))
##
## > xscale <- pin[1]/xdelta
##
## > yscale <- pin[2]/ydelta
##
## > scale <- min(xscale, yscale)
##
## > xadd <- 0.5*(pin[1]/scale - xdelta)
##
## > yadd <- 0.5*(pin[2]/scale - ydelta)
##
## > plot(numeric(0), numeric(0),
## + xlim = range(x)+c(-1,1)*xadd, ylim = range(y)+c(-1,1)*yadd,
## + type = "n", ann = FALSE)
##
## > usr <- par("usr")
##
## > rect(usr[1], usr[3], usr[2], usr[4], col="green3")
##
## > contour(x, y, volcano, levels = lev, col="yellow", lty="solid", add=TRUE)
##
## > box()
##
## > title("A Topographic Map of Maunga Whau", font= 4)
##
## > title(xlab = "Meters North", ylab = "Meters West", font= 3)
##
## > mtext("10 Meter Contour Spacing", side=3, line=0.35, outer=FALSE,
## + at = mean(par("usr")[1:2]), cex=0.7, font=3)
##
## > ## Conditioning plots
## >
## > par(bg="cornsilk")
##
## > coplot(lat ~ long | depth, data = quakes, pch = 21, bg = "green3")
##
## > par(opar)
demo(image)
##
##
## demo(image)
## ---- ~~~~~
##
## > # Copyright (C) 1997-2009 The R Core Team
## >
## > require(datasets)
##
## > require(grDevices); require(graphics)
##
## > x <- 10*(1:nrow(volcano)); x.at <- seq(100, 800, by=100)
##
## > y <- 10*(1:ncol(volcano)); y.at <- seq(100, 600, by=100)
##
## > # Using Terrain Colors
## >
## > image(x, y, volcano, col=terrain.colors(100),axes=FALSE)
##
## > contour(x, y, volcano, levels=seq(90, 200, by=5), add=TRUE, col="brown")
##
## > axis(1, at=x.at)
##
## > axis(2, at=y.at)
##
## > box()
##
## > title(main="Maunga Whau Volcano", sub = "col=terrain.colors(100)", font.main=4)
##
## > # Using Heat Colors
## >
## > image(x, y, volcano, col=heat.colors(100), axes=FALSE)
##
## > contour(x, y, volcano, levels=seq(90, 200, by=5), add=TRUE, col="brown")
##
## > axis(1, at=x.at)
##
## > axis(2, at=y.at)
##
## > box()
##
## > title(main="Maunga Whau Volcano", sub = "col=heat.colors(100)", font.main=4)
##
## > # Using Gray Scale
## >
## > image(x, y, volcano, col=gray(100:200/200), axes=FALSE)
##
## > contour(x, y, volcano, levels=seq(90, 200, by=5), add=TRUE, col="black")
##
## > axis(1, at=x.at)
##
## > axis(2, at=y.at)
##
## > box()
##
## > title(main="Maunga Whau Volcano \n col=gray(100:200/200)", font.main=4)
##
## > ## Filled Contours are even nicer sometimes :
## > example(filled.contour)
##
## flld.c> require("grDevices") # for colours
##
## flld.c> filled.contour(volcano, asp = 1) # simple
##
## flld.c> x <- 10*1:nrow(volcano)
##
## flld.c> y <- 10*1:ncol(volcano)
##
## flld.c> filled.contour(x, y, volcano,
## flld.c+ color.palette = function(n) hcl.colors(n, "terrain"),
## flld.c+ plot.title = title(main = "The Topography of Maunga Whau",
## flld.c+ xlab = "Meters North", ylab = "Meters West"),
## flld.c+ plot.axes = { axis(1, seq(100, 800, by = 100))
## flld.c+ axis(2, seq(100, 600, by = 100)) },
## flld.c+ key.title = title(main = "Height\n(meters)"),
## flld.c+ key.axes = axis(4, seq(90, 190, by = 10))) # maybe also asp = 1
##
## flld.c> mtext(paste("filled.contour(.) from", R.version.string),
## flld.c+ side = 1, line = 4, adj = 1, cex = .66)
##
## flld.c> # Annotating a filled contour plot
## flld.c> a <- expand.grid(1:20, 1:20)
##
## flld.c> b <- matrix(a[,1] + a[,2], 20)
##
## flld.c> filled.contour(x = 1:20, y = 1:20, z = b,
## flld.c+ plot.axes = { axis(1); axis(2); points(10, 10) })
##
## flld.c> ## Persian Rug Art:
## flld.c> x <- y <- seq(-4*pi, 4*pi, length.out = 27)
##
## flld.c> r <- sqrt(outer(x^2, y^2, `+`))
##
## flld.c> filled.contour(cos(r^2)*exp(-r/(2*pi)), axes = FALSE)
##
## flld.c> ## rather, the key *should* be labeled:
## flld.c> filled.contour(cos(r^2)*exp(-r/(2*pi)), frame.plot = FALSE,
## flld.c+ plot.axes = {})
demo(persp)
##
##
## demo(persp)
## ---- ~~~~~
##
## > ### Demos for persp() plots -- things not in example(persp)
## > ### -------------------------
## >
## > require(datasets)
##
## > require(grDevices); require(graphics)
##
## > ## (1) The Obligatory Mathematical surface.
## > ## Rotated sinc function.
## >
## > x <- seq(-10, 10, length.out = 50)
##
## > y <- x
##
## > rotsinc <- function(x,y)
## + {
## + sinc <- function(x) { y <- sin(x)/x ; y[is.na(y)] <- 1; y }
## + 10 * sinc( sqrt(x^2+y^2) )
## + }
##
## > sinc.exp <- expression(z == Sinc(sqrt(x^2 + y^2)))
##
## > z <- outer(x, y, rotsinc)
##
## > oldpar <- par(bg = "white")
##
## > persp(x, y, z, theta = 30, phi = 30, expand = 0.5, col = "lightblue")
##
## > title(sub=".")## work around persp+plotmath bug
##
## > title(main = sinc.exp)
##
## > persp(x, y, z, theta = 30, phi = 30, expand = 0.5, col = "lightblue",
## + ltheta = 120, shade = 0.75, ticktype = "detailed",
## + xlab = "X", ylab = "Y", zlab = "Z")
##
## > title(sub=".")## work around persp+plotmath bug
##
## > title(main = sinc.exp)
##
## > ## (2) Visualizing a simple DEM model
## >
## > z <- 2 * volcano # Exaggerate the relief
##
## > x <- 10 * (1:nrow(z)) # 10 meter spacing (S to N)
##
## > y <- 10 * (1:ncol(z)) # 10 meter spacing (E to W)
##
## > persp(x, y, z, theta = 120, phi = 15, scale = FALSE, axes = FALSE)
##
## > ## (3) Now something more complex
## > ## We border the surface, to make it more "slice like"
## > ## and color the top and sides of the surface differently.
## >
## > z0 <- min(z) - 20
##
## > z <- rbind(z0, cbind(z0, z, z0), z0)
##
## > x <- c(min(x) - 1e-10, x, max(x) + 1e-10)
##
## > y <- c(min(y) - 1e-10, y, max(y) + 1e-10)
##
## > fill <- matrix("green3", nrow = nrow(z)-1, ncol = ncol(z)-1)
##
## > fill[ , i2 <- c(1,ncol(fill))] <- "gray"
##
## > fill[i1 <- c(1,nrow(fill)) , ] <- "gray"
##
## > par(bg = "lightblue")
##
## > persp(x, y, z, theta = 120, phi = 15, col = fill, scale = FALSE, axes = FALSE)
##
## > title(main = "Maunga Whau\nOne of 50 Volcanoes in the Auckland Region.",
## + font.main = 4)
##
## > par(bg = "slategray")
##
## > persp(x, y, z, theta = 135, phi = 30, col = fill, scale = FALSE,
## + ltheta = -120, lphi = 15, shade = 0.65, axes = FALSE)
##
## > ## Don't draw the grid lines : border = NA
## > persp(x, y, z, theta = 135, phi = 30, col = "green3", scale = FALSE,
## + ltheta = -120, shade = 0.75, border = NA, box = FALSE)
##
## > ## `color gradient in the soil' :
## > fcol <- fill ; fcol[] <- terrain.colors(nrow(fcol))
##
## > persp(x, y, z, theta = 135, phi = 30, col = fcol, scale = FALSE,
## + ltheta = -120, shade = 0.3, border = NA, box = FALSE)
##
## > ## `image like' colors on top :
## > fcol <- fill
##
## > zi <- volcano[ -1,-1] + volcano[ -1,-61] +
## + volcano[-87,-1] + volcano[-87,-61] ## / 4
##
## > fcol[-i1,-i2] <-
## + terrain.colors(20)[cut(zi,
## + stats::quantile(zi, seq(0,1, length.out = 21)),
## + include.lowest = TRUE)]
##
## > persp(x, y, 2*z, theta = 110, phi = 40, col = fcol, scale = FALSE,
## + ltheta = -120, shade = 0.4, border = NA, box = FALSE)
##
## > ## reset par():
## > par(oldpar)
Un aspecto importante en el uso de RStudio enfocado en el análisis de datos sociales es el manejo de base de datos. Esto puede referir tanto a bases que quien efectúa el análisis haya construido como a bases de datos de estudio sociales realizados por otros.
En este manual de apoyo docente se utilizará una base de datos secundaria. Tanto para enseñar los procedimientos de importación, validación y modificación de datos, como para los posteriores aspectos de análisis estadístico descriptivo.
En los siguientes ejemplos se trabajará con la base de datos de la Encuesta Nacional de Opinión Pública del Centro de Estudios Públicos (Encuesta CEP, de aquí en más) correspondiente a su versión 81 cuyo trabajo en terreno se efectuó entre Septiembre y Octubre de 2017. Esta encuesta busca caracterizar las actitudes y opiniones, políticas, sociales y económicas de la población chilena, destacando las necesidades, principales preocupaciones y preferencias de todos los habitantes del territorio nacional. Es una de las fuentes de información más importantes para estudiar la opinión pública en Chile, respecto a temas coyunturales (CEP 2017).
Toda la información de esta encuesta se puede encontrar en el apartado Encuesta CEP de la página institucional del centro de estudios mencionado. Desde este sitio en línea el usuario podrá descargar las bases de datos históricas de esta encuesta, junto con sus manuales de uso metodológico.
Descarga de una base de datos de interés
Una primera acción a realizar es la descarga de la base de datos indicada al computador. Esto puede realizarse accediendo a la página web indicada del CEP. En la botonera superior de esta página se encuentra un botón de acceso a la Encuesta CEP (fecha de consulta: 10 de abril 2019).
Imagen 1: Página del CEP
Dentro de esa página se debe buscar el apartado Encuestas anteriores para encontrar en el repositorio la Encuesta CEP Septiembre-Octubre 2017.
Imagen 2: Repositorio de bases de datos históricas de la encuesta CEP
Allí, oprimiendo el botón Base de datos se podrán descargar los archivos de interés.
Luego de clickear el enlace de descarga se debe guardar descargar archivo comprimido como el que se ve en la siguiente imagen (con la extensión “.rar”). Para poder descomprimir el archivo original basta con tener instalado el programa WinRar (puede descargarse desde esta página http://winrar.es/ ). Haciendo clic derecho sobre el archivo descargado, y seleccionando la opción Extraer aquí se obtendrá el archivo original (de extensión “.sav”“)
Imagen 3: Archivos resultantes luego de descomprimir el archivo descargado
Como se ve en la imagen 3 el archivo resultante es un archivo de extensión .sav, lo que indica que se trata de un archivo para abrir y trabajar en SPSS. De aquí en más se usará este archivo para desarrollar ejemplos de análisis estadístico. Además de la base de datos, se descarga el manual de uso de la encuesta en formato PDF.
¿Cómo abrir bases de datos desde formato SPSS?
¿Por qué puede ser importante manejar herramientas que permitan que desde RStudio interactuemos con diferentes formatos de bases de datos? Imagine lo siguiente: durante la formación universitaria usted ha aprendido a usar RStudio como herramienta de análisis estadístico. Pero sucede que al incorporarse a un centro de estudios sociales y observa que el resto de los analistas trabaja fundamentalmente con otro software de análisis estadístico: SPSS. Si usted quiere lograr dialogar con tales especialistas - que muy probablemente no estarán dispuestos a aprender a usar RStudio - deberá lograr al menos abrir bases de datos en formato SPSS, para así poder hacer los análisis que se le encomiende.
Para ello se utiliza un paquete específico llamado haven que cuenta con funciones para abrir bases de datos desde diferentes formatos. Para eso primero se descarga e instala el paquete en el computador:
install.packages("haven")
La función específica a utilizar es \(read_spss\) cuya ejecución es muy simple. Luego de la función se indica el nombre del archivo a leer entre paréntesis y entre comillas. Como se observa a continuación, esta función sirve para leer la base de datos de la Encuesta CEP descargada en el anterior apartado. Es importante recordar que para usar una fuente de datos en análisis futuros se debe guardar como un “objeto” en el entorno de trabajo. Para eso, la ejecución de la función debe asignarse a un objeto nuevo o preexistente mediante la función de asignación. En este caso se asigna a un nuevo objeto llamado “CEP” (CEP <- lectura base de datos)
library(haven) #Debemos asegurarnos de que el paquete a ejecutar,
## Warning: package 'haven' was built under R version 4.2.3
#está cargado en la sesión de Rtudio.
CEP <- read_spss("Base.sav") # Nombre del archivo entre comillas.
Como se observa en la imagen a continuación, luego de tal operación ya se cuenta con un nuevo objeto (CEP) en el entorno de trabajo, que almacena la base de datos de la Encuesta CEP.
Imagen 4: Base de datos CEP cargada como objeto en entorno de trabajo
Con la base de datos cargada en el entorno de trabajo podemos explorar sus principales características. Como se observa en la imagen 5 existe un botón (flecha azul indicando hacia abajo) que permite desplegar la estructura interna de la base de datos existente en el entorno de trabajo.
Imagen 5: Información de los vectores que componen la base de datos
Así se logra explorar rápidamente las principales características de la estructura interna de la base de datos cargada, es decir, los atributos de las diferentes variables en su interior. Entre estas características - a las que se puede acceder con el comando \(attributes\) - se encuentran las etiquetas, entre otros elementos de formato, que generalmente se observan en la vista de variables en SPSS. Si bien en R no se cuenta con una visualización directa de tales elementos, estos se almacenan como atributos en las variables importadas.
También es posible usar el botón con forma de planilla ubicado a la derecha del objeto en el entorno de trabajo para visualizar la base de datos como una planilla y así poder inspeccionarla visualmente. Esto permitirá explorar rápidamente su estructura y contenidos: observar los nombres de las variables, el tipo de valores que almacenan, etc.
Imagen 6.6: Visualización de base CEP como planilla
¿Cómo abrir bases de datos desde diferentes formatos de Microsoft Excel?
Otro formato relevante para el trabajo con base de datos es el ecosistema integrado por los diferentes tipos de planillas de cálculo vinculadas al software Microsoft Excel. Muchas veces la digitación de encuestas se efectúa en softwares de estas características, por lo que un formato primario para almacenar bases de datos es, sin duda, las planillas de cálculo. Es más, muchas veces diversos fenómenos sociales son registrados por la actividad humana en este tipo de archivos.
Si bien el formato planilla de cálculo generalmente se asocia a Microsoft Excel, debido a la masividad de sus productos, refiere a un formato más general. Hoy en día, también se asocia a aplicaciones en línea como Hojas de Cálculo de Google u otros software de funcionalidad de oficina como Calc en sus versiones de Libre Office y Open Office.
A continuación se indican dos modalidades para importar bases de datos en formato de hoja de cálculo a R.
La primera opción es trabajar directamente con un archivo con formato para Microsoft Excel 2007 o superior. Se trata de archivos con una extensión .xlsx que tienen un formato de libro, es decir, pueden soportar en su interior a más de una hoja de trabajo. Para efectos caso de este manual se ha convertido la base de datos de la Encuesta CEP trabajada en el apartado anterior a tal formato. Este procedimiento es sencillo si se tiene instalado SPSS. Para efectos de este tutorial el archivo se puede descargar desde el siguiente enlace (indicado en la presentación de este documento), y que lleva por nombre CEP_sep-oct_2017.xlsx*.
Al abrir el archivo desde el explorador de archivos se observa que tiene dos hojas. Una primera almacena el registro de las respuestas con un número de identificación por caso y la fecha de respuesta de la encuesta. La segunda hoja es la base de datos propiamente tal y es la que interesa cargar en el entorno de trabajo de R.
Entonces, ¿cómo cargar una base de datos desde este formato al entorno de trabajo en R? En esta instancia se usará el paquete readxl que previamente se debe descargar e instalar en el computador. Específicamente se usará la función read_excel de este paquete en su versión más simple - sin argumentos adicionales - indicando solamente el nombre del archivo a leer. Su resultado se guardará en un nuevo objeto llamado CEP_excel.
install.packages("readxl") #Descarga e instalación del paquete
library(readxl) #Cargar paquete en sesión de trabajo de R
## Warning: package 'readxl' was built under R version 4.2.3
CEP_excel <- read_excel("Base.xls") #Leer libro excel
Al observar la planilla cargada en el entorno de trabajo se verá que no es la que interesa para desarrollar análisis estadísticos. Por defecto la función lee la primera hoja de del libro de trabajo Excel, por lo que en este caso cargó la planilla con los datos de identificación de cada respondiente, siendo que interesa la lectura de la segunda hoja del libro de trabajo que contiene la base de datos propiamente tal.
Imagen 9: Hoja de identificación de respondientes cargada como base de datos
Para solucionar este problema se repetirá la operación pero utilizando un argumento extra en la función. Mediante el argumento sheet = se le indica al programa la posición o nombre de la hoja que interesa leer en el interior del libro de trabajo. En este caso se indicará que interesa leer la hoja ubicada la posición “2” del libro, o la hoja de nombre “DATOS” (que es lo mismo para efectos de este manual). A continuación se sobrescribe el objeto creado en el entorno de trabajo, actualizando la función con esta información.
CEP_excel<- read_excel("Base.xls", sheet= 1)
#indica posición de hoja en el libro de trabajo.
CEP_excel <-read_excel("Base.xls", sheet = "Base")
#indica nombre de hoja en libro de trabajo.
Observando el entorno de trabajo se verá que el objeto CEP_excel se ha actualizado y ahora presenta 220 variables.
Imagen 10: Hoja con respuestas a encuesta cargada como base de datos
Ahora bien, la función read_excel considera a la primera línea de datos como los nombres de las variables de forma automática. Hay veces en que en la primera fila no se encuentra el nombre de las variables habiendo primero otro tipo de información. Por eso resulta importante indicarle a la función desde qué fila comenzar a leer los datos. En el caso de que la información relevante comience en la fila 2, siendo tal fila la que contiene el nombre de las variables se agrega a estos argumentos la opción skip = 1 para que el software salte u omita la primera fila y comience la lectura de los datos en la fila 2.
CEP_excel <- read_excel("Base.xls", sheet = 1, skip = 1)
#indica posición de la hoja en el libro de trabajo (Desde donde empieza).
Una segunda forma de importar a R archivos del tipo hoja de cálculo es trabajar con el formato *CSV (comma separated values) o archivo de valores delimitados por comas. Se trata de un tipo de archivo más simple que un libro de Microsoft Excel, en la medida que puede contener sólo una - y no varias - hoja de trabajo. Si bien este tipo de archivos pueden encontrarse al descargar una base de datos secundaria, para efectos de este manual de apoyo docente se creará un archivo de este tipo a partir del archivo originalmente almacenado en formato Excel (.xlsx).
Para esto, como se observa en las imágenes 11 y 12, basta con que desde el archivo excel en cuestión, se utilice a la opción Guardar como y seleccionando el formato señalado. El programa debería advertir al usuario que el formato seleccionado no soporta múltiples hojas de trabajo (imagen 12), por lo que guardará sólo la hoja activa. Asegurando que la planilla que interesa guardar es la hoja que está seleccionada se hace clic sobre la opción Aceptar.
Imagen 11: Guardar planilla de libro de Excel en formato CSV
Imagen 12: Advertencia al guardar en formato CSV
El archivo resultante puede ser leído como planilla de Microsoft Excel. Sin embargo, para entender su estructura interna primero se abrirá con la aplicación Bloc de Notas mediante la opción Abrir con….
Imagen 13: Abrir archivo CSV con Bloc de Notas
Imagen 14: Estructura interna del archivo CSV
Se observa en la imagen 14 que el archivo presenta una tabulación del tipo matriz de datos - en este caso, la base de datos que contiene las respuestas a una encuesta social - cuyos valores ( cada variable) están separados por el signo punto y coma (;). Este signo delimita cada columna de casos.
Resulta importante considerar esta estructura pues en la notación anglosajona que subyace al lenguaje original de R (el inglés) se usa la coma para separar valores y el punto para denotar valores decimales. En el caso de la notación de habla hispana se emplea la coma para denotar decimales y el punto y coma para separar las observaciones.
Este “detalle” importa pues las funciones básicas para leer archivos CSV viene configuradas por defecto para entender a las comas como separador de los casos. Esto se muestra en el siguiente ejemplo.
CEP_csv <- read.csv("Base.csv")
En este caso, la ejecución del comando implica un error de ejecución.
Imagen 15: Mensaje de error al leer archivo CSV
El programa avisa que la cantidad de columnas que resultan de la lectura de la planilla es mayor a la cantidad de nombres de variables, por lo que no logra leer el archivo. Esto sucede debido a que las comas presentes en los casos, que denotan valores decimales, el software las ha entendido como separador de casos.
Para solucionar tal problema se indican dos opciones.
1. La primera es ocupar la misma función, pero agregando un argumento sep = mediante el cual se indica qué signo debe considerar para separar los valores.
2. La segunda opción es usar una variación de la función original (read.csv2), configurada para que considere las comas como notación de decimales y el punto y coma como separador de valores.
CEP_csv <- read.csv("Base.csv", sep = ";")
CEP_csv2 <- read.csv2("Base.csv")
Imagen 16: Bases de datos Excel y CSV con la misma estructura de casos y variables
El resultado muestra que se han leído tres bases de datos coincidentes en términos de estructura, por lo que se ha logrado llegar al mismo resultado usando funciones diferentes.
Para el desarrollo de todos los ejemplos posteriores de este manual se considerará una de las bases de datos leídas. Para ello se “limpiará” por primera vez el entorno de trabajo dejando solamente la base de datos nombrada como CEP_csv.
remove(CEP_csv2, CEP_excel)
Construir una base de datos sólo con variables de interés: exploración de bases de datos y recodificación de variables
Lograr cargar una base de datos a la sesión de R es un paso inicial que permite disponer de un conjunto de datos “en bruto” que, muy probablemente se deberán configurar para lograr efectuar los análisis de interés. Se sugiere trabajar solamente con aquellas variables a analizar y crear una nueva base de datos sólo con la información de interés, sin editar la fuente original de datos. Para los siguientes ejercicios se seleccionarán siete variables desde la base de datos de la Encuesta CEP ya mencionada. En la siguiente tabla se indica una descripción de la variable, su nombre en la base original y el nuevo nombre a asignar en la nueva base de datos.
Para construir una nueva base de datos solo con las variables especificadas será preciso efectuar varias operaciones, que se detallan en este apartado. Para este tipo de manipulación de bases de datos y variables existen diferentes herramientas: algunas forman parte de las funcionalidades básicas de R, pero otras provienen del paquete dplyr. En la siguiente tabla se resumen las principales características de las funciones a utilizar.
Así, de forma específica para la selección de las variables ya señaladas se utilizará la función select del paquete dplyr, como se observa a continuación:
library(dplyr) #Cargar paquete, si no está cargado desde antes.
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
CEP <- select(CEP_excel , pond = POND, sexo = SEXO, region = REGION, edad = DS_P2_EXACTA,satisfaccion_vida = SV_1, satisfaccion_chilenos = SV_2, eval_econ = MB_P2)
#Se indica base de datos, el nombre de variable a crear y los datos que la compondrán.
CEP #Visualización de la base
## # A tibble: 1,424 × 7
## pond sexo region edad satisfaccion_vida satisfaccion_chilenos eval_econ
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1.34 2 13 18 8 3 1
## 2 1.27 2 1 57 10 5 2
## 3 0.605 2 14 25 10 10 4
## 4 1.03 2 13 37 8 5 3
## 5 0.675 2 14 50 5 5 2
## 6 0.292 2 8 60 9 5 2
## 7 0.694 2 9 66 9 5 4
## 8 1.34 2 13 19 6 8 2
## 9 0.787 2 7 34 6 7 3
## 10 1.03 2 13 39 10 10 3
## # … with 1,414 more rows
Luego de ejecutar ese comando se habrá creado un nuevo objeto llamado CEP (del tipo base de datos) en el entorno de trabajo. Esta base de datos tiene la misma cantidad de casos (1.424) que la base de datos original (CEP_csv) pero presenta solamente las siete variables seleccionadas y renombradas mediante el comando select. Para corroborar el resultado de la construcción de esta nueva base de datos se la visualizará en formato planilla. Como se observa a continuación, esta base de datos contiene solamente las siete variables de interés.
Imagen 17: Base de datos con selección de variables
Imagen 18: Visualización base de datos construida
Una vez hecha esta operación de selección de variables se habrá creado una nueva base de datos que contiene solo aquellas variables que resultan de interés para los análisis. Antes de proseguir conviene guardar tal objeto como una base de datos de formato R. Para esto se usa el comando save para guardar un objeto de formato (o extensión) .RData.
save(CEP, file = "seleccion_CEP.RData")
Este comando creará un nuevo archivo de formato RData en la carpeta de trabajo. Este archivo puede usarse para enviar esta base de datos específica a un equipo de trabajo en el contexto de una investigación colectiva, o sencillamente contar con un archivo que ya contenga la base de datos solamente con las variables de interés.
Imagen 19: Nueva base de datos guardada en carpeta de trabajo
El siguiente paso es explorar las características de tales variables y configurarlas para poder ejecutar los análisis estadísticos que precisemos para nuestro proceso de investigación. Se sugiere dejar en cero el espacio de trabajo utilizando algunas de las formas ya indicadas para hacerlo.
rm(list = ls())
Para continuar con este manual se sugiere cargar a la sesión de R (desde la carpeta de trabajo) la base construida solo con las variables de interés y guardada en formato RData (seleccion_CEP.RData).
load("seleccion_CEP.RData")
Una primera opción para conocer las características de la base de datos es explorar los nombres de sus variables. Mediante el comando names aplicado sobre el objeto CEP se obtiene una lista con los nombres de cada columna de la base de datos. Esto resulta de utilidad para los procedimientos de manejo de datos pues no siempre se puede determinar con certeza el nombre de un objeto leyendo el encabezado de la columna. Al visualizarlos como se observa en el siguiente resultado se puede estar seguro del nombre exacto determinando si hay espacios en blanco antes o después de las letras que impliquen un nombre diferente al leído de manera directa.
Ejercicio 6.10
names(CEP)
## [1] "pond" "sexo" "region"
## [4] "edad" "satisfaccion_vida" "satisfaccion_chilenos"
## [7] "eval_econ"
Otro comando de utilidad para conocer características generales de la base de datos es la función dim. Este comando, como se observa en el ejercicio 11, permite conocer la dimensión de la base de datos que se está explorando. El resultado de esta función arroja dos números: el primer número indica la cantidad de filas de la base de datos, mientras que el segundo número indica la cantidad de columnas. Las bases de datos de estudios sociales están construidas de manera que cada fila representa un caso y cada columna una variable, por lo que la dimensión de una base de datos indicará la cantidad de casos y variables (en ese orden).
dim(CEP)
## [1] 1424 7
Finalmente, para conocer las características generales de las variables en la base de datos, y comenzar a evaluar que tipo de recodificaciones se deben realizar, se utiliza la función summary para obtener estadísticos descriptivos de cada una de las variables.
summary(CEP)
## pond sexo region edad
## Min. :0.1557 Min. :1.000 Min. : 1.00 Min. :18.00
## 1st Qu.:0.6419 1st Qu.:1.000 1st Qu.: 6.00 1st Qu.:36.00
## Median :0.8470 Median :2.000 Median : 9.00 Median :50.00
## Mean :1.0000 Mean :1.612 Mean : 9.16 Mean :49.87
## 3rd Qu.:1.2251 3rd Qu.:2.000 3rd Qu.:13.00 3rd Qu.:64.00
## Max. :8.6120 Max. :2.000 Max. :15.00 Max. :97.00
## satisfaccion_vida satisfaccion_chilenos eval_econ
## Min. : 1.000 Min. : 1.000 Min. :1.000
## 1st Qu.: 6.000 1st Qu.: 4.000 1st Qu.:2.000
## Median : 7.000 Median : 5.000 Median :3.000
## Mean : 7.924 Mean : 8.879 Mean :2.821
## 3rd Qu.: 9.000 3rd Qu.: 6.000 3rd Qu.:3.000
## Max. :99.000 Max. :99.000 Max. :9.000
¿Por qué resultan de interés estos resultados? Considerando los estadísticos descriptivos construidos se observa que algunas variables contienen códigos que refieren a casos perdidos. Es el caso de las variables satisfaccion_vida, satisfaccion_chilenos y eval_econ. Como se observa en la imagen 6.20, estas variables presentan valores 9 o 99 como máximos, siendo que se trata de variables que presentan un rango menor de valores posibles: la variable satisfaccion_vida presenta valores 99 como máximos cuando es una escala de 1 a 10; lo mismo ocurre con la variable satisfaccion_chilenos; finalmente, algo similar ocurre con la variable eval_econ, que presenta valores 9 como máximos, cuando es una escala del 1 al 5.
Tal información constituye una primera alerta sobre los ajustes a hacer sobre los datos y por tanto son un valioso insumo para el proceso de recodificación de variables.
Para avanzar en los análisis se efectuarán algunas recodificaciones. Primero que todo se observan las características de cada variable a recodificar. Para eso se utiliza las funciones table para crear una tabla de frecuencias absoluta y observar los valores de la variable, y class para conocer de qué tipo de objeto se trata (en este caso, la variable sexo).
table(CEP$sexo)
##
## 1 2
## 553 871
class(CEP$sexo)
## [1] "numeric"
El primer resultado muestra una variable que presenta sólo valores 1 y 2, a la vez que el segundo resultado informa que se trata de un vector de tipo integer.
Luego de conocer estas características se cuenta con información suficiente para modificar la variable e incorporarla al análisis. Para recodificar se usará el comando mutate del paquete dplyr. Este comando permite transformar una variable guardando el resultado de tal operación en una nueva variable dentro de la misma base de datos. Resulta útil toda vez que permite editar variables sin perder la información (la variable) original.
Como se observa en las siguientes líneas de comando, el resultado de mutate se asigna a la base de datos. Luego de la función mutate, los argumentos son:
i) el conjunto de datos del que provienen las variables (CEP, nuestra base de datos),
ii) el nombre de la nueva variable a crear (sexo_chr) seguido de un igual y la operación de transformación de la variable original; en este caso, se emplea el comando recode del paquete dplyr cuyos argumentos son: la variable que se quiere transformar (CEP$sexo) y luego las equivalencias de cada categoría de respuesta. En este caso se indica que “1” se asigna a “hombre” y que “2” se asigna a “mujer”. Vale la pena enfatizar que cada categoría se expresa entre comillas y que cada equivalencia se separa de la siguiente mediante una coma.
CEP <- mutate(CEP, sexo_chr = recode(CEP$sexo, "1" = "hombre", "2" = "mujer"))
table(CEP$sexo_chr)
##
## hombre mujer
## 553 871
class(CEP$sexo_chr)
## [1] "character"
Al ejecutar tal recodificación y solicitar una tabla de frecuencias simple así como la información sobre el tipo de variable creada, ahora la tabla contiene las categorías “hombre” y “mujer”, haciendo más fácil su lectura, mientras que la variable recodificada es del tipo “character”, al contener ahora una codificación basada en letras. Si bien es un ejemplo muy simple, permite entender una primera utilidad de la recodificación de variables, que tiene que ver con facilitar nuestro acercamiento a los datos.
Adicionalmente, se aplicará una segunda recodificación sobre esta variable para convertirla en un vector de tipo factor, aplicando etiquetas para las dos categorías de respuesta (“Hombre” y “Mujer”). ¿Para qué resulta de utilidad esto si ya contamos con resultados entendibles?: pues bien, en el ejercicio anterior al recodificar a valores de tipo character, la información numérica almacenada en la variable CEP$sexo se perdió al ser reemplazada por letras. Para determinados análisis, precisaremos contar con ambos niveles de información: los valores numéricos asociados a las respuestas, pero también las etiquetas de tales valores, que indican lo que representan los números almacenados en la base de datos en términos analíticos.
CEP <- mutate(CEP, sexo_factor = factor(CEP$sexo, labels =c("Hombre", "Mujer")))
Como se observa en el ejercicio 14, el comando factor permite incorporar las etiquetas a cada código de respuesta. En este caso, al ingresar los valores Hombre y Mujer, en el argumento labels, lo que estamos haciendo es asociar las etiquetas indicadas a los valores 1 y 2 en que está codificada la variable.
Ahora se recodificará la variable región indicando que sólo habrán dos categorías: “otras regiones” y “región Metropolitana”. Como en el caso anterior, primero conviene observar las características generales de la variable: interesa saber cuantos casos están en la categoría Región Metropolitana para así asegurar que luego de la transformación de la variable, la estructura de casos siga el mismo patrón. Primero exploraremos las características generales de la variable de interés.
table(CEP$region) #Observar características de la variables
##
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
## 24 57 24 52 150 82 94 192 98 69 5 17 501 39 20
class(CEP$region)
## [1] "numeric"
En base a estos primeros resultados se observa que se trata de una variable cuyo nivel de medición es nominal. La categoría “13”, que equivale a la región Metropolitana, presenta 501 casos. También sabemos que la variable está configurada como un vector de tipo integer. Así, usando el comando mutate se recodificará la variable en una nueva denominada region_factor; sin embargo, al interior del comando que se utilizaría por defecto, introduciremos una variante.
Si introducimos una función como argumento dentro una función determinada, R siempre aplica las funciones que corresponden al paquete de la función principal inicialmente indicada (en este caso, es el paquete dplyr). Para el caso del comando recode, el provisto por este paquete no resulta de utilidad para recodificar variables con gran cantidad de categorías, pues no permite realizar una recodificación por tramos de información. Por ello, antes de la función recode - que implicaría utilizar aquella incorporada en dplyr - incorporamos el argumento car::, lo que fuerza a que se ejecute el comando recode incluido en un paquete distinto llamado car, para así poder recodificar indicando tramos de datos.
Para el caso de la siguiente recodificación se indica el mismo valor para los códigos del 1 al 12, y del 13 al 14, mientras que se asigna un valor diferente para el valor 13 (correspondiente a la región Metropolitana). Primero se recodifica en una nueva variable indicando los tramos de recodificación ya explicados. Posteriormente se sobrescribe la variable asignándole el resultado de la operación de convertirla en una variable de tipo factor, aplicando etiquetas para las dos categorías de respuesta ahora existentes.
#Transformar a una variable distinta, categorías "Otras regiones" y "Región Metropolitana".
CEP <- mutate(CEP, region_factor = car::recode(CEP$region, "1:12 = 1; 13 = 2; 14:15 = 1"))
#Sobreescribir variable con resultado de convertir a factor incorporando etiquetas
CEP$region_factor <- factor(CEP$region_factor,
labels = c("Otras regiones", "Región Metropolitana"))
table(CEP$region_factor) #Se sigue manteniendo la estuctura de casos.
##
## Otras regiones Región Metropolitana
## 923 501
class(CEP$region_factor) #Cambia el formato del objeto.
## [1] "factor"
El resultado muestra que, por un lado, se mantuvo la estructura de casos, con 501 casos en la región Metropolitana y 923 casos en otras regiones. Eso indica que la recodificación se hizo de manera adecuada. Además, la tabla muestra etiquetas y al solicitar el formato de la variable indica que es un vector de tipo factor. En síntesis, la recodificación se efectuó de manera óptima.
Ahora se recodificarán las variables satisfacción con la vida propia y percepción de satisfación que los chilenos en general sienten con su propia vida. Primero se exploran ambas variables.
#Explorar variable "satisfacción con la propia vida"
class(CEP$satisfaccion_vida)
## [1] "numeric"
table(CEP$satisfaccion_vida)
##
## 1 2 3 4 5 6 7 8 9 10 88 99
## 20 10 30 63 176 169 256 244 142 304 4 6
summary(CEP$satisfaccion_vida)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.000 6.000 7.000 7.924 9.000 99.000
El análisis exploratorio de esta primera variable permite identificar que la variable está compuesta por valores discretos entre 1 y 10, aunque se observa también valores 88 y 99 en la distribución que - debido a la lectura del cuestionario y la ficha técnica de la encuesta - se sabe que denotan las categorías “no sabe”y “no responde.” Es posible afirmar entonces que el nivel de medición de esta variables de tipo intervalar, pues incluye más de 7 categorías de respuesta.
#Explorar variable "percepción de la satisfacción que los chilenos tienen con su vida"
class(CEP$satisfaccion_chilenos)
## [1] "numeric"
table(CEP$satisfaccion_chilenos)
##
## 1 2 3 4 5 6 7 8 9 10 88 99
## 33 21 97 216 469 241 160 56 24 47 52 8
summary(CEP$satisfaccion_chilenos)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.000 4.000 5.000 8.879 6.000 99.000
Algo similar ocurre con esta segunda variable. Al igual que la anterior, presenta valores 88 y 99 como casos válidos, lo que altera la distribución de casos ceñida a la escala de 1 a 10. Tal como en la variable anterior, también es posible afirmar - dado que presenta más de 7 categorías - que el nivel de medición de esta variable es de intervalo.
Por tanto, se efectuará una operación de transformación de estas variables, asignando el valor lógico NA a aquellos valores que representan casos perdidos (sin respuesta). Este valor lógico servirá pues R lo reconoce como un valor que no puede utilizarse para operaciones matemáticas, a la vez que no representa un valor alfanumérico contable, por lo que no alterará el tipo de vector (en este caso, numérico de tipo integer). Entonces se ejecutan las siguientes funciones, para cada variable:
CEP$satisfaccion_vida[CEP$satisfaccion_vida==88]<- NA #Asignar NA a casos con valor 88
CEP$satisfaccion_vida[CEP$satisfaccion_vida==99]<- NA #Asignar NA a casos con valor 99
CEP$satisfaccion_chilenos[CEP$satisfaccion_chilenos==88]<- NA
CEP$satisfaccion_chilenos[CEP$satisfaccion_chilenos==99]<- NA
Habiendo ejecutado tales instrucciones, R habrá cambiado los valores 88 y 99 de cada variable por el valor lógico NA. Para evaluar si la transformación de datos resultó en el sentido esperado, a continuación se solicita una tabla de frecuencias y estadísticos de resumen para cada variable.
table(CEP$satisfaccion_vida) #Ver resultado de codificación
##
## 1 2 3 4 5 6 7 8 9 10
## 20 10 30 63 176 169 256 244 142 304
summary(CEP$satisfaccion_vida)
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 1.000 6.000 7.000 7.311 9.000 10.000 10
table(CEP$satisfaccion_chilenos)
##
## 1 2 3 4 5 6 7 8 9 10
## 33 21 97 216 469 241 160 56 24 47
summary(CEP$satisfaccion_chilenos)
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 1.000 4.000 5.000 5.334 6.000 10.000 60
Como se observa en los resultados, si luego de la recodificación se solicitan tablas de frecuencias para ambas variables, estas ya no mostrarán los casos perdidos (88 y 99). Y si se solicitan estadísticos descriptivos, se observará que ahora los valores mínomo y máximo de la variable en la base de datos, coinciden con los valores mínimos y máximos del rango de la variable discernible desde el cuestionario y ficha técnica del estudio, indicándose además la cantidad de casos NA que han resultado excluidos de los análisis.
Ahora se explorará la variable que mide la percepción de la persona entrevistada respecto a la situación económica del país. Como se observa en el ejercicio 17, se utiliza el comando class, table y summary para efectuar una aproximación exploratoria a la variable.
#Explorar variable "evaluación de la economía"
class(CEP$eval_econ)
## [1] "numeric"
table(CEP$eval_econ)
##
## 1 2 3 4 5 8 9
## 74 397 730 204 5 8 6
summary(CEP$eval_econ)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.000 2.000 3.000 2.821 3.000 9.000
Como se observa en los resultados, la variable es una escala de acuerdo compuesta por valores del 1 al 5, aunque también se observan valores 88 y 99 en la distribución que gracias a la información disponible en el cuestionario y ficha técnica de este estudio, sabemos que representan a a las categorías “no sabe” y “no contesta” ( casos perdidos ). Así, es posible afirmar que se trata de una variable de nivel de medición ordinal.
Para efectos de simplificar las respuestas para su análisis, esta variable se recodificará en tres categorías de evaluación (negativa, neutra y positiva). Asimismo, en el ejercicio 17 también se indicará una segunda forma de asignar casos perdidos (valores NA) esta vez desde una recodificación vía comando mutate.
la variable está compuesta por valores discretos entre 1 y 10, aunque se observa también valores 88 y 99 en la distribución que - debido a la lectura del cuestionario y la ficha técnica de la encuesta - se sabe que denotan las categorías “no sabe”y “no responde”. Es posible afirmar entonces que el nivel de medición de esta variables de tipo intervalar, pues incluye más de 7 categorías de respuesta.
#Recodificar variable a 3 tramos: positiva, neutra, negativa
#Valores perdidos se asignan en la misma codificación,
#con argumento else ("todos los demás valores")
CEP <- mutate(CEP, eval_econ_factor =
car::recode(CEP$eval_econ,
"1:2 = 1; 3 = 2; 4:5 = 3; else = NA"))
CEP$eval_econ_factor <- factor(CEP$eval_econ_factor,
labels = c("Positiva", "Neutra", "Negativa"))
table(CEP$eval_econ_factor)
##
## Positiva Neutra Negativa
## 471 730 209
summary(CEP$eval_econ_factor)
## Positiva Neutra Negativa NA's
## 471 730 209 14
Habiendo efectuado estas transformaciones sobre las variables originales, ya se cuenta con variables de los diferentes tipos para efectuar análisis de estadística univariada: variables cuyo nivel de medición es nominal, ordinal, intervalar y de razón. En el siguiente capítulo se verá como ejecutar análisis estadísticos univariados, en sus variantes de alcance muestral y poblacional.
Referencias
CEP. 2017. “Manual Del Usuario Encuesta CEP Nº 81. Estudio Nacional de Opinión Pública Nº 51 – Tercera Serie, Septiembre-Octubre 2017.” Centro de Estudios Públicos. https://www.cepchile.cl/estudio-nacional-de-opinion-publica-septiembre-octubre-2017/cep/2017-10-25/105022.html.
https://rpubs.com/aafernandez1976
Habiendo construido y configurado una base de datos ad hoc para nuestros análisis, en el presente capítulo se presenta un modo de calcular frecuencias, medidas de tendencia central, de posición y de dispersión a nivel de la estimación puntual y del parámetro poblacional.
Para efectos de este manual distinguiremos aquellos estadísticos que se calculan desde una muestra (probabilística o no probabilística) de aquellos que se obtienen vía registros administrativos o de censos. Específicamente, en el caso de estadísticos que se estiman a partir de una muestra probabilística avanzaremos desde la estimación puntual al cálculo del parámetro poblacional (con su respectivo intervalo de confianza, nivel de confianza y error de estimación). Esto último nos permitirá estimar los valores que la estimación puntual alcanza en la población (que representa la muestra) y si existen diferencias estadísticamente significativas a nivel del parámetro entre dos grupos (prueba de hipótesis).
El presente capitulo se ordena como sigue: en el primer apartado se indicará cómo realizar la estimación puntual de diferentes estadísticos descriptivos usando R (medidas de tendencia central, distribuciones de frecuencia y medidas de posición, medidas de dispersión y forma de distribuciones) y su exportación a planillas de cálculo para su presentación en reportes de investigación académica o profesional; en el segundo apartado se presentará la forma de realizar inferencia estadística, construyendo intervalos de confianza de parámetros, a partir de estimaciones puntuales, y calculando indicadores insesgados a partir de factores de expansión.
Medidas de tendencia central: cálculo de la media, mediana y moda
El cálculo de la media es sencillo en cuanto a los códigos y aritmética necesarios. En las dos siguientes líneas de comandos se observa el cálculo de una media simple y luego de una media recortada. En ambos casos se utiliza el comando mean, indicando la variable sobre la cual se desea hacer el cálculo, acompañado del argumento na.rm = TRUE para excluir del cálculo a los casos dados como perdidos (valores NA que ya fueron codificados al preparar las variables). Para el caso de la media recortada, se agrega el argumento trim, que permite indicar la proporción de casos que se eliminan en cada extremo de la distribución. Así, se obtiene la estimación puntual de ambos estadísticos.
mean(CEP$satisfaccion_vida, na.rm = TRUE)
## [1] 7.311174
mean(CEP$satisfaccion_vida, na.rm = TRUE, trim = 0.025)
## [1] 7.390625
Como se observa en el ejercicio 1 la media aritmética arroja un valor de 7,31 en una una escala de 1 a 10, mientras la media recortada al 5% presenta un valor de 7,39 (ambos valores se redondean al segundo decimal). Como fue anticipado, la interpretación de estos resultados (y sus diferencias) se explorarán en el siguiente apartado.
Para el caso de la mediana, se trata de un comando similar. Se utiliza el comando median, indicando la variable sobre la cual se desea hacer el cálculo, acompañado del argumento na.rm = TRUE para excluir del cálculo a los casos dados como perdidos (valores NA que ya fueron codificados al preparar las variables). De esa forma se obtiene la estimación puntual de este estadístico.
median(CEP$satisfaccion_vida, na.rm = TRUE)
## [1] 7
Como se ve en el ejercicio 2, el valor de la mediana es 7, lo que indica que la mitad de los casos indica tener una satsfacción con la propia vida menor a 7, considerando una escala de 1 a 10.
Para el cálculo de la moda, debe haberse instalado previamente el paquete modeest y haberlo cargado en la sesión de trabajo; esto permitirá contar con una función que calculará de forma automática el o los valores más frecuentes de la distribución, herramienta que no existe en los paquetes básicos de R. Así, como se ve en el ejercicio 3, se utiliza el comando mfv para calcular la moda de la variable edad. El resultado de este comando indica el valor, o los valores, con más frecuencia dentro de la distribución de casos de la variable señalada.
#Ejecutar previamente install.packages("modeest")
# library(modeest)
mode <- function(x) {
return(as.numeric(names(which.max(table(CEP$edad)))))
}
moda<-mode(x)
#Indica el o los valores con más frecuencia
En este caso, la moda de la variable edad es 50 años, lo que indica que el valor que más se repite en la distribución de casos de tal variable corresponde a tal respuesta.
Medidas de posición: cálculo de frecuencias absolutas y relativas, cuantiles
El cálculo de tablas de frecuencias absolutas para una variable se efectúa mediante el comando table, indicando como argumento de la función la variable sobre la cual se ejecuta el cálculo. En este sub apartado, trabajaremos con la variable ordinal de nuestra base de datos que refiere a una evaluación de la economía (CEP$eval_econ_factor) por parte de la persona encuestada.
tabla <- table(CEP$eval_econ_factor)
tabla
##
## Positiva Neutra Negativa
## 471 730 209
De forma muy sucinta, los resultados del ejercicio 4 muestran un importante predominio de opiniones “neutras” con 730 casos, seguida por una concentración de 471 casos en respuestas “positivas” y 209 casos en respuestas “negativas”.
Para facilitar la comparación con otros datos, estos resultados pueden expresarse como frecuencias relativas. Para el cálculo de frecuencias relativas, se usa el objeto tabla de frecuencias simples ya construido (tabla): como se ve en el ejercicio 5, sobre tal objeto, se ejecuta la función prop.table que construye una nueva tabla en la que cada celda es la división simple entre la cantidad de casos de la categoría y el total de casos.
#Frecuencias relativas (proporciones)
prop.table(tabla)
##
## Positiva Neutra Negativa
## 0.3340426 0.5177305 0.1482270
De tal forma, los resultados del ejercicio 5 confirman lo ya observado con frecuencias absolutas. En este caso, los valores se presentan como proporciones de una unidad (1). Así, se observa que la categoría “neutra” concentra la mayor cantidad de casos con aproximadamente la mitad de los mismos (0,52), mientras que la categoría “positiva” concentra un tercio de las respuestas (0,33).
Ahora bien, para facilitar aún más la lectura y divulgación de estos resultados resultará de utilidad convertir estas frecuencias relativas en porcentajes. Como se observa en el ejercicio 7.6, para construir una tabla de porcentajes se multiplica el resultado de calcular la tabla de proporciones por 100; si a ello además se le agrega la función round, es posible configurar un resultado redondeado a dos decimales.
#Frecuencias relativas (porcentajes)
prop.table(tabla)*100
##
## Positiva Neutra Negativa
## 33.40426 51.77305 14.82270
round((prop.table(tabla)*100),2)
##
## Positiva Neutra Negativa
## 33.40 51.77 14.82
De tal modo, los resultados del ejercicio 6 permiten afirmar que la mayor cantidad de respuestas se concentran en una evaluación “neutra” con el 51,8% del total, seguido por evaluaciones positivas que representan un 33,4% del total de respuestas, y finalmente las respuestas que reflejan evaluaciones negativas, que alcanzan un 14,8% de los casos.
Una última configuración que resulta de utilidad para el examen de distribuciones de frecuencias es la modalidad de cantidades acumuladas. Este tipo de tablas se construye a partir de tablas de frecuencias absolutas y relativas (sea en su modalidad proporcional o porcentual). En el ejercicio 7 se aplica la función cumsum a las diferentes modalidades de tablas de frecuencias construidas a partir del objeto tabla, que almacena una tabla de frecuencias absolutas. Esta función (cumsum) permite construir tablas de frecuencias acumuladas, a partir del cálculo de cada distribución de frecuencias ya construidas en los ejercicios anteriores.
#Frecuencias absolutas acumuladas
cumsum(tabla)
## Positiva Neutra Negativa
## 471 1201 1410
#Frecuencias relativas acumuladas
cumsum(prop.table(tabla))
## Positiva Neutra Negativa
## 0.3340426 0.8517730 1.0000000
#Porcentaje acumulado redondado en dos decimales
round(cumsum(prop.table(tabla)*100),2)
## Positiva Neutra Negativa
## 33.40 85.18 100.00
Así los resultados del ejercicio 7, permiten observar que una amplia mayoría de los casos presenta una evaluación neutra o positiva de la economía nacional con el 85,2% de los casos (en términos relativos), lo que equivale a 1.201 casos en términos absolutos. Sólo un 14,8% de las respuestas se concentran en una evaluación negativa de la economía.
Una última forma de describir la concentración de casos de una distribución es mediante los cuantiles. Tales estadísticos permiten describir la concentración de casos según cualquier posición relativa de los datos en la distribución, que resulte de interés. Para el cálculo de cuantiles la función quantile permite calcular los casos equivalentes a diferentes proporciones de la distribución (ejercicio 8). El primer argumento de la función es la variable a considerar, en este caso CEP$atisfaccion_chilenos, que refiere a la evaluación que la persona cree los chilenos tiene de su vida, en una escala de 1 a 10; luego, mediante el argumento prob se indica en formato vector los cuantiles a calcular (expresados como proporción). Si existen casos perdidos (codificados como NA) se debe indicar el argumento na.rm = TRUE.
quantile(CEP$satisfaccion_chilenos, prob = c(0.25, 0.5, 0.75), na.rm = TRUE)
## 25% 50% 75%
## 4 5 6
De tal forma, en el ejercicio 8 se solicita al software el cálculo del cuartil 1 (percentil 25 o caso que corta el 25% de la distribución), el cuartil 2 (percentil 50, mediana o caso que corta el 50% de la distribución) y el cuartil 3 (percentil 75, o caso que corta el 75% de la distribución). Los resultados obtenidos en el ejercicio 7.6 permiten concluir que solamente el 25% superior de los casos cree que los chilenos tienen una evaluación de su vida mayor o igual a 6, mientras que la mitad de los casos cree que los chilenos tienen una evaluación de su vida igual o menor a 5, en una escala del 1 al 10.
Medidas de dispersión: rango, varianza, desviación estándar y coeficiente de variación
En relación a las medidas de dispersión se partirá por cálculo del rango. Los valores mínimo y máximo de la distribución pueden calcularse de manera simultánea con la función range, indicando como argumentos la variable de interés y adicionando también el argumento na.rm = TRUE en el caso de que hubieran sido codificados como NA los valores perdidos. Como se observa en el ejercicio 9, los valores mínimo y máximo pueden calcularse de forma independiente, con las funciones min y max respectivamente, que siguen la misma lógica que la función range. En este caso se trabajará sobre la variable CEP$edad.
range(CEP$edad, na.rm = TRUE)
## [1] 18 97
min(CEP$edad, na.rm = TRUE)
## [1] 18
max(CEP$edad, na.rm = TRUE)
## [1] 97
max(CEP$edad, na.rm = TRUE) - min(CEP$edad, na.rm = TRUE)
## [1] 79
Luego, para conocer el rango de los valores se debe restar el valor máximo al mínimo. Los resultados del ejercicio 9 indican que, para el caso de la variable edad el rango es de 79 años, siendo su valor mínimo 18 años y su valor máximo 97 años.
Otras dos medidas de uso generalizado para describir la dispersión de una variable son la varianza y la desviación estándar. El cálculo de ambas medidas, sigue la misma lógica en lo que refiere a la estructura de la función utilizada. Respectivamente, las funciones utilizadas son var para varianza y sd, para desviación estándar. Como se observa en el ejercicio 10 también debe incluirse el argumento na.rm = TRUE en el caso de que hubieran sido codificados como NA los valores perdidos. En este ejercicio se trabaja sobre la variable CEP$atisfaccion_chilenos, que refiere a la evaluación que la persona cree los chilenos tiene de su vida, en una escala de 1 a 10.
var(CEP$satisfaccion_chilenos, na.rm = TRUE)
## [1] 3.017771
sd(CEP$satisfaccion_chilenos, na.rm = TRUE)
## [1] 1.737173
Atendiendo a los resultados del ejercicio 10 puede observarse que la dispersión de la variable que representa la evaluación que cada entrevistado/a cree que los chilenos tienen sobre su propia vida es de 3,02 unidades atendiendo a la varianza y de 1,73 unidades atendiendo a la desviación estándar.
Otra medida que debemos aprender a calcular para nuestros análisis es el coeficiente de variación. Como se observa en el siguiente ejercicio (11), el coeficiente de variación puede calcularse de manera sencilla, ejecutando una división simple entre la desviación estándar y la media de la variable de interés. Sin embargo, también existe la función coefficient.variation (del paquete FinCal) que permite calcular este estadístico indicando como argumento la media y desviación estándar de la variable de interés. En este caso tales valores se incluyen indicando la función para calcularlos: de tal modo se evitan errores humanos al momento de especificar los valores, ya sea por tipeo del valor o por su aproximación decimal. En este ejercicio trabajamos con la variable edad.
sd(CEP$edad)/mean(CEP$edad)
## [1] 0.3566407
#Asegurarse de ejecutar previamente el comando "install.packages("FinCal")"
library(FinCal)
## Warning: package 'FinCal' was built under R version 4.2.3
coefficient.variation(sd=sd(CEP$edad), avg = mean(CEP$edad))
## [1] 0.3566407
De este modo, los resultados del ejercicio 11 muestran que la variable edad presenta una dispersión del 35,67% según el valor obtenido en el coeficiente de variación.
Forma de una distribución: simetría, curtosis y normalidad
Un último grupo de estadísticos que resulta de importancia para el análisis estadístico de variables refiere a aquellos que dan cuenta de la forma general que asume la distribución de una variable. En específico, mostraremos en este apartado el cálculo e interpretación de estadísticos que dan cuenta de la simetría, curtosis y normalidad de una variable.
En el siguiente ejercicio (12) se muestra cómo calcular los coeficientes de simetría y curtosis mediante las funciones skew y kurtosi presentes en el paquete psych.
library(psych)
## Warning: package 'psych' was built under R version 4.2.3
##
## Attaching package: 'psych'
## The following objects are masked from 'package:FinCal':
##
## geometric.mean, harmonic.mean
skew(CEP$satisfaccion_vida)
## [1] -0.5335721
kurtosi(CEP$satisfaccion_vida)
## [1] -0.1461258
Como se observa en los resultados del ejercicio 12, los coeficientes de simetría y curtosis para la variable satisfaccion_vida son de -0,53 y -0.15 unidades de forma respectiva. Dado que tales valores están expresados en la unidad de medida de la variable en cuestión, no son útiles para comparar variables de diferentes unidades de medición. Es por eso que en el ejercicio 7.13 se calculan los coeficientes de simetría y curtosis estandarizados: para ello, se divide la simetría calculada, por la raíz cuadrada del valor 6 (este número es una constante en la fórmula) dividido en la cantidad de casos de la variable (en este caso, 1.401).
skew(CEP$satisfaccion_vida)/sqrt(6/1401)
## [1] -8.153359
kurtosi(CEP$satisfaccion_vida)/sqrt(6/1401)
## [1] -2.232906
Un criterio general para determinar si los coeficientes de simetría y curtosis reflejan una variable semejante a una distribución normal es que ambos valores se encuentren entre -2 y 2. Como puede observarse en los resultados, ambos coeficientes escapan a tal rango por lo que se observa una distribución poco similar a una normal.
Finalmente, la forma más certera para evaluar si la distribución de datos de una variable se comporta según los parámetros que asume una distribución normal es aplicando un test estadístico específicamente orientado a tal evaluación. Para ello, se diferencia el test de Shapiro Wilk y el de Kolmogorov Smirnov. El primero se adecua a muestras pequeñas (menores a 50 casos), mientras que el segundo sirve para muestras de entre 50 y 1.000 casos.
El ejercicio 14 muestra la aplicación de tales pruebas sobre la variable edad. En específico, se utilizan las funciones shapiro.test y ks.test del paquete básico de R(stats). Su ejecución es bastante sencilla pues sólo requieren como argumento la variable sobre la cual se aplicará la prueba.
#Prueba de Shapiro Wilk (muestras pequeñas)
shapiro.test(CEP$edad)
##
## Shapiro-Wilk normality test
##
## data: CEP$edad
## W = 0.97588, p-value = 1.037e-14
#Prueba Kolmogorov Smirnov (muestras grandes)
ks.test(CEP$edad, "pnorm", mean(CEP$edad, na.rm=T), sd(CEP$edad,na.rm=T))
## Warning in ks.test.default(CEP$edad, "pnorm", mean(CEP$edad, na.rm = T), : ties
## should not be present for the Kolmogorov-Smirnov test
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: CEP$edad
## D = 0.056969, p-value = 0.0001935
## alternative hypothesis: two-sided
Para interpretar uno u otro estadístico, debemos hacerlo atendiendo a la cantidad de casos de la variable en cuestión. En este caso, la variable edad presenta 1424 casos. Esto provoca que ninguno de los resultados sea confiable pues las pruebas tienden a fallar cuando no se cumple el supuesto del margen de casos tolerado.
No obstante, vale la pena recordar que para determinar normalidad univariada se busca un valor p mayor a 0,05 para no rechazar la hipótesis nula de igualdad entre la distribución normal teórica y la distribución empírica de datos que estamos evaluado. Si los resultados fueran válidos en términos estadísticos, en ambas situaciones deberíamos rechazar la normalidad pues se aprueba la hipótesis alternativa dado que ambas pruebas presentan valores p menores a 0,05.
Tablas para informes e interpretación de resultados
En el presente apartado se presenta una forma para poder exportar resultados construidos en R a un formato compatible con softwares de tipo planilla de datos como Excel de Microsoft Office o Calc de Libre u Open Office. Es importante recalcar que aquí se enseña una forma simple para exportar resultados, que no se apoya en paquetes especializados en tales procedimientos. Una aproximación general a formas más avanzadas de construcción de reportes de resultados con formato, se presenta en el capítulo 9 de este manual.
En términos generales, el método de tratamiento de datos que presentaremos a continuación supone tres operaciones: en primer lugar, los resultados calculados deben existir en forma de matriz de datos en el entorno de trabajo; en segundo término, tales objetos deben imprimirse a un archivo de tipo CSV (para lo cual emplearemos la función write.csv2), el cual podrá ser abierto con alguno de los softwares tipo planilla de cálculo ya mencionados; así, en tercer lugar, tales datos podrán manipularse y editarse, ajustando su formato según lo que se desee, siendo fácilmente exportables hacia procesadores de texto.
Distribuciones de frecuencias
En primera instancia se guardarán como objetos las tablas de frecuencias ya existentes en nuestro entorno de trabajo de R (creadas en apartados anteriores de este capítulo). En el ejercicio 15 se observa la creación de tres objetos que contendrán, de forma respectiva, una tabla de frecuencias absolutas, una tabla de frecuencias relativas (porcentajes) y una tabla de frecuencia relativa acumulada (porcentajes) para la variable evaluación de la economía chilena (eval_econ_factor).
f <- table(CEP$eval_econ_factor)
f_porc <- round((prop.table(tabla)*100),2)
f_porc_acum <- round(cumsum(prop.table(tabla)*100),2)
Luego, con la función write.csv2 tales tablas pueden imprimirse, de forma individual, a archivos de tipo CSV. Como se observa en la continuación del ejercicio 15, al ejecutar esta función se indica el objeto a imprimir como planilla y luego con el argumento file se indica (entre comillas) el nombre del archivo a crear con su respectiva extensión (.csv en este caso).
write.csv2(f, file = "Tabla 1.csv")
write.csv2(f_porc, file= "Tabla 2.csv")
write.csv2(f_porc_acum, file= "Tabla 3.csv")
El resultado de la anterior operación será la creación de tres archivos de tipo CSV en la ubicación definida como carpeta de trabajo para la sesión de R. La aparición de tales archivos en la carpeta de trabajo que hemos definido para nuestra sesión de R se observan en la imagen a continuación.
Imagen 1: Exportación de resultados a archivos CSV
A partir de estos archivos - como se observa en la siguiente imagen - es posible unificar y editar las tablas en un formato adecuado para su inserción en reportes de investigación académica o profesional. No se profundiza aquí en el manejo de datos en softwares tipo planilla de cálculo, pero las simples operaciones copiar y pegar permiten consolidar los diferentes resultados en una sola tabla. Las tres primeras tablas de la imagen 2 son los datos impresos por R en archivos CSV independientes; la cuarta tabla es una tabla configurada manualmente mediante Microsoft Excel, lista para ser copiada y pegada en cualquier procesador de texto.
Imagen 2: Tabla de frecuencias en archivo CSV
Ahora bien, es importante recalcar cómo se debe realizar la interpretación de una tabla de resultados como la presentada (a continuación se incorpora en un formato adecuado para este documento interactivo).
Reiterando lo señalado para los resultados del ejercicio 7, se observa que una amplia mayoría de los casos presenta una evaluación neutra o positiva de la economía nacional con el 85,2% de los casos (en términos relativos), lo que equivale a 1201 casos en términos absolutos. Dentro de ello destaca que un tercio (33,4% o 471 casos) de las personas declaran una opinión neutra. Finalmente, es posible señalar que sólo un 14,8% de las respuestas (209 casos) se concentran en una evaluación negativa de la economía.
Estadísticos descriptivos<
A continuación, se explicarán dos formas para calcular de manera agregada diversos estadísticos descriptivos de interés, para la posterior construcción de una tabla exportable a formato planilla que contenga tales valores.
Una primera forma es utilizando la función summary, que se ejecuta con la variable de interés como principal argumento (en este caso, la variable CEP$satisfaccion_vida). A diferencia del cálculo de estadísticos de forma individual, en este caso no es necesario indicar mediante argumento el tratamiento de los casos codificados como NA pues esta función reconoce de manera automática tales valores lógicos y los excluye del análisis.
Mediante este simple comando (ejercicio 16) se obtiene una tabla de estadísticos de resumen: valores mínimo y máximo, primer cuartil, mediana, media, tercer cuartil y cantidad de valores codificados como perdidos (valores NA).
summary(CEP$satisfaccion_vida)
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 1.000 6.000 7.000 7.311 9.000 10.000 10
Ahora bien, para exportar estos resultados a una planilla de cálculo se deberán ejecutar diversas operaciones que se muestran en el siguiente ejemplo (17). Primero, los resultados del comando summary se guardan como un objeto específico (al cual se le asigna el nombre de descriptivos).
#Guardar summary como objeto
descriptivos <- summary(CEP$satisfaccion_vida)
Luego, si se le aplican las funciones names y as.numeric a tal objeto (descriptivos) se obtiene un vector que contiene los encabezados de los resultados y otro vector que contiene su valores numéricos: en la continuación del ejercicio 17 se observa el resultado de aplicar ambas funciones sobre el objeto mencionado, que muestran en primer lugar los encabezados de la tabla y luego sus valores.
Finalmente, usando la función as.data.frame se configura como matriz de datos el resultado de unir como filas - mediante la función rbind (unir como filas, o row bind en inglés) el encabezado de los estadísticos descriptivos y sus valores numéricos. Tal objeto se nombra como descr_sat_vida y se configura como un objeto de tipo data.frame.
#Ver nombres y valores del objeto
names(descriptivos)
## [1] "Min." "1st Qu." "Median" "Mean" "3rd Qu." "Max." "NA's"
as.numeric(descriptivos)
## [1] 1.000000 6.000000 7.000000 7.311174 9.000000 10.000000 10.000000
#Configurar como matriz de datos
descr_sat_vida <-as.data.frame(rbind(names(descriptivos), as.numeric(descriptivos)))
View(descr_sat_vida)
El resultado del ejercicio (como se observa en la siguiente imagen) muestra cómo ya hemos configurado el resultado en formato matriz de datos.
Comando summary configurado como matriz de datos
Así, a partir de esta matriz se puede crear un nuevo archivo tipo planilla de cálculo con estos resultados. Que quedará guardado en la carpeta de trabajo con el nombre de “Tabla 4.csv”.
#Exportar matriz a archivo CSV
write.csv2(descr_sat_vida, file = "Tabla 4.csv")
Como se observa en la siguiente imagen, la planilla superior es la exportada desde R, mientras la inferior es la que resulta luego de una edición simple enfocada en editar los encabezados de la tabla y editar el formato de la misma destacando todos los bordes de la tabla y centrando sus valores (entre otras ediciones).
Resultados de summary exportados a archivo CSV
Ahora bien, la función summary no calcula todos los estadísticos que pueden ser de interés. Para construir una tabla completamente personalizada es necesario calcular cada valor que sea de interés y disponerlo en una matriz de datos para así poder guardar una tabla en formato planilla de cálculo.
Como se observa en el \(ejercicio 18\), cada valor construido se concatena como un sólo vector \((descriptivos_satvida)\) de valores numéricos, que será una fila de la planilla a construir (la que contendrá los valores de cada estadístico). Adicionalmente, se construye un vector llamado nombres con los encabezados de cada estadístico a escribir en la tabla, que se utilizará como fila de títulos de la planilla de cálculo a exportar.
#Cálculo simple de estadíticos descriptivos
min <- min(CEP$satisfaccion_vida, na.rm = TRUE)
q1 <- quantile(CEP$satisfaccion_vida, probs = 0.25, na.rm = TRUE)
media <- mean.default(CEP$satisfaccion_vida, na.rm = TRUE)
media_rec <- mean.default(CEP$satisfaccion_vida, trim = 0.025, na.rm = TRUE)
mediana <- median.default(CEP$satisfaccion_vida, na.rm = TRUE)
mode <- function(x) {
return(as.numeric(names(which.max(table(CEP$satisfaccion_vida)))))
}
moda<-mode(x)
var <- var(CEP$satisfaccion_vida, na.rm = TRUE)
desvest <- sd(CEP$satisfaccion_vida, na.rm = TRUE)
q3 <- quantile(CEP$satisfaccion_vida, probs = 0.75, na.rm = TRUE)
max <- max(CEP$satisfaccion_vida, na.rm = TRUE)
s <- skew(CEP$satisfaccion_vida)
c <- kurtosi(CEP$satisfaccion_vida)
#Valores de estadísticos como vector
descriptivos_satvida <- as.numeric(c(min, q1,media, media_rec, mediana, moda,var, desvest, q3, max, s, c))
#Encabezados de cada estadístico como un vector
nombres <- c("Mínimo", "Q1", "Media", "Media recortada", "Mediana", "Moda","Varianza", "Desviación Estándar", "Q3", "Máximo", "Simetría", "Curtosis")
Con esta información ya es posible configurar una planilla de estadísticos descriptivos completamente personalizada. Para ello, como se muestra en la continuación del ejercicio 18, se configura una matriz de datos (función as.data.frame) a partir de los vectores de encabezados de los estadísticos (nombres) y los valores de tales coeficientes (descriptivos_satvida), que se almacena como un objeto de nombre \(descr2\). Finalmente, esta matriz de datos se exporta a una planilla de cálculo mediante la función \(write.csv2\).
descr2 <- as.data.frame(rbind(nombres,descriptivos_satvida))
write.csv2(descr2, file = "Tabla 5.csv")
Todo esto resulta en una planilla como la que se observa en la siguiente imagen, editable para la construcción de reportes formales.
Tabla de estadísticos muestrales personalizada y exportada a CSV
En base a una edición muy simple en software tipo planilla de cálculo, se puede llegar a darle un formato como el que se observa en la siguiente imagen, útil para su presentación en reportes formales de resultados.
Tabla de estadísticos con edición de formato, lista para ser incorporada a procesador de texto
Ahora bien, es importante recalcar cómo se debe realizar la interpretación de una tabla de resultados como la presentada (a continuación se incorpora en un formato adecuado para este documento interactivo).
En la tabla de estadísticos descriptivos presentada destacan diferentes elementos. En primera instancia, los valores mínimo (1) y máximo (10) indican que el rango de las respuestas observadas abarca todas las respuestas posibles de la pregunta del cuestionario.
En segunda instancia, el valor de la media (7,31) es bastante similar al de la mediana, lo que indica que la mitad de las personas de la distribución le asigna la nota 7 a la satisfacción con su propia vida, en la escala de 1 a 10 ya señalada.
En tercera instancia, se puede señalar - respecto a la dispersión de los datos - que la desviación estándar es de 2,11 unidades.
Finalmente, sobre la forma de la distribución se destacan dos elementos: dado que el coeficiente de simetría es negativo (-0,53), se está ante una variable de tipo asimétrica hacia la izquierda de la distribución de casos, en la cual los casos tienden a concentrase hacia la derecha de la media, es decir, en las puntuaciones más altas de la distribución de respuestas; dado que el coeficiente de curtosis es negativo (-0,15), estamos ante una distribución platicúrtica, en la que existe una menor concentración de casos en torno a la media y presenta una forma más “achatada”.
En el apartado anterior se revisaron procedimientos para estimar de forma puntual estadísticos univariados provenientes de una muestra probabilística. En el presente apartado se explicará como calcular intervalos de confianza para parámetros a partir de muestras probabilísticas. Concretamente se indicarán los procedimientos para calcular intervalos de confianza para medias y proporciones. En primera instancia se presenta los procedimientos más simples que presenta el lenguaje R, adecuados para diseños de muestras del tipo aleatorio simple; en segundo término, se expondrá como calcular estimadores insesgados a partir del uso de factores de expansión, y cómo estimar intervalos de confianza apropiados ante la utilización de diseños muestrales complejos.
Cálculo de intervalos de confianza para proporciones
Para calcular el intervalo de confianza de una proporción es de utilidad la función \(exactci\) del paquete \(PropCIs\) (recordar instalarlo vía \(install.packages("PropCIs")\) de manera previa).
La información que necesita esta función es la cantidad de casos que coinciden con la condición de interés (la categoría de la variable cuya frecuencia relativa interesa) y la cantidad que representa a la totalidad de casos de la muestra (n muestral). En este caso, a partir de la encuesta CEP se evaluará la hipótesis de si la mayoría de las y los chilenos declara tener una opinión negativa de la situación económica del país.
Como se observa a continuación en el ejercicio 19, solicitando una tabla de frecuencias es posible identificar la cantidad de casos probables que tienen una percepción negativa de la situación económica del país (730). Por otra parte, solicitando la cantidad de filas de la base de datos (usando la función nrow) se conoce el n total de casos (1.424). Con ambos valores es posible aplicar la función exactci; un tercer argumento a indicar es el nivel de confianza (conf.level), que se expresa en términos numéricos y proporcionales (para este caso, un 95% de confianza se expresa como 0.95).
library(PropCIs)
table(CEP$eval_econ)
##
## 1 2 3 4 5 8 9
## 74 397 730 204 5 8 6
nrow(CEP)
## [1] 1424
exactci(x = 730, n = 1424, conf.level = 0.95)
##
##
##
## data:
##
## 95 percent confidence interval:
## 0.4863248 0.5389039
Así, el resultado de esta función es el límite superior e inferior del intervalo de confianza, construido a partir de la probabilidad ya indicada (proporción de personas que declaran tener una percepción positiva de la situación económica del país).
Con este resultado es posible concluir que, con un 95% de confianza, no es posible afirmar que la proporción de personas que declaran tener una percepción positiva de la situación económica del país sea más que la mitad de población nacional, debido a que el parámetro poblacional se sitúa entre el 48,6% y 53,9% de los casos.
Una segunda variante para el cálculo de parámetros poblacionales es la estimación del intervalo de confianza de una media. Para calcular el intervalo de confianza de una media es útil la función ci.mean del paquete Publish (recordar instalarlo vía install.packages(“Publish”) de manera previa). Su uso se detalla en el ejercicio 20.
La información que necesita esta función solamente es la variable a utilizar. En el caso que se observa a continuación, se calcula el intervalo de confianza para la media de la variable nivel de satisfacción con la propia vida.
El resultado es bastante sencillo e indica el valor del estadístico muestral (el valor de la media bajo la etiqueta mean) a la vez que indica el límite superior e inferior del intervalo de confianza para la media poblacional, o parámetro, bajo la etiqueta CI-95%; esto último indica el nivel de confianza utilizado para la construcción del intervalo.
library(Publish)
## Warning: package 'Publish' was built under R version 4.2.3
## Loading required package: prodlim
## Warning: package 'prodlim' was built under R version 4.2.3
ci.mean(CEP$satisfaccion_vida) #Nivel de confianza por defecto.
## mean CI-95%
## 7.31 [7.20;7.42]
Como se observa en los resultados del ejercicio 20, al ejecutar el comando con su configuración por defecto, este utiliza un 95% de confianza. Con este valor se puede afirmar con un 95% de confianza que la media poblacional, es decir el promedio de satisfacción que los chilenos declaran tener respecto a su propia vida, se encuentra entre los valores 7,2 y 7,42, en una escala de 1 a 10.
Ahora bien, si al comando básico se le agrega el argumento alpha es posible definir el valor complementario al nivel de confianza, es decir la proporción de error aceptable para la estimación. En la continuación del ejercicio 7.20 se define el valor 0,2 o 20%.
ci.mean(CEP$satisfaccion_vida, alpha = 0.2) #Definición manual del nivel de confianza.
## mean CI-80%
## 7.31 [7.24;7.38]
Lo realizado en la continuación del ejercicio 20 implica que, observando los resultados, se puede afirmar con un nivel de confianza del 80% que el promedio de satisfacción que los chilenos declaran tener respecto a su propia vida, se encuentra entre los valores \(7,24\) y \(7,38\), en una escala de 1 a 10.
Un uso muy frecuente para este tipo de cálculos es evaluar si la diferencia entre una misma media para dos grupos es diferente, de manera estadísticamente significativa. Supondremos que se busca calcular el mismo indicador anterior pero efectuando el cálculo de intervalo de confianza para medias, diferenciando según hombres y mujeres.
Como se observa en el código a continuación (ejercicio 21), para indicar a la función cual es la variable de clasificación para construir los grupos, se debe indicar la variable de interés seguida a continuación de una virgulilla (mismo signo usado en la ñ, \(~\)) que la separa de la variable de clasificación. Luego se debe indicar el conjunto de datos de donde provienen tales variables con el argumento data. Como antes, si no se le indica el valor del \[alpha\], la función asume por defecto un cálculo con un 95% de confianza.
ci.mean(satisfaccion_vida~sexo_factor, data=CEP)
## sexo_factor mean CI-95%
## Hombre 7.46 [7.29;7.62]
## Mujer 7.22 [7.07;7.37]
Como se observa en los resultados, las medias muestrales permitirán afirmar que hombres y mujeres presentan una evaluación levemente diferente en relación a la satisfacción que declaran tener con su propia vida. Específicamente, la media de esta variables es de 7,46 para los hombres y de 7,22 para las mujeres. Sin embargo, al observar el comportamiento de tales mediciones a nivel poblacional, es posible afirmar que las variables no difieren de manera estadísticamente significativa. Para un nivel de confianza del 95%, la media de esta variable para los hombres se sitúa entre los valores 7,29 y 7,62 para los hombres y entre los valores de 7,07 y 7,37 para las mujeres. Los valores indican que los intervalos se intersectan, es decir, existe una elevada probabilidad de que estos parámetros sean muy similares e incluso iguales. Por lo tanto, no es posible afirmar que tales valores sean diferentes, de manera estadísticamente significativa, pues existe una alta posibilidad de que sean iguales.
Coeficientes de expansión para estimación de parámetros poblacionales
Las investigaciones sociales basadas en un diseño muestral probabilístico (aleatorio simple, estratificado, por conglomerado o diseños complejos) deben utilizar un ponderador en la estimación de variables de interés para alcanzar validez sobre la población objetivo. Lo anterior se relaciona con las probabilidades de selección de las distintas unidades de muestreo y da cuenta del número de personas de la población que representa cada individuo de la muestra. Este ponderador es conocido como factor de expansión.
Los factores de expansión buscan incorporar en las estimaciones puntuales surgidas desde la muestra, las características y el peso que tienen esos atributos en la población. De esa forma, se producen estimadores insesgados, es decir, que están centrados en el parámetro población. Así, denominaremos “sesgo de un estimador” a la diferencia entre la esperanza (o valor esperado) del estimador y el valor observado del parámetro a estimar. Decimos que un estimador es “insesgado” o “centrado” cuando su sesgo es nulo: o sea, cuando el valor esperado es igual al valor observado.
INE. 2019a. “Base de Datos Encuesta Nacional de Empleo. Trimestre Enero-Febrero-Marzo de 2019.” Archivo en formato CSV. Santiago de Chile: Instituto Nacional de Estadísticas. https://www.ine.cl/estadisticas/laborales/ene.
INE. 2019b. “Población Total de 15 Años Y Más Por Situación En La Fuerza de Trabajo, Nivel Nacional Y Regional, Ambos Sexos.” Santiago de Chile: Instituto Nacional de Estadísticas. https://www.ine.cl/estadisticas/laborales/ene.
Lumley, Thomas. 2019. “Survey Analysis in R.” Survey Package. http://r-survey.r-forge.r-project.org/survey/index.html.
En este capítulo se presentan funciones para la creación de gráficos con una sola variable cuantitativa.
Esta función permite crear el gráfico llamado de tallo y hoja. Este gráfico fue propuesto por Tukey (1977) y a pesar de no ser un gráfico para presentación definitiva, se utiliza a la vez que el analista recoge la información para ver rápidamente la distribución de los datos.
¿Qué muestra este gráfico?*
El centro de la distribución.
La forma general de la distribución:
Ventajas del gráfico:
Muy fácil de realizar y puede hacerse a mano.
Fácil de entender.
Desventajas del gráfico:
El gráfico es tosco y no sirve para presentaciones definitivas.
Funciona cuando el número de observaciones no es muy grande.
No permite comparar claramente diferentes poblaciones
Como ilustración vamos a crear el gráfico de tallo y hoja para la variable altura (cm) de un grupo de estudiantes de la universidad. Primero se leerán los datos disponibles en github y luego se usará la función \(stem\) para obtener el gráfico. A continuación el código usado.
url <- 'https://tinyurl.com/k55nnlu'
datos <- read.table(file=url, header=T)
stem(datos$altura)
##
## The decimal point is 1 digit(s) to the right of the |
##
## 14 | 7
## 15 | 3
## 15 | 679
## 16 | 0001
## 16 | 68888
## 17 | 001334
## 17 | 5678899
## 18 | 000033
## 18 | 88
## 19 | 1
De este gráfico sencillo se puede ver que la variable presenta una mayor frecuencia para alturas ente 170 y 179 cm y que no tiene una distribución simétrica.
La función boxplot sirve para crear un diagrama de cajas y bigote para una variable cuantitativa. La estructura de la función boxplot con los argumentos más comunes de uso se muestran a continuación.
Los argumentos de la función boxplot son:
x: vector numérico con los datos para crear el boxplot.
formula: fórmula con la estructura x ~ g para indicar que las observaciones en el vector x van a ser agrupadas de acuerdo a los niveles del factor g.
data: marco de datos con las variables.
subset: un vector opcional para especificar un subconjunto de observaciones a ser usadas en el proceso de ajuste.
na.action: una función la cual indica lo que debería pasar cuando
los datos contienen
NA’s’’.
range: valor numérico que indica la extensión de los bigotes. Si es positivo, los bigotes se extenderán hasta el punto más extremo de tal manera que el bigote no supere veces el rango intercuatílico (\(IQ\)). Un valor de cero hace que los bigotes se extiendan hasta los datos extremos.
width: un vector con los anchos relativos de las cajas.
varwidth: Si es TRUE, las cajas son dibujadas con anchos proporcionales a las raíces cuadradas del número de observaciones en los grupos.
notch: si es TRUE, una cuña es dibujada a cada lado de las cajas. Cuando las cuñas de dos gráficos de caja no se traslapan, entonces las medianas son significativamente diferentes a un nivel del 5%.
names: un con las etiquetas a ser impresas debajo de cada boxplot.
plot: si es TRUE (por defecto) entonces se produce el gráfico, de lo contrario, se producen los resúmenes de los boxplots.
col: vector con los colores a usar en el cuerpo de las cajas.
log: para indicar si las coordenadas x o y o serán graficadas en escala logarítmica.
…: otros parámetros gráficos que pueden ser pasados como argumentos para el boxplot.
Como ilustración vamos a crear dos boxplot para la variable altura (cm) de un grupo de estudiantes de la universidad, el primer boxplot será vertical y el segundo horizontal. Primero se leerán los datos disponibles en github y luego se usará la función boxplot para obtener ambos gráfico. A continuación el código usado.
url <- 'https://tinyurl.com/k55nnlu'
datos <- read.table(file=url, header=T)
par(mfrow=c(1, 2))
boxplot(x=datos$altura, ylab='Altura (cm)')
boxplot(x=datos$altura, xlab='Altura (cm)', horizontal=TRUE)
En la Figura se presentan los boxplots obtenidos con las instrucciones anteriores. El segundo y tercer boxplot son el mismo, lo único que se modificó fueron los nombres o etiquetas a colocar debajo de cada boxplot por medio del argumento names y la orientación.
Es posible crear boxplots comparativos usando 1 o 2 variables cualitativas. A continuación se construyen dos boxplots para la variable precio de apartamentos usados en Medellín. En el primer boxtplot diferencia por la variable balcón (no, si) y en el segundo se diferencia por los cruces de las variables parqueadero y ubicación (en Laureles y Poblado).
A continuación se muestra el código necesario para construir los boxplots. En el primero se usa la fórmula precio ~ balcon para crear el boxplot del precio diferenciando por los dos niveles de la variable balcón. En el segundo se usa la fórmula precio ~ ubicacion parqueadero* pero se limitan las ubicaciones a sólo dos, Laureles y Aburrá sur, por esa razón se usa el parámetro subset para incluir la restricción. Se agregó también drop=TRUE para que en el segundo boxplot no aparezcan las otras ubicaciones.
url <- 'https://tinyurl.com/hwhb769'
datos <- read.table(url, header=T)
par(mfrow=c(1, 2))
boxplot(precio ~ balcon, data=datos,
col=c('lightblue', 'pink'),
xlab='¿Tiene balcón?', main='(a)',
ylab='Precio (millones $)')
boxplot(precio ~ ubicacion * parqueadero,
data=datos, drop=TRUE, col=c('lightblue', 'pink'),
subset=ubicacion %in% c('laureles', 'aburra sur'),
xlab='Ubicación y parqueadero',
ylab='Precio (millones $)', main='(b)')
En la Figura se muestran los boxplots. Se le recomienda al lector construir nuevamente el segundo boxplot pero eliminando drop=TRUE para que vea el efecto que tiene sobre el dibujo.
La función hist sirve para crear el histograma a una variable cuantitativa. Como argumentos esta función recibe un vector con los datos y opcionalmente puede ingresarse como argumento adicional el número de intervalos a ser graficados o en su defecto el número de intervalos se determina con la fórmula de Sturges.
Los programas de computador usualmente construyen los histogramas automáticamente, sin embargo, un buen programa debe permitirnos elegir el número de intervalos del histograma. Si usted posee un programa que no le permite hacer cambios, cambie de programa.
La estructura de la función hist con los argumentos más comunes de uso se muestran a continuación.
x: vector numérico de valores para construir el histograma.
breaks: puede ser un número entero que indica el número aproximado de clases o un vector cuyos elementos indican los límites de los intervalos.
freq: argumento lógico; si se especifica como TRUE, el histograma presentará frecuencias absolutas o conteo de datos para cada intervalo; si se especifica como FALSE el histograma presentar las frecuencias relativas (en porcentaje). Por defecto, este argumento toma el valor de TRUE siempre y cuando los intervalos sean de igual ancho.
include.lowest: argumento lógico; si se especifica como TRUE, un x[i] igual a los equal a un valor breaks se incluirá en la primera barra, si el argumento right = TRUE, o en la última en caso contrario.
right: argumento lógico; si es TRUE, los intervalos son abiertos a la izquierda y cerrados a la derecha \((a,b]\). Para la primera clase o intervalo si include.lowest=TRUE el valor más pequeño de los datos será incluido en éste. Si es FALSE los intervalos serán de la forma \((a,b]\) y el argumento include.lowest=TRUE tendrá el significado de incluir elmás alto’’.
col: para definir el color de las barras. Por defecto, NULL produce barras sin fondo.
border: para definir el color de los bordes de las barras.
plot: argumento lógico. Por defecto es TRUE, y el resultado es el gráfico del histograma; si se especifica como FALSE el resultado es una lista de conteos por cada intervalo.
labels: argumento lógico o carácter. Si se especifica como TRUE coloca etiquetas arriba de cada barra.
…: parámetros gráficos adicionales a title y axis.
Vamos a construir varios histogramas para los tiempos de 180 corredores de la media maratón de CONAVI realizada hace algunos años. A continuación se muestra la forma de ingresar los 180 datos.
maraton <- c(
10253, 10302, 10307, 10309, 10349, 10353, 10409, 10442, 10447,
10452, 10504, 10517, 10530, 10540, 10549, 10549, 10606, 10612,
10646, 10648, 10655, 10707, 10726, 10731, 10737, 10743, 10808,
10833, 10843, 10920, 10938, 10949, 10954, 10956, 10958, 11004,
11009, 11024, 11037, 11045, 11046, 11049, 11104, 11127, 11205,
11207, 11215, 11226, 11233, 11239, 11307, 11330, 11342, 11351,
11405, 11413, 11438, 11453, 11500, 11501, 11502, 11503, 11527,
11544, 11549, 11559, 11612, 11617, 11635, 11655, 11731, 11735,
11746, 11800, 11814, 11828, 11832, 11841, 11909, 11926, 11937,
11940, 11947, 11952, 12005, 12044, 12113, 12209, 12230, 12258,
12309, 12327, 12341, 12413, 12433, 12440, 12447, 12530, 12600,
12617, 12640, 12700, 12706, 12727, 12840, 12851, 12851, 12937,
13019, 13040, 13110, 13114, 13122, 13155, 13205, 13210, 13220,
13228, 13307, 13316, 13335, 13420, 13425, 13435, 13435, 13448,
13456, 13536, 13608, 13612, 13620, 13646, 13705, 13730, 13730,
13730, 13747, 13810, 13850, 13854, 13901, 13905, 13907, 13912,
13920, 14000, 14010, 14025, 14152, 14208, 14230, 14344, 14400,
14455, 14509, 14552, 14652, 15009, 15026, 15242, 15406, 15409,
15528, 15549, 15644, 15758, 15837, 15916, 15926, 15948, 20055,
20416, 20520, 20600, 20732, 20748, 20916, 21149, 21714, 23256)
Los datos están codificados como por seis números en el formato hmmss, donde h corresponde a las horas, mm a los minutos y ss a los segundos necesarios para completar la maratón. Antes de construir los histogramas es necesario convertir los tiempos anteriores almacenados en maraton a horas, para esto se utiliza el siguiente código.
horas <- maraton %/% 10000
min <- (maraton - horas * 10000) %/% 100
seg <- maraton - horas * 10000 - min * 100
Tiempos <- horas + min / 60 + seg / 3600
A continuación se muestra el código para construir cuatro histogramas con 2, 4, 8 y 16 intervalos para los tiempos a partir de la variable Tiempos.
par(mfrow=c(2,2))
hist(x=Tiempos, breaks=2, main="", xlab="Tiempo (horas)",
ylab="Frecuencias", las=1)
mtext("(A)", side=1, line=4, font=1)
hist(x=Tiempos, breaks=4, main="", xlab="Tiempo (horas)",
ylab="Frecuencias", las=1)
mtext("(B)", side=1, line=4, font=1)
hist(x=Tiempos, breaks=8, main="", xlab="Tiempo (horas)",
ylab="Frecuencias")
mtext("(C)", side=1, line=4, font=1)
hist(x=Tiempos, breaks=16, main="", xlab="Tiempo (horas)",
ylab="Frecuencias")
mtext("(D)", side=1, line=4, font=1)
En la Figura se presentan los cuatro histogramas. El histograma C, con 8 barras, muestra más claramente la asimetría (este es el que la mayoría de los programas produce por defecto, ya que la regla de Sturges para este conjunto de datos aproxima a 8 barras). Si consideramos más barras por ejemplo 16, como tenemos en D, se refina más la información y empezamos a notar multimodalidad. En el código anterior se incluyó las = 1 para conseguir que los número del eje Y queden escritos de forma horizontal.
A continuación vamos a construir cuatro histogramas: el primero con dos intervalos intervalos y puntos de corte dados por el mínimo, la mediana y el máximo; el segundo con tres intervalos y puntos de corte dados por el mínimo, cuartiles 1, 2, 3 y máximo; el cuarto con diez intervalos y puntos de corte dados por los deciles; y el último con veinte intervalos y puntos de corte dados por cuantiles 5, 10, , 95. En el código mostrado a continuación se presenta la creación de los puntos de corte y los cuatro histogramas.
puntos1 <- c(quantile(Tiempos, probs=c(0, 0.5, 1)))
puntos2 <- c(quantile(Tiempos, probs=c(0, 0.25, 0.5, 0.75, 1)))
puntos3 <- c(quantile(Tiempos, probs=seq(0, 1, by=0.1)))
puntos4 <- c(quantile(Tiempos, probs=seq(0, 1, by=0.05)))
par(mfrow=c(2, 2))
hist(Tiempos, breaks=puntos1, freq=FALSE, ylim=c(0,2), labels=TRUE,
main="", ylab="Densidad")
mtext("(A)", side=1, line=4, font=1)
hist(Tiempos, breaks=puntos2, freq=FALSE, ylim=c(0,2), labels=TRUE,
main="", ylab="Densidad")
mtext("(B)", side=1, line=4, font=1)
hist(Tiempos, breaks=puntos3, freq=FALSE, ylim=c(0,2),
main="", ylab="Densidad")
mtext("(C)", side=1, line=4, font=1)
hist(Tiempos, breaks=puntos4, freq=FALSE, ylim=c(0,2),
main="", ylab="Densidad")
mtext("(D)", side=1, line=4, font=1)
Figure: Histogramas para el tiempo en la media maratón de CONAVI. A: histograma con dos intervalos, B: histograma con cuatro intervalos, C: histograma con diez intervalos, C: histograma con veinte intervalos.
Nota: En estos histogramas, las alturas corresponden a las intensidades (frec. relativa/long. intervalo). Por tanto, el área de cada rectángulo da cuenta de las frecuencias relativas. Para el caso (A) ambos intervalos tienen igual área y cada uno contiene 50% de los datos. esto puede verificarse así:
Los gráficos cuantil cuantil (quantile-quantile plot) son una ayuda para explorar si un conjunto de datos o muestra proviene de una población con cierta distribución.
La función qqnorm sirve para explorar la normalidad de una muestra mientras que la función qqplot es de propósito más general, sirve para crear el gráfico cuantil cuantil para cualquier distribución.
La estructura de las funciones con los argumentos más comunes de uso se muestran a continuación.
La función __qqnorm_ sólo necesita que se le ingrese el vector con la muestra por medio del parámetro y, la función qqplot necesita de la muestra en el parámetro y y que se ingrese en el parámetro x los cuantiles de la población candidata.
Existe otra función útil y es qqline, esta función sirve para agregar una línea de referencia al gráfico cuantil cuantil obtenido con qqnorm.
Ejemplo Simular 30 observaciones de una distribución \(N(?=10,?=4)\) y construir el gráfico cuantil cuantil.
El código para simular la muestra y crear el gráfico cuantil cuantil se muestra a continuación.
muestra <- rnorm(n=30, mean=10, sd=4)
par(mfrow=c(1, 2))
qqnorm(y=muestra)
qqline(y=muestra)
qqnorm(y=muestra, main='', ylab='Cuantiles muestrales',
xlab='Cuantiles teóricos', las=1)
qqline(y=muestra, col='blue', lwd=2, lty=2)
En la izquierda de la Figura está el gráfico cuantil cuantil sin editar, en la derecha se encuentra el gráfico luego de modificar los nombres de los ejes, grosor y color de la línea de referencia.
Simular 100 observaciones de una distribución \[Weibull(1,1)\] y construir dos gráficos cuantil cuantil, el primero tomando como referencia los cuantiles de una \[N(0,1)\] y el segundo tomando los cuantiles de la \[Weibull(1,1)\].
El código para simular la muestra y crear los gráficos cuantil cuantil se muestra a continuación.
n <- 100
muestra <- rweibull(n=n, shape=1, scale=1)
par(mfrow=c(1, 2))
qqplot(y=muestra, x=qnorm(ppoints(n)))
qqplot(y=muestra, x=qweibull(ppoints(n), shape=1, scale=1))
En la Figura están los gráficos cuantil cuantil solicitados. Del pánel izquierdo de la figura vemos que los puntos NO están alineados, esto indica que la muestra no proviene de la distribución \(N(0,1)\), esto es un resultado esperado ya que sabíamos que la muestra no fue generada de una normal. En el pánel derecho de la misma figura vemos que los puntos SI están alineados, esto indica que la muestra generada puede provenir de una población \[Weibull(1,1)\]. Los nombres de los ejes en la Figura anterior pueden ser editados para presentar una figura con mejor apariencia.
Los gráficos de densidad son muy útiles porque permiten ver el(los) intervalo(s) donde una variable cuantitativa puede ocurrir con mayor probabilidad.
La función density crea la información de la densidad y la función plot dibuja la densidad.
La estructura de la función density con los argumentos más comunes de uso se muestra a continuación.
Los argumentos de la función density son:
x: vector con los datos para los cuales se quiere la densidad.
bw: ancho de banda.
kernel: núcleo de suavización a usar, los posibles valores son gaussian, rectangular, triangular, epanechnikov, biweight, cosine o optcosine, el valor por defecto es gaussian.
na.rm: valor lógico, si es TRUE se eliminan los valores con NA para construir la densidad, el valor por defecto es FALSE.
Simular mil observaciones de una \[N(0,1)\], aplicar la función density al vector y explorar el contenido de la salida.
Primero se generan las observaciones y se almacenan en el objeto y, luego se aplica la función density y el resultado se guarda en el objeto res, para explorar lo que almacena res se usa la función names. A continuación el código utilizado.
y <- rnorm(n=1000)
res <- density(y)
names(res)
## [1] "x" "y" "bw" "n" "call" "data.name"
## [7] "has.na"
De la salida anterior se observa que la lista res tiene 7 elementos, los dos primeros son los vectores con las coordenadas para dibujar la densidad, los restantes elementos con información adicional.
Con los datos generados en el ejemplo anterior construir la densidad para varios núcleo y para varios valores de ancho de banda.
En el siguiente código se construyen 4 densidades para diferentes núcleos.
par(mfrow=c(2, 2))
plot(density(y, kernel='gaussian'))
plot(density(y, kernel='triangular'))
plot(density(y, kernel='cosine'))
plot(density(y, kernel='rectangular'))
En la Figura se muestran las densidades para 4 elecciones del núcleo. En la práctica se usa el núcleo que está por defecto (gaussian) ya que el objetivo de una densidad es ver la zonas donde es más probable encontrar observaciones de la variable.
En el siguiente código se construyen 4 densidades para diferentes anchos de banda.
par(mfrow=c(2, 2))
plot(density(y, bw=0.1))
plot(density(y, bw=0.2241)) # bw obtenido antes
plot(density(y, bw=0.5))
plot(density(y, bw=1))
En la Figura se muestran las densidades para 4 elecciones del parámetro ancho de banda bw, el valor de 0.2241 fue el valor calculado automáticamente por R y fue obtenido de la Figura anterior, los otros valores fueron elegidos arbitrariamente para ver los cambios en la densidad. El usar un ancho de banda pequeño la densidad queda muy rugosa y usar un valor muy grande la suaviza, se recomienda usar el valor automático.
Construir un gráfico de densidad para la variable peso corporal de la base de datos medidas_cuerpo, luego construir la densidad para la misma variable pero diferenciando por sexo.
url <- 'https://raw.githubusercontent.com/fhernanb/datos/master/medidas_cuerpo'
datos <- read.table(file=url, header=T)
par(mfrow=c(1, 2))
plot(density(datos$peso), main='Distribución del peso corporal',
xlab='Peso corporal (kg)', ylab='Densidad', lwd=4)
den.hom <- with(datos, density(peso[sexo == 'Hombre']))
den.muj <- with(datos, density(peso[sexo == 'Mujer']))
plot(den.hom, xlim=c(20, 120),
main='Distribución del peso corporal por género', ylab='Densidad',
xlab='Peso corporal (kg)', lwd=4, col='blue')
lines(den.muj, lwd=4, col='red')
legend('topright', legend=c('Hombres', 'Mujeres'), bty='n',
lwd=3, col=c('blue', 'red'))
En el panel izquierdo de la Figura se muestra la densidad para la variable peso, de esta figura se observa que tiene dos sectores de mayor densidad, alrededor de 50 kg y alrededor de 80 kg. En el panel izquierdo están la densidades del peso corporal para hombres y mujeres, aquí se observa claramente la diferencia entre los pesos de hombres y mujeres.
En este capítulo se presentan funciones para la creación de gráficos que involucran varias variables cuantitativas.
Los gráficos de dispersión son muy útiles porque permiten ver la relación que existe entre dos variables cuantitativas, la función plot permite crear este tipo de gráficos. La estructura de la función plot con los argumentos más usuales se muestra a continuación:
Los argumentos de la función plot son:
‘p’ para obtener puntos, esta es la opción por defecto.
‘l’ para obtener líneas.
‘b’ para obtener los puntos y líneas que unen los puntos.
‘c’ para obtener sólo las líneas y dejando los espacios donde estaban los puntos obtenidos con la opción ‘b’.
‘o’ para obtener los puntos y lineas superpuestas.
‘h’ para obtener líneas verticales desde el origen hasta el
valor
$yi$ de cada punto, similar a un histograma.
‘s’ para obtener escalones.
‘S’ similar al anterior.
‘n’ para que no dibuje.
Crear 16 parejas de puntos tales que \(x = − 5 , − 4 , ... , 9 , 10\) con \(y = − 10 + ( x − 3 )^2\) , dibujar los nueve diagramas de dispersión de \(y\) contra \(x\) usando todos los valores posibles para el parámetro type.
A continuación se muestra el código para crear las 16 parejas de \(x\) e \(y\) . Los nueve diagramas de dispersión se observan en la Figura , de esta figura se observa claramente el efecto que tiene el parámetro type en la construcción del diagrama de dispersión.
x <- -5:10
y <- -10 + (x-3)^2
par(mfrow=c(3, 3))
plot(x=x, y=y, type='p', main="con type='p'")
plot(x=x, y=y, type='l', main="con type='l'")
plot(x=x, y=y, type='b', main="con type='b'")
plot(x=x, y=y, type='c', main="con type='c'")
plot(x=x, y=y, type='o', main="con type='o'")
plot(x=x, y=y, type='h', main="con type='h'")
plot(x=x, y=y, type='s', main="con type='s'")
plot(x=x, y=y, type='S', main="con type='S'")
plot(x=x, y=y, type='n', main="con type='n'")
Efecto del parámetro type en la función
plot.
Como ilustración vamos a crear dos diagramas de dispersión entre el precio de apartamentos usados en la ciudad de Medellín y el área de los apartamentos. El primer diagrama es un gráfico básico mientras que el segundo es un diagrama mejorado.
El código necesario para cargar la base de datos y construir los diagramas de dispersión se muestra a continuación.
url <- 'https://tinyurl.com/hwhb769'
datos <- read.table(file=url, header=T)
par(mfrow=c(1, 2))
plot(x=datos$mt2, y=datos$precio)
plot(x=datos$mt2, y=datos$precio, pch='.', xlab='Área del apartamento (m2)', ylab='Precio (millones de pesos)')
Figura: Diagrama de dispersión del precio del apartamento versus área del apartamento. A la izquierda el diagrama de dispersión sin editar y a la derecha el diagrama de dispersión mejorado.
En la Figura se presentan los diagramas de dispersión entre precio y área de los apartamentos. En el panel de la izquierda está el diagrama básico, de este diagrama se observa claramente que a medida que los apartamentos tienen mayor área el precio promedio y la variabilidad del precio aumentan. En la parte inferior izquierda de este diagrama se tiene una zona de alta densidad de puntos y por esa razón se observa una mancha en el diagrama. Para la construcción del diagrama de dispersión mostrado en el panel derecho se usó el parámetro pch=‘.’ con el objetivo de obtener pequeños puntos que representen cada apartamento y que no se traslapen debido a que se tienen 694 observaciones en la base de datos.
La función sunflowerplot sirve para crear gráficos de dispersión en los cuales hay parejas repetidas que se superponen y que por lo tanto no se podrían apreciar.
Si una pareja está una sola vez se representará por un punto; si la pareja se repite dos veces se representará por un punto y dos rayitas rojas (pétalos); si la pareja se repite tres veces se representará por un punto y tres rayitas rojas (pétalos) y así sucesivamente.
La estructura de la función sunflowerplot con los argumentos más usuales se muestra a continuación:
Los argumentos de la función plot son:
x: vector numérico con las coordenadas del eje horizontal.
y: vector numérico con las coordenadas del eje vertical.
pch: valor o vector numérico con el tipo de punto a usar, por defecto pch=1. Para conocer los diferentes símbolos que se pueden obtener con el parámetro pch.
seg.col: color de los pétalos.
seg.lwd: ancho de los pétalos.
…: otros parámetros gráficos que pueden ser pasados como argumentos para plot.
A continuación se presenta un ejemplo en los cuales se usa la función sunflowerplot para crear gráficos de dispersión.
Suponga que se tiene el punto \((1, 1)\) una sola vez, el punto \((2, 3)\) repetido cuatro veces, el punto \((3, 2)\) repetido seis veces y el punto \((4, 5)\) repetido dos veces. El objetivo es crear dos gráficos de dispersión con los trece puntos, uno con la función plot y el otro con la función sunflowerplot para comparar los resultados.
A continuación se muestra el código con la creación de los vectores y los dos gráficos de dispersión con plot y sunflowerplot.
x <- c(1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 4, 4)
y <- c(1, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 5, 5)
par(mfrow=c(1, 2))
plot(x, y, main='a) Usando plot')
sunflowerplot(x, y, seg.col='blue', seg.lwd=2, cex=1.5, main='b) Usando sunflowerplot')
En la Figura se muestran los dos gráficos de dispersión para los datos. En el panel de la izquierda está el diagrama obtenido con la función plot pero sólo se observan 4 puntos cuando en realidad eran 13, esto se debe a que hay puntos repetidos y éstos quedan unos sobre otros. En el panel de la derecha está el diagrama obtenido con la función sunflowerplot, en este diagrama se observan cuatro objetos: el punto con coordenadas (1, 1) se representó con un punto ya que él no se repite, el punto (2, 3) está representado por un punto y cuatro pétalos de color azul ya que él se repite cuatro veces, y los puntos (3, 2) y (4, 5) se represetan con seis y dos pétalos porque esas son las veces que ellos se repiten.
Como ilustración vamos a crear un diagrama de dispersión con la función sunflowerplot para las variables número de alcobas y número de baños de la base de datos de apartamentos usados.
El código necesario para cargar la base de datos y construir el diagrama de dispersión se muestra a continuación. El diagrama de dispersión se construyó excluyendo la información del apartamento 594 porque éste aparece con 14 alcobas.
url <- 'https://tinyurl.com/hwhb769'
datos <- read.table(file=url, header=T)
datos <- datos[-594, ] # Sacando el apto 594
sunflowerplot(x=datos$alcobas, y=datos$banos)
En la Figura se muestra el diagrama de dispersión entre número de baños
y número de alcobas. De este diagrama se observa que la mayor parte de
los apartamentos de la base de datos tienen 2 o 3 alcobas con 2 o 3
baños, se nota que sólo 8 apartamentos tiene 5 alcobas y que sólo un
apatamento tiene seis baños.
La función symbols sirve para construir diagramas de dispersión en dos dimensiones incluyendo información adicional de variables cuantitativas.
Como ilustración vamos a crear un diagrama de dispersión entre el precio de apartamentos usados en la ciudad de Medellín y el área de los apartamentos pero incluyendo otras variables.
url <- 'https://tinyurl.com/hwhb769'
datos <- read.table(url, header=T)
subdat <- subset(datos, ubicacion == 'centro')
par(mfrow=c(1, 2))
with(subdat, symbols(x=mt2, y=precio, circles=alcobas, las=1, inches=0.4, fg='red', main='Radio = N° alcobas'))
with(subdat, symbols(x=mt2, y=precio, squares=alcobas, las=1, inches=0.4, fg='dodgerblue4', bg='pink', main='Lado = N° alcobas'))
Diagrama de dispersión con los símbolos circle y squares para incluir
más variables.
par(mfrow=c(1, 2))
with(subdat, symbols(x=mt2, y=precio, rectangles=cbind(estrato, alcobas), las=1, inches=0.4, fg='chartreuse4', bg='yellow', main='Ancho = estrato y Alto = N° alcobas'))
with(subdat, symbols(x=mt2, y=precio, stars=cbind(estrato, alcobas, banos), las=1, inches=0.4, fg='purple1', bg='tomato', main='Estrato, alcobas y n° baños'))
Diagrama de dispersión con los símbolos rectangles y stars para
incluir más variables.
Las matrices de dispersión obtenidas con la función pairs proporcionan un método simple de presentar las relaciones entre pares de variables cuantitativas y son la versión múltiple de la función plot. Este gráfico consiste en una matriz donde cada entrada presenta un gráfico de dispersión sencillo. Un inconveniente es que si tenemos muchas variables el tamaño de cada entrada se reduce demasiado impidiendo ver con claridad las relaciones entre los pares de variables. La celda \(( i , j )\) de una matriz de dispersión contiene el gráfico de dispersión de la columna \(i\) versus la columna \(j\) de la matriz de datos.
La estructura de la función pairs con los argumentos más usuales se muestra a continuación:
Los argumentos de la función pairs son:
x: matriz o marco de datos con la información de las variables cuantitativas a incluir en la matriz de dispersión.
labels: vector opcional con los nombres a colocar en la diagonal, por defecto se usan los nombres de columnas del objeto x.
panel: función usual de la forma function(x,y,…) a ser usada para determinar el contenido de los páneles. Por defecto es points, indicando que se graficarán los puntos de los pares de variables. Es posible utilizar aquí otras funciones diseñadas por el usuario.
…: Indica que es posible agregar otros parámetros gráficos, tales como pch y col, con los cuales puede especificarse un vector de símbolos y colores a ser usados en los scatterplots.
lower.panel, upper.panel: función usual de la forma function(x,y,…) para definir lo que se desea dibujar en los paneles abajo y arriba de la diagonal.
diag.panel: función usual de la forma function(x,y,…) para definir lo que se desea dibujar en la diagonal.
text.panel: Es opcional. Permite que una función: function(x, y, labels, cex, font, …) sea aplicada a los paneles de la diagonal.
label.pos: Para especificar la posición
y de las etiquetas en el text panel.
cex.labels, font.labels: Parámetros para la expansión de caracteres y fuentes a usar en las etiquetas de las variables.
row1attop: Parámetro lógico con el cual se especifica si el gráfico para especificar si el diseño lucirá como una matriz con su primera fila en la parte superior o como un gráfico con fila uno en la base. Por defecto es lo primero.
Dibujar una matriz de dispersión para las variables precio, área, número de alcobas y número de baños de la base de datos sobre apartamentos en Medellín.
A continuación se muestra el código usado para crear el gráfico solicitado. El objeto datos corresponde a la base de datos completa mientras que datos.num es el marco de datos sólo con las variables de interés precio, área, número de alcobas y número de baños.
url <- 'https://tinyurl.com/hwhb769'
datos <- read.table(file=url, header=T)
datos.num <- datos[, c('precio', 'mt2', 'alcobas', 'banos')]
pairs(datos.num)
Figura: Matriz de dispersión para las variables precio, área, número de alcobas y número de baños de la base de datos sobre apartamentos en Medellín.
Volver a construir la Figura anterior editando los nombres de las variables, usando cruces azules en lugar de puntos, en escala logaritmica, con marcas horizontales en el eje vertical y eliminando los diagramas de dispersión abajo de la diagonal.
Para obtener la nueva matriz de dispersión con los cambios solicitados se usa el siguiente código. En la Figura se presenta la nueva matriz de dispersión.
pairs(datos.num, lower.panel=NULL, cex.labels=1.5, log='xy',
main='Matriz de dispersión', las=1,
labels=c('Precio', 'Área', 'Num alcobas', 'Num baños'),
pch=3, cex=0.6, col='dodgerblue2')
Matriz de dispersión modificando los parámetros adicionales de la
función pairs.
Ejemplo Construir una matriz de dispersión con las variables precio,
área y avaluo para apartamentos que cumplan la condición
\(100m^2<area<130m^2\).
Adicionalmente, se deben diferenciar los apartamentos sin parqueadero
con color rojo y los apartamentos con parqueadero con color verde.
Para crear una matriz de dispersión se puede tambien usar la base de datos original llamada datos que contiene todas las variables y usar una fórmula con la ayuda del operador ~ para indicar las variables de interés. La fórmula NO debe contener nada del lado izquierdo mientras que en el lado derecho se colocan todas las variables a considerar en la matriz de dispersión, por esta razón es que en el códido mostrado abajo se inicia con la instrucción \(~ precio + mt2 + avaluo\). Para incluir condiciones se usa el parámetro subset de la siguiente manera: \(subset=mt2 > 100 & mt2 < 130\) . A continuación el código completo para construir la matriz de dispersión solicitada.
col1 <- ifelse(datos$parqueadero == 'no', 'red', 'green3')
pairs(~ precio + mt2 + avaluo, data=datos,
lower.panel=NULL, col=col1,
subset=mt2 > 100 & mt2 < 130, pch=19, cex=0.8,
main="Matriz de dispersión para aptos con
100 < área < 130 mt2")
En la Figura se presenta la matriz de dispersión solicitada, los puntos
rojos representan los apartamento sin parqueadero mientras que los
puntos verdes son los apartamento que si tienen parqueadero.
¿Es posible agregar una leyenda a una matriz de dispersión?
Claro que es posible, se construye la matriz de dispersión y se deja en el lienzo del dibujo un espacio para colocar la leyenda. A continuación se muestra un ejemplo disponible en Stackoverflow.
pairs(iris[1:4], main="Anderson's Iris Data -- 3 species",
pch=21, bg=c("red", "green3", "blue")[iris$Species],
oma=c(4, 4, 6, 12))
par(xpd=TRUE)
legend(0.85, 0.7, as.vector(unique(iris$Species)), bty='n',
fill=c("red", "green3", "blue"))
Figura : Matriz de dispersión con leyenda.
¿Es posible modificar el contenido de los páneles de una matriz de dispersión?
Claro que es posible, para hacer esto se definen funciones que hagan lo que se desea ver tanto en la diagonal como arriba y abajo de la misma.
Como ejemplo vamos a construir una matriz de dispersión que cumpla:
sobre la diagonal un diagrama de dispersión para las variables involucradas y la recta de regresión ajustada,
en la diagonal un histograma para la variable,
debajo de la diagonal el coeficiente de correlación entre las variables involucradas y usando un tamaño de fuente proporcional a la fuerza de correlación.
Para obtener esta matriz de dispersión especial se definen a
continuación las funciones panel.reg,
panel.hist y panel.cor, a continuación el
código utilizado. Luego se usa la función pairs y se indica
qué función debe actuar en cada uno de los parámetros
upper.panel, diag.panel y
lower.panel.
# Función para dibujar los puntos y agregar la recta de regresión
panel.reg <- function (x, y)
{
points(x, y, pch=20)
abline(lm(y ~ x), lwd=2, col='dodgerblue2')
}
# Función para crear el histograma
panel.hist <- function(x, ...)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(usr[1:2], 0, 1.5) )
h <- hist(x, plot = FALSE)
breaks <- h$breaks; nB <- length(breaks)
y <- h$counts; y <- y/max(y)
rect(breaks[-nB], 0, breaks[-1], y, col="dodgerblue2", ...)
}
# Función para obtener la correlación
panel.cor <- function(x, y, digits=2, prefix="", cex.cor)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(0, 1, 0, 1))
r <- abs(cor(x, y))
txt <- format(c(r, 0.123456789), digits=digits)[1]
txt <- paste(prefix, txt, sep="")
if(missing(cex.cor)) cex <- 0.8/strwidth(txt)
text(0.5, 0.5, txt, cex = cex * r)
}
pairs(datos.num,
upper.panel = panel.reg,
diag.panel = panel.hist,
lower.panel = panel.cor)
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
## Warning in par(usr): argument 1 does not name a graphical parameter
En la Figura se presenta la matriz de dispersión con las modificaciones
en cada uno de los páneles. Cualquier usuario puede modificar las
funciones
panel.reg, panel.hist y panel.cor para
personalizar la apariencia de los contenidos.
La función panel.smooth está disponible en R para que el usuario pueda incluir arriba o abajo de la diagonal un diagrama de dispersión con una línea resultado de un ajuste suavizado. Abajo se muestra el código de cómo incluir la función panel.smooth y en la Figura se muestra gráfico obtenido.
pairs(datos.num,
upper.panel = panel.reg,
diag.panel = panel.hist,
lower.panel = panel.smooth)
Matriz de dispersión usando la función panel.smooth.
Al construir un gráfico con pairs se recomienda no incluir demasiadas variables porque el gráfico se satura y no se pueden apreciar los patrones con facilidad. Se recomienda máximo 10 variables.
La función persp dibuja superfices en tres dimensiones y es posible rotar la superficie para obtener una perpectiva apropiada. La estructura de la función persp con los argumentos más usuales se muestra a continuación:
Los argumentos de la función plot son:
x: vector numérico con los valores de x donde fue evaluada la función o superficie.
y: vector numérico con los valores de y donde fue evaluada la función o superficie.
z: matriz que contiene las alturas z de la supercifie para cada combinación de x e y.
main: vector numérico con las coordenadas del eje vertical.
sub: vector numérico con las coordenadas del eje vertical.
theta, phi: ángulo para la visión de la superficie, theta para la dirección azimutal y phi para latitud. Ver Figura 3.14 para una ilustración de los ángulos.
r: distancia entre el centro de la caja de dibujo al punto de vista.
col: color de la superficie.
border: color para el borde de la superficie.
box: valor lógico para indicar si se quiere dibujar la caja que contiene la superficie, por defecto es TRUE.
axes: valor lógico para indicar si se desean marcas en los ejes y nombres de los ejes, por defecto es TRUE. Si box=‘FALSE’ no aparecen marcas ni nombres de los ejes.
expand: factor de expansión aplicado a los valores en el eje z.
ticktype: tipo de marcas a colocar en los ejes, simple no dibuja nada y detailed coloca números a los ejes.
nticks: número aproximado de marcas en los ejes.
Ilustración de los ángulos theta y phi para la función persp.
Dibujar la superficie asociada a la función \(f(x,y)=sen(x^2+y^2)\) para \(−2≤x≤2\) y \(−2≤y≤2\) . Usar 4 combinaciones de los parámetros \(theta\) y \(phi\) para obtener un buen punto de vista de la superficie.
Lo primero que se debe hacer es crear la función \(f(x,y)\) la cual se va a llamar fun. Luego se definen los vectores x e y tomando por ejemplo 25 puntos equiespaciados en el intervalo \([−2,2]\). Luego se usa la función outer para crear la rejilla o matriz que contiene los valores de \(f(x,y)\) para cada combinación de x e y, los resultados se almacenan en el objeto z. Por último se dibujan 4 perspectivas de la función variando los parámetros theta y phi de la función persp. A continuación el código utilizado.
fun <- function(x, y) sin(x^2 + y^2)
x <- seq(from=-2, to=2, length.out=25)
y <- seq(from=-2, to=2, length.out=25)
z <- outer(x, y, fun)
par(mfrow=c(2, 2), mar=c(1, 1, 2, 1))
persp(x, y, z, zlim=c(-1, 1.5), theta=0, phi=0, col='aquamarine',
main='(A) theta=0, phi=0')
persp(x, y, z, zlim=c(-1, 1.5), theta=15, phi=15, col='lightpink',
main='(B) theta=15, phi=15')
persp(x, y, z, zlim=c(-1, 1.5), theta=45, phi=30, col='yellow1',
main='(c) theta=45, phi=30')
persp(x, y, z, zlim=c(-1, 1.5), theta=60, phi=50, col='lightblue',
main='(D) theta=60, phi=50')
En la Figura se presentan las 4 perspectivas de la función f(x,y)=sen(x2+y2). De los 4 páneles se nota que (C) y (D) muestran mejor la superficie de interés.
Al aumentar el valor del parámetro length.out en la creación de los vectores x e y se obtendrá una rejilla más tupida, se recomienda modificar este valor para obtener una superficie apropiada.
Dibujar la superficie de una distribución normal bivariada con vector de medias \(μ=(5,12)^⊤\) , varianzas unitarias y covarianza con valor de -0.8. Explorar el efecto de los parámetros ticktype, nticks, expand, axes y box.
Primero se define el vector de medias y la matriz de varianzas y covarianzas, luego se carga el paquete mvtnorm que contiene la función dmvnorm que calcula la densidad dado el vector de medias y la matriz de varianzas y covarianzas. Se construye la función fun y se vectoriza para luego obtener las alturas de la superficie con la ayuda de outer. Por último se dibujan tres perspectivas diferentes para la densidad modificando los parámetros ticktype, nticks, expand, axes y box, a continuación el código usado.
media <- c(5, 12)
varianza <- matrix(c(1, -0.8, -0.8, 1), ncol=2)
require(mvtnorm)
## Loading required package: mvtnorm
fun <- function(x, y) dmvnorm(c(x, y), mean=media, sigma=varianza)
fun <- Vectorize(fun)
x <- seq(from=2, to=8, length.out=30)
y <- seq(from=9, to=15, length.out=30)
z <- outer(x, y, fun)
par(mfrow=c(1, 3), mar=c(1, 1, 2, 1))
persp(x, y, z, theta=30, phi=30, ticktype = "detailed", nticks=4)
persp(x, y, z, theta=30, phi=30, col='salmon1', expand=0.5, axes=FALSE)
persp(x, y, z, theta=30, phi=30, col='springgreen1', expand=0.2, box=FALSE)
En la Figura se presentan las 3 perspectivas para la densidad. Note los efectos que ticktype, nticks, expand, axes y box tienen sobre los dibujos de las perspectivas.
La función contour dibuja gráficos contornos. La estructura de la función contour con los argumentos más usuales se muestra a continuación:
Los argumentos de la función son:
x, y: vectores numéricos en los cuales se evaluó la función de interés para construir el objeto z. Ambos vectores deben estar ordenados.
z: matriz con las alturas de la función de interés, por lo general creada con la función outer.
xlim, ylim, zlim: límites de los ejes x, y e z respectivamente.
nlevels: número aproximado de niveles o cortes en la superficie a representar.
col: color a usar en las líneas de contornos.
La función contour tiene otros parámetros adicionales que el lector puede consultar en la ayuda usando help(contour).
Generar una muestra aleatoria de 50 observaciones de una distribución
normal con parámetros
\(μ=170\) y \(σ^2=25\) . Dibujar un gráfico de contornos
para la superficie de log-verosimilitud.
La muestra aleatoria se genera con el siguiente código.
y <- rnorm(n=50, mean=170, sd=5) # sd es desviación
y
## [1] 171.0481 171.6045 170.6749 165.0498 174.3570 170.7016 176.4470 165.0395
## [9] 165.7232 169.2311 169.7462 174.2024 172.7362 164.5975 174.7319 164.1579
## [17] 167.2196 157.2810 169.2308 175.3657 177.2424 172.5461 164.8191 172.7631
## [25] 168.9983 174.1448 174.9266 161.9388 180.9789 171.7199 166.8648 174.4082
## [33] 166.8582 174.4845 171.6742 167.6304 169.5858 168.5463 168.9841 164.0178
## [41] 168.4886 168.0983 172.7877 165.0556 176.3748 165.5551 162.7987 175.7504
## [49] 179.2729 171.9723
Para dibujar los contornos solicitados se debe primero construir la función de log-verosimilitud llamada ll. A continuación el código para crear ll, mayores detalles de cómo construir funciones de log-verosimilitud se pueden consultar en Hernández and Usuga (2018).
ll <- function(a, b) sum(dnorm(x=y, mean=a, sd=b, log=TRUE))
ll <- Vectorize(ll) # Para vectorizar la función
Una vez construída la función ll se deben construir los vectores con las coordenadas horizontal y vertical donde se evalua la función ll. En el código mostrado abajo se tienen dos vectores xx e yy obtenidos como secuencias desde el menor valor hasta el mayor valor para cada uno de los parámetros \(μ\) y \(σ\) de la distribución normal, el valor by=0.5 indica el tamaño de paso de la secuencia. Luego se construye la matriz zz usando la función outer evaluando ll en xx e yy. Por último la función contour se aplica sobre los elementos xx, yy e zz. En la Figura se muestra el gráfico de contornos con aproximadamente 50 niveles.
xx <- seq(from=160, to=180, by=0.5)
yy <- seq(from=3, to=7, by=0.5)
zz <- outer(X=xx, Y=yy, ll)
contour(x=xx, y=yy, z=zz, nlevels=50,
col=gray(0.3), lwd=2, lty='solid',
xlab=expression(mu), ylab=expression(sigma))
Gráfico de contornos para la función de log-verosimilitud para el ejemplo sobre normal.
La función filled.contour dibuja gráficos contornos pero usando una paleta de colores, este tipo de representación se denomina gráfico de nivel. La estructura de la función filled.contour con los argumentos más usuales se muestra a continuación:
Los argumentos de la función son:
La función filled.contour tiene otros parámetros adicionales que el lector puede consultar en la ayuda usando help(filled.contour).
Abajo se muestra el código para crear un gráfico de nivel usando los datos del ejemplo anterior. Observe que las figuras anterior y la siguiente son bastante similares.
filled.contour(x=xx, y=yy, z=zz, nlevels=20,
xlab=expression(mu), ylab=expression(sigma),
color = topo.colors)
Gráfico de nivel para la función de log-verosimilitud para el ejemplo sobre normal.
¿Es posible crear un gráfico de contornos y de niveles en la misma figura?
Claro que si, abajo se muestra el código para crear primero el gráfico de niveles (en colores) y luego superponer el gráficos de contornos (las lineas).
filled.contour(x=xx, y=yy, z=zz, nlevels=20,
color.palette = topo.colors,
plot.axes=contour(xx, yy, zz, nlevels=20, add=TRUE))
Gráfico de nivel y contornos para la función de log-verosimilitud para
el ejemplo sobre normal.
La función image dibuja un gráfico de calor similar al obtenido con la función filled.contour. La estructura de la función image con los argumentos más usuales se muestra a continuación:
Los argumentos de la función son:
x, y: vectores numéricos en los cuales se evaluó la función de interés para construir el objeto z. Ambos vectores deben estar ordenados.
z: matriz con las alturas de la función de interés, por lo general creada con la función outer.
Para la muestra aleatoria obtenida en el ejemplo anterior, dibujar un gráfico con image para la superficie de log-verosimilitud.
Usando los objetos xx, yy e zz creados en el ejemplo anterior se puede construir el gráfico solicitado, a continuación el código utilizado.
image(x=xx, y=yy, z=zz,
xlab=expression(mu), ylab=expression(sigma))
Gráfico para la función de log-verosimilitud para el ejemplo sobre
normal.
La función kde2d pertenece al paquete MASS y es utilizada para crear densidades para dos variables cuantitativas. La estructura de la función kde2d con los argumentos más usuales se muestra a continuación:
Los argumentos de la función son:
x: vector con la variable para el eje X.
y: vector con la variable para el eje Y..
h: vector con los anchos de banda en las direcciones X e Y.
n: número de puntos para construir la rejilla.
lims: límites del rectángulo de datos a considerar, debe ser un vector de la forma c(xl, xu, yl, yu). Este parámetro por defecto es c(range(x), range(y)).
La base de datos medidas del cuerpo cuenta con 6 variables registradas a un grupo de 36 estudiantes de la universidad, dos de esas variables son la altura y el peso corporal. Se desea construir un gráfico de densidad bivariada para altura y peso.
El código mostrado a continuación hace la lectura de la base de datos y luego se construyen dos densidades, la primera con n=5 y la segunda con n=50, esto para ver el efecto del parámetro n.
url <- 'https://raw.githubusercontent.com/fhernanb/datos/master/medidas_cuerpo'
datos <- read.table(file=url, header=T)
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
require(MASS) # Se debe cargar este paquete
f1 <- kde2d(x=datos$peso, y=datos$altura, n=5)
f2 <- kde2d(x=datos$peso, y=datos$altura, n=50)
En el código mostrado a continuación se dibujan las dos densidades usando un gráfico de calor usando la función image.
par(mfrow=c(1, 2))
image(f1, xlab='Peso', ylab='Estatura', main='n=5')
image(f2, xlab='Peso', ylab='Estatura', main='n=50')
Gráfico de densidad bivariada para el peso corporal y la estatura de
un grupo de estudiantes. A la izquierda la densidad con n=5
y a la derecha con n=50. Figure 3.21: Gráfico de densidad
bivariada para el peso corporal y la estatura de un grupo de
estudiantes. A la izquierda la densidad con n=5 y a la derecha con
n=50.
La función interaction.plot dibuja gráficos de interacción. La estructura de la función interaction.plot con los argumentos más usuales se muestra a continuación:
Los argumentos de la función son:
response: vector numérico con la variable respuesta.
x.factor: factor 1 a ubicar en el eje horizontal.
trace.factor: factor 2 para diferenciar las líneas.
fun: función a aplicar para a response para cada combinación de
x.factor y trace.factor.
legend: valor lógico para incluir o no leyenda.
trace.label: nombre a colocar en la leyenda.
La función interaction.plot tiene otros parámetros adicionales que el lector puede consultar en la ayuda usando help(interaction.plot).
Se realizó un experimento para determinar cómo influye el material de la batería y la temperatura del medio ambiente sobre la duración en horas de la batería. Se desea construir un gráfico de interacción entre Temperatura y Material para ver el efecto sobre la duración promedio de las baterías. Los datos y el código para generar el gráfico solicitado se muestran a continuación.
horas <- c(130, 155, 74, 180, 150, 188, 159, 126, 138, 110, 168,
160, 34, 40, 80, 75, 136, 122, 106, 115, 174, 120, 150,
139, 20, 70, 82, 58, 25, 70, 58, 45, 96, 104, 82, 60)
temperatura <- rep(c(15, 70, 125), each=12)
material <- rep(1:3, each=4, times=3)
interaction.plot(x.factor=temperatura, trace.factor=material,
response=horas, trace.label='Material',
xlab='Temperatura',
ylab='Duración promedio (horas)',
col=c('dodgerblue3', 'chartreuse4', 'salmon3'),
fun=mean, lwd=2, las=1, fixed=T)
Gráfico de interacción entre Temperatura y Material sobre la duración promedio de las baterías.
Los gráficos de espagueti son usados para representar la evolución de una variable medida para un grupo de sujetos en diferentes momentos del tiempo. La función interaction.plot se puede usar para obtener este tipo de gráficos, a continuación un ejemplo.
El ejemplo aquí presentado fue tomado de este enlace. El objetivo es crear un gráfico de espagueti para mostrar la evolución de la variable tolerancia a través del tiempo para cada uno de los 16 individuos estudiados. El código para descargar la base de datos y construir el gráfico se muestran a continuación.
dt <- read.table("https://stats.idre.ucla.edu/stat/r/faq/tolpp.csv",
sep=",", header=T)
library(gplots)
## Warning: package 'gplots' was built under R version 4.2.3
##
## Attaching package: 'gplots'
## The following object is masked from 'package:stats':
##
## lowess
require(gplots) # Paquete especial para crear una paleta
palette(rich.colors(16)) # de 16 colores diferentes y suaves
interaction.plot(response=dt$tolerance,
x.factor=dt$time, col=1:16, lwd=2,
trace.factor=dt$id, las=1, lty=1,
xlab="Tiempo", ylab="Tolerancia", legend=F)
Gráfico de espagueti para ver la evolución de la variable tolerancia.
En este capítulo se presentan funciones para la creación de gráficos para variables cualitativas.
Los gráficos de barras son útiles para representar las frecuencias absolutas o relativas asociadas a los niveles de una variable cualitativa y la función barplot se usa para obtener un gráfico de barras. La estructura de la función barplot con los argumentos más comunes de uso se muestra a continuación.
Los argumentos de la función barplot son:
height: vector o matriz con la información de las frecuencias absolutas o relativas.
beside: valor lógico para indicar si las barras deben estar al lado cuando la información ingresada es una matriz.
horiz: valor lógico para indicar si el diagrama de barras debe ser horizontal, por defecto es FALSE.
La función barplot tiene otros parámetros que pueden ser consultados en la ayuda de la función por medio de la instrucción ?barplot.
Suponga que queremos construir un diagrama de barras para las frecuencias relativas de la variable estrato socioeconómico del apartamento de la base de datos sobre apartamentos usados en Medellín.
A continuación se muestra el código necesario para cargar la base de datos aptos2015. Antes de construir el diagrama de barras solicitado es necesario construir la tabla de frecuencias para la variable estrato, para esto se usa la función table y los resultados se almacenan en el objeto tabla1 que contiene las frecuencias absolutas. Para obtener las frecuencias relativas se usa luego la función prop.table sobre el objeto tabla1.
url <- 'https://tinyurl.com/hwhb769'
datos <- read.table(file=url, header=T)
tabla1 <- table(datos$estrato)
tabla1 <- prop.table(tabla1)
tabla1
##
## 2 3 4 5 6
## 0.01152738 0.23198847 0.19884726 0.20893372 0.34870317
Una vez se tiene el objeto con la información de las frecuencias relativas se puede dibujar el diagrama de barras usando el siguiente código.
barplot(tabla1, xlab='Estrato socioeconómico',
ylab='Frecuencia relativa', las=1)
En la Figura se presenta el diagrama de barras solicitado. Se observa que hay pocos apartamentos (1.15%) pertenecientes al estrato dos, los estratos tres, cuatro y cinco aportan porcentajes similares a la base de datos y que el estrato 6 es el que más apartamentos aporta a la base de datos, 34.87%.
Algunas veces se acostumbra a colocar las frecuencias relativas sobre la parte superior de las barras para facilitar la lectura. A continuación se presenta el código para replicar la Figura con las frecuencias relativas. Lo primero que se hace es dibujar el diagrama de barras y almacenar la información de él en el objeto xx para luego poder usar la ubicación de cada una de las barras. Note que se agregó también ylim=c(0, 0.45) para conseguir una ampliación del eje vertical, esto para lograr que se vea el número sobre la barra del estrato 6. Luego se usa la función text para incluir un texto en las coordenadas x=xx y y=tabla1, el parámetro pos=3 coloca el texto en la parte superior de las coordenadas y el parámetro label sirve para indicar lo que se desea escribir en las coordenadas indicadas, en este caso son las frecuencias relativas almacenadas en tabla1.
xx <- barplot(tabla1, ylim=c(0, 0.45), col=gray(0.9),
xlab='Estrato socioeconómico',
ylab='Frecuencia relativa', las=1)
text(x=xx, y=tabla1, pos=3, cex=0.8, col="red",
label=round(tabla1, 4))
En la Figura se muestra el diagrama de barras modificado. Note que si no se hubiese usado ylim=c(0, 0.45) al dibujar el diagrama, la marca sobre la última barra no se podría ver.
Suponga que queremos construir un diagrama de barras para comparar la variable presencia de parqueadero con el estrato socioeconómico en la base de datos sobre apartamentos usados en Medellín.
La función barplot también puede ser usada para representar una tabla de frecuencia con dos variables. Para obtener la tabla de frecuencia para relacionar parqueadero con estrato se usa el siguiente código.
tabla2 <- table(datos$parqueadero, datos$estrato)
tabla2
##
## 2 3 4 5 6
## no 5 88 24 8 1
## si 3 73 114 137 241
El anterior resultado es la tabla de contingencia entre las variables parqueadero y estrato, de esta tabla vemos que para estratos superiores el número de apartamentos que si tienen parqueadero es mayor que los apartamentos sin parqueadero. La tabla anterior se puede representar gráficamente usando el siguiente código.
barplot(tabla2)
En la Figura se muestra el gráfico de barras sin editar, el color negro representa la frecuencia de los apartamentos sin parqueadero (no) y el color gris representa los apartamentos con parqueadero (si), las barras están una encima de la otra y la comparación no es tan clara como debería. Para mejorar la comparación se usa el parámetro besides=TRUE, a continuación el código utilizado.
barplot(tabla2, beside=TRUE)
En la Figura está el gráfico de barras obtenido agregando besides=TRUE para que las barras se ubiquen una junto a la otra. Este gráfico se puede mejorar aún más colocando una leyenda para identificar las barras, nombrando los ejes y usando otros colores, a continuación el código utilizado.
barplot(tabla2, beside = TRUE, las=1,
xlab='Estrato', ylab='Frecuencia',
col = c("lightblue", "mistyrose"),
ylim = c(0, 250))
legend('topleft', legend=rownames(tabla2), bty='n',
fill=c("lightblue", "mistyrose"))
En la Figura se observa el gráfico de barras solicitado, se observa claramente que en los estratos 4, 5 y 6 predominan los aparatamentos con parqueadero.
Es posible construir una tabla de contingencia de frecuencia relativa para ver cómo es el comportamiento de tener o no parquedadero dentro de cada estrato, el siguiente código construye la tabla3 con la información necesaria. La función prop.table permite obtener la tabla de frecuencias relativas a partir de una tabla de frecuencias absolutas, el parámetro margin sirve para indicar si las frecuencias relativas se deben obtener por fila (margin=1) o por columnas (margin=2).
tabla3 <- prop.table(tabla2, margin=2)
tabla3
##
## 2 3 4 5 6
## no 0.625000000 0.546583851 0.173913043 0.055172414 0.004132231
## si 0.375000000 0.453416149 0.826086957 0.944827586 0.995867769
De la anterior tabla se ve que el porcentaje de apartamentos con parqueadero supera enormemente el los apartamentos sin parqueadero para los estratos 6, 5 y 4. El código para generar un gráfico asociado a la tabla3 se muestra a continuación.
barplot(tabla3,
beside = TRUE, las=1,
xlab='Estrato', ylab='Frecuencia relativa',
col = c("lightblue", "mistyrose"),
ylim = c(0, 1))
legend('topleft', legend=rownames(tabla2), bty='n',
fill=c("lightblue", "mistyrose"))
Relación entre la presencia de parqueadero y el estrato socioeconómico.
En R es posible construir gráficos de pastel para representar una tabla de frecuencia relativa o absoluta, sin embargo este tipo de gráficos no es recomendable y en la ayuda de la función se hace la siguiente advertencia:
La estructura de la función pie con los argumentos más comunes de uso se muestra a continuación.
Los argumentos de la función pie son:
x: vector con elementos no negativos que representan las frecuencias de los niveles de la variable cualitativa.
labels: vector con los nombres a colocar en cada parte del pastel, por defecto se usan los nombres del vector x.
Dibujar un gráfico de pastel para las frecuencias relativas de la variable estrato socioeconómico del apartamento de la base de datos sobre apartamentos usados en Medellín.
La tabla1 construída en el primer ejemplo de barplot se utiliza para construir el gráfico solicitado. Abajo el código necesario para construir el gráfico.
nombres <- paste('Estrato ', 2:6)
pie(x=tabla1, labels=nombres,
main='Gráfico de pastel NO recomendado!!!')
La Figura presenta el gráfico de pastel construído con la instrucción anterior.
Los gráficos de puntos son útiles para representar tablas de frecuencia (de 1 o 2 vías) o tablas de resumen en relación a una o dos variables. La estructura de la función dotchart se muestra a continuación.
Los argumentos de la función dotchart son:
x: vector o matriz con la información de las frecuencias o medida de resumen a representar. Si x es una matriz las columnas representarán agrupaciones.
labels: vector con los nombres a usar para los puntos, por defecto toma los nombres de las filas de la matriz x.
groups: vector con los nombres a usar para los grupos, por defecto toma los nombres de las columnas de la matriz x.
cex: tamaño de los nombres a mostrar en los ejes.
pt.cex: tamaño del punto.
pch: tipo de punto a usar.
color: tipo de color usar para los puntos.
lcolor: color para la línea asociada a cada punto.
…: otros parámetros gráficos que pueden ser pasados como argumentos.
###Ejemplo
Suponga que queremos explorar el rendimiento de del combustible (mpg) para los 32 autos de la base de datos mtcars usando un diagrama de puntos.
A continuación se muestra el código necesario para crear el gráfico de puntos deseado.
mtcars <- mtcars[order(mtcars$mpg), ] # Para re-ordenar la base
dotchart(mtcars$mpg, labels=row.names(mtcars),
cex=0.6, xlab="mpg")
En la Figura se tiene el gráfico de puntos para la variable mpg, de esta figura se observa fácilmente que el Toyota Corolla tiene el mayor rendimiento de combustible y que el Cadillac Fleetwood el menor.
Suponga que se tiene una tabla de contingencia con la información del número de hombres y mujeres que sufren de una enfermedad rara en cuatro ciudades importantes del mundo, a continuación la matriz x con la información recolectada.
x <- matrix(c(4, 6, 30, 18, 7, 13, 35, 20),
ncol=4, byrow=TRUE)
colnames(x) <- c('Madrid', 'Londres', 'Paris', 'Miami')
rownames(x) <- c('Hombre', 'Mujer')
x
## Madrid Londres Paris Miami
## Hombre 4 6 30 18
## Mujer 7 13 35 20
En la salida anterior se presenta la matriz x, de esta matriz se observa que en Madrid hay 11 personas, 4 hombres y 7 mujeres que sufren de la enfermedad, las demás columnas se interpretan de forma similar.
Para construir un gráfico de puntos con el objetivo de presentar la información de la matriz x se utiliza el código mostrado abajo. En la Figura se presenta el gráfico de puntos obtenido y de esta figura se nota claramente que en París es donde hay más personas que sufren de la enfermedad.
dotchart(x=x)
El gráfico de puntos presentado en la Figura se puede mejorar usando los otros argumentos disponibles en la función dotchart, a continuación el código y en la Figura siguiente el resultado.
dotchart(x=x,
pt.cex=2, pch=c(8, 21), color=c('blue', 'red'),
lcolor='black',
xlab='Número de personas')
Gráfico de puntos mejorado para una tabla de contingencia de 2
variables.
Use funciones o procedimientos (varias líneas) de R para responder cada una de las siguientes preguntas.
Todas las preguntas siguientes están relacionadas con la base de datos sobre apartamentos usados en la ciudad de Medellín.
Construya un diagrama de barras para representar las frecuencias ABSOLUTAS de la variable ubicación.
Vuelva a construir el mismo diagrama de barras anterior pero de forma horizontal y agregando números de color azul para indicar las frecuencias.
Construya una tabla de dos vías para las variables ubicación y parqueadero.
Construya una tabla de frecuencias relativas para ver cómo se comporta la variable parqueadero dentro de cada ubicación.
Dibuje un diagrama de barras para la tabla de frecuencias del punto anterior.