逻辑回归

一、逻辑回归

    逻辑斯蒂是统计学习中的经典分类算法!

二、虚拟变量

    定义：虚拟变量又称虚设变量、名义变量或哑变量，用以反映质的属性的一个人工变量,是量化了的质变量，通常取值为0或1。引入哑变量可使线形回归模型变得更复杂，但对问题描述更简明，一个方程能达到俩个方程的作用，而且接近现实。例如，反映文程度的虚拟变量可取为：1:本科学历；0：非本科学历

    一般地，在虚拟变量的设置中：基础类型、肯定类型取值为1；比较类型，否定类型取值为0。

    设置原则：

   在模型中引入多个虚拟变量时，虚拟变量的个数应按下列原则确定：

1）如果有m种互斥的属性类型，在模型中引入（m-1）个虚拟变量，否则会导致多重共线性。称作虚拟变量陷阱。　例如，性别有2个互斥的属性，引用2-1=1个虚拟变量；再如，文化程度分小学、初中、高中、大学、研究生5类，引用4个虚拟变量。

2）关于定型变量中哪个取0哪个取1是任意的，不影响检验结果。

3）若定型变量取值为0，所对应的类别称为基础类别。

4)对于多于两个类别的定型变量可采用设一个虚拟变量，而对于不同类别采取赋值不同的方法处理。

 

三、逻辑回归实例（包含R代码和Python代码）

    该数据集是共包含4个变量，其中admit为因变量，gre、gpa、prestige为自变量，且prestige有1,2,3,4四个等级，需要设置虚拟变量。

（a）R代码

data<-read.csv("C:/Users/ymz/Downloads/binary.csv",header=TRUE)
head(data)
names(data)<-c("admit","gre","gpa","prestige")#其中数据集中prestige名是后来更换，原本是rank
data$prestige<-as.factor(data$prestige)
library(nnet)
dummy<-class.ind(data$prestige)#设置虚拟变量
head(dummy)#得到虚拟变量向量
dummy<-as.data.frame(dummy)#将格式转化成数据框
names(dummy)<-c('prestige_1','prestige_2','prestige_3','prestige_4')#对数据框每一列重命名
head(dummy)
#ndata<-merge(x = data, y = dummy, by = NULL)#合并数据框（但是好像顺序打乱？数据量数目不对）
#head(ndata)
ndata<-cbind(data,dummy)#合并数据框，按列合并，数据数目对应相同
head(ndata)
newdata1<-ndata[,c(-4,-5)]#去掉因变量和一个为基准的虚拟变量
head(newdata1)
#glm.fit=glm(data$admit~data$gre+data$gpa+dummy$prestige_2+dummy$prestige_3+dummy$prestige_4,family=binomial(link="logit"))#注意，虚拟变量取N-1个进入模型，未进入模型的作基准。
glm.fit=glm(admit~.,data=newdata1,family=binomial(link="logit"))
summary(glm.fit)#查看参数估计结果
predict <- predict(glm.fit,type='response',newdata=newdata1)
real <- newdata$admit
res <- data.frame(real,predict =ifelse(predict>0.5,'1','0'))
table(real,predict =ifelse(predict>0.5,'1','0'))#计算正确率　

 代码结果显示：

由代码显示结果知：正确率为（254+30）/(254+19+97+30)=0.71.

 

  （b）python代码

import numpy as np
import pandas as pd
import statsmodels.api as sm
import pylab as pl
from pandas.core import datetools
from scipy.stats import chisquare


df=pd.read_csv("C:/Users/ymz/Downloads/binary.csv")
print(df.head())
df.columns = ["admit", "gre", "gpa", "prestige"] #对数据集中的rank重命名为prestige
#print(df.columns)
print(df.describe())#与R中summary相似,描述数据的情况，比如每一列数数目情况、均值、标准差、分位数等。
print(df.std())#输出每一列的标准差

# 频率表，表示prestige与admin的值相应的数量关系
print(pd.crosstab(df['admit'], df['prestige'], rownames=['admit']))
# plot all of the columns
# df.hist()
# pl.show()
#虚拟变量，也叫哑变量，可用来表示分类变量、非数量因素可能产生的影响。在计量经济学模型，需要经常考虑属性因素的影响。例如，职业、文化程度、季节等属性因素往往很难直接度量它们的大小。只能给出它们的“Yes—D=1”或”No—D=0”，或者它们的程度或等级。为了反映属性因素和提高模型的精度，必须将属性因素“量化”。通过构造0-1型的人工变量来量化属性因素。

#pandas提供了一系列分类变量的控制。我们可以用get_dummies来将”prestige”一列虚拟化。

#get_dummies为每个指定的列创建了新的带二分类预测变量的DataFrame，在本例中，prestige有四个级别：1，2，3以及4（1代表最有声望），prestige作为分类变量更加合适。当调用get_dummies时，会产生四列的dataframe，每一列表示四个级别中的一个。
dummy_ranks = pd.get_dummies(df['prestige'], prefix='prestige')#创建虚拟变量（哑变量）
print (dummy_ranks.head())
# 为逻辑回归创建所需的data frame
# 除admit、gre、gpa外，加入了上面常见的虚拟变量（注意，引入的虚拟变量列数应为虚拟变量总列数减1，减去的1列作为基准）
cols_to_keep = ['admit', 'gre', 'gpa']
data = df[cols_to_keep].join(dummy_ranks.ix[:, 'prestige_2':])#将自变量数据整合到一起
print(data.head())
data['intercept']=1.0#注意，python中需要加入截距项常数项，而R中不需要
#这样，数据原本的prestige属性就被prestige_x代替了，例如原本的数值为2，则prestige_2为1，prestige_1、prestige_3、prestige_4都为0。

#将新的虚拟变量加入到了原始的数据集中后，就不再需要原来的prestige列了。在此要强调一点，生成m个虚拟变量后，只要引入m-1个虚拟变量到数据集中，未引入的一个是作为基准对比的。

#最后，还需加上常数intercept，statemodels实现的逻辑回归需要显式指定

#实际上完成逻辑回归是相当简单的，首先指定要预测变量的列，接着指定模型用于做预测的列，剩下的就由算法包去完成了。

#本例中要预测的是admin列，使用到gre、gpa和虚拟变量prestige_2、prestige_3、prestige_4。prestige_1作为基准，所以排除掉，以防止多元共线性(multicollinearity)和引入分类变量的所有虚拟变量值所导致的陷阱(dummy variable trap)。

# 指定作为训练变量的列，不含目标列`admit`
train_cols = data.columns[1:]#第一列到最后一列，去掉第零列
# Index([gre, gpa, prestige_2, prestige_3, prestige_4], dtype=object)

logit = sm.Logit(data['admit'], data[train_cols])#建立模型

# 拟合模型
result = logit.fit()

import copy

combos = copy.deepcopy(data)

# 数据中的列要跟预测时用到的列一致
predict_cols = combos.columns[1:]

# 预测集也要添加intercept变量
combos['intercept'] = 1.0
# 进行预测，并将预测评分存入 predict 列中
combos['predict'] = result.predict(combos[predict_cols])
# 预测完成后，predict 的值是介于 [0, 1] 间的概率值
# 我们可以根据需要，提取预测结果
# 例如，假定 predict > 0.5，则表示会被录取
# 在这边我们检验一下上述选取结果的精确度
total = 0
hit = 0
for value in combos.values:
    # 预测分数 predict, 是数据中的最后一列
    predict = value[-1]
    # 实际录取结果
    admit = int(value[0])

    # 假定预测概率大于0.5则表示预测被录取
    if predict > 0.5:
        total += 1
        # 表示预测命中为1
        if admit == 1:
            hit += 1



# 输出结果
print('Total: %d, Hit: %d, Precision: %.2f' % (total, hit, 100.0 * hit / total))
# Total: 49, Hit: 30, Precision: 61.22
print(result.summary2())#查看逻辑回归详细结果

#在最后进行模型解释时，所有类别哑变量的回归系数，均表示该哑变量与参照相比之后对因变量的影响。
print(result.conf_int())#查看系数的置信区间
print(np.exp(result.params))

这计算的是Precsion ，预测为1中正确的数和预测为1的数的比例。