1

我有一个 shapefile“网格”，它基本上是美国和加拿大的网格图。每个网格都有自己的特定 ID (grid_code)。我还有一个 csv 文件，其中包含与各个网格相关的鸟类物种的相对丰度值。由于每个物种都有其独特的分布，因此物种表仅包含所有可用网格代码的子集。

我正在尝试创建一个主题图，显示每个物种的相对丰度。为此，我首先使用 leftjoin 函数将 csv 文件加入到网格 shapefile 的属性表中，其中 grid_code 是公共列。在为 tmap 指定调色板后，我正在为每个物种的相对丰度创建一个 png 文件。

当我手动输入 i 的值时，以下代码可以正常工作；但是，当我运行整个代码时，它会在 leftjoin 函数之后停止，并为每个物种创建一个空白 png 文件。我已经尝试了几件事，但没有解决这个问题。这只是代码放置的问题，还是我的 for 循环有什么问题？

谢谢。

# reading in csv file with relative abundance info
rel_ab <- read.csv("ra_test.csv", h=T, stringsAsFactors = F)

# specifying color palette for tmap
colfunc <- colorRampPalette(c("lightpink", "red"))

# running for loop code where I'm reading and reprojecting the grid
# shapefile each time, then joining it to the species table, and finally 
# plotting the tmap

# only want to work with one species at a time, so first:
ra <- unique(rel_ab$species)

for (i in ra){
   species <- rel_ab[which(rel_ab$species == i),]
   grid <- readOGR(dsn = "grid_clip", "grid") # reading grid shapefile
   grid83<- spTransform(grid, CRS("+init=epsg:5070")) # transorming shapefile
   grid83@data <- left_join(grid83@data, species) 

# grid83 now contains a column that includes all the grid codes PLUS a 
# column for relative abundance values (V2) for species i. Again, because a
# species is not found in every part of US/CANADA, there are a lot of V2
# rows that are NA.

# Here is what grid83@data looks like:

head(grid83@data)

AREA PERIMETER TEMP_ TEMP_ID GRID_CODE BCR STATE STATE_CODE STRATA_ KMSQ  X V2 V3 species
1 463738112  87471.12  6548   92624      6981   5   ALK          3     S94  461 NA NA NA    <NA>
2 463749632  87485.10  6551   92625      6982   5   ALK          3     S94  461 NA NA NA    <NA>
3 463766720  87499.32  6554   92626      6983   4    YT         93     S68  461 NA NA NA    <NA>

   png(width = 1000, height = 1000, filename = paste("RA_maps/",i,".png", sep=""))

   tm_shape(grid83)+
    tm_fill("V2", title = "BBS Summer Distribution", textNA = "None counted", 
  style = "fixed", breaks = c(Inf, 100, 30, 10, 3, 1, 0.05), colorNA = "white", palette = colfunc(6))+
    tm_shape(uscan83)+ 
    tm_borders(lwd = 1.5)+
   tm_shape(gray83)+
    tm_fill(col="grey85")+
    tm_borders()+
   tm_layout(
     inner.margins = c(.02, .02, .02, .02),
     legend.title.size = 1.8,
     legend.text.size = 1,
     legend.bg.color = "white",
     legend.bg.alpha = 0,
     legend.width = 0.22)

   dev.off()}
4

1 回答 1

2

print如果你在 for 循环中运行它，你必须明确地绘制 tmap 图。或者，save_tmap包中有一个可能有用的。请参阅以下代码块，其中还包含一些其他 tmap 函数。由于我没有形状和数据，我无法自己测试。

for (i in ra){
   species <- rel_ab[which(rel_ab$species == i),]
   grid83 <- read_shape("grid_clip/grid.shp", projection = 5070)
   grid83 <- append_data(grid83, species) 

   m <- tm_shape(grid83)+
    tm_fill("V2", title = "BBS Summer Distribution", textNA = "None counted", 
  style = "fixed", breaks = c(Inf, 100, 30, 10, 3, 1, 0.05), colorNA = "white", palette = colfunc(6))+
    tm_shape(uscan83)+ 
    tm_borders(lwd = 1.5)+
   tm_shape(gray83)+
    tm_fill(col="grey85")+
    tm_borders()+
   tm_layout(
     inner.margins = c(.02, .02, .02, .02),
     legend.title.size = 1.8,
     legend.text.size = 1,
     legend.bg.color = "white",
     legend.bg.alpha = 0,
     legend.width = 0.22)

   save_tmap(m, width = 1000, height = 1000, units="px", filename = paste("RA_maps/",i,".png", sep=""))
}
于 2016-02-08T13:14:15.510 回答