Análisis de diamantes

library(tidyverse)
library(modelr)
options(na.action = na.warn) # nos indicará si falta algun dato

Modelo: los diamantes más malos son los más caros

Información del dataset “diamonds”

Para este caso trataremos con el dataset “diamonds”. Información general:

Format: A data frame with 53940 rows and 10 variables. A dataset containing the prices and other attributes of almost 54,000 diamonds. The variables are as follows:

• price: price in US dollars ($326–$18,823)

• carat: weight of the diamond (0.2–5.01

• cut: quality of the cut (Fair, Good, Very Good, Premium, Ideal)

• color: diamond colour, from J (worst) to D (best)

• clarity: a measurement of how clear the diamond is (I1 (worst), SI2, SI1, VS2, VS1, VVS2, VVS1, IF (best))

• x: length in mm (0–10.74)

• y: width in mm (0–58.9)

• z: depth in mm (0–31.8)

• depth: total depth percentage = z / mean(x, y) = 2 * z / (x + y) (43–79)

Business question

¿Por qué los diamantes de baja calidad son más caros?

# El precio de los diamantes de peor corte son de media más caros
diamonds %>%
  ggplot(aes(cut, price)) +
    geom_boxplot()

# Los diamantes con peor color tienen un precio muy elevado (color J es el peor)
diamonds %>%
  ggplot(aes(color, price)) +
    geom_boxplot()

# La peor claridad parece ser el más caro
diamonds %>%
  ggplot(aes(clarity, price)) +
    geom_boxplot()

Hipótesis y análisis gráfico

Hipótesis: puede que el peso del diamante sea el que define el precio:

diamonds %>%
  ggplot(aes(carat, price)) +
    geom_hex(bins = 50) # separados por color según densidad de muestras

Visualizando la variable predictora

Este gráfico nos muestra que la variable carat parece tener una relación muy clara con el precio. Para desarrollar nuestra toería procesaremos los datos, fijándonos solo en los “normales”, y dejando de lado los que se salen de lo “normal”, los que estan más alla de 2’5. Por lo tando filtramos los outliers:

diamonds2 <- diamonds %>%
  filter(carat <= 2.5) %>% # quitamos una parte muy pequeña de la información, menos del 3.5% del dataset
  mutate(lprice = log2(price), lcarat = log2(carat))# como vemos esta tendencia de crecer de forma exponencial con el kilate, vamos a transformar los datos con el logaritmo en base 2 para buscar así la relación lineal

diamonds2 %>%
  ggplot(aes(lcarat, lprice)) +
    geom_hex(bins = 50)

Con esta tranformación hemos visto que hay que hay una relación lineal.

NOTA: la tranformación logarítmica es muy útil cuando vemos que la tendencia crece muy rápido tanto en X como en Y, porque este crecimiento acelerado tiende a tranformarse en un patrón lineal, ayudándonos a crear un modelo que sea mas fácil trabajar (un modelo lineal).

Creación del modelo

Ahora podemos empezar a ver si podemos ajustar un modelo con la función lm():

mod_diamonds <- lm(lprice ~ lcarat, data = diamonds2)

# Que nos dice el modelo?
grid <- diamonds2 %>%
  data_grid(carat = seq_range(carat, 30)) %>% # hacmos 30 divisiones
  mutate(lcarat = log2(carat)) %>%
  add_predictions(mod_diamonds, "lprice") %>%
  mutate(price = 2 ^ lprice)
grid

Entendamos la parrilla que hemos hecho:

1. Tenemos carat en su valor original

2. lcarat como la transformación del logaritmo de carat

3. lprice es el la predicción que ha hecho nuestro modelo

4. price es la predicción final después de elevar a la 2 el la predicción del precio para volver al precio original en dólares

Visualicación del modelo

diamonds2 %>%
  ggplot(aes(carat, price)) +
    geom_hex(bins = 50) +
    geom_line(data = grid, color = "red", size = 1)

Este modelo nos dice que los diamantes que son más grandes, son mas baratos de los que realmente uno espera, dado que el máximo de precio no sobrepasa nunca los $20000.

Comprobación de los errores de los modelos

# Añadimos una columna con los residuos
diamonds2 <- diamonds2 %>%
  add_residuals(mod_diamonds, "lresid")

diamonds2 %>%
  ggplot(aes(lcarat, lresid)) +
    geom_hex(bins = 50)

Este modelo nos muestra que más o menos ha eliminado este patrón lineal. Todo está más o menos concentrado en torno al centro. Si que hay un patrón en la derecha del gráfico, pero el color es un azul fuerte, indicando que casi no hay observaciones.

Probando los residuos con diferentes variables

Vamos a ver si el residuo se distribuye de igual forma entre todas las categorías:

El precio de los diamantes de peor corte són de media más caros

diamonds2 %>%
  ggplot(aes(cut, lresid)) +
    geom_boxplot()

Los diamantes con peor color tienen un precio muy elevado (color J es el peor)

diamonds2 %>%
  ggplot(aes(color, lresid)) +
    geom_boxplot()

La peor claridad parece ser el más caro

diamonds2 %>%
  ggplot(aes(clarity, lresid)) +
    geom_boxplot()

Todos estan más o menos alrededor de 0. Ahora podemos ver que la claridad del diamante si que afecta a precio.

Con el análisi de resiudos queremos explicar la parte del modelo que no queda explicada, y en que escala estan. En caso que sea negativo, podemos concluir que la clase por ejemplo, de claridad, es más barata que la que nos predice el modelo, y a la inversa.

Modelo2: mejora del modelo anterior

Vamos a incluir color, corte y claridad para hacer explícita el efecto de estas tres variables catgóricas en relación al precio del diamante.

mod_diamonds2 <- lm(lprice ~ lcarat + color + cut + clarity,
                    data = diamonds2)

Creación parrilla para las predicciones

grid <- diamonds2 %>%
  data_grid(cut, .model = mod_diamonds2) %>% # parrilla de usando la variable cut a partir del modelo de mod_diamonds2 para hacer más senzilla la injesta de datos
  add_predictions(mod_diamonds2) # le añado la predicción
grid

Visaulización de los resultados para todas las variables predictoras

ggplot(grid, aes(cut, pred)) +
  geom_point() +
  ggtitle("Modelo respecto el corte")

Tenemos un precio que va en función del corte del diamante, a mejor corte, más precio.

diamonds2 %>%
  data_grid(color, .model = mod_diamonds2) %>%
  add_predictions(mod_diamonds2) %>%
  ggplot(aes(color, pred)) +
    geom_point() +
    ggtitle("Modelo respecto al color")

Se ve claramente como baja el precio respecto al color.

diamonds2 %>%
  data_grid(clarity, .model = mod_diamonds2) %>%
  add_predictions(mod_diamonds2) %>%
  ggplot(aes(clarity, pred)) +
    geom_point() +
    ggtitle("Modelo respecto la claridad")

Sube el precio en relación a la claridad.

Evaluando los residuos

diamonds2 <- diamonds2 %>%
  add_residuals(mod_diamonds2, "lresid2")

ggplot(diamonds2, aes(lcarat, lresid2)) +
  geom_hex(bins = 50)

La gran mayoría de valores residules estan en torno a 0, excepto algunos que son muy elevados con una diferencia de 2 puntos. Es decir, que el precio debería ser 4 veces superior al predicho, y a la inversa con los valores negativos.

diamonds2 %>%
  filter(abs(lresid2) >1) %>% # filtramos por los residuos mayores que 1 para evaluar los errores
  add_predictions(mod_diamonds2) %>% # añadimos las predicciones
  mutate(pred = round(2^pred)) %>% # tranformación para ver los precios reales
  select(price, pred, cut, carat, cut, color, clarity, x:z) %>% # seleccionames las variables que nos interesan
  arrange(price) # ordenamos por el precio

Solucionar la pregunta de si es un problema que el modelo no se adapte a este caso. En el caso del diamante “Fair” con una predicción de un precio de 2644, puede que fuera una buena opción comprarlo, dado que el precio esta por debajo de la predicción.

Conclusión y resultado

Con este modelo hemos explicado que no es verdad que los diamantes más malos, más mal cortados y más feos sean los mas baratos, sino que la variable que molesta en la interpretación de los datos es el tamaño de los diamantes, “carat”.

Algunas representaciones gráficas para entender los resultados:

Plot1

ggplot(diamonds2, aes(carat, lresid2)) +
  geom_hex() +
  facet_grid( ~ color)

Plot2

ggplot(diamonds2, aes(carat, lresid2)) +
  geom_hex() +
  facet_grid( ~ cut)

Plot3

ggplot(diamonds2, aes(carat, lresid2)) +
  geom_hex() +
  facet_grid( ~ clarity)