判别分析——逻辑斯蒂回归

对股票市场的数据进行分析，预测数据的走向

#股票市场数据
library(ISLR)

#ISLR库的股票市场数据的数值
names(Smarket)

#汇总数据
summary(Smarket)

#判断相关性散点图
pairs(Smarket)

#计算所有定量预测变量两两之间相关系数的矩阵
cor(Smarket[,-9])

#Volume与Year强相关
attach(Smarket)
plot(Volume)

接下来预测Direction(up or down)

#逻辑斯蒂回归,属于二项分布族的类别
glm.fit <- glm(Direction~Lag1+Lag2+Lag3+Lag4+Lag5+Volume,data=Smarket,family=binomial)
summary(glm.fit)
#p值大于0.05，所以没有充分的理由表明Direction和Lag1之间有确切的联系

#提取模型系数
coef(glm.fit)

#提取模型系数及其推断情况
summary(glm.fit)$coef

#提取p值
summary(glm.fit)$coef[,4]

#predict(type="response")预测市场上涨的概率P(Y=1|X),即对数发生比.可看帮助文档。
如果没有提供数据集给predict函数，默认使用训练数据的拟合逻辑斯蒂回归模型的概率
glm.probs <- predict(glm.fit,type="response")
glm.probs[1:4] #这里给出前四个预测概率

#创建一个哑变量,factor
contrasts(Direction)

#根据预测市场上涨概率比0.5大还是小，判断预测的结果
glm.pre<- rep("Down",1250)
glm.pre[glm.probs>.5] <- "Up"

#产生混淆矩阵来判断有多少观测被正确或错误分类
table(glm.pre,Direction)
        Direction
glm.pre Down  Up
  Down  145 141
  Up    457 507
  
#计算正确预测的比例
正确预测的比例：
  (145+507)/1250
mean(glm.pre==Direction)
[1] 0.5216  

总结：以上数据为用训练拟合模型的数据来做预测，测试错误率为1-52.16%=47.8%，
而训练数据集偏向于低估测试错误率。

接下来改变预测数据集，用2001-2004年的数据做训练，2005年数据做预测

#测试数据集
test <- Smarket[Year>2004,]
dim(test)
Direction.test <- test$Direction

#根据subset参数对2005年之前的数据拟合一个逻辑斯蒂回归模型
glm.fit <- glm(Direction~Lag1+Lag2+Lag3+Lag4+Lag5+Volume,data=Smarket,family = binomial,subset=(Year<2005))
summary(glm.fit1)
#subset是一个逻辑值，用于挑选观测数据集的子集，即训练数据集
#得到测试集2005年每一天的股票上涨预测概率

glm.probs <- predict(glm.fit,type="response",test) 
table(glm.probs,Direction.test)

#glm.probs是一个数值向量，将它转化成因子类型，可以用列联表统计
glm.pred <- rep("Down",252)
glm.pred[glm.probs > .5] <- "Up"
table(glm.pred,Direction.test)
      Direction.test
glm.pred Down Up
  Down   77 97
  Up     34 44

#预测正确率
  mean(glm.pred==Direction.test)
  [1] 0.4801587

总结：测试错误率更高了，而在实际生活中，人们一般很少用前几日的投资回报率来预测未来市场的表现。

优化模型：由逻辑斯蒂回归模型的p值，Lag1和Lag2两个预测变量的p值稍小，所以去除与响应变量无关的预测变量，
只用Lag1和Lag2两个预测变量。

glm.fit <- glm(Direction~Lag1+Lag2,data=Smarket,family=binomial,subset = (Year<2005))
  summary(glm.fit)
  
  glm.probs <- predict(glm.fit,test,type="response")
  glm.pred <- rep("Down",252)
  glm.pred[glm.probs > 0.5] <- "Up"
  table(glm.pred,Direction.test)
  
         Direction.test
  glm.pred Down  Up
    Down   35  35
    Up     76 106
    
    #预测正确率
    mean(glm.pred==Direction.test)
    
      106/(106+76)#预测市场上涨的正确率58%
      35/(35+35)#预测市场下跌的正确率50%