Processing math: 100%

Lecture Content

Note: Document prepared for Spatial Socioeconometric Modeling and materials build on Gimond (2017) and Evangelista and Beskow (2018)

1 Differences between GIS and Spatial Analysis

1.1 Types of graphs (see Gimond (2017) for more details)

1.1.1 Reference maps

  • Examples include USGS maps, hiking maps, road maps
    • Used to navigate landscapes or identify locations of points-of-interest.
      Hiking map example, taken from Alltrails.com [Alltrails.com](https://www.alltrails.com/explore/us/pennsylvania/philadelphia?b_tl_lat=40.04416780544861&b_tl_lng=-75.24702782941932&b_br_lat=40.008575806662606&b_br_lng=-75.13587108344186)

      Hiking map example, taken from Alltrails.com Alltrails.com

1.1.2 Presentation maps

  • Used in journals and in outlets such as NY times and Wall Street
  • Designed to convey a very specific narrative of the author’s choosing
The New York Times, source <https://twitter.com/i/events/1319396419556544513>

The New York Times, source https://twitter.com/i/events/1319396419556544513

1.1.3 Statistical maps

  • Raw data is manipulated (merged, joined, cleanned) in such a way to tease out patterns otherwise not discernible in they original form
  • Sometimes benefit from being explored outside of a spatial context (regression tables, for example).

GIS deals more with the first form (reference maps), we focus more on presentation and statistical mapping here.

2 Geospatial Point Density (Evangelista and Beskow (2018))

Example taken from Vizual-Statistix

(kernels <- eval(formals(density.default)$kernel))
## [1] "gaussian"     "epanechnikov" "rectangular"  "triangular"   "biweight"    
## [6] "cosine"       "optcosine"
plot (density(0, bw = 1), xlab = "",
      main = "R's density() kernels with bw = 1")
for(i in 2:length(kernels))
   lines(density(0, bw = 1, kernel =  kernels[i]), col = i)
legend(1.5,.4, legend = kernels, col = seq(kernels),
       lty = 1, cex = .8, y.intersp = 1)
Example Kernel density

Example Kernel density

For example, Gimond (2017) shows the following example:

Region and point occurrence, see @gimond2017intro

Region and point occurrence, see Gimond (2017)

Basic (unweighted) Kernel point occurrence, see @gimond2017intro

Basic (unweighted) Kernel point occurrence, see Gimond (2017)

However, kernel weights as a function of distance from events (points) can be applied

Weighted Kernel point occurence, see @gimond2017intro

Weighted Kernel point occurence, see Gimond (2017)

Nonetheless, even cells without occurrence, may have a non-zero value. + Referred to as the square mesh approach*

2.1 Alternative to square mesh

  • The pointdensity() algorithm returns values only for event locations
    • This improves clarity on both the location and density of event behavior.

2.2 Statistical approach

  • Point density requires two parameters
    • g, a grid size that represents the fraction of latitude and longitude degree separation in the grid or binning system that will be used for aggregation, and
    • r, a positive integer value that represents the radius in grid steps for the neighborhood of interest, which yields a complexity of O(n(2r)2).
  • Typically, (2r)2≪G.
    • Figure 1 compared 800,000 crime locations across a mesh that measured 800 x 750 = 600,000 = G
    • The mesh of points were separated approximately 0.1 km.
    • Figure 4 resulted from pointdensity() with a radius of 10 steps (1 kilometer total) and the same grid size of 0.1 km,
    • The local neighborhood was (2r)2=400 points (or 400≪600,000)
    • This reduced computational complexity achieved by the pointdensity() algorithm largely results from the advantage of assuming and managing an unbounded spatial extent
  • Common disadvantages of density approaches
    • Include unnecessary computations
    • Always require pre-existing knowledge regarding the spatial extent (see grid examples above)
      • For example, what is our criterion to establishing those quadrants of interest?
  • Since the pointdensity() algorithm processes data without the requirement for a spatial extent, it is possible to efficiently process data that includes significant spatial separation

3 Methodology

3.1 The pointdensity() algorithm

  • applies elementary geometry with a geospatial projection in order to build the point density measurements.
Longitudinal tolerance, see @evangelista2018geospatial

Longitudinal tolerance, see Evangelista and Beskow (2018)

  • The center of the circle represents the grid point (any point in the Cartesian plane) closest to an event or binned collection of events

  • r is the radius distance selected (also the radius of the sphere)

  • The density of every grid point within r will increase by the amount of density associated with the center grid point

  • The algorithm iterates north and south within r (hence the degrees are not only 45o) along the latitudes of the grid (recall latitudes are the horizontal lines) to calculate the tolerance t of longitudinal lines that are within r

    • The previous figure shows how the neighborhood radius, r, is projected across a square grid 45o from a grid point

  • To find t, the algorithm relies on a spherical Pythagorean theorem

    • Theorem applied to right triangles on spheres (see Veljan (2000))
    • Spherical triangle is any any 3-sided region enclosed by sides that are arcs of great circles.
      • a great circle on a sphere is any circle whose center coincides with the center of the sphere
      • The regular Pythagorean theorem is a special case of the spherical version, as r goes to infinity, the spherical geometry becomes more and more like regular planar geometry.

  • The binning tolerance ensures that only grid points within r increase in density

  • Additionally, binning reduces the computing cost significantly

    • By grouping occurrences into classes, binning maps n points to G grid points renders G≪n
      Example with three points, see @evangelista2018geospatial

      Example with three points, see Evangelista and Beskow (2018)

  • Three grid points

    • r is four grid steps

  • The figure shows the final density count for every point within the neighborhood radius of the three events

3.2 Current implementation

  • pointdensity() algorithm reduces a list of n points to a list of m grid points on a mesh (or network).
  • The original list of n points, sorted by their latitude and longitude, reduces to m points once every point is rounded to the nearest point on the mesh.
  • Since only areas with points (occurrences) exist, for each pf the m points on the network, an initial density of one or greater exists.

3.2.1 Density calculation function

  • calc_density()
    • starts with enumeration of all latitudes within the radius of the event point
    • For each latitude, there is a longitudinal distance that represents the tolerance of the radius.
      • This tolerance is computed using the spherical Pythagorean theorem
  • After some matrix algebra, we have that:
    • All points outside of the radius have a negative value.
    • Once all elements in the matrix are set to zero and all positive ones are set to 1, we can multiply this binary matrix by the original matrix (Hadamard product)
Example Hadamard product, see [Hadamard Product (Element-wise Multiplication section here)](https://medium.com/linear-algebra/part-14-dot-and-hadamard-product-b7e0723b9133)

Example Hadamard product, see Hadamard Product (Element-wise Multiplication section here)

  • The resulting matrix contains the longitudinal value when there is an event within r
  • If there is temporal element, we can add these non-zero instances as a count

4 Application

We will rely on the following

# install.packages("pointdensityP")
library(ggmap)
library(KernSmooth)
library(pointdensityP)
## Warning: package 'pointdensityP' was built under R version 4.0.3
Comparison kernel VS pointdensity

Comparison kernel VS pointdensity

Data available here

SD<-read.table("incidents-5y.csv", sep = ",", header = TRUE)

4.1 ggmap version

The approach in the paper uses a now for profit service, we will stick with freeware

bbox <- c(left = -117.45, bottom = 32.52, right = -116.66, top = 33.39) #SD
sd <- get_stamenmap(bbox, zoom = 12, maptype = "toner-lite")
map_base<-sd

map_base<-ggmap(sd, extent = 'device', darken = c(.01, "black")) +
   theme(legend.key.size = grid::unit(1.2,"lines"),
         legend.title = element_text(size = 16, face = "bold"),
         legend.text = element_text(size = 14))
  • After this we are ready
SD_density <- pointdensity(df = SD, lat_col = "lat", lon_col = "lon", date_col = "date", grid_size = 0.1, radius = 1)
## 
## The radius was adjusted to  1.0008 km in order to accomodate the grid size
## 
## 
## The grid size is  0.001  measured in degrees
## 
## There are  47096  unique grids that require  17001656  measurements...
## 
## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |                                                                      |   1%
  |                                                                            
  |=                                                                     |   1%
  |                                                                            
  |=                                                                     |   2%
  |                                                                            
  |==                                                                    |   2%
  |                                                                            
  |==                                                                    |   3%
  |                                                                            
  |===                                                                   |   4%
  |                                                                            
  |===                                                                   |   5%
  |                                                                            
  |====                                                                  |   5%
  |                                                                            
  |====                                                                  |   6%
  |                                                                            
  |=====                                                                 |   7%
  |                                                                            
  |=====                                                                 |   8%
  |                                                                            
  |======                                                                |   8%
  |                                                                            
  |======                                                                |   9%
  |                                                                            
  |=======                                                               |   9%
  |                                                                            
  |=======                                                               |  10%
  |                                                                            
  |=======                                                               |  11%
  |                                                                            
  |========                                                              |  11%
  |                                                                            
  |========                                                              |  12%
  |                                                                            
  |=========                                                             |  12%
  |                                                                            
  |=========                                                             |  13%
  |                                                                            
  |==========                                                            |  14%
  |                                                                            
  |==========                                                            |  15%
  |                                                                            
  |===========                                                           |  15%
  |                                                                            
  |===========                                                           |  16%
  |                                                                            
  |============                                                          |  17%
  |                                                                            
  |============                                                          |  18%
  |                                                                            
  |=============                                                         |  18%
  |                                                                            
  |=============                                                         |  19%
  |                                                                            
  |==============                                                        |  19%
  |                                                                            
  |==============                                                        |  20%
  |                                                                            
  |==============                                                        |  21%
  |                                                                            
  |===============                                                       |  21%
  |                                                                            
  |===============                                                       |  22%
  |                                                                            
  |================                                                      |  22%
  |                                                                            
  |================                                                      |  23%
  |                                                                            
  |================                                                      |  24%
  |                                                                            
  |=================                                                     |  24%
  |                                                                            
  |=================                                                     |  25%
  |                                                                            
  |==================                                                    |  25%
  |                                                                            
  |==================                                                    |  26%
  |                                                                            
  |===================                                                   |  27%
  |                                                                            
  |===================                                                   |  28%
  |                                                                            
  |====================                                                  |  28%
  |                                                                            
  |====================                                                  |  29%
  |                                                                            
  |=====================                                                 |  29%
  |                                                                            
  |=====================                                                 |  30%
  |                                                                            
  |=====================                                                 |  31%
  |                                                                            
  |======================                                                |  31%
  |                                                                            
  |======================                                                |  32%
  |                                                                            
  |=======================                                               |  32%
  |                                                                            
  |=======================                                               |  33%
  |                                                                            
  |=======================                                               |  34%
  |                                                                            
  |========================                                              |  34%
  |                                                                            
  |========================                                              |  35%
  |                                                                            
  |=========================                                             |  35%
  |                                                                            
  |=========================                                             |  36%
  |                                                                            
  |==========================                                            |  37%
  |                                                                            
  |==========================                                            |  38%
  |                                                                            
  |===========================                                           |  38%
  |                                                                            
  |===========================                                           |  39%
  |                                                                            
  |============================                                          |  39%
  |                                                                            
  |============================                                          |  40%
  |                                                                            
  |============================                                          |  41%
  |                                                                            
  |=============================                                         |  41%
  |                                                                            
  |=============================                                         |  42%
  |                                                                            
  |==============================                                        |  42%
  |                                                                            
  |==============================                                        |  43%
  |                                                                            
  |==============================                                        |  44%
  |                                                                            
  |===============================                                       |  44%
  |                                                                            
  |===============================                                       |  45%
  |                                                                            
  |================================                                      |  45%
  |                                                                            
  |================================                                      |  46%
  |                                                                            
  |=================================                                     |  47%
  |                                                                            
  |=================================                                     |  48%
  |                                                                            
  |==================================                                    |  48%
  |                                                                            
  |==================================                                    |  49%
  |                                                                            
  |===================================                                   |  49%
  |                                                                            
  |===================================                                   |  50%
  |                                                                            
  |===================================                                   |  51%
  |                                                                            
  |====================================                                  |  51%
  |                                                                            
  |====================================                                  |  52%
  |                                                                            
  |=====================================                                 |  52%
  |                                                                            
  |=====================================                                 |  53%
  |                                                                            
  |=====================================                                 |  54%
  |                                                                            
  |======================================                                |  54%
  |                                                                            
  |======================================                                |  55%
  |                                                                            
  |=======================================                               |  55%
  |                                                                            
  |=======================================                               |  56%
  |                                                                            
  |========================================                              |  56%
  |                                                                            
  |========================================                              |  57%
  |                                                                            
  |========================================                              |  58%
  |                                                                            
  |=========================================                             |  58%
  |                                                                            
  |=========================================                             |  59%
  |                                                                            
  |==========================================                            |  59%
  |                                                                            
  |==========================================                            |  60%
  |                                                                            
  |==========================================                            |  61%
  |                                                                            
  |===========================================                           |  61%
  |                                                                            
  |===========================================                           |  62%
  |                                                                            
  |============================================                          |  62%
  |                                                                            
  |============================================                          |  63%
  |                                                                            
  |=============================================                         |  64%
  |                                                                            
  |=============================================                         |  65%
  |                                                                            
  |==============================================                        |  65%
  |                                                                            
  |==============================================                        |  66%
  |                                                                            
  |===============================================                       |  66%
  |                                                                            
  |===============================================                       |  67%
  |                                                                            
  |===============================================                       |  68%
  |                                                                            
  |================================================                      |  68%
  |                                                                            
  |================================================                      |  69%
  |                                                                            
  |=================================================                     |  69%
  |                                                                            
  |=================================================                     |  70%
  |                                                                            
  |=================================================                     |  71%
  |                                                                            
  |==================================================                    |  71%
  |                                                                            
  |==================================================                    |  72%
  |                                                                            
  |===================================================                   |  72%
  |                                                                            
  |===================================================                   |  73%
  |                                                                            
  |====================================================                  |  74%
  |                                                                            
  |====================================================                  |  75%
  |                                                                            
  |=====================================================                 |  75%
  |                                                                            
  |=====================================================                 |  76%
  |                                                                            
  |======================================================                |  76%
  |                                                                            
  |======================================================                |  77%
  |                                                                            
  |======================================================                |  78%
  |                                                                            
  |=======================================================               |  78%
  |                                                                            
  |=======================================================               |  79%
  |                                                                            
  |========================================================              |  79%
  |                                                                            
  |========================================================              |  80%
  |                                                                            
  |========================================================              |  81%
  |                                                                            
  |=========================================================             |  81%
  |                                                                            
  |=========================================================             |  82%
  |                                                                            
  |==========================================================            |  82%
  |                                                                            
  |==========================================================            |  83%
  |                                                                            
  |===========================================================           |  84%
  |                                                                            
  |===========================================================           |  85%
  |                                                                            
  |============================================================          |  85%
  |                                                                            
  |============================================================          |  86%
  |                                                                            
  |=============================================================         |  87%
  |                                                                            
  |=============================================================         |  88%
  |                                                                            
  |==============================================================        |  88%
  |                                                                            
  |==============================================================        |  89%
  |                                                                            
  |===============================================================       |  89%
  |                                                                            
  |===============================================================       |  90%
  |                                                                            
  |===============================================================       |  91%
  |                                                                            
  |================================================================      |  91%
  |                                                                            
  |================================================================      |  92%
  |                                                                            
  |=================================================================     |  92%
  |                                                                            
  |=================================================================     |  93%
  |                                                                            
  |==================================================================    |  94%
  |                                                                            
  |==================================================================    |  95%
  |                                                                            
  |===================================================================   |  95%
  |                                                                            
  |===================================================================   |  96%
  |                                                                            
  |====================================================================  |  97%
  |                                                                            
  |====================================================================  |  98%
  |                                                                            
  |===================================================================== |  98%
  |                                                                            
  |===================================================================== |  99%
  |                                                                            
  |======================================================================|  99%
  |                                                                            
  |======================================================================| 100%done...
# SD_density$count[SD_density$count>10000] <- 10000 
## creates discriminating scale
# png("SD_pointdensity_test.png", width = 1000, height = 1000, units = "px")
map_base + geom_point(aes(x = lon, y = lat, colour = count), shape = 16, size = 0.5, data = SD_density) + scale_colour_gradient(low = "green", high = "red") + labs(color = "density count\n2km radius\n") + theme(legend.position = c(0.1, 0.2), legend.background = element_rect(fill = NA), legend.key.size = unit(.5, "cm"),  legend.text = element_text(size = 8), legend.title = element_text(size = 10))
## Warning: Removed 74030 rows containing missing values (geom_point).
Replicating Figure 4 in paper

Replicating Figure 4 in paper 

# dev.off()

4.2 t-map version

library(classInt)
library(sf)
library(spdep)
library(tigris)
options(tigris_use_cache = TRUE)
library(acs)
library(stringr)
library(tmaptools)
library(tmap)
library(plyr)
library(viridis)

zip<-zctas(starts_with = c("90", "91", "92", "93", "94", "95", "96"), class="sp")
## Warning in proj4string(obj): CRS object has comment, which is lost in output
zip$lat<-as.numeric(as.character(zip$INTPTLAT10))
zip$lon<-as.numeric(as.character(zip$INTPTLON10))
zip<-zip[zip$lon > -117.45 &
                  zip$lon < -116.66 &
                  zip$lat > 32.52 &
                  zip$lat < 33.39, ]

projcrs <- "+proj=longlat +datum=NAD83 +no_defs"
crimespoints <- st_as_sf(x = SD_density,
           coords = c("lon", "lat"),
           crs = projcrs)

library(classInt)
class <- classIntervals(SD_density$count, 9, style = "quantile")    
class
## style: quantile
##      [0,129)    [129,748)   [748,1268)  [1268,1760)  [1760,2364)  [2364,3161) 
##        88532        88725        88562        88760        88251        88997 
##  [3161,4646)  [4646,7873) [7873,31212] 
##        88784        88668        88699
tmap_mode("plot")
tm_shape(zip) +
    tm_polygons() +
    tm_layout(bg.color = "grey51",
    title = "Crime distribution over five years",
    title.position = c("right", "top"), title.size = 1.1, title.color = "white",
    legend.position = c("left", "bottom"), legend.text.size = 0.85,
    legend.width = 0.25, legend.text.color = "white", legend.title.color="white") +
    # tm_credits("Data source: ACS\n*Missing values contain ZCTAs with no information",
    # position = c(0.25, 0.02), size = .75, col="white")+
    tm_borders(col=rgb(31,31,31,max=250,250/3)) +
    tm_shape(crimespoints) + 
    tm_dots(col="count", palette = "inferno", size=.025, title = "Crimes", 
            style = "quantile", n = 9, legend.show = FALSE) +
    tm_add_legend('fill', 
    col = viridis::viridis(9, alpha = 1, begin = 0, end = 1, direction = 1, option = "B"),
    border.col = "grey40",
    size = 1, labels = c("[0,129)", "[129,748)",   "[748,1268)", "[1268,1760)",  "[1760,2364)",  
                         "[2364,3161)",  "[3161,4646)",  "[4646,7873)", "[7873,31212]"),
    title="Crimes")
## Warning in sp::proj4string(obj): CRS object has comment, which is lost in output
## Warning: One tm layer group has duplicated layer types, which are omitted. To
## draw multiple layers of the same type, use multiple layer groups (i.e. specify
## tm_shape prior to each of them).
Replicating Figure 4 in paper with `tmap`

Replicating Figure 4 in paper with tmap

References

Evangelista, Paul F, and David Beskow. 2018. “Geospatial Point Density.” R Journal 10 (2). https://journal.r-project.org/archive/2018/RJ-2018-061/index.html.

Gimond, Manuel. 2017. “Intro to Gis and Spatial Analysis.” https://mgimond.github.io/Spatial/index.html.

Veljan, Darko. 2000. “The 2500-Year-Old Pythagorean Theorem.” Mathematics Magazine 73 (4): 259–72.