Meta Polys in Massachussetts

Meta Polys in Sample Areas

This shows prototype generation for “meta polys” for two sample commuting zones in MA, Springfield and Boston czs. I compound most of the division type we used, turn them all into cleaned lines, and then get polygonal subdivisions of the commuting zones. These sample measures are generated with: - Places - School Districts - County lines - Highways (all limited-access) - Rails

For water, I had previously used a shapefile that is too-low resolution for generating alongside these. I could use the census tigris water, which is ultra high detail, but haven’t included yet. (The very high detail also makes it somewhat slow to work with, even over just one cz)

Quick output plots and descriptives are included in code, with interactive maps at the end.

suppressMessages(library(sf))
suppressMessages(library(dplyr))
suppressMessages(library(lwgeom))
suppressMessages(library(mapview))
# aggregated shapefiles and did some cleaning; not showing code in markdown.
load("prepped-for-sample-areas.RData")

# generate for springfield ------------------------------------------------

spf.divs <- st_intersection( all.divs,
                             st_buffer(tmp.czs[2, ]
                                       , 100) )

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

spf.divs <- st_collection_extract(spf.divs
                                  , "LINESTRING")


spf.polys <- polygonal.div(   spf
                              ,spf.divs
                              ,return.sf = T
                              ,min.size = 1e4 # m^2 (very small)
                              ,min.population.count = 10 # 30 # even 30 will trim out some suburban enclaves in this example..
                              )

## generating cz - 20800

## Registered S3 method overwritten by 'geojsonlint':
##   method         from 
##   print.location dplyr

## Warning in class(input) == c("sf", "tbl_df", "tbl", "data.frame"): longer object
## length is not a multiple of shorter object length

spf.polys$area <- st_area(spf.polys$geometry)
plot(spf.polys['id'], main = "Springfield polygons")

# summaries of sizes/interpolated populations for each sub-polygon
quantile(spf.polys$population, seq(0,1, .1)  )

##          0%         10%         20%         30%         40%         50% 
##    10.12414    29.50247    62.44928   177.60450   417.02814   694.85167 
##         60%         70%         80%         90%        100% 
##  1047.86730  1654.97523  2464.66699  4776.07698 93713.27637

quantile(spf.polys$area, seq(0,1, .1)  )

## Units: [m^2]
##           0%          10%          20%          30%          40%          50% 
##     35310.49    244195.05    585306.23   1297815.81   1907586.49   3407101.95 
##          60%          70%          80%          90%         100% 
##   6095125.50  11919564.32  22548328.08  41925536.55 210025475.74

tibble(spf.polys)[,c("population", "pop.perc", "area")] %>% summary()

##    population          pop.perc              area          
##  Min.   :   10.12   Min.   :1.461e-05   Min.   :    35310  
##  1st Qu.:  104.87   1st Qu.:1.513e-04   1st Qu.:   922499  
##  Median :  694.85   Median :1.003e-03   Median :  3407102  
##  Mean   : 2085.27   Mean   :3.009e-03   Mean   : 14648249  
##  3rd Qu.: 1941.52   3rd Qu.:2.802e-03   3rd Qu.: 16781705  
##  Max.   :93713.28   Max.   :1.352e-01   Max.   :210025476

max(spf.polys$id) # total polygonal division (given parameters)

## [1] 331

spf.divs <- rmapshaper::ms_explode(spf.divs)

# generate for boston -----------------------------------------------------

boston.divs <- st_intersection( all.divs,
                             st_buffer(tmp.czs[1,]
                                       , 100) )

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

boston.divs <- st_collection_extract(boston.divs
                                     , "LINESTRING")

boston.polys <- polygonal.div(   boston
                                ,boston.divs
                                ,return.sf = T
                                ,min.size = 1e4
                                ,min.population.count = 10  )

## generating cz - 20500

## Warning in class(input) == c("sf", "tbl_df", "tbl", "data.frame"): longer object
## length is not a multiple of shorter object length

boston.polys$area <- st_area(boston.polys$geometry)
plot(boston.polys['id'], main="Boston polygons")

# summaries of sizes/interpolated populations for each sub-polygon
quantile(boston.polys$population, seq(0,1, .1)  )

##           0%          10%          20%          30%          40%          50% 
##     10.40092     54.90010    150.65654    319.05581    591.98870   1012.91529 
##          60%          70%          80%          90%         100% 
##   1717.14275   2781.19595   4346.26002   8270.24110 155253.25151

quantile(boston.polys$area, seq(0,1, .1)  )

## Units: [m^2]
##           0%          10%          20%          30%          40%          50% 
## 1.009143e+04 9.347098e+04 3.166344e+05 7.928928e+05 1.550009e+06 2.721366e+06 
##          60%          70%          80%          90%         100% 
## 4.414009e+06 6.570189e+06 1.080726e+07 2.121161e+07 1.934246e+09

tibble(boston.polys)[,c("population", "pop.perc", "area")] %>% summary()

##    population          pop.perc              area          
##  Min.   :    10.4   Min.   :2.020e-06   Min.   :1.009e+04  
##  1st Qu.:   237.5   1st Qu.:4.613e-05   1st Qu.:5.538e+05  
##  Median :  1012.9   Median :1.967e-04   Median :2.721e+06  
##  Mean   :  3256.7   Mean   :6.326e-04   Mean   :9.218e+06  
##  3rd Qu.:  3409.7   3rd Qu.:6.623e-04   3rd Qu.:8.632e+06  
##  Max.   :155253.2   Max.   :3.016e-02   Max.   :1.934e+09

max(boston.polys$id) # total polygonal division (given parameters)

## [1] 1578

boston.divs <- rmapshaper::ms_explode(boston.divs)

Meta Polys in Massachussetts

kira

10/20/2020

Meta Polys in Sample Areas

Mapview