Census data - why is it important? - brief history? - what questions can be answered with census data? - why do we need to perform spatial analysis on census data - how variable are socioeconomic patterns in the Victoria Capital Regional District (CRD)?
Spatial Autocorrelation determines how different neighbours are from each other, by measuring how similar observations are to other observations that are close to each other. In geography, Tobler’s law states that closer things are generally more similar to each other than to things further away. - describe global vs local methods of analysis - how can these methods be applied to census data? - entice reader with this! - just give general background to see if this tutorial will be of any use
How to perform spatial analysis of census variables to assess their spatial autocorrelation.
This tutorial will describe how to use R to perform a spatial autocorrelation analysis on multiple census variables, by using both global and local methods to create outputs that will help to visualize the socioeconomic variability of the population in the CRD. Following these instructions, and using the R code chunks, should reproduce the output results that are included.
The median income variable (MdInc) from the Victoria Capital Regional District (CRD) 2016 Census Data will be analyzed to determine if either of them has any spatial autocorrelation with themselves, over space.
The following R code can be used as is to create a full spatial autocorrelation for one variable, and produce maps and tables as outputs. As one of the benefits of R is creating reproducable results, if another variable has been chosen to analyze, the same code can be reused, and run again after only making a few changes in the variable names, and the titles produced for output images files. Note for Different Variable: After each of the following code segments, a note will be added to identify the code segments needed to be changed as necessary for a different variable.
ToCII.1. Set the working directory to where the data will be located and the output images will be stored. This is how R will know how to the locate the data the code needs to run.
dir <- "/Working/Data"
setwd(dir)
getwd()
In order to use R to do spatial autocorrelation analysis, the following R packages are needed to run the code. Each package contains a bundle of open-source shareable code, documentation and data, each of which is a tool that performs different functions.
If any, or all, of these R packages have not been previously installed, the function install.packages("packageName") is used to install each separate package. This step will only need to be done once. If more than one package needs to be installed, this can be done with the c() function to combine all the packages into a list. Even though packages are already on your computer, they may need updating, and R may prompt with a message to restart R prior to the install - restarting may, or may not be necessary.
install.packages(c("plyr", "dplyr", "spdep", "GISTools", "raster", "maptools", "rgdal", "tmap", "BAMMtools", "shinyjs"))
After the R packages have been installed, their Libraries need to be called in order for the R code to access the package-specific functions. The libraries will need to be recalled for this code if the R console has been closed.
ToClibrary("plyr")
library("dplyr")
library("spdep")
library("GISTools")
library("raster")
library("maptools")
library("rgdal")
library("tmap")
library("BAMMtools")
library("shinyjs")
IV.1. Download the most recent census shape and data files from Statistics Canada (2016), and save all the data files in the computer's working directory, as specified in section II.
IV.2. Read in the data from the census tract boundaries shapefile, using the rgdal package (don't include the .shp file extension), and save the data in a object, which is a spatial polygon dataframe, named tracts.
tracts <- readOGR(dsn = ".", layer = "Vic_Census")
## OGR data source with driver: ESRI Shapefile
## Source: "/Users/geographydepartment/Documents/Geog418/Assign3-Geog418/Data/Working/Data", layer: "Vic_Census"
## with 78 features
## It has 2 fields
tracts
## class : SpatialPolygonsDataFrame
## features : 78
## extent : 3903547, 3974156, 1910542, 1962451 (xmin, xmax, ymin, ymax)
## crs : +proj=lcc +lat_1=49 +lat_2=77 +lat_0=63.390675 +lon_0=-91.86666666666666 +x_0=6200000 +y_0=3000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
## variables : 2
## names : CTUID, GUID
## min values : 9350001.00, 9350001
## max values : 9350180.06, 9350180.06
IV.3. Read the data from the 2016 census data .csv file (include the .csv file extension), assigning the newly created data frame object to a variable named census.16.
census.16 <- read.csv("CensusTractData.csv")
IV.4. Summary(census.16) provides summary statistics of the attributes for each of the variable columns of the Census Tract Data, including the minimum, median, mean, maximum, and 1st and 3rd Quartiles. The summary also records the number of missing NA values in each column of data, which will be removed in step IV.6.
summary(census.16)
## some summary stats have been ommitted ...
## GUID CUID Pop PopCh
## Min. :9350000 Min. : 0.0 Min. : 94 Min. :-60.200
## 1st Qu.:9350102 1st Qu.:101.5 1st Qu.: 3479 1st Qu.: 1.725
## Median :9350128 Median :128.0 Median : 4689 Median : 4.100
## Mean :9350110 Mean :109.6 Mean : 9430 Mean : 4.599
## 3rd Qu.:9350154 3rd Qu.:154.0 3rd Qu.: 6156 3rd Qu.: 7.275
## Max. :9350180 Max. :180.1 Max. :367770 Max. : 44.200
## NA's :1 NA's :1
##
## Divorced Widowed MdInc MDAftTaxInc
## Min. : 0.0 Min. : 0.0 Min. :12384 Min. :12192
## 1st Qu.: 196.2 1st Qu.: 135.0 1st Qu.:34056 1st Qu.:30667
## Median : 302.5 Median : 215.0 Median :37481 Median :33380
## Mean : 671.4 Mean : 476.2 Mean :37557 Mean :33388
## 3rd Qu.: 473.8 3rd Qu.: 328.8 3rd Qu.:42240 3rd Qu.:37328
## Max. :26195.0 Max. :18565.0 Max. :57557 Max. :48939
## NA's :1 NA's :1 NA's :4 NA's :4
##
## MVisMin Owner Renter Band
## Min. : 0.00 Min. : 20.0 Min. : 0.0 Min. : 0.000
## 1st Qu.: 10.00 1st Qu.: 937.5 1st Qu.: 308.8 1st Qu.: 0.000
## Median : 20.00 Median : 1277.5 Median : 602.5 Median : 0.000
## Mean : 38.78 Mean : 2612.4 Mean : 1557.2 Mean : 2.821
## 3rd Qu.: 28.75 3rd Qu.: 1690.0 3rd Qu.: 1071.2 3rd Qu.: 0.000
## Max. :1505.00 Max. :101885.0 Max. :60745.0 Max. :105.000
## NA's :1 NA's :1 NA's :1 NA's :1
IV.5. Merge the two data frames census tract boundaries with the census data, using the common column GUID. Assign this merged data to the variable name crd.data.
crd.data <- merge(tracts, census.16, by = "GUID")
IV.6.i) Remove any missing 'na' values from the merged data, as there is no numerical value associated with this row of data, and thus no calculations can be done on any data that doesn't contain a value. The function is formatted to start with an ! explanation mark, which means 'not', thus !is.na() is used to specify which data in crd.data is not na. Once the dataset has the rows containing a na value removed, the remaining data is reassigned back to the crd.data data frame.
IV.6.ii) To specify a column to use from the dataset, a $ dollar sign is placed between the data set name and the column name.
crd.data <- crd.data[!is.na(crd.data$MdInc),]
IV.7.i) Notes for Different Variable: The na values will also be needed to removed from each different variable. To choose a different column from the census data, replace the column name after the $ crd/data$NewColumnName.
IV.7.ii) Notes for Different Variable: In order to compare the Renter variable between different polygons, the variable needs to be 'normalized', by dividing each count of Renters in a polygon by the population of that polygon, and creating a new column with the % of population that is a renter in any polygon. These values can now be compared across different polygons.
# crd.data$RenterNorm <- crd.data$Renter / crd.data$Pop
# class(crd.data$RenterNorm)
IV.8. To understand more about the data, use class() to see that crd.data is a "SpatialPolygonsDataFrame", and use summary() to produce a table of the data frame summary statistics, including the column names of each different variable.
class(crd.data)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"
summary(crd.data)
## Object of class SpatialPolygonsDataFrame
## Coordinates:
## min max
## x 3903547 3974156
## y 1910542 1962451
## Is projected: TRUE
## proj4string :
## [+proj=lcc +lat_1=49 +lat_2=77 +lat_0=63.390675
## +lon_0=-91.86666666666666 +x_0=6200000 +y_0=3000000 +ellps=GRS80
## +towgs84=0,0,0,0,0,0,0 +units=m +no_defs]
##
## Data attributes:
##
## some summary stats have been ommitted ...
## GUID CTUID CUID Pop
## Min. :9350001 9350001.00: 1 Min. : 1.0 Min. : 332
## 1st Qu.:9350101 9350002.00: 1 1st Qu.:101.2 1st Qu.: 3538
## Median :9350126 9350003.01: 1 Median :126.5 Median : 4744
## Mean :9350108 9350003.02: 1 Mean :108.4 Mean : 4964
## 3rd Qu.:9350153 9350004.00: 1 3rd Qu.:152.8 3rd Qu.: 6156
## Max. :9350180 9350005.00: 1 Max. :180.1 Max. :10393
## (Other) :68
##
## Separated Divorced Widowed MdInc
## Min. : 5.0 Min. : 10.0 Min. : 5.0 Min. :12384
## 1st Qu.: 60.0 1st Qu.: 200.0 1st Qu.:138.8 1st Qu.:33878
## Median : 95.0 Median : 310.0 Median :215.0 Median :37563
## Mean :117.0 Mean : 353.6 Mean :250.9 Mean :37558
## 3rd Qu.:158.8 3rd Qu.: 473.8 3rd Qu.:328.8 3rd Qu.:42250
## Max. :380.0 Max. :1065.0 Max. :670.0 Max. :57557
##
## Vminority MVisMin Owner Renter
## Min. : 0.000 Min. : 0.00 Min. : 65.0 Min. : 15.0
## 1st Qu.: 0.000 1st Qu.:10.00 1st Qu.: 966.2 1st Qu.: 325.0
## Median : 0.000 Median :20.00 Median :1287.5 Median : 620.0
## Mean : 5.608 Mean :20.54 Mean :1375.5 Mean : 820.3
## 3rd Qu.:10.000 3rd Qu.:28.75 3rd Qu.:1690.0 3rd Qu.:1071.2
## Max. :35.000 Max. :80.00 Max. :3255.0 Max. :3945.0
##
ToC
A map will be created using the tmap package.
V.1. Select a colour palette to create the chloropleth map using tmaptools::palette_explorer(), after first removing the comment #. In the following R code, the # hashtag symbol is used to make the line of code a 'comment', so that the code will not be run. Removing the # from the beginning of the line, will allow the code to be run, and leaving the # will ensure that this code won't be run if not necessary.
To generate the colour code needed for the map palette, check the Code Generator button for tmap layer function code, then copy/paste the code found under the colour palette. Note, the palette explorer window must be closed in order for R to continue processing the code, which is why it is initially commented out.
# tmaptools::palette_explorer()
V.2. The shape of the polygons comes from the merged crd.data file.
V.2.i) The polygons are classified by colour based on the values in the MdInc column (col="MdInc").
V.2.ii) The title is used for the classification legend.
V.2.iii) The classification style Jenks uses natural breaks to set up the best classification of variable values, based on the data itself. n=6 chooses 6 different classification breaks.
V.2.iv) Add a scale bar and compass in specific positions. The first value is RIGHT / LEFT and the second value is TOP/BOTTOM (note these positions must be capitalized).
V.2.v) The layout function adds a title to the map, in a specific postion. To make room for north arrow and scale bar, place inner.margins around map and image boundary (bottom, left, top, right).
V.3. Note for Different Variable: For the new variable, the crd.data will reflect the new data column as created in Step IV.7 after removing the na values from the new variable column. Change the following values tm_polygons(col = "NewColumnName", title= "New Title"); tm_layout(title= "New Title"). If comparing results from different census variables, each would need to use the same classification scheme and number of breaks.
map_PopDens <- tm_shape(crd.data) +
tm_polygons(col = "MdInc",
title = "Median $ Income",
style = "jenks",
palette = "viridis", n = 6) +
tm_scale_bar(width = 0.22, position = c("RIGHT", "BOTTOM")) +
tm_compass(position = c("RIGHT", "TOP")) +
tm_layout(title = "CRD 2016 Census Tracts Median $ Income", title.position = c("LEFT", "TOP"), inner.margins = c(.08, .03, .08, .03))
V.4. To print the plot in R window, run the variable name map_PopDens. In order to save the map, the file object png is called using a filename to save, call the plot, then close the graphic device connection with dev.off() in order to return back to the R code console.
map_PopDens
V.5. In order to to save the file, and then close the plot window the function dev.off() must be used.
V.5.i) Note for Different Variable: Change the following values name of png file to be saved "newName.png".
# create png of Median Income
png("map_PopDens_MdInc.png")
map_PopDens
dev.off()
V.5.ii) For Mac OSX users, 'quartz_off_screen' is a graphic device driver. Calling dev.off() again should return user to the RStudio graphic device or 'null device'.
## quartz_off_screen
## 2
V.6. Interpreting the output: The finished map shows the 2016 median income range for the CRD census tracts. The darkest blue colours correspond to the lowest median income range of $12,384 - $12,704, and the yellow colour the highest range of $47,872 - $57,557. Some polygons of similar income ranges appear to grouped closer together in the $37,931-$42,261 and $42,262-$47,872 ranges.
Neighbourhoods are defined by weighing schemes to determine which polygons with adjoining boundaries are considered neighbours. The Rook selects neighbouring polygons with shared boundarie in only 4 directions: top, bottom, and both sides, giving these neighbours a value of 1, and all other shared boundary neighbours a value of 0. The Queen Neighbour weighing scheme selects the same neighbours as the Rook scheme, as well as all the diagonally adjoining polygons, giving them all a value of 1. - describe output - what does it show?
VI.1. Using the spdep package,
VI.1.i) poly2nb builds a list of neighbours from polygon regions whose boundaries are continguous: https://www.rdocumentation.org/packages/spdep/versions/1.1-3/topics/poly2nb
VI.1.ii) nb2lines stores spatial coordinates as vectors (needed for weights): https://www.rdocumentation.org/packages/spdep/versions/1.1-3/topics/nb2lines
VI.1.iii) queen - FALSE defines the neighbourhood as Rook https://cran.r-project.org/web/packages/spdep/vignettes/nb.pdf
VI.2. Note for Different Variable: the crd.data will reflect the new data column as created in Step IV.7 after removing the na values from the new variable column. The tm_layout value for title will need to be changed for the new variable.
VI.3. Defining Queen Neighbours
crd.nb <- poly2nb(crd.data)
crd.net <- nb2lines(crd.nb,coords=coordinates(crd.data))
VI.3.i) Map Queen Neighbours (default weight scheme).
map_crd.net <- tm_shape(crd.data) + tm_borders(col='lightgrey') +
tm_shape(crd.net) + tm_lines(col='red') +
# add scale bar - first value is TOP, second value is BOTTOM
tm_scale_bar(width = 0.22, position = c("RIGHT", "BOTTOM")) +
# add compass
tm_compass(position = c("RIGHT", "TOP")) +
tm_layout(title = "CRD 2016 Census Tracts Median $ Income\nQueen Neighbours", title.position = c("LEFT", "TOP"), inner.margins = c(.08, .03, .08, .03))
VI.4. Defining Rook Neighbours
crd.nb2 <- poly2nb(crd.data, queen = FALSE)
crd.net2 <- nb2lines(crd.nb2,coords=coordinates(crd.data))
VI.4.i) Mapping Rook Neighbours
map_crd.net2 <- tm_shape(crd.data) + tm_borders(col='lightgrey') +
tm_shape(crd.net) + tm_lines(col='red', lwd = 2) +
tm_shape(crd.net2) + tm_lines(col='blue', lwd = 2) +
# add scale bar - first value is TOP, second value is BOTTOM
tm_scale_bar(width = 0.22, position = c("RIGHT", "BOTTOM")) +
# add compass
tm_compass(position = c("RIGHT", "TOP")) +
tm_layout(title = "CRD 2016 Census Tracts Median $ Income\nRook (blue) & Queen (blue & red) Neighbours", title.position = c("LEFT", "TOP"), inner.margins = c(.08, .03, .08, .03))
VI.5. Create a png image file for the Queen Neighbourhoods weighting (crd.net), and the Rooks Neighbourhoods (crd.net2).
VI.6. Note for Different Variable: change the file name to reflect the new variable when saving the file.
map_crd.net
map_crd.net2
# create png of Median Income CRD net
png("map_crd.net_MdInc.png")
map_crd.net
dev.off()
# create png of Median Income CRD net 2
png("map_crd.net2_MdInc.png")
map_crd.net2
dev.off()
VI.7. Don't be concerned about the warning that current projection of shape crd.net is unknown.
## Warning: Current projection of shape crd.net unknown and cannot be
## determined.
VI.8. Create Neighbourhood map for Queen's weight
VI.8. Create Neighbourhood map for Rook's weight plus Queens weight
VI.9. Interpreting the output of Rook & Queen Neighbours:
ToCVII.1. Using the spdep package, nb2listw uses the neighbours list previously created and adds weight: https://www.rdocumentation.org/packages/spdep/versions/1.1-3/topics/nb2listw
crd.lw <- nb2listw(crd.nb, zero.policy = TRUE, style = "W")
print.listw(crd.lw, zero.policy = TRUE)
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 74
## Number of nonzero links: 424
## Percentage nonzero weights: 7.742878
## Average number of links: 5.72973
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 74 5476 74 29.46931 322.5924
VIII.1. The lagged means, lag.listw, is calculated using the spdep package.
VIII.2. A map is created using the tmap package.
VIII.3. Note for Different Variable: The column name needs to be changed to the different variable name, and the name of the new column created should reflect the new variable name. This new column name is used by tm_polygon to colour the classification scheme. The title for the legend and tm_layout will also need to be changed to the new variable name.
# Calculate Lag Means
crd.data$IncLagMeans = lag.listw(crd.lw, crd.data$MdInc, zero.policy = TRUE)
map_LagMean <- tm_shape(crd.data) +
tm_polygons(col = "IncLagMeans",
title = "Median $ Income\nLagged Means",
style = "jenks",
palette = "viridis", n = 6) +
# add scale bar - first value is TOP, second value is BOTTOM
tm_scale_bar(width = 0.22, position = c("RIGHT", "BOTTOM")) +
# add compass
tm_compass(position = c("RIGHT", "TOP")) +
tm_layout(title = "CRD 2016 Census Tracts Median $ Income\nLagged Means", title.position = c("LEFT", "TOP"), inner.margins = c(.08, .03, .08, .03))
Interpreting the output of the Lagged Mean Map:
Note for Different Variable: Change the file name to reflect the new variable when saving the file.
map_LagMean
# create png of map LagMean Median Income CRD net 2
png("map_LagMean_MdInc.png")
map_LagMean
dev.off()
Global Moran's I Test creates one I value for every polygon as if they were all one polygon.
IX.1. First a calculation will be done to test the random expected range of -0.8193252 to 1.0821054.
# Global Moran's I test
mi <- moran.test(crd.data$MdInc, crd.lw, zero.policy = TRUE)
mi
##
## Moran I test under randomisation
##
## data: crd.data$MdInc
## weights: crd.lw
##
## Moran I statistic standard deviate = 2.2665, p-value = 0.01171
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.143170949 -0.013698630 0.004790316
moran.range <- function(lw) {
wmat <- listw2mat(lw)
return(range(eigen((wmat + t(wmat))/2)$values))
}
moran.range(crd.lw)
## [1] -0.8193252 1.0821054
IX.2. Calculate the z value, by calculating the values of the Global Moran's I, the estimated value of Global Moran's I, and the Variable. The z value is Global Moran's I less the Estimated Global Moran's I, divided by the square root of the variance.
## calculate z value
mI <- mi$estimate[[1]]
eI <- mi$estimate[[2]]
var <- mi$estimate[[3]]
z <- (mI - eI) / var^(1/2)
IX.3. Add the data variables to create a table of values
# Add Stats data objects to data.frame ########
data.for.table1 = data.frame(mI, eI, var, z)
data.for.table1
## mI eI var z
## 1 0.1431709 -0.01369863 0.004790316 2.266505
table1 <- data.for.table1
IX.4. Interpreting the output: mI > eI a positive clustered autocorrelation, and z > 1.96 for a greater than 95% confidence that the Global Moran's I value is significantly greater than would
IX.5. Note for Different Variable: change the file name to reflect the new variable.
# write csv file with Global Moran's results for I, E, Var, Z, P
write.csv(table1, "MoransGlobal_MdInc.csv", row.names = FALSE)
ToC
The Local Moran's i test is used as a local indicator of spatial association (Lisa), by calculating separate i values for every polygon in order to test how spatially autocorrelated each polygon is,
- how to perform test with R - describe output - what does it show? If the local moran's i value is positive, the polygon will have a similar value to its neighbours, and be considered clustered. If the local moran's i value is negative, the polygon will have a different value to its neighbours, and be considered dispersed.
X.1. The lisa test creates 5 different variables for local moran's i, the local moran's i, the Expected i, the Variance of i, the Z value and the p value. These values are assigned to new columns and given new column names.
X.2. Note for Different Variable: Change the new variable column name after the $ e.g. crd.data$NewColumnName
# Lisa Test
lisa.test <- localmoran(crd.data$MdInc, crd.lw)
crd.data$Ii <- lisa.test[,1]
crd.data$E.Ii<- lisa.test[,2]
crd.data$Var.Ii<- lisa.test[,3]
crd.data$Z.Ii<- lisa.test[,4]
crd.data$P<- lisa.test[,5]
X.3. Calculate significance values
crd.data$IiLessEIi <- crd.data$Ii - crd.data$E.Ii
# if z value is less than 95% significance level, value is FALSE
crd.data$Z.Ho <- crd.data$Z.Ii < 1.96
# if Ii < EIi, value is negative
crd.data$Ii.Neg <- crd.data$IiLessEIi < 0
X.4. Add the data values to a data frame to create a new table of values for the Local Moran's Ii.
# make a data frome from the Lisa Test Data
data.for.table2 = data.frame(crd.data$Ii, crd.data$E.Ii, crd.data$Var.Ii, crd.data$IiLessEIi, crd.data$Ii.Neg, crd.data$Z.Ii, crd.data$Z.Ho, crd.data$P)
table2 <- data.for.table2
table2
## crd.data.Ii crd.data.E.Ii crd.data.Var.Ii crd.data.IiLessEIi
## 1 0.2744151305 -0.01369863 0.14693811 0.2881137606
## 2 0.0557309861 -0.01369863 0.22656398 0.0694296162
## 3 -0.2349448907 -0.01369863 0.14693811 -0.2212462606
## 4 0.3326591723 -0.01369863 0.22656398 0.3463578025
## 5 0.1990156088 -0.01369863 0.17878846 0.2127142389
## 6 -0.0338955771 -0.01369863 0.14693811 -0.0201969469
## 7 0.0128136506 -0.01369863 0.17878846 0.0265122808
## 8 0.7919282302 -0.01369863 0.17878846 0.8056268603
## 9 -0.0464760893 -0.01369863 0.17878846 -0.0327774592
## 10 0.4657065327 -0.01369863 0.12418786 0.4794051628
## 11 0.0067879072 -0.01369863 0.10712518 0.0204865373
## 12 0.2947025748 -0.01369863 0.09385420 0.3084012049
## 13 -0.0923032394 -0.01369863 0.46544159 -0.0786046092
## 14 -0.0285935895 -0.01369863 0.30618985 -0.0148949593
## 15 0.4126788509 -0.01369863 0.09385420 0.4263774810
## 16 -0.2187477999 -0.01369863 0.10712518 -0.2050491697
## 17 0.4433240328 -0.01369863 0.12418786 0.4570226629
## 18 1.2205874502 -0.01369863 0.30618985 1.2342860804
## 19 0.4616826274 -0.01369863 0.17878846 0.4753812575
## 20 0.2638695873 -0.01369863 0.12418786 0.2775682174
## 21 2.1078481863 -0.01369863 0.22656398 2.1215468164
## 22 -0.0124306146 -0.01369863 0.17878846 0.0012680156
## 23 0.0980896658 -0.01369863 0.08323742 0.1117882959
## 24 -0.0231942136 -0.01369863 0.22656398 -0.0094955835
## 25 0.1047856106 -0.01369863 0.12418786 0.1184842408
## 26 0.0246750827 -0.01369863 0.10712518 0.0383737128
## 27 0.1900350829 -0.01369863 0.14693811 0.2037337130
## 28 0.0019320024 -0.01369863 0.08323742 0.0156306325
## 29 0.9284834191 -0.01369863 0.17878846 0.9421820492
## 30 0.2032721811 -0.01369863 0.17878846 0.2169708112
## 31 0.1259530137 -0.01369863 0.12418786 0.1396516438
## 32 0.4535832667 -0.01369863 0.17878846 0.4672818968
## 33 0.0760594107 -0.01369863 0.14693811 0.0897580408
## 34 0.6622247407 -0.01369863 0.30618985 0.6759233709
## 35 0.2429889826 -0.01369863 0.17878846 0.2566876128
## 36 0.0356433844 -0.01369863 0.17878846 0.0493420145
## 37 -0.0004154355 -0.01369863 0.17878846 0.0132831947
## 38 -0.0193496543 -0.01369863 0.12418786 -0.0056510242
## 39 0.0074249521 -0.01369863 0.17878846 0.0211235823
## 40 0.1856654215 -0.01369863 0.22656398 0.1993640517
## 41 0.0595949852 -0.01369863 0.12418786 0.0732936153
## 42 -0.2359649540 -0.01369863 0.12418786 -0.2222663239
## 43 -0.0579399062 -0.01369863 0.14693811 -0.0442412761
## 44 0.1041228751 -0.01369863 0.30618985 0.1178215053
## 45 0.1458425508 -0.01369863 0.14693811 0.1595411809
## 46 -0.1574332175 -0.01369863 0.17878846 -0.1437345874
## 47 0.0536021511 -0.01369863 0.22656398 0.0673007812
## 48 0.6232149028 -0.01369863 0.14693811 0.6369135329
## 49 0.5484882380 -0.01369863 0.10712518 0.5621868682
## 50 -0.0440773555 -0.01369863 0.02923030 -0.0303787253
## 51 -1.9458534186 -0.01369863 0.30618985 -1.9321547885
## 52 -0.3484095342 -0.01369863 0.14693811 -0.3347109041
## 53 -0.0252503561 -0.01369863 0.12418786 -0.0115517259
## 54 1.0843288469 -0.01369863 0.14693811 1.0980274770
## 55 0.2836125755 -0.01369863 0.10712518 0.2973112057
## 56 0.0384246510 -0.01369863 0.22656398 0.0521232811
## 57 0.0025467591 -0.01369863 0.14693811 0.0162453892
## 58 -0.0346348845 -0.01369863 0.22656398 -0.0209362543
## 59 -3.7243153785 -0.01369863 0.94319681 -3.7106167484
## 60 0.0117550365 -0.01369863 0.09385420 0.0254536667
## 61 0.2614551227 -0.01369863 0.22656398 0.2751537528
## 62 0.1399328018 -0.01369863 0.46544159 0.1536314319
## 63 0.4384667710 -0.01369863 0.46544159 0.4521654011
## 64 -0.3612599654 -0.01369863 0.94319681 -0.3475613352
## 65 0.2206921127 -0.01369863 0.14693811 0.2343907428
## 66 0.4853407163 -0.01369863 0.14693811 0.4990393464
## 67 0.3210323517 -0.01369863 0.22656398 0.3347309818
## 68 -0.0802313849 -0.01369863 0.14693811 -0.0665327548
## 69 0.3159098132 -0.01369863 0.12418786 0.3296084434
## 70 0.6415313555 -0.01369863 0.12418786 0.6552299857
## 71 0.4949972138 -0.01369863 0.17878846 0.5086958439
## 72 -0.0146686830 -0.01369863 0.30618985 -0.0009700528
## 73 0.2927142780 -0.01369863 0.09385420 0.3064129082
## 74 1.0828675345 -0.01369863 0.30618985 1.0965661646
## crd.data.Ii.Neg crd.data.Z.Ii crd.data.Z.Ho crd.data.P
## 1 FALSE 0.751617311 TRUE 2.261406e-01
## 2 FALSE 0.145864407 TRUE 4.420142e-01
## 3 TRUE -0.577176596 TRUE 7.180899e-01
## 4 FALSE 0.727661741 TRUE 2.334103e-01
## 5 FALSE 0.503068148 TRUE 3.074582e-01
## 6 TRUE -0.052688823 TRUE 5.210101e-01
## 7 FALSE 0.062701416 TRUE 4.750021e-01
## 8 FALSE 1.905303636 TRUE 2.837031e-02
## 9 TRUE -0.077518533 TRUE 5.308945e-01
## 10 FALSE 1.360389058 TRUE 8.685342e-02
## 11 FALSE 0.062592570 TRUE 4.750455e-01
## 12 FALSE 1.006674841 TRUE 1.570455e-01
## 13 TRUE -0.115216696 TRUE 5.458633e-01
## 14 TRUE -0.026918070 TRUE 5.107375e-01
## 15 FALSE 1.391769798 TRUE 8.199606e-02
## 16 TRUE -0.626487255 TRUE 7.345023e-01
## 17 FALSE 1.296875124 TRUE 9.733708e-02
## 18 FALSE 2.230593514 FALSE 1.285403e-02
## 19 FALSE 1.124274380 TRUE 1.304483e-01
## 20 FALSE 0.787644346 TRUE 2.154524e-01
## 21 FALSE 4.457149340 FALSE 4.152837e-06
## 22 FALSE 0.002998851 TRUE 4.988036e-01
## 23 FALSE 0.387469072 TRUE 3.492045e-01
## 24 TRUE -0.019949234 TRUE 5.079581e-01
## 25 FALSE 0.336218041 TRUE 3.683532e-01
## 26 FALSE 0.117243303 TRUE 4.533336e-01
## 27 FALSE 0.531490704 TRUE 2.975394e-01
## 28 FALSE 0.054177288 TRUE 4.783970e-01
## 29 FALSE 2.228255999 FALSE 1.293172e-02
## 30 FALSE 0.513134921 TRUE 3.039285e-01
## 31 FALSE 0.396283943 TRUE 3.459478e-01
## 32 FALSE 1.105119431 TRUE 1.345539e-01
## 33 FALSE 0.234156457 TRUE 4.074318e-01
## 34 FALSE 1.221524176 TRUE 1.109438e-01
## 35 FALSE 0.607064965 TRUE 2.719039e-01
## 36 FALSE 0.116693626 TRUE 4.535514e-01
## 37 FALSE 0.031414691 TRUE 4.874694e-01
## 38 TRUE -0.016035688 TRUE 5.063970e-01
## 39 FALSE 0.049957170 TRUE 4.800783e-01
## 40 FALSE 0.418843150 TRUE 3.376654e-01
## 41 FALSE 0.207982392 TRUE 4.176214e-01
## 42 TRUE -0.630716351 TRUE 7.358870e-01
## 43 TRUE -0.115414512 TRUE 5.459417e-01
## 44 FALSE 0.212926233 TRUE 4.156923e-01
## 45 FALSE 0.416203354 TRUE 3.386306e-01
## 46 TRUE -0.339931605 TRUE 6.330460e-01
## 47 FALSE 0.141391946 TRUE 4.437802e-01
## 48 FALSE 1.661549369 TRUE 4.830158e-02
## 49 FALSE 1.717650982 TRUE 4.293015e-02
## 50 TRUE -0.177685867 TRUE 5.705152e-01
## 51 TRUE -3.491777155 TRUE 9.997601e-01
## 52 TRUE -0.873177696 TRUE 8.087169e-01
## 53 TRUE -0.032779876 TRUE 5.130749e-01
## 54 FALSE 2.864481232 FALSE 2.088465e-03
## 55 FALSE 0.908375690 TRUE 1.818399e-01
## 56 FALSE 0.109505596 TRUE 4.564007e-01
## 57 FALSE 0.042380189 TRUE 4.830978e-01
## 58 TRUE -0.043984894 TRUE 5.175418e-01
## 59 TRUE -3.820717593 TRUE 9.999335e-01
## 60 FALSE 0.083085168 TRUE 4.668919e-01
## 61 FALSE 0.578069434 TRUE 2.816086e-01
## 62 FALSE 0.225189161 TRUE 4.109161e-01
## 63 FALSE 0.662772886 TRUE 2.537380e-01
## 64 TRUE -0.357874121 TRUE 6.397812e-01
## 65 FALSE 0.611467288 TRUE 2.704451e-01
## 66 FALSE 1.301869827 TRUE 9.648044e-02
## 67 FALSE 0.703235000 TRUE 2.409547e-01
## 68 TRUE -0.173567448 TRUE 5.688973e-01
## 69 FALSE 0.935316835 TRUE 1.748125e-01
## 70 FALSE 1.859320199 TRUE 3.149088e-02
## 71 FALSE 1.203063216 TRUE 1.144759e-01
## 72 TRUE -0.001753073 TRUE 5.006994e-01
## 73 FALSE 1.000184696 TRUE 1.586106e-01
## 74 FALSE 1.981707007 FALSE 2.375602e-02
head(table2)
## crd.data.Ii crd.data.E.Ii crd.data.Var.Ii crd.data.IiLessEIi
## 1 0.27441513 -0.01369863 0.1469381 0.28811376
## 2 0.05573099 -0.01369863 0.2265640 0.06942962
## 3 -0.23494489 -0.01369863 0.1469381 -0.22124626
## 4 0.33265917 -0.01369863 0.2265640 0.34635780
## 5 0.19901561 -0.01369863 0.1787885 0.21271424
## 6 -0.03389558 -0.01369863 0.1469381 -0.02019695
## crd.data.Ii.Neg crd.data.Z.Ii crd.data.Z.Ho crd.data.P
## 1 FALSE 0.75161731 TRUE 0.2261406
## 2 FALSE 0.14586441 TRUE 0.4420142
## 3 TRUE -0.57717660 TRUE 0.7180899
## 4 FALSE 0.72766174 TRUE 0.2334103
## 5 FALSE 0.50306815 TRUE 0.3074582
## 6 TRUE -0.05268882 TRUE 0.5210101
X.5. Note for Different Variable: change the file name to reflect the new variable.
# write csv file with Lisa Test results for I, E, Var, Z, P
write.csv(table2, "LisaTest_MdInc.csv", row.names = FALSE)
ToC
XI.1. Note for Different Variable: Change the following values tm_polygons(col = "NewColumnName", title= "New Title"); tm_layout(title= "New Title"). For the new variable, the crd.data will reflect the new data column as created in Step IV.6 after removing the na values from the new variable column.
# Map Lisa
map_LISA <- tm_shape(crd.data) +
tm_polygons(col = "Ii",
title = "Local Moran's i",
style = "jenks",
palette = "viridis", n = 6) +
# add scale bar - first value is TOP, second value is BOTTOM
tm_scale_bar(width = 0.22, position = c("RIGHT", "BOTTOM")) +
# add compass
tm_compass(position = c("RIGHT", "TOP")) +
tm_layout(title = "CRD 2016 Census Tracts Median $ Income\nLocal Moran's i", title.position = c("LEFT", "TOP"), inner.margins = c(.08, .03, .08, .03))
XI.2. Note for Different Variable: change the file name to reflect the new variable.
map_LISA
# create png of Map LISA MdInc
png("map_LISA_MdInc.png")
map_LISA
dev.off()
ToC
XII.1. Note for Different Variable: When running the code for the second variable, the new column name will need to replace thw item crd.data$MdInc, the xlab="New Variable Label Name", the ylab = "New Variable Label Name", main = New VariableName title.
# Scatter Plot
moran.plot(crd.data$MdInc, crd.lw, zero.policy=NULL, spChk=NULL, labels=NULL, xlab="Median $ Income",
ylab="Spatially Lagged Mean Median $ Income", quiet=NULL, main ="Local Moran's i Plot for CRD 2016 Median $ Income")
XII.2. Note for Different Variable: change the saved file name to reflect the new variable. When running the code for the second variable, the new column name will need to replace thw item crd.data$MdInc, the xlab="New Variable Label Name", the ylab = "New Variable Label Name", main = New VariableName title.
# create png of Moran scatter plot
png("moran.plot_MdInc.png")
moran.plot(crd.data$MdInc, crd.lw, zero.policy=NULL, spChk=NULL, labels=NULL, xlab="Median $ Income",
ylab="Spatially Lagged Median $ Income", quiet=NULL, main ="Local Moran's i Plot for CRD 2016 Median $ Income")
dev.off()
https://www.rdocumentation.org/packages/plyr/versions/1.8.4 ↳
https://github.com/rstudio/cheatsheets/blob/master/data-transformation.pdf ↳
https://cran.r-project.org/web/packages/spdep/index.html ↳
https://cran.r-project.org/web/packages/GISTools/GISTools.pdf ↳
https://cran.r-project.org/web/packages/raster/raster.pdf ↳
https://cran.r-project.org/web/packages/maptools/maptools.pdf ↳
https://cran.r-project.org/web/packages/BAMMtools/BAMMtools.pdf ↳
Built with R v.3.6.1 with R Markdown in RStudio v.1.1.463.
Moran's Global I Test creates one i value for every polygon as if they were all one polygon.
# Moran's I test
mi <- moran.test(crd.data$RenterNorm, crd.lw, zero.policy = TRUE)
mi
##
## Moran I test under randomisation
##
## data: crd.data$RenterNorm
## weights: crd.lw
##
## Moran I statistic standard deviate = 8.6868, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.595415952 -0.013698630 0.004916714
moran.range <- function(lw) {
wmat <- listw2mat(lw)
return(range(eigen((wmat + t(wmat))/2)$values))
}
moran.range(crd.lw)
## [1] -0.8193252 1.0821054
## calculate z value
mI <- mi$estimate[[1]]
eI <- mi$estimate[[2]]
var <- mi$estimate[[3]]
z <- (mI - eI) / var^(1/2)
# Add Stats data objects to data.frame ########
data.for.table3 = data.frame(mI, eI, var, z)
data.for.table3
## mI eI var z
## 1 0.595416 -0.01369863 0.004916714 8.686834
table3 <- data.for.table3
# write csv file with Global Moran's results for I, E, Var, Z, P
write.csv(table3, "MoransGlobal_RenterNorm.csv", row.names = FALSE)
Separate i values are created for every single polygon.
# Lisa Test
lisa.test <- localmoran(crd.data$RenterNorm, crd.lw)
crd.data$Ii <- lisa.test[,1]
crd.data$E.Ii<- lisa.test[,2]
crd.data$Var.Ii<- lisa.test[,3]
crd.data$Z.Ii<- lisa.test[,4]
crd.data$P<- lisa.test[,5]
crd.data$IiLessEIi <- crd.data$Ii - crd.data$E.Ii
# if z value is less than 95% significance level, value is FALSE
crd.data$Z.Ho <- crd.data$Z.Ii < 1.96
crd.data$Z.Ho
## [1] TRUE TRUE TRUE TRUE TRUE FALSE FALSE FALSE FALSE TRUE FALSE
## [12] TRUE TRUE TRUE FALSE FALSE FALSE TRUE TRUE TRUE TRUE FALSE
## [23] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE
## [34] TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [45] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE
## [56] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [67] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
crd.data$Ii.Neg <- crd.data$IiLessEIi < 0
crd.data$Ii.Neg
## [1] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE TRUE FALSE
## [12] FALSE FALSE FALSE FALSE FALSE FALSE FALSE TRUE TRUE FALSE FALSE
## [23] FALSE TRUE FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE FALSE
## [34] TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [45] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [56] FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE FALSE FALSE FALSE
## [67] FALSE TRUE FALSE FALSE FALSE FALSE FALSE FALSE
# make a data frome from the Lisa Test Data
data.for.table4 = data.frame(crd.data$Ii, crd.data$E.Ii, crd.data$Var.Ii, crd.data$IiLessEIi, crd.data$Ii.Neg, crd.data$Z.Ii, crd.data$Z.Ho, crd.data$P)
table4 <- data.for.table4
table4
## crd.data.Ii crd.data.E.Ii crd.data.Var.Ii crd.data.IiLessEIi
## 1 0.233821626 -0.01369863 0.15045712 0.247520256
## 2 0.032943097 -0.01369863 0.23218304 0.046641727
## 3 0.549070589 -0.01369863 0.15045712 0.562769219
## 4 0.503692530 -0.01369863 0.23218304 0.517391160
## 5 0.060472308 -0.01369863 0.18314749 0.074170938
## 6 1.556366078 -0.01369863 0.15045712 1.570064708
## 7 3.391433333 -0.01369863 0.18314749 3.405131963
## 8 1.553256264 -0.01369863 0.18314749 1.566954894
## 9 1.140736503 -0.01369863 0.18314749 1.154435133
## 10 -0.071622535 -0.01369863 0.12710685 -0.057923905
## 11 1.333661488 -0.01369863 0.10959416 1.347360119
## 12 0.360821928 -0.01369863 0.09597317 0.374520558
## 13 0.205149532 -0.01369863 0.47736081 0.218848163
## 14 0.316206463 -0.01369863 0.31390896 0.329905094
## 15 3.643380638 -0.01369863 0.09597317 3.657079268
## 16 1.379862981 -0.01369863 0.10959416 1.393561611
## 17 1.839610501 -0.01369863 0.12710685 1.853309131
## 18 0.359894929 -0.01369863 0.31390896 0.373593560
## 19 -0.045775803 -0.01369863 0.18314749 -0.032077173
## 20 -0.018460853 -0.01369863 0.12710685 -0.004762223
## 21 0.640176636 -0.01369863 0.23218304 0.653875266
## 22 1.897336263 -0.01369863 0.18314749 1.911034893
## 23 0.098397443 -0.01369863 0.08507638 0.112096073
## 24 -0.026062633 -0.01369863 0.23218304 -0.012364003
## 25 0.147756323 -0.01369863 0.12710685 0.161454953
## 26 0.431022202 -0.01369863 0.10959416 0.444720832
## 27 0.356600204 -0.01369863 0.15045712 0.370298834
## 28 0.304032189 -0.01369863 0.08507638 0.317730819
## 29 0.577861864 -0.01369863 0.18314749 0.591560494
## 30 -0.168322406 -0.01369863 0.18314749 -0.154623776
## 31 -0.001621923 -0.01369863 0.12710685 0.012076707
## 32 4.474925356 -0.01369863 0.18314749 4.488623986
## 33 0.064613726 -0.01369863 0.15045712 0.078312356
## 34 -0.904832954 -0.01369863 0.31390896 -0.891134324
## 35 0.011019807 -0.01369863 0.18314749 0.024718437
## 36 3.129918431 -0.01369863 0.18314749 3.143617061
## 37 0.136728114 -0.01369863 0.18314749 0.150426744
## 38 0.199580888 -0.01369863 0.12710685 0.213279518
## 39 0.097877099 -0.01369863 0.18314749 0.111575729
## 40 0.529096109 -0.01369863 0.23218304 0.542794739
## 41 0.236191920 -0.01369863 0.12710685 0.249890550
## 42 0.226356428 -0.01369863 0.12710685 0.240055058
## 43 0.604421772 -0.01369863 0.15045712 0.618120402
## 44 0.442216620 -0.01369863 0.31390896 0.455915250
## 45 0.147383189 -0.01369863 0.15045712 0.161081819
## 46 0.713133316 -0.01369863 0.18314749 0.726831946
## 47 0.387150645 -0.01369863 0.23218304 0.400849276
## 48 0.612030700 -0.01369863 0.15045712 0.625729331
## 49 0.156866528 -0.01369863 0.10959416 0.170565159
## 50 0.234052890 -0.01369863 0.02964488 0.247751521
## 51 0.662880165 -0.01369863 0.31390896 0.676578795
## 52 0.630290771 -0.01369863 0.15045712 0.643989401
## 53 0.223824443 -0.01369863 0.12710685 0.237523073
## 54 1.838131703 -0.01369863 0.15045712 1.851830333
## 55 0.305859518 -0.01369863 0.10959416 0.319558149
## 56 0.387356307 -0.01369863 0.23218304 0.401054937
## 57 0.099487976 -0.01369863 0.15045712 0.113186606
## 58 0.046748357 -0.01369863 0.23218304 0.060446987
## 59 0.998762311 -0.01369863 0.96771634 1.012460942
## 60 0.043229883 -0.01369863 0.09597317 0.056928513
## 61 -0.252114651 -0.01369863 0.23218304 -0.238416021
## 62 0.355471073 -0.01369863 0.47736081 0.369169703
## 63 0.645739004 -0.01369863 0.47736081 0.659437634
## 64 0.290854679 -0.01369863 0.96771634 0.304553309
## 65 0.299887489 -0.01369863 0.15045712 0.313586119
## 66 0.601525828 -0.01369863 0.15045712 0.615224458
## 67 0.811880946 -0.01369863 0.23218304 0.825579576
## 68 -0.037825907 -0.01369863 0.15045712 -0.024127277
## 69 0.151580204 -0.01369863 0.12710685 0.165278834
## 70 0.262426181 -0.01369863 0.12710685 0.276124811
## 71 0.584495330 -0.01369863 0.18314749 0.598193960
## 72 0.157706771 -0.01369863 0.31390896 0.171405401
## 73 0.104039684 -0.01369863 0.09597317 0.117738314
## 74 0.768114071 -0.01369863 0.31390896 0.781812701
## crd.data.Ii.Neg crd.data.Z.Ii crd.data.Z.Ho crd.data.P
## 1 FALSE 0.63812297 TRUE 2.616968e-01
## 2 FALSE 0.09679644 TRUE 4.614440e-01
## 3 FALSE 1.45085485 TRUE 7.341014e-02
## 4 FALSE 1.07375144 TRUE 1.414671e-01
## 5 FALSE 0.17331386 TRUE 4.312024e-01
## 6 FALSE 4.04772671 FALSE 2.585874e-05
## 7 FALSE 7.95670883 FALSE 8.833791e-16
## 8 FALSE 3.66147450 FALSE 1.253839e-04
## 9 FALSE 2.69754721 FALSE 3.492619e-03
## 10 TRUE -0.16247006 TRUE 5.645321e-01
## 11 FALSE 4.06995858 FALSE 2.351075e-05
## 12 FALSE 1.20892885 TRUE 1.133451e-01
## 13 FALSE 0.31675212 TRUE 3.757159e-01
## 14 FALSE 0.58882625 TRUE 2.779889e-01
## 15 FALSE 11.80482233 FALSE 1.842843e-32
## 16 FALSE 4.20951901 FALSE 1.279575e-05
## 17 FALSE 5.19832443 FALSE 1.005465e-07
## 18 FALSE 0.66680297 TRUE 2.524490e-01
## 19 TRUE -0.07495414 TRUE 5.298744e-01
## 20 TRUE -0.01335750 TRUE 5.053287e-01
## 21 FALSE 1.35699942 TRUE 8.739070e-02
## 22 FALSE 4.46547986 FALSE 3.994474e-06
## 23 FALSE 0.38431373 TRUE 3.503730e-01
## 24 TRUE -0.02565924 TRUE 5.102354e-01
## 25 FALSE 0.45286305 TRUE 3.253237e-01
## 26 FALSE 1.34336421 TRUE 8.957703e-02
## 27 FALSE 0.95465395 TRUE 1.698764e-01
## 28 FALSE 1.08931840 TRUE 1.380067e-01
## 29 FALSE 1.38228846 TRUE 8.344157e-02
## 30 TRUE -0.36130652 TRUE 6.410648e-01
## 31 FALSE 0.03387381 TRUE 4.864889e-01
## 32 FALSE 10.48848459 FALSE 4.879043e-26
## 33 FALSE 0.20189424 TRUE 4.199997e-01
## 34 TRUE -1.59052798 TRUE 9.441421e-01
## 35 FALSE 0.05775911 TRUE 4.769703e-01
## 36 FALSE 7.34563180 FALSE 1.023947e-13
## 37 FALSE 0.35149939 TRUE 3.626069e-01
## 38 FALSE 0.59822515 TRUE 2.748449e-01
## 39 FALSE 0.26071694 TRUE 3.971554e-01
## 40 FALSE 1.12647195 TRUE 1.299829e-01
## 41 FALSE 0.70091499 TRUE 2.416780e-01
## 42 FALSE 0.67332754 TRUE 2.503695e-01
## 43 FALSE 1.59355372 TRUE 5.551801e-02
## 44 FALSE 0.81373362 TRUE 2.078988e-01
## 45 FALSE 0.41527918 TRUE 3.389688e-01
## 46 FALSE 1.69837476 TRUE 4.471853e-02
## 47 FALSE 0.83188991 TRUE 2.027355e-01
## 48 FALSE 1.61317003 TRUE 5.335379e-02
## 49 FALSE 0.51522464 TRUE 3.031980e-01
## 50 FALSE 1.43893589 TRUE 7.508434e-02
## 51 FALSE 1.20758171 TRUE 1.136041e-01
## 52 FALSE 1.66024565 TRUE 4.843252e-02
## 53 FALSE 0.66622560 TRUE 2.526335e-01
## 54 FALSE 4.77413642 FALSE 9.024001e-07
## 55 FALSE 0.96528642 TRUE 1.672007e-01
## 56 FALSE 0.83231672 TRUE 2.026151e-01
## 57 FALSE 0.29180227 TRUE 3.852189e-01
## 58 FALSE 0.12544675 TRUE 4.500849e-01
## 59 FALSE 1.02921058 TRUE 1.516904e-01
## 60 FALSE 0.18376167 TRUE 4.271002e-01
## 61 TRUE -0.49478917 TRUE 6.896255e-01
## 62 FALSE 0.53432153 TRUE 2.965595e-01
## 63 FALSE 0.95444377 TRUE 1.699295e-01
## 64 FALSE 0.30959168 TRUE 3.784357e-01
## 65 FALSE 0.80844497 TRUE 2.094172e-01
## 66 FALSE 1.58608780 TRUE 5.635970e-02
## 67 FALSE 1.71334055 TRUE 4.332495e-02
## 68 TRUE -0.06220166 TRUE 5.247989e-01
## 69 FALSE 0.46358861 TRUE 3.214713e-01
## 70 FALSE 0.77449915 TRUE 2.193178e-01
## 71 FALSE 1.39778875 TRUE 8.108826e-02
## 72 FALSE 0.30593041 TRUE 3.798288e-01
## 73 FALSE 0.38005189 TRUE 3.519534e-01
## 74 FALSE 1.39540689 TRUE 8.144659e-02
head(table4)
## crd.data.Ii crd.data.E.Ii crd.data.Var.Ii crd.data.IiLessEIi
## 1 0.23382163 -0.01369863 0.1504571 0.24752026
## 2 0.03294310 -0.01369863 0.2321830 0.04664173
## 3 0.54907059 -0.01369863 0.1504571 0.56276922
## 4 0.50369253 -0.01369863 0.2321830 0.51739116
## 5 0.06047231 -0.01369863 0.1831475 0.07417094
## 6 1.55636608 -0.01369863 0.1504571 1.57006471
## crd.data.Ii.Neg crd.data.Z.Ii crd.data.Z.Ho crd.data.P
## 1 FALSE 0.63812297 TRUE 2.616968e-01
## 2 FALSE 0.09679644 TRUE 4.614440e-01
## 3 FALSE 1.45085485 TRUE 7.341014e-02
## 4 FALSE 1.07375144 TRUE 1.414671e-01
## 5 FALSE 0.17331386 TRUE 4.312024e-01
## 6 FALSE 4.04772671 FALSE 2.585874e-05
# write csv file with Lisa Test results for I, E, Var, Z, P
write.csv(table4, "LisaTest_RenterNorm.csv", row.names = FALSE)
Using the tmap package
# Map Lisa
map_LISA <- tm_shape(crd.data) +
tm_polygons(col = "Ii",
title = "Local Moran's i",
style = "jenks",
palette = "viridis", n = 6) +
# add scale bar - first value is TOP, second value is BOTTOM
tm_scale_bar(width = 0.22, position = c("RIGHT", "BOTTOM")) +
# add compass
tm_compass(position = c("RIGHT", "TOP")) +
tm_layout(title = "CRD 2016 Census Tracts Renter % of Population\nLocal Moran's i", title.position = c("LEFT", "TOP"), inner.margins = c(.08, .03, .08, .03))
map_LISA
# create png of Map LISA Renter Normalized
png("map_LISA_RenterNorm.png")
map_LISA
dev.off()
When running the code for the second variable, the new column name will need to replace this item crd.data$Renter.
# Scatter Plot
moran.plot(crd.data$RenterNorm, crd.lw, zero.policy=NULL, spChk=NULL, labels=NULL, xlab="Renter % of Population",
ylab="Spatially Lagged Renter % of Population", quiet=NULL, main ="Local Moran's i Plot for CRD 2016 Renter % of Population")
# create png of Moran scatter plot
png("moran.plot_RenterNorm.png")
moran.plot(crd.data$RenterNorm, crd.lw, zero.policy=NULL, spChk=NULL, labels=NULL, xlab="Renter % of Population",
ylab="Spatially Lagged Renter % of Population", quiet=NULL, main ="Local Moran's i Plot for CRD 2016 Renter % of Population")
dev.off()