#第1步：讀入資料
setwd("C:/Users/IBM/Desktop/logsistic建模/3.2建模")#設定當前的工作目錄，重要！
safe<-read.csv("model2.csv",header=T)
str(safe) #看錶頭
colnames(safe)#看列名
attach(safe)
cor(safe)   #變數間相關係數陣

#第2步：使用glm函式計算.
glm.safe<-glm(is_run~reg_age+reg_life+heart_num+max_period
+min_period+median25_period+median50_period+median75_period
+lost_num_m30+avg_income_fund+avg_amt_fund+avg_nums_fund
+avg_amt_p2p+avg_nums_p2p+log_num+avg_income_p2p
,data=safe,family=binomial)#logistic模型
summary(glm.safe)

#第3步：逐步迴歸
glm.safe1<-step(glm.safe)
summary(glm.safe1)#逐步迴歸

#第4步：模型解讀
coef(glm.safe1)#用以解釋模型
exp(glm.safe1$coefficients) #解釋Odds比與x的關係
exp(confint(glm.safe1)) #迴歸係數的置信區間
xp05<-0/glm.safe1$coefficients[]
xp05 #求使得pi為0.5的x
ratio05<-glm.safe1$coefficients[]*0.25
ratio05 #pi為0.5處的pi關於x的變化率

#第5步：模型評價
#5.1計算決定係數
R2cox<-1-exp((glm.safe1$deviance-glm.safe1$null.deviance)/length(is_run))

#5.2計算Cox-Snell擬合優度
cat("Cox-Snell R2=",R2cox,"\n")
R2nag<-R2cox/(1-exp((-glm.safe1$null.deviance)/length(is_run)))

#5.3計算Nagelkerke擬合優度
cat("Nagelkerke R2=",R2nag,"\n")

#5.4殘差分析
plot(residuals(glm.safe1))

#5.5異常值診斷
library(car)
influencePlot(glm.safe1)

#第6步：分類表
fitt.pi<-fitted(glm.safe1)#同predict(glm.safe1,data.frame(x2=x2),type="resp")
ypred<-1*(fitt.pi>0.5) #1*邏輯變數就變成了0和1變數
ypred
length(ypred)
n<-table(safe$is_run,ypred)
n
Percantage<-c(n[1,1]/sum(n[1,]),n[2,2]/sum(n[2,]))
Percantage
rbind(Category,Percantage)