2020年9月3日木曜日

野生生物の生息適地と分布モデリング

午前
・標本整理

午後
・論文書き
・授業準備


野生生物の生息適地と分布モデリング: Rプログラムによる実践

外来種の侵入予測や希少種の生息適地の推定、さらには分布制限要因の解析など、GISを用いた分布推定は様々な場面で利用されている。

そして、地図上に出力されたカッコいい推定地図を見ていると、作ってみたくなるのが性だろう、、、そんな貴方に、必携の一冊。

ずっと購入しなければと思っていたが、実際に手に入れると、その文量と内容の豊富さでナカナカ読み終わらなかった。

内容は多岐にわたり、分布推定の基礎知識となるニッチ概念の解説から、GISデータの取り扱い方法、空間統計やモデリングの解説、そして、様々な推定方法 (距離ベース、回帰ベース、分類ベース、最大エントロピーetc) による実践例が紹介されている。

また、モデリングの際にエンドユーザーを悩ませる、モデルやデータの評価方法についても解説がある。

そして、最も重要なこととして、タイトルにもある通り、これらの解析は全てRで行われており、そのコードが書かれている!

私は以前、Rを用いて日本における外来性ワラジムシ類の分布推定をしたことがあるが、勉強を始めると様々な解析方法が開発されていることが分かり、何がどう違くて、どれを使えば良いのか分からず非常に混乱した。

そのような悩みはこの1冊で解消される!

内容は非常に充実しており、、、充実し過ぎており、一度目を通しただけで完全に理解することは、私には不可能であった。実際に解析する際に辞書的に使うのが現実的か。

内容が豊富なだけに、Rも分布推定も全く手を出したことがない人が、どれくらい理解できるかやや心配だが(大学生向けらしいが)、Rコードを実行しながら読み進めることで、実体験を通して理解できるように思われる。

とりあえず、やってみたい人は、19章に解析事例があるので、これを試してみると良いかも。

また、第I部の第3章は、分布推定とは関係なく、生物の分布規定要因を勉強するのに非常に良い資料だと思う。

ちなみに、分布推定・生息適地推定のモデリングでは、生態的ニッチモデル(ENM)、種分布モデル(SDM)、生息地分布モデル(HDM)なども、様々な名称がある。これらは違いがあり、厳密には使い分けられるそうだが、この本では、生息適地モデル(HSM)に統一されている。

もくじ
1.本書の概要
[I 生息適地モデリングのが概要、理論、前提]
2.生息適地モデリングの手順
3.何が種分布を決めるのか?
4.ニッチをモデル化する—概念と実際のデータ
5.生息適地モデルの仮定
[II データの編集と収集、サンプリングデザイン、空間スケール]
6.生息適地を予測する環境データの選択と編集
7.モデリング対象種のデータ収集と設計
8.生態学的スケーリング:「空間」「時間」「主題」の解像度と規模
[III モデリングのアプローチとモデル校正]
9.エンベロープアプローチと距離ベースアプローチ
10.回帰ベースアプローチ
11.分類アプローチと機械学習システム
12.ブースティングとバギングによるアプローチ
13.最大エントロピー
14.アンサンブルモデリングとモデル平均化
[IV モデルの評価:過誤と不確実性]
15.モデル精度の測定—どの指標を用いるべきか?
16.モデル性能の評価—どのデータを用いるべきか?
[V 空間的な予測と時間的な予測]
17.時間・空間上へのモデルの投影
[VI 本書で用いたデータとツールおよび発展的な事例]
18.本書の分析で用いたデータとツール
19.Biomod2モデリングパッケージを用いた分析事例
[VII 将来展望] 20.生息適地モデルの将来

このデータでやってみるか。

2020年9月2日水曜日

思い出す

午前
・標本整理

午後
・論文書き
・授業準備

久しぶりにサソリモドキ 論文に手を出す。間を空けすぎたので、思い出すだけでも結構苦戦。あと少しなのだが、そのあと少しが終わらない。

2020年8月31日月曜日

データが無い期間

午前
・標本整理

午後
・修論手伝い(データ解析)
・学会業務

修論のデータ解析。データが無い期間を0個として表現する方法でつまづいた、、、未だに、つまづいている。

これといった用事はなかったのに、午後は時折届くメールへの対応でバタバタした。

2020年8月28日金曜日

お返しします

午前
・データ解析

午後
・データ解析

消耗品が切れてしまったので、一日中データ解析。数ヶ月待たせてしまったが、共同研究者にやっとdraft reportができた。

ここ数日苦戦したのが、haplotype networkの描画。結局、解決はしていないけど。

library(pegas) を使った方法だと、プロットのサイズの制御ができなかった。1枚だけなら良いのだが、今回は6つのサブグループがあり、グループ間でサイズの統一ができず断念。

library(haplotypes) では、プロットの配置の際に、シミュレーション?が行われているようで、出力するたびに微妙に配置が異なる。

完全な制御はできていないが、
coo1 <- plot(p, net=1, mode="fruchtermanreingold",,,,)
plot(p, net=1, coord=coo1, label.pos=5,,,,)
とすると、大きくズレることはない気がする?気のせいかも。

とはいえ、どうにかこうにか、haplotipe network図を作成し、
系統樹を作成し、
地図に投影できた。ゴチャゴチャ、、、でも、しばらくは考えたくない。一旦、お返しします。

2020年8月27日木曜日

どハマり

午前
・データ解析
・統計ゼミ

午後
・個人ゼミ
・データ解析

昨日で解決したと思ったHaplotype netwokの描画で、どハマり中、、、。

2020年8月26日水曜日

ハプロタイプ;シーケンス名変更;「.fasta」を出力

午前
・標本整理

午後
・授業準備
・データ解析

かなり苦戦したが、どうにかできるようになった。

library(haplotypes)を用いてハプロタイプごとにシーケンスをまとめ、さらに、シーケンス名を変更したのち、「.fasta」を出力する方法。
#library(haplotypes) #「::」で指定する
#library(ape) #「::」で指定する
x <- haplotypes::read.fas("???.fasta")
h <- haplotypes::haplotype(x) #to class"haplotype"
#h@nhap #ハプロタイプ数の確認
h_dna <- haplotypes::as.dna(h) #to class"DNA"
names(h_dna) <- c(1:34) #シーケンス名の変更;class"DNA"で使う
h_dnabin <- haplotypes::as.DNAbin(h_dna) #to class"DNAbin"
ape::write.dna(h_dnabin, file = 'test.fasta', format = 'fasta')
library()は使わず、「::」で指定した方が良い?詳しく検証していないが、2つのlibraryで同じ関数名が使われている?

関数が「class」に依存するので、順番を間違えるとハマる。「class:haplotype」→「class:DNA」→「class:DNAbin」で「ape::write.dna」が実行できる。

また、「names()」は「class:DNA」で使う。

2020年8月25日火曜日

packages"haplotypes"

午前
・標本整理

午後
・データ解析
・手伝い
・会議

変異が大きいため、pegas::haplotype、pages::haploNetを使った方法では、ハプロタイプネットワークがうまく描けなかった。
replot()を使って、頑張ればどうにかなりそうだけど、あまり現実できではない。

で、昨夜、packages"haplotypes"というのを見つけた。

ハプロタイプネットワークの老舗ソフトTCSでも利用されている最節約法によるサブネットワーク化?を計算することができる。pagesでもできるかも?
install.packages("haplotypes")
library(haplotypes)
#マニュアルでは「.fas」となっているが、うまくいかなかった。
x <- read.fas("***.fasta")
h <- haplotypes::haplotype(x)
p <- parsimnet(x, prob=.95)
rownames(p@d[[1]])<-c(1:8)
plot(p,
   net=1, #network番号の指定
   label.pos=5,
   vertex.cex=c(rep(5, p@nhap[[1]]), #円の大きさの指定
		rep(1, c(nrow(p@d[[1]])-p@nhap[[1]]))),
   vertex.col=c(rep("magenta", p@nhap[[1]]), #の色の指定
   		rep("deepskyblue3", c(nrow(p@d[[1]])-p@nhap[[1]]))))
パッケージのマニュアルでは、「p@nhap」となっているが、「p@nhap[[1]]」と「[[1]]」の指定が必要。この数字はネットワーク番号に合わせて変更する。

例えば、2つ目のネットワークの緑色で出力する場合は下記となる。
rownames(p@d[[2]])<-c(7:11)
plot(p,
   net=2,
   label.pos=5,
   vertex.cex=c(rep(5, p@nhap[[2]]),
		rep(1, c(nrow(p@d[[2]])-p@nhap[[2]]))),
   vertex.col=c(rep("green", p@nhap[[2]]),
   		rep("deepskyblue3", c(nrow(p@d[[2]])-p@nhap[[2]]))))