とりあえずグラフでも書いてみる

　空間統計は、あまりにもディープな世界であることが判明したため、見なかったことに。
  む～奥が深いですな。

　とりあえず、追加的な解析。まずはランダム効果入れて乱数計算。正答率がどれくらいか見てみます。
　非常にR的なスクリプト。ベクトルがさも、一つの数値かのように扱えるのは便利極まりないです。

##glmmの結果の検証
##ランダム要因を入れて正答率計算
####################################

setwd("d:/【研究】葉枯れ")
dat<-read.csv("gis_output.csv",header=T)
v.output <- NULL
##dat$KARE.binは1,0の従属変数
base <- dat$KARE.bin

for(i in 1:100){
rel <- 1/(1+exp(-(0.1165*dat$DBH-5.327+rnorm(1,0,0.000024))))
diff <- base-rel
no.diff <- diff[0.5>diff&-0.5<diff]
output <- length(no.diff)/length(base)*100
v.output <- c(v.output,output)
}

次は、ラグプロット入れたグラフづくり。　ラグプロット入ると「お、Rだ」的な雰囲気になりますな。 後、最近ようやくESSの使い方もはまってきたんですが、もう無茶便利。dabbrevのabbrevでサクサクっと字書いて、ちゃっと実行。Rのコンソールで書く気起きなくなります。まだ、使われてない方は是非。

for(i in 1:100){
 par(new=T,mai=c(1,1,1,1))

  plot(dat$DBH,dat$H,ylim=c(0,1),xlim=c(0,80),cex.axis=2,cex.lab=2,lty=10,xlab="胸高直径(cm)",ylab="確率",type="n")
  line <- function(x) {1/(1+exp(-(0.1165*x-5.327+rnorm(1,0,0.000024))))}
  x <- seq(0,80,length=100)
  lines(x,line(x),col="yellow",lwd=3)
 }

##0,1のラグプロットを書いてみる
##jitter()は重なるバーをずらして表示する関数
abs <- dat[dat$KARE.bin==0,]
par(new=T,mai=c(1,1,1,1))
rug(jitter(abs$DBH, amount = .01),lwd=1)

pre <- dat[dat$KARE.bin==1,]
par(new=T,mai=c(1,1,1,1))
rug(jitter(pre$DBH, amount = .01),lwd=1,side=3)

　出来たグラフはこんな感じ。

　ま～今回の全てはこのグラフに集約されるわけで、なんていうか・・・・ショボイ。　ランダム効果の意味ほとんどないし。
　つーことはあれかい、glm()でも別に良いってこと？いいんすglmmやってみたかっただけなんで。