Memo
QuantumGIS
A. 地図の準備
A-1. 日本の基盤地図（海岸線や道路，市町村境界など）
1．国土地理院HP「基盤地図情報：http://www.gsi.go.jp/kiban/」からダウンロード
＊利用には登録必要
2．基盤地図（GML形式）をダウンロード
＊必要な県を選ぶ
＊必要な市町村を選ぶ
＊ZIP形式で保存される
3．ダウンロードしたデータをShapeファイルに変換．
＊ダウンロードしたデータはKML形式なので，Shape形式に変換する必要がある．
＊変換はコンバートソフト「基盤地図情報ビューア・コンバーター」（FDGV.exe）を使う．
＊「基盤地図情報ビューア・コンバーター」（FDGV.exe）は「http://fgd.gsi.go.jp/download/」ダウンロード．
4．基盤地図情報ビューア・コンバーター」（FDGV.exe）を起動
＊「ファイル」-「新規プロジェクト」
＊プロジェクトのタイトルをつける．例えば，「基盤地図情報沖縄」
＊ダウンロードしたデータを選択する．KML形式でもZIPファイルでも良い．
＊「コンバート」-「シェープファイルへ出力」
＊もし，標高データをダウンロードした際は，「コンバート」-「標高メッシュをシェープファイルへ出力」を選択する．
5．ShapeファイルをQGISで読み込む
B-1. 国外の基盤地図（海岸線や道路，市町村境界など）
1．地球地図HPから全球データや国ごとのデータをダウンロード
http://www.gsi.go.jp/kankyochiri/globalmap.html
＊利用には地球地図国際運営委員会(ISCGM)に登録必要
2．各国版のダウンロード
＊「National data」から必要な国を探す．
＊Oceaniaでデータがあるのは，「New Zealand」・「Papua New Guinea」・「Samoa」のみ
＊Asiaはそこそこ揃っている「Philippines」もある
3．ダウンロードしたShapeファイルをQGISに載せる
＊．地球地図（日本）版もダウンロードできる （こちらは登録必要なし）
http://www.gsi.go.jp/kankyochiri/gm_jpn.html
＊ベクタとして「行政界・水系・人口集中域・交通」がある．（第2版2011年公開））
＊ラスタとして「標高・土地被覆・土地利用・植生」がある．（第1.1 版2006年公開）
B. GPSで取った調査地点をプロットする
1．GPSデータ（ID，緯度，経度）をcsv形式で保存．
2．「レイヤ」－「デリミテッドテキストレイヤの追加」でcsvファイルを選択．
3．「OK」
4．プロットされる．
5．データをshapeファイルに保存する（保存しないと終了したら消えてしまう）
レイヤを選択して右クリック－「名前をつけて保存」
C. 空中写真を幾何補正して載せる
1．ラスタ画像の挿入
「プラグイン」－「Georeference」－「ファイルを開く」－「ラスタを開く」
2．座標をあわせる
ラスタ画像上で「ポイントの追加」，基盤地図上で「マップキャンパスより」を選び，
交差点などを探して座標合わせをする．
しかし，陸上の地図しかないので，海も含めて合わせる際は写真の6点の座標を調べておき，直接X座標（経度XXX.XXXXXX度）とY座標（XX.XXXXXX度）を入力する．
3．自動幾何補正
「ファイル」－「ジオリファレンシングの開始」－「変換タイプ」－「出力ラスタ」で任のファイル名を決める
4．作成したレイヤを載せる
「レイヤ」－「ラスタライヤの追加」
＊「いちとり式教材レシピ」参照
D. 複数の図をレイヤ―する
1．座標系・投影法が異なっているとうまく重ならない
2．「シェープファイルを選び右クリック」-「レイヤCRS」で現状を確認する
3． 基本はWGS84（EPSG4326）になっているはず
4．「WGS84」になっていない場合は，画面右下の地球儀みたいなボタン「CRSステータス」をクリック
5．「オンザフライCRSを有効にする」にチェックを入れて，「WGS84」を選択する
6．おそらく，きれいに重なるはず．
E. イラストレーターでトレースした海岸線を載せる
1．PlugX-Shapeをインストールする．（3万円くらいする）
2．トレースした海岸線を準備する（イラストレーターで作成）．
3．「ウィンドウ」-「PlugX-Shape 4.0」を選ぶ
4．「座標マーカーレイヤ―作成」を選ぶ．4点か3点か選べる．
5．「テキストオブジェクト」（Tの文字のやつ」を使って座標系を入力する．
例 30度40分50秒，135度45分50秒ならば，30/40/50&135/45/50
＊文字を書く際は「ドラッグ」して入力エリアをつくるのではなく，
「クリック」して書き出すことがポイント．
6．「シェープファイル書き出し」レイヤ―の選択，シェープタイプ（点とか線とか）を選んで次へ．座標の変換をする．
7．シェープファイルを入れた3つのファイルが出来ている．
8．QGIS上で「ベクタレイヤ」の追加で，XXX.shpを選択する．
＊詳細はSHAPE4.0CSユーザーズガイドに出ている．
F. 数値の大小をシンボルサイズの大小として表現
1．「レイヤ」-「プロパティ」-「機能拡張」-「サイズ縮尺フィールド」
2．示したいデータが入っているラベルを選択
3．右上部の「大きさ」を適切な数値に調整する
GPS map 62sc (Garmin) & Google Earth
A. Google EarthウェイポイントをGarmin GPSにインポート
1．Google Earth上でウェイポイントを作成する．
2．作成したら保存する．右クリックで「名前を付けて場所を保存」「kml」形式を選択．
3．ソフトウェア「カシミール3D」を立ち上げる．
4．カシミール3Dにkmlデータをインポートする．
＊「ファイル」-「GPS各種ファイルからの追加」で地図上にプロットされる．
5．GPSとPCをUSBケーブルでつなぐ．
6．カシミール3DからGPSへインポートする
＊「通信」-「GPSへアップロード予約」必要なウェイポイントを選択する
＊「Garminドライブ名」をGPSをつないでいるドライブにする．
アップロードを押すとGPSへインポートされる．
B. Garmin GPSデータをGoogle Earthにインポート
1．Google Earthを立ち上げる．
2．GPSとPCをUSBケーブルでつなぐ．
3．Google Earth上の「ツール」-「GPS」-「インポート」
＊ガーミンを選択すること．「ウェイポイント」や「ルート」などが選択出来る．
4．GPSに入力されているすべてのポイントがインポートされる．
C. Garmin GPSデータをカシミールにインポート
1．カシミールを立ち上げる．
2．GPSとPCをUSBケーブルでつなぐ．
3．カシミールの「通信」-「GPSからダウンロード」-「ウェイポイント」で
必要なデータを選択
＊「Garminドライブ名」をGPSとつないでいるドライブにする．
4．「カシミールに保存」を押す．
D. GPSデータをQGISにインポート
1．まずC．GPSデータをカシミールにインポートする．
2．カシミールの「編集」-「GPSデータ編集」-「ウェイポイント」
3．必要なデータをすべて選択して「右クリック」-「ファイルへ書き出し」-「○○.txt形式」で保存する．
＊緯度経度形式が「DEG」であることを確認してファイル名を付ける．
4．保存されたtxtファイルは「,カンマ形式」なので，EXCELを使って，必要なデータを切り出す．
5．ID，緯度，経度としてcsv形式で保存．
E. カシミールデータをGoogle Earthにインポート
1．カシミールを立ち上げる．
2．「ファイル」-「GPS各種ファイルに書き出す」-「ウェイポイント」
3．必要なデータを選択して右クリック-「ファイルへの書き出し」
＊ファイルの種類が「KML」であることを確認してファイル名を付ける．
4． Google Earthを立ち上げる．
5．「ファイル」-「開く」でKMLファイルを選択する．
6．KMLデータがGoogle Earth上にプロットされる．
F. カシミールのデータをCSVにインポート
1．ファイル-GPS各種ファイルに書き出す
2．必要なデータを選択（ウェイポイント，ルート，トラック）
3．ファイル-GPSファイルツール-選択データのCSV書き出し
G. カシミールに多量の緯度経度データをインポート
1． 緯度経度などの情報をcsv形式に保存する．
以下の情報が必要である．
種別,名称,緯度,経度,標高,日付,時刻(UTC),測地系,方位,仰角,区間,アイコン,GPSでの名前,説明
入力例
W1,TEST,35.37340,139.37401,6.0,2000-02-23,05:36:20,WGS84,-9999.0,-9999.0,1,909012,XXX,
2．カシミールを立ち上げる．
3．「ツール」-「GPSファイルツール」-「CSV形式から読み込み」
＊dd.dddd形式かdd.mmss形式か選ぶ．
＊もし，「ツール」-「GPSファイルツール」-「CSV形式から読み込み」が無かったら
プラグインをインストールする．
ファイルは「GPSファイルツール プラグイン Ver 1.2.1」を入れた（2014年12月29日）．
http://www.kashmir3d.com/kash/kashget.html#plugin_gfil
R command
A. 棒グラフの作成1
A-2. 棒グラフの作成2
A-3. 棒グラフの作成3（有義波高と有義波周期）
A-4 折れ線グラフの作成（地形断面図）
B. 2標本の検定（Wilcoxon rank sum test ウィルコクソン順位和検定）
B-2. 2標本の検定（スチューデントのt test）
C. 関連2標本の検定（Wilconxon matched-pairs sigined-ranked test ウィルコクソン符号順位検定）
D. 無相関検定（Peasonの積率相関係数の無相関検定，Spearmanの順位相関係数の無相関検定）
E. 時系列解析・スペクトル解析
F. クラスター解析
G. 計量NMDS解析
H. 重回帰分析
I. 正規分布の検定
J. 3群以上の比較（対応がない，ノンパラメトリック）Kruskal-Wallis検定
K.多重比較（対応がない，ノンパラメトリック）Steel-Dwass検定
L 波浪解析
A. 棒グラフの作成1
1．データの保存
ローカルディスク(:C)の「R」フォルダに保存
各地点5本のトランセクトで2012年と2013年のデータがある．
＊データは算術平均値（coralcovermean.csv）と標準偏差（coralcoversd.csv）を使用する．
2．データの読み込み
> coralcovermean<-read.csv ("C:/R/coralcovermean.csv",header=TRUE,row.names=1)
3．データの確認
> coralcovermean
4．棒グラフの作成
> barplot(coralcovermean)
もし，以下が表示されたら．
--------------------------------------------------------------------
> barplot(coralcovermean)
以下にエラー barplot.default(coralcovermean) :
'height' はベクトルか行列でなければなりません.
データの型を確認する．
> class(coralcovermean)
"data.frame" と出たら，"matrix"に変換する必要がある．
変換方法
> coralcovermean<-as.matrix(coralcovermean)
--------------------------------------------------------------------
軸を見やすいように変更する．
> barplot(coralcovermean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3)
・ylab="Coral cover (%)" ラベルをY軸にラベルを付ける
・ylim=c(0,100) 0から100までに変更する
・beside=TRUE 2010年と2013年のデータをサイトごとに横に並べる
・las=3 X軸のラベルを縦置きにする
標準偏差を付けて表示させる．
いったん，算術平均値をMEANとして置きかえる．そのMEANに対して標準偏差の値を+，-して表示させる．
・angle = 90 バーの角度
・length = 0.03 バーの長さ
> MEAN<-barplot(coralcovermean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3)
>arrows(MEAN, coralcovermean - coralcoversd, MEAN, coralcovermean + coralcoversd, angle = 90, length = 0.03)
>arrows(MEAN, coralcovermean + coralcoversd, MEAN, coralcovermean - coralcoversd, angle = 90, length = 0.03)
A-2. 棒グラフの作成2
1つのサイトのデータのみを標準偏差付きの棒グラフにする（多くのデータが1枚の図あったら見にくい）
例：1列目の「siteA.3m」を図化する．
> siteA.3mmean<-coralcovermean[,1]
> siteA.3msd<-coralcoversd[,1]
> siteA.3mfig<-barplot(siteA.3mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteA.3m")
> arrows(siteA.3mfig, siteA.3mmean - siteA.3msd, siteA.3mfig, siteA.3mmean + siteA.3msd, angle = 90, length = 0.03)
> arrows(siteA.3mfig, siteA.3mmean + siteA.3msd, siteA.3mfig, siteA.3mmean - siteA.3msd, angle = 90, length = 0.03)
------------------------------------------------------------
1枚の図に複数のサイトのデータを載せる
6つのサイトのデータのみを標準偏差付きの棒グラフにする（1枚の図に6つ載せる）
例：1-6列目の「siteA.3m」，「siteA.10m」,「siteB.3m」，「siteB.10m」,「siteC.3m」，「siteC.10m」を図化する．
画面を分割するコマンド
例 3行2列で6つのスペースをまず作る．
par(mfrow=c(3,2))
そのうえで，6つのデータを順番に打ち込む
---------------------------------------------------------
siteA.3mmean<-coralcovermean[,1]
siteA.3msd<-coralcoversd[,1]
siteA.3mfig<-barplot(siteA.3mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteA.3m")
arrows(siteA.3mfig, siteA.3mmean - siteA.3msd, siteA.3mfig, siteA.3mmean + siteA.3msd, angle = 90, length = 0.03)
arrows(siteA.3mfig, siteA.3mmean + siteA.3msd, siteA.3mfig, siteA.3mmean - siteA.3msd, angle = 90, length = 0.03)
siteA.10mmean<-coralcovermean[,2]
siteA.3msd<-coralcoversd[,2]
siteA.10mfig<-barplot(siteA.10mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteA.10m")
arrows(siteA.10mfig, siteA.10mmean - siteA.3msd, siteA.10mfig, siteA.10mmean + siteA.3msd, angle = 90, length = 0.03)
arrows(siteA.10mfig, siteA.10mmean + siteA.3msd, siteA.10mfig, siteA.10mmean - siteA.3msd, angle = 90, length = 0.03)
siteB.3mmean<-coralcovermean[,3]
siteB.3msd<-coralcoversd[,3]
siteB.3mfig<-barplot(siteB.3mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteB.3m")
arrows(siteB.3mfig, siteB.3mmean - siteB.3msd, siteB.3mfig, siteB.3mmean + siteB.3msd, angle = 90, length = 0.03)
arrows(siteB.3mfig, siteB.3mmean + siteB.3msd, siteB.3mfig, siteB.3mmean - siteB.3msd, angle = 90, length = 0.03)
siteB.10mmean<-coralcovermean[,4]
siteB.10msd<-coralcoversd[,4]
siteB.10mfig<-barplot(siteB.10mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteB.10m")
arrows(siteB.10mfig, siteB.10mmean - siteB.10msd, siteB.10mfig, siteB.10mmean + siteB.10msd, angle = 90, length = 0.03)
arrows(siteB.10mfig, siteB.10mmean + siteB.10msd, siteB.10mfig, siteB.10mmean - siteB.10msd, angle = 90, length = 0.03)
siteC.3mmean<-coralcovermean[,5]
siteC.3msd<-coralcoversd[,5]
siteC.3mfig<-barplot(siteC.3mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteC.3m")
arrows(siteC.3mfig, siteC.3mmean - siteC.3msd, siteC.3mfig, siteC.3mmean + siteC.3msd, angle = 90, length = 0.03)
arrows(siteC.3mfig, siteC.3mmean + siteC.3msd, siteC.3mfig, siteC.3mmean - siteC.3msd, angle = 90, length = 0.03)
siteC.10mmean<-coralcovermean[,6]
siteC.10msd<-coralcoversd[,6]
siteC.10mfig<-barplot(siteC.10mmean,ylab="Coral cover (%)",ylim=c(0,100),beside=TRUE,las=3,xlab="siteC.10m")
arrows(siteC.10mfig, siteC.10mmean - siteC.10msd, siteC.10mfig, siteC.10mmean + siteC.10msd, angle = 90, length = 0.03)
arrows(siteC.10mfig, siteC.10mmean + siteC.10msd, siteC.10mfig, siteC.10mmean -siteC.10msd, angle = 90, length = 0.03)
---------------------------------------------------------
A-3. 棒グラフの作成3（有義波高と有義波周期）
1．データの保存
ローカルディスク(:C)の「R」フォルダにcsv形式で保存（ファイル名：waveishigaki.csv）
＊1996から2013までの1日2回のデータ 1996/8/1/0,1996/8/1/12・・・
＊年の区切りをわかりやすくするために各年の最初にX1 20m，25sのダミーデータを入れている
2．データの読み込み
waveishigaki<-read.csv ("C:/R/waveishigaki.csv",header=T, row.names=1)
3．行と列の並び替え
ishigakiwave<-t(waveishigaki)
4．ラベルの関連付け
height<-ishigakiwave[1,]
period<-ishigakiwave[2,]
5．棒グラフの作成
barplot(height)
barplot(period)
6．軸を見やすいように変更する．
barplot(height,main="Hs: Ishigaki O point", xlab= "Year",ylab="Hs (m)",ylim=c(0,20),cex.main=2,cex.lab=1,cex.names=0.5,cex.axis=1)
barplot(period,main="Ts: Ishigaki O point", xlab= "Year",ylab="Ts (s)",ylim=c(0,20),cex.main=2,cex.lab=1,cex.names=0.5,cex.axis=1)
・main="Ishigaki O point: Hs" タイトルを付ける
・xlab= "Year" X軸にラベルを付ける
・ylab="Hs (m)" Y軸にラベルを付ける
・ylim=c(0,20) Y軸の数値を0から20に設定する
・cex.main=2 タイトルのフォントサイズを2に設定する
・cex.lab=1 軸のラベルのフォントサイズを1に設定する
・cex.names=0.5 X軸のフォントサイズを0.5に設定する
・cex.axis=1 Y軸のフォントサイズを1に設定する
-----------------------------------------------------
ヒストグラムの作成（有義波高と有義波周期）
1．データの保存
ローカルディスク(:C)の「R」フォルダにcsv形式で保存（ファイル名：waveishigaki2.csv）
＊1996から2013までの1日2回のデータ 1996/8/1/0,1996/8/1/12・・・
2．データの読み込み
waveishigaki2<-read.csv ("C:/R/waveishigaki2.csv",header=T, row.names=1)
3．行と列の並び替え
ishigakiwave2<-t(waveishigaki2)
4．ラベルの関連付け
height<-ishigakiwave2[1,]
period<-ishigakiwave2[2,]
5．ヒストグラムの作成
hist(height)
hist(period)
折れ線グラフの作成（Wave set up）
6．見やすいように変更する．
hist(height,main="Hs: Ishigaki O point", xlab= "Hs (m)",ylab="Frequency",ylim=c(0,2000),cex.main=2,cex.lab=1,cex.axis=1,breaks=seq(0,20,0.2),xaxs="i",yaxs="i")
hist(period,main="Ts: Ishigaki O point", xlab= "Hs (m)",ylab="Frequency",ylim=c(0,3000),cex.main=2,cex.lab=1,cex.axis=1,breaks=seq(0,20,0.2),xaxs="i",yaxs="i")
・main="Ishigaki O point: Hs" タイトルを付ける
・xlab= "Hs (m)" X軸にラベルを付ける
・ylab="Frequency" Y軸にラベルを付ける
・ylim=c(0,2000) Y軸の数値を0から2000に設定する
・cex.main=2 タイトルのフォントサイズを2に設定する
・cex.lab=1 軸のラベルのフォントサイズを1に設定する
・cex.axis=1 Y軸のフォントサイズを1に設定する
・breaks=seq(0,20,0.2) 階級幅を0から20まで0.2間隔で区切る
・xaxs="i",yaxs="i" X軸とY軸を最小値で交差させる
・xaxp=c(0,2000,500) X軸の0から2000を500ごとに目盛を打つ
＊ height とPeriodの平均値と中央値を求める
算術平均値
meanheight<-mean(height)
meanperiod<-mean(period)
中央値
medianheight<-median(height)
medianperiod<-median(period)
------------------------------------------
ヒストグラムの作成（Skewness）
1．データの保存
ローカルディスク(:C)の「R」フォルダにcsv形式で保存（ファイル名：akaishiskewness.csv）
＊degradeとhealthy（64ケース）の海岸のwave skewness
2．データの読み込み
akaishiskewness<-read.csv ("C:/R/akaishiskewness.csv",header=T, row.names=1)
3．行と列の並び替え
skewnessakaishi<-t(akaishiskewness)
4．ラベルの関連付け
skdegraded<-skewnessakaishi[1,]
skhealthy<-skewnessakaishi[2,]
5．ヒストグラムの作成
hist(skdegraded)
hist(skhealthy)
6．見やすいように変更する．
hist(skdegraded,xlab= "Skewness(-)",ylab="Frequency",ylim=c(0,50),cex.main=2,cex.lab=1,cex.axis=1,breaks=seq(-1,1,0.2),xaxs="i",yaxs="i")
hist(skhealthy,xlab= "Skewness(-)",ylab="Frequency",ylim=c(0,50),cex.main=2,cex.lab=1,cex.axis=1,breaks=seq(-1,1,0.2),xaxs="i",yaxs="i")
・main="Ishigaki O point: Hs" タイトルを付ける
・xlab= "Hs (m)" X軸にラベルを付ける
・ylab="Frequency" Y軸にラベルを付ける
・ylim=c(0,2000) Y軸の数値を0から2000に設定する
・cex.main=2 タイトルのフォントサイズを2に設定する
・cex.lab=1 軸のラベルのフォントサイズを1に設定する
・cex.axis=1 Y軸のフォントサイズを1に設定する
・breaks=seq(0,20,0.2) 階級幅を0から20まで0.2間隔で区切る
・xaxs="i",yaxs="i" X軸とY軸を最小値で交差させる
・xaxp=c(0,2000,500) X軸の0から2000を500ごとに目盛を打つ
＊ height とPeriodの平均値と中央値を求める
算術平均値
meanheight<-mean(height)
meanperiod<-mean(period)
中央値
medianheight<-median(height)
medianperiod<-median(period)
A-4 折れ線グラフの作成（地形断面図）
1．データの保存
Rのある以下に保存する
C:\Users\chuki\Dropbox\My PC (DESKTOP-5MLLL02)\Documents\R
Data1----（1行目）、Length（2行目）とDepth（3行目）のデータがある．
2．データの読み込み
shirahocross<-read.csv ("R/shirahocross.csv",header=TRUE,row.names=1)
3．データ型を確認する
class(shirahocross)
[1] "data.frame"
＊matrixに変換する必要がある
4．matrixへの変換
shirahocross2<-as.matrix(shirahocross)
5．ラベルの関連ずけ
length<-shirahocross2[1,]
depth<-shirahocross2[2,]
6．折れ線グラフの作成（広範囲）
plot(length,depth,type="b")
＊type="b"でラベルつき
＊type="l"で線のみ
7．グラフの調整
plot(length,depth,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,2100))
＊X軸：0～2100 xlim=c(0,2100)
＊Y軸：10～-60 ylim=c(10,60)
＊X軸ラベル：xlab="Length [m]"
＊Y軸ラベル：ylab="Water depth/ elevation [m]"
＊X軸、Y軸の交点の設定：axes = FALSE)
8．Y軸、Y軸のラベル設定
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=2100)
-----------------------------------------------------------
9.6．折れ線グラフの作成（詳細）
plot(length,depth,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-6,5),xlim=c(0,1600))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=1600)
10. 複数のグラフの作図
par(mfrow=c(2,1))
plot(length,depth,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,2100))
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=2100)
text(500,-50,"Shiraho reef")
plot(length,depth,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-5,5),xlim=c(0,1600))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=1600)
text(500,-4,"Shiraho reef")
11. 他の4地点のデータ入力
ikeicross<-read.csv ("R/ikeicross.csv",header=TRUE,row.names=1)
ikeicross2<-as.matrix(ikeicross)
lengthike<-ikeicross2[1,]
depthike<-ikeicross2[2,]
oharacross<-read.csv ("R/oharacross.csv",header=TRUE,row.names=1)
oharacross2<-as.matrix(oharacross)
lengthoha<-oharacross2[1,]
depthoha<-oharacross2[2,]
tomoricross<-read.csv ("R/tomoricross.csv",header=TRUE,row.names=1)
tomoricross2<-as.matrix(tomoricross)
lengthtom<-tomoricross2[1,]
depthtom<-tomoricross2[2,]
aridacross<-read.csv ("R/aridacross.csv",header=TRUE,row.names=1)
aridacross2<-as.matrix(aridacross)
lengthari<-aridacross2[1,]
depthari<-aridacross2[2,]
12．多地点の地形断面図の作図（広範囲＆詳細）
par(mfrow=c(2,1))
plot(length,depth,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,2100))
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=2100)
text(500,-50,"Shiraho reef")
plot(length,depth,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-5,5),xlim=c(0,1600))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=1600)
text(500,-4,"Shiraho reef")
par(mfrow=c(2,1))
plot(lengthike,depthike,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,1200))
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=1200)
text(500,-50,"Ikei reef")
plot(lengthike,depthike,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-5,5),xlim=c(0,500))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=500)
text(450,-4,"Ikei reef")
par(mfrow=c(2,1))
plot(lengthoha,depthoha,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,1600))
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=1600)
text(500,-50,"Ohara reef")
plot(lengthoha,depthoha,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-5,5),xlim=c(0,900))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=900)
text(600,-4,"Ohara reef")
par(mfrow=c(2,1))
plot(lengthtom,depthtom,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,4200))
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=4200)
text(500,-50,"Tomori reef")
plot(lengthtom,depthtom,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-5,5),xlim=c(0,700))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=700)
text(400,-4,"Tomori reef")
par(mfrow=c(2,1))
plot(lengthari,depthari,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-60,10),xlim=c(0,4000))
axis(1, pos=-60)
axis(2, pos=0)
axis(4, pos=4000)
text(500,-50,"Arida site")
plot(lengthari,depthari,type="l",xlab="Length [m]",ylab="Water depth/ elevation [m]",axes = FALSE, ylim=c(-5,5),xlim=c(0,500))
axis(1, pos=-5)
axis(2, pos=0)
axis(4, pos=500)
text(300,-4,"Arida site")
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
B. 2標本の検定（Wilcoxon rank sum test ウィルコクソン順位和検定）
1．データをcsv形式で用意する．
2．データの読み込み
deheal<-read.csv ("C:/R/degradedhealthy.csv",header=TRUE,row.names=1)
3．データの確認
deheal
4．ラベルの関連付け
degraded1<-deheal[1,]
healthy1<-deheal[2,]
5．matrixデータに変更
degraded<-as.matrix(degraded1)
healthy<-as.matrix(healthy1)
6．パッケージ"exactRankTests"をインストールする．
パッケージを読み込む
library(exactRankTests)
*既にローカルに保存済み．
必要ならば，
install.packages("exactRankTests")
7．Wilcoxon rank sum testを実施
wilcox.exact(degraded, healthy,paired=F)
*paired=Fを設定することでウィルコクソン順位和検定（Wilcoxon rank sum test）（ノンパラメトリック）
もし，paired=Tを設定するとウィルコクソン符号順位検定（Wilcoxon signed rant test）になる
8．結果
##############結果###################
Asymptotic Wilcoxon rank sum test
data: degraded and healthy
W = 2514, p-value = 0.02634
alternative hypothesis: true mu is not equal to 0
#####################################
H0帰無仮説=差なし
H1対立仮説=差あり（有意水準5%とすれば: p-value < 0.05）
＊＊0.02634は<0.05なので，degradedとhealthyには差がある
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
B-2. 2標本の検定（スチューデントのt test）
石垣島明石データ波高11-15m、波周期14-18sのデータ（Case4-256, 1136-1150,1286-1300）
deheal<-read.csv ("C:/R/newtest.csv",header=TRUE,row.names=1)
deheal
degraded<-deheal[1,]
healthy<-deheal[2,]
degraded<-as.matrix(degraded)
healthy<-as.matrix(healthy)
shapiro.test(degraded)
shapiro.test(healthy)
wilcox.exact(degraded, healthy,paired=F)
--------正規性の検定--------
> shapiro.test(degraded)
Shapiro-Wilk normality test
data: degraded
W = 0.96056, p-value = 0.01547
#p>0.01で帰無仮説（正規分布）が保留．→正規分布である
> shapiro.test(healthy)
Shapiro-Wilk normality test
data: healthy
W = 0.97795, p-value = 0.1882
#p>0.01で帰無仮説（正規分布）が保留．→正規分布である
--------等分散性の検定--------
> var.test(degraded,healthy)
F test to compare two variances
data: degraded and healthy
F = 1.6764, num df = 78, denom df = 78, p-value = 0.02372
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
1.071992 2.621451
sample estimates:
ratio of variances
1.676358
#p>0.01で帰無仮説（等分散）が保留．→等分散である
--------スチューデントのt検定--------
> t.test(degraded,healthy,var.equal=T,paired=F)
Two Sample t-test
data: degraded and healthy
t = 3.3676, df = 156, p-value = 0.0009552
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
0.05338073 0.20484712
sample estimates:
mean of x mean of y
2.235316 2.106203
#p<0.01で帰無仮説（2群間の平均に差が無い）が棄却された．→2群に差がある。
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
C. 関連2標本の検定（Wilconxon matched-pairs sigined-ranked test ウィルコクソン符号順位検定）
パッケージ"exactRankTests"をインストールする．
ツールバーからサイトを選びインストールする．
中身の確認
> help(package=exactRankTests)
パッケージを読み込む
> library(exactRankTests)
データの準備
> coralcover<-read.csv ("C:/R/coralcover.csv",header=TRUE)
検定（標本数が5つなので片側(右側)検定した）
x：siteAの3mの2010年
y：siteAの3mの2013年
alternative="g"：右側検定（右側だと”g”，左側だと"ｌ"）
paired=T：符号順位検定（paired=Fにすると順位和検定）
x<-coralcover[,1]
y<-coralcover[,5]
wilcox.exact(x,y,paired=T,alternative="g")
---------------------------
Exact Wilcoxon signed rank test
data: x and y
V = 15, p-value = 0.03125
alternative hypothesis: true mu is greater than 0
---------------------------
p値が0.03125なので，有意水準5%とすればデータXとデータYには差がある．
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
D. 無相関検定（Peasonの積率相関係数の無相関検定，Spearmanの順位相関係数の無相関検定）
1．データをcsv形式で用意する．例：14feb23瀬底コユビ直径データ.xls
2．データの読み込み
coralsize<-read.csv ("C:/R/coralsize.csv",header=TRUE,row.names=1)
3．データの確認
coralsize
4．ラベルの関連付け
diameter1<-coralsize[1,]
weight1<-coralsize[2,]
5．matrixデータに変更
diameter<-as.matrix(diameter1)
weight<-as.matrix(weight1)
6.Peasonの積率相関係数の無相関検定(パラメトリック法)
cor.test(diameter, weight, method="pearson")
----------結果--------------------------------------
Pearson's product-moment correlation
data: diameter and weight
t = 5.9516, df = 13, p-value = 4.812e-05
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.6105137 0.9509222
sample estimates:
cor
0.8552901
----------------------------------------------------
H0帰無仮説=相関はない
H1対立仮説=相関はある（p-value < 0.01）
＊＊p-value = 4.812e-05 ＜0.01なので，diameterとweightには相関がある
7.Spearmanの順位相関係数の無相関検定(ノンパラメトリック法)
cor.test(diameter, weight, method="spearman")
----------結果--------------------------------------
Spearman's rank correlation rho
data: diameter and weight
S = 57.7528, p-value = 5.907e-06
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
0.89687
警告メッセージ：
In cor.test.default(diameter, weight, method = "spearman") :
タイのため正確な p 値を計算することができません
----------------------------------------------------
H0帰無仮説=相関はない
H1対立仮説=相関はある（p-value < 0.01）
＊＊p-value = 5.907e-06 ＜0.01なので，diameterとweightには相関がある
------------------------------------------------------------------------------------------
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
E. 時系列解析・スペクトル解析
1．時系列データの描画（時系列，スペクトル）
*白保の3600sの水位変化（ファイル名：specsh39.csv）
2．データの読み込み
specsh39<-read.csv("specsh39.csv",header=TRUE)
3．データの確認
specsh39
4．ラベルの関連付け
ANS<-specsh39[,2]
＊時間と入射の水位だけ
reef<-specsh39[,3]
＊時間と礁地の水位だけ
5．データを時系列に変換
ts.ANS <- ts(ANS)
ts.reef <- ts(reef)
6.時系列プロット
plot(ts.ANS,xlab="Time [s]",ylab="Water level[m]",las=1)
＊規則波が確認できる．
plot(ts.reef,xlab="Time [s]",ylab="Water level[m]",las=1)
**いろいろな波があることを確認できる
7．パワースペクトル密度図
spectrum(ANS,method="ar")
ar：パラメトリックの自己回帰法
spectrum(ANS,method="pgram")
pgram：ノンパラメトリックのピリオドグラム法
＊0.05Hz (20s)入射波をスペクトル図で確認できる．
spectrum(reef,method="ar")
spectrum(reef,method="pgram")
*礁池では明瞭なピークを確認できない．
---------------------
1．時系列データの描画（時系列:Melekeok）
*現在1800s間の水位変化（ファイル名：melekeokwaterlevel.csv）
2．データの読み込み
melekeokwaterlevel<-read.csv("C:/R/melekeokwaterlevel.csv",header=TRUE,row.names=1)
3．データの確認
melekeokwaterlevel
4．データを時系列に変換
melewater<-as.ts(melekeokwaterlevel)
5．ラベルの関連付け
pre1850mHTSSsetup<-melewater[1,]
pre1850mHTSS<-melewater[2,]
pre1850mHTsetup<-melewater[3,]
pre1850mHT<-melewater[4,]
pre1850mMSLSSsetup<-melewater[5,]
pre1850mMSLSS<-melewater[6,]
6.時系列プロット
ts.plot(pre1850mHTSSsetup)
6-1．必要に応じてラベルと軸の調整
ts.plot(pre1850mHTSS,xlab="Time (s)", ylab="Water level (m)",ylim=c(0,4),main="pre 8.7m13sMSL+1.8m (HT+SS) water level")
ts.plot(pre1850mMSLSS,xlab="Time (s)", ylab="Water level (m)",ylim=c(0,4),main="pre 8.7m13sMSL+1.8m (HT+SS) water level")
7.一枚の図に複数の図を描画
par(mfrow=c(2,2),xaxs="i",yaxs="i")
ts.plot(pre1850mHTSSsetup,xlab="Time (s)", ylab="Water level (m)",ylim=c(-1,4),main="pre 8.7m13sMSL+1.8m (HT+SS) wave setup")
ts.plot(pre1850mHTSS,xlab="Time (s)", ylab="Water level (m)",ylim=c(-1,4),main="pre 8.7m13sMSL+1.8m (HT+SS) water level")
ts.plot(pre1850mHTsetup,xlab="Time (s)", ylab="Water level (m)",ylim=c(-1,4),main="pre 8.7m13sMSL+0.8m (HT) wave setup")
ts.plot(pre1850mHT,xlab="Time (s)", ylab="Water level (m)",ylim=c(-1,4),main="pre 8.7m13sMSL+0.8m (HT) water level")
---------------------
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
F. クラスター解析
1．データ読み込み
allcoral<-read.csv("C:/R/allcoral.csv",header=T, row.names=1)
#ラベル付でデータを入力
2．行と列の並び替え
coralall<-t(allcoral)
3．調査地ごとの距離関数が計算（この場合はユークリッド距離）
coralall.d<-dist(coralall)
4．Bray-Curis距離を扱えるパッケージの”vegan”をインストールする
install.packages("vegan")
5．veganがインストールされているのか確認
library(vegan)
6．veganにある関するvegdistでbray（Bray-Curtis）を選択
coralall.d<-vegdist(coralall,"bray")
7．作図
plot(hclust(vegdist(coralall,"bray"),"average"),hang=-1)
#以上のまとめ
#brayで距離を算出し，average（郡平均法）でクラスターを作成し，図をプロットする．
#枝の高さを揃えるために"hang=-1"を入力．
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
G. 計量NMDS解析
1．データ読み込み
takanaka<-read.csv("C:/R/takanaka.csv",header=T,row.names=1)
#ラベル付でデータを入力（宝島と中之島データ）
2．行と列の並び替え
nakataka<-t(takanaka)
#行と列の並び替え
3．調査地ごとの距離関数が計算（この場合はユークリッド距離）
nakataka.d<-dist(nakataka)
4．Bray-Curis距離を扱えるパッケージの”vegan”をインストールする
install.packages("vegan")
5．veganがインストールされているのか確認
library(vegan)
6．veganにある関するvegdistでbray（Bray-Curtis）を選択
nakataka.d<-vegdist(nakataka,"bray")
7．cmdscale関数
nakataka.cmd<-cmdscale(nakataka.d,k=2,eig=FALSE)
8．作図
plot(nakataka.cmd)
9．ラベルを表示
text(nakataka.cmd,row.names(nakataka.cmd))
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
H. 重回帰分析
1．データ読み込み
multi <- read.table("C:/R/multi.txt", header = TRUE)
2．データ関連付け
waveheight <- multi$waveheight
waveperiod <- multi$waveperiod
waterlevel <- multi$waterlevel
reefgrowth <- multi$reefgrowth
reduction <- multi$reduction
#waveheigh, waveperiod, waterlevel, reefgrowthを独立変数，reductionを従属変数とする
3．ヒストグラムでデータを確認
par(mfrow = c(1, 5))
hist(waveheight)
hist(waveperiod)
hist(waterlevel)
hist(reefgrowth)
hist(reduction)
4．単位が異なるので標準化しておく
waveheightstd <- scale(waveheight)
waveperiodstd <- scale(waveperiod)
waterlevelstd <- scale(waterlevel)
reefgrowthstd <- scale(reefgrowth)
5．重回帰分析
result3 <- lm(reduction ~ waveheightstd + waveperiodstd + waterlevelstd + reefgrowthstd)
summary(result3)
---------------------
Call:
lm(formula = reduction ~ waveheightstd + waveperiodstd + waterlevelstd +
reefgrowthstd)
Residuals:
Min 1Q Median 3Q Max
-10.9253 -2.6242 0.3345 2.7298 7.7218
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 65.6550 0.1530 429.079 < 2e-16 ***
waveheightstd 1.1374 0.1533 7.419 5.13e-13 ***
waveperiodstd -2.6266 0.1532 -17.143 < 2e-16 ***
waterlevelstd -3.8343 0.1568 -24.453 < 2e-16 ***
reefgrowthstd 2.3959 0.1567 15.288 < 2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 3.428 on 497 degrees of freedom
Multiple R-squared: 0.6765, Adjusted R-squared: 0.6739
F-statistic: 259.9 on 4 and 497 DF, p-value: < 2.2e-16
----------------------------
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
I. 正規分布の検定
1．データ読み込み
givekeep2 <- read.table("C:/R/givekeep2.txt", header = TRUE)
2．データ関連付け
give <- givekeep2$give
keep <- givekeep2$keep
3．検定
ks.test(give, keep)
4．結果
Two-sample Kolmogorov-Smirnov test
data: give and keep
D = 0.2863, p-value = 2.974e-09
alternative hypothesis: two-sided
警告メッセージ：
In ks.test(give, keep) :
タイがあるため、正しい p 値を計算することが
#p<0.01で帰無仮説（正規分布）が棄却された．
#正規分布では無い（独立2群）----Mann-Whitney U検定（Wilcoxon 順位和検定 Wilcoxon rank sum test)
#正規分布では無い（関連2群）----Wilcoxon符号順位検定（Wilcoxon signed rank test)
#Mann-Whitney検定（Wilcoxon 順位和検定）------> wilcox.test(give,keep)
Wilcoxon rank sum test with continuity correction
data: give and keep
W = 20405, p-value = 9e-11
alternative hypothesis: true location shift is not equal to 0
#Wilcoxon符号順位検定(両側検定)-----> wilcox.exact(give,keep, paired=T)
Asymptotic Wilcoxon signed rank test
data: give and keep
V = 1931.5, p-value < 2.2e-16
alternative hypothesis: true mu is not equal to 0
#Wilcoxon符号順位検定（片側検定）-----> wilcox.exact(give,keep, paired=T,alternative="l")
Asymptotic Wilcoxon signed rank test
data: give and keep
V = 1931.5, p-value < 2.2e-16
alternative hypothesis: true mu is less than 0
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
J. 3群以上の比較（対応がない，ノンパラメトリック）Kruskal-Wallis検定
1．データ読み込み
evelyn <- read.table("C:/R/evelyn.txt", header = TRUE)
2．データの確認
evelyn
#データはアンスタック形式
3．スタック形式データに変換
evelyn1 <- stack(evelyn)
4．ラベルをわかりやすく変更
colnames(evelyn1) <- c("photo", "group")
5．Kruskal-Wallis検定
kruskal.test(evelyn1$photo ~ evelyn1$group)
6．結果
Kruskal-Wallis rank sum test
data: evelyn1$photo by evelyn1$group
Kruskal-Wallis chi-squared = 9.374, df = 3, p-value = 0.02471
p<0.05水準でgroup間に有意差がある
7.箱ひげ図で確認
boxplot(photo~group,data=evelyn1, xlab="Origin Sewer",ylab="photo from June to Oct")
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
K.多重比較（対応がない，ノンパラメトリック）Steel-Dwass検定
1．データ読み込み
evelyn <- read.table("C:/R/evelyn.txt", header = TRUE)
2．データの確認
evelyn
#データはアンスタック形式
3．スタック形式データに変換
evelyn1 <- stack(evelyn)
4．ラベルをわかりやすく変更
colnames(evelyn1) <- c("photo", "group")
5．Steel-Dwass検定
http://aoki2.si.gunma-u.ac.jp/R/Steel-Dwass.htmlにて詳細
#インストール必要あり
#groupを数値しておく必要があるため，as.numeric()としている
Steel.Dwass(evelyn1$photo, as.numeric(evelyn1$group))
6．結果
t p
1:2 1.9012827 0.22743028
1:3 0.2816715 0.99220791
1:4 0.4929252 0.96069350
2:3 2.0684478 0.16349605
2:4 2.8564280 0.02225976
3:4 1.4117977 0.49188417
#2:4とはmed：contのことである．
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
L. 波浪解析
CADMAS-SURFのデータの波浪解析
1．"data.tran"の不要データを削除してtxt形式で保存（例：casea001.txt）
2. textファイルをＲにてcsv形式で保存
casea001<-read.table ("R/CADMAS-SURF/Arida/casea001.txt")
colnames(casea010)<-c("Time", "ANS","m3600m","m3700m","m3820m","m3870m","m3920m","m3930m")
＊ヘッダの追加（元々、data.tranにヘッダが無いので、csvファイルに追加しておく。）
*例 A列（1秒ごと） B列（1000mの水位）C列（1100mの水位）D列（1200mの水位）
*A列に時間（Time），B列以降は水位データ（ラベルは便宜上"m1000m"としておくと読み込みエラーが起きない）。
write.table(casea001,file="R/CADMAS-SURF/Arida/casea001.csv",quote=F, sep=",", row.names=FALSE, col.names=FALSE)
＊データは同じフォルダに保存されている。
3. パッケージのインストール
install.packages("oceanwaves")
4. パッケージの読み込み
library(oceanwaves)
5．データの読み込み
casea001csv<-read.csv ("R/CADMAS-SURF/Arida/casea001.csv",header=TRUE,row.names=1)
6.必要な期間を抽出して名前を付ける
casea001e<-casea001csv[1801:3600,]
＊例えば計算後半の1800s分のデータとして1801-3600を切り出して"casea001e"としてラベルをつける。
7.波浪解析
casea001eANS<-waveStatsZC(casea001e$ANS, Fs = 1)
*1秒ごとのデータなので、サンプリング周波数は1Hzとしておく（Fs＝1）
8．データの保存
write.table(casea001eANS,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=TRUE)
*出力データから "" を省くために、quote=F をつける
*カンマ区切り sep="," としておく
*col.names=TRUE としておき、ラベルを付ける
*上書きをOKとするため、append=T としておく。
9.値の出力
6項目が出力される：Hsig Hmean H10 Hmax Tmean Tsig．連続処理
10.連続処理
波浪解析は測定地点ごとに値が出るため、一測線上の複数地点のデータを一つのファイルで出力できたほうが良い。
casea001eANS<-waveStatsZC(casea001e$ANS, Fs = 1)
write.table(casea001eANS,file="casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=TRUE)
casea001em3600m<-waveStatsZC(casea001e$m3600m, Fs = 1)
write.table(casea001em3600m,file="casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3700m<-waveStatsZC(casea001e$m3700m, Fs = 1)
write.table(casea001em3700m,file="casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3820m<-waveStatsZC(casea001e$m3820m, Fs = 1)
write.table(casea001em3820m,file="casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3870m<-waveStatsZC(casea001e$m3870m, Fs = 1)
write.table(casea001em3870m,file="casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3920m<-waveStatsZC(casea001e$m3920m, Fs = 1)
write.table(casea001em3920m,file="casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
＊２つ目のwrite.table以降はヘッダが不要のため、col.names=FALSE としている。
11. 測定地点ごとのラベルを追加する。
casea001out<-read.csv ("R/CADMAS-SURF/Arida/casea001out.csv", header=FALSE)
rownames(casea001out)<-c("casea001","ANS","m3600m","m3700m","m3820m","m3870m","m3920m")
write.table(casea001out,file="R/CADMAS-SURF/Arida/casea001outout.csv",quote=F, sep=",",)
12．まとめると、1つのtxtファイルから波浪計算出力コード：
casea001<-read.table ("R/CADMAS-SURF/Arida/casea001.txt")
colnames(casea001)<-c("Time", "ANS","m3600m","m3700m","m3820m","m3870m","m3920m","m3930m")
write.table(casea001,file="R/CADMAS-SURF/Arida/casea001.csv",quote=F, sep=",", row.names=FALSE, col.names=TRUE)
casea001csv<-read.csv ("R/CADMAS-SURF/Arida/casea001.csv",header=TRUE,row.names=1)
casea001e<-casea001csv[1801:3600,]
casea001eANS<-waveStatsZC(casea001e$ANS, Fs = 1)
write.table(casea001eANS,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=TRUE)
casea001em3600m<-waveStatsZC(casea001e$m3600m, Fs = 1)
write.table(casea001em3600m,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3700m<-waveStatsZC(casea001e$m3700m, Fs = 1)
write.table(casea001em3700m,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3820m<-waveStatsZC(casea001e$m3820m, Fs = 1)
write.table(casea001em3820m,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3870m<-waveStatsZC(casea001e$m3870m, Fs = 1)
write.table(casea001em3870m,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001em3920m<-waveStatsZC(casea001e$m3920m, Fs = 1)
write.table(casea001em3920m,file="R/CADMAS-SURF/Arida/casea001out.csv",quote=F, sep=",", append=T, row.names=FALSE, col.names=FALSE)
casea001out<-read.csv ("R/CADMAS-SURF/Arida/casea001out.csv", header=FALSE)
rownames(casea001out)<-c("casea001","ANS","m3600m","m3700m","m3820m","m3870m","m3920m")
write.table(casea001out,file="R/CADMAS-SURF/Arida/outcasea001out.csv",quote=F, sep=",",)
13．複数ファイルから必要なセルを切り出して出力
casea001010<-cbind(casea001out[5:7,1],casea002out[5:7,1],casea003out[5:7,1],casea004out[5:7,1],casea005out[5:7,1],casea006out[5:7,1],casea007out[5:7,1],casea008out[5:7,1],casea009out[5:7,1],casea010out[5:7,1])
*10個のファイルについて、必要な3つのデータ（"m3820m","m3870m","m3920m"についてのHs）を切り出して、横に並べる。
*ファイル名は手入力ではなく、Excelの関数concatを用いる。
14．わかりやすいように、ケースごとのラベルを追加する。
colnames(casea001010)<-c("casea001out","casea002out","casea003out","casea004out","casea005out","casea006out","casea007out","casea008out","casea009out","casea010out")
15. 行列を入れ替えて、一旦、出力する。
outcasea001010<-t(casea001010)
write.table(outcasea001010,file="R/CADMAS-SURF/Arida/outcasea001010.csv",quote=F, sep=",",append=T,row.names=TRUE, col.names=FALSE)
16．一旦、出力したファイルを読みだして、わかりやすいように測定地点ごとのラベルを追加して最終出力する。
outoutcasea001010<-read.csv ("R/CADMAS-SURF/Arida/outcasea001010.csv", header=FALSE)
colnames(outoutcasea001010)<-c("casea001010","m3820m","m3870m","m3920mshore")
write.table(outoutcasea001010,file="R/CADMAS-SURF/Arida/outoutcasea001010.csv",quote=F, sep=",",row.names=FALSE, col.names=TRUE)
17．まとめると、以下になる。
casea001010<-cbind(casea001out[5:7,1],casea002out[5:7,1],casea003out[5:7,1],casea004out[5:7,1],casea005out[5:7,1],casea006out[5:7,1],casea007out[5:7,1],casea008out[5:7,1],casea009out[5:7,1],casea010out[5:7,1])
colnames(casea001010)<-c("casea001out","casea002out","casea003out","casea004out","casea005out","casea006out","casea007out","casea008out","casea009out","casea010out")
outcasea001010<-t(casea001010)
write.table(outcasea001010,file="R/CADMAS-SURF/Arida/outcasea001010.csv",quote=F, sep=",",append=T,row.names=TRUE, col.names=FALSE)
outoutcasea001010<-read.csv ("R/CADMAS-SURF/Arida/outcasea001010.csv", header=FALSE)
colnames(outoutcasea001010)<-c("casea001010","m3820m","m3870m","m3920mshore")
write.table(outoutcasea001010,file="R/CADMAS-SURF/Arida/outoutcasea001010.csv",quote=F, sep=",",row.names=FALSE, col.names=TRUE)
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
その他
------------------------------------------------
1．フォントの確認
> windowsFonts()
$serif
[1] "TT Times New Roman"
$sans
[1] "TT Arial"
$mono
[1] "TT Courier New"
2．フォントサイズの確認
> par("ps")
[1] 12
3．フォントサイズの変更（20ポイントに変更）
> par(ps=20)
4．余白の設定
# 下・左・上・右の順で余白を設定
> par(oma = c(3, 0, 0, 0))
下に3インチの余白が出来る．
------------------------------------------------
4．任意のデータの抽出
例：1行目の1番目のデータを抽出
> coralcover[1,1]
＊「35.5」が出力される
5．各サイトにおける5トランセクトの算術平均値の算出
例：「siteA.2010.3m」と「siteA.2013.3m」を計算
> meansiteA.2010.3m<-mean(coralcover$siteA.2010.3m)
> meansiteA.2010.3m
＊「37.222」が出力される
> meansiteA.2013.3m<-mean(coralcover$siteA.2013.3m)
> meansiteA.2013.3m
＊「1.608」が出力される
> meansiteA.2013.3m<-mean(coralcover$siteA.2013.3m)
> meansiteA.2013.3m
＊「1.608」が出力される
6．各サイトにおける5トランセクトの標準偏差の算出
例：「siteA.2010.3m」と「siteA.2013.3m」を計算
> sdsiteA.2010.3m<-sd(coralcover$siteA.2010.3m)
> sdsiteA.2010.3m
＊「10.81175」が出力される
> sdsiteA.2013.3m<-sd(coralcover$siteA.2013.3m)
> sdsiteA.2013.3m
＊「1.941216」が出力される
TSMaster8
A. 水深データの準備
1．A列に時間，B列に水深データを入力して，CSV形式にて保存．
＊1行目はヘッダ：A1が「Time」，B1が「水深[m]」
B. ソフト「TSMaster8」を開く
1．ファイルを開く
＊「ファイル」-「読み込み」-「エクセル/カンマ区切りファイル」
C. 読み込み設定
1．「計測開始時刻」に適切な時間を入力する
2．「計測時間間隔」に入力する．0.5s間隔の測定なら，500 msec．
3．Readの1列目はTimeなので出力する必要はない．「Yes」を「No」に変更する．
「右クリック」-「読み込み解除」でNoになる
4．「先頭スキップ行数」 1を選択する．1列目が見えなくなる．
5．Aの係数を「1」から「0.001」に変更する．これはデータを四捨五入して計算を速める．
6．「表のアップデート」を押す
7．「データの読み込み」を押す
8．時系列ファイルの保存
「計測時間間隔」と「ファイルの保存先」を入力して，「保存」する
9．軸の設定をする
X軸：時刻に設定する．手動にて入力する
Y軸：Depth[m]を入力．適切なメモリ幅を入力．
D. 波形解析
1．「解析」-「波形処理」-緑の「→」を押す
2．波形処理データが見える．
＊計算方法で「ゼロアップクロス」もしくは「ゼロダウンクロス」が選べる
振幅が波高となる．周期もわかる．
GMT
A. 距離投影
計測したGPSデータを測線（始点：Kayangel-e0m，終点：Kayangel-e300m）に投影して始点Kayangel-e0mからの距離を算出する．
1．以下のようなCSVファイルを用意する．経度（x），緯度（y），何かの値（z）
ここでは，4点のデータ（GPS453，GPS454，GPS455，Reef edge）を投影する．
ファイル名はkaya.csv（短いファイル名がなぜか良い）
2．GMTで以下のコマンドを入力する．
project kaya.csv -C134.720772/8.083801 -E134.723466/8.084214 –Hi -Q > kayangel.csv
パラメータの意味は以下の通り
-C で投影の中心を指定する．ここではKayangel-e0mの座標を指定．
-E で投影の端を指定する．ここではKayangel-e300mの座標を指定．
-Hi は入力ファイルにヘッダ行があることを指示する．
-Q はz,p,qの値をkm単位とする.
3．kayangel.csvが出力される．
このファイルはタブ区切りなので，メモ帳などを利用して，再度Excelに入力する．
*文字列の引用符を「なし」とするとちゃんと読み込める．
4．D行の値が始点（Kayangel-e0m）からの距離である．単位はkm．必要に応じて距離mのカラムを用意する．
メモ PDFのフォント埋め込み方法
（1）Word --> PDFとして保存
（2）Adobeで保存したPDFを開く
（3）「ファイル」-「印刷」-「プロパティ」-「PDF設定」を「プレス品質」にして「OK」
（4）印刷-「OK」