R学习日志（二）-- 逻辑回归与模型选择：逻辑回归、AIC/BIC、ROC曲线 - 好文

Task3 逻辑回归

数据读取

a<-read.csv(file.choose())
edit(a)
summary(a)

数据分组

a1=a[a$ROW<640,-1]
a2=a[a$ROW>639,-1]

#根据内外样本比例对总体数据进行分组并命名，可随机打乱次序后分组

方差分析

options(digits = 9)

#调整显示小数位数

glm0.a1=glm(VIO~1,family = binomial(link=logit),data=a1)
glm1.a1=glm(VIO~SEX+MAR+CHI+REV+BAC+MAS+DOC+FOR+SIX+UPS+FRI+CON,
family = binomial(link=logit),data=a1)

anova(glm0.a1,glm1.a1)

#建立逻辑回归模型，进行方差分析（其中glm0.a1是空模型: ;

而glm1.a1是全模型）

模型整体显著性水平

1-pchisq(130.4875,12)

#130.4875为anova中的Deviance，12为自由度之差

输出回归结果

summary(glm1.a)

中心化、标准化

scale(b$估计值,center = T, scale = T)
scale(c$估计值,center = F, scale = T)
#变量系数差异太大可以把连续变量标准化中心化缩小差距

连续变量的排序柱状图（正负柱状图）

install.packages("ggpolt2",destdir = "D:/R_projects")
library(ggplot2)

x <-c("收入水平","微博好友数","消费观念")
y <-c(0.00001388,-0.0248,-1.4140) #设置x,y变量
df <- data.frame(x = x,y = y) # 判断y是否为正值
df <- transform(df,judge = ifelse(y>0,"YES","NO")) # 去除图例用theme()主题函数
ggplot(df,aes(x=reorder(x,y),y = y,fill = judge))+
geom_bar(stat = "identity",col="darkred")+
theme(legend.position= "")+
scale_fill_manual(values = c("olivedrab","grey"))+
ylab("标准化系数估计值")+xlab("连续变量")#reoder排序

Task4 模型选择与精度评价

AIC/BIC选择变量

AIC(glm1.a) #AIC
AIC(glm1.a,k=log(length(a[,1]))) #BIC

logit.aic=step(glm1.a,trace=0)
summary(logit.aic) #经过AIC选择后的变量系数输出表

n=length(a1[,1])
logit.bic=step(glm1.a,k=log(n),trace=0)
summary(logit.bic) #经过BIC选择后的变量系数输出表

（AIC、BIC都认可的变量是最为重要的变量）

预测结果及P(Y=1)的概率

pre<-predict(glm1.a,data=a)
p<-exp(pre)/(1+exp(pre))

绘制P(Y=1)的直方图

b1<-read.csv(file.choose())
breaks2<-seq(0,1,0.02)
hist(b1$p,breaks = breaks2, freq=TRUE,
include.lowest = TRUE,right = TRUE,
density =1000, col = "olivedrab", border = "white",
main="(Y=1)的概率直方图",xlab="概率区间",ylab="频率")

ROC曲线

install.packages("pROC",destdir = "D:/R_projects")
library(pROC)

#b1表中必须有0-1变量的outcome(good or bad)作为真实情况,与预测概率p做比较绘制ROC曲线

data(b1)
rocobj1<-roc(b1$outcome,b1$p)
auc(rocobj1) #得出AUC的值，越大模型越精确
plot(rocobj1,print.auc=TRUE, auc.polygon=TRUE, grid=c(0.1, 0.2),
grid.col=c("olivedrab", "grey"), max.auc.polygon=TRUE,
auc.polygon.col="olivedrab3", print.thres=TRUE,
main="ROC曲线",xlab="特异度",ylab="灵敏度") #ROC曲线越往左上角突出模型拟合越好

热门工具换一换