This example demonstrates the process of calculating segregation and
diversity metrics for the user-defined scale. In the
raceland
package, scale is defined by the size
parameter expressed as a number of cells. RL grids have a resolution of
30 m. Thus the scale will be the multiplication of 30. For example, size
= 10 means 10*30m = 300m.
More technical details about calculating racial diversity and
segregation metric for a specific scale is included in the
raceland
tutorial https://cran.r-project.org/web/packages/raceland/vignettes/raceland-intro2.html
The inputs to calculate racial diversity and segregation metrics for the user-defined scale are RL grid racial ID, and RL grid population density from NRGD2020. This example uses two layers included in the Atlanta dataset. Atlanta dataset can be downloaded from http://www.socscape.edu.pl/index.php?id=nrgd)
The example below demonstrates how to calculate metrics for the scale
of 4.5 km (size = 150). The results are returned as a table showing a
value of segregation/diversity in each tile. Additionally, using
create_grid()
function from the raceland
package we can generate spatial object and join the metrics for further
visualization.
library(terra)
library(raceland)
#set working directory to the folder containing data using function setwd("")
##setwd("")
#read RL grid racial ID to R using terra package
= rast("Atlanta_dataset/RL_grids/rl_msa.tif")
rl #read RL grid population density to R using terra package
= rast("Atlanta_dataset/RL_grids/rd_msa.tif")
rd
#calculate metrics using calculate_metrics from raceland package. Metrics are calculated for each subtile of the specified size.
#Parameter threshold defines how many NA cells can be in the data. In this case, the metrics will be calculated if up to 70% of NA cells are in a given subtile. The metrics are not calculated for subtiles with more than 70% NA cells, and the NA is returned.
= calculate_metrics(rl, rd, neighbourhood = 4, fun = "mean", size = 150, threshold = 0.7)
metr
#metrics can be joined with the spatial object to prepare visualization.
#create_grid from the raceland package creates spatial object
= create_grid(rl, size = 150)
grid_sf = dplyr::left_join(grid_sf, metr, by = c("row", "col")) attr_grid
The code below allows to plot racial diversity and segregation map at the specific scale showing the spatial distribution of segregation/diversity within a city.
plot(attr_grid["ent"], breaks = seq(0, 2.6, by = 0.25), key.pos = 1,
pal = rev(hcl.colors(length(seq(0, 2.6, by = 0.25)) - 1, palette = "RdBu")),
bty = "n", main = "Racial diversity (Entropy)")
plot(attr_grid["mutinf"], breaks = seq(0, 1, by = 0.1), key.pos = 1,
pal = rev(hcl.colors(length(seq(0, 1, by = 0.1)) - 1, palette = "RdBu")),
bty = "n", main = "Racial segregation (Mutual information)")
The code below calculates metrics for a range of user-defined scales. In this example, the scale is defined as 30, 60, 90, and 120 cells, corresponding to the scale of 900m, 1.8km, 2.7km, and 3.6km. The results are saved in the “results” directory created in the working directory. The results directory contains two subdirectories: the shp subdirectory will contain ESRI Shapefiles that can be used for the visualization, the final subdirectory will contain text file metrics.csv with the segregation and diversity indices for each scale, and the segregation_profile.pdf.
library(terra)
library(raceland)
library(sf)
library(ggplot2)
#define scale: size 30 cells -> 900 m, 60 cells -> 1.8km, 90 cells -> 2.7km, 120 cells -> 3.6km
= c(30, 60, 90, 120)
scales
#set working directory to the folder containing data using function setwd("")
##setwd("")
#read RL grid racial ID to R using terra package
= rast("Atlanta_dataset/RL_grids/rl_msa.tif")
rl #read RL grid population density to R using terra package
= rast("Atlanta_dataset/RL_grids/rd_msa.tif")
rd
# create directory with the results, directory "result" will be created as sub-directory in the working directory
= getwd()
pf dir.create(file.path(pf, "results"), showWarnings = FALSE)
dir.create(file.path(pf, "results", "final"), showWarnings = FALSE)
dir.create(file.path(pf, "results", "shp"), showWarnings = FALSE)
= data.frame(size = integer())
complete_smr_df for (size in scales) {
#######################CALCULATE METRICS FOR SPECIFIED SCALE######################
#calculate metrics using calculate_metrics() from raceland package. Metrics are calculate for each subtile of the specified size.
= calculate_metrics(rl, rd,
metr_df neighbourhood = 4, fun = "mean",
size = size, shift = NULL, threshold = 0.7)
#based on entropy calculate Hill's number
$hill = 2^metr_df$ent
metr_df#calculate NMI - normalized mutual information (mutual information divided by entropy)
$NMI = metr_df$mutinf/metr_df$ent
metr_df
# Summarize metrics
## number of subtiles at a given scale
= metr_df[!is.na(metr_df$ent), ]
sel = nrow(sel)
nsubtiles
#Calculate an average value of metrics from all subtiles.
= c(scale = (size*30)/1000, n = nsubtiles, E = mean(metr_df$ent, na.rm = TRUE), MI = mean(metr_df$mutinf, na.rm = TRUE), HILL = mean(metr_df$hill, na.rm=TRUE), NMI = mean(metr_df$NMI, na.rm=TRUE))
complete_smr
#complete_smr_df is a data.frame with all metrics calculated for specified scales.
= rbind(complete_smr_df, complete_smr)
complete_smr_df
####################### CREATE SPATIAL OBJECT WITH METRICS ######################
# Create spatial object with metrics that will be save as ESRI Shapefile and can be used to show how the segregation/diversity change within a city at the specific scale.
# Create spatial object
= create_grid(rl, size = size)
grid_sf
# Join metrics to the grid
= dplyr::left_join(grid_sf, metr_df, by = c("row", "col"))
attr_grid = attr_grid[!is.na(attr_grid$ent),]
sel_grid
# Save the grid as shapefile
st_write(attr_grid,
file.path("results", "shp",
paste("metr_stat_", size, ".shp", sep = "")), delete_dsn = TRUE)
#stop size loop }
## Deleting source `results/shp/metr_stat_30.shp' using driver `ESRI Shapefile'
## Writing layer `metr_stat_30' to data source
## `results/shp/metr_stat_30.shp' using driver `ESRI Shapefile'
## Writing 49725 features with 9 fields and geometry type Polygon.
## Deleting source `results/shp/metr_stat_60.shp' using driver `ESRI Shapefile'
## Writing layer `metr_stat_60' to data source
## `results/shp/metr_stat_60.shp' using driver `ESRI Shapefile'
## Writing 12543 features with 9 fields and geometry type Polygon.
## Deleting source `results/shp/metr_stat_90.shp' using driver `ESRI Shapefile'
## Writing layer `metr_stat_90' to data source
## `results/shp/metr_stat_90.shp' using driver `ESRI Shapefile'
## Writing 5550 features with 9 fields and geometry type Polygon.
## Deleting source `results/shp/metr_stat_120.shp' using driver `ESRI Shapefile'
## Writing layer `metr_stat_120' to data source
## `results/shp/metr_stat_120.shp' using driver `ESRI Shapefile'
## Writing 3192 features with 9 fields and geometry type Polygon.
colnames(complete_smr_df) = c("scale", "n", "E", "MI", "HILL", "NMI")
write.csv(complete_smr_df, file.path("results", "final", "metrics.csv"))
####################### CREATE SEGREGATION PROFILE ######################
#plot segregation profile
= ggplot(complete_smr_df, aes(x = scale, y = MI)) +
p geom_line(linewidth= 0.7) +
scale_x_continuous(breaks = complete_smr_df$scale) +
labs(x = "scale (km)", y = "Mutual information") +
ylim(0, max(complete_smr_df$MI)) +
theme_bw()
p
#save segregation profile
ggsave(file.path("results", "final", "segregation_profile.pdf"))