<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/introducing-the-kernelheaping-package-ii/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://www.inwt-statistics.com/read-blog/introducing-the-kernelheaping-package-ii.html"> INWT-Blog-RBloggers</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<pre>In the first part of <a rel="nofollow" href="https://www.inwt-statistics.com/read-blog/introducing-the-kernelheaping-package-512.html" target="_blank">Introducing the Kernelheaping Package</a> I showed how to compute and plot kernel density estimates on rounded or interval censored data using the Kernelheaping package. Now, let’s make a big leap forward to the 2-dimensional case. Interval censoring can be generalised to rectangles or alternatively even arbitrary shapes. That may include counties, zip codes, electoral districts or administrative districts. Standard area-level mapping methods such as chloropleth maps suffer from very different area sizes or odd area shapes which can greatly distort the visual impression. The Kernelheaping package provides a way to convert these area-level data to a smooth point estimate. For the German capital city of Berlin, for example, there exists an open data initiative, where data on e.g. demographics is publicly available.  We first load a dataset on the Berlin population, which can be downloaded from: <a rel="nofollow" href="https://www.statistik-berlin-brandenburg.de/opendata/EWR201512E_Matrix.csv" target="_blank">https://www.statistik-berlin-brandenburg.de/opendata/EWR201512E_Matrix.csv</a> ```r library(dplyr) library(fields) library(ggplot2) library(Kernelheaping) library(maptools) library(RColorBrewer) library(rgdal) gpclibPermit() ``` ```r data &lt;- read.csv2(&quot;EWR201512E_Matrix.csv&quot;) ``` This dataset contains the numbers of inhabitants in certain age groups for each administrative districts. Afterwards, we load a shapefile with these administrative districts, available from:  <a rel="nofollow" href="https://www.statistik-berlin-brandenburg.de/opendata/RBS_OD_LOR_2015_12.zip" target="_blank">https://www.statistik-berlin-brandenburg.de/opendata/RBS_OD_LOR_2015_12.zip</a> ```r berlin &lt;- readOGR(&quot;RBS_OD_LOR_2015_12/RBS_OD_LOR_2015_12.shp&quot;) ``` ```r berlin &lt;- spTransform(berlin, CRS(&quot;+proj=longlat +datum=WGS84&quot;)) ``` Now, we form our input data set, which contains the area/polygon centers and the variable of interest, whose density should be calculated. In this case we' like to calculate the spatial density of people between 65 and 80 years of age: ```r dataIn &lt;- lapply(berlin@polygons, function(x) x@labpt) %>%    do.call(rbind, .) %>%    cbind(data$E_E65U80) ``` In the next step we calculate the bivariate kernel density with the “dshapebivr” function (this may take some minutes) using the prepared data and the shape file: ```r est &lt;- dshapebivr(data = dataIn, burnin = 5, samples = 10,                    adaptive = FALSE, shapefile = berlin,                    gridsize = 325, boundary = TRUE) ``` To plot the map with &quot;ggplot2&quot;, we need to perform some additional data preparation steps: ```r berlin@data$id &lt;- as.character(berlin@data$PLR) berlin@data$E_E65U80 &lt;- data$E_E65U80 berlinPoints &lt;- fortify(berlin, region = &quot;id&quot;) berlinDf &lt;- left_join(berlinPoints, berlin@data, by = &quot;id&quot;) kData &lt;- data.frame(expand.grid(long = est$Mestimates$eval.points[[1]],                                 lat = est$Mestimates$eval.points[[2]]),                     Density = as.vector(est$Mestimates$estimate)) %>%    filter(Density > 0) ``` Now, we are able to plot the density together with the administrative districts: ```r ggplot(kData) +   geom_raster(aes(long, lat, fill = Density)) +    ggtitle(&quot;Bivariate density of Inhabitants between 65 and 80 years&quot;) +   scale_fill_gradientn(colours = c(&quot;#FFFFFF&quot;, &quot;#5c87c2&quot;, &quot;#19224e&quot;)) +   geom_path(color = &quot;#000000&quot;, data = berlinDf, aes(long, lat, group = group)) +   coord_quickmap() ```</pre>
<div>
<figure> <img width="320" alt="" src="https://i1.wp.com/www.inwt-statistics.com/assets/images/b/Kernelheaping_2_BivariateDensity-3890c83b.png?resize=320%2C244&ssl=1" data-recalc-dims="1" height="244"> </figure>
</div>
<pre>This map gives a much better overall impression of the distribution of older people than a simple choropleth map:  ```r ggplot(berlinDf) +   geom_polygon(aes(x = long, y = lat, fill = E_E65U80, group = id)) +    ggtitle(&quot;Number of Inhabitants between 65 and 80 years by district&quot;) +   scale_fill_gradientn(colours = c(&quot;#FFFFFF&quot;, &quot;#5c87c2&quot;, &quot;#19224e&quot;), &quot;n&quot;) +   geom_path(color = &quot;#000000&quot;, data = berlinDf, aes(long, lat, group = group)) +   coord_quickmap() ```</pre>
<div>
<figure> <img width="320" alt="" src="https://i2.wp.com/www.inwt-statistics.com/assets/images/d/Kernelheaping_2_choropleth-8aaca67d.png?resize=320%2C244&ssl=1" data-recalc-dims="1" height="244"> </figure>
</div>
<pre>Often, as the case with Berlin we may have large uninhabited areas such as woods or lakes. Furthermore, we may like to compute the proportion of older people compared to the overall population in a spatial setting. The third part of this series shows how you can compute boundary corrected and smooth proportion estimates with the Kernelheaping package.</pre>
<div>
<p>Further parts of the article series <strong>Introducing the Kernelheaping Package</strong>:</p>
<ul>
<li><a rel="nofollow" href="https://www.inwt-statistics.com/read-blog/introducing-the-kernelheaping-package-512.html" target="_blank">Part 1: Univariate Kernel Density Estimation for Heaped Data</a></li>
<li><strong>Part 2: Bivariate Kernel Density Estimation for Heaped Data</strong></li>
</ul></div>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="https://www.inwt-statistics.com/read-blog/introducing-the-kernelheaping-package-ii.html"> INWT-Blog-RBloggers</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img width="1" alt="" src="http://feeds.feedburner.com/~r/RBloggers/~4/zIMkmhkYSkE" height="1"></a></div>