連載
» 2010年10月04日 00時00分 公開

実践! Rで学ぶ統計解析の基礎(6):夏の異常気象をオープン・データで確認 (2/2)

[柏野雄太,バクフー株式会社]
前のページへ 1|2       

夏の異常気象:NASA GISSによるグローバル表面温度データ

 さて、ここで視点をグローバルに変えて地球全体を考えてみましょう。今年の夏は地球全体をで見た場合どこでどのくらい「アノーマリー」だったのか、ということを、世界地図に空間プロットして確認してみましょう。

 Rを用いてデータを地図上に空間プロットするには、Geographic Infomation System(GIS)関係の機能を実装した地図表示パッケージや空間分析パッケージが必要になります。空間プロットについては、地図情報が絡むことで普通のプロットより複雑になるせいか、まだデ・ファクト・スタンダードといえるパッケージは出ていないようです。そこで、今回はいろいろなパッケージを組み合わせて空間プロットを実現します。具体的には、spパッケージを中心に、世界地図のポリゴンはmaptoolsパッケージにあるものを利用し、実際の空間プロットにはfieldsにあるas.image.plotを利用したいと思います。

 地球の地表面温度と海表面温度のアノーマリーを定期的に解析しているのが、先ほどもでてきたコロンビア大学にあるNASAのGoddard Institute for Space Studies(GISS)のSurface Temperature Analysisチームです。ちなみに、このNASAのGISSは、地球温暖化理論の父といわれるJames Hansenが長く所長を勤めていることからも分かるように、地球温暖化研究の中心的役割を果たしている機関の1つです。もともとは金星や水星の惑星探査データの解析や惑星大気のモデリングというNASAの王道的な業務をやっていたのですが、1970年代後半に入り惑星探査計画が下火になってリストラされる運命だったところを、地球温暖化研究という「ビジネスモデル」の転換をしてサバイバルしたそうです。

 地球の地表面温度と海表面温度のアノーマリーデータはGISSのSurface Temperature Analysisの以下のサイトから手に入れることができます。

Global Maps from GHCN Data

 今回は“Mean Period”を6、7、8月の3カ月データで取り、“Smoothing Radius”を250kmにしたデータを利用します。データはここに上げておきました。

 念のために記述しますが、このデータがどのようにして加工されたか、という詳細な過程を知りたい方は、データ加工用のソースコード(Fortran、C、Pythonおよびシェルスクリプト)と説明がオープンになっていますので、それを参照ください。

 さて、必要なパッケージをインストールし、データを読み込むところまで話を進めましょう。パッケージのインストールに指定するCRANサイトは筑波大学を選ぶとよいでしょう。

> install.packages("sp")
> install.packages("maptools")
> install.packages("fields")
> gtemp <- read.csv("http://spreadsheets.google.com/pub?key=0AlBuJgqcP5f3dHVCYm5SZUY4cDN1VFFMeWNKQ3IxbFE&hl=en&output=csv", header = T, skip = 1, sep = ',')
> names(gtemp) <- c("i", "j", "lon", "lat", "anom")
> head(gtemp)
  i j  lon lat anom
1 1 1 -179 -89 0.29
2 2 1 -177 -89 0.29
3 3 1 -175 -89 0.29
4 4 1 -173 -89 0.29
5 5 1 -171 -89 0.29
6 6 1 -169 -89 0.29
> gtemp$anom[grep(9999.00, gtemp$anom, fixed=T)] <- NA
> summary(gtemp)
       i                j             lon              lat           anom          
 Min.   :  1.00   Min.   : 1.0   Min.   :-179.0   Min.   :-89   Min.   :  -5.7000  
 1st Qu.: 45.75   1st Qu.:23.0   1st Qu.: -89.5   1st Qu.:-45   1st Qu.:   0.0200  
 Median : 90.50   Median :45.5   Median :   0.0   Median :  0   Median :   0.4600  
 Mean   : 90.50   Mean   :45.5   Mean   :   0.0   Mean   :  0   Mean   :   0.5536  
 3rd Qu.:135.25   3rd Qu.:68.0   3rd Qu.:  89.5   3rd Qu.: 45   3rd Qu.:   0.9800  
 Max.   :180.00   Max.   :90.0   Max.   : 179.0   Max.   : 89   Max.   :   6.0400  
                                                                NA's   :3729.0000  

 上記のコードにより、データはgtempという名前のデータフレームに格納されました。i、jはただのインデックスで、lonlatがそれぞれ経度、緯度という地表面の座標を表します。anomが今年の夏の地表面温度/海表面温度のアノーマリーです。上記のコードの最後から2番目でやっていることには少し注意が必要です。GISSデータに含まれている“9999.00”というデータは、データが入手できないという意味ですが、このままだとRでは正しい数値として認識されてしまいます。そこで、これをRにおける“Not Avairable”の意であるNAに置き換えました。

 ここまでデータを用意できたら、それを地図データに展開するべくspパッケージで空間プロット用のデータフレームに加工していきます。具体的には以下の手順です。

> library(sp)
> coordinates(gtemp) = c("lon", "lat")
> summary(gtemp)
Object of class SpatialPointsDataFrame
Coordinates:
     min max
lon -179 179
lat  -89  89
Is projected: NA 
proj4string : [NA]
Number of points: 16200
Data attributes:
       i                j             anom          
 Min.   :  1.00   Min.   : 1.0   Min.   :  -5.7000  
 1st Qu.: 45.75   1st Qu.:23.0   1st Qu.:   0.0200  
 Median : 90.50   Median :45.5   Median :   0.4600  
 Mean   : 90.50   Mean   :45.5   Mean   :   0.5536  
 3rd Qu.:135.25   3rd Qu.:68.0   3rd Qu.:   0.9800  
 Max.   :180.00   Max.   :90.0   Max.   :   6.0400  
                                 NA's   :3729.0000  
> gridded(gtemp) <- TRUE
> summary(gtemp)
Object of class SpatialPixelsDataFrame
Coordinates:
     min max
lon -180 180
lat  -90  90
Is projected: NA 
proj4string : [NA]
Number of points: 16200
Grid attributes:
    cellcentre.offset cellsize cells.dim
lon              -179        2       180
lat               -89        2        90
Data attributes:
       i                j             anom          
 Min.   :  1.00   Min.   : 1.0   Min.   :  -5.7000  
 1st Qu.: 45.75   1st Qu.:23.0   1st Qu.:   0.0200  
 Median : 90.50   Median :45.5   Median :   0.4600  
 Mean   : 90.50   Mean   :45.5   Mean   :   0.5536  
 3rd Qu.:135.25   3rd Qu.:68.0   3rd Qu.:   0.9800  
 Max.   :180.00   Max.   :90.0   Max.   :   6.0400  
                                 NA's   :3729.0000  

 上のコードですが、まず通常のデータフレームにcoodinatesで座標を指定することでSpatialPointsDataFrameに変換しました。次にgriddedを利用して全体がグリッドデータであることを明示してSpatialPixelsDataFrameにしました。

 このSpatialPixelsDataFrameを実際のプロットイメージにするのが、as.image.SpatialGridDataFrameで、これを利用すると最後に経度、緯度、そして表面温度アノーマリーの値のリストになります。このリストがプロットイメージです。

> anomimage <- as.image.SpatialGridDataFrame(gtemp["anom"])
> summary(anomimage)
  Length Class  Mode   
x   180  -none- numeric
y    90  -none- numeric
z 16200  -none- numeric
> class(anomimage)
[1] "list"

 さて、いよいよこのプロットイメージを空間プロットします。空間プロットにはspパッケージのimage関数を利用し、世界地図はmaptoolswrld_simlを利用しましょう。

> image(anomimage)
> library(maptools)
> data(wrld_simpl)
> plot(wrld_simpl, add = T)
プロットイメージを空間プロット

 これを見ると、確かにプロットはされていますが、同色のカラースキームが利用されているので、値の違いがこれではよく分かりませんし、カラーバーが表示されていないので値の大きさが把握できません。これを改善するために、fieldsパッケージのimg.plotを利用します。

> library(fields)
> image.plot(anomimage, nlevel=24)
> plot(wrld_simpl, add = T)
> grid(col = "grey", lty = 1)
カラースキームの調整とカラーバーを追加

 image.plotnlevelは値を色に割り振る際の分解能を指定します。

 上図なら、地球全体でこの夏に何が起きていたかがかなり分かります。この夏は緯度50度線の付近で温度の過度な上昇というアノーマリーが発生していました。ロシアやウクライナ、カザフスタン、モンゴルでは例年の平均気温より5度から6度も高いという異常さを表しています。また日本などの極東やヨーロッパ全域、そしてアメリカ東部も平均気温が2度以上も高くなっています。

 結果の詳細はGISSのSurface Temperature Analysisチームの以下の報告にありますので、興味がある方は参照ください。

 この中緯度付近で起きた世界的な熱波は、直接的には地球規模の温暖化とは関係がなく、偏西風のブロッキング現象により引き起こされたと説明する研究者もいます。

 そして地球規模の温暖化が、そのブロッキング現象をより引き起こしやすくするという地球シミュレーターを利用した数値実験の結果もあるようです。

 この夏の異常な温度上昇という現象がどのようなメカニズムで引き起こされたかということは、今後の気候学、気象学の研究を待つことになりますが、このように直近のオープン・データを利用し、少し手間をかけるだけで、気候や気象の状況をモニタリングするということさえ個人でできるようになりました。本当に良い時代になったものだと思います。

次回について

 次回もまた、オープンソースのRを利用して、オープンデータにアクセスし、オープンアイディアの活用したいと思います。それでは、また2週間後に会いましょう。

Index

夏の異常気象をオープン・データで確認

Page1
今回の前口上
夏の異常気象:東京の気温データ

Page2
夏の異常気象:NASA GISSによるグローバル表面温度データ
次回について


前のページへ 1|2       

Copyright © ITmedia, Inc. All Rights Reserved.

RSSについて

アイティメディアIDについて

メールマガジン登録

@ITのメールマガジンは、 もちろん、すべて無料です。ぜひメールマガジンをご購読ください。