前回の[気象庁の観測データを使う (1)] [気象観測データのグラフ作成と基本的な統計処理 (Rによる) (1)]に つづいて、複数の観測点の気象データをつかって現象の空間分布をみることにする。
この教材をつくりはじめた2020年12月17日に、北関東の大雪のニュースがあった。 「藤原」という地名を聞いた。 (東武線の新藤原のあたり(これは栃木県)かと思ったが、そうではなかった。) 群馬県みなかみ町の藤原だ。ここには気象庁のアメダスの観測所がある。 わたしはこの事例を専門的に追いかけてはいないが、 報道や気象庁ウェブサイトで知ったことを紹介しておきたい。
気象庁ウェブサイトの「防災情報」の中の「気象情報」のページには、 気象庁やその地方気象台などが発表した、注意すべき気象状況の情報についての文章が収録されている。 全国のものと、府県別のものがある。
群馬県のところを見ると、前橋地方気象台が出した「大雪に関する群馬県気象情報」がある。 2020年12月20日夜現在最新のものは、12月20日10時59分発表の「第6号」である。 この「号」の番号はときどき1からつけなおされている。 12月16日17時00分発表の(前回の)「第6号」に、つぎのような部分がある。 (箇条書きの部分の数値は、かっこの外が12時間降雪量、かっこの中が積雪深である。 なお、わたしはこの「降雪量」を「降雪深」とよぶことにしている。 別記事[雪について観測されている量]も見てほしい。)
［雪の実況］ １６日１６時現在の１２時間降雪量（積雪の深さ）は みなかみ町藤原 ３２センチ（１６３センチ） みなかみ町幸知 ２２センチ（ ７６センチ） となっています（アメダス速報値）。また、みなかみ町藤原では１２月１５日の積雪差日合計が１１７センチとなり、１９８９年の統計開始以来の極値を更新しました。（これまでの極値は１９９１年１月６日の８７センチ）
なお、この「みなかみ町幸知」はアメダスの地点「みなかみ」と思われる。 ただし、前橋地方気象台のウェブサイト (https://www.jma-net.go.jp/maebashi/) のお知らせ記事によれば、 アメダス「みなかみ」は2020年11月26日に 幸知から湯原に移転しているので、たぶん正しくは「みなかみ町湯原」での観測値だろう。 (なお、下の教材中でふれる気象庁のアメダス地点表のファイル(12月19日にダウンロードしたもの)では、みなかみ観測点の所在地は湯原となっている。)
また、2020年12月18日の毎日新聞の報道 (有料記事だがhttps://mainichi.jp/articles/20201218/ddl/k10/040/079000c?pid=14509)にはつぎのようにあった。
前橋地方気象台によると、７２時間降雪量は、みなかみ町藤原で同日午後４時までに２１９センチに達し、統計開始以来最多を更新した。 ... 同町みなかみでも７２時間降雪量が１５９センチとなり観測史上１位を記録。同町藤原の積雪量は１７日午後５時現在で２メートル超となった。
なお、今回のこの地点の積雪深は新記録ではない。 藤原の積雪深の極値は 2006年1月26日の 301 cm、 みなかみは 2006年1月28日の 275 cm である (気象庁の「過去の気象データ検索」サイトで検索した)。
そこで、時空間分布を見る対象を、この大雪の事例をふくめるように、 2020年12月の、群馬県北部を中心とする約 150 km四方の領域にふくまれる地点の、 気象庁のアメダスの観測データとした。
観測値の空間分布を図示するためには、観測地点の位置を知る必要がある。 地球上の位置は、ふつう、緯度・経度 (および 高さ) で表現される。
前回の教材でとりあげた気象庁ウェブサイトの「過去の気象データ ダウンロード」では、 地点を選択する際に、その地点の緯度・経度がポップアップの形で表示される。 しかし、ダウンロードされたCSVファイルには緯度・経度の情報ははいっていない。 ポップアップを書きとるのはうまくない方法だ。 気象庁の観測地点の緯度・経度をふくむ一覧表をつかうべきだ。
気象庁ウェブサイトの 「知識・解説」の「気象」のうち「気象に関する観測と予報の技術の解説」のうちに「気象の観測」というリンクがある。 そのページの表題は「気象観測・気象衛星」となっている。 そのページの「地上・地域気象観測」のところに、 「地上気象観測所について」と「地域気象観測システム(アメダス)について」がある。 気象台などの(いわゆる)気象官署での観測は「地上気象観測所について」のほうに書かれているのだが、 ほとんどの気象官署はアメダスの観測点にもなっているので、 ここではアメダスのほうにしぼってみることにする。
「地域気象観測システム(アメダス)について」 http://www.jma.go.jp/jma/kishou/know/amedas/kaisetsu.html のページには、アメダスではどんな観測をしているかの解説とともに、 「観測所一覧」がある。 このうち、CSV形式の「地域気象観測所一覧 [ZIP圧縮形式]」のファイル http://www.jma.go.jp/jma/kishou/know/amedas/ame_master.zip をつかう。 2020年12月19日にダウンロードしたものを、つぎのところにコピーしておく。 ファイル形式は zip で、展開すると ame_master.csv というCSVファイル (日本語文字コードはShift_JIS)がはいっている。
なお、PDFファイルの内容は、このファイルと、もうひとつの「地域気象観測所一覧（雪） 」との内容を そのままあわせたもののようだ。
このファイルに収録された地点は、現在のアメダスの観測点である。 過去に観測していたが現在は観測していない地点のことや、観測所の移動経歴などを知りたいときは、 同じページからリンクされている「アメダス地点情報履歴ファイル」を見る必要がある。
アメダスの地点には5けたのID番号がついている。この番号は国際地点番号ではない。
気象庁の ame_master.csv は、そのまま表計算ソフトウェアで表として見ることもできるし、 R で read.csv で読むこともできる。 しかし、わたしがやりたい作業にとってはちょっと不便な形をしているので、あらかじめ組みかえておくことにする。
組みかえたおもな点は、緯度・経度の単位である。 気象庁の地点表の緯度・経度は、それぞれ「度」と「分」の列にわかれている。 「度」のほうには整数値がはいり、「分」のほうには小数点下1けたまでの値がはいっている。 緯度・経度の表現は、60進法の度・分・秒 (秒よりこまかいところが必要ならば秒の10進小数)か、 度とその10進小数か、がふつうだと思うが、気象庁の地点表はどちらでもないのだ。 わたしは、「度」の単位であつかいたいので、「度」の列の値に「分」の列の値の60分の1をたして換算する。
【なお、日本国内の地点ならばすべて北緯・東経なので簡単なのだが、 世界の地点をあつかうときは南緯・西経のあつかいをきめておく必要がある。 南緯・西経を負にするならば、度の列だけでなく分の列の値も符号をかえる対象になる。】
組みかえのもうひとつは、地名のカタカナ表記である。 観測所名の表記は漢字が基本となっているが、それだけでは読みかたがわからない。 気象庁の地点表にはカタカナ表記がある。 ただし、1バイトのカタカナの文字コード (漢字の半分の幅で表示されることが多かったので「半角カナ」とよばれることが多い)が つかわれている。 その体系では、濁点・半濁点はそれがつくカナ文字とは独立した文字のようにあつかわれている。 いまでは1バイトカタカナには対応していないソフトウェアが多いので、 漢字をふくむ日本語文字コードのカタカナ(いわゆる「全角カナ」)に変換しておきたい。 ネット検索でみつけた horio さんのブログに書かれた R プログラム をつかった。 (わたしはむかし同様なプログラムを Awk 言語で書いたことがあるが、 そのときの作業かぎりのもので、汎用性はもたせなかったと思う。)
積雪深の観測地点は別の表になっているが、ame_master.csv の「備考1」の列に 積雪深の観測点としての地点IDがはいっている。 例外的に積雪深だけの観測点があるかもしれないが、積雪深を観測するアメダスの地点のほとんどはこれでわかる。 積雪深を観測していない地点は「－」(ダッシュ)がはいっているが、NA にかえておくことにする。
ame_master.csv には「所在地」の列がある。住所にあたる地名 (番地まではない) が書かれている。 それにくわえて、気象官署でもある地点のばあいは、官署名がはいっている。 (くぎり文字として「全角」の空白がはいっている。) 官署名だけ拾い出しておくことにする。
ame_master.csv には IDも地名も同じ地点が複数行ふくまれていることがある。 観測機器がちょっと離れたところにわかれて置かれているのだ。 緯度・経度・高さがちがうといってもわずかなので、さきに出てきたほうだけを残すことにする。
以上の方針で組みかえをした。
観測点の空間分布を見るには、地図上に表現したほうがよい。 しかし、わたしは (地図は別のソフトウェアでかいているので） R で地図をかく応用パッケージのつかいかたを指導することが (まだ) できない。
ここでは、横軸に経度、縦軸に緯度をとった x-y 平面上の位置として、観測地点の位置を表示することにする。 これも一種の地図投影方法とみなすことができ、そのばあいは「正距円筒図法」とよばれる。 【地図投影方法については、2000年ごろ地理情報学・地図学の教材としてつくった [地図投影]のページを見てほしい。】実際の地球では、一定の経度間隔が、高緯度ほど短い距離間隔に対応する。 正距円筒図法の地図上の距離はそうなっていないことに注意する必要がある。 緯度と経度を同じ尺度で画面に表示すると、赤道に近いところ以外では、現実よりも横長になる。 つぎのプログラム例では、東西の経度幅 2度と、南北の緯度幅 1.5 度を画面で同じ長さにした。 これは、現実よりもやや縦長になっている。
アメダスの地点のうち 北緯 36～37.5度、東経 138.5～140.5度の範囲にあるものの位置と 地点名を、amedas_stn.csv のデータによって、R のグラフィックス機能で表示した。
観測している気象要素を、記号の形と色の両方でしめした。
地理的背景として、海岸線、地形、県境などをいっしょにしめすことがのぞましいのだが、 残念ながら準備ができていない。 地図帳、または、国土地理院の「地理院地図」( https://maps.gsi.go.jp/ ) などを併用して見てほしい。
ここでつかった地点表ファイルは日本全国をふくんでいるので、 これとちがった地域についての図にかえることは、 プログラムの緯度経度範囲 (xlimit, ylimit の値として与えている) を変えればすぐできる。 (ただし、図の下の凡例をしめす位置も、あわせて変える必要がある。) しかし、このあとしめす気象データは、上の図の緯度経度範囲だけダウンロードしたものなので、 他の地域についてやるには、ダウンロードからやりなおす必要がある。
空間的にも時間的にも、広く見たい、こまかく見たいという欲求もあり、データ処理能力がかぎられているという事情もある。 ここでは、空間については、上記の図と同じ 北緯 36～37.5度、東経 138.5～140.5度の範囲にふくまれるアメダス観測点(気象官署もふくむ)にかぎることにした。 時間については、2020年12月1日から19日までの日値をとりあげることにした。 (いまも大雪の可能性のある状態がつづいているので、リアルタイムにデータを更新できるとよいのだが、 わたしの作業手順では、急にデータをさしかえるのはむずかしい。)
前回と同じく、気象庁の「過去の気象データ ダウンロード」のサイトから CSV 形式のデータを取得した。 一度にダウンロードできるデータの量の制約はあるが、時間方向の件数がすくなければ観測地点数をふやすことができる。 地点選択では都道府県別の地図の画面で地点を指定する必要があるが、複数の都道府県の地点を指定していっしょにダウンロードすることができる。 この方法で、降水量、積雪深、降雪深 [注]、日照時間、日平均気温の、 それぞれ対象とした緯度経度範囲にふくまれる地点を指定してデータをダウンロードした。 地点数は降水量がいちばんおおく、積雪深と降雪深の地点、日照時間と気温の地点は、それぞれ、降水量の地点の部分集合にあたる。
ダウンロードしたファイル(ファイル名を data.csv からかえた)を教材ウェブサイトに置く。気象庁サイトが混雑しないように、実習にはこれをつかってみてほしい。 (もしあたらしい時期のデータをくわえたければ、ダウンロードしなおすのもよいが、 そのときは、下にのべる対応表づくりもやりなおす必要が生じるだろう。)
[2020-12-22 補足] CSVファイルの読みとばすべき行数は、日最大積雪深のファイルでは 6行、その他の4つのファイルでは 5行である。
作業を始めてみて、このダウンロードデータの形式に由来するむずかしさにぶつかった。 観測地点と、その地点のデータがふくまれるCSVファイルの列との対応を、機械的につけることがむずかしいのだ。
対策として、ダウンロードしたCSVファイルを見て、手作業で、地点名と、CSVファイルでその地点がはいっている列番号との対応表を つくっておくことにした。 ただし、「表を見る」作業は表計算ソフトウェアでやるのが便利で、そこでは列は A～Z, AA, AB, ... とアルファベットで表示されている。 R での処理では列の番号が整数でほしいのだが、表計算の列記号から番号への変換は機械にやらせることにする。 その変換のついでに、地点表 amedas_stn.csv を読んで、地点の緯度・経度も列番号とならべて表にしておくことにする。
ここまで準備すれば、ダウンロードしたデータファイルと、key 表の2つを読んで、 日ごとの観測値の空間分布を図化することができる。
図には、観測点の位置に記号(まる)をかき、その上側に数字で観測値をしめすことにした。 ただし、観測値が欠けているばあいは「NA」と表示することにした。 また、観測が不完全なばあい(「品質情報」が「8」でないばあい) は、観測値を表示する文字の色を赤にした。
図形出力先は、画面と、画像ファイル(png)とのどちらでもできるようにしたが、 ここにリンクしたプログラムは、画像ファイルに出すものにしてある。 画面に出したいときの変更点は、「png(...)」と「dev.off()」の行をコメントにする(行の頭に「#」をいれる)ことと、 コメントにしてある「par(ask=TRUE)」をいかす(行の頭にいれてある「#?#」などをけす)ことである。
プログラムでは19日ぶんの図をつづけてつくるが、ここではそのうちある1日ぶんの図をのせておく。 [2020-12-22 日最大積雪深の図をさしかえました。さしかえ前は、入力データの読みとばすべき行数がちがっていたので、1日ずれたものでした。] (降水量が「NA」となっている地点は藤原で、雪があまりに多くて雨量計が正常に動作しなかったのだ。)
データのあるすべての観測点についてそれぞれ、時間(ここでは日単位なので日付)を横軸、 観測値を縦軸にとったグラフをつくるプログラムを作成した。
降雪深、最大積雪深、日照時間、降水量は、数値のゼロに絶対的な意味があるので、 ゼロから観測値までの棒と、観測値の位置での 〇(まる)の両方で表示することにした。 ただし、観測が不完全な(「品質情報」が 8 でない)ばあいは、まるをかかないことにした。
気温は、℃単位のゼロは約束上の原点にすぎず絶対的な意味がないので、 棒はつかわず、二重のまるで表示することにした。 ただし、観測が不完全なばあいは、大きなまるをかかないことにした。
正方形の1画面の上下に同じ横軸のグラフを2つ配置することにした。
日降雪深と日最大積雪深、日照時間と日平均気温は、同じ地点の両要素を同じ画面に配置する。 日降水量は、画面表示のばあいは、2地点ずつ表示する。 ただし、画像ファイルにしまうばあいは、1地点ごとに地点名をつけたファイルにしたいので、 画面の上半分だけに図がある形になる。図形出力先は、画面と、画像ファイル(png)とのどちらでもできるようにしたが、 ここにリンクしたプログラムは、画像ファイルに出すものにしてある。 画面に出したいときの変更点は、「png(...)」と「dev.off()」の行をコメントにする(行の頭に「#」をいれる)ことと、 コメントにしてある「par(ask=TRUE)」をいかす(行の頭にいれてある「#?#」などをけす)ことである。 降水量のグラフについては、1画面に2つの図を入れる指定「par(mfrow=c(2,1))」が2か所にあるが、 画面に出したいときは「par(ask=TRUE)」につづくほうを、画像ファイルに出したいときは「png(...)」につづくほうを 生かす (他方をコメントにする)。
次のプログラム例では、画像ファイル名には、地点名が日本語文字(漢字・ひらがな)ではいるようになっている。 実行する計算機の OS によっては、このままでは画像ファイルがうまくできないかもしれない。
[2020-12-27] 画像ファイル名に、地点名がASCIIコードのローマ字ではいるように修正したプログラム例を用意した。 こちらならばどの OS でも動作すると思う (わたしが確認したのは Windows 10 日本語版だけだが)。 【[2021-01-14] ただし、地名の頭文字だけ大文字になっていることに注意。 HTMLで参照されるファイル名は大文字・小文字の区別があり、まちがえると参照されない。 OSによっては大文字・小文字だけがちがうファイルが区別されないこともある。】
プログラムではデータのある全地点の図をつづけてつくるが、ここではそのうち1地点の図をのせておく。 藤原の降水量にはデータ欠損や不完全な観測が見られる。 【[2021-01-14] (画像ファイル名の大文字・小文字をまちがえたせいで) 日降雪深と日最大積雪深の図が表示されていなかったのを修正しました。】
