Найти максимальную точку с в каждом многоугольнике для набора многоугольников R

Question

Я уверен, что на этот вопрос ответили в другом месте, но я не смог придумать его с помощью поиска.

У меня есть точки, представляющие города в стране, а также население для каждого города. У меня также есть полигональный файл округов. Я хочу найти местоположение самого большого города в каждом графстве.

Как это может быть сделано?

Вот некоторые данные

структура (список (Страна = c ("нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас") ",
"нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас", "нас" "," us "), City = c(" cabarrus", "cox store", "cal-vel", "briarwood townhouses", "barker heights", "davie "
Перекресток "," Деревня Краб-Пойнт "," Азалия "," Честерфилд "," Чарльзмонт "," Коннор "," Клевер Гарден "," Корриер Хайтс "," Каллисонс "," Крест-Вью-акрс "," Клегг "," Ханаан park", "chantilly", "belgrade", "brices crossroads", "bluff", "butner", "bottom", "bandy", "bostian heights"), AccentCity = c("Cabarrus", "Cox Store"), "Cal-Vel", "Briarwood Townhouses", "Barker Heights", "Davie Crossroads", "Crab Point Village", "Azalea", "Chesterfield", "Charlesmont", "Connor", "Clover Garden", "Corriher Heights", "Callisons", "Crestview Acres", "Clegg", "Canaan Park", "Chantilly", "Belgrade", "Brices Crossroads", "Bluff", "Butner", "Bottom", "Bandy", "Bostian Heights"), Region = c("NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC", "NC"), население = c(NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, A_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_, NA_integer_), широта = (35,2369444, 35,275, 351, 351), 36,429), 36,429 (29), 35, 351, 8, 35, 8, 35, 8, 35, 8, 8, 35, 8, 8, 35, 8, 8, 35, 8, 8, 8, 35, 8, 8, 35, 8, 8, 35, 8, 8, 29, 8, 8, 35, 8, 29, 28, 8, 35, 28, 8, 35, 28, 8, 12, 12, 12, 12, 35, 5, 35, 35, 8, 12, 35, 8, 12, 13, 35, 8, 12, 13, 35, 8, 12, 35, 28, 35, 28, 8, 12, 12, 12, 25, 35, 37, 351, 351, ПБ 35.7419444, 36.1883333, 35.5605556, 35.0841667, 35.0213889, 35.8655556, 36.2761111, 36.3016667, 34.88, 34.8186111, 35.8377778, 36.1319444, 36.4747222, 35.6419442, 35.75.588, 284, 282, 684, 2882, 2882, 2882, 2882, 284, 2882, 2882, 2882, 2882, 2882, 2882, 2882, 2882, 2882, 2882, 541, США, 35.744444, 381, 288, 2818, 2818, США, 35.7419444. -80,535, -76,7305556, -82,4713889, -81,6611111, -81,5127778, -78,1486111, -79,4630556, -80,635, -76,7255556, -80,5427778, -78,8497222, -79,7852778, -76,17,178,87,87,88,88, 87,88,88,88 и 681,87,87,87,88,1755,3681,85,85,1641,35,1641,3641,05,1755,35,1645,05,15,160% составляют -77678771, 87771, 8677, -7, 858, -76, -7, 8700, -76, -77, -77, -77, -77, -77, -77, -77,., -80.7741667, -81.09, -80.9294444)),.Names = c ("Страна", "Город", "AccentCity", "Регион", "Население", "Широта", "Долгота"), row.names = с (544L, 889L, 551L, 434L, 190L, 975L, 894L, 147L, 717L, 700L, 831L, 773L, 862L, 559L, 915L, 753L, 584L, 695L, 262L, 437L, 372L, 537L, 406L, 178L, 02L), класс = "data.frame")

и некоторый код для чтения в Северной Каролине

xx <- readShapePoly(system.file("shapes/sids.shp", package="maptools")[1],
                IDvar="FIPSNO", proj4string=CRS("+proj=longlat +ellps=clrk66"))

plot(xx)

Я хочу найти город с максимальным населением в каждом округе. извините, у меня нет воспроизводимого примера. Если бы я это сделал, у меня был бы ответ!

Answer 1

Короткий ответ: вы должны использовать gContains(...) в упаковке rgeos,

Вот длинный ответ.

В приведенном ниже коде мы берем шейп-файл с высоким разрешением для округов Северной Каролины из базы данных GADM и геокодированный набор данных для городов Северной Каролины из базы данных геологической службы США. Последний уже имеет информацию округа, но мы игнорируем это. Затем мы сопоставляем города с соответствующими округами, используя gContains(...)добавьте эту информацию во фрейм данных городов и определите самый большой город в каждом округе с помощью пакета data.table. Большая часть работы состоит из 4 строк кода ближе к концу.

library(raster)   # for getData(...);   you may not need this
library(foreign)  # for read.dbf(...);  you may not need this
library(rgeos)    # for gContains(...); loads package sp as well

setwd("< directory for downloaded data >")
# get North Carolina Counties shapefile from GADM database
USA         <- getData("GADM",country="USA",level=2)   # level 2 is counties
NC.counties <- USA[USA$NAME_1=="North Carolina",]      # North Carolina Counties
# get North Carolina Cities data from USGS database
url <- "http://dds.cr.usgs.gov/pub/data/nationalatlas/citiesx010g_shp_nt00962.tar.gz"
download.file(url,"cities.tar.gz")
untar("cities.tar.gz")
data      <- read.dbf("citiesx010g.dbf",as.is=TRUE)
NC.data   <- data[data$STATE=="NC",c("NAME","COUNTY","LATITUDE","LONGITUDE","POP_2010")]
## --- evverything up to here is just to set up the example

# convert cities data.frame to SpatialPointsDataFrame
NC.cities <- SpatialPointsDataFrame(NC.data[,c("LONGITUDE","LATITUDE")],
                                    data=NC.data,
                                    proj4string=CRS(proj4string(NC.counties)))
# map cities to counties
city.cnty   <- gContains(NC.counties,NC.cities,byid=TRUE)
# add county information to cities data
NC.data$county <- apply(city.cnty,1,function(cnty)ifelse(any(cnty),NC.counties@data[cnty,]$NAME_2,NA))
# identify largest city in each county
library(data.table)
result <- setDT(NC.data)[,.SD[which.max(POP_2010)],by="county"]
head(result)
#      county             NAME   COUNTY LATITUDE LONGITUDE POP_2010
# 1:  Jackson        Cullowhee  Jackson 35.31371 -83.17653     6228
# 2:   Graham     Robbinsville   Graham 35.32287 -83.80740      620
# 3:   Wilkes North Wilkesboro   Wilkes 36.15847 -81.14758     4245
# 4:    Rowan        Salisbury    Rowan 35.67097 -80.47423    33662
# 5: Buncombe        Asheville Buncombe 35.60095 -82.55402    83393
# 6:    Wayne        Goldsboro    Wayne 35.38488 -77.99277    36437

Рабочей лошадкой здесь является строка:

city.cnty   <- gContains(NC.counties,NC.cities,byid=TRUE)

Это сравнивает каждую точку в SpatialPointsDataFrame NC.Cities для каждого полигона в SpatialPolygonsDataFrame NC.counties и возвращает логическую матрицу, в которой строки представляют города, а столбцы представляют округа, а [i,j] элемент TRUE если город i в округе j, FALSE иначе. Мы обрабатываем матрицу построчно в следующем выражении:

NC.data$county <- apply(city.cnty,1,function(cnty)ifelse(any(cnty),NC.counties@data[cnty,]$NAME_2,NA))

используя каждую строку подряд для индексации таблицы атрибутов для NC.counties чтобы извлечь название округа.

Данные, которые вы указали в своем вопросе, имеют некоторые проблемы, которые тем не менее поучительны. Во-первых, шейп-файл NC в maptools Пакет имеет относительно низкое разрешение. В частности, это означает, что некоторые из прибрежных островов полностью отсутствуют, поэтому любой город на одном из этих островов не будет привязан к графству. У вас могут быть те же проблемы с вашими реальными данными, так что следите за ними.

Во-вторых, сравнивая COUNTY столбец в исходном наборе данных USGS с county В колонке, которую мы добавили, есть 3 (из 865) округов, которые не согласны. Оказывается, в этих случаях база данных USGS была неправильной (или устаревшей). У вас может быть та же проблема, так что следите за этим тоже.

В-третьих, еще три города не отображались ни в одном графстве. Все они были прибрежными городами и, вероятно, отражали небольшие неточности в шейп-файле Северной Каролины. У тебя ночью тоже есть эта проблема.