2016-09-01 14 views
4

私は、1900年以来フロリダと交差したハリケーン/熱帯低気圧を示すマップの小さな倍数グリッドをプロットしようとしています。私はいくつかの空間クエリーを使って大西洋嵐のデータベースをサブセット化しましたこのプロジェクトのために。forループを使って小さなループをプロットするR

フロリダ州のポリゴン、いくつかの連続した州、フロリダのいくつかの大都市、もちろんオケエホビ湖の上に、私の限られた数のハリケーントラックの行shapefileをプロットしています。

library(maptools) 
library(gdata) 
library(RColorBrewer) 

setwd("~/hurricanes") 

# read shapefiles 
florida <- readShapePoly("florida.shp") 
south <- readShapePoly("south.shp") 
hurricanes <- readShapeLines("hurricanes-florida.shp") 

cities <- readShapePoints("cities.shp") 
lakes <- readShapePoly("lakes.shp") 

# miami, orlando and tallahassee (in FL) 
cities <- subset(cities, ST == "FL") 

# don't need ALL the 'canes 
hurricanes1900 <- subset(hurricanes, Season >= 1900) 

mycolors <- brewer.pal(5, "YlOrRd") 

pdf(file = "hurricanemaps.pdf", ,width=8,height=20) 
    par(mfrow=c(15,5), mar=c(1,1,1,1)) 
    for(i in 1:nrow(hurricanes1900)) 
     { 
     plot(south, col="#e6e6e6", border = "#999999") 
     plot(florida, col="#999999", border = "#999999", add = TRUE) 
     plot(lakes, col="#ffffff", border = "#999999", add = TRUE) 
     plot(cities, pch=10, cex=.1,col="#000000", bg="#e38d2c", lwd=1, add = TRUE) 
     plot(hurricanes1900[i,], col = mycolors[cut(hurricanes$MAX_Wind_W, breaks = 5)], 
      lwd=3, add = TRUE); title(hurricanes1900$Title[i]) 
    } 
dev.off() 

三大問題私は遭遇しています:

1)ループが私に各嵐のマップを与えている、ここで簡単なコードです。私は、毎年(たとえ嵐のない年であっても)、その年のすべての嵐、好ましくはラベル付きのフロリダ/南地図をグリッドに作成することを望んでいます。

2)ループの各特定の行にあるものだけでなく、すべての嵐の中で風速の色を設定したいと思います。そうすれば、強い嵐(1992年のアンドリューのような)は、唯一の年の嵐であっても、より暗く表示されます。おそらく、私はそれに応じてスタイルを付けることができるカテゴリ(H1、H2など)変数を記録してこれを解決することができます。

3)私がNo.1を見つけ出すことができると仮定して、私は各嵐の経路上でレンダリングするラベルを得るのに問題があります。 maptoolsのドキュメントは素晴らしいものではありません。とにかく

、ここでの出力は(タイトルはシェープファイル内の2つのフィールドを連結したものです)これまでのところです:

enter image description here

本当の問題は、あなたの助けを事前に1号のおかげです。

+0

小倍数は、ファセットと呼ばれる 'ggplot'、というパッケージにも実装されています。こちらをご覧ください:https://www.r-bloggers.com/world-map-panel-plots-with-ggplot2-2-0-ggalt/ – Chris

答えて

0

再現性のあるデータがないため、このデモンストレーションのデータを収集しました。次回からSOユーザーのための最小限の再現可能なデータを提供してください。それはあなたがより多くの助けを受けるのを助けるでしょう。

あなたが達成したいことは、ggplot2で行うことができます。毎年マップを作成する場合は、facet_wrap()を使用できます。風に基づいて色を追加する場合は、パスを描画するときにaes()で行うことができます。ハリケーンの名前を追加したい場合は、ggrepelパッケージを使用して、注釈を簡単に提供できます。滑らかなパスを描きたい場合は、さらにデータ処理が必要です。

library(stringi) 
library(tibble) 
library(raster) 
library(ggplot2) 
library(ggthemes) 
library(ggrepel) 
library(RColorBrewer) 
library(data.table) 

# Get some data. Credit to hmbrmstr for a few lines in the following code. 

mylist <- c("http://weather.unisys.com/hurricane/atlantic/2007H/BARRY/track.dat", 
      "http://weather.unisys.com/hurricane/atlantic/2007H/TEN/track.dat", 
      "http://weather.unisys.com/hurricane/atlantic/2006H/ERNESTO/track.dat", 
      "http://weather.unisys.com/hurricane/atlantic/2006H/ALBERTO/track.dat") 

temp <- rbindlist(
      lapply(mylist, function(x){ 
       foo <- readLines(x) 
       foo <- read.table(textConnection(gsub("TROPICAL ", "TROPICAL_", 
            foo[3:length(foo)])), 
            header=TRUE, stringsAsFactors=FALSE) 

       year <- stri_extract_first(str = x, regex = "[0-9]+") 
       name <- stri_extract_first(str = x, regex = "[A-Z]{2,}") 

       cbind(foo, year, name) 

      } 
     )) 


### Add a fake row for 2017 
temp <- temp %>% 
     add_row(ADV = NA, LAT = NA, LON = NA, TIME = NA, WIND = NA, 
       PR = NA, STAT = NA, year = 2017, name = NA) 

### Prepare a map 
usa <- getData('GADM', country = "usa", level = 1) 
mymap <- subset(usa, NAME_1 %in% c("Florida", "Arkansas", "Louisiana", 
            "Alabama", "Georgia", "Tennessee", 
            "Mississippi", 
            "North Carolina", "South Carolina")) 
mymap <- fortify(mymap) 


# Create a data.table for labeling hurricanes later. 
label <- temp[, .SD[1], by = name][complete.cases(name)] 

g <- ggplot() + 
    geom_map(data = mymap, map = mymap, 
       aes(x = long, y = lat, group = group, map_id = id), 
        color = "black", size = 0.2, fill = "white") + 
    geom_path(data = temp, aes(x = LON, y = LAT, group = name, color = WIND), size = 1) + 
    scale_color_gradientn(colours = rev(brewer.pal(5, "Spectral")), name = "Wind (mph)") + 
    facet_wrap(~ year) + 
    coord_map() + 
    theme_map() + 
    geom_text_repel(data = label, 
        aes(x = LON, y = LAT, label = name), 
        size = 2, 
        force = 1, 
        max.iter = 2e3, 
        nudge_x = 1, 
        nudge_y = -1) + 
    theme(legend.position = "right") 

enter image description here

関連する問題