2014-03-31 64 views
1

我一直在試圖用shapefile連接數據框並繪製結果。我試圖按照@ jlhoward對this question的回答中提出的方法使用。我有national dataset of vaccination rates by post code。我試圖將它與澳大利亞統計局郵政編碼ESRI shapefile合併,並根據郵政編碼按照其他問題繪製結果。將數據框附加到shapefile並繪製它

這是我當前的嘗試坐:

library(rgdal) 
library(maptools) 
library(ggplot2) 
library(plyr) 
setwd("~/Google Drive/R/PC_Shapes") 
vac.data <- read.csv(file = "Postcode2013.csv", header=TRUE, sep=" ", na.string="NA", dec=".", strip.white=TRUE) 
postcode <- readOGR("POA06aAUST_region.shp", layer="POA06aAUST_region") 
[email protected]$id <- rownames([email protected]) 
postcode.df <- fortify(postcode) 
postcode.df <- join(postcode.df, [email protected], by="id") 
postcode.df <- merge(postcode.df, vac.data, all=TRUE) 
ggp <- ggplot(data=postcode.df, aes(x=long, y=lat, group=group)) 
ggp <- ggp + geom_polygon(aes(fill=LEVEL))   
ggp <- ggp + geom_path(color="grey", linestyle=2) 
ggp <- ggp + coord_equal() 
ggp <- ggp + scale_fill_gradient(low = "#ffffcc", high = "#ff4444", space = "Lab", na.value = "grey50", guide = "colourbar") 
ggp <- ggp + labs(title="Vaccination Rates: Australia") 
print(ggp) 

我想我的問題在於以下兩行中,我知道我需要分配by.x =和/或by.y =:但是我不斷收到我不清楚它們來自哪裏的錯誤。我不知道我想要在這裏實現...

postcode.df <- join(postcode.df, [email protected], by="id") 
postcode.df <- merge(postcode.df, vac.data, all=TRUE) 

我shape文件在這一點上結束了超過5,500,000的觀察和R開始奮鬥。

其值得注意的是在ABS shapefile中有一些我沒有數據的郵編。我不知道如何排除它們。他們可能是一個問題。在之前的嘗試中,我嘗試了這種方法:

library("sp","rgdal","plyr") 
setwd("~/Google Drive/R/PC_Shapes") 
ogrListLayers("POA06aAUST_region.shp") 
postcode <- readOGR("POA06aAUST_region.shp", layer="POA06aAUST_region") 
vacs <- read.csv("~/Google Drive/R/PC_Shapes/Postcode2013.csv") 
PNI <- melt(vacs, id=c("Postcode","Percent.not.fully.immunised")) 
postcode$POA_2006 %in% PNI$Postcode 
postcode$POA_2006[which(!postcode$POA_2006 %in% PNI$Postcode)] 
levels(postcode$POA_2006[which(!postcode$POA_2006 %in% PNI$Postcode)]) 

如果任何人有任何想法我摔倒了,我非常感謝任何提示。如果這是一個明顯的問題,我對R很感興趣。

+0

無法獲取「國家郵政編碼疫苗接種率數據集」。 – hrbrmstr

+0

@hrbrmstr道歉。它現在應該可以訪問。 – vengefulsealion

回答

2

很多東西在這裏錯了。 read.csv行... sep =「,」,而不是「」。 要確保你合併在正確的列上。使用head(df)來查看你想要合併的df的第一對夫婦行,或者使用str(df)來查看一堆關於它的信息。

祝你好運。

library(rgdal) 
library(maptools) 
library(ggplot2) 
library(plyr) 
gpclibPermit() 

vac.data <- read.csv(file = "Postcode2013.csv", header=TRUE, sep=",", na.string="NA", dec=".", strip.white=TRUE) 
postcode <- readOGR("POA06aAUST_region.shp", layer="POA06aAUST_region") 
# took too long to fortify on whole data set 
postcode <- postcode[1:50,] 
[email protected]$id <- rownames([email protected]) 
pts <- fortify(postcode,region="id") 
postcode.df <- merge(pts,postcode,by="id", stringsAsFactors=F) 
postcode.df$id <- as.numeric(postcode.df$id) 
postcode.df2 <- merge(postcode.df, vac.data, by.x="POA_2006", by.y="PC_2006") 
postcode.df2 <- postcode.df2[order(postcode.df2$id,postcode.df2$order),] 

ggplot() + geom_polygon(aes(x=long,y=lat, group=group, 
          fill=Percent.not.fully.immunised), 
         data=postcode.df2) 
+0

感謝Cory的反饋,非常感謝。 – vengefulsealion