Использование Rcartogram на объекте SpatialPolygonsDataFrame

Question

Я пытаюсь сделать то же самое, заданное в этом вопросе, Cartogram + choropleth map in R, но начиная с SpatialPolygonsDataFrame и надеясь получить тот же тип объекта.

Я мог бы сохранить объект в виде шейп-файла, использовать scapetoad, снова открыть его и преобразовать обратно, но я предпочел бы все это внутри R, чтобы процедура была полностью воспроизводимой, и чтобы я мог кодировать десятки вариантов автоматически ,

Я раскодировал код Rcartogram на github и добавил свои усилия до сих пор here.

По существу, эта демонстрация создает SpatialGrid по карте, просматривает плотность населения в каждой точке сетки и преобразует ее в матрицу плотности в формате, который необходим для работы cartogram(). Все идет нормально.

Но как интерполировать исходные точки карты на основе вывода cartogram()?

Здесь есть две проблемы. Первый заключается в том, чтобы получить карту и сетку в одни и те же единицы, чтобы обеспечить интерполяцию. Второй - доступ к каждой точке каждого многоугольника, интерполирование и сохранение всех в правильном порядке.

Сетка находится в единицах сетки, а карта находится в проецируемых единицах (в случае примера длинной латы). Либо сетка должна быть проецирована в длинную, либо карту в единицы сетки. Моя мысль состоит в том, чтобы сделать поддельную CRS и использовать ее вместе с функцией spTransform() в package(rgdal), так как она обрабатывает каждую точку объекта с минимальной суматохой.

Доступ к каждой точке затруднен, потому что они несколько слоев вниз в объект SpPDF: объект> полигоны> Полигоны> линии> координаты, которые я думаю. Любые идеи, как получить доступ к ним, сохраняя при этом структуру общей карты нетронутой?

Answer 1

Эта проблема может быть решена с помощью пакета getcartr, который можно загрузить по адресу Chris Brunsdon's GitHub, так же красиво излагается в сообщении this.

Функция quick.carto делает именно то, что вы хотите - принимает входной сигнал SpatialPolygonsDataFrame и имеет выход SpatialPolygonsDataFrame.

воспроизводящих суть, например, в блоге здесь, в случае, если ссылка идет мертвым, со своим собственным стилем смешивается в & опечаток исправлено:

(Shapefile; World Bank population data)

library(getcartr) 
library(maptools) 
library(data.table) 

world <- readShapePoly("TM_WORLD_BORDERS-0.3.shp") 
#I use data.table, see blog post if you want a base approach; 
# data.table wonks may be struck by the following step as seeming odd; 
# see here: http://stackoverflow.com/questions/32380338 
# and here: https://github.com/Rdatatable/data.table/issues/1310 
# for some background on what's going on. 
world@data <- setDT(world@data) 

world.pop <- fread("sp.pop.totl_Indicator_en_csv_v2.csv", 
        select = c("Country Code", "2013"), 
        col.names = c("ISO3", "pop")) 

world@data[world.pop, Population := as.numeric(i.pop), on = "ISO3"] 

#calling quick.carto has internal calls to the 
# necessary functions from Rcartogram 
world.carto <- quick.carto(world, world$Population, blur = 0) 

#plotting with a color scale 
x <- world@data[!is.na(Population), log10(Population)] 
ramp <- colorRampPalette(c("navy", "deepskyblue"))(21L) 
xseq <- seq(from = min(x), to = max(x), length.out = 21L) 
#annoying to deal with NAs... 
cols <- ramp[sapply(x, function(y) 
    if (length(z <- which.min(abs(xseq - y)))) z else NA)] 

plot(world.carto, col = cols, 
    main = paste0("Cartogram of the World's", 
        " Population by Country (2013)"))