0.1 Nota de los autores:

Para la realización de este proyecto se han seguido las etapas: carga de datos, análisis exploratorio, creación de variables nuevas, selección de variables, construcción de modelos y evaluación de los mismos.

El proceso de evolución entre estas fases durante el proyecto ha sido iterativo, tal y como es, generalmente, en proyectos que involucran ciencia de datos. Esto quiere decir que en determinados momentos del proyecto se ha tenido que volver a atrás desde un punto avanzado, con tal de obtener mejores resultados (e.g volver a realizar análisis exploratorio o plantear una nueva selección de características tras obtener resultados de un modelo).

En esta memoria, sin embargo, se van a reflejar todos los pasos seguidos en cada etapa como si estos hubiesen ocurrido de forma secuencial, con tal de hacerle al lector la lectura más clara. El transcurso real ha sido mucho más complejo.

1 Importación y comprensión del problema

Primero, se importan los ficheros de datos disponibles mediante la función read_delim del paquete readr

dmodel <- read_delim("Modelar_UH2019.txt",delim = "|", escape_double = FALSE, trim_ws = TRUE)
destim <- read_delim("Estimar_UH2019.txt",delim = "|", escape_double = FALSE, trim_ws = TRUE)

Los datos de modelaje tiene un draframe con 53 variables y 9958 filas, siendo una de ellas la variable objetivo, TARGET. Cada fila del dataframe se corresponde con los datos referidos a un inmueble. Para los datos de estimación disponemos de 1104 observaciones.

Se muestra la estructura del dataframe de datos para poder observar las variables de las que se dispone y el tipo de datos que componen cada variable

glimpse(dmodel)
Observations: 9,958
Variables: 53
$ HY_id                           <dbl> 6028499, 6028500, 6028502, 6028512, 6028513, 6028514, 6...
$ HY_cod_postal                   <chr> "18151", "29915", "03740", "30592", "30592", "30592", "...
$ HY_provincia                    <chr> "Granada", "Málaga", "Alicante", "Murcia", "Murcia", "M...
$ HY_descripcion                  <chr> "Plaza de Garaje ubicada en la planta sótano de un edif...
$ HY_distribucion                 <chr> NA, NA, "VIVIENDA-PISO", NA, NA, "Se distribuye en reci...
$ HY_tipo                         <chr> "Garaje", "Casa de pueblo", "Piso", "Garaje", "Garaje",...
$ HY_antiguedad                   <dbl> NA, NA, 2008, NA, NA, NA, NA, 2009, NA, NA, NA, NA, NA,...
$ HY_metros_utiles                <dbl> NA, NA, 92.58, 11.81, 14.29, 48.10, 48.05, 931.17, 11.8...
$ HY_metros_totales               <dbl> 35.00, 151.98, 102.78, 12.46, 15.43, 54.24, 55.90, 956....
$ HY_num_banos                    <dbl> 0, 2, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 2, 2, 2, 1...
$ HY_cert_energ                   <chr> NA, NA, NA, NA, NA, "F", "F", NA, NA, NA, NA, NA, NA, N...
$ HY_num_terrazas                 <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0...
$ HY_ascensor                     <dbl> 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0...
$ HY_trastero                     <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0...
$ HY_num_garajes                  <dbl> 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0...
$ HY_precio                       <dbl> 12300, 115000, 88000, 5800, 7400, 55200, 53700, 1010000...
$ HY_precio_anterior              <dbl> NA, 125000, 93500, NA, NA, 61300, 59600, 1230000, 6200,...
$ IDEA_area                       <dbl> 670.87, 2216.24, 2040.53, 905.44, 905.44, 4116.56, 4116...
$ IDEA_poblacion                  <dbl> 13372, 1507, 5810, 7004, 7004, 1008, 1008, 13767, 844, ...
$ IDEA_densidad                   <dbl> 19.93, 0.68, 2.85, 7.74, 7.74, 0.24, 0.24, 5.18, 13.57,...
$ IDEA_pc_1960                    <dbl> 0.00, 0.12, 0.15, 0.07, 0.07, 0.05, 0.05, 0.03, 0.00, 0...
$ IDEA_pc_1960_69                 <dbl> 0.00, 0.12, 0.15, 0.07, 0.07, 0.05, 0.05, 0.03, 0.00, 0...
$ IDEA_pc_1970_79                 <dbl> 0.00, 0.18, 0.12, 0.04, 0.04, 0.03, 0.03, 0.11, 0.00, 0...
$ IDEA_pc_1980_89                 <dbl> 0.27, 0.18, 0.10, 0.04, 0.04, 0.01, 0.01, 0.13, 0.32, 0...
$ IDEA_pc_1990_99                 <dbl> 0.31, 0.18, 0.11, 0.17, 0.17, 0.03, 0.03, 0.26, 0.00, 0...
$ IDEA_pc_2000_10                 <dbl> 0.41, 0.22, 0.36, 0.62, 0.62, 0.84, 0.84, 0.43, 0.68, 0...
$ IDEA_pc_comercio                <dbl> 0.01, 0.02, 0.03, 0.01, 0.01, 0.00, 0.00, 0.01, 0.00, 0...
$ IDEA_pc_industria               <dbl> 0.00, 0.00, 0.05, 0.01, 0.01, 0.00, 0.00, 0.00, 0.00, 0...
$ IDEA_pc_oficina                 <dbl> 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0...
$ IDEA_pc_otros                   <dbl> 0.10, 0.06, 0.17, 0.06, 0.06, 0.05, 0.05, 0.06, 0.02, 0...
$ IDEA_pc_residencial             <dbl> 0.55, 0.57, 0.59, 0.45, 0.45, 0.52, 0.52, 0.52, 0.12, 0...
$ IDEA_pc_trast_parking           <dbl> 0.34, 0.35, 0.17, 0.47, 0.47, 0.42, 0.42, 0.40, 0.85, 0...
$ IDEA_ind_tienda                 <dbl> 0.20, 0.00, 0.76, 0.51, 0.51, 0.33, 0.33, 0.18, 1.00, 0...
$ IDEA_ind_turismo                <dbl> 0.60, 0.00, 0.08, 0.06, 0.06, 0.00, 0.00, 0.36, 0.00, 1...
$ IDEA_ind_alimentacion           <dbl> 0.20, 0.00, 0.15, 0.43, 0.43, 0.67, 0.67, 0.45, 0.00, 0...
$ IDEA_ind_riqueza                <dbl> 0.20, 0.17, 0.16, 0.31, 0.31, 0.11, 0.11, 0.16, 0.16, 0...
$ IDEA_rent_alquiler              <dbl> 4.21, NA, 6.52, 4.39, 4.39, 5.69, 5.69, 5.78, NA, 4.54,...
$ IDEA_ind_elasticidad            <dbl> 2.0, NA, NA, NA, NA, NA, NA, 4.5, NA, 2.5, 2.5, 2.5, 2....
$ IDEA_ind_liquidez               <dbl> 1, NA, NA, NA, NA, NA, NA, 0, NA, 0, 0, 0, 0, 0, NA, 0,...
$ IDEA_unitprice_sale_residential <dbl> 1043.31, 708.18, 1080.75, 1133.81, 1133.81, 1023.28, 10...
$ IDEA_price_sale_residential     <dbl> 163841.89, 121250.00, 144013.80, 229375.00, 229375.00, ...
$ IDEA_stock_sale_residential     <dbl> 3417, 33, 393, 654, 654, 504, 504, 1086, 21, 552, 552, ...
$ IDEA_demand_sale_residential    <dbl> 48.34, 0.00, 49.51, 49.95, 49.95, 47.87, 47.87, 49.63, ...
$ IDEA_unitprice_rent_residential <dbl> 5.09, 4.08, 4.43, 5.03, 5.03, 5.02, 5.02, 4.99, NA, 4.5...
$ IDEA_price_rent_residential     <dbl> 626.38, 700.00, 525.00, 863.17, 863.17, 433.86, 433.86,...
$ IDEA_stock_rent_residential     <dbl> 237, 6, 15, 168, 168, 24, 24, 93, NA, 27, 27, 27, 27, 2...
$ IDEA_demand_rent_residential    <dbl> 51.01, 50.69, 51.82, 50.74, 50.74, 51.69, 51.69, 51.06,...
$ GA_page_views                   <dbl> 71, 684, 651, 92, 3, 318, 39, 168, 7, 182, 30, 21, 10, ...
$ GA_mean_bounce                  <dbl> 7.14, 16.36, 15.12, 15.82, 0.00, 20.41, 4.55, 0.50, 33....
$ GA_exit_rate                    <dbl> 10.65, 5.25, 11.21, 14.55, 0.00, 11.87, 18.07, 7.26, 50...
$ GA_quincena_ini                 <dbl> 2, 1, 1, 2, 2, 1, 2, 1, 6, 2, 4, 2, 6, 1, 1, 1, 1, 1, 1...
$ GA_quincena_ult                 <dbl> 46, 47, 47, 47, 39, 47, 46, 46, 45, 47, 47, 47, 47, 37,...
$ TARGET                          <dbl> 32.16, 30.71, 63.48, 108.27, 7.50, 41.06, 43.17, 80.40,...

En estos datos se tienen tres grandes grupos de variables. Las que comienzan con HY_ son aquellas relacionadas expresamente con el inmueble, como su antigüedad, el número de baños o el precio. Las que comienzan por IDEA_ son aquellas variables de mercado de la zona en la que está localizado el inmueble, como la densidad de población de la zona, el porcentaje de uso residencial, el indicador de niveles de restaurantes, el precio medio del \(m^2\) en la zona, etc. El último grupo, aquellas que comienzan por GA_ están relacionadas con la visualización del inmueble en la web. Algunos ejemplos son la quincena de aparición del inmueble en la web o hábitos del usuario como el porcentaje de sesiones en las que el usuario no ha interactuado con la página web antes de cerrarla. Finalmente contamos con la variable output TARGET, que nos dice la duración media, en segundos, de la visita web al inmueble.

El objetivo del problema es predecir el tiempo medio de visita en la página web de un inmueble (TARGET) en base a las variables que han sido facilitadas en los ficheros. La métrica para la cual se desea optimizar esta predicción es el error absoluto de la mediana ( median absolute error), con lo cual lo que interesará es que el modelo ajuste bien a los valores centrales y no será necesario que ajuste bien a los valores extremos.

2 Preprocesado de datos

Ahora que se conoce el problema y las características básicas de los datos (tipo y dimensiones de los datos) de los que se dispone, se procede a estudiarlos más a fondo para enriquecerlos generando nuevas variables que se considerán útiles para la resolución del problema.

Cabe señalar que este preprocesado también se realizará con los datos de estimación, para que los modelos predictores que se utilicen sean aplicables a ambos datasets.

2.1 Datos de tipo texto

Los datos de tipo texto libre no son adecuados para los modelos de predicción habituales, ya que no se pueden categorizar a priori. Se extraen, pues, de los campos que incluyen textos largos (HY_descripcion y HY_distribucion) otras variables que sí pueden ser procesadas por el modelo. Se añaden tres por cada campo, una que indique si existe o no, otra que indique el número de carácteres y otra el número de palabras.

HY_distribucion_exist: El campo distribución es o no es faltante HY_distribucion_length: Número de caracteres de la distribución HY_distribucion_words: Número de palabras en la distribución HY_descripcion_exist: El campo descripción es o no es faltante HY_descripcion_length: Número de caracteres de la descripción HY_descripcion_words: Número de palabras en la descripción

En el caso de que falte el campo de descripción o distribución, las variables correspondientes tomarán el valor de 0 en vez del típico NA, ya que la interpretación práctica es la misma tanto si no hay descripción como si esta tiene una longitud de 0 y puede ser más beneficioso para el modelo.

dmodel$HY_distribucion_length<-sapply(dmodel$HY_distribucion, function (x) ifelse(is.na(x),0,nchar(x)))
dmodel$HY_distribucion_exist<-sapply(dmodel$HY_distribucion,function(x)ifelse(is.na(x),0,1))
dmodel$HY_descripcion_length<-sapply(dmodel$HY_descripcion, function (x) ifelse(is.na(x),0,nchar(x)))
dmodel$HY_descripcion_exist<-sapply(dmodel$HY_descripcion,function(x)ifelse(is.na(x),0,1))
dmodel$HY_descripcion_words<-sapply(dmodel$HY_descripcion,function(x) ifelse(is.na(x),0,length(unlist(strsplit(x," ")))))
dmodel$HY_distribucion_words<-sapply(dmodel$HY_distribucion,function(x) ifelse(is.na(x),0,length(unlist(strsplit(x," ")))))

destim$HY_distribucion_length<-sapply(destim$HY_distribucion, function (x) ifelse(is.na(x),0,nchar(x)))
destim$HY_distribucion_exist<-sapply(destim$HY_distribucion,function(x)ifelse(is.na(x),0,1))
destim$HY_descripcion_length<-sapply(destim$HY_descripcion, function (x) ifelse(is.na(x),0,nchar(x)))
destim$HY_descripcion_exist<-sapply(destim$HY_descripcion,function(x)ifelse(is.na(x),0,1))
destim$HY_descripcion_words<-sapply(destim$HY_descripcion,function(x) ifelse(is.na(x),0,length(unlist(strsplit(x," ")))))
destim$HY_distribucion_words<-sapply(destim$HY_distribucion,function(x) ifelse(is.na(x),0,length(unlist(strsplit(x," ")))))

2.2 Datos de imágenes

Se dispone también de imágenes de los inmuebles (en concreto 47382 imágenes asociadas a los inmuebles de los datasets de modelado y de estimación), otro tipo de datos que tampoco permiten el uso de modelos de regresión sencillos. Se decidió, al igual que con los campos de texto, transformar estas imágenes en otras variables que pueden ser introducidas al modelo. La primera variable que se va a crear será el número de imágenes que hay para cada inmueble, y también se añadirá el valor medio de la intensidad de los píxeles de las imágenes de cada inmueble (lo cual puede dar una idea de la luminosidad de las imágenes). También se añadirá la media de la primera imagen y el número de dimensiones de la misma, ya que se piensa que las primeras impresiones son las más importantes. Finalmente, se añadirá una variable que nos dirá el grado de similitud media entre las imágenes de cada inmueble y otra que nos indique cuántas imágenes de estas superan un umbral de similitud, por lo que serían consideradas la misma foto. En resumen, creamos:

IM_n: Número de imágenes del inmueble IM_mean: Media de todas las imágenes del inmueble IM_mean_1: Media de la primera imagen del inmueble IM_size_1: Número de píxeles de la primera imagen del inmueble IM_sim_mean: Similitud media de todas las imágenes del inmueble. NA si sólo hay una imagen. IM_n_rep: Número de imágenes repetidas. 0 si sólo hay una imagen (no puede estar repetida)

#Se establece la ruta al directorio donde están las imágenes que se van a usar.
path = "./imagenes_inmuebles_haya/"

#Se obtiene un vector con la ruta de cada imagen que se va a emplear
files <- dir(path, full.names = TRUE)
files<-files[!files=="./imagenes_inmuebles_haya/6031545__posifoto3__1K9lJvlehV.jpg"] #ESTA SE ELIMINA PORQUE ESTÁ CORRUPTA

#Se crea dataframe auxiliar, la primera columna contiene la ID de un inmueble a y la segunda columna contiene el número de imágenes en la página de dicho inmueble.
images_por_inmueble<-data.frame(HY_id=str_extract(files,"[0-9]{7}"))%>%group_by(HY_id)%>%summarise(IM_n=n())


#Se crea dataframe auxiliar, la primera columna contiene la ID de una vivienda y la segunda columna contiene la media de sus imágenes.
media_imagenes<-data.frame(HY_id=str_extract(files,"[0-9]{7}"),IM_mean=sapply(files,function(x) mean(readJPEG(x))))%>%group_by(HY_id)%>%summarise(IM_mean=mean(media))

#Vamos a obtener la media y las dimensiones de la primera foto
dimensiones<-function(x){
  #INPUT: Paths,
     if(str_detect(x,"posifoto1")==T){ #Si el path inclye posifoto1 será la primera foto
         a<-readJPEG(x)
         return(dim(a)[1]*dim(a)[2])
     }else{0}}
imagenes_1<-data.frame(HY_id=str_extract(files,"[0-9]{7}"),
                        Num_foto=str_extract(files,"posifoto[:digit:]{1,}"),
                        IM_mean_1=sapply(files,function(x){ifelse(str_detect(x,"posifoto1")==T,
                                                              mean(readJPEG(x)),
                                                              0)}),
                        IM_size_1=sapply(files,function(x){dimensiones(x)}))

#Creamos una función que nos devuelve la similitud entre las fotos
calculadoraSimilitudes<-function(paths){
  #INPUT: Paths a las imágenes de cad inmueble
  #OUTPUT: Un vector con las similitudes entre las fotos
  #LEEMOS LAS IMAGENES Y LAS REDIMENSIONAMOS
  ims<-lapply(paths,function(x) resize(readJPEG(x),w=350,h=300))
  #COMPROBAMOS QUE TENGAN TODAS LAS DIMENSIONES y si no añadimos nosotros
  for (im in 1:length(ims)){
    if(length(ims[[im]])<315000){
      ims[[im]]=abind(ims[[im]],ims[[im]],ims[[im]],along=3)
    }
  }
  #CALCULAMOS LAS SIMILITUDES
  sim=c()
  if(length(ims)>1){
  for(i in 1:(length(ims)-1)){
    for (j in (i+1):length(ims)){
      sim=c(sim,similarity(ims[[i]],ims[[j]]))
    }
  }
  }else{sim=NA} #SI NO HAY MÁS DE UNA IMAGEN, NO HAY NADA A LO QUE PARECERSE Y PONEMOS NA
  return(sim)
}
#Creamos dataframe auxiliar
imagenes<-data.frame(paths=files,HY_id=str_extract(files,"[0-9]{7}"),IM_sim_mean=rep(NA,length(files)),IM_n_rep=rep(NA,length(files)))
for (id in unique(imagenes$HY_id)){
  sim<-calculadoraSimilitudes(as.character(unlist(imagenes%>%filter(HY_id==id)%>%select(paths))))
  imagenes[imagenes$HY_id==id,"IM_sim_mean"]=mean(sim)
  imagenes[imagenes$HY_id==id,"IM_n_rep"]=ifelse(sum(is.na(sim)),0,sum(sim>1.2))
}
imagenes=imagenes%>%select(-paths)
#INCORPORAMOS TODAS LAS VARIABLES A LOS DATASETS
images_por_inmueble$HY_id<-as.numeric(as.character(images_por_inmueble$HY_id))
#Añadimos al dataset original el número de imágenes por inmueble
dmodel<-left_join(dmodel,images_por_inmueble,by="HY_id")
dmodel$IM_n[is.na(dmodel$IM_n)]=0 #Si no hay disponibles imágenes, ponemos que hay 0
destim<-left_join(destim,images_por_inmueble,by="HY_id")
destim$IM_n[is.na(destim$IM_n)]=0#Si no hay disponibles imágenes, ponemos que hay 0


media_imagenes$HY_id<-as.numeric(as.character(media_imagenes$HY_id))
#Añadimos al dataset original el valor medio de las imágenes por inmueble
dmodel<-left_join(dmodel,media_imagenes,by="HY_id")
#Añadimos al dataset estim estas  variables
destim<-left_join(destim,media_imagenes,by="HY_id")


imagenes_1_filt<-filter(imagenes_1,Num_foto=="posifoto1")%>%select(-Num_foto)
imagenes_1_filt$HY_id<-as.numeric(as.character(imagenes_1_filt$HY_id))
dmodel<-left_join(dmodel,imagenes_1_filt,by="HY_id")
destim<-left_join(destim,imagenes_1_filt,by="HY_id")

imagenes$HY_id<-as.numeric(as.character(imagenes$HY_id))
imagenes=unique(imagenes)
dmodel<-left_join(dmodel,imagenes,by="HY_id")
dmodel$IM_n_rep[is.na(dmodel$IM_n_rep)]=0 #Si hay 0 imágenes, hay 0 imágenes repetidas
destim<-left_join(destim,imagenes,by="HY_id")
destim$IM_n_rep[is.na(destim$IM_n_rep)]=0 #Si hay 0 imágenes, hay 0 imágenes repetidas


#Eliminamos la variables auxiliares
rm(images_por_inmueble)
rm(media_imagenes)
rm(imagenes_1_filt,imagenes_1)
rm(imagenes)
rm(dimensiones)
rm(files,path)

2.3 Datos faltantes

El siguiente paso es estudiar los datos faltantes que hay en el dataset dmodel.

miss_var_summary(dmodel)
# A tibble: 65 x 3
   variable                        n_miss pct_miss
   <chr>                            <int>    <dbl>
 1 HY_cert_energ                     8300   83.4  
 2 HY_distribucion                   6745   67.7  
 3 HY_antiguedad                     5484   55.1  
 4 IDEA_ind_elasticidad              5106   51.3  
 5 IDEA_ind_liquidez                 5106   51.3  
 6 HY_descripcion                    3942   39.6  
 7 HY_precio_anterior                3158   31.7  
 8 IDEA_rent_alquiler                3062   30.7  
 9 IDEA_unitprice_rent_residential   2998   30.1  
10 IDEA_price_rent_residential       2998   30.1  
11 IDEA_stock_rent_residential       2998   30.1  
12 HY_metros_utiles                  2886   29.0  
13 IM_sim_mean                       2727   27.4  
14 IDEA_ind_tienda                   2717   27.3  
15 IDEA_ind_turismo                  2717   27.3  
16 IDEA_ind_alimentacion             2717   27.3  
17 IDEA_pc_1960                      2710   27.2  
18 IDEA_pc_1960_69                   2710   27.2  
19 IDEA_pc_1970_79                   2710   27.2  
20 IDEA_pc_1980_89                   2710   27.2  
21 IDEA_pc_1990_99                   2710   27.2  
22 IDEA_pc_2000_10                   2710   27.2  
23 IDEA_pc_comercio                  2710   27.2  
24 IDEA_pc_industria                 2710   27.2  
25 IDEA_pc_oficina                   2710   27.2  
26 IDEA_pc_otros                     2710   27.2  
27 IDEA_pc_residencial               2710   27.2  
28 IDEA_pc_trast_parking             2710   27.2  
29 IDEA_poblacion                    2635   26.5  
30 IDEA_densidad                     2635   26.5  
31 IDEA_unitprice_sale_residential   2635   26.5  
32 IDEA_price_sale_residential       2635   26.5  
33 IDEA_stock_sale_residential       2635   26.5  
34 IDEA_demand_sale_residential      2635   26.5  
35 IDEA_demand_rent_residential      2635   26.5  
36 IDEA_area                         2631   26.4  
37 IDEA_ind_riqueza                  2631   26.4  
38 IM_mean_1                         2121   21.3  
39 IM_size_1                         2121   21.3  
40 IM_mean                           2116   21.2  
41 HY_provincia                        42    0.422
42 HY_metros_totales                   34    0.341
43 HY_id                                0    0    
44 HY_cod_postal                        0    0    
45 HY_tipo                              0    0    
46 HY_num_banos                         0    0    
47 HY_num_terrazas                      0    0    
48 HY_ascensor                          0    0    
49 HY_trastero                          0    0    
50 HY_num_garajes                       0    0    
51 HY_precio                            0    0    
52 GA_page_views                        0    0    
53 GA_mean_bounce                       0    0    
54 GA_exit_rate                         0    0    
55 GA_quincena_ini                      0    0    
56 GA_quincena_ult                      0    0    
57 TARGET                               0    0    
58 HY_distribucion_length               0    0    
59 HY_distribucion_exist                0    0    
60 HY_descripcion_length                0    0    
61 HY_descripcion_exist                 0    0    
62 HY_descripcion_words                 0    0    
63 HY_distribucion_words                0    0    
64 IM_n                                 0    0    
65 IM_n_rep                             0    0    
vis_miss(dmodel)

Con este gráfico se aprecia a simple vista como las variables IDEA_ van a resultar conflictivas, ya que tienen un elevado número de valores faltantes. Además, hay otras variables del grupo HY_ que también lo van a ser, como el certificado energético, que alcanza el 83% de datos faltantes, la distribución y la antigüedad. Por otro lado, se observa cómo hay ciertas variables sin datos faltantes como aquellas que se corresponden con información sobre la página web y la mayoría del grupo HY_.

Además, se puede ver esta misma distribución de valores faltantes en destim:

vis_miss(destim)

con lo cual, las transformaciones que se apliquen a dmodel serán también aplicadas a destim.

Se sustituye la variable HY_cer_energ, que es la que más datos faltantes tiene, por una variable dummy que indique simplemente si está presente o no en cada inmuble.

*HY_cert_energ_exist: Existe o no el certificado energético

dmodel$HY_cert_energ_exist=sapply(dmodel$HY_cert_energ,function(x) ifelse(is.na(x),0,1))
destim$HY_cert_energ_exist=sapply(destim$HY_cert_energ,function(x) ifelse(is.na(x),0,1))

Para intentar enriquecer los datos de las variables de tipo IDEA_, se realiza un promedio de los valores de estas variables por códigos postales, ya que en principio parece lógico que las características zonales vayan asociadas con el código postal. Luego se imputan los valores faltantes dependiendo del código postal. Si no se dispone información de un código postal en concreto, se mantienen los valores faltantes.

#Se juntan los dataframes de modelización y estimaciones para obtener todos los valores posibles de las variables de IDEA.
dmerge<-rbind(dmodel%>%select(-TARGET),destim)

a<-unlist(dmerge%>%filter(is.na(IDEA_rent_alquiler))%>%group_by(HY_cod_postal)%>%summarise())

b<-unlist(dmerge%>%filter(!is.na(IDEA_rent_alquiler))%>%group_by(HY_cod_postal)%>%summarise())

#Se ve cuántos códigos postales que tienen filas con valores NA tienen también filas sin ellos para poder obtener la información económinca de la zona.
c<-a[a%in%b]

#Se obtiene la media de los indicadores económicos de un mismo código postal, valor que se considerará representativo de dicho código postal y se empleará para sustituirlo en aquella filas con datos faltantes.
dmerge_mean<-dmerge%>%filter(HY_cod_postal%in%c)%>%
  select(HY_cod_postal,starts_with("IDEA"))%>%
  group_by(HY_cod_postal)%>%
  summarise_all(funs(mean),na.rm=T)

#Se incluyen en los dataframes originales los valores obtenidos para cada código postal para aquellas filas en las que no hay datos.
for (i in 1:nrow(dmodel)) {
  if (dmodel$HY_cod_postal[i]%in%dmerge_mean$HY_cod_postal) {
    dmodel[i,18:47]<-dmerge_mean[dmerge_mean$HY_cod_postal==dmodel$HY_cod_postal[i],2:31]
  }
}

for (i in 1:nrow(destim)) {
  if (destim$HY_cod_postal[i]%in%dmerge_mean$HY_cod_postal) {
    destim[i,18:47]<-dmerge_mean[dmerge_mean$HY_cod_postal==destim$HY_cod_postal[i],2:31]
  }
}

#Aquellos datos con un valor NaN debido a la operación de la media realizada antes se convierten de nuevo a tipo NA.
for (k in 18:47) {
  dmodel[is.nan(unlist(dmodel[,k])),k]<-NA
}

for (k in 18:4) {
  destim[is.nan(unlist(destim[,k])),k]<-NA
}

rm(a,b,c,i,k,dmerge,dmerge_mean)

3 Análisis exploratorio de datos

A continuación llega una de las partes más complicadas a la hora de generar un modelo: el análisis exploratorio de datos. Se desea estudiar estadísticamente la distribución de los valores de cada una de las variables para ver posibles valores faltantes, valores atípicos, relaciones entre variables…

Se obtiene un gráfico boxplot de aquellas variables numéricas para comprobar la presencia de outliers, pero para simplificar este report se omiten. En ellos se puede apreciar como casi todas las variables muestran la presencia de outliers. Se va a estudiar, por ejemplo, el caso del número de baños de los inmuebles. Tenemos inmuebles con las siguientes cantidades de baños:

unique(dmodel$HY_num_banos)
[1]  0  2  1  3  4  6 99  5 20

Se ve cómo hay casas en las que hay 20 o 99 baños, esto da una idea de que al introducirse a mano se introducen errores, por lo que habrá que filtrar los datos disponibles.

Hay un inmueble con una antigüedad de 2088 en dmodel, por lo que se procede a sustituirlo por un NA.

dmodel$HY_antiguedad[dmodel$HY_antiguedad==2088]=NA

Se puede ver cómo afecta la cantidad de carácteres de la descripción a la duración media de las visitas a los inmuebles.

dmodel%>%filter(HY_descripcion_exist==1)%>%
  ggplot(aes(x=HY_descripcion_length,y=TARGET))+
  stat_density_2d(aes(fill = ..level..), geom = "polygon", colour="darkgreen")+
  scale_fill_distiller(palette=11, direction=1)

Se puede ver claramente que conforme más corta es la descripción, menos dura la visita, en la mayoría de los casos.

Obteniendo los estadísticos de las variables numéricas presentes en el dataset se puede estudiar también su distribución.

datos=unlist(apply(dmodel[,sapply(dmodel, is.numeric)],2,function(x) list(mean(x,na.rm=T),median(x,na.rm=T),min(x,na.rm=T),max(x,na.rm=T),sd(x,na.rm=T),IQR(x,na.rm=T),quantile(x,p=0.05,names=F,na.rm=T),quantile(x,p=0.25,names=F,na.rm=T),quantile(x,p=0.75,names=F,na.rm=T),quantile(x,p=0.95,names=F,na.rm=T))))
metrics=as.data.frame(matrix(datos,10,ncol(dmodel[,sapply(dmodel, is.numeric)]))
                   ,stringsAsFactors = F
                   ,row.names=c('Media','Mediana','Mínimo','Máximo','Desviación típica','Intervalo intercuartil','Percentil 5','Percentil 25','Percentil 75','Percentil 95'))
names(metrics)<-names(dmodel[,sapply(dmodel, is.numeric)])
metrics=round(metrics,3)

metrics
                           HY_id HY_antiguedad HY_metros_utiles HY_metros_totales HY_num_banos
Media                  6272841.7      1997.286          471.424           543.392        0.820
Mediana                6080821.0      2007.000           63.000            81.000        1.000
Mínimo                 6028499.0      1877.000            0.000             0.000        0.000
Máximo                 7017665.0      2016.000      1820000.000       1820000.000       99.000
Desviación típica       342104.2        21.049        21920.308         19219.100        2.913
Intervalo intercuartil  489208.2        17.000           75.715            83.280        1.000
Percentil 5            6030907.8      1960.000            3.750            14.000        0.000
Percentil 25           6039047.5      1992.000           14.285            31.700        0.000
Percentil 75           6528255.8      2009.000           90.000           114.980        1.000
Percentil 95           6950459.6      2011.000          243.455           450.000        2.000
                       HY_num_terrazas HY_ascensor HY_trastero HY_num_garajes  HY_precio
Media                            0.020       0.270       0.062          0.097    76724.3
Mediana                          0.000       0.000       0.000          0.000    59700.0
Mínimo                           0.000       0.000       0.000          0.000        0.0
Máximo                          10.000       1.000       1.000          2.000 10400000.0
Desviación típica                0.185       0.444       0.241          0.297   163436.1
Intervalo intercuartil           0.000       1.000       0.000          0.000    81200.0
Percentil 5                      0.000       0.000       0.000          0.000     4600.0
Percentil 25                     0.000       0.000       0.000          0.000    10500.0
Percentil 75                     0.000       1.000       0.000          0.000    91700.0
Percentil 95                     0.000       1.000       1.000          1.000   201000.0
                       HY_precio_anterior IDEA_area IDEA_poblacion IDEA_densidad IDEA_pc_1960
Media                            89392.54  5233.819      10150.689        33.752        0.110
Mediana                          64100.00  1599.820       8229.000         7.740        0.104
Mínimo                               1.00     6.840          1.000         0.010        0.000
Máximo                        30833100.00 95349.800      54589.000       526.030        0.480
Desviación típica               431414.61 10395.859       7726.604        58.544        0.083
Intervalo intercuartil           84800.00  4329.856       9145.357        44.204        0.117
Percentil 5                       5100.00    90.532       1234.000         0.300        0.000
Percentil 25                     12300.00   577.134       4518.500         1.910        0.043
Percentil 75                     97100.00  4906.990      13663.857        46.114        0.160
Percentil 95                    227000.00 23921.585      25869.732       152.063        0.270
                       IDEA_pc_1960_69 IDEA_pc_1970_79 IDEA_pc_1980_89 IDEA_pc_1990_99
Media                            0.110           0.158           0.136           0.151
Mediana                          0.104           0.143           0.130           0.140
Mínimo                           0.000           0.000           0.000           0.000
Máximo                           0.480           0.990           0.750           0.800
Desviación típica                0.083           0.096           0.062           0.064
Intervalo intercuartil           0.117           0.112           0.062           0.078
Percentil 5                      0.000           0.030           0.054           0.070
Percentil 25                     0.043           0.090           0.098           0.110
Percentil 75                     0.160           0.202           0.160           0.188
Percentil 95                     0.270           0.329           0.223           0.260
                       IDEA_pc_2000_10 IDEA_pc_comercio IDEA_pc_industria IDEA_pc_oficina
Media                            0.334            0.022             0.008           0.007
Mediana                          0.300            0.020             0.004           0.007
Mínimo                           0.000            0.000             0.000           0.000
Máximo                           1.000            0.170             0.230           0.120
Desviación típica                0.160            0.014             0.012           0.008
Intervalo intercuartil           0.225            0.018             0.010           0.010
Percentil 5                      0.117            0.002             0.000           0.000
Percentil 25                     0.214            0.010             0.000           0.000
Percentil 75                     0.439            0.028             0.010           0.010
Percentil 95                     0.633            0.049             0.029           0.020
                       IDEA_pc_otros IDEA_pc_residencial IDEA_pc_trast_parking IDEA_ind_tienda
Media                          0.054               0.500                 0.408           0.328
Mediana                        0.040               0.500                 0.418           0.310
Mínimo                         0.000               0.120                 0.020           0.000
Máximo                         0.430               0.910                 0.680           1.000
Desviación típica              0.042               0.070                 0.076           0.240
Intervalo intercuartil         0.044               0.098                 0.090           0.360
Percentil 5                    0.015               0.396                 0.270           0.000
Percentil 25                   0.024               0.452                 0.368           0.140
Percentil 75                   0.068               0.550                 0.458           0.500
Percentil 95                   0.146               0.610                 0.520           0.711
                       IDEA_ind_turismo IDEA_ind_alimentacion IDEA_ind_riqueza IDEA_rent_alquiler
Media                             0.325                 0.261            0.163              5.382
Mediana                           0.242                 0.240            0.145              5.210
Mínimo                            0.000                 0.000            0.030              2.500
Máximo                            1.000                 1.000            1.000             12.000
Desviación típica                 0.290                 0.212            0.066              1.179
Intervalo intercuartil            0.439                 0.300            0.060              1.276
Percentil 5                       0.000                 0.000            0.100              3.920
Percentil 25                      0.101                 0.100            0.121              4.629
Percentil 75                      0.540                 0.400            0.181              5.905
Percentil 95                      1.000                 0.669            0.290              7.230
                       IDEA_ind_elasticidad IDEA_ind_liquidez IDEA_unitprice_sale_residential
Media                                 3.087             0.018                        1007.614
Mediana                               3.000             0.000                         895.269
Mínimo                                1.000            -1.000                         273.830
Máximo                                5.000             1.000                        6363.400
Desviación típica                     1.053             0.424                         400.447
Intervalo intercuartil                1.734             0.000                         359.440
Percentil 5                           1.103            -0.950                         592.610
Percentil 25                          2.206             0.000                         773.380
Percentil 75                          3.940             0.000                        1132.820
Percentil 95                          5.000             1.000                        1740.417
                       IDEA_price_sale_residential IDEA_stock_sale_residential
Media                                    143922.74                    1122.042
Mediana                                  124710.10                    1005.000
Mínimo                                    28575.00                       3.000
Máximo                                  1360885.71                    4539.000
Desviación típica                         73898.65                     777.097
Intervalo intercuartil                    51321.25                     944.484
Percentil 5                               82937.33                     166.250
Percentil 25                             104383.69                     540.069
Percentil 75                             155704.94                    1484.553
Percentil 95                             280075.38                    2388.750
                       IDEA_demand_sale_residential IDEA_unitprice_rent_residential
Media                                        48.836                           5.098
Mediana                                      49.133                           4.920
Mínimo                                        0.000                           1.000
Máximo                                       52.380                          18.750
Desviación típica                             3.725                           1.628
Intervalo intercuartil                        0.878                           1.751
Percentil 5                                  47.580                           2.970
Percentil 25                                 48.712                           4.111
Percentil 75                                 49.590                           5.862
Percentil 95                                 50.370                           8.120
                       IDEA_price_rent_residential IDEA_stock_rent_residential
Media                                      592.797                     112.633
Mediana                                    530.280                      63.375
Mínimo                                     200.000                       3.000
Máximo                                    3023.800                    5685.000
Desviación típica                          263.841                     219.663
Intervalo intercuartil                     225.860                      99.000
Percentil 5                                350.000                       6.000
Percentil 25                               447.635                      30.000
Percentil 75                               673.495                     129.000
Percentil 95                               997.500                     339.000
                       IDEA_demand_rent_residential GA_page_views GA_mean_bounce GA_exit_rate
Media                                        48.238       190.965         18.716       22.210
Mediana                                      51.350        30.000         14.290       16.960
Mínimo                                        0.000         0.000          0.000        0.000
Máximo                                       54.220     13698.000        100.000      100.000
Desviación típica                            11.641       482.313         19.984       20.995
Intervalo intercuartil                        1.089       162.000         29.170       21.388
Percentil 5                                   5.418         1.000          0.000        0.000
Percentil 25                                 50.700         6.000          0.000        8.612
Percentil 75                                 51.789       168.000         29.170       30.000
Percentil 95                                 52.460       902.150         55.267       65.509
                       GA_quincena_ini GA_quincena_ult   TARGET HY_distribucion_length
Media                           14.421          42.547   74.156                 43.038
Mediana                          9.000          46.000   56.275                  0.000
Mínimo                           1.000           1.000    0.000                  0.000
Máximo                          47.000          47.000 1770.000                500.000
Desviación típica               14.091           8.404   89.739                 75.924
Intervalo intercuartil          24.000           5.000   61.853                 83.000
Percentil 5                      1.000          22.000    4.500                  0.000
Percentil 25                     2.000          42.000   28.655                  0.000
Percentil 75                    26.000          47.000   90.508                 83.000
Percentil 95                    41.000          47.000  194.179                207.150
                       HY_distribucion_exist HY_descripcion_length HY_descripcion_exist
Media                                  0.323               431.821                0.604
Mediana                                0.000               172.000                1.000
Mínimo                                 0.000                 0.000                0.000
Máximo                                 1.000              3000.000                1.000
Desviación típica                      0.468               562.364                0.489
Intervalo intercuartil                 1.000               715.000                1.000
Percentil 5                            0.000                 0.000                0.000
Percentil 25                           0.000                 0.000                0.000
Percentil 75                           1.000               715.000                1.000
Percentil 95                           1.000              1640.000                1.000
                       HY_descripcion_words HY_distribucion_words   IM_n IM_mean IM_mean_1
Media                                67.338                 6.637  4.282   0.479     0.483
Mediana                              30.000                 0.000  4.000   0.476     0.480
Mínimo                                0.000                 0.000  0.000   0.161     0.118
Máximo                              505.000                95.000 10.000   0.920     0.995
Desviación típica                    86.901                12.028  3.450   0.063     0.110
Intervalo intercuartil              112.000                12.000  6.000   0.076     0.118
Percentil 5                           0.000                 0.000  0.000   0.385     0.313
Percentil 25                          0.000                 0.000  1.000   0.441     0.421
Percentil 75                        112.000                12.000  7.000   0.517     0.539
Percentil 95                        249.150                32.000 10.000   0.581     0.636
                       IM_size_1 IM_sim_mean IM_n_rep HY_cert_energ_exist
Media                   299548.4       1.035    0.145               0.166
Mediana                 220820.0       1.027    0.000               0.000
Mínimo                   15000.0       1.001    0.000               0.000
Máximo                 2073600.0       1.509   12.000               1.000
Desviación típica       300038.5       0.044    0.554               0.373
Intervalo intercuartil  174594.0       0.014    0.000               0.000
Percentil 5              65852.0       1.012    0.000               0.000
Percentil 25            148950.0       1.021    0.000               0.000
Percentil 75            323544.0       1.035    0.000               0.000
Percentil 95            874800.0       1.072    1.000               1.000

En el caso de TARGET, la variable de salida, se observa cómo hay una desviación entre el valor de la media y la mediana de la variable target lo que da una idea de la distribución de los valores. Se procede a representar la variable TARGET en forma de histograma.

3.1 Estudio de la variable TARGET

ggplot(data=dmodel,aes(TARGET))+
  geom_histogram(color="red",fill='darkred',bins=40)+
  ylab("segundos")+
  ggtitle("Duración media de la visita a una página de un inmueble")

Se ve que la mayor parte de los valores se encuentran en el rango de 0 a 250 segundos de media, por lo que se va a representar de nuevo el histograma centrándose en esa zona.

ggplot(data=dmodel%>%filter(TARGET<250),aes(TARGET))+
  geom_histogram(color="red",fill='darkred')+
  ylab("segundos")+
  ggtitle("Duración media de la visita a una página de un inmueble")

Se ve cómo la distibución está centrada en valores bajos, de alrededor de un minuto, pero tiene una cola hacia la derecha que hace que el valor de la media se desplace respecto a la mediana y sea mayor.

3.2 Estudio de la variable TARGET respecto a los tipos de inmueble

Se puede ver también los distintos tipos de inmueble, sus frecuencias y sus tiempos medios de visitas:

dmodel%>%group_by(HY_tipo)%>%summarise(Cantidad=n(), 'Cantidad relativa' = n()/9958, `Tiempo medio de visita (s)`=mean(TARGET))%>%arrange(desc(Cantidad))
# A tibble: 20 x 4
   HY_tipo                      Cantidad `Cantidad relativa` `Tiempo medio de visita (s)`
   <chr>                           <int>               <dbl>                        <dbl>
 1 Piso                             4182            0.420                            73.9
 2 Garaje                           2761            0.277                            63.3
 3 Local                             827            0.0830                           85.6
 4 Casa de pueblo                    557            0.0559                           75.6
 5 Chalet adosado                    433            0.0435                           85.1
 6 Trastero                          282            0.0283                           70.2
 7 Chalet independiente              209            0.0210                           91.8
 8 Solar                             198            0.0199                           80.7
 9 Dúplex                            184            0.0185                           83.6
10 Nave adosada                      130            0.0131                          112. 
11 Oficina                           103            0.0103                          110. 
12 Nave aislada                       38            0.00382                         108. 
13 Edificio turístico                 17            0.00171                          95.6
14 Otro                               12            0.00121                          64.6
15 Chalet pareado                      9            0.000904                         87.1
16 Ático                               8            0.000803                         75.3
17 Almacén                             3            0.000301                         71.7
18 Edificio comercial                  2            0.000201                         99.6
19 Suelo rústico no urbanizable        2            0.000201                         89.2
20 Aparcamiento                        1            0.000100                         24.8

La mayoría son pisos o garajes.

Si se quitasen aquellos registros con un TARGET superior a 300 (considerado outlier), se perderán sólo el 1.9783089% de los datos, por lo que se continuará teniendo una distribución de tipos de inmueble similar.

dmodel%>%filter(TARGET<300)%>%group_by(HY_tipo)%>%summarise(Cantidad=n(), 'Cantidad relativa' = n()/nrow(dmodel%>%filter(TARGET<300)), `Tiempo medio de visita (s)`=mean(TARGET),`Desviación típica del tiempo (s)`=sd(TARGET))%>%arrange(desc(Cantidad))
# A tibble: 20 x 5
   HY_tipo             Cantidad `Cantidad relativ~ `Tiempo medio de visi~ `Desviación típica del ti~
   <chr>                  <int>              <dbl>                  <dbl>                      <dbl>
 1 Piso                    4119           0.422                      67.2                       48.8
 2 Garaje                  2691           0.276                      50.2                       55.1
 3 Local                    813           0.0833                     79.0                       45.4
 4 Casa de pueblo           551           0.0564                     72.3                       46.7
 5 Chalet adosado           426           0.0436                     78.6                       55.7
 6 Trastero                 273           0.0280                     54.6                       50.6
 7 Chalet independien~      208           0.0213                     90.7                       50.1
 8 Solar                    195           0.0200                     77.0                       53.3
 9 Dúplex                   180           0.0184                     76.2                       48.1
10 Nave adosada             122           0.0125                     90.2                       54.5
11 Oficina                   94           0.00963                    59.6                       54.0
12 Nave aislada              36           0.00369                    94.7                       54.3
13 Edificio turístico        17           0.00174                    95.6                       58.9
14 Otro                      12           0.00123                    64.6                       80.9
15 Ático                      8           0.000820                   75.3                       35.5
16 Chalet pareado             8           0.000820                   47.7                       53.9
17 Almacén                    3           0.000307                   71.7                       24.7
18 Edificio comercial         2           0.000205                   99.6                       30.0
19 Suelo rústico no u~        2           0.000205                   89.2                       43.7
20 Aparcamiento               1           0.000102                   24.8                      NaN  

Se puede ver que más de dos tercios de los inmuebles son pisos o garajes y tienen un tiempo medio de visita de alrededor de un minuto. Los inmuebles cuyas visitas duran más son las naves y las oficinas, aunque no representen ni un 2% del dataset.

Se va a observar también a continuación la distribución de los valores de la variable TARGET en función del tipo de inmueble que sea, esto se va a realizar mediante la representación de varios boxplots

ggplot(dmodel%>%filter(TARGET<250)%>%group_by(HY_tipo))+
  geom_boxplot(aes(x=HY_tipo,y=TARGET,col=HY_tipo))+theme(legend.position="none",axis.text.x = element_text(angle = 45, hjust = 1))

Intuyendo agrupaciones en los tipos de inmueble, creamos HY_tipo_simplif, con la que distinguimos tres grupos de inmuebles. HY_tipo_simplif2 tendrá 5 subcategorías que corresponderán a las que se realizan en la página web de Haya.

HY_tipo_simplif: crea 3 grupos distintos dependiendo del tipo de inmueble. En el subgrupo de Viviendas quedarán los inmuebles Piso, Ático, Chalet pareado, Dúplex, Chalet adosado, Casa de pueblo, Chalet independiente; en el subgrupo de Garajes habrá Aparcamiento, Trastero y Garaje; finalmente, en Resto quedarán todos los inmuebles que no pertenezcan a los tipos anteriormente mencionados. —HY_tipo_simplif2: tendrá 5 subcategorías que se corresponden a las que se realizan en la página web de Haya. En el subgrupo de Viviendas quedan los inmuebles Piso, Ático, Chalet pareado, Dúplex, Chalet adosado, Casa de pueblo, Chalet independiente; en el subgrupo de Garajes estarán Aparcamiento, Trastero, Garaje y Otro; en Locales quedan Local,Oficina,Almacén,Edificio turístico, Edificio comercial; en la categoría Naves se incluyen Nave adosada,Nave aislada; finalmente en el grupo Terrenos se recogen los inmuebles de tipo Solar,Suelo rústico no urbanizable.

dmodel$HY_tipo_simplif<-rep(0,nrow(dmodel))
for (i in 1:nrow(dmodel)) {
  if(dmodel$HY_tipo[i]%in%c("Piso", "Ático", "Chalet pareado", "Dúplex", "Chalet adosado", "Casa de pueblo", "Chalet independiente")){dmodel$HY_tipo_simplif[i]<-"Vivienda"}
  if(dmodel$HY_tipo[i]%in%c("Aparcamiento","Trastero","Garaje")){dmodel$HY_tipo_simplif[i]<-"Garajes"}
  if(!dmodel$HY_tipo[i]%in%c("Piso", "Ático", "Chalet pareado", "Dúplex", "Chalet adosado", "Casa de pueblo", "Chalet independiente","Aparcamiento","Trastero","Garaje")){dmodel$HY_tipo_simplif[i]<-"Resto"}
}
destim$HY_tipo_simplif<-rep(0,nrow(destim))
for (i in 1:nrow(destim)) {
  if(destim$HY_tipo[i]%in%c("Piso", "Ático", "Chalet pareado", "Dúplex", "Chalet adosado", "Casa de pueblo", "Chalet independiente")){destim$HY_tipo_simplif[i]<-"Vivienda"}
  if(destim$HY_tipo[i]%in%c("Aparcamiento","Trastero","Garaje")){destim$HY_tipo_simplif[i]<-"Garajes"}
  if(!destim$HY_tipo[i]%in%c("Piso", "Ático", "Chalet pareado", "Dúplex", "Chalet adosado", "Casa de pueblo", "Chalet independiente","Aparcamiento","Trastero","Garaje")){destim$HY_tipo_simplif[i]<-"Resto"}
}

dmodel$HY_tipo_simplif2<-rep(0,nrow(dmodel))
for (i in 1:nrow(dmodel)) {
  if(dmodel$HY_tipo[i]%in%c("Piso", "Ático", "Chalet pareado", "Dúplex", "Chalet adosado", "Casa de pueblo", "Chalet independiente")){dmodel$HY_tipo_simplif2[i]<-"Vivienda"}
  if(dmodel$HY_tipo[i]%in%c("Aparcamiento","Trastero","Garaje","Otro")){dmodel$HY_tipo_simplif2[i]<-"Garajes"}
  if(dmodel$HY_tipo[i]%in%c("Local","Oficina","Almacén","Edificio turístico", "Edificio comercial")){dmodel$HY_tipo_simplif2[i]<-"Locales"}
  if(dmodel$HY_tipo[i]%in%c("Nave adosada","Nave aislada")){dmodel$HY_tipo_simplif2[i]<-"Naves"}
  if(dmodel$HY_tipo[i]%in%c("Solar","Suelo rústico no urbanizable")){dmodel$HY_tipo_simplif2[i]<-"Terrenos"}
}
destim$HY_tipo_simplif2<-rep(0,nrow(destim))
for (i in 1:nrow(destim)) {
  if(destim$HY_tipo[i]%in%c("Piso", "Ático", "Chalet pareado", "Dúplex", "Chalet adosado", "Casa de pueblo", "Chalet independiente")){destim$HY_tipo_simplif2[i]<-"Vivienda"}
  if(destim$HY_tipo[i]%in%c("Aparcamiento","Trastero","Garaje")){destim$HY_tipo_simplif2[i]<-"Garajes"}
  if(destim$HY_tipo[i]%in%c("Local","Oficina","Almacén","Edificio turístico", "Edificio comercial")){destim$HY_tipo_simplif2[i]<-"Locales"}
  if(destim$HY_tipo[i]%in%c("Nave adosada","Nave aislada")){destim$HY_tipo_simplif2[i]<-"Naves"}
  if(destim$HY_tipo[i]%in%c("Solar","Suelo rústico no urbanizable")){destim$HY_tipo_simplif2[i]<-"Terrenos"}
}

Se pueden ver las estadísticas por subgrupo

dmodel%>%filter(TARGET<300)%>%group_by(HY_tipo_simplif)%>%summarise(Cantidad=n(), 'Cantidad relativa' = n()/nrow(dmodel%>%filter(TARGET<300)), `Tiempo medio de visita (s)`=mean(TARGET),`Desviación típica del tiempo (s)`=sd(TARGET))%>%arrange(desc(Cantidad))
# A tibble: 3 x 5
  HY_tipo_simplif Cantidad `Cantidad relativa` `Tiempo medio de visita~ `Desviación típica del tiem~
  <chr>              <int>               <dbl>                    <dbl>                        <dbl>
1 Vivienda            5500               0.563                     69.8                         49.4
2 Garajes             2965               0.304                     50.6                         54.7
3 Resto               1296               0.133                     78.9                         49.4
dmodel%>%filter(TARGET<300)%>%group_by(HY_tipo_simplif2)%>%summarise(Cantidad=n(), 'Cantidad relativa' = n()/nrow(dmodel%>%filter(TARGET<300)), `Tiempo medio de visita (s)`=mean(TARGET),`Desviación típica del tiempo (s)`=sd(TARGET))%>%arrange(desc(Cantidad))
# A tibble: 5 x 5
  HY_tipo_simplif2 Cantidad `Cantidad relativ~ `Tiempo medio de visita~ `Desviación típica del tiem~
  <chr>               <int>              <dbl>                    <dbl>                        <dbl>
1 Vivienda             5500             0.563                      69.8                         49.4
2 Garajes              2977             0.305                      50.7                         54.8
3 Locales               929             0.0952                     77.3                         46.9
4 Terrenos              197             0.0202                     77.1                         53.1
5 Naves                 158             0.0162                     91.3                         54.3

Como se ha podido intuir la presencia de valores extremos de TARGET, vamos a estudiar en más profundidad los cuantiles de esta variable dependiendeo del tipo de inmueble.

dmodel %>% group_by(HY_tipo_simplif) %>%
     summarise(`75%`=quantile(TARGET, probs=0.75),
               `80%`=quantile(TARGET, probs=0.8),
               `85%`=quantile(TARGET, probs=0.85),
               `90%`=quantile(TARGET, probs=0.9),
               `95%`=quantile(TARGET, probs=0.95),
               avg=mean(TARGET),
               n=n(),
               n_rel=n()/nrow(dmodel))
# A tibble: 3 x 9
  HY_tipo_simplif `75%` `80%` `85%` `90%` `95%`   avg     n n_rel
  <chr>           <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl>
1 Garajes          71.0  86.7  107.  141.  219.  63.9  3044 0.306
2 Resto           107.  120.   133.  152.  201.  90.0  1332 0.134
3 Vivienda         91.1 101.   116.  140.  185.  76.0  5582 0.561
dmodel %>% group_by(HY_tipo_simplif2) %>%
     summarise(`75%`=quantile(TARGET, probs=0.75),
               `80%`=quantile(TARGET, probs=0.8),
               `85%`=quantile(TARGET, probs=0.85),
               `90%`=quantile(TARGET, probs=0.9),
               `95%`=quantile(TARGET, probs=0.95),
               avg=mean(TARGET),
               n=n(),
               n_rel=n()/nrow(dmodel))
# A tibble: 5 x 9
  HY_tipo_simplif2 `75%` `80%` `85%` `90%` `95%`   avg     n  n_rel
  <chr>            <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <int>  <dbl>
1 Garajes           71.1  86.8  107.  141.  219.  63.9  3056 0.307 
2 Locales          104.  116.   129.  146.  187.  88.5   952 0.0956
3 Naves            128.  139.   151.  196.  302. 111.    168 0.0169
4 Terrenos         110.  121.   138.  156.  188.  80.7   200 0.0201
5 Vivienda          91.1 101.   116.  140.  185.  76.0  5582 0.561 

Se puede ver que los cuantiles tienen unos valores de TARGET bastante distintos entre los tipos de inmueble, por lo que puede ser aconsejable en el conjunto de entrenamiento no simplemente filtrar TARGET, si no filtrar determinados valores dependiendo del tipo de inmueble.

ggplot(data=dmodel,aes(TARGET),alpha=0.2)+
  geom_histogram()+
  geom_vline(xintercept=quantile(unlist(dmodel%>%select(TARGET)),0.8), colour="red") +
  geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Garajes")%>%select(TARGET)),0.8), label="\n80%", y=3000), angle=90,colour="red", text=element_text(size=11)) +
  geom_vline(xintercept=quantile(unlist(dmodel%>%select(TARGET)),0.85), colour="blue") +
  geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Garajes")%>%select(TARGET)),0.85), label="\n85%", y=2000), angle=90,colour="blue", text=element_text(size=11)) +
  geom_vline(xintercept=quantile(unlist(dmodel%>%select(TARGET)),0.9), colour="green") +
  geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Garajes")%>%select(TARGET)),0.9), label="\n90%", y=1000), angle=90,colour="green", text=element_text(size=11)) +
  geom_vline(xintercept=quantile(unlist(dmodel%>%select(TARGET)),0.95), colour="black") +
  geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Garajes")%>%select(TARGET)),0.95), label="\n95%", y=500), angle=90,colour="black", text=element_text(size=11)) +
  ggtitle("Duración media de la visita a una página de un inmueble")
ggplot(data=dmodel%>%filter(TARGET<500)%>%group_by(HY_tipo_simplif),aes(TARGET,fill=HY_tipo_simplif))+
  geom_histogram(alpha=0.8)+
    geom_vline(xintercept=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Vivienda")%>%select(TARGET)),0.8), colour="black") +
    geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Vivienda")%>%select(TARGET)),0.8), label="80% Vivienda", y=1000), angle=90,colour="black", text=element_text(size=11)) +
  geom_vline(xintercept=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Garajes")%>%select(TARGET)),0.8), colour="black") +
  geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Garajes")%>%select(TARGET)),0.8), label="80% Garajes\n", y=1000), angle=90,colour="black", text=element_text(size=11)) +
  geom_vline(xintercept=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Resto")%>%select(TARGET)),0.8), colour="black") +
  geom_text(aes(x=quantile(unlist(dmodel%>%filter(HY_tipo_simplif=="Resto")%>%select(TARGET)),0.8), label="\n80% Resto", y=1000), angle=90,colour="black", text=element_text(size=11)) +
  ggtitle("Duración media de la visita a una página de un inmueble según tipo")

3.3 Análisis en función de la existencia de descripción

Se va a comprobar ahora la distribución de la variable TARGET que se quiere modelizar en función de la existencia o no de descripción en la página del inmueble, ya que en principio cuanto más larga (más detallada) sea la descripción, más largas serán las visitas. Lo primero que se hace es ver un resumen de los estadísticos más importantes, para después pasar a representar la variable empleando un boxplot y posteriormente un violin plot.

Se hace esto filtrando para valores de TARGET menores de 300, ya que la existencia de outliers hace más difícil el análisis gráfico posterior y resulta más conveniente de momento ignorar dichos outliers.

Primero se obtienen las distribuciones de estos datos.

dmodel_NOdesc<-dmodel%>%filter(HY_descripcion_exist==0,TARGET<300)
dmodel_SIdesc<-dmodel%>%filter(HY_descripcion_exist==1,TARGET<300)
cat("En los inmuebles con descripción tienen los siguientes estadísticos la variable TARGET:\n")
En los inmuebles con descripción tienen los siguientes estadísticos la variable TARGET:
summary(dmodel_SIdesc$TARGET)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   0.00   30.88   56.49   65.08   86.71  299.67 
cat("En los inmuebles sin descripción tienen los siguientes estadísticos la variable TARGET:\n")
En los inmuebles sin descripción tienen los siguientes estadísticos la variable TARGET:
summary(dmodel_NOdesc$TARGET)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   0.00   23.98   52.76   65.28   90.06  298.00 

A continuación, se representan gráficamente dichas distribuciones. Se observa cómo los resultados obtenidos son muy similares en los inmuebles con y sin descripción, sin embargo la mediana y el valor del primer cuartil son ligeramente menores en aquellos casos sin descripción. Observando la forma del violin plot representado se ve cómo la parte inferior de aquellos casos sin descripción es mayor que la de los casos con descripción, sugiriendo que en esa situación hay una mayor cantidad de gente que abre la página e inmediatamente, al ver que no hay descripción la cierra.

ggplot(dmodel%>%mutate(Tag=as.factor(ifelse(is.na(HY_descripcion),"Sin descripción","Con descripción"))), aes(x=Tag, y=TARGET, fill=Tag)) + 
  geom_violin(alpha=0.8)+
  ylim(c(0,250))+
  ggtitle("variación de TARGET según la presencia de descripción en el inmueble")

Sin embargo, no parece que sea una variable que influya tanto como se esperaría en el tiempo medio de visita.

3.4 Análisis en función de las imágenes

Como una de las hipótesis iniciales también es que a mayor número de fotos mayor tiempo en la página, se van a realizar una serie de gráficos para intentar visualizar la existencia de esta relación.

ggplot(dmodel, aes(x=as.factor(IM_n), y=TARGET, col=as.factor(IM_n))) + geom_violin(alpha=0.6)+theme(legend.position = "none")+
  ggtitle("TARGET según el número de imágnees")

Esta visualización iría a priori en contra de las hipótesis planteadas anteriormente, pues los inmuebles con más de 5 imágenes no llegan a alcanzar valores tan altos como las de menos de 5. Esto puede deberse a que hay un mayor número de inmuebles con pocas imágenes, por lo que la dispersión es más grande en este último caso. Se repite la visualización retirando los valores considerados extremos de TARGET.

ggplot(dmodel%>%filter(TARGET<200), aes(x=as.factor(IM_n), y=TARGET, col=as.factor(IM_n))) + geom_violin(alpha=0.6)+theme(legend.position = "none")+
  ggtitle("TARGET según el número de imágenes, para TARGETs no extremos")

En esta visualización sí que se puede ver como la media de TARGET es mayor conforme mayor es el número de imágenes. Quitar aquellos inmuebles con un TARGETsuperior a 350 implicaría eliminar sólo el -97.5137578.

Se comprueba ahora si influye el hecho de que existan fotos repetidas en el tiempo de visita.

dmodel%>%filter(TARGET<350,IM_n_rep<=5)%>%group_by(IM_n_rep)%>%summarize(`Cantidad`=n())
# A tibble: 6 x 2
  IM_n_rep Cantidad
     <dbl>    <int>
1        0     8893
2        1      622
3        2      159
4        3       89
5        4       35
6        5        6
ggplot(dmodel%>%filter(TARGET<350,IM_n_rep<=5), aes(x=as.factor(IM_n_rep), y=TARGET, col=as.factor(IM_n_rep))) + geom_violin(alpha=0.6)+theme(legend.position = "none")+
  ggtitle("TARGET según el número de imágnees repetidas")

En esta visualización se puede ver que cuantas más fotos repetidas hay, menor es el tiempo de visualización en general. La media de TARGET puede ser mayor en los inmuebles con más imágenes repetidas, pero su dispersión es menor, haciendo el tiempo de visita de estos inmuebles más consistente. Al haber más inmuebles con 0 imágenes repetidas, su dispersión aumenta y es más fácil encontrar valores distintos de TARGET a pesar del número de imágenes repetidas.

Se visualiza ahora si influye la similitud media entre las imágenes del inmueble

ggplot(dmodel%>%filter(TARGET<350), aes(x=IM_sim_mean, y=TARGET, col=as.factor(IM_n))) + geom_point(alpha=0.6)

En este gráfico se aprecia como en la mayoria de los casos hay pocas imágenes similares, lo cual es bueno para los clientes. Los inmuebles con imágenes con una similitud media mayor son aquellos con 2 o 3 imágenes. Sin embargo, el valor de TARGET parece estar bien distribuido para todas las similitudes. Puede que una combinación de la variable IM_n e IM_sim_mean sea una buena predictora, ya que la distribución de TARGET según el número de imágenes no es tan clara.

dmodel%>%group_by(IM_n)%>%summarize(n=n(),sim_mean=mean(IM_sim_mean),n_rep=mean(IM_n_rep),mean=mean(IM_mean,na.rm=T))
# A tibble: 11 x 5
    IM_n     n sim_mean  n_rep    mean
   <dbl> <int>    <dbl>  <dbl>   <dbl>
 1     0  2116    NA    0      NaN    
 2     1   611    NA    0        0.494
 3     2   828     1.07 0.165    0.496
 4     3  1161     1.04 0.106    0.483
 5     4   956     1.03 0.100    0.470
 6     5   864     1.03 0.0972   0.474
 7     6   605     1.03 0.180    0.471
 8     7   590     1.03 0.269    0.471
 9     8   454     1.03 0.233    0.468
10     9   380     1.03 0.261    0.473
11    10  1393     1.03 0.381    0.480

La mayoría de los inmuebles tiene 0, 3 o 4 imágenes. La similitud media y el número de imágenes repetidas aumenta con el número de imagenes.

3.5 Inmuebles repetidos

Tras el análisis previo se aprecia el hecho de que existen un gran número de inmuebles repetidos exactamente, con la descripción, la provincia, la fecha de subida y eliminación del inmueble de la página web, pero con distintos valores de TARGET. Se agrupan los datos por estas variables y se obtienen algunos de sus estadísticos para ver este hecho:

descripRepTodo<-dmodel%>%filter(HY_descripcion_exist==1)%>%
  group_by(HY_descripcion,HY_cod_postal,GA_quincena_ini,GA_quincena_ult)%>%
  summarise(Repetidos=n(),VisitasMedia=mean(GA_page_views),VisitasSD=sd(GA_page_views),MediaTarget=mean(TARGET),MedianaTarget=median(TARGET),DesvEstandarTarget=sd(TARGET))%>% filter(Repetidos>1)%>%arrange(desc(Repetidos))
head(descripRepTodo)
# A tibble: 6 x 10
# Groups:   HY_descripcion, HY_cod_postal, GA_quincena_ini [6]
  HY_descripcion HY_cod_postal GA_quincena_ini GA_quincena_ult Repetidos VisitasMedia VisitasSD
  <chr>          <chr>                   <dbl>           <dbl>     <int>        <dbl>     <dbl>
1 Plaza de gara~ 46183                      34              34        37         1.57     1.19 
2 Piso de dos d~ 04825                       1              46        15        19.2     13.3  
3 "Garaje ubica~ 29200                      20              20        14         1.21     0.426
4 Plaza de gara~ 04007                       2              47        13        48.8     26.9  
5 Plaza de gara~ 12593                      10              43        12         8.08     3.48 
6 Plaza de Gara~ 12130                      21              37        12         4.08     1.24 
# ... with 3 more variables: MediaTarget <dbl>, MedianaTarget <dbl>, DesvEstandarTarget <dbl>

Se tienen un total de 1015 cuando en realidad debería haber 328.

Se realiza la misma agrupación que antes sin las fechas de subida/bajada ya que puede ser que los propietarios suban y quiten un inmueble periódicamente para situarse mejor en el recomendador de la plataforma que gestiona los inmuebles:

descripRep<-dmodel%>%filter(HY_descripcion_exist==1)%>%
  group_by(HY_descripcion,HY_cod_postal)%>%
  summarise(Repetidos=n(),VisitasMedia=mean(GA_page_views),VisitasSD=sd(GA_page_views),MediaTarget=mean(TARGET),MedianaTarget=median(TARGET),DesvEstandarTarget=sd(TARGET))%>% filter(Repetidos>1)%>%arrange(desc(Repetidos))
head(descripRep)
# A tibble: 6 x 8
# Groups:   HY_descripcion [6]
  HY_descripcion HY_cod_postal Repetidos VisitasMedia VisitasSD MediaTarget MedianaTarget
  <chr>          <chr>             <int>        <dbl>     <dbl>       <dbl>         <dbl>
1 Plaza de gara~ 46183                76         3.68      6.19        88.6         42.8 
2 Promoción de ~ 47015                71         8.03      8.98        37.3         19.3 
3 Plaza de gara~ 04720                56         7.66      8.54        31.7         13.4 
4 Plaza de gara~ 03440                52         6.42      4.14        20.5          8.82
5 "Promoción de~ 26004                44        11        11.9         52.8         36.8 
6 "Garaje ubica~ 29200                39         1.74      1.65       116.          22   
# ... with 1 more variable: DesvEstandarTarget <dbl>

Se tienen un total de 3002 cuando en realidad debería haber 465.

Para tener en cuenta la repetición de inmbuebles, se añaden dos variables con el número de veces que están repetidos los inmuebles. Como estos inmuebles también pueden estar repetidos en destim, se utilian también los datos de este dataset.

#Creamos dataframes auxiliares
dtotalRep<-rbind(dmodel%>%select(HY_id,HY_descripcion_exist,HY_descripcion,HY_cod_postal,GA_quincena_ini,GA_quincena_ult,HY_metros_totales),destim%>%select(HY_id,HY_descripcion_exist,HY_descripcion,HY_cod_postal,GA_quincena_ini,GA_quincena_ult,HY_metros_totales))

descripRep<-dtotalRep%>%filter(HY_descripcion_exist==1)%>%group_by(HY_descripcion,HY_cod_postal,HY_metros_totales)%>% summarise(Repetidos=n())%>%arrange(desc(Repetidos))
descripTodoRep<-dtotalRep%>%filter(HY_descripcion_exist==1)%>%group_by(HY_descripcion,HY_cod_postal,HY_metros_totales,GA_quincena_ini,GA_quincena_ult)%>% summarise(Repetidos=n())%>%arrange(desc(Repetidos))

dmodel$HY_num_veces_rep<-rep(0,nrow(dmodel))
dmodel$HY_num_veces_todo_rep<-rep(0,nrow(dmodel))

for (x in 1:nrow(dmodel)){
    dmodel$HY_num_veces_rep[x]=ifelse(is.na(dmodel$HY_descripcion[x]),NA,ifelse(dmodel$HY_descripcion[x]%in% descripRep$HY_descripcion,
                                      descripRep$Repetidos[Reduce(intersect,list(which(descripRep$HY_descripcion==dmodel$HY_descripcion[x]) , which(descripRep$HY_cod_postal==dmodel$HY_cod_postal[x]) , which(descripRep$HY_metros_totales==dmodel$HY_metros_totales[x]) ))],
                                      0))
    dmodel$HY_num_veces_todo_rep[x]=ifelse(is.na(dmodel$HY_descripcion[x]),NA,ifelse(dmodel$HY_descripcion[x]%in% descripTodoRep$HY_descripcion,
                                           descripTodoRep$Repetidos[Reduce(intersect,list(which(descripTodoRep$HY_descripcion==dmodel$HY_descripcion[x]) , which(descripTodoRep$HY_cod_postal==dmodel$HY_cod_postal[x]) , which(descripTodoRep$HY_metros_totales==dmodel$HY_metros_totales[x]) , which(descripTodoRep$GA_quincena_ini==dmodel$GA_quincena_ini[x]) , which(descripTodoRep$GA_quincena_ult==dmodel$GA_quincena_ult[x])))],
                                           0))
}

destim$HY_num_veces_rep<-rep(0,nrow(destim))
destim$HY_num_veces_todo_rep<-rep(0,nrow(destim))

for (x in 1:nrow(destim)){
    destim$HY_num_veces_rep[x]=ifelse(is.na(destim$HY_descripcion[x]),
                                     NA,
                                     ifelse(destim$HY_descripcion[x]%in% descripRep$HY_descripcion,
                                            descripRep$Repetidos[Reduce(intersect,list(which(descripRep$HY_descripcion==destim$HY_descripcion[x]), which(descripRep$HY_cod_postal==destim$HY_cod_postal[x]), which(descripRep$HY_metros_totales==destim$HY_metros_totales[x])))],
                                            0))
    destim$HY_num_veces_todo_rep[x]=ifelse(is.na(destim$HY_descripcion[x]),
                                          NA,
                                          ifelse(destim$HY_descripcion[x]%in% descripTodoRep$HY_descripcion,
                                                 descripTodoRep$Repetidos[Reduce(intersect,list(which(descripTodoRep$HY_descripcion==destim$HY_descripcion[x]), which(descripTodoRep$HY_cod_postal==destim$HY_cod_postal[x]), which(descripTodoRep$HY_metros_totales==destim$HY_metros_totales[x]), which(descripTodoRep$GA_quincena_ini==destim$GA_quincena_ini[x]), which(descripTodoRep$GA_quincena_ult==destim$GA_quincena_ult[x])))],
                                                 0))
}

# destim$HY_num_veces_todo_rep[is.na(destim$HY_num_veces_todo_rep)]=0
# destim$HY_num_veces_rep[is.na(destim$HY_num_veces_rep)]=0
# dmodel$HY_num_veces_todo_rep[is.na(dmodel$HY_num_veces_todo_rep)]=0
# dmodel$HY_num_veces_rep[is.na(dmodel$HY_num_veces_rep)]=0

rm(dtotalRep)
rm(descripRep)
rm(descripTodoRep)

Se añaden, pues, las siguientes variables:

*HY_num_veces_rep: Número de veces que aparece un inmueble con la misma descripción en el mismo código postal con los mismos metros totales. NA si no hay descripción.

*HY_num_veces_todo_rep: Número de veces en las que aparece, además, en la misma semana y se elimina en la misma semana. NA si no hay descripción.

Se ve la relación de estas variables con TARGET:

ggplot(dmodel%>%filter(TARGET<350), aes(x=HY_num_veces_todo_rep, y=TARGET)) + geom_point(alpha=0.6)
ggplot(dmodel%>%filter(TARGET<350), aes(x=HY_num_veces_rep, y=TARGET)) + geom_point(alpha=0.6)

Se puede ver que, como en las visualizaciones anteriores, aquellos que tienen un número menor de repeticiones tienen valores de TARGET más elevados, pero se puede deber a que hay más inmuebles sin repeticiones. Los inmuebles repetidos tienen visitas de, por lo general,un minuto/minuto y medio.

3.6 Estudio de las variables de zona

Se tiene una variable HY_provincia que nos dice la provincia en la que se haya el inmueble.

ggplot(dmodel%>%group_by(HY_provincia)%>%filter(TARGET<150),aes(y=TARGET,x=HY_provincia,col=HY_provincia))+geom_violin()+theme(axis.text.x = element_text(angle = 90, hjust = 1),legend.position="none")
dmodel%>%group_by(HY_provincia)%>%filter(TARGET<150)%>%summarize(Cantidad=n(),`Cantidad relativa`=n()/nrow(dmodel),`TARGET medio`=mean(TARGET,na.rm=T),`TARGET sd`=sd(TARGET,na.rm=T))
# A tibble: 43 x 5
   HY_provincia Cantidad `Cantidad relativa` `TARGET medio` `TARGET sd`
   <chr>           <int>               <dbl>          <dbl>       <dbl>
 1 <NA>               40            0.00402            48.0       28.8 
 2 A Coruña           15            0.00151            49.4       42.0 
 3 Albacete           37            0.00372            56.7       27.5 
 4 Alicante          719            0.0722             61.0       37.4 
 5 Almería          1602            0.161              47.4       36.0 
 6 Asturias            1            0.000100           63.0       NA   
 7 Avila               9            0.000904           70.9       38.1 
 8 Badajoz             1            0.000100           96.6       NA   
 9 Baleares           72            0.00723            64.8       31.7 
10 Barcelona         147            0.0148             71.5       33.6 
11 Burgos              2            0.000201           39.7       20.8 
12 Cáceres             1            0.000100          144.        NA   
13 Cádiz              64            0.00643            65.9       32.7 
14 Castellón        1219            0.122              55.8       38.4 
15 Ciudad Real         4            0.000402           27.2       16.8 
16 Córdoba            28            0.00281            48.4       42.9 
17 Cuenca             19            0.00191            45.8       27.3 
18 Gerona             36            0.00362            67.3       41.7 
19 Granada           355            0.0356             43.6       35.1 
20 Guadalajara        19            0.00191            34.7       25.6 
21 Huelva              9            0.000904           79.2       24.7 
22 Huesca              3            0.000301           59.4       11.8 
23 Jaén                1            0.000100           53.6       NA   
24 La Rioja           42            0.00422            39.6       30.3 
25 Las Palmas        108            0.0108             65.1       38.9 
26 León               17            0.00171            61.5       28.9 
27 Lérida              5            0.000502           77.3       20.0 
28 Madrid             73            0.00733            61.1       36.3 
29 Málaga            479            0.0481             47.4       35.0 
30 Murcia           1386            0.139              57.5       33.4 
31 Navarra             1            0.000100           65.4       NA   
32 Palencia           17            0.00171            56.2       47.7 
33 Pontevedra          5            0.000502           45.9       47.2 
34 Segovia             1            0.000100           92.8       NA   
35 Sevilla             9            0.000904           62.1       46.7 
36 Soria              32            0.00321            35.4       36.4 
37 Tarragona         269            0.0270             59.5       34.8 
38 Tenerife            2            0.000201           76.0        9.00
39 Teruel             20            0.00201            56.0       30.9 
40 Toledo             14            0.00141            53.6       21.5 
41 Valencia         1899            0.191              59.8       35.4 
42 Valladolid        281            0.0282             38.4       33.9 
43 Zaragoza            7            0.000703           90.7       37.2 

Viendo la gráfica se podría pensar que el hecho de tener el inmueble en determinadas provincias influye bastante en el TARGET, pero viendo la tabla de estadísticas se puede ver que esas provincias con una desviación típica de TARGET pequeña como Ciudad Real o Huesca son aquellas con un número menor de inmuebles, por lo que no resultaría fiable incluir esta variable como predictora al modelo. Se cuenta, pues, con que la información que pueda aportar la zona al modelo la incluyan las variables IDEA_.

Sin embargo, sí parece importante reflejar de alguna forma en el modelo la facilidad que hay para encontrar inmuebles en una provincia, ya que se intuye que a mayor oferta, el tiempo de visita estará más repartido entre los inmuebles de una provincia. Se añade, pues, una variable que indique cuántas veces aparece una provincia. Se hará lo mismo con los códigos postales. Para que estas dos variables sean lo más ricas posibles, se imputan las provincias faltantes de las que se dispone el código postal. Se añade ahora: HY_provincia_frec: Número de veces que aparece una provincia HY_cod_postal_frec: Número de veces que aparece un código postal

#Imptutamos NAs de provincia con codigo postal
dmodel$HY_cod_postal[is.na(dmodel$HY_provincia)]
 [1] "04630" "04120" "30157" "30157" "04630" "30110" "30100" "04120" "04869" "04867" "03206" "30100"
[13] "30150" "30110" "03698" "03330" "04610" "04639" "03600" "03450" "04271" "03150" "03160" "03185"
[25] "04860" "04250" "07300" "03201" "04720" "04720" "04639" "04600" "30110" "04620" "04120" "06460"
[37] "30110" "30120" "08930" "29560" "04738" "04869"
dmodel$HY_provincia[is.na(dmodel$HY_provincia)]=c(rep("Almería",2),rep("Murcia",2),"Almería",rep("Murcia",2),rep("Almería",3),"Alicante",rep("Murcia",3),rep("Alicante",2),rep("Almería",2),rep("Alicante",2),"Almería",rep("Alicante",3),rep("Almería",2),"Baleares",rep("Alicante",1),rep("Almería",4),"Murcia","Almería","Almería","Badajoz",rep("Murcia",2),"Barcelona","Málaga",rep("Almería",2))

destim$HY_provincia[is.na(destim$HY_provincia)]=c("Almería","Málaga")

#Calculamos frecuencias aparición provincias y códigos postales
provincia_frec=rbind(dmodel%>%select(-TARGET),destim)%>%group_by(HY_provincia)%>%summarize(HY_provincia_frec=n())
cod_postal_frec=rbind(dmodel%>%select(-TARGET),destim)%>%group_by(HY_cod_postal)%>%summarize(HY_cod_postal_frec=n())

#Añadimos
dmodel<-left_join(dmodel,provincia_frec,by="HY_provincia")
dmodel<-left_join(dmodel,cod_postal_frec,by="HY_cod_postal")

destim<-left_join(destim,provincia_frec,by="HY_provincia")
destim<-left_join(destim,cod_postal_frec,by="HY_cod_postal")

rm(cod_postal_frec)
rm(provincia_frec)

3.7 Nuevas variables

Por último, de los análisis realizados al principio de esta sección y en base a la lógica se decide añadir las siguientes variables nuevas: GA_quincena_ini_year: Quincena del año en la que se colgó el anuncio en la página web GA_quincena_ult_year: Quincena del año en la que se quitó el anuncio en la página web GA_duracion: Número de quincenas que estuvo colgado el inmueble en la página web HY_precio_metro_cuadrado: Precio por metro cuadrado de inmueble HY_dif_precio: Diferencia en valor absoluto entre el precio anterior y el actual. NA si no hay precio anterior. HY_dif_precio_prc: Diferencia en porcentaje entre el precio anterior y el actual. NA si no hay precio anterior. HY_antiguedad_real: Antigüedad real del inmueble, calculado como 2018-HY_antiguedad GA_years: Número de años que ha estado el inmueble en la página web. 1,2 o over si ha estado durante todo el tiempo del que diponemos datos.

Además, buscando reforzar la teoría de que a mayor número de imágenes y mayor número de fotos, mayor tiempo en el inmueble, se crean estas dos variables nuevas:

HY_descripcion_words_2= HY_descripcion_words^2 IM_n_2 = IM_n^2

dmodel<-dmodel%>%mutate(GA_quincena_ini_year=ifelse(GA_quincena_ini>24,GA_quincena_ini-24,GA_quincena_ini),
                        GA_quincena_ult_year=ifelse(GA_quincena_ult>24,GA_quincena_ult-24,GA_quincena_ult),
                        GA_duracion=GA_quincena_ult-GA_quincena_ini,
                        HY_antiguedad_real=2018-HY_antiguedad,
                        HY_precio_metro_cuadrado=HY_precio/HY_metros_totales,
                        HY_dif_precio=ifelse(is.na(HY_precio-HY_precio_anterior),NA,HY_precio-HY_precio_anterior),
                        HY_dif_precio_prc=ifelse(is.na(HY_precio-HY_precio_anterior),NA,(HY_precio-HY_precio_anterior)/HY_precio_anterior),
                        IM_n_2=IM_n^2,
                        HY_descripcion_words_2=HY_descripcion_words^2)
destim<-destim%>%mutate(GA_quincena_ini_year=ifelse(GA_quincena_ini>24,GA_quincena_ini-24,GA_quincena_ini),
                        GA_quincena_ult_year=ifelse(GA_quincena_ult>24,GA_quincena_ult-24,GA_quincena_ult),
                        GA_duracion=GA_quincena_ult-GA_quincena_ini,
                        HY_antiguedad_real=2018-HY_antiguedad,
                        HY_precio_metro_cuadrado=HY_precio/HY_metros_totales,
                        HY_dif_precio=ifelse(is.na(HY_precio-HY_precio_anterior),NA,HY_precio-HY_precio_anterior),
                        HY_dif_precio_prc=ifelse(is.na(HY_precio-HY_precio_anterior),NA,(HY_precio-HY_precio_anterior)/HY_precio_anterior),
                        IM_n_2=IM_n^2,
                        HY_descripcion_words_2=HY_descripcion_words^2)


dmodel$HY_precio_metro_cuadrado[is.nan(dmodel$HY_precio_metro_cuadrado)]=NA
destim$HY_precio_metro_cuadrado[is.nan(destim$HY_precio_metro_cuadrado)]=NA

dmodel$GA_years<-ifelse(dmodel$GA_duracion<24,"1",ifelse(dmodel$GA_quincena_ini==1&dmodel$GA_quincena_ult==47,"over",2))
destim$GA_years<-ifelse(destim$GA_duracion<24,"1",ifelse(destim$GA_quincena_ini==1&destim$GA_quincena_ult==47,"over",2))

Además, también se realizó una exploración de las variables para ver su distribución utilizando una aplicación de Shiny que se adjunta, por lo que se aplican logaritmos a algunas de ellas. Se sumará 1 a las variables antes de logaritmizarlas para obtener valores reales. IDEA_pc_otros,IDEA_ind_tienda,IDEA_ind_alimentacion son ratios entre 0 y 1, por lo que se van a multiplicar por 100 para evitar dispersarlas erroneamente.

—HY_metros_utiles_log. —HY_metros_totales_log. —HY_precio_log. —IDEA_area_log. —IDEA_poblacion_log. —IDEA_densidad_log. —IDEA_rent_alquiler_log. —IDEA_unitprice_sale_residential_log. —IDEA_unitprice_rent_residential_log. —IDEA_price_sale_residential_log. —IDEA_price_rent_residential_log. —GA_page_views_log. —GA_mean_bounce_log. —GA_exit_rate_log. —HY_distribucion_length_log. —HY_distribucion_words_log. —HY_descripcion_length_log. —HY_descripcion_words_log. —HY_antiguedad_real_log. —IDEA_pc_otros_log. —IDEA_ind_tienda_log. —IDEA_ind_alimentacion_log. —HY_dif_precio_log.

LOGARITMIZADOR3000<-function(dataset, variables){
  #INPUT: variables es un vector de variables (char) del dataset que se quieren logaritmizar
  #OUTPUT: dataset logaritmizado
  for (variable in variables){
    dataset[paste0(variable,'_log')]=log(dataset[[variable]]+1)
  }
  dataset$HY_dif_precio_log=sapply(dataset$HY_dif_precio,function(x) ifelse(x<0,log(-x+1),log(x+1)))
  return(dataset)
}

variables<-c("HY_metros_utiles","HY_metros_totales","HY_precio","IDEA_area","IDEA_poblacion","IDEA_densidad","IDEA_rent_alquiler","IDEA_unitprice_sale_residential","IDEA_unitprice_rent_residential","IDEA_price_sale_residential","IDEA_price_rent_residential","GA_page_views","GA_mean_bounce","GA_exit_rate","HY_distribucion_length","HY_distribucion_words","HY_descripcion_length","HY_descripcion_words","HY_antiguedad_real","IDEA_pc_otros","IDEA_ind_tienda","IDEA_ind_alimentacion","HY_dif_precio")
dmodel<-mutate(dmodel,IDEA_pc_otros=100*IDEA_pc_otros,IDEA_ind_tienda=100*IDEA_ind_tienda,IDEA_ind_alimentacion=100*IDEA_ind_alimentacion)
destim<-mutate(destim,IDEA_pc_otros=100*IDEA_pc_otros,IDEA_ind_tienda=100*IDEA_ind_tienda,IDEA_ind_alimentacion=100*IDEA_ind_alimentacion)

dmodel<-LOGARITMIZADOR3000(dmodel,variables)
destim<-LOGARITMIZADOR3000(destim,variables)

Se añade también una columna (HY_num_habs) que sea el número de habitaciones del inmueble. Este será extraído mediante expresiones regulares del campo descripción o de distribución. La función numHabs() se encuentra en el fichero adjunto ExpresionesRegularesNumHabs.R

source('ExpresionesRegularesNumHabs.R')
dmodel<-numHabs(dmodel)
destim<-numHabs(destim)

3.8 Análisis de las variables GA_

Las variables GA_ eran las que estaban relacionadas directamente con datos de la interacción de los usuarios con la página web, con lo cual a simple vista parece que serán variables muy representativas para poder estimar el tiempo medio de visita de un inmueble en la página web. Por eso se les dedica una sección a parte para estudiar sus distribuciones.

3.8.1 Análisis en función de GA_quincena_ini_year

dmodel%>%ggplot()+geom_histogram(aes(x=GA_quincena_ult_year))
dmodel%>%ggplot()+geom_histogram(aes(x=GA_quincena_ini_year))

Se podría decir que durante el periodo estudiado la tendencia fue subir la mayoría de inmuebles a principios de año y dar de baja la mayoría de inmuebles a finales del periodo. Esto se puede deber a que los pisos que ya se encontraban en la página web cuando empezó el periodo de extracción de datos se puso que su quincena inicial era 1, y a que todos los que quedaban cuando se terminó el periodo se etiquetaron con el mismo número.

3.8.2 Análisis en función de GA_mean_bounce

También podría estudiarse la distribución de la variable GA_mean_bounce que es el porcentaje de veces que el usuario abandonó la página sin interactuar con ella.

a<-data.frame(Tag="Con descripción",GA_mean_bounce=dmodel_SIdesc$GA_mean_bounce)
b<-data.frame(Tag="Sin descripción",GA_mean_bounce=dmodel_NOdesc$GA_mean_bounce)
datos_bounce<-rbind(a,b)

ggplot(datos_bounce, aes(x=Tag, y=GA_mean_bounce, fill=Tag)) + geom_violin()

Se ve cómo el porcentage de gente que abandona la página sin interactuar con ella es mayor cuando el inmueble no cuenta con descripción respecto a cuando sí tiene. En aquellos inmuebles que tienen descripción se ve cómo los usuarios interaccionan en mucha mayor medida con la página respecto a aquellos que no tienen descripción.

ggplot(dmodel,aes(x=TARGET,y=GA_mean_bounce))+geom_violin()

Se puede ver dos distribuciones más o menos claras, diferenciándose para cuango GA_mean_bouncetoma un valor menor a 10, que nos pueden permitir conseguir un filtro en los datos.

p1<-ggplot(dmodel%>%filter(GA_mean_bounce<10),aes(x=TARGET,y=GA_mean_bounce))+geom_violin()
p2<-ggplot(dmodel%>%filter(GA_mean_bounce>=10),aes(x=TARGET,y=GA_mean_bounce))+geom_violin()
multiplot(p1,p2,cols=2)

Se puede ver que realmente el cambio no está en 10, sino en aproximadamente 0.

dmodel%>%group_by(GA_mean_bounce==0)%>%summarize(Cantidad=n(),`TARGET medio`=mean(TARGET),`TARGET sd`=sd(TARGET))
# A tibble: 2 x 4
  `GA_mean_bounce == 0` Cantidad `TARGET medio` `TARGET sd`
  <lgl>                    <int>          <dbl>       <dbl>
1 FALSE                     6917           75.7        60.2
2 TRUE                      3041           70.6       135. 
dmodel%>%group_by(GA_mean_bounce==0)%>%filter(TARGET<150)%>%summarize(Cantidad=n(),`TARGET medio`=mean(TARGET),`TARGET sd`=sd(TARGET))
# A tibble: 2 x 4
  `GA_mean_bounce == 0` Cantidad `TARGET medio` `TARGET sd`
  <lgl>                    <int>          <dbl>       <dbl>
1 FALSE                     6363           62.9        34.1
2 TRUE                      2707           36.2        34.4

3.8.3 Análisis en función de GA_years

ggplot(dmodel%>%filter(TARGET<150))+
  geom_violin(aes(x=GA_years,y=TARGET,color=GA_years))+theme(legend.position="none")+ggtitle("Influencia del número de años en TARGET")

Se puede ver que influye bastante el número de años que ha estado el inmueble en la página web, alargándose las distribuciones para los inmuebles con más tiempo en la misma, aumentando el tiempo medio de visita.

3.9 Correlaciones

Se va a estudiar ahora las posibles correlaciones lineales entre cada una de las variables para ello se van a emplear matrices de correlaciones, que dan la variación entre las variables relacionándolas 1 a 1. Se van a obtener 4 matrices distintas, la primera conteniendo la información relativa al inmueble (HY_), la segunda con la información económica de la zona (IDEA_) y la tercera con la información de la página web (GA_) y la cuarta con las imágenes (IM_). Además en todas ellas se incluye la variable TARGET que es la que se quiere modelizar y con la que se busca principalmente establecer las relaciones.

corrplot(cor(dmodel%>%select(HY_descripcion_length,HY_antiguedad:HY_precio_anterior,TARGET,-HY_cert_energ),use = "complete.obs"),method="ellipse",type="upper",tl.col="black",main="HY Variables")
corrplot(cor(dmodel%>%select(starts_with('IDEA'),TARGET),use = "complete.obs"),method="ellipse",type="upper",tl.cex = 0.5,tl.col="black",main="IDEA Variables")
corrplot(cor(dmodel%>%select(starts_with('GA'),-GA_years,TARGET),use = "complete.obs"),method="ellipse",type="upper",tl.col="black",main="GA Variables")
corrplot(cor(dmodel%>%select(starts_with('IM'),TARGET),use = "complete.obs"),method="ellipse",type="upper",tl.col="black",main="IM Variables")

Se ve que no hay apenas relaciones lineales con la variable TARGET. A continuación van a volver a relizarse estas mismas matrices de correlaciones, pero en esta ocasión se filtrarán aquellos valores muy grandes de la variable TARGET para poder eliminar los outliers que no interesa estudiar.

A continuación se filtra para valores menores de 200 y observaciones con los datos completos.

corrplot(cor(dmodel[,sapply(dmodel,is.numeric)]%>%filter(TARGET<200),use = "complete.obs"),method="ellipse",type="upper",tl.col="black",main="HY Variables")
rm(temp)
corrplot(cor(dmodel%>%select(starts_with('IDEA'),TARGET)%>%filter(TARGET<200),use = "complete.obs"),method="ellipse",type="upper",tl.cex = 0.5,tl.col="black",main="IDEA Variables")
corrplot(cor(dmodel%>%select(starts_with('GA'),-GA_years,TARGET)%>%filter(TARGET<200),use = "complete.obs"),method="ellipse",type="upper",tl.col="black",main="GA Variables")
corrplot(cor(dmodel%>%select(starts_with('IM'),TARGET)%>%filter(TARGET<200),use = "complete.obs"),method="ellipse",type="upper",tl.col="black",main="IM Variables")

Se observa cómo cuanto menor sea el valor máximo de TARGET considerado y, por lo tanto, mayor sea el número de outliers excluido, las correlaciones entre la variable output y el resto de variables van aumentando, aunque en ningún caso este valor del coeficiente de correlación llega a ser significativo, sino que se encuentra en un rango cercano a 0.4. Después de visualizar estas relaciones, probamos a realizar modelos lineales, pero dado que ninguna variable tiene un alta correalción con TARGET, no se consiguieron buenos resultados.

3.10 Filtros

Mientras se hacía el análisis exploratorio, se han observado de diversas incongrugencias que se van a “reparar” sustituyéndolas por NAs. Algunas de estas son, como ya se ha mencionado antes, los inmubles con más de 20 baños o una antigüedad superior a 2018. Por otro lado, los inmuebles con un número de metros totales de 999,99999,9999 o 0 parecen sospechosos, por los que se sustituirán por NA.

dmodel2<-dmodel
dmodel2$HY_antiguedad<-sapply(dmodel2$HY_antiguedad,function(x) ifelse(!is.na(x)& x>2018,NA,x))
dmodel2$HY_metros_totales<-sapply(dmodel2$HY_metros_totales,function(x) ifelse(!is.na(x)& x%in%c(99999,9999,999,0),NA,x))
dmodel2$HY_num_banos<-sapply(dmodel2$HY_num_banos,function(x) ifelse(!is.na(x) & x>=20,NA,x))
dmodel2$HY_num_terrazas<-sapply(dmodel2$HY_num_terrazas,function(x) ifelse(!is.na(x) & x>=10,NA,x))
dmodel2$IDEA_demand_rent_residential<-sapply(dmodel2$IDEA_demand_rent_residential,function(x) ifelse(!is.na(x) & x==0,NA,x))
dmodel2$IDEA_demand_sale_residential<-sapply(dmodel2$IDEA_demand_sale_residential,function(x) ifelse(!is.na(x) & x==0,NA,x))

destim2<-destim
destim2$HY_antiguedad<-sapply(destim2$HY_antiguedad,function(x) ifelse(!is.na(x)& x>2018,NA,x))
destim2$HY_metros_totales<-sapply(destim2$HY_metros_totales,function(x) ifelse(!is.na(x)& x%in%c(99999,9999,999,0),NA,x))
destim2$HY_num_banos<-sapply(destim2$HY_num_banos,function(x) ifelse(!is.na(x) & x>=20,NA,x))
destim2$HY_num_terrazas<-sapply(destim2$HY_num_terrazas,function(x) ifelse(!is.na(x) & x>=10,NA,x))
destim2$IDEA_demand_rent_residential<-sapply(destim2$IDEA_demand_rent_residential,function(x) ifelse(!is.na(x) & x==0,NA,x))
destim2$IDEA_demand_sale_residential<-sapply(destim2$IDEA_demand_sale_residential,function(x) ifelse(!is.na(x) & x==0,NA,x))

Antes de proceder a la producción de modelos, se definen los conjuntos de datos y funciones que se usarán en la siguiente etapa. - dmodelC: dmodel, quitando variables con bastantes NAs - Función completos: Dado un dataset, devuelve sólo las observaciones completas - Función quitarOutliers: Dado un dataset y nombres de variables del mismo, quita los outliers de la variable que se quiera - Función quitarXSup: Dado un dataset y nombres de variables del mismo, quita las observaciones que en la VARIABLE que se quiera tengan un valor superior al percentil X - Función imputacionNAs: Dado un dataset, se imputan los valores faltantes del mismo con la moda si son categóricos o con la media si son numéricos.

##Creación de datasets
dmodelC<-dmodel%>%select(-HY_id,-HY_descripcion,-HY_cod_postal,-HY_distribucion,-HY_cert_energ,-HY_antiguedad,-HY_precio_anterior,-IDEA_ind_elasticidad,-IDEA_ind_liquidez,-HY_metros_utiles) 
dmodelC$HY_provincia<-as.factor(dmodelC$HY_provincia)
dmodelC$HY_tipo<-as.factor(dmodelC$HY_tipo)
dmodelC$HY_tipo_simplif<-as.factor(dmodelC$HY_tipo_simplif)
dmodelC$HY_tipo_simplif2<-as.factor(dmodelC$HY_tipo_simplif2)
dmodelC$GA_years<-as.factor(dmodelC$GA_years)

## FUNCIONES DE INTERES
completos<-function(dataset){dataset[complete.cases(dataset),]}

quitarOutliers<-function(dataset,variable){dataset[!dataset[,variable]%in%boxplot.stats(dataset[,variable])$out,]}

quitarXSup<-function(dataset,variable,x){dataset[dataset[,variable]<quantile(unlist(dataset[,variable]),x),]}

imputacionNAs<-function(dataset){ dataset<-dataset %>% 
  mutate_if(is.numeric, 
            .funs = funs(
              ifelse(is.na(.), 
                     median(., na.rm = TRUE),
                     .))) %>%
  mutate_if(is.character, 
            .funs = funs(
              ifelse(is.na(.), 
                     getmode(.),
                     .)))
return(dataset)}

4 Machine Learning

Con todo el conocimiento y tratamiento que se ha llevado a cabo se procede ya a construir modelos predictores utilizando la librería mlr.

4.1 Elección y entrenamiento de modelos, fase inicial

Se crea la función corredorModelos a la que se le pasa el dataset que se quiere probar y los nombres de los modelos a probar. Con esta función lo que se quiere obtener es el mejor modelo dado un dataset. Tiene incluida una selección de características con método Wrapper para escoger un buen subset de características para cada modelo. Se define la función corredorModelos.

corredorModelos<-function(dataset,learners,crossval=T,norm=F){
  #crossval=T hará cross validation, mientras que crossval=F hará un único entrenamiento
  #norm=T hará que se normalicen las variables numéricas
  #OUTPUT: Performance
  if (norm){dataset<-normalizeFeatures(dataset,target="TARGET",method="standardize")}
  dataset<-createDummyFeatures(dataset,target="TARGET")
  task<-makeRegrTask(data=dataset,target="TARGET")
  
  cv=makeResampleDesc("CV",iters=5) 
  #ctrl=makeFeatSelControlRandom(maxit=10)
  ctrl=makeFeatSelControlSequential(method="sfs",max.features =20) 
  rdesc=makeResampleDesc("Holdout")
  perf=list()

  for (learner in learners){
    print(paste('learner:',learner))
    sfeats = selectFeatures(learner = learner, task = task, resampling = cv,
    control = ctrl, show.info =TRUE)
    task=subsetTask(task,features=sfeats$x)
    print(paste('features:',sfeats$x))
    if(crossval==TRUE){
      res=resample(learner,task,cv,medae,show.info = F)
      perf[learner]=res$aggr
    }else{
      ho<-makeResampleInstance("Holdout",task)
      tsk.train<-subsetTask(task,ho$train.inds[[1]])
      tsk.test<-subsetTask(task,ho$test.inds[[1]])
      model=train(learner,tsk.train)
      pred=predict(model,tsk.test)
      perf[learner]=performance(pred=pred,measure=medae)
    }
    print(paste('learner',learner,'performance',perf[learner]))
  }
  return(perf)
}

Para saber qué modelos se puede usar para el problema se utiliza la función auxiliar de mlr:

task<-makeRegrTask(data=completos(dmodelC),target="TARGET") 
listLearners(task) #para ver los learners aptos para nuestra tarea
             class                                                                    name
1 regr.bartMachine                                      Bayesian Additive Regression Trees
2       regr.bcart                                                           Bayesian CART
3        regr.brnn                Bayesian regularization for feed-forward neural networks
4        regr.btgp                                         Bayesian Treed Gaussian Process
5     regr.btgpllm Bayesian Treed Gaussian Process with jumps to the Limiting Linear Model
6        regr.btlm                                             Bayesian Treed Linear Model
   short.name     package type installed numerics factors ordered missings weights  prob oneclass
1 bartmachine bartMachine regr      TRUE     TRUE    TRUE   FALSE     TRUE   FALSE FALSE    FALSE
2       bcart         tgp regr      TRUE     TRUE    TRUE   FALSE    FALSE   FALSE FALSE    FALSE
3        brnn        brnn regr      TRUE     TRUE    TRUE   FALSE    FALSE   FALSE FALSE    FALSE
4        btgp         tgp regr      TRUE     TRUE    TRUE   FALSE    FALSE   FALSE FALSE    FALSE
5     btgpllm         tgp regr      TRUE     TRUE    TRUE   FALSE    FALSE   FALSE FALSE    FALSE
6        btlm         tgp regr      TRUE     TRUE    TRUE   FALSE    FALSE   FALSE FALSE    FALSE
  twoclass multiclass class.weights featimp oobpreds functionals single.functional    se lcens
1    FALSE      FALSE         FALSE   FALSE    FALSE       FALSE             FALSE FALSE FALSE
2    FALSE      FALSE         FALSE   FALSE    FALSE       FALSE             FALSE  TRUE FALSE
3    FALSE      FALSE         FALSE   FALSE    FALSE       FALSE             FALSE FALSE FALSE
4    FALSE      FALSE         FALSE   FALSE    FALSE       FALSE             FALSE  TRUE FALSE
5    FALSE      FALSE         FALSE   FALSE    FALSE       FALSE             FALSE  TRUE FALSE
6    FALSE      FALSE         FALSE   FALSE    FALSE       FALSE             FALSE  TRUE FALSE
  rcens icens
1 FALSE FALSE
2 FALSE FALSE
3 FALSE FALSE
4 FALSE FALSE
5 FALSE FALSE
6 FALSE FALSE
... (#rows: 41, #cols: 24)

En un principio se filtran los datos según algunos de los criterios comentados antes y algunos criterios nuevos, creando el que se pensaba que sería un buen dataset tras el análisis exploratorio de datos: -Un inmueble no puede tener el valor 0 en la variable HY_metros_totales ya que no tiene sentido. -Un inmueble del tipo de los que se tienen en los datos no tendrá más de 20 baños. -El precio de un inmueble no puede ser 0 ya que los resultados de la variable TARGET que produzca serán erróneos. -Si el 100% de la gente que ha visitado un inmueble ha salido sin interactuar con la página, dichos inmuebles no generarán valores de la variable TARGET que aporten conocimiento. -Si el tiempo medio de visita a un inmueble es 0 quiere decir que se debe a un error o alguna situación semejante por lo que tampoco se tendrán en cuenta los registros en los que se de esto.

Como todos los modelos no son capaces de trabajar con valores faltantes, se van a utilizar la función completos() para construir el dataset que se manda a corredorModelos(). De esta forma quedan con pocos inmuebles, pero sirve para realizar una primera aproximación a los distintos modelos. También se trabajó con imputacionNAs(), pero producía peores resultados. Como el dataset tiene bastantes valores faltantes, parece lógico pensar que se preferirán modelos que sean capaces de trabajar con ellos como random forest, xgboost o cubist.

dataset<-dmodelC%>%
  dplyr::select(-HY_provincia,-HY_num_terrazas,-IDEA_unitprice_rent_residential,IDEA_price_rent_residential,-IDEA_stock_rent_residential,-HY_ascensor)%>%
  filter(HY_metros_totales!=0,!HY_num_banos%in%c(20,99),HY_precio!=0,GA_page_views>0,GA_mean_bounce!=100,GA_exit_rate!=100,TARGET!=0)

dataset<-as.data.frame(dataset)
dataset<-quitarOutliers(dataset,"GA_mean_bounce")
dataset<-quitarOutliers(dataset,"TARGET")
#dataset<-imputacionNAs(dataset)
dataset<-completos(dataset) 
dataset<-convierteAFactor(dataset)
learners=c('regr.bcart','regr.kknn','regr.svm','regr.nnet','regr.randomForest','regr.xgboost')
performance<-corredorModelos(dataset,learners,crossval=T,norm=F) 

Se prueba lo mismo tanto con norm=Fcomo norm=T en varias etapas del proceso, pero funcionaban mejor los modelos cuando las variables no estaban normalizadas. Se vio que SVM producía los mejores resultados, por lo que se continuó trabajando con él. Más adelante, sin embargo, comprobamos el funcionamiento del modelo Cubist y, ya que funciona mejor que SVM, se apusta por este modelo. Cubist es un modelo basado en reglas de decisión que en cada hoja aplica un modelo lineal a los datos, por lo que puede trabajar con variables categóricas, numéricas y valores faltantes. Si además se dummifican las variables categóricas, aparte de para el proceso de separación del árbol, también pueden ser usadas en los modelos lineales.

Se fueron probando selecciones aleatorias de características como se puede ver en el siguiente código para conseguir unos buenos resultados. Así se consiguió la selección de características que produjo los resultados de la preeentrega de la fase local. En este segmento de código ya se puede ver como se implementa la predicción corregida, por la que en los casos en los que el inmueble tenga GA_page_views=0 o GA_mean_bounce=100 se asignará directamente un 0 al TARGET predicho. Además, también se empieza a filtar los Outliers de TARGET.

dmodelC_selec<-dmodelC%>%
  dplyr::select(-HY_provincia,-HY_num_terrazas,-IDEA_unitprice_rent_residential,IDEA_price_rent_residential,-IDEA_stock_rent_residential,-HY_ascensor)%>%
  filter(HY_metros_totales!=0,!HY_num_banos%in%c(20,99),HY_precio!=0,GA_page_views>0,GA_mean_bounce!=100,GA_exit_rate!=100,TARGET!=0)
dmodelC_selec<-select(dmodelC_selec,-HY_tipo)
dataset<-as.data.frame(dmodelC_selec)
dataset<-convierteAFactor(dataset)
dataset<-createDummyFeatures(dataset,target="TARGET")
row_indices<-sample(1:nrow(dataset),size=0.8*nrow(dataset))
train_set<-dataset[row_indices,]
test_set<-dataset[-row_indices,]
test_set<-completos(test_set)
#Preprocesado extra para los datos de train
train_set<-quitarOutliers(train_set,"GA_mean_bounce")
train_set<-quitarOutliers(train_set,"TARGET")
train_set<-completos(train_set)

task<-makeRegrTask(data=train_set,target="TARGET")
cv=makeResampleDesc("CV",iters=5) 
ctrl=makeFeatSelControlRandom(maxit=100) ## cambiar a secuencial
rdesc=makeResampleDesc("Holdout")

sfeats = selectFeatures(learner = "regr.cubist", task = task, resampling = cv,
control = ctrl, show.info = TRUE,measure=medae)
task=subsetTask(task,features=sfeats$x)
print(paste('features:',sfeats$x))
model=train("regr.cubist",task)
pred=predict(model,task)
perf=performance(pred=pred,measure=medae)

test_set_model<-completos(test_set)
test_set_model<-test_set_model%>%select(unlist(sfeats$x))
pred=predict(model,newdata=test_set_model)
perf_test=measureMEDAE(test_set$TARGET,pred$data$response)
pred_corregida=pred$data$response
for (i in 1:nrow(test_set)){
  pred_corregida[i]<-ifelse((test_set$GA_exit_rate[i]==100|test_set$GA_page_views[i]==0),0,pred_corregida[i])
}
perf_corregido=measureMEDAE(test_set$TARGET,pred_corregida)
print(paste('perf_train:',perf,', perf_test:',perf_test,', perf_corregido:',perf_corregido))

Las características seleccionadas fueron: [1] “HY_num_garajes” “HY_precio”
[3] “IDEA_densidad” “IDEA_pc_1960_69”
[5] “IDEA_pc_2000_10” “IDEA_pc_comercio”
[7] “IDEA_pc_industria” “IDEA_ind_riqueza”
[9] “IDEA_rent_alquiler” “IDEA_demand_sale_residential” [11] “IDEA_price_rent_residential” “GA_page_views”
[13] “GA_mean_bounce” “GA_exit_rate”
[15] “GA_quincena_ini” “GA_quincena_ult”
[17] “HY_distribucion_exist” “HY_descripcion_length”
[19] “HY_descripcion_exist” “HY_precio_metro_cuadrado”
[21] “GA_duracion”

Sin embargo, este método de selección de características no era el óptimo, ya que depende mucho de la inicialización. Lo ideal habría sido correr una búsqueda exhaustiva, pero dado que requiere mucho poder computacional, nos decidimos por aplicar métodos Filter.

4.2 Selección de características con métodos Filter

Para todo buen modelo predictor un paso fundamental es la selección de características para reducir la complejidad del modelo intentando no perder acierto en la predicción. Se utilizaron métodos Filter de la librería mlr en base a distintas métricas. En este punto, además, tras el uso del Shiny comentado anteriormente, se apreció que quizás sería mucho mejor intentar predecir el log(TARGET) en vez de TARGET dado su distribución. En este punto del proyecto, se empieza a trabajar en paralelo con modelos para predecir TARGET y para predecir log(TARGET), obteniendo tanto los mejores parámetros como carácterísticas para ambos posibles modelos. Además, dado que el modelo escogido puede trabajar con valores faltantes, se trabaja con dmodel2 en vez de dmodel, ya que parece un dataset más preciso.

Los métodos escogidos para ver la importancia de las características son los siguientes: - cforest.importance: Se escoge este método porque está basado en árboles de decisión y nuestro modelo es, en gran parte, un árbol de decisión. - gain.ratio: Devuelve la independencia/dependencia entre cada característica y el TARGET utilizando el test de independencia estadísitca Chi-cuadrado - information.gain: Método basado en la entropía que devuelve la información que se gana utilizando una característica respecto a TARGET.

set.seed(13)
#Selección de características para log(TARGET)
A=dmodel2%>%select(-HY_id,-HY_descripcion,-HY_distribucion,-HY_provincia,-HY_cod_postal)
A$HY_tipo<-as.factor(A$HY_tipo)
A$HY_tipo_simplif<-as.factor(A$HY_tipo_simplif)
A$HY_tipo_simplif2<-as.factor(A$HY_tipo_simplif2)
A$HY_cert_energ<-as.factor(A$HY_cert_energ)
A$GA_years<-as.factor(A$GA_years)
A=A[sample(1:nrow(A),4000),]
A$TARGET=log(A$TARGET+1)

#Selección de caracetrísticas para TARGET con variables dummy
A2=createDummyFeatures(A)
dmodel.Task=makeRegrTask(id=deparse(substitute(A2)),A2,target='TARGET' )
fv = generateFilterValuesData(dmodel.Task,method =c("cforest.importance","gain.ratio","information.gain"))
plotFilterValues(fv) + ggpubr::theme_pubr()
impVariablesTodasDummy<-as.data.frame(fv$data)
save(impVariablesTodasDummy,file="impVariablesTodasDummy.RData")

dmodel.Task=makeRegrTask(id=deparse(substitute(A)),A,target='TARGET' )
fv = generateFilterValuesData(dmodel.Task,method =c("cforest.importance","gain.ratio","information.gain"))
plotFilterValues(fv) + ggpubr::theme_pubr()
impVariablesTodas<-as.data.frame(fv$data)
save(impVariablesTodas,file="impVariablesTodas.RData")

Tras ver los resultados de los distintos métodos, se plantea como características importantes: - GA: GA_exit_rate, GA_page_views_log, GA_duracion, GA_mean_bounce - HY: HY_dif_precio_log, HY_distribucion_words_log, HY_metros_totales_log, HY_provincia_frec, HY_cod_postal_frec, HY_descripcion_words_2, HY_num_veces_rep, HY_precio_metro_cuadrado - IM: IM_n, IM_n_2, IM_sim_mean - IDEA: IDEA_demand_sale_residential, IDEA_unitprice_sale_residential_log, IDEA_unitprice_rent_residential, IDEA_price_sale_residential, IDEA_price_rent_residential, IDEA_ind_turismo, IDEA_pc_otros_log, IDEA_pc_2000_10, IDEA_area, IDEA_pc_Comercio

4.3 Tuning

Buscando mejores resultados, se desea saber cuáles son los mejores hiperparámetros para el modelo. El modelo puede ajustarse con los siguientes parámetros: - committees: Este es el número de modelos que se construyen para obtener la predicción final - rules: Número de reglas máximas que construye cada modelo, pero pueden ser menos. - extrapolation: Como Cubist usa modelos lineales, este parámetro nos indica cuánto (%) se han de ajustar las predicciones dadas por este modelo al conjunto de entrenamiento - unbiased: Indica si las reglas que construyen los modelos son unbiased o no.

features<-names(dmodel2%>%select(starts_with('GA'),-ends_with("_log"),HY_precio_anterior,HY_dif_precio_prc,IM_n,IM_sim_mean,HY_num_veces_todo_rep,HY_descripcion_words,HY_num_habs,IDEA_densidad,HY_metros_totales,IDEA_price_rent_residential,IDEA_price_sale_residential,IDEA_ind_turismo,IDEA_pc_2000_10,HY_tipo_simplif2,HY_tipo_simplif,HY_distribucion_length))
dmodel2$HY_provincia<-as.factor(dmodel2$HY_provincia)
dmodel2$HY_tipo<-as.factor(dmodel2$HY_tipo)
dmodel2$HY_tipo_simplif<-as.factor(dmodel2$HY_tipo_simplif)
dmodel2$HY_tipo_simplif2<-as.factor(dmodel2$HY_tipo_simplif2)
dmodel2$HY_cert_energ_exist<-as.factor(dmodel2$HY_cert_energ_exist)
dmodel2$HY_descripcion_exist<-as.factor(dmodel2$HY_descripcion_exist)
dmodel2$HY_distribucion_exist<-as.factor(dmodel2$HY_distribucion_exist)

dataset<-as.data.frame(dmodel2%>%mutate(TARGET_log=log(TARGET+1))%>%select(features,TARGET))
dataset<-createDummyFeatures(dataset)

task=makeRegrTask(id=deparse(substitute(dataset)),data=dataset,target="TARGET ")
ps=makeParamSet(
  makeDiscreteParam("committees",values=c(1,25,50,75,100)),  
  makeLogicalParam(id="unbiased"),
  makeDiscreteParam("rules",values=c(10,100,200,300,400)),
  makeDiscreteParam("extrapolation",values=c(0,25,50,75,100))
)

ctrl=makeTuneControlGrid()

rdesc=makeResampleDesc("CV",iters=5)

res=tuneParams("regr.cubist",
               task=task,
               resampling=rdesc,
               par.set=ps,
               control=ctrl,
               measures = list(medae,setAggregation(medae,test.sd))
               )

El resultado de este código se puede encontrar adjunto. Como varios hiperparámetros producian resultados muy similares, se fueron probando entre los que provocaban menos error, teniendo siempre en cuenta que la función tuneParams() no permite filtrar el conjunto de entrenamiento o corregir las predicciones. Finalmente, los hiperparámetros escogidos son: - committees: 100 - rules: 100 - extrapolation: 100, aunque se obtienen unos resultados muy similares con el valor 0 - unbiased: TRUE

4.4 Entrenamientos finales

Se comprueba la selección de características y parámetros final mediante un CV, después del cual se realizará el entrenamiento final con todos los datos. Se probó a filtrar los datos de train y se filtraron los inmuebles con valores de TARGETcon un valor superior a un percentiles. Entre los percentiles, el percentil 80 resultó ser el que produjo mejores resultados, eliminando aproximadamente el 20 de los inmuebles de training. En un primer momento se filtraron todos aquellos superiores al percentil 80, pero después se comprobó que el percentil 80 no era el mimso para los distintos tipos de inmueblo (en HY_tipo_simplif), por lo que se filtra el valor del percentil para las tres categorias diferentes de tipo inmuebles simplificado. Esto elimina nuevamente aproximadamente el 20% de los datos de training. Más adelante, se comprueba que estos filtros se pueden ajustar para producir mejores resultados de la siguiente forma: Se filtran del conjunto de training los inmuebles de tipo Garaje con un target superior a 75 segundos, eliminando el 25% de los garajes; se filtran los inmuebles de tipo Resto con un target superior a 133, eliminando el 15% y los inmuebles de tipo Vivienda con target superior a 101, eliminando el 20% de ellos (estos valores se obtuvieron con los cuantiles calculados anteriormente). Todo esto siguiendo la distribución de TARGET según HY_tipo_simplif. Esto se hizo así para que la distribución de la variable TARGET estuviera más concentrada para cada tipo de inmueble y así el modelo se ajustará más a los datos centrales, con lo que se obtiene un error mediano absoluto mejor.

Además, como se ha comentado anteriormente, también se eliminan del conjunto de training los datos que tengan TARGET=0, ya que son los que tienen GA_mean_bounce=100 o GA_exit_rate=100 y a estos se les asigna directamente el valor de 0 en la predición.

features<-c("GA_exit_rate",
            "GA_page_views_log",
            "GA_duracion",
            "IM_n",
            "IM_n_2",
            "HY_dif_precio_log",
            "HY_distribucion_words_log",
            "HY_metros_totales_log",
            "HY_precio_log",
            "HY_provincia_frec",
            "HY_cod_postal_frec",
            "IDEA_demand_sale_residential",
            "HY_descripcion_words_2",
            "IM_sim_mean",
            "HY_num_veces_rep",
            "IDEA_unitprice_sale_residential_log",
            "IDEA_unitprice_rent_residential",
            "IDEA_price_sale_residential",
            "IDEA_price_rent_residential",
            "HY_precio_metro_cuadrado",
            "IDEA_ind_turismo",
            "IDEA_pc_otros",
            "IDEA_pc_2000_10",
            "IDEA_area",
            "IDEA_pc_comercio"
            )


dmodel2$HY_tipo<-as.factor(dmodel2$HY_tipo)
dmodel2$HY_tipo_simplif<-as.factor(dmodel2$HY_tipo_simplif)
dmodel2$HY_tipo_simplif2<-as.factor(dmodel2$HY_tipo_simplif2)
dmodel2$HY_cert_energ_exist<-as.factor(dmodel2$HY_cert_energ_exist)
dmodel2$HY_descripcion_exist<-as.factor(dmodel2$HY_descripcion_exist)
dmodel2$HY_distribucion_exist<-as.factor(dmodel2$HY_distribucion_exist)

dataset<-as.data.frame(dmodel2%>%select(-HY_id,-HY_cod_postal, -HY_descripcion, -HY_distribucion,HY_cert_energ,-HY_tipo,-HY_provincia))
dataset<-createDummyFeatures(dataset)
a=Sys.time()
for (i in 616:620) {
print(i)
set.seed(i)

y=log(dmodel2$TARGET+1)
folds<-createFolds(y,k=5)
resultados<-c()
for(k in 1:5){
    testIndexes <-folds[[k]]
    testData <- dataset[testIndexes, ]
    trainData <- dataset[-testIndexes, ]
    # trainData<-quitarXSup(trainData,"TARGET",0.8)
    # trainData<-trainData%>%filter(TARGET!=0)
    #trainData<-quitarOutliers(trainData,"GA_mean_bounce")
    trainData<-trainData%>%filter(TARGET!=0,
                                   HY_tipo_simplif.Garajes==1&TARGET<75 | #86.7
                                   HY_tipo_simplif.Resto==1&TARGET<133 |
                                   HY_tipo_simplif.Vivienda==1&TARGET<101)
    #train_target<-trainData$TARGET_log
    train_target<-log(trainData$TARGET)
    
    
    model<-Cubist::cubist(x=trainData%>%select(features,starts_with("HY_tipo_simplif"),starts_with("GA_years")),
                          y=train_target, 
                          committees=100, 
                          control=Cubist::cubistControl(unbiased=TRUE, rules=100, seed=11,extrapolation=100) )
    
    #Se ve el funcionamiento en los datos de train
    pred=predict(model, trainData%>%select(features,starts_with("HY_tipo_simplif"),starts_with("GA_years")))
    perf=measureMEDAE(trainData$TARGET,exp(pred))
    #perf=measureMEDAE(trainData$TARGET,pred)

    #Se ve cómo se ajusta a los datos de test
    test_set_model<-testData%>%select(features,starts_with("HY_tipo_simplif"),starts_with("GA_years"))
    pred=predict(model,newdata=test_set_model)
    pred=exp(pred)
    perf_test=measureMEDAE(testData$TARGET,pred)
    pred_corregida=pred
    for (i in 1:nrow(testData)){
      pred_corregida[i]<-ifelse((testData$GA_exit_rate[i]==100|testData$GA_page_views[i]==0),0,pred_corregida[i])
    }
    perf_corregido=measureMEDAE(testData$TARGET,pred_corregida)
    resultados<-append(resultados,perf_corregido)
    print(paste('perf train',perf,'perf test',perf_test,'perf test corregida', perf_corregido,'fold', k))
}
print(paste("El medae medio en test corregido es:",mean(resultados)))
}
[1] 616
[1] "perf train 11.288105554235 perf test 16.2884973493037 perf test corregida 16.2884973493037 fold 1"
[1] "perf train 11.4109419888093 perf test 16.0092195780928 perf test corregida 16.0092195780928 fold 2"
[1] "perf train 11.5567554664554 perf test 15.5589222223542 perf test corregida 15.547237705653 fold 3"
[1] "perf train 11.2936240115384 perf test 16.6897511707747 perf test corregida 16.6421315174515 fold 4"
[1] "perf train 11.5134970808001 perf test 15.3864604719718 perf test corregida 15.3619795833265 fold 5"
[1] "El medae medio en test corregido es: 15.9698131467655"
[1] 617
[1] "perf train 11.3562836624596 perf test 16.9529812298567 perf test corregida 16.9466110780554 fold 1"
[1] "perf train 11.6113660925459 perf test 15.3269083288688 perf test corregida 15.3269083288688 fold 2"
[1] "perf train 11.4811375891383 perf test 15.8074874022873 perf test corregida 15.7809438445369 fold 3"
[1] "perf train 11.3656512265716 perf test 16.1703473813458 perf test corregida 16.1703473813458 fold 4"
[1] "perf train 11.5444290966408 perf test 16.5849986618242 perf test corregida 16.5820694574029 fold 5"
[1] "El medae medio en test corregido es: 16.161376018042"
[1] 618
[1] "perf train 11.31595534416 perf test 16.2003502254183 perf test corregida 16.1667036355369 fold 1"
[1] "perf train 11.4682369621512 perf test 16.5040551270734 perf test corregida 16.5040551270734 fold 2"
[1] "perf train 11.5194547927907 perf test 15.864064443461 perf test corregida 15.864064443461 fold 3"
[1] "perf train 11.4220120881968 perf test 16.0934092295536 perf test corregida 15.9856686238449 fold 4"
[1] "perf train 11.4892889098532 perf test 15.8874897107846 perf test corregida 15.8874897107846 fold 5"
[1] "El medae medio en test corregido es: 16.0815963081402"
[1] 619
[1] "perf train 11.5233509726897 perf test 15.6154689667447 perf test corregida 15.6154689667447 fold 1"
[1] "perf train 11.4055621060394 perf test 15.6595771173747 perf test corregida 15.601399545587 fold 2"
[1] "perf train 11.5684715669727 perf test 16.5959202418992 perf test corregida 16.59067276558 fold 3"
[1] "perf train 11.2944341197542 perf test 16.6486594255529 perf test corregida 16.6242084787683 fold 4"
[1] "perf train 11.2804109653015 perf test 16.2409418054047 perf test corregida 16.2409418054047 fold 5"
[1] "El medae medio en test corregido es: 16.1345383124169"
[1] 620
[1] "perf train 11.5651667947233 perf test 16.4229305683908 perf test corregida 16.4229305683908 fold 1"
[1] "perf train 11.3488129975313 perf test 16.2363453056877 perf test corregida 16.2170947052449 fold 2"
[1] "perf train 11.3149674034122 perf test 16.7607517798815 perf test corregida 16.7607517798815 fold 3"
[1] "perf train 11.4120880479901 perf test 15.9627772638626 perf test corregida 15.8634987887803 fold 4"
[1] "perf train 11.6000435687965 perf test 15.5585664968895 perf test corregida 15.5494911990796 fold 5"
[1] "El medae medio en test corregido es: 16.1627534082754"
Sys.time()-a
Time difference of 18.71244 mins
ggplot()+geom_point(aes(x=pred_corregida,y=testData$TARGET),color="coral")+geom_line(aes(x=c(0,150),y=c(0,150)),color="green")
ggplot()+geom_violin(aes(x=0,y=pred_corregida))

Ahora que ya se está en un margen de error aceptable, con un error mínimo de 15.95 en CV y un error medio en 5 semillas de 16.06, se procede a entrenar el modelo definitivo con todos los datos, apartir del cual se obtendrá las predicciones para destim. Con esto se busca construir un modelo entrenado con un mayor número de datos, por lo que tendrá más información para enfrentarse a los datos a predecir.

set.seed(616)
#Se elige el dataset que se va a usar
dataset_model<-select(dataset,features,TARGET,starts_with("HY_tipo_simplif"),starts_with("GA_years"))

dataset_model<-dataset_model%>%filter(TARGET!=0,
                                   (HY_tipo_simplif.Garajes==1&TARGET<75) |
                                   (HY_tipo_simplif.Resto==1&TARGET<133) |
                                   (HY_tipo_simplif.Vivienda==1&TARGET<101))

dataset_model$TARGET_log<-log(dataset_model$TARGET)



model<-Cubist::cubist(x=dataset_model%>%select(-TARGET,-TARGET_log),y=dataset_model$TARGET_log, committees=100,
                  control=Cubist::cubistControl(unbiased=TRUE, rules=100, seed=11) )

#Se ve el funcionamiento en los datos de train
pred=predict(model, dataset_model)
pred=exp(pred)
perf=measureMEDAE(dataset_model$TARGET,pred)
print(paste('Performance en train:',perf))
[1] "Performance en train: 11.4384843206968"
#Se ve el funcionamiento en los datos de train sin filtrar
pred=predict(model, dataset)
pred=exp(pred)
perf=measureMEDAE(dataset$TARGET,pred)
print(paste('Performance en train sin filtrar:',perf))
[1] "Performance en train sin filtrar: 15.206168152184"
#Se ve el funcionamiento tras corregir las predicciones
pred_corregida=pred
for (i in 1:nrow(dataset)){
  pred_corregida[i]<-ifelse((dataset$GA_exit_rate[i]==100|dataset$GA_page_views[i]==0),0,pred_corregida[i])
}
perf=measureMEDAE(dataset$TARGET,pred_corregida)
print(paste('Performance en train sin filtrar,tras correccion:',perf))
[1] "Performance en train sin filtrar,tras correccion: 15.1900071017578"
summary(pred_corregida)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   0.00   22.50   43.27   40.87   57.11  132.93 
ggplot()+geom_point(aes(x=pred_corregida,y=dataset$TARGET),color="coral")+geom_line(aes(x=c(0,150),y=c(0,150)),color="green")+ggtitle("Diferencias entre TARGET predicho y el real")
ggplot()+geom_violin(aes(x=0,y=pred_corregida))+ggtitle("Distribución del TARGET predicho en train")

Se muestra la estructura del modelo

#summary(model)

Evaluation on training data (7693 cases):

Average  |error|          0.5321641
Relative |error|               0.77
Correlation coefficient        0.65


Attribute usage:
  Conds  Model

   97%    91%    GA_page_views_log
   49%    74%    GA_exit_rate
   20%    66%    HY_cod_postal_frec
   15%    59%    GA_duracion
   11%    61%    HY_provincia_frec
   10%    54%    HY_metros_totales_log
    7%    21%    IDEA_price_rent_residential
    5%    48%    HY_precio_log
    4%    33%    HY_descripcion_words_2
    4%    16%    IDEA_demand_sale_residential
    4%    16%    IDEA_pc_otros
    4%    18%    IM_sim_mean
    4%    52%    IDEA_unitprice_sale_residential_log
    4%    16%    IDEA_unitprice_rent_residential
    4%    24%    HY_precio_metro_cuadrado
    3%    14%    IDEA_ind_turismo
    3%    33%    HY_num_veces_rep
    2%    36%    HY_tipo_simplif.Resto
    2%    35%    IDEA_pc_comercio
    2%    69%    HY_tipo_simplif.Garajes
    2%    23%    IDEA_area
    2%    16%    IDEA_pc_2000_10
    1%    62%    IM_n_2
    1%    60%    IM_n
    1%    24%    HY_dif_precio_log
          26%    IDEA_price_sale_residential
          36%    HY_distribucion_words_log
          24%    HY_tipo_simplif.Vivienda
          10%    HY_tipo_simplif2.Locales
           4%    HY_tipo_simplif2.Garajes

En esta última parte podemos ver las variables que ha utilizado para realizar las reglas de decisión y las que ha usado para hacerla regresión lineal. Se puede ver que GA_page_views_log es una variable muy usada, junto con GA_exit_rate, GA_duracion, HY_cod_postal_frec y HY_provincia_frec. Las tres primeras nos dan información sobre el inmueble en la web y cómo interaccionan los usuarios con él, y las últimas nos aporta información de la oferta que hay en la zona donde se ubica.

En la etapa de regresión toman importancia, de las variables IDEA, el precio de alquiler y el de venta de inmuebles residenciales de la zona. También es importante el número de imágenes que hay y la similitud media que hay entre ellas, al igual que la cantidad de palabras de la descripción y la distribución. En general, se prefieren las variables HY_tipo_simplif que las HY_tipo_simplif2, que son las disponibles en la página web.

4.5 Obtención de resultados para destim

Producción de resultados en destim.

destim2$HY_provincia<-as.factor(destim2$HY_provincia)
destim2$HY_tipo<-as.factor(destim2$HY_tipo)
destim2$HY_tipo_simplif<-as.factor(destim2$HY_tipo_simplif)
destim2$HY_tipo_simplif2<-as.factor(destim2$HY_tipo_simplif2)
destim2$HY_cert_energ_exist<-as.factor(destim2$HY_cert_energ_exist)
destim2$HY_descripcion_exist<-as.factor(destim2$HY_descripcion_exist)
destim2$HY_distribucion_exist<-as.factor(destim2$HY_distribucion_exist)

dataset<-as.data.frame(destim2%>%select(-HY_id,-HY_cod_postal, -HY_descripcion, -HY_distribucion,HY_cert_energ,-HY_tipo))
dataset<-createDummyFeatures(dataset)

destim_model<-select(destim2,features,starts_with("GA_years"),starts_with("HY_tipo_simplif"))
destim_model<-createDummyFeatures(destim_model)
pred=predict(model,newdata=destim_model)
pred_corregida=exp(pred)

for (i in 1:nrow(destim_model)){
  pred_corregida[i]<-ifelse((destim_model$GA_exit_rate[i]==100|destim_model$GA_page_views[i]==0),0,pred_corregida[i])
}

ggplot()+geom_violin(aes(x=0,y=pred_corregida))+ggtitle("Distribución de las predicciones finales")
output<-as.data.frame(cbind(destim2$HY_id,pred_corregida))
names(output)<-c("HY_id","TM_Est")
write.table(output,file="Mongo.txt",dec=".",sep = "|",row.names=F)

summary(output$TM_Est)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   0.00   23.25   42.24   40.86   56.77   98.88 

Estas predicciones son similares a las de entregas anteriores.

5 Conclusiones

En este proyecto se ha invertido una gran parte de los recursos temporales en la realización de análisis exploratorio de datos (EDA), convirtiéndose en la gran apuesta, el cual ha permitido descubrir muchas relaciones entre los datos y TARGET. Se han creado un gran número de variables nuevas a partir de las iniciales, por lo que el mayor reto ha sido realizar una selección de características buena. Para conseguirla, aparte de apoyarse en los métodos Filter, se ha ido avanzando en base al conocimiento que se ha obtenido respecto al problema a través de las fases de la competición con el EDA.

Respecto a la elección y construcción de modelos, se probaron modelos ensemble y modelos separados por tipos. Sin embargo, el mejor modelo conseguido hasta la fecha ha sido el modelo simple que presentamos en este documento. El éxito del mismo seguramente se debe a su reducida complejidad, ya que trabaja primero con reglas decisiones y luego afina el resultado con modelos lineales.