Code not running - Optional Lab #6

[Niagara1000]

Hi Dr. @lecy ,

I am trying to run the bolded code below:

library( geojsonio )   # read shapefiles
library( sp )          # work with shapefiles
library( sf )          # work with shapefiles - simple features format
library( tmap )        # theme maps
library( dplyr )       # data wrangling
library( pander )      # nice tables 

github <- "https://raw.githubusercontent.com/benbalter/dc-maps/master/maps/2006-traffic-volume.geojson"
**traffic <- geojson_read(x="https://raw.githubusercontent.com/benbalter/dc-maps/master/maps/2006-traffic-volume.geojson", method="external", what="sp")**
sp <- sf::st_read(github)

plot( traffic, col="steelblue" )

But, the traffic variable is not getting created. it returns the following error:

Error in st_sf(x, ..., agr = agr, sf_column_name = sf_column_name) : no simple features geometry column present

When I click on "Show Traceback" at the right of the output, I get this:

What should I do to fix this?

Thank you!

[Niagara1000]

Dr. @lecy ,

I found how to fix the above error. I had to remove the method argument from the function.

My next question is:

In the instructions for Lab 6, it says " Instructions: Create an R script that will convert all US Metro Area shapefiles into Dorling cartograms, one new shapefile for each metro area. "

Q1. Where do I find all the US Metro Area Shapefiles? Q2. What is the meaning of " .. cartograms, one new shapefile for .. " ? If we are to convert the shapefiles into cartograms, then why do we need one new shapefile for each metro area?

Thank you!

[lecy]

So this worked for you?

github <- "https://raw.githubusercontent.com/benbalter/dc-maps/master/maps/2006-traffic-volume.geojson"
traffic <- geojson_read(x="https://raw.githubusercontent.com/benbalter/dc-maps/master/maps/2006-traffic-volume.geojson",  what="sp")
sp <- sf::st_read(github)

plot( traffic, col="steelblue" )

I recall that they did a big update the geojsonio package recently, so they may have deprecated that argument.

The rest of the code is from Lab-02 in CPP 529. If you have not taken that course yet it might be harder to follow, but the get_acs() function downloads both the population for the MSA tracts as well as the shapefiles.

There is conversion between "simple features" format (sf) and the R spatial format (sp) in order for the dorling code to work.

library( geojsonio )   # read shapefiles
library( sp )          # work with shapefiles
library( sf )          # work with shapefiles - simple features format
library( tmap )        # theme maps
library( dplyr )       # data wrangling
library( pander )      # nice tables 

crosswalk <- "https://raw.githubusercontent.com/DS4PS/cpp-529-master/master/data/cbsatocountycrosswalk.csv"
crosswalk <- read.csv( crosswalk, stringsAsFactors=F, colClasses="character" )

# search for citie names by strings, use the ^ anchor for "begins with" 
grep( "^MIN", crosswalk$msaname, value=TRUE ) 

# select all FIPS for Minneapolis
these.minneapolis <- crosswalk$msaname == "MINNEAPOLIS-ST. PAUL, MN-WI"
these.fips <- crosswalk$fipscounty[ these.minneapolis ]
these.fips <- na.omit( these.fips )

state.fips <- substr( these.fips, 1, 2 )
county.fips <- substr( these.fips, 3, 5 )

dat <- data.frame( name="MINNEAPOLIS-ST. PAUL, MN-WI",
                   state.fips, county.fips, fips=these.fips )               

library( tidycensus )

# census_api_key("YOUR KEY GOES HERE")
# key <- "abc123"
# census_api_key( key )

# Minneapolis metro area spans two states - 
# Minnesota = 27
# Wisconsin = 55

msp.pop1 <-
get_acs( geography = "tract", variables = "B01003_001",
         state = "27", county = county.fips[state.fips=="27"], geometry = TRUE ) %>% 
         select( GEOID, estimate ) %>%
         rename( POP=estimate )

msp.pop2 <-
get_acs( geography = "tract", variables = "B01003_001",
         state = "55", county = county.fips[state.fips=="55"], geometry = TRUE ) %>% 
         select( GEOID, estimate ) %>%
         rename( POP=estimate )

msp.pop <- rbind( msp.pop1, msp.pop2 )

plot( msp.pop )

# convert sf map object to an sp version
msp.sp <- as_Spatial( msp )
class( msp.sp )

# project map and remove empty tracts
msp.sp <- spTransform( msp.sp, CRS("+init=epsg:3395"))
msp.sp <- msp.sp[ msp.sp$POP != 0 & (! is.na( msp.sp$POP )) , ]

# convert census tract polygons to dorling cartogram
# no idea why k=0.03 works, but it does - default is k=5
msp.sp$pop.w <- msp.sp$POP / 9000 # max(msp.sp$POP)   # standardizes it to max of 1.5
msp_dorling <- cartogram_dorling( x=msp.sp, weight="pop.w", k=0.05 )
plot( msp_dorling )