Rayshader.En el presente código se repasa la lógica detrás del proceso de análisis mediante bootstrap y el graficado de coeficientes en 3D.
Para la creación de la base de Datos, descargaremos las bases de la librerÃa HistData, dentro de la cual se ubica la base de GaltonFamilies, la cual es una medición de las alturas de padres e hijos, como dice a continuación:
This data set lists the individual observations for 934 children in 205 families on which Galton (1886) based his cross-tabulation.
In addition to the question of the relation between heights of parents and their offspring, for which this data is mainly famous, Galton had another purpose which the data in this form allows to address: Does marriage selection indicate a relationship between the heights of husbands and wives, a topic he called assortative mating? Keen p. 297-298 provides a brief discussion of this topic.
# Removemos objetos previos
rm(list = ls(all = TRUE)) #clear workspace
# Librerias
library(HistData)
# Leemos los datos.
data("GaltonFamilies")
# Almacenamos la base en un objeto de nombre _d_.
d <- GaltonFamilies
# Generamos una var. dummy para cuando el individuo
# es masculino o femenino
d$H <- 1 # Caso Hombres
d$H[d$gender == "female"] <- 0 # Caso Mujeres
# GENERAMOS UNA FUNCION PARA HACER BOOTSTRAP #
my_boots <- function(){
# Sampleo
samp <- sample(1:nrow(d), size = nrow(d), replace = TRUE)
# Estimamos los parametros del modelo...
m <- glm(childHeight ~ H, data = d[samp,])
# Quedandonos con los coeficientes de m, y el sigma
c(coef(m), sigma = sigma(m))
}
# Guardamos las replicaciones de bootstrap de los parametros cada vez que hago
# una replicación de bootstrap
boots <- replicate(50000, my_boots())
boots <- t(boots)
head(boots)
## (Intercept) H sigma
## [1,] 64.18696 5.048456 2.438772
## [2,] 64.05541 5.271368 2.402322
## [3,] 64.19382 5.251919 2.432279
## [4,] 64.04742 5.317001 2.464453
## [5,] 64.13260 5.145109 2.479575
## [6,] 64.10022 5.160624 2.474592
A la base de coeficientes resultante del Bootstrap, le sacamos estadÃsticas descriptivas.
# Covarianza de boots
cov(boots)
## (Intercept) H sigma
## (Intercept) 0.0122066889 -0.0122232537 0.0001191759
## H -0.0122232537 0.0264409926 -0.0002895936
## sigma 0.0001191759 -0.0002895936 0.0038144836
# Varianza de boots
sqrt(diag(cov(boots)))
## (Intercept) H sigma
## 0.11048389 0.16260687 0.06176151
# Sacamos la desv est a la boots/alpha
sd(boots[,1])
## [1] 0.1104839
# Sacamos la desv est a la boots/beta
sd(boots[,2])
## [1] 0.1626069
# Sacamos los coeficientes de la muestra completa
colMeans(boots)
## (Intercept) H sigma
## 64.103877 5.129569 2.493916
Llego el momento del retrato hablado.
Vemos la distribución del alpha y beta.
# Distribución boostrap de alpha
plot(density(boots[,1], adjust = 4))
# Distribución bootstrap de beta
plot(density(boots[,2], adjust = 4))
La distribucion es multivariada porque las variables covarian conjuntamente. (Nota de la clase)
plot(boots[,1:2])
En la gráfica se muestra la relación entre los coeficientes obtenida en el proceso de bootstrap.
# GRAFICA EN 3D DE LINEAS.
library(hexbin)
# Otro plot de esas variables, pero con otra funcion
bin <- hexbin(boots[,1:2], xbins = 50)
plot(bin)
library(MASS)
den3d <- kde2d(boots[,1], boots[,2])
persp(den3d, box = FALSE)
# GRAFICA EN 3D - RAYSHADER #
# Libreria para hacer renders en 3D
library(rayshader)
# Libreria para hacer graficas y usar el simbolo %>%
library(tidyverse)
# Elaboramos la gráfica de hexagonos.
ggbins <- ggplot(boots %>% as.data.frame()) +
geom_hex(aes(x = boots[,1], y = boots[,2]))
# Convertimos una grafica ggplot en 2D a una en 3D
plot_gg(ggbins, # Objeto ggplot
width = 5, # Ancho del ggplot
height = 5, # Alto del ggplot
multicore = TRUE, # Arg. para utilizar todos los nucleos de la compu en el procesamiento
scale = 250, # Multiplicador para escala vertical del 3D
zoom = 0.7, # Zoom de la ventana
theta = 10, # Angulo horizontal de presentacion
phi = 30, # Angulo vertical de presentacion
windowsize = c(800, 800))
render_snapshot(clear = TRUE)