ポリゴン空間結合のxy点を新しい単位に再コード化する方法は？

Question

私は、データフレームからxyポイントをマップし、それらを新しいポリゴンに関連するマーカーで再コード化して、そのユニットにリンクすることができます。私の理解から、これには空間結合が必要です。データフレームの変換プロセスに関するものの

R convert zipcode or lat/long to county

https://gis.stackexchange.com/questions/137621/join-spatial-point-data-to-polygons-in-r

が、私は何かが欠けていると私は、私はエラーを作っていますどこか分からない。私はこれらのリンクを追ってきました。私は、例としてオープンデータAPIの迅速で汚れた再現可能なコードを貼り付けています。私が使っている形状ファイルは、ここで見つけるコミュニティ地区NYCファイルです：

http://www1.nyc.gov/site/planning/data-maps/open-data/districts-download-metadata.page

コードは以下の通りです：

この場合、再現可能なサンプルが開いているデータから、これらの点を取り、それらをマッピングすることができますようニューヨーク市のコミュニティ地区ファイル。助けを歓迎します！

#libraries-------------------------- 

     library(ggplot2) 
     library(ggmap) 
     library(sp) 
     library(jsonlite) 
     library(RJSONIO) 

#call api data-------------------------- 

     df = fromJSON("https://data.cityofnewyork.us/resource/24t3-xqyv.json?$query= SELECT Lat, Long_") 
     df = df[1:10000] 
     df <- data.frame(t(matrix(unlist(df),  nrow=length(unlist(df[1]))))) 
     names(df)[names(df) == 'X2'] = 'x' 
     names(df)[names(df) == 'X1'] = 'y' 
     df = df[, c("x", "y")] 
     df$x = as.numeric(as.character(df$x)) 
     df$y = as.numeric(as.character(df$y)) 
     df$x = round(df$x, 4) 
       df$y = round(df$y, 4) 
       df$x[df$x < -74.2] = NA 
       df$y[df$y < 40.5] = NA 
       df = na.omit(df) 

    #map data---------------------------- 




    cd = readOGR("nycd.shp", layer = "nycd") 
      cd = spTransform(cd, CRS("+proj=longlat +datum=WGS84")) 
      cd_f = fortify(cd) 


#map data 
      nyc = ggplot() + 
      geom_polygon(aes(x=long, 
       y=lat, group=group), fill='grey', 
      size=.2,color='black', data=cd_f, alpha=1) 


      nyc + geom_point(aes(x = x, y = y), data = df, size = 1) 


#spatial join--------------------------- 


# Convert pointsDF to a SpatialPoints object 
        pointsSP <- SpatialPoints(df, 
         proj4string=CRS("+proj=longlat +datum=wgs84")) 

       Error in CRS("+proj=longlat +datum=wgs84") : 
        unknown elliptical parameter name 

        # Use 'over' to get _indices_ of the Polygons object containing each point 
        indices <- over(pointsSP, cd) 

Answer 1

正しいCRSを指定する必要があります。私はその後spatialreference.org（WGS 1984）

# Convert pointsDF to a SpatialPoints object 
pointsSP <- SpatialPoints(df, 
          proj4string=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")) 

でそれを見つけ、あなたはcdのCRSにポイントを再投影する必要があります。いつものよう

pointsSP <- spTransform(pointsSP, proj4string(cd)) 

# Use 'over' to get _indices_ of the Polygons object containing each point 
indices <- over(pointsSP, cd) 
head(indices) 
# BoroCD Shape_Leng Shape_Area 
# 1 407 139168.79 328355685 
# 2 407 139168.79 328355685 
# 3 407 139168.79 328355685 
# 4 407 139168.79 328355685 
# 5 407 139168.79 328355685 
# 6 105 35287.62 43796717