1

Import the corona_cologne dataset into your R session.
There are several corona_cologne files in the ../data/ folder, you would only need the one with the .shp extension.

2

Use the same command as before, but instead of defining the path to the file location plug in the following string:

"https://geoportal.stadt-koeln.de/arcgis/rest/services/Politik_und_Verwaltung/covid_stadtteile/MapServer/1/query?where=id+is+not+null&text=&objectIds=&time=&geometry=&geometryType=esriGeometryEnvelope&inSR=&spatialRel=esriSpatialRelIntersects&distance=&units=esriSRUnit_Foot&relationParam=&outFields=*&returnGeometry=true&returnTrueCurves=false&maxAllowableOffset=&geometryPrecision=&outSR=4326&havingClause=&returnIdsOnly=false&returnCountOnly=false&orderByFields=&groupByFieldsForStatistics=&outStatistics=&returnZ=false&returnM=false&gdbVersion=&historicMoment=&returnDistinctValues=false&resultOffset=&resultRecordCount=&returnExtentOnly=false&datumTransformation=&parameterValues=&rangeValues=&quantizationParameters=&featureEncoding=esriDefault&f=geojson"

Have a look at the results. What do you think has happened?
If you don’t like long strings occupying a lot of horizontal screen space, consider cutting it in pieces. I prefer using the paste()/paste0() function or functions from additional packages such as glue::glue() from the glue package.

3

Transform the CRS of the dataset to EPSG:3035 using the sf::st_transform() function.
You only need the four digit number for defining the CRS.

4

Load both the immigrants_cologne and the inhabitants_cologne raster datasets in your R session.
You would need the stars::read_stars() function.

5

Aggregate the number of immigrants in each district using sum as function. Does it work?
It’s the same function as on slide XXX.

Let’s see what is happening here by comparing the exact CRS strings:

sf::st_crs(corona_cologne)
## Coordinate Reference System:
##   User input: EPSG:3035 
##   wkt:
## PROJCRS["ETRS89-extended / LAEA Europe",
##     BASEGEOGCRS["ETRS89",
##         DATUM["European Terrestrial Reference System 1989",
##             ELLIPSOID["GRS 1980",6378137,298.257222101,
##                 LENGTHUNIT["metre",1]]],
##         PRIMEM["Greenwich",0,
##             ANGLEUNIT["degree",0.0174532925199433]],
##         ID["EPSG",4258]],
##     CONVERSION["Europe Equal Area 2001",
##         METHOD["Lambert Azimuthal Equal Area",
##             ID["EPSG",9820]],
##         PARAMETER["Latitude of natural origin",52,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8801]],
##         PARAMETER["Longitude of natural origin",10,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8802]],
##         PARAMETER["False easting",4321000,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8806]],
##         PARAMETER["False northing",3210000,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8807]]],
##     CS[Cartesian,2],
##         AXIS["northing (Y)",north,
##             ORDER[1],
##             LENGTHUNIT["metre",1]],
##         AXIS["easting (X)",east,
##             ORDER[2],
##             LENGTHUNIT["metre",1]],
##     USAGE[
##         SCOPE["Statistical analysis."],
##         AREA["Europe - European Union (EU) countries and candidates. Europe - onshore and offshore: Albania; Andorra; Austria; Belgium; Bosnia and Herzegovina; Bulgaria; Croatia; Cyprus; Czechia; Denmark; Estonia; Faroe Islands; Finland; France; Germany; Gibraltar; Greece; Hungary; Iceland; Ireland; Italy; Kosovo; Latvia; Liechtenstein; Lithuania; Luxembourg; Malta; Monaco; Montenegro; Netherlands; North Macedonia; Norway including Svalbard and Jan Mayen; Poland; Portugal including Madeira and Azores; Romania; San Marino; Serbia; Slovakia; Slovenia; Spain including Canary Islands; Sweden; Switzerland; Turkey; United Kingdom (UK) including Channel Islands and Isle of Man; Vatican City State."],
##         BBOX[24.6,-35.58,84.17,44.83]],
##     ID["EPSG",3035]]
sf::st_crs(immigrants_cologne)
## Coordinate Reference System:
##   User input: ETRS89-extended / LAEA Europe (with axis order normalized for visualization) 
##   wkt:
## PROJCRS["ETRS89-extended / LAEA Europe (with axis order normalized for visualization)",
##     BASEGEOGCRS["ETRS89",
##         DATUM["European Terrestrial Reference System 1989",
##             ELLIPSOID["GRS 1980",6378137,298.257222101004,
##                 LENGTHUNIT["metre",1]]],
##         PRIMEM["Greenwich",0,
##             ANGLEUNIT["degree",0.0174532925199433]],
##         ID["EPSG",4258]],
##     CONVERSION["Lambert Azimuthal Equal Area",
##         METHOD["Lambert Azimuthal Equal Area",
##             ID["EPSG",9820]],
##         PARAMETER["Latitude of natural origin",52,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8801]],
##         PARAMETER["Longitude of natural origin",10,
##             ANGLEUNIT["degree",0.0174532925199433],
##             ID["EPSG",8802]],
##         PARAMETER["False easting",4321000,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8806]],
##         PARAMETER["False northing",3210000,
##             LENGTHUNIT["metre",1],
##             ID["EPSG",8807]]],
##     CS[Cartesian,2],
##         AXIS["(E)",east,
##             ORDER[1],
##             LENGTHUNIT["metre",1]],
##         AXIS["(N)",north,
##             ORDER[2],
##             LENGTHUNIT["metre",1]]]

Indeed, they are a bit different. We should change that.

6

Transform the raster data (again) to EPSG:3035.
stars and sf objects talk fluently to each other. You could either use sf::st_transform or stars::st_transform_proj() for this task. Try it out!

7

Aggregate the number of immigrants and inhabitants in each district using sum as function. Does it work now?
It’s the same function as on slide XXX.

8

Create an immigrant share variable from the two aggregated layers and add them to the corona_cologne dataset. Plot it if you will.

You can either do it on-the-fly using dplyr::mutate() as on the slides (no. 33), or create first a new object and then add it using base-R (corona_cologne$immigrant_share <- ...). You could even combine the approaches.

Either way, what might be a bit cumbersome is that you have to extract the first list element in each of the raster layers (e.g., immigrants_count[[1]])