第一次痛苦地尝试做一张空间地图

Question

我正在努力让我的第一张地图工作。我已经阅读了我能找到的所有文档，但我无法将它们全部放在一起以在地图上查看我的数据。

这是我到目前为止所做的。 1. 我创建了一个非常基本的数据表，其中包含 3 个观察值和 5 个变量作为一个非常简单的起点。

str(Datawithlatlongnotvector)

    'data.frame':   3 obs. of  5 variables:
     $ Client: Factor w/ 3 levels "Jan","Piet","Susan": 2 1 3
     $ Sales : int  100 1000 15000
     $ Lat   : num  26.2 33.9 23.9
     $ Lon   : num  28 18.4 29.4
     $ Area  : Factor w/ 3 levels "Gauteng","Limpopo",..: 1 3 2

（该地区是南非的省份，也是我下载的SHP文件，见下文）

1. 我下载了一张南非地图并将所有 3 个文件（.dbf、shp 和 shx）文件放在同一目录中 - 之前的错误，但我从另一个用户的问题中找到了答案。 http://www.mapmakerdata.co.uk.s3-website-eu-west-1.amazonaws.com/library/stacks/Africa/South%20Africa/index.htm 并选择简单底图。

2. 我创建的地图如下：

   SAMap <- readOGR(dsn = ".", layer = "SOU-level_1")

我可以用 plot(SAMap) 绘制显示省份的国家地图

4. 我也可以绘制数据

   plot(datawithlatlong)

5. 我看到了如何制作 SpatialPointsData 框架的说明，我做到了：

   coordinates(Datawithlatlong) = ~Lat + Lon

6. 我不知道如何将它们组合在一起并执行以下操作： 在地图上用不同的颜色显示数据（100,1000 和 15000），即 1 到 500 之间是一种颜色，501 到 10 000 之间是一种颜色，10 000 以上是一种颜色。

Answer 1

也许可以尝试ggplot2 的一些功能，例如：

   map = ggplot(df, aes(long, lat, fill = Sales_cat)) +  scale_fill_brewer(type = "seq", palette = "Oranges", name = "Sales") + geom_polygon() 

使用scale_fill_brewer，您可以在地图上用颜色表示比例。您应该创建一个factor 变量来表示根据销售范围的类别（“Sales_cat”）。在任何情况下，形状文件都必须转换为data.frame。

Answer 2

尝试将“SAMap”作为国家 shapefile，将“datawithlatlong”作为转换为 SpatialPointDataFrame 的数据：

library(maptools)
library(classInt)
library(RColorBrewer)

# Prepare colour pallete  
plotclr <- brewer.pal(3,"PuRd")
class<-classIntervals(datawithlatlong@data$sales, n=3, style="fixed", fixedBreaks=c(0, 500,1000,10000)) # you can adjust the intervals here 
colcode <- findColours(class, plotclr)

# Plot country map
plot(SAMap,xlim=c(16, 38.0), ylim=c(-46,-23))# plot your polygon shapefile with appropriate xlim and ylim (extent)

# Plot dataframe convereted to SPDF (in your step 5)
plot(datawithlatlong, col=colcode, add=T,pch=19)

# Creating the legend
legend(16.2, -42, legend=names(attr(colcode, "table")), fill=attr(colcode, "palette"), cex=0.6, bty="n") # adjust the x and y for fixing appropriate location for the legend

Answer 3

我生成了一个更大的数据集，因为我认为只有 3 个点很难看出事情是如何运作的。

library(rgdal)
library(tmap)
library(ggmap)
library(randomNames)

#I downloaded the shapefile with the administrative area polygons
map <- readOGR(dsn = ".", layer = "SOU")

#the coordinate system is not part of the loaded object hence I added this information
proj4string(map) <- CRS("+init=epsg:4326")

# Some sample data with random client names and random region 
ADM2   <- sample(map@data$ADM2, replace = TRUE, 50)
name   <- randomNames(50)
sales  <- sample(0:5000, 50)

clientData <- data.frame(id = 1:50, name, region = as.character(ADM2), sales,
                         stringsAsFactors = FALSE)

#In order to add the geoinformation for each client I used the awesome 
#function `ggmap::geocode` which takes a character string as input an 
#provides the lon and lat for the region, city ...

geoinfo <- geocode(clientData$region, messaging = FALSE)

# Use this information to build a Point layer
clientData_point <- SpatialPointsDataFrame(geoinfo, data = clientData)
proj4string(clientData_point) <- CRS("+init=epsg:4326")

现在我希望回答这个问题的部分：

# Adding all sales which occured in one region
# If there are 3 clients in one region, the sales of the three are 
# summed up and returned in a new layer
sales_map <- aggregate(x = clientData_point[ ,4], by = map, FUN = sum)

# Building a map using the `tmap` package`
tm_shape(sales_map) + tm_polygons(col = "sales")

编辑：

这是一个ggplot2 解决方案，因为您似乎想坚持下去。

首先，对于 ggplot，您必须将您的 SpatialPolygonDataFrame 转换为普通的 data.frame。幸运的是，broom::tidy() 会自动完成这项工作。 其次，您的Lat 值缺少-。我添加了它。 第三，我重命名了您的对象以减少输入。

point_layer<- structure(list(Client = structure(c(2L, 1L, 3L),
                            .Label = c("Jan", "Piet", "Susan"), 
                             class = "factor"), 
                             Sales = c(100, 1000, 15000 ), 
                             Lat = c(-26.2041, -33.9249, -23.8962), 
                             Lon = c(28.0473, 18.4241, 29.4486), 
                             Area = structure(c(1L, 3L, 2L), 
                            .Label = c("Gauteng", "Limpopo", "Western Cape"), 
                             class = "factor"), 
                             Sale_range = structure(c(1L, 2L, 4L), 
                            .Label = c("(1,500]", "(500,2e+03]", "(2e+03,5e+03]", "(5e+03,5e+04]"), 
                            class = "factor")), 
                            .Names = c("Client", "Sales", "Lat", "Lon", "Area", "Sale_range"), 
                             row.names = c(NA, -3L), class = "data.frame") 

point_layer$Sale_range <- cut(point_layer$Sales, c(1,500.0,2000.0,5000.0,50000.0 )) 

library(broom)
library(ggplot2)

ggplot_map <- tidy(map)

ggplot() + geom_polygon(ggplot_map, mapping = aes(x = long, y = lat, group = group),
                        fill = "grey65", color = "black") +
           geom_point(point_layer, mapping = aes(x = Lon, y = Lat, col = Sale_range)) +
  scale_colour_brewer(type = "seq", palette = "Oranges", direction = 1)