Assignment 4
#Research Question: How do variables such as poverty, education (or lack thereof), and hardship have an effect on the number of crimes committed in a certain area, specifically Chicago in 2015?
library(sf) ## Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUElibrary(dplyr) ##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(spData) ## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`library(ggplot2)
library(ggthemes)
library(spdep)
library(spatialreg)## Loading required package: Matrix##
## Attaching package: 'spatialreg'## The following objects are masked from 'package:spdep':
##
## get.ClusterOption, get.coresOption, get.mcOption,
## get.VerboseOption, get.ZeroPolicyOption, set.ClusterOption,
## set.coresOption, set.mcOption, set.VerboseOption,
## set.ZeroPolicyOptionlibrary(GWmodel)## Loading required package: robustbase## Loading required package: sp## Loading required package: Rcpp## Welcome to GWmodel version 2.3-2.library(tidyr)##
## Attaching package: 'tidyr'## The following objects are masked from 'package:Matrix':
##
## expand, pack, unpacksetwd("/Users/victoriamalaspina/Documents/DIDA 370/chicago_airbnb")
chi_data <- st_read("airbnb_Chicago 2015.shp")## Reading layer `airbnb_Chicago 2015' from data source
## `/Users/victoriamalaspina/Documents/DIDA 370/chicago_airbnb/airbnb_Chicago 2015.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 77 features and 20 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -87.94011 ymin: 41.64454 xmax: -87.52414 ymax: 42.02304
## Geodetic CRS: WGS 84ggplot()+
geom_sf(data = chi_data, aes(fill=as.numeric(poverty)))+
scale_fill_steps(
name = "Poverty",
low = "lightsteelblue1",
high = "tomato1",
n.breaks = 5,
show.limits = T)+
theme_void()
#Variables w/test
model <- chi_data %>%
mutate(unemp = unemployed/population,
educ = without_hs/population,
har = harship_in/population,
p_crime = num_crimes/population * 100) %>%
select(c(poverty, unemp, educ, har, p_crime, geometry)) %>%
filter(!poverty %in% "NA") %>%
mutate(poverty = as.numeric(poverty))
summary(lm1 <- lm(poverty ~ unemp + educ + har, p_crime, data = model))##
## Call:
## lm(formula = poverty ~ unemp + educ + har, data = model, subset = p_crime)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.1267 -2.5660 -0.0969 4.6320 17.2217
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.583 1.089 20.747 < 2e-16 ***
## unemp -37176.334 3780.551 -9.834 5.18e-15 ***
## educ -37535.979 3801.333 -9.874 4.35e-15 ***
## har 26339.074 2184.978 12.055 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.539 on 73 degrees of freedom
## Multiple R-squared: 0.7411, Adjusted R-squared: 0.7305
## F-statistic: 69.65 on 3 and 73 DF, p-value: < 2.2e-16#spatautocorr w resid w moran I
chi_list <- model %>%
poly2nb(st_geometry(model)) %>%
nb2listw(zero.policy = TRUE)
lm.morantest(lm1, chi_list)##
## Global Moran I for regression residuals
##
## data:
## model: lm(formula = poverty ~ unemp + educ + har, data = model, subset
## = p_crime)
## weights: chi_list
##
## Moran I statistic standard deviate = 0.9194, p-value = 0.1789
## alternative hypothesis: greater
## sample estimates:
## Observed Moran I Expectation Variance
## 0.050022038 -0.015875033 0.005137208#The Moran I test above indicates that there is not a significance in the spatial autocorrelation of the residuals. The p-value of the variables is measured at 0.1789, which is not a significant value, which is typically measured around 0.05.
lm_error <- errorsarlm(poverty ~ unemp + educ + har + p_crime,
data = model,
listw = chi_list,
zero.policy = TRUE,
na.action = na.omit)
summary(lm_error)##
## Call:errorsarlm(formula = poverty ~ unemp + educ + har + p_crime,
## data = model, listw = chi_list, na.action = na.omit, zero.policy = TRUE)
##
## Residuals:
## Min 1Q Median 3Q Max
## -17.79143 -3.56286 -0.68924 3.53254 23.29828
##
## Type: error
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 9.8014e+00 2.3331e+00 4.2009 2.658e-05
## unemp -3.9772e+03 1.7774e+03 -2.2376 0.02524
## educ 7.9326e+01 2.7520e+03 0.0288 0.97700
## har 1.7535e+03 9.4369e+02 1.8581 0.06316
## p_crime 4.7688e-01 7.4344e-02 6.4144 1.414e-10
##
## Lambda: 0.55255, LR test value: 11.612, p-value: 0.0006554
## Asymptotic standard error: 0.11806
## z-value: 4.6803, p-value: 2.8644e-06
## Wald statistic: 21.905, p-value: 2.8644e-06
##
## Log likelihood: -251.402 for error model
## ML residual variance (sigma squared): 37.252, (sigma: 6.1035)
## Number of observations: 77
## Number of parameters estimated: 7
## AIC: 516.8, (AIC for lm: 526.42)lm_lag <- lagsarlm(poverty ~ unemp + educ + har + p_crime,
data = model,
listw = chi_list,
zero.policy = TRUE,
na.action = na.omit)
summary(lm_lag)##
## Call:lagsarlm(formula = poverty ~ unemp + educ + har + p_crime, data = model,
## listw = chi_list, na.action = na.omit, zero.policy = TRUE)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.83602 -4.28968 -0.76901 3.25398 22.34944
##
## Type: lag
## Coefficients: (asymptotic standard errors)
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.0476e+00 2.0941e+00 0.9778 0.32818
## unemp -3.3381e+03 1.8713e+03 -1.7839 0.07445
## educ 1.7033e+03 2.4398e+03 0.6981 0.48509
## har 1.1788e+03 9.5267e+02 1.2374 0.21595
## p_crime 4.0962e-01 6.7541e-02 6.0647 1.322e-09
##
## Rho: 0.4061, LR test value: 14.366, p-value: 0.00015047
## Asymptotic standard error: 0.09733
## z-value: 4.1724, p-value: 3.0141e-05
## Wald statistic: 17.409, p-value: 3.0141e-05
##
## Log likelihood: -250.0247 for lag model
## ML residual variance (sigma squared): 37.301, (sigma: 6.1075)
## Number of observations: 77
## Number of parameters estimated: 7
## AIC: 514.05, (AIC for lm: 526.42)
## LM test for residual autocorrelation
## test value: 0.31283, p-value: 0.57595#Error model and lag models both have insignificant p-values. Additionally, their AIC values are similar, though still not much better than the AIC for the lagrange multiplier.
seResiduals <- rep(0, length(model$poverty))
resIndex <- lm_error$residuals %>% names() %>% as.integer();
seResiduals[resIndex] <- lm_error$residuals
chi_list %>%
moran.test(seResiduals, ., zero.policy = TRUE)##
## Moran I test under randomisation
##
## data: seResiduals
## weights: .
##
## Moran I statistic standard deviate = 0.33901, p-value = 0.3673
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.010792044 -0.013157895 0.004991022seResiduals <- rep(0, length(model$poverty))
resIndex <- lm_lag$residuals %>% names() %>% as.integer();
seResiduals[resIndex] <- lm_lag$residuals
chi_list %>%
moran.test(seResiduals, ., zero.policy = TRUE) ##
## Moran I test under randomisation
##
## data: seResiduals
## weights: .
##
## Moran I statistic standard deviate = 0.60985, p-value = 0.271
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.030190758 -0.013157895 0.005052419#Again, these tests show similar p-values which are not very significant, and the lag model is only very slightly more significant than the error model.
LM <- lm.LMtests(lm1, chi_list, test = "all")## Please update scripts to use lm.RStests in place of lm.LMtestsLM##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = poverty ~ unemp + educ + har, data = model, subset
## = p_crime)
## test weights: listw
##
## RSerr = 0.45207, df = 1, p-value = 0.5014
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = poverty ~ unemp + educ + har, data = model, subset
## = p_crime)
## test weights: listw
##
## RSlag = 1.1491, df = 1, p-value = 0.2837
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = poverty ~ unemp + educ + har, data = model, subset
## = p_crime)
## test weights: listw
##
## adjRSerr = 0.0018964, df = 1, p-value = 0.9653
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = poverty ~ unemp + educ + har, data = model, subset
## = p_crime)
## test weights: listw
##
## adjRSlag = 0.69894, df = 1, p-value = 0.4031
##
##
## Rao's score (a.k.a Lagrange multiplier) diagnostics for spatial
## dependence
##
## data:
## model: lm(formula = poverty ~ unemp + educ + har, data = model, subset
## = p_crime)
## test weights: listw
##
## SARMA = 1.151, df = 2, p-value = 0.5624#Based on these scores, it is clear to see that while they are not very significant, the lag model shows the most significant spatial dependence compared to the other models. Based on these results, I would conclude that there is not much spatial dependence on these variables owing to the p-values shown within the data. While some are slightly more significant than others, they all reflect scores that are not very significant thus implying that the variables (poverty, unemployment, education, and hardship) are not as correlated as I previously thought.