Doel
In deze blog wil ik laten zien hoe je tegenwoordig geografische kaarten in R kunt maken. Het pakket sf (Simple Features) is de moderne standaard daarvoor en is de natuurlijke opvolger van sp. Dit sf pakket kan goed met het databewerkingspakket tidyverse werken en met onderdelen daarvan zoals dplyr en het visualisatiepakket ggplot. Een pakket waar je, zeker in combinatie met sf, goed kaarten mee kunt maken is tmap (statische en interactieve kaarten). Ik had mij bij het maken van deze blog de taak gesteld om te laten zien hoe je Nederland geografisch af kunt beelden. Dat laat ik hieronder zien en daar gebruik ik die moderne R-pakketten voor. Ik kwam er bij het maken ook achter dat je van de open data van het CBS ook makkelijk kaarten van Nederland kunt maken. Dat laat ik je hieronder ook zien. In een latere fase kom ik hier nog eens op terug.
Pakketten
Eerst maar eens enkele moderne R-pakketten binnenhalen.
Geographic Data Science (GDS) en Geographic Information Science (GISc)
In een mooie blog laat Katie Jolly het verschil zien tussen Geografische Datawetenschap en Geografische Informatiewetenschap Blog van Katie Jolly. Bij R en Python praat je meer over geografische datawetenschap of geocomputation, zoals Lovelace e.a. dat noemen Geocomputatie. Daarnaast heb je Geografische Informatie Systemen zoals QGIS en ArcMap as GISc (vaak GIS-programma’s genoemd). Een van de grote verschillen, en daar heeft ze gelijk in, is de de reproduceerbaarheid van Geografische Datawetenschap ten opzichte van GIS-programma’s. Met het volgende overzicht, dat zijn ontleent aan Lovelace e.a., maakt zij dat duidelijk: 
Deze tabel komt uit Robin Lovelace’s blog-post Can Geographic Data Save the World?.
Met sf in combinatie met tidyverse en enkele specifieke pakketten kom je een heel eind tegenwoordig en kun je dat doen wat dure GIS-programma’s ook kunnen en soms nog meer.
R’s ruimtelijke ecosysteem
Er is tegenwoordig een hele variëteit aan ruimtelijke data-analyse pakketten beschikbaar.
- sp: Classes and Methods for Spatial Data
- sf: Simple Features for R (bouwt voort op sp)
- spdep: Spatial Dependence: Weighting Schemes, Statistics, and Models
- lwgeom: Binding to the liblwgeom library
- …
De nieuwste en de beste op dit moment is sf. Over ruimtelijk analyseren is een hele serie tutorials beschikbaar op internet, ook al lopen ze in de tijd wat achter, is te vinden op de r-spatial site. Daarnaast is er een hele serie specieke pakketten voor data-visualisatie waar je mee kunt werken en die de gebruiker veel flexibiliteit bieden. Terwijl sf zelf heel veel kan, werkt het heel goed in combinatie met deze twee pakketten.
De laatste maand heb ik wat gewerkt met sf in combinatie met ggplot en tmap. Ik wilde mij dit wat beter eigen maken en het gebruiken voor het maken van kaarten van Nederland en delen van Nederland. Hieronder vind je achtereenvolgens wat informatie over sf en ggplot en daarna iets over tamp.
Structure of sf data
Een figuur uit Geocomputation with R, zie ook mijn vorig blog Geocomputation with R.
Simple Features is een hierarchisch data model dat een behoorlijk breed palet aan geometrische types representeert. Het komt erop neer dat simple features een dataframe is met in elke rij ruimtelijke gegevens (een bevolkingsgegeven, een punt, een stad, …) met een list-kolom met coördinaten waarmee de geografische vorm gemaakt kan worden. Hieronder leer je hoe het werkt.
Projecties
Het is goed om te weten dat projecties een 3D oppervlakte naar een 2D-oppervlakte omvormen. Verschillende projecties laten de geografische vorm er heel anders uit zien.
Open data file binnen tmap
Laten we eerst eens een databestand van Nederland binnenhalen dat in het pakket tmap opgeslagen zit. Dit is een databestand met de steden van Nederland en enkele gegevens van deze steden. Onderstaande werkt dus alleen als je de library(tmap) hebt geopend.
Laten we met sf eens kijken wat we in huis hebben gehaald. We tonen enkele plot mogelijkheden van sf. Eerst maar eens de grove kaart van ons land.
Welke data zitten er in het bestand, in combinatie met de geografische kaders? In het databestand zitten de codes, de namen, de provincies, de populatie, populatie_man, populatie_vrouw, populaties over verschillende leeftijdsgroepen etc.
Stel dat we maximaal 15 variabelen willen afdrukken.
Of alleen maar een, de eerste.
We kunnen de gegevens nu ook goed vanuit sf met ggplot afdrukken. Stel bijvoorbeeld dat we de populatie van 0-14 jaar (een variabele die in het databestand zit) willen laten zien en hoe dat percentage in de verschillende steden in Nederland eruit ziet. We breken de percentages op in van 14 tot 40 procent (je kunt zien dat de percentages daar tussen liggen).
Omdat het met tidyverse, en dus met dplyr werkt, kunnen we ook de standaard databewerkingscodes gebruiken. Dat maakt het werken met geografische data een stuk makkelijker. Als we alleen Zuid-Holland willen laten zien, maken we het databestand Zuid-Holland.
Dan kunnen we hier de eerste variabele laten zien en dan zien we de gemeenten van deze provincie afgebeeld.
En hetzelfde als hierboven. Hoe zit het met de jonge bevolking in de ZuidHollandse steden? We zien het hieronder met inzet van het pakket ggplot.
Overstappen naar tmap
Daar waar sf goed is voor het binnenhalen en bewerken van de data op allerlei manieren, daar is tmap heel goed in het maken van de kaarten op een eenvoudige manier. tmapis gemaakt om met grote flexibiliteit kaarten te kunnen maken. Het heeft dezelfde gelaagde structuur als ggplot. Je vindt op de website een document om er makkelijk mee te kunnen beginnentmap: get started. En de ontwikkelaar (Martijn Tennekes van het CBS ) heeft er een inzichtelijk artikel over geschreven in Journal of Statistical Software artike.
Stel dat je er een percentage bevolking aan wilt toevoegen en dit percentage per provincie wilt afbeelden. Met de volgende code doe je dat.
Stel dat we twee grafieken naast elkaar willen zetten. Dat doe je zo.
Met tmap kun je ook en net zo makkelijk interactieve kaarten maken.
Open data via CBS
Ik kwam er ook achter dat je kaarten ook kunt maken via de open-data mogelijkheden van het CBS (Het kan, maar het vraagt nog wel wat oefening de komnende tijd). Open hiervoor het pakket cbsodataR
Zoek vervolgens op welke data beschikbaar zijn:
[1] "WijkenEnBuurten"
[2] ""
[3] "Gemeentenaam_1"
[4] "SoortRegio_2"
[5] "Codering_3"
[6] "IndelingswijzigingWijkenEnBuurten_4"
[7] ""
[8] "AantalInwoners_5"
[9] ""
[10] "Mannen_6"
[11] "Vrouwen_7"
[12] ""
[13] "k_0Tot15Jaar_8"
[14] "k_15Tot25Jaar_9"
[15] "k_25Tot45Jaar_10"
[16] "k_45Tot65Jaar_11"
[17] "k_65JaarOfOuder_12"
[18] ""
[19] "Ongehuwd_13"
[20] "Gehuwd_14"
[21] "Gescheiden_15"
[22] "Verweduwd_16"
[23] ""
[24] "WestersTotaal_17"
[25] ""
[26] "NietWestersTotaal_18"
[27] "Marokko_19"
[28] "NederlandseAntillenEnAruba_20"
[29] "Suriname_21"
[30] "Turkije_22"
[31] "OverigNietWesters_23"
[32] ""
[33] "GeboorteTotaal_24"
[34] "GeboorteRelatief_25"
[35] "SterfteTotaal_26"
[36] "SterfteRelatief_27"
[37] ""
[38] "HuishoudensTotaal_28"
[39] "Eenpersoonshuishoudens_29"
[40] "HuishoudensZonderKinderen_30"
[41] "HuishoudensMetKinderen_31"
[42] "GemiddeldeHuishoudensgrootte_32"
[43] "Bevolkingsdichtheid_33"
[44] ""
[45] "Woningvoorraad_34"
[46] "GemiddeldeWoningwaarde_35"
[47] ""
[48] "PercentageEengezinswoning_36"
[49] "PercentageMeergezinswoning_37"
[50] ""
[51] "PercentageBewoond_38"
[52] "PercentageOnbewoond_39"
[53] ""
[54] "Koopwoningen_40"
[55] ""
[56] "HuurwoningenTotaal_41"
[57] "InBezitWoningcorporatie_42"
[58] "InBezitOverigeVerhuurders_43"
[59] "EigendomOnbekend_44"
[60] ""
[61] "BouwjaarVoor2000_45"
[62] "BouwjaarVanaf2000_46"
[63] ""
[64] ""
[65] "GemiddeldElektriciteitsverbruikTotaal_47"
[66] ""
[67] "Appartement_48"
[68] "Tussenwoning_49"
[69] "Hoekwoning_50"
[70] "TweeOnderEenKapWoning_51"
[71] "VrijstaandeWoning_52"
[72] ""
[73] "Huurwoning_53"
[74] "EigenWoning_54"
[75] ""
[76] "GemiddeldAardgasverbruikTotaal_55"
[77] ""
[78] "Appartement_56"
[79] "Tussenwoning_57"
[80] "Hoekwoning_58"
[81] "TweeOnderEenKapWoning_59"
[82] "VrijstaandeWoning_60"
[83] ""
[84] "Huurwoning_61"
[85] "EigenWoning_62"
[86] "PercentageWoningenMetStadsverwarming_63"
[87] ""
[88] ""
[89] "AantalInkomensontvangers_64"
[90] "GemiddeldInkomenPerInkomensontvanger_65"
[91] "GemiddeldInkomenPerInwoner_66"
[92] "k_40PersonenMetLaagsteInkomen_67"
[93] "k_20PersonenMetHoogsteInkomen_68"
[94] "Actieven1575Jaar_69"
[95] ""
[96] "k_40HuishoudensMetLaagsteInkomen_70"
[97] "k_20HuishoudensMetHoogsteInkomen_71"
[98] "HuishoudensMetEenLaagInkomen_72"
[99] "HuishOnderOfRondSociaalMinimum_73"
[100] ""
[101] "PersonenPerSoortUitkeringBijstand_74"
[102] "PersonenPerSoortUitkeringAO_75"
[103] "PersonenPerSoortUitkeringWW_76"
[104] "PersonenPerSoortUitkeringAOW_77"
[105] ""
[106] "BedrijfsvestigingenTotaal_78"
[107] ""
[108] "ALandbouwBosbouwEnVisserij_79"
[109] "BFNijverheidEnEnergie_80"
[110] "GIHandelEnHoreca_81"
[111] "HJVervoerInformatieEnCommunicatie_82"
[112] "KLFinancieleDienstenOnroerendGoed_83"
[113] "MNZakelijkeDienstverlening_84"
[114] "RUCultuurRecreatieOverigeDiensten_85"
[115] ""
[116] ""
[117] "PersonenautoSTotaal_86"
[118] ""
[119] "PersonenautoSJongerDan6Jaar_87"
[120] "PersonenautoS6JaarEnOuder_88"
[121] ""
[122] "PersonenautoSBrandstofBenzine_89"
[123] "PersonenautoSOverigeBrandstof_90"
[124] "PersonenautoSPerHuishouden_91"
[125] "PersonenautoSNaarOppervlakte_92"
[126] "Motorfietsen_93"
[127] ""
[128] "AfstandTotHuisartsenpraktijk_94"
[129] "AfstandTotGroteSupermarkt_95"
[130] "AfstandTotKinderdagverblijf_96"
[131] ""
[132] "AfstandTotSchool_97"
[133] "ScholenBinnen3Km_98"
[134] ""
[135] "OppervlakteTotaal_99"
[136] "OppervlakteLand_100"
[137] "OppervlakteWater_101"
[138] ""
[139] "MeestVoorkomendePostcode_102"
[140] "Dekkingspercentage_103"
[141] ""
[142] "MateVanStedelijkheid_104"
[143] "Omgevingsadressendichtheid_105"
[144] ""
[145] "TotaalDiefstalUitWoningSchuurED_106"
[146] "VernielingMisdrijfTegenOpenbareOrde_107"
[147] "GeweldsEnSeksueleMisdrijven_108" Gebruik de WijkenenBuurten-data en de GeboorteRelatief_25 en verwijder spaties uit regiocodes.
Haal de kaart met gemeentegrenzen op van PDOK
Reading layer `OGRGeoJSON' from data source `https://geodata.nationaalgeoregister.nl/cbsgebiedsindelingen/wfs?request=GetFeature&service=WFS&version=2.0.0&typeName=cbs_gemeente_2017_gegeneraliseerd&outputFormat=json' using driver `GeoJSON'
Simple feature collection with 388 features and 5 fields
geometry type: MULTIPOLYGON
dimension: XY
bbox: xmin: 13565.4 ymin: 306846.9 xmax: 277992.8 ymax: 619291
epsg (SRID): 28992
proj4string: +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +towgs84=565.2369,50.0087,465.658,-0.406857,0.350733,-1.87035,4.0812 +units=m +no_defsKoppel CBS-data aan geodata met regiocodes
Maak een thematische kaart
Ik zal kijken of ik de komende maanden een eenvoudige tutorial over dit onderwerp kan maken. Als ik deze klaar heb, kom ik terug op dit onderwerp.