Analisi esplorativa del Mercato Immobiliare del Texas

Nel dataset presentato possiamo trovare questi tipi di variabili:

Abbiamo come variabili del tempo: year, month, e months_inventory, che possiamo utilizzare per avere una rappresentazione più specifica nel tempo delle azioni immobiliari delle varie città.

Parto dal creare un summary delle variabili un pò più generali, quindi year, month e sales, a mio avviso.

library(ggplot2)
library(dplyr)
## 
## Caricamento pacchetto: 'dplyr'
## I seguenti oggetti sono mascherati da 'package:stats':
## 
##     filter, lag
## I seguenti oggetti sono mascherati da 'package:base':
## 
##     intersect, setdiff, setequal, union
options(scipen = 999)
df = read.csv("realestate_texas.csv")
year_sum = summary(df$year)
sales_month = summary(df$month)
sales_sum = summary(df$sales)
print("Here's every helpful summary:")
## [1] "Here's every helpful summary:"
print("Year")
## [1] "Year"
print(year_sum)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    2010    2011    2012    2012    2013    2014
print("------------------------------------------")
## [1] "------------------------------------------"
print("Month")
## [1] "Month"
print(sales_month)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1.00    3.75    6.50    6.50    9.25   12.00
print("------------------------------------------")
## [1] "------------------------------------------"
print("Sales")
## [1] "Sales"
print(sales_sum)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    79.0   127.0   175.5   192.3   247.0   423.0

In seguito definisco delle funzioni per calcolare automaticamente index di posizione, varianza, forma e gini. E poi, calcoliamo la varianza più alta:

position.idx = function(x){
  mean = mean(x)
  min = min(x)
  max = max(x)
  median = median(x)
  quant = quantile(x)
  geom_mean = prod(x) ** (1 / length(x))
  armonic_mean = 1 / (sum(1 / x) / length(x))
  
  return(c(mean, min, max, median, quant, geom_mean, armonic_mean))
}


variance.idx = function(x){
  iqr = IQR(x)
  va = var(x)
  ds = sd(x)
  cv = ds / mean(x) * 100
  return(c(iqr, va, ds, cv))
}

gini.idx = function(x){
  ni = table(x)
  fi = ni / length(x)
  fi2 = fi ** 2
  J = length(ni)
  
  gini = 1 - sum(fi2)
  gini.norm = gini / ((J - 1) / J)
  return(gini.norm)
}

forma.idx = function(x){
  mu = mean(x)
  sigma = sd(x)
  n = length(x)
  m3 = sum((x - mu) ** 3) / n
  asim.index = m3 / sigma ** 3
  
  m4 = sum((x - mu) ** 4) / n
  kurtosis.index = m4 / sigma ** 4 - 3
  
  return(c(asim.index, kurtosis.index))
}

# Position index
pos_sales = position.idx(df$sales)
pos_volume = position.idx(df$volume)
pos_listings = position.idx(df$listings)

# variance index
variance_sales = variance.idx(df$sales)
variance_volume = variance.idx(df$volume)
variance_listings = variance.idx(df$listings)

print("Sales Variance")
## [1] "Sales Variance"
print(variance_sales[2])
## [1] 6344.3
print("------------------------------------------")
## [1] "------------------------------------------"
print("Volume Variance")
## [1] "Volume Variance"
print(variance_volume[2]) # Varianza positiva maggiore
## [1] 277.2707
print("------------------------------------------")
## [1] "------------------------------------------"
print("Listings Variance")
## [1] "Listings Variance"
print(variance_listings[2])
## [1] 566569
# Forma Index
forma_sales = forma.idx(df$sales)
forma_volume = forma.idx(df$volume)
forma_median_price = forma.idx(df$median_price)
forma_listings = forma.idx(df$listings)

print("------------------------------------------")
## [1] "------------------------------------------"
print("Quantile Sales")
## [1] "Quantile Sales"
quantile(df$sales)
##    0%   25%   50%   75%  100% 
##  79.0 127.0 175.5 247.0 423.0
print("------------------------------------------")
## [1] "------------------------------------------"
print("Quantile Volume")
## [1] "Quantile Volume"
quantile(df$volume)
##      0%     25%     50%     75%    100% 
##  8.1660 17.6595 27.0625 40.8930 83.5470
print("------------------------------------------")
## [1] "------------------------------------------"
print("Quantile Months Inventory")
## [1] "Quantile Months Inventory"
quantile(df$months_inventory)
##    0%   25%   50%   75%  100% 
##  3.40  7.80  8.95 10.95 14.90
# Gini idx
gini_months_inv = gini.idx(df$months_inventory)

Ora cerchiamo di capire quale sia la variabile con la più alta variabilità. Per fortuna, non è difficile da capire, avendo già calcolato la varianza nella funzione variance.idx, e possiamo notare come la variabile del Listings sia la più alta.

Ora bisogna capire quale sia la distribuzione più asimmetrica, anche questo è facile da calcolare grazie alla funzione forma.idx.

print("Distribuzione Sales")
## [1] "Distribuzione Sales"
print(forma_sales[1])
## [1] 0.7136206
print("------------------------------------------")
## [1] "------------------------------------------"
print("Distribuzione Volume")
## [1] "Distribuzione Volume"
print(forma_volume[1]) #  Distribuzione asimmetrica positiva
## [1] 0.8792182
print("------------------------------------------")
## [1] "------------------------------------------"
print("Distribuzione Prezzo Mediano")
## [1] "Distribuzione Prezzo Mediano"
print(forma_median_price[1]) # Distr asimmetrica negativa
## [1] -0.3622768
print("------------------------------------------")
## [1] "------------------------------------------"
print("Distribuzione Listings")
## [1] "Distribuzione Listings"
print(forma_listings[1])
## [1] 0.6454431

Da questo capiamo che la variabile del volume delle case ha una distribuzione asimmetrica positiva, sicuramente derivante dal fatto che la media è più alta della mediana e della moda, in più, visto che il risultato è il più alto sopra lo zero tra le variabili calcolate, possiamo definire sia quella più variabile.

Al contrario notiamo che invece, il prezzo medio ha una distribuzione asimmetrica negativa, probabilmente risultante dal fatto che il prezzo medio più alto è nettamente superiore alla media e mediana.

In seguito ho deciso di fare una distribuzione di frequenza per la variabile median_price, visto che non ne vale la pena di calcolare gli indici di posizione.

df$median_classes = cut(df$median_price, breaks = 5)
ni_median = table(df$median_classes)
fi_median = table(df$median_classes) / dim(df)[1]
Ni_median = cumsum(ni_median)
Fi_median = Ni_median / dim(df)[1]
median_classes_table = cbind(ni_median, fi_median, Ni_median, Fi_median)
median_classes_table
##                     ni_median fi_median Ni_median Fi_median
## (7.37e+04,9.5e+04]         18 0.0750000        18 0.0750000
## (9.5e+04,1.16e+05]         40 0.1666667        58 0.2416667
## (1.16e+05,1.38e+05]        73 0.3041667       131 0.5458333
## (1.38e+05,1.59e+05]        84 0.3500000       215 0.8958333
## (1.59e+05,1.8e+05]         25 0.1041667       240 1.0000000

Passaggio 4

Per questo passaggio selezioniamo la variabile quantitativa sales, e creiamo come con median price, una distribuzione di frequenza, ed in più calcoliamo l’indice di gini per questa variabile.

df$sales_classes = cut(df$sales, breaks = 5)

ni = table(df$sales_classes)
fi = table(df$sales_classes) / dim(df)[1]
Ni = cumsum(ni)
Fi = Ni / dim(df)[1]

sales_classes_table = table(df$sales_classes)
sales_classes_table
## 
## (78.7,148]  (148,217]  (217,285]  (285,354]  (354,423] 
##         84         77         41         27         11
ggplot(data = df)+
  geom_bar(aes(x = sales_classes, fill = city),
           stats = "count",
           col = "black")
## Warning in geom_bar(aes(x = sales_classes, fill = city), stats = "count", :
## Ignoring unknown parameters: `stats`

gini_idx_sales = df %>% group_by(df$city) %>% summarise(gini_sales = gini.idx(df$sales))


print(gini_idx_sales)
## # A tibble: 4 × 2
##   `df$city`             gini_sales
##   <chr>                      <dbl>
## 1 Beaumont                   0.998
## 2 Bryan-College Station      0.998
## 3 Tyler                      0.998
## 4 Wichita Falls              0.998

Dal grafico notiamo come la variabile delle vendite della classe più bassa è molto più alta delle altre, sopratutto per la città di Beaumont mentre nelle ultime due classi, per le zone con più vendita, il range di città è diminuito a solamente Bryan-College Station e Tyler.

Possiamo notare finalmente che l’indice di gini della variabilità sales è per tutte le città uguale, ed è anche estremamente vicino all’1, questo vuol dire che mostra un indice di eterogeneità massima e quindi equidistribuzione tra tutte le città.

Calcolo delle Probabilità

Per fare questo calcolo abbiamo bisogno di solamente una formula, ovvero il numero del numero di osservazioni per città ed il numero totale di osservazioni, quindi, 60 osservazioni moltiplicato su 4 città che risulta 240 osservazioni totali (oppure, più semplicemente, ottenere la length del dataset).

Ottenuti questi due dobbiamo semplicemente dividerli, e visto che le 4 città hanno lo stesso numero di osservazioni in questo dataset, possiamo concludere che la città Beaumount ha un 25% di probabilità di venir pescata randomicamente, così come tutte le altre città.

sample_city = sample(df$city, 1000000, replace = TRUE)

ggplot()+
  geom_histogram(aes(x=sample_city,
                     y = after_stat(count/sum(count))),
                 stat = "count",
                 col = "black",
                 fill = c("#A04248", "#754447", "#4B3839", "#332D2E"))+
  theme_light()+
  ylab("Probability")+
  xlab("City")+
  scale_y_continuous(breaks = seq(0, 0.5, 0.01))
## Warning in geom_histogram(aes(x = sample_city, y =
## after_stat(count/sum(count))), : Ignoring unknown parameters: `binwidth`,
## `bins`, and `pad`

table = table(df$city)
n = length(df$city)
probs = table[1] / n
print("Probability of it being \"Beaumont\" is: 0.25, or 25%")
## [1] "Probability of it being \"Beaumont\" is: 0.25, or 25%"

Ora rifacciamo lo stesso calcolo per capire quanta probabilità ci sia che peschiamo il mese di luglio:

sample_month = sample(df$month, 1000000, replace = TRUE)
ggplot()+
  geom_histogram(aes(x=sample_month,
                     y = after_stat(count/sum(count))),
                 stat = "count",
                 col = "black",
                 fill = "#754447")+
  labs(x = df$month)+
  theme_light()+
  ylab("Probability")+
  xlab("Month")+
  scale_y_continuous(breaks = seq(0, 0.5, 0.01))+
  scale_x_continuous(breaks = seq(0, 12, 1))
## Warning in geom_histogram(aes(x = sample_month, y =
## after_stat(count/sum(count))), : Ignoring unknown parameters: `binwidth`,
## `bins`, and `pad`

table_month = table(df$month)
n_month = length(df$month)
probs_month = table_month[1] / n_month

print("Probability of it being July is: 0.08333333 or 8.3%")
## [1] "Probability of it being July is: 0.08333333 or 8.3%"

Possiamo notare che in questo caso i vari mesi non hanno esattamente la stessa probabilità di uscire da una pescata randomica, specificatamente nel mese di luglio risulta che abbiamo un 8.3% di riuscita.

Infine, dobbiamo calcolare la probabilità che esca il mese di dicembre specificatamente del 2012:

sample_year = sample(df$year, 1000000, replace = TRUE)
ggplot()+
  geom_histogram(aes(x=sample_year,
                     y = after_stat(count/sum(count))),
                 stat = "count",
                 col = "black",
                 fill = "#754447")+
  labs(x = df$month)+
  theme_light()+
  ylab("Probability")+
  xlab("Year")+
  scale_y_continuous(breaks = seq(0, 0.5, 0.01))
## Warning in geom_histogram(aes(x = sample_year, y =
## after_stat(count/sum(count))), : Ignoring unknown parameters: `binwidth`,
## `bins`, and `pad`

table_year = table(df$year)
n_year = length(df$year)
probs_year = table_year[1] / n_year
probs = probs_month * probs_year
print("Probability of it being December of 2012 is 0.016 or 1.66%")
## [1] "Probability of it being December of 2012 is 0.016 or 1.66%"

Per calcolare la probabilità in questo caso, basta calcolare come abbiamo già visto la probabilità singola per l’anno, ed in seguito moltiplicarla per la probabilità del mese, così da avere una combinazione delle due.

Creazione di variabili Prezzo Medio e dell’efficacia degli annunci

Queste due operazioni possono essere interessanti per capire molto meglio secondo quali criteri specifici una persona sia tendente ad acquistare una specifica casa in una zona specifica, dipendentemente dal prezzo, dalla location, e sicuramente dal volume della casa.

Per poter calcolare efficacemente il prezzo medio possiamo semplicemente partire dalla mediana, e nonostante sia una distribuzione normale, possiamo essere specifici e calcolare il logaritmo della mediana, calcolare da esso la deviazione standard, moltiplicarla per 0.5, elevare il valore risultate per 2 ed infine moltiplicarlo per il prezzo mediano.

In questo caso notiamo che non c’è molta differenza tra la media e la mediana, questo perchè, i dati del dataset utilizzato tendono ad una distribuzione normale.

Ora che abbiamo la media possiamo creare una nuova colonna del dataset.

mean = df$median_price * exp(0.5 * sd(log(df$median_price))**2)
df$mean_price = mean

Per il calcolo dell’efficienza è un pò più facile, dovendo semplicemente dividere le vendite per il numero totale degli annunci, e moltiplicarlo per 100 per ottenere una percentuale.

df$listings
##   [1] 1533 1586 1689 1708 1771 1803 1857 1830 1829 1779 1742 1646 1677 1691 1762
##  [16] 1767 1832 1845 1822 1789 1785 1722 1687 1596 1647 1666 1730 1735 1765 1724
##  [31] 1749 1683 1704 1671 1652 1570 1581 1620 1567 1696 1659 1675 1708 1675 1681
##  [46] 1655 1624 1534 1575 1636 1539 1604 1620 1672 1657 1617 1501 1575 1544 1500
##  [61] 1298 1439 1577 1984 1613 1588 1646 1599 1549 1530 1510 1416 1480 1562 1723
##  [76] 1833 1840 1758 1662 1581 1536 1492 1445 1362 1486 1634 1806 1834 1793 1734
##  [91] 1669 1518 1498 1458 1442 1442 1486 1599 1750 1680 1581 1462 1385 1385 1201
## [106] 1155 1132 1057 1199 1218 1261 1271 1212 1152 1041 1016 1022 1031  973  882
## [121] 2727 2763 2729 3014 3175 3294 3272 3267 3296 3156 3042 2878 2852 2938 3101
## [136] 3196 3266 3256 3263 3218 3094 3056 2876 2720 2811 2857 2868 2940 2981 3041
## [151] 3072 3042 2953 2897 2830 2633 2658 2666 2788 2920 2946 2986 2998 2953 2917
## [166] 2852 2701 2500 2609 2625 2737 2778 2744 2855 2875 2791 2696 2602 2460 2272
## [181]  908  915  946  904  914  972  993 1022 1028 1005  968  938  955  950  968
## [196]  996 1052 1030 1029 1004 1005  963  902  844  859  861  887  907  914  934
## [211]  941  933  931  907  877  801  854  850  868  895  900  923  844  830  812
## [226]  796  777  743  746  774  838  852  899  961  941  973  940  905  870  821
df$sales
##   [1]  83 108 182 200 202 189 164 174 124 150 150 148 108 108 146 166 143 177
##  [19] 163 160 127 155 124 151 110 135 175 176 197 173 182 218 182 193 162 160
##  [37] 159 140 170 198 246 232 206 273 213 188 177 212 148 186 182 208 246 254
##  [55] 212 262 224 260 180 202  89 107 176 233 282 286 189 196 122 100 101 130
##  [73]  94 101 162 186 238 284 244 196 131 114 124 135 115 124 149 193 294 293
##  [91] 292 296 149 163 159 134 141 125 187 267 341 357 402 328 186 164 166 190
## [109] 152 171 275 303 353 377 403 298 204 218 169 200 160 181 250 316 282 262
## [127] 255 238 220 202 155 209 143 181 241 253 271 313 278 295 258 227 208 198
## [145] 169 225 271 253 289 302 322 317 281 289 210 234 197 227 298 289 326 335
## [163] 369 357 287 272 253 239 238 244 282 323 388 423 371 347 361 369 300 332
## [181]  89  91 147 167 165 129 104 130 132 113  97 117  80  79 119 111 128 111
## [199] 127 135 113  98  93  81 105  90 116 125 102 130 132 123  95  97 119 115
## [217]  79 101 159 124 144 121 150 149 128 114  94  92  89  93 102 143 140 123
## [235] 150 137 110 112  96 109
efficiency = (df$sales / df$listings) * 100
df$efficiency = efficiency

ggplot(data = df)+
  geom_density(aes(x = efficiency, col = city))+
  theme_minimal()+
  xlab("Efficiency in %")+
  ylab("Density")

Analisi condizionata

Utilizziamo la libreria dplyr per concatenare la summarise della media delle vendite e del volume alla raggruppazione in base alla città.

df %>% 
  group_by(city,year,month) %>% 
  summarise(mean_sales = mean(sales), 
            mean_volume = mean(volume))
## `summarise()` has grouped output by 'city', 'year'. You can override using the
## `.groups` argument.
## # A tibble: 240 × 5
## # Groups:   city, year [20]
##    city      year month mean_sales mean_volume
##    <chr>    <int> <int>      <dbl>       <dbl>
##  1 Beaumont  2010     1         83        14.2
##  2 Beaumont  2010     2        108        17.7
##  3 Beaumont  2010     3        182        28.7
##  4 Beaumont  2010     4        200        26.8
##  5 Beaumont  2010     5        202        28.8
##  6 Beaumont  2010     6        189        27.2
##  7 Beaumont  2010     7        164        22.7
##  8 Beaumont  2010     8        174        25.2
##  9 Beaumont  2010     9        124        17.2
## 10 Beaumont  2010    10        150        23.9
## # ℹ 230 more rows

Visualizzazioni con ggplot2

Per quest’ultimo passaggio ho deciso di rappresentare i quattro plot a mio parere, più significativi:

ggplot(data = df)+
  geom_line(aes(x = month, y = sales, color = city),
                 stat = "identity")+
  scale_x_continuous(breaks = seq(0, 12, 1))+
  ylab("Sales")+
  xlab("Month")+
  facet_wrap(~ year)+
  theme_minimal()

Nel primo abbiamo una rappresentazione divisa per anno che calcola le vendite ogni mese a seconda della città. Così possiamo capire quali tipi di abitazioni sono state più vendute in quale periodo, notando che nel 2014 c’è stato il numero di vendite più alto generalmente, e che tendenzialmente in tutti gli anni, abbiamo il picco di vendite che si aggirano dal mese di Maggio fino a Luglio.

ggplot(data = df)+
  geom_histogram(aes(x = month, y = volume, fill = city),
                 stat = "identity")+
  scale_x_continuous(breaks = seq(0, 12, 1))+
  ylab("Volume")+
  xlab("Month")+
  facet_wrap(~ year)+
  theme_minimal()
## Warning in geom_histogram(aes(x = month, y = volume, fill = city), stat =
## "identity"): Ignoring unknown parameters: `binwidth`, `bins`, and `pad`

In questo secondo grafico abbiamo mantenuto la divisione tra anni e città, con lunghezza temporale divisa in mesi, però calcolando in quali mesi ci sono state le abitazioni con i volume più alti.

Mettendo a confronto il primo ed il secondo grafico, notiamo che le vendite più alte sono state nei periodi dove ci sono stati gli appartamenti volume più alto.

ggplot(data = df)+
  geom_boxplot(aes(x = median_price, fill = city))+
  xlab("Prezzo mediano")+
  theme_minimal()

In quest’ultimo grafico abbiamo un boxplot diviso per città basandosi sul prezzo medio, notiamo che il prezzo medio è più alto sulla città di Bryan-College Station, però si ritrova sul secondo quartile, mentre Wichita Falls si ritrova sull’ultimo quartile in alto, però con la media dei prezzi più bassa.

ggplot(data = df)+
  geom_histogram(aes(x = month, y = sales, fill = city),
                 stat = "identity",
                 position = "stack")+
  scale_x_continuous(breaks = seq(0, 12, 1))+
  ylab("Sales")+
  xlab("Month")+
  facet_wrap(~ year)+
  theme_minimal()
## Warning in geom_histogram(aes(x = month, y = sales, fill = city), stat =
## "identity", : Ignoring unknown parameters: `binwidth`, `bins`, and `pad`

ggplot(data = df)+
  geom_histogram(aes(x = month, y = sales, fill = city),
                 stat = "identity",
                 position = "fill")+
  scale_x_continuous(breaks = seq(0, 12, 1))+
  ylab("Sales")+
  xlab("Month")+
  facet_wrap(~ year)+
  theme_minimal()
## Warning in geom_histogram(aes(x = month, y = sales, fill = city), stat =
## "identity", : Ignoring unknown parameters: `binwidth`, `bins`, and `pad`

Qui possiamo vedere il confronto tra le vendite totale per mese in divisione tra le città nel corso dei vari anni in un grafico normalizzato in confronto alla controparte non normalizzata.

Conclusione

Per concludere, in questa analisi del dataset, abbiamo scoperto che in media, le persone hanno preferito comprare immobili più grandi, nonostante il prezzo, nei periodi centrali dell’anno, ovvero tra Maggio e Luglio, e specificatamente che siano nelle città di Beaumont e Bryan College Station.

Abbiamo scoperto che gli annunci che hanno portato a più vendite si possono ritrovare nelle città di Tyler e Whichita Falls, che riportano la percentuale più alta di vendite effettuate nel totale.

Infine, scopriamo con i dati in questo dataset che se dovessimo cercare un immobile in un periodo randomico dell’anno, avremmo una distribuzione di possibilità molto equilibrata in generale, per ogni città.