Projekt Ekonometria przestrzenna
library(htmltools)
library(sp)
library(spdep)
library(sf)
library(lmtest)
library(whitestrap)
library(spatialreg)
library(tmap)
library(tseries)
library(dplyr)
library(readxl)
library(knitr)Wstęp
W ostatnich latach rynek nieruchomości w Polsce przechodzi dynamiczne zmiany, co stanowi istotny temat badań zarówno dla ekonomistów, urbanistów, jak i potencjalnych inwestorów. Ceny mieszkań są kluczowym wskaźnikiem gospodarczym, wpływającym na dostępność lokali mieszkalnych oraz na decyzje zakupowe konsumentów, a także decyzje inwestorów. W obszarach silnie zurbanizowanych oraz w ich pobliżu, ceny mieszkań są zazwyczaj wyższe. W niniejszym projekcie badawczym skupiamy się na analizie przestrzennej zmienności cen za 1 metr kwadratowy mieszkań w Polsce, z podziałem na poszczególne powiaty. Ceny mieszkań nie kształtują się równomiernie na terenie całego kraju.
Celem projektu jest zidentyfikowanie i zrozumienie czynników wpływających na różnice w cenach mieszkań w różnych częściach kraju. Zrozumienie tych różnic ma kluczowe znaczenie dla kształtowania polityki mieszkaniowej, planowania przestrzennego oraz dla deweloperów i inwestorów szukających optymalnych lokalizacji dla swoich przedsięwzięć.
Generowanie danych (jedynie przy pierwszym włączeniu kodu)
1.Wczytanie danych (jeżeli dane zostały już wcześniej wygenerowane)
Opis zmiennych
Zmienna objaśniana
CENA_1M2 - Średnia cena za 1 m2 lokali mieszkalnych sprzedanych w ramach transakcji rynkowych (zł)
Zmienne objaśniające
BEZR_REJ - Bezrobotni zarejestrowani pozostający bez pracy dłużej niż
1 rok [osoby]
ZAS_MIESZK - Zasoby mieszkaniowe
MALZENSTWA - Małżeństwa na 1000 ludności wg lokalizacji
UR_REJ - Stopa bezrobocia rejestrowanego
L_TRANS - Liczba transakcji kupna/sprzedaży lokali mieszkalnych
MIESZK_ODD - Mieszkania oddane do użytkowania
POZWOLENIA - Pozwolenia wydane na budowę mieszkań i zgłoszenia budowy z
projektem budowlanym
WYNAGR - Przeciętne miesięczne wynagrodzenia brutto
WSK_DROGI - Drogi - wskaźniki [km/100km2]
L_ZGON - Liczba zgonów
TER_ZIEL - Tereny zieleni - liczba parków spacerowo-wypoczynkowych
L_PRZYST - Liczba czynnych przystanków autobusowych i tramwajowych
GEST_ZAL - Gęstość zaludnienia (liczba osób na 1 km ^2)
2.Dodanie sąsiedztwa I rzędu
3.Utworzenie mapy sąsiedztwa
coords <- coordinates(as(powiaty, "Spatial"))
plot(st_geometry(powiaty), border="grey", col = "light blue",asp=1.4)
plot(queen1, coords, add=TRUE)
title("Mapa sąsiedztwa powiatów w Polsce")4.Utworzenie heatmapy
tmap_options(check.and.fix = TRUE)
tm_shape(powiaty) +
tm_fill("CENA_1M", title = "Cena za metr kwadratowy", palette = "YlOrRd",breaks = c(0,2000,3000,4000,5000,6000,7000,Inf), legend.show = TRUE) +
tm_layout(main.title = "Mapa cieplna cen za metr mieszkań")5. Test Morana - Statystyki lokalne i globalne.
Badamy stopień intensywności danej cechy w obiektach przestrzennych.
##
## Moran I test under randomisation
##
## data: powiaty$CENA_1M
## weights: Wqueen1
##
## Moran I statistic standard deviate = 9.987, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.337563495 -0.002638522 0.001160397
## Ii E.Ii Var.Ii Z.Ii
## Min. :-2.26662 Min. :-1.447e-02 Min. :0.00000 Min. :-2.7168
## 1st Qu.:-0.01517 1st Qu.:-3.335e-03 1st Qu.:0.01853 1st Qu.:-0.1091
## Median : 0.11300 Median :-1.143e-03 Median :0.08764 Median : 0.6517
## Mean : 0.33756 Mean :-2.639e-03 Mean :0.31999 Mean : 0.6007
## 3rd Qu.: 0.46025 3rd Qu.:-2.658e-04 3rd Qu.:0.30376 3rd Qu.: 1.3415
## Max. : 5.02424 Max. :-1.000e-09 Max. :4.97258 Max. : 4.9706
## Pr(z != E(Ii))
## Min. :0.0000007
## 1st Qu.:0.1415508
## Median :0.3526931
## Mean :0.4052572
## 3rd Qu.:0.6506720
## Max. :0.9988542
I. ćwiartka -> zależności High-High, III. Ćwiartka -> zależności Low-Low Statystyka morana = 0.338 P-value=2.2e-16 < alfa = 0,05, więc odrzucamy H0 na korzyść hipotezy alternatywnej, mówiącej o występowaniu zależności przestrzennych wśród cen mieszkań pomiędzy powiatami.
6. Dodanie logarytmów i utworzenie modelu OLS - liniowego modelu potęgowego
Model ze zmiennymi: Gestość zaludnienia, wynagrodzenia, zasoby
mieszkaniowe(liczba mieszkań) Eliminacji zmiennych dokonaliśmy na
podstawie prób i testów modeli w programie GRETL.
powiaty <- powiaty %>%
mutate(WYNAGR = ifelse(JPT_KOD_J == "2467", 5345.45, WYNAGR))
OLS <- lm(log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) + log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data=powiaty)
OLS##
## Call:
## lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) + log(powiaty$WYNAGR) +
## log(powiaty$ZAS_MIE), data = powiaty)
##
## Coefficients:
## (Intercept) log(powiaty$GEST_ZA) log(powiaty$WYNAGR)
## 4.38485 0.03895 0.30340
## log(powiaty$ZAS_MIE)
## 0.12006
7. Test Morana - Wybór modelu regresji przestrzennej
##
## Global Moran I for regression residuals
##
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## weights: Wqueen1
##
## Moran I statistic standard deviate = 8.0868, p-value = 3.061e-16
## alternative hypothesis: greater
## sample estimates:
## Observed Moran I Expectation Variance
## 0.269991597 -0.004434827 0.001151578
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## test weights: listw
##
## RSerr = 62.006, df = 1, p-value = 3.442e-15
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## test weights: listw
##
## RSlag = 53.905, df = 1, p-value = 2.104e-13
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## test weights: listw
##
## adjRSerr = 8.9288, df = 1, p-value = 0.002807
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## test weights: listw
##
## adjRSlag = 0.82762, df = 1, p-value = 0.363
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## test weights: listw
##
## SARMA = 62.834, df = 2, p-value = 2.265e-14
columbus.lagrange <- lm.LMtests(OLS, listw=Wqueen1, test=c("LMerr","RLMerr","LMlag","RLMlag","SARMA"))
summary(columbus.lagrange)## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
## data:
## model: lm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty)
## test weights: listw
##
## statistic parameter p.value
## RSerr 62.00621 1 3.442e-15 ***
## adjRSerr 8.92881 1 0.002807 **
## RSlag 53.90502 1 2.104e-13 ***
## adjRSlag 0.82762 1 0.362962
## SARMA 62.83383 2 2.265e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Test Morana - H0: brak zależności przestrzennej p-value blieskie 0
< alfa -> występuje zależność przestrzenna
I para:
RSerr - dotyczy modelu SEM ;
H0 - wybieramy OLS. Dodanie przestrzennego błędu nic nie wnosi ;
H1: wybieramy SEM, wprowadzamy reg. przestrzenna <-
RSlag - dotyczy modelu SAR ;
H0 - wybieramy OLS. Dodanie przestrzennego błędu nic nie wnosi;
H1: wybieramy SAR, wprowadzamy reg. przestrzenna W obu odrzucamy H0,
<-
żeby wybrać model, patrzymy na drugą parę.
II para:
adjRSerr - p-value<alfa, odrzucamy H0
adjRSlag - p-value>alfa, brak podstaw do odrzucenia H0.
Na podstawie drugiej miary wybieramy model SEM.
SARMA - połaczenie modeli, bardzo małe p-value to jest najlepsza opcja.
8.1. Durbin - model wspólnego czynnika
DURBIN <- lagsarlm(log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) + log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data=powiaty, listw=Wqueen1, type="mixed");
summary(DURBIN)##
## Call:lagsarlm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty,
## listw = Wqueen1, type = "mixed")
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4260211 -0.1028220 0.0050804 0.0959232 0.6214819
##
## Type: mixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.256654 1.278964 0.9826 0.325826
## log(powiaty$GEST_ZA) 0.057553 0.010088 5.7051 1.163e-08
## log(powiaty$WYNAGR) 0.207686 0.101205 2.0521 0.040157
## log(powiaty$ZAS_MIE) 0.088895 0.017999 4.9389 7.855e-07
## lag.log(powiaty$GEST_ZA) -0.088120 0.018034 -4.8864 1.027e-06
## lag.log(powiaty$WYNAGR) 0.016740 0.168098 0.0996 0.920671
## lag.log(powiaty$ZAS_MIE) 0.093449 0.035966 2.5983 0.009369
##
## Rho: 0.41411, LR test value: 43.477, p-value: 4.29e-11
## Asymptotic standard error: 0.062287
## z-value: 6.6484, p-value: 2.9623e-11
## Wald statistic: 44.202, p-value: 2.9623e-11
##
## Log likelihood: 144.1586 for mixed model
## ML residual variance (sigma squared): 0.026469, (sigma: 0.16269)
## Number of observations: 380
## Number of parameters estimated: 9
## AIC: -270.32, (AIC for lm: -228.84)
## LM test for residual autocorrelation
## test value: 13.595, p-value: 0.00022676
SCM2 <- lagsarlm(log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) + log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data=powiaty, listw=Wqueen1, Durbin = ~log(powiaty$ZAS_MIE)+
log(powiaty$GEST_ZA));
summary(SCM2)##
## Call:lagsarlm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty,
## listw = Wqueen1, Durbin = ~log(powiaty$ZAS_MIE) + log(powiaty$GEST_ZA))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4260160 -0.1022567 0.0055942 0.0957840 0.6223171
##
## Type: mixed
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.3534254 0.8323450 1.6260 0.103941
## log(powiaty$GEST_ZA) 0.0573568 0.0098622 5.8158 6.034e-09
## log(powiaty$WYNAGR) 0.2120192 0.0913373 2.3213 0.020272
## log(powiaty$ZAS_MIE) 0.0888738 0.0179977 4.9381 7.890e-07
## lag.log(powiaty$ZAS_MIE) 0.0935208 0.0359411 2.6021 0.009267
## lag.log(powiaty$GEST_ZA) -0.0874753 0.0168519 -5.1908 2.094e-07
##
## Rho: 0.41506, LR test value: 44.877, p-value: 2.098e-11
## Asymptotic standard error: 0.061656
## z-value: 6.7319, p-value: 1.6748e-11
## Wald statistic: 45.318, p-value: 1.6748e-11
##
## Log likelihood: 144.1537 for mixed model
## ML residual variance (sigma squared): 0.026465, (sigma: 0.16268)
## Number of observations: 380
## Number of parameters estimated: 8
## AIC: -272.31, (AIC for lm: -229.43)
## LM test for residual autocorrelation
## test value: 10.785, p-value: 0.001023
8.2. SARMA - model mieszany SAR+SMA
SARMA<- sacsarlm(log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) + log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data=powiaty, listw=Wqueen1)
summary(SARMA)##
## Call:sacsarlm(formula = log(powiaty$CENA_1M) ~ log(powiaty$GEST_ZA) +
## log(powiaty$WYNAGR) + log(powiaty$ZAS_MIE), data = powiaty,
## listw = Wqueen1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4167510 -0.0942371 0.0034049 0.0909543 0.5358800
##
## Type: sac
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 7.9420341 1.3074681 6.0744 1.245e-09
## log(powiaty$GEST_ZA) 0.0657181 0.0092633 7.0945 1.298e-12
## log(powiaty$WYNAGR) 0.2690093 0.0960455 2.8009 0.005097
## log(powiaty$ZAS_MIE) 0.0733567 0.0167245 4.3862 1.154e-05
##
## Rho: -0.34495
## Asymptotic standard error: 0.12021
## z-value: -2.8696, p-value: 0.00411
## Lambda: 0.7139
## Asymptotic standard error: 0.070467
## z-value: 10.131, p-value: < 2.22e-16
##
## LR test value: 58.888, p-value: 1.6309e-13
##
## Log likelihood: 138.4168 for sac model
## ML residual variance (sigma squared): 0.024461, (sigma: 0.1564)
## Number of observations: 380
## Number of parameters estimated: 7
## AIC: -262.83, (AIC for lm: -207.95)
9. Decyzja
Model DURBINA ma najmniejsze kryt. informacyjne AKAIKE i dlatego wybieramy go, pomimo, że model SARMA też ma dobre statystyki.
10. Model Durbina - interpretacje
Współczynniki:
log(powiaty$GEST_ZA): 0.057 (p < 0.0001)
Wzrost gęstości zaludnienia o 1% powoduje wzrost ceny mieszkań o około 0.057%. Ten współczynnik jest dodatni i bardzo istotny, co wskazuje na silny wpływ gestacji(zagęsczenia ludności) na ceny mieszkań.
log(powiaty$WYNAGR): 0.212 (p = 0.020)
Wzrost wynagrodzeń o 1% prowadzi do wzrostu ceny mieszkań o około 0.212%. Współczynnik ten jest pozytywny i istotny, co oznacza, że wyższe wynagrodzenia są związane z wyższymi cenami mieszkań.
log(powiaty$ZAS_MIE): 0.089 (p < 0.0001)
Wzrost zasobów mieszkaniowych o 1% powoduje wzrost ceny mieszkań o około 0.089%. Ten współczynnik jest dodatni i bardzo istotny, co wskazuje, że większa dostępność mieszkań (zasobów mieskzaniowych) przyczynia się do wzrostu cen mieszkań.
Współczynniki Opóźnień Przestrzennych.
lag.log(powiaty$ZAS_MIE): 0.094 (p = 0.009)
Przestrzenne opóźnienie zasobów mieszkaniowych o 1% powoduje wzrost ceny mieszkań o około 0.094%. Współczynnik ten jest dodatni i istotny, co sugeruje, że ceny mieszkań są pozytywnie skorelowane z liczbą zasobów mieszkaniowych w sąsiednich powiatach.
lag.log(powiaty$GEST_ZA): -0.087 (p < 0.0001)
Przestrzenne opóźnienie gęstości zaludnienia o 1% powoduje spadek ceny mieszkań o około 0.087%. Ten współczynnik jest negatywny i bardzo istotny, co wskazuje, że wyższe zagęszczenie ludności w sąsiednich powiatach ma negatywny wpływ na ceny mieszkań.