各ポリゴン内の観測のコリオプレットをプロットする方法は？私は、次のポリゴンを用いchoroplethマップをプロットしたい

Question

：

library(sp) 
library(sf) 
library(spatialEco) 
library(tigris) 

br_tracts <- tracts(state = 'LA', county = 'East Baton Rouge', cb = T, year = 2016) 
plot(br_tracts) 

そして、私は、各多角形内に入る観測点の数に色をマップしたいです。以下は私のサンプルデータです。 dtSpatialを作成するには、チグリスパッケージをアップロードする必要があります。

dtSpatial <- new("SpatialPointsDataFrame" 
    , data = structure(list(parenttype = c("garbage", "garbage", "garbage", 
"garbage", "garbage", "garbage", "garbage", "garbage", "garbage", 
"garbage", "garbage", "garbage", "garbage", "garbage", "garbage", 
"garbage", "garbage", "garbage", "garbage", "garbage")), .Names = "parenttype", row.names = c(NA, 
-20L), class = c("tbl_df", "tbl", "data.frame")) 
    , coords.nrs = numeric(0) 
    , coords = structure(c(-91.043777, -91.026382, 0, -91.027748, 0, -91.08049, 
-91.047173, -91.172501, -91.162384, -91.139465, -91.087585, -91.152748, 
-91.163086, -91.185814, -91.135101, 0, -91.105972, 0, -91.168846, 
-91.041435, 30.452148, 30.447191, 0, 30.412008, 0, 30.415289, 
30.420155, 30.430065, 30.478041, 30.460482, 30.429127, 30.469275, 
30.436682, 30.420218, 30.453447, 0, 30.431898, 0, 30.466148, 
30.416723), .Dim = c(20L, 2L), .Dimnames = list(NULL, c("longitude", 
"latitude"))) 
    , bbox = structure(c(-91.185814, 0, 0, 30.478041), .Dim = c(2L, 2L), .Dimnames = list(
    c("longitude", "latitude"), c("min", "max"))) 
    , proj4string = new("CRS" 
    , projargs = NA_character_ 
) 
) 

与えられた観測がポリゴンで表示され、私が試したかどうかを確認するには、次の

pts.poly <- point.in.polygon(point.x = dtSpatial$longitude, 
          point.y = dtSpatial$latitude, 
          pol.x = fortify(br_tracts)$long, 
          pol.y = fortify(br_tracts)$lat) 

> pts.poly 
[1] 1 1 0 0 0 0 0 1 1 0 1 0 0 0 0 0 1 0 0 0 

をしかし、どのように私は、各ポリゴン内の観測値の数に色をマッピングしていますか？

Answer 1

ここでは1つの方法を示します。あなたのデータを見て、最初にしたいのは、dtSpatialに投影を追加することです。次に、各ポリゴンに存在するデータポイントの数を数えたいとします。これは、GISToolsパッケージのpoly.counts()を使用して実現できます。出力はベクトルです。この関数はポリゴンIDを使用してデータポイントをカウントしているため、IDとして名前が表示されます。 stack()を使用してデータフレームに変換します。 br_tractsをggplotのデータフレームに変換します。最後に、マップを描きます。私はgeom_cartogram()を使用しました。私はこの関数を2回使用しています。私は最初にポリゴンを描画しています。そして、私はcountを使って色を塗りつぶしています。私はこの仕事をするmap_idと一致しています。注：geom_map()、geom_polygon()を使用できます。 coord_map()およびtheme_map()はオプションです。

library(tigris) 
library(GISTools) 
library(RColorBrewer) 
library(ggalt) 
library(ggthemes) 

# Add projection to dtSpatial, which is identical to projection of br_tracts 
proj4string(dtSpatial) <- CRS(proj4string(br_tracts)) 

# How many data poiints stay in each polygon. Polygon ID is the name of the vector. 

count <- poly.counts(pts = dtSpatial, polys = br_tracts) 
count <- stack(count) 

mymap <- fortify(br_tracts) 

ggplot() + 
geom_cartogram(data = mymap, aes(x = long, y = lat, map_id = id), 
       map = mymap) + 
geom_cartogram(data = count, aes(fill = values, map_id = ind), 
       map = mymap, color = "black", size = 0.3) + 
scale_fill_gradientn(colours = rev(brewer.pal(10, "Spectral"))) + 
coord_map() + 
theme_map()