# ロジスティック回帰分析をする方法

beta0 <- -3	# 切片を設定
beta1 <- 0.6	# 係数を設定
grad <- seq(0, 10, 0.1)	# 環境傾度として1から10まで0.1区切りのものを設定
fact <- 1 / (1 + exp(-beta0 - beta1 * grad))	# 理論上の出現確率を計算
pa <- rbinom(length(fact), 1, prob=fact)	# 出現確率をもとに出現の有無を乱数で発生
logit <- glm(pa ~ 1 + grad, family = binomial(logit))	# 種の有無をもとにロジスティック回帰
curve( 1 / ( 1 + exp (-logit[[1]][1] - logit[[1]][2] * x)), xlim=c(0,10), ylim=c(0,1))	# 種の有無から推定した出現確率
rug(side=3, jitter(grad[pa==1]))	# グラフの上部に出現地点を表示
rug(side=1, jitter(grad[pa==0]))	# グラフの下部に非出現地点を表示 
par(new=T)	# グラフを上書きに設定
curve( 1 / ( 1 + exp (3 - 0.6 * x)), xlim=c(0,10), ylim=c(0,1), col=2)	# 理論上の出現確率を表示
par(new=F)	# 上書きの設定を解除