目标:向 Map 添加密度
经过测试的包:ggplot,ggmap,sp,rgdal
根据rdgal包,我有一个特定州的 Map .
library(rgdal)
us<-getData('GADM', country='USA', level=2)
WV <-subset(us, NAME_1=="West Virginia")
plot(WV)
我想在 Map 上绘制一些处方数据 . 基于geocode()方法,在以下帖子中指定:R convert zipcode or lat/long to county,我将县名转换为一组纬度和经度坐标 .
library(ggmap)
WV.opioid <- opioid.dat[opioid.dat$State.Name == "West Virginia",]
WV.coordinates <- geocode(as.character(WV.opioid$County.Name))
WV.opioid <- cbind(WV.opioid, WV.coordinates)
WV.opioid <- WV.opioid[which(!is.na(WV.opioid$lat)),]
rownames(WV.opioid) <- NULL
我试图在 Map 上绘制一列数据集:
library(ggplot2)
#opioid data from following link: https://www.datafiles.samhsa.gov/study-series/national-survey-drug-use-and-health-nsduh-nid13517
opioid.map <- ggplot(WV.opioid, aes(x = lon, y = lat,
fill = cut_number(Diff.Opioid.Prescription.Rates,
5))) +
geom_polygon(color = "gray10", size = 0.2) +
coord_equal() +
viridis::scale_fill_viridis(discrete = TRUE) +
labs(title = "Test",
fill = NULL) +
theme_void() +
theme(legend.position = "bottom",
panel.background = element_rect(fill = NA, colour = "#cccccc"))
opioid.map
我不确定如何修复这个情节,并希望得到帮助 . 我还尝试使用rdgal包并将我的数据添加到SpatialPoints对象但尚未成功 .
1 回答
如果你没有纠正modifiable areal unit problem,那么它真的是正确的 .
如果阿片类药物数据是按邮政编码,则
tigris
包中有该分辨率的shapefile .我认为主要问题在于你在第二个代码块中进行的转换 . 一旦你不再这样做,绘制处方率的氯含量应该相当简单 .