00.イントロ
＜３年生へお願い＞
●ハザードマップ学習．
前回も出しましたが，引き続きお願いします．伊勢市で配られているものや，どうゆうものがあるか，あるは皆さんにどんなアイディアがあるかなど特にテーマを絞る必要はありません．なるべく早急に調べてください．忙しいとは思いますがよろしくお願いします．
＜皆さんへお願い＞
●１２月７日（金）の防災訓練について
１２月７日（金）に大紀町にて防災訓練があります．授業などがある方は仕方ありませんが参加できる方は参加していただくようお願いします．
今日は，前回作成したラスターデータに属性情報を与え，それをRで表示してみます．
なお，「mlterm」は画面がちらつくため使用しないほうがいいかもしれません．
01.ファイルの転送
まず，前回作成した作業用のディレクトリに移動します．
cd R_test/Isebay
GRASSを立ち上げます．（Location：Isebay，Mapset：ユーザー名）
まず，GRASSのラスターファイルのリストをとります．
g.list rast
※今後，特に指定しませんがファイルを作成した後や削除した後など必要に応じてリストをとってください．
前回作成した「g5…」のファイルは今後必要ないため削除します．
g.remove rast=g5035,g5036,g5135,g5136,g5235,g5236
モニターを立ち上げ「geology」の属性情報などを見てみましょう．
d.mon x0
d.rast geology
d.what.rast
クリックしてみてください．すると，番号（カテゴリー番号）は入っていますが，その他の情報は入っておらず番号が何を示しているのかわかりません．
そこで，詳しい情報をラスターマップに与えます．
まず，「windows」でLegend（属性情報）について書かれたファイルを確認します．
「マイコンピュータ」→「Nドライブ」→「Geology」→「Legend」→「LEGED1J.HTM」（日本語） or 「LEGEND1E.HTM」（英語）
この年代については，各自勉強してください．（
参考）
このファイルは「htm」形式です．このままではGRASSに入れられません．
そこで，同じフォルダの「category」を見てください．（テキストファイルなので「メモ帳」か何かで見てください．）これは，GRASSで扱える形式になっています．
「category」を「eosmai」の今作業しているディレクトリに移動しましょう．
GRASSを実行しているプロンプトで自分のいる位置を確認してください．
pwd
リストをとります．
ls -al
※こちらも今後ページに表記がない場合でも適宜実行してください．
ファイルを転送させるには「ftp:File Transfer Protocol」という機能を使います．
先ほどは，「windows」から「category」を探しにいきましたが，今度は「Linux」の中から探しにいきます．
※「eosmer」にユーザーがない方は以下の作業はできません！！
lftp -u ユーザ名 192.168.2.1
「category」があるディレクトリに移動します．
cd home/earth/Geology
ls -al
cd Legend
「category」を転送させます．
get category
※「9670bytes transferred」と出ればOKです．
「lftp」から抜けます．
exit
転送できなかったっ人（ユーザーがなかった人）は，転送できた人からコピーしましょう．
cp /home/takehiro/R_test/Isebay/category .
※ちゃんとタブを使いましょう．
リストを取ります．
ls
「category」を見てみましょう．
less category
03.「awk」を使う
まず，「category」のファイルの行数を調べましょう．
wc category
すると173行ということがわかります．
ところで，「category」の５列目以降の説明は長すぎるため今回は省きたいと思います．手作業で５列目を消してもいいので173行は多すぎます．
そこで，「awk」という機能を使い，必要な列（１～４列目まで）を抜き出します．（
awkについて）
gawk '{print $1,$2,$3,$4}' category
すると，抜き出した結果が表示されると思います．
今回は，この抜き出した内容を使いたいのでファイルとして保存します．このLinuxの機能を「リダイレクション」といいます．
>：リダイレクション（意味：結果を「>」より右のファイルに書き込みますよー．）
gawk '{print $1,$2,$3,$4}' category > cat
「less」で見てみましょう．
less cat
04.ラスターマップに属性情報を与える
ラスターマップの属性情報を変更するコマンドに「r.reclass」がありました．今回はカテゴリーなどを直接入力するのではなく，先ほど作成した「cat」というファイルを用います．ここでも「リダイレクション」機能を使います．
<：（意味：「<」より右のファイルの中身を参考に左のコマンドを実行しますよー．）
r.reclass in=geology out=geology_cat < cat
g.list rast
地図を表示してみましょう．
d.rast geology_cat
ちょっと色が変ですが表示できましたか？
では，どんな情報が入っているか見てみましょう．
d.what.rast
ちゃんと説明も入っていましたか？
先ほどの変な色を昔の「geology」と同じカラーパレットにして，きれいな図に戻します．
まず，「geology_cat」に入っている情報を見ます．
r.info geology_cat
「600」までありますが，使っているのはこのうちの一部です．
カラーパレットを変更します．
r.colors geology_cat rast=geology
表示させて見ましょう．
d.rast geology_cat
05.統計情報の表示
では，どのような地質の地域が多いのか調べてみましょう．
この地図に関する統計情報を表示させて見ます．
r.report help
Paramtersのunitは対象とするカテゴリが全体に対してどれだけ占めているかを表示する単位を指定するものです．
k：k㎡を表示
p：%を表示
では，「geology_cat」の統計情報を見てみましょう．
.report geology_cat u=k,p | less
ここで，「|」は「パイプ」いうLinuxの機能で，「r.report」の結果を「less」に受け渡すという働きがあります．つまり，「r.report」の結果を「less」を使ってみることが出来るのです．これにより，画面をスクロールして結果を見ることが出来ます．
では，表示している結果の中で一番多いデータは何ですか？
undefin:37が一番多いと思います．しかし，これは海（海はデータが入っていません）です．陸地で一番占有面積が大きいのはk2-k2:8.58です．
06.MASKを作る
今回の解析では海は対象としません．
そこで，「MASK」を作成します．昔やっているのでコマンドだけ表示します．
r.reclass in=dem out=sea_land
&-999 thru -3 = 0 Sea
&-2 thru 2187 = 1 Land
end
※陸地を「-2」と指定しているのは地盤沈下している地域があるためです．
ちなみに間違えた人はやり直すとき以下のように入力してください．
r.reclass in=dem out=sea_land --o
すでにあるファイルをもう一回作るときは--o：overwrite
地図が出来ているか確認してみましょう．
g.list rast
d.rast sea_land
値を確認してみてください．
d.what.rast
MASKを作ります．
r.mask help
r.mask in=sea_land m=1
「MASK」というファイルがあるか確認してください．
g.list rast
もう一度「sea_land」を表示してみてください．
d.rast sea_land
海が真っ白けになりましたか？（NULL:データが入っていないことです）
もう一度，統計情報を表示させて見ましょう．
r.report geology_cat u=k,p -Nh
※-N：NULLのところは無視．
※-h：ヘッダーの削除．（hをつけないものと比較してみてください．）
※オプション（フラグ）はどこに書いても良いです．
「05.」で表示したものと割合などが違う値を示していると思います．
本来は，この結果を使って解析したり，論文にこの結果を載せたりします．そこで，ファイルとして保存する方法を学んでおきましょう．
r.report geology_cat u=k,p -Nh > geo_rep
「ls」や「less」でファイルを確かめてみてください．
07.「R」を使う
では，前回と同様にGRASSからRを使って見ましょう．
R
前回の内容が保存されているか見てみましょう．
ls()
すると「"G" "ise"」と表示されたと思います．
前回の続きなので「ise」を選択します．
ise
上矢印で下のコマンドを探してください．
library(spgrass6)
G <- gmeta6()
ちゃんとRが覚えていてくれているんですね．
下のコマンドも同様に探し「geology」を「geology_cat」に直してください．
ise <- rast.get6(c("geology_cat","dem"), cat = c(TRUE,FALSE), ignore.stderr = TRUE)
「ise」という関数の中に「dem」と「geology_cat」が入りました．
図を表示させて見ましょう．
image(ise,attr=1)
※なぜか表示できません！！
現在GRASSからRに入っています．つまり，ベースにGRASSが走っているわけです．そのため，RのプロンプトのままGRASSを使うことができます．
system("r.stats -cl geology_cat")
GRASSのコマンドをRの中から使う場合は，GRASSのコマンドを「system("○○○")」（○○○がGRASSのコマンド）と入力してあげればいいということです．
地層と標高の関係を見ましょう．表で見てもいいのですが，グラフの方がわかりやすく表示することが出来ます．
縦軸：標高（demの値：-1～2187）
横軸：地質（カテゴリー名）
→どの地層がどれくらいの標高に位置しているか．
今，iseという関数の中にgeology_catとDEMが入っています．
グラフを書いてみましょう．
ちなみに，「boxplot」とは「箱ひげ図」のことです．
boxplot(ise$dem ~ise$geology_cat)
この図は，それぞれの地層がどの標高に分布しているかを示した図です．
長方形の中心が正規分布の中央値を示し，長方形の上下端が中央値から25%づつのところを示しています．（
参考）
色をつけてみましょう．
boxplot(ise$dem ~ise$geology_cat,col="green")
横向きの図を描いてみましょう．
boxplot(ise$dem ~ise$geology_cat,col="green",horizontal=T)
もちろん，オプションはどこにつけても良いです．
今日はここまでです．
今日やったコマンドを保存して終了します．
q()→y
exit
今日は終わりです．
お疲れ様でした．
最終更新：2008年11月20日 19:26
