0

我想用内布拉斯加州的县级作物产量数据为 STAT 类做一些空间统计分析。为此，我需要每个县的地理质心的经度和纬度。有人知道如何在 R 中做到这一点吗？我知道它可以在 ArcGIS 中完成，但我现在无法访问它。

4

4 回答 4

3

您没有提供任何从哪里获得 shapefile 的详细信息，但我从这里获得了一个，您可以gCentroid从这里使用rgeos：

library(rgdal)
library(sp)
library(rgeos)

nebraska <- readOGR("CountyBoundsUTM/", "CountyUTM")

gCentroid(nebraska, byid=TRUE)

## SpatialPoints:
##           x       y
## 0  721768.5 4636738
## 1  430938.8 4524651
## 2  698036.4 4566570
## 3  370970.6 4641340
## ...
## 89 623301.6 4603228
## 90 618883.0 4486931
## 91 439295.3 4582756
## 92 493680.8 4522680
## Coordinate Reference System (CRS) arguments: +proj=utm +zone=14 +datum=NAD83
## +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0
于 2014-10-27T03:11:14.887 回答
0

您也可以使用 提取SpatialPolygons*对象coordinates的质心，但不会像SpatialPoints使用 一样返回质心rgeos::gCentroid。

例如：

library(rgdal)
download.file('http://dds.cr.usgs.gov/pub/data/nationalatlas/countyp020_nt00009.tar.gz', 
              f <- tempfile()) # ~ 4.5 Mb
untar(f, exdir=tempdir())
counties <- readOGR(tempdir(), 'countyp020')
xy <- coordinates(counties)
head(xy)

#        [,1]     [,2]
# 0 -153.3905 69.30193
# 1 -156.0582 71.33094
# 2 -155.6695 71.24763
# 3 -155.5164 71.23148
# 4 -155.1846 71.18189
# 5 -155.6126 71.00725

请注意，正如@Spacedman 在评论中指出的那样，多边形应首先投影到平面坐标系。

于 2014-10-27T07:24:04.587 回答
0

sf包的最新版本（我认为是从版本 1 开始）使用 Google 的S2 库进行球面几何计算。优点是质心计算不仅仅是平面的。相关方法是st_centroid()。具有重要空间范围的区域的示例：

library(rnaturalearth)
# ne_countries() returns 'sp'-type data by default
nc <- ne_countries(continent = "Asia", returnclass = "sf")
  
library(sf)
# long-lat data in WGS84
st_crs(nc)
# use st_geometry() to plot only the polygons and not the associated data
plot(st_geometry(nc), axes = T)
plot(st_centroid(st_geometry(nc)), pch = "+", col = "red", add = T)

# 'sf' integrates nicely with 'ggplot2':
library(ggplot2)
ggplot(nc) + geom_sf() +
    geom_sf(aes(geometry = st_centroid(st_geometry(nc))), colour = "red")
于 2022-02-18T10:01:16.450 回答
-2

您可以使用包中的get_map()函数将ggplot2包中的美国县地图数据提取maps到数据框中。然后，您可以按县（或您想用来定义地理中心的任何方法）计算纬度/经度列范围的中点。

于 2014-10-27T02:31:10.513 回答