Generación de un gráfico en 3D mediante la librería 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.

Proceso de generación de la información.

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

Coeficientes promedio.

# 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)

Gráfica 2D.

plot(boots[,1:2])

Gráfica en 3D de líneas.

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

# 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)