Tipo ~ Precio

En primer lugar estudiaremos gráficamente como afecta el tipo de inmueble en el precio del mismo.

ggplot(datos,aes(Tipo,Precio))+
  geom_boxplot(aes(fill=Tipo))

Como vemos, el hecho de añadir Chalets en nuestro análisis descriptivo daña la visualización del precio. Como evidenciamos y además es obvio, el precio de los chalets en la ciudad de Elche es muy superior a cualquier otro tipo de vivienda, alcanzando estos inmuebles unos precios excesivamente altos. También observamos como existe un piso cuyo valos supera el 1.5M de Euros. Para poder obtener otra interpretación de esta relación, eliminaremos los chalets y los inmuebles cuyo valor supera el 1.0M de Euros.

datos=datos[-which(datos$Tipo=="Chalet"),]
datos=datos[-which(datos$Precio>=1000000),]

ggplot(datos,aes(Tipo,Precio))+
  geom_boxplot(aes(fill=Tipo))+
  scale_y_continuous(labels = label_comma())

Ahora vemos un gráfico mucho más uniforme. En este gráfico vemos qu los precios más bajos los encontramos en los estuidos y las plantas bajas. Seguramente los estudios tengan esos precios tan bajos por el número de metros cuadrados, habitaciones y aseos. Las plantas bajas por otra parte tienen un precio más bajo por el hecho de estar a pie de calle y ser menos valorado que pisos en más altura. Entre todos los inmuebles (excepto los Chalets evidentemente), lo más caro es una casa, entendiendo esta como un bungalow o un adosado.

Barrio ~ Precio

Estudiaremos en segundo lugar la relación del precio con el barrio en el que se encuentra.

ggplot(datos,aes(Barrio,Precio))+
  geom_boxplot(aes(fill=Barrio))+
  scale_y_continuous(labels = label_comma())+
  theme(axis.text.x=element_text(angle=90))

Encontramos por un lado obviedades que podemos intuir antes del análisis por lógica y es que los barrios donde el precio de la vivienda es más bajo es en aquellos donde el nivel socio económico es más bajo: Carrús, Palmeral y algunas zonas de Altabix. Vemos observaciones irrelevantes para nuestro estudio, como puede ser UMH, donde se ha ofrecido es ubicación con un solo inmueble, podríamos pensar en dejar El Palmleral, pero para ofrecer más uniformidad vamos a eliminarlo.

datos=datos[-which(datos$Barrio=="UMH"),]
datos=datos[-which(datos$Barrio=="El Palmeral"),]


ggplot(datos,aes(Barrio,Precio))+
  geom_boxplot(aes(fill=Barrio))+
  scale_y_continuous(labels = label_comma())+
  theme(axis.text.x=element_text(angle=90))

Aunque pueda parecer que el barrio más caro es Nou Altabix, yo considero el Centro como el más caro, porque aunque el precio de la vivienda sea generalmente más bajo, los pisos buenos son más caros que en Nou Altabix, llegando algunos hasta los 750.000€. Ahora bien, también es cierto que en la zona Centro es la más variable de todas, pues a pesar de tener pisos excesivamente caros, también hay algunos realmente bajos, llegando inclusio al precio más bajo que encontramos en barrios como Carrús o Altabix.

También vemos como el Raval es una zona realmente característica, pues el precio medio de la vivienda es prácticamente el percentil 25, pero observamos que hay muchos inmuebles que superan los 500.000€ con mucha facilidad. Esto indica que existen inmuebles realmente “lujosos” en esa zona, pero que la mayoría de ellos tienen un precio muy normal.

Piscina/Garaje ~ Precio

A continuación, veremos como afecta el hecho de tener piscina o garaje en las viviendas.

piscina=ggplot(datos,aes(Piscina,Precio))+
  geom_boxplot(aes(fill=Piscina))+
  scale_y_continuous(labels = label_comma())

garaje=garaje=ggplot(datos,aes(Garaje,Precio))+
  geom_boxplot(aes(fill=Garaje))+
  scale_y_continuous(labels = label_comma())

grid.arrange(piscina,garaje,nrow=1)

Viendo ambos gráficos, vemos claramente como el hecho de tener piscina es mucho más relevante para el precio. Este dato hay que cogerlo con pinzas, porque no estamos etudiando aisladamente la influencia de la piscina, por lo que aplicando lógica con el apartado anterior vemos lo siguiente. Las viviendas que disponen de piscina son aquellas que se construyen de obra nueva y en las zonas más actuales de Elche: Nou Altabix, Travalón y Ciudad Universitaria. Estas zonas son las más caras, por lo que el hecho de tener pisicna hace que suba el precio, no solo por tener piscina, sino por la ubicación de dichas viviendas entre otros factores.

Metros cuadrados ~ Precio (concicionado al Tipo)

Para estudiar la relación del precio con los metros cuadrados del inmuebles hemos añadido la variable del tipo de inmueble que vamos a tratar. Hemos hecho esto, porque si estudiamos aisladamente los metros cuadrados con el precio no vamos a obtener ningún tipo de conclusión, pues no es lo mismo tener 150m2 en el Pla, en Carrús, en Ciudad Jardín o en el Centro.

Para esta relación hemos realizadeo un gráfico de dispersión animado, donde podemos seleccionar y desseleccionar el tipo de inmueble para la visualización clickando sobre la leyenda.

a=ggplot(datos, aes(Metros2,Precio,color=Tipo))+
  geom_point()+
  stat_smooth(method="lm",
              formula=y~x,
              geom="smooth")+
  theme_bw()
  
ggplotly(a)

Las relaciones más importantes las obtenemos en los dúplex, donde el precio de la vivienda aumenta considerablemente con el aumento de metros cuadrados. Siendo sin ninguna duda la vivienda con más metros la más cara.También encontramos esta relación en las casas hasta los 380m2. A partir de ahí no existe relación.

Como vemos, no existen relaciones aparentemente claras entre el tipo de inmueble y el precio de la vivienda, por lo que vamos a a probar ahora en función del barrio en el que se encuentra.

Habitaciones/Aseos ~ Precio (Condicionado al barrio)

De nuevo podemos seleccionar y deseleccionar el barrio clickando en la leyenda

b=ggplot(datos, aes(Metros2,Precio,color=Barrio))+
  geom_point()+
  stat_smooth(method="lm",
              formula=y~x,
              geom="smooth")+
  theme_bw()
  
ggplotly(b)

Los barrios en los que sí que aumenta el precio de la vivienda con el aumento de los metros cuadrados son: Altabix, Centro, Ciudad Jardín, El Raval, Pont Nou. Esto sucede por la cotización que existe en estas zonas. Altabix, Centro, Raval y Pont Nou son zonas con una muy buena localización y barrios antiguos donde todo está ya construido, por lo que el metro cuadrado está muy cotizado, y los inmuebles con muchos metros valen mucho dinero. Ciudad Jardín por otra parte son casas con terrenos, y son parte de una zona lujosa de Elche, por lo que de nuevo el precio de cada metro es muy caro.

Explicación algoritmo

Para la estimación de la proporción hemos cogido un número n de observaciones al azar de la base de datos que compondrá el conjunto de entrenamiento (train). Dentro de este train vamos a aplicar la estimación de la proporción y el respectivo intervalo de confianza, pero para hacerlo más significativo vamos a proceder a aplicar Monte Carlo. Dentro de este conjunto de entrenamiento tomamos n muestras aleatorias de tamaño 0.1*n y almacenamos los datos relativos al límite superior, inferior y central en tres bases de datos distintas. Tras acabar el bucle aplicamos la media de los n límites superiores, inferiores y centrales obtenidos y ya tenemos la estimación de la proporción.

Para probar cómo de bien funciona respecto a los demás datos, cogemos todos los datos que no han sido utilizados en el conjunto de entrenamiento (train), a los cuales vamos a llamar test. Comparamos ahora cuál es la proporción de días en los que se produce lluvia en el periodo test y del periodo train y vemos cómo de fiable es nuestra estimación.

Estimación de la proporción de días de lluvia en Elche

Para los días de lluvia en la ciudad de Elche, hemos utilizado un mínimo de 0.2mm de lluvia al día para considerarse lluvia. Además hemos escogido un tamaño n igual a 1000 observaciones para el conjunto de entrenamiento y un alpha de 0.05 para hacer una estimación del intervalo de confianza de 95%.

set.seed(5)
estimacion_prop = function(criterio,n,alpha){
  datos$riesgo=ifelse(datos$lluvia>=criterio,1,0)
  train1=sample(nrow(datos),n)
  inferior=c()
  superior=c()
  media=c()
  for (i in 1:n){
    azar=sample(train1,n*0.1)
    train=datos[azar,]
    
    N=nrow(datos)
    n2=n*0.1
    y=length(which(train$riesgo==1))
    
    prop=y/n2
    q=1-prop
    
    z=qnorm(1-alpha/2, mean=0,sd=1,lower.tail=T)
    
    semi=z*sqrt((prop*q/(n2-1))*((N-n2)/N)+1/(2*n2));semi
    lower=prop-semi
    upper=prop+semi
    
    inferior[i]=lower
    superior[i]=upper
    media[i]=prop

  }
  vacios1=c(seq(1:nrow(datos)))
  vacios=data.frame(vacios1)
  vacios$x=rep(0)
  vacios$x[intersect(vacios1,train1)]=1
  test=datos[which(vacios$x==0),]
  prop_real=length(which(test$riesgo==1))/nrow(test)*100
cat("La estimación de la proporción de lluvia ofrece un intervalo de confianza
definido por:\n Limite Inferior:",mean(inferior)*100,"%\n Estimacion Central:",mean(media)*100,"%\n Limite superior:",mean(superior)*100,"%")

cat("\n\nPara este ejercicio, hemos utilizado",n,"datos para el conjunto train, por lo que hemos dejado el resto sin utilizar. De los datos no utilizados para el problema obtenemos una proporción de lluvia de:",prop_real)
}


criterio=0.2
n=1000
alpha=0.05
estimacion_prop(criterio,n,alpha)

La estimación de la proporción de lluvia ofrece un intervalo de confianza
definido por:
 Limite Inferior: 0.5019036 %
 Estimacion Central: 16.037 %
 Limite superior: 31.5721 %

Para este ejercicio, hemos utilizado 1000 datos para el conjunto train, por lo que hemos dejado el resto sin utilizar. De los datos no utilizados para el problema obtenemos una proporción de lluvia de: 16.45455

Del resultado obtenido podemos obtener las siguientes conclusiones. En primer lugar, que el tamaño muestral del conjunto train es suficientemente bueno, ya que la predicción muestral es muy ajustada a lo que sucede en el conjunto test. También me gustaría destacar de estos resultados que cuando el límite inferior es menor a 0, este no tiene sentido para este tipo de datos, pues es lo mismo -7.43% (por ejemplo) que 0%, indicando ambos que el límite inferior es que no llueve.

Para estas estimaciones utilizamos la aleatoriedad con la función “sample”; me gustaría destacar que esta aleatoriedad se puede controlar con la función set.seed(), es decir, fijando una semilla de aleatoriedad. En mi caso he elegido la semilla 5, que ofrece los valores que se pueden observar en el resultado.

Si multiplicamos la proporción por el número de días que tiene el año, es decir 16.037/100*365, obtendremos el número de días que llueve al año en la ciudad de Elche, que es 58.53, según el criterio podríamos decir que es o 58 o 59, con un intervalo al 95% que tiene un límite inferior de 2 días y un límite superior de 115 días.

Si hacemos una rápida búsqueda en internet, llegaremos a la siguiente página https://www.meteoelche.com/hysnoaa.php?select=2021.TXT&type=anual, donde podemos ver el número de días que llovió en Elche el año 2021. Este dato es exactamente 59 días, por lo que podemos decir que hemos obtenido una estimación prácticamente perfecta.

Estimación de la proporción de días de lluvia en Elche con toda la base de datos

Cambiando los criterios del tamaño del conjunto n podemos dar con la estimación del conjunto completo de la base de datos, para ello n = longitud de la base de datos.

set.seed(5)
estimacion_prop = function(criterio,n,alpha){
  datos$riesgo=ifelse(datos$lluvia>=criterio,1,0)
  train1=sample(nrow(datos),n)
  inferior=c()
  superior=c()
  media=c()
  for (i in 1:n){
    azar=sample(train1,n*0.1)
    train=datos[azar,]
    
    N=nrow(datos)
    n2=n*0.1
    y=length(which(train$riesgo==1))
    
    prop=y/n2
    q=1-prop
    
    z=qnorm(1-alpha/2, mean=0,sd=1,lower.tail=T)
    
    semi=z*sqrt((prop*q/(n2-1))*((N-n2)/N)+1/(2*n2));semi
    lower=prop-semi
    upper=prop+semi
    
    inferior[i]=lower
    superior[i]=upper
    media[i]=prop

  }
  vacios1=c(seq(1:nrow(datos)))
  vacios=data.frame(vacios1)
  vacios$x=rep(0)
  vacios$x[intersect(vacios1,train1)]=1
  test=datos[which(vacios$x==0),]
  prop_real=length(which(test$riesgo==1))/nrow(test)*100
cat("La estimación de la proporción de lluvia ofrece un intervalo de confianza
definido por:\n Limite Inferior:",mean(inferior)*100,"%\n Estimacion Central:",mean(media)*100,"%\n Limite superior:",mean(superior)*100,"%")
}


criterio=0.2
n=nrow(datos)
alpha=0.05
estimacion_prop(criterio,n,alpha)

La estimación de la proporción de lluvia ofrece un intervalo de confianza
definido por:
 Limite Inferior: 5.560401 %
 Estimacion Central: 16.23016 %
 Limite superior: 26.89992 %

El conjunto de la base de datos nos ofrece una estimación del 16.23016%, que si la multiplicamos por 365 días de nuevo, obtenemos 59.24 días. De nuevo vemos como esta apreciación vuelve a ser prácticamente perfecta con el año 2021.

Si nos fijamos en el intervalo de confianza vemos como con un tamaño muestral más grande ajustamos mejor el intervalo de confianza y ahora es: - Límite inferior: 20 días - Límite superior: 98 días

Seguramente 2021 haya sido una gran coincidencia, y es que la estimación del clima puede llevar implícitos algunos factores que se escapan de las estimaciones muestrales y poblacionales. Por ejemplo este año hemos sufrido un momento de sequía total, y hemos visto llover tan solo 30 días este año en la ciudad. Aún con unas precipitaciones extremadamente bajas, el intervalo de confianza ha garantizado las precipitaciones para las que estaba dispuesto.

Tras observar todos los años desde 2013 (año a partir del cual dejamos de tener datos en nuestra base de datos) hasta 2023 (último año con registros), se observa como todos los años el número de días lluviosos por año ha entrado sin ningún tipo de problemas en el intervalo de confianza.

No solo obtenemos resultados interesantes cambiando el tamaño de la muestra del periodo train, sino también cambiando el criterio. Si ahora queremos calcular cuántos días al año no se asistirá a clase porque se decretará estado naranja de alerta climática, tan solo tenemos que ajustar las precipitaciones a 40 mm.

set.seed(3)
estimacion_prop = function(criterio,n,alpha){
  datos$riesgo=ifelse(datos$lluvia>=criterio,1,0)
  train1=sample(nrow(datos),n)
  inferior=c()
  superior=c()
  media=c()
  for (i in 1:n){
    azar=sample(train1,n*0.1)
    train=datos[azar,]
    
    N=nrow(datos)
    n2=n*0.1
    y=length(which(train$riesgo==1))
    
    prop=y/n2
    q=1-prop
    
    z=qnorm(1-alpha/2, mean=0,sd=1,lower.tail=T)
    
    semi=z*sqrt((prop*q/(n2-1))*((N-n2)/N)+1/(2*n2));semi
    lower=prop-semi
    upper=prop+semi
    
    inferior[i]=lower
    superior[i]=upper
    media[i]=prop

  }
  vacios1=c(seq(1:nrow(datos)))
  vacios=data.frame(vacios1)
  vacios$x=rep(0)
  vacios$x[intersect(vacios1,train1)]=1
  test=datos[which(vacios$x==0),]
  prop_real=length(which(test$riesgo==1))/nrow(test)*100
cat("La estimación de la proporción de lluvia ofrece un intervalo de confianza
definido por:\n Limite Inferior:",mean(inferior)*100,"%\n Estimacion Central:",mean(media)*100,"%\n Limite superior:",mean(superior)*100,"%")

cat("\n\nPara este ejercicio, hemos utilizado",n,"datos para el conjunto train, por lo que hemos dejado el resto sin utilizar. De los datos no utilizados para el problema obtenemos una proporción de lluvia de:",prop_real)
}


criterio=40
n=1000
alpha=0.05
estimacion_prop(criterio,n,alpha)

La estimación de la proporción de lluvia ofrece un intervalo de confianza
definido por:
 Limite Inferior: -13.77652 %
 Estimacion Central: 0.095 %
 Limite superior: 13.96652 %

Para este ejercicio, hemos utilizado 1000 datos para el conjunto train, por lo que hemos dejado el resto sin utilizar. De los datos no utilizados para el problema obtenemos una proporción de lluvia de: 0.09090909

Como vemos, se estima una proporción de 0.095% de los días. Esto suma un total de aproximadamente 1 día al año, con un límite inferior de 0 días y un límite superior de 50 (exageradamente alto a simple vista).

Si nos vamos a la emeroteca y buscamos en google “suspensión clases Elche” y entramos en el apartado de noticias, veremos que tan solo el 22 de Mayo de 2023 (1 día en todo el año) se suspendieron las clases, por lo que de nuevo se ha ajustado perfectamente la estimación.

Estimación del tamaño de muestra óptimo para el estudio realizado

Una de las principales disyuntivas es, cómo de grande debe ser la muestra que tomemos para que los resultados que obtengamos sean representativos. Este problema es realmente importante, pues si el coste de obtener cada observación es muy alto, conocer de antemano el número de observaciones que tenemos que estudiar supondrá una gran ventaja para la entidad que requiera los resultados. Para ello aplicamos la siguiente función:

estim_muestra=function(alpha,d,n){
set.seed(2)
datos$riesgo=ifelse(datos$lluvia>=0.1,1,0)
train1=sample(nrow(datos),n)
train=datos[train1,]
  
N=nrow(datos)
y=length(which(train$riesgo==1));y
  
prop=y/N
q=1-prop

z=qnorm(1-alpha/2, mean=0, sd=1, lower.tail=T)
num=(prop*q*z**2+d**2)*N;num
den=(d**2*N+prop*q*z**2);den
nopt=ceiling(num/den)
return(nopt)
}
alpha=0.05
d=0.02
n=1000
estim_muestra(alpha,d,n)

[1] 563

Para nuestro caso, donde partimos de un tamaño poblacional de 2100 observaciones, y suponemos tener una muestra de tan solo 1000, vemos que para no cometer un error mayor a 2 puntos porcentuales, con una muestra de 563 observaciones sería suficiente para garantizar un 95% de confianza.