10-5 统计模型如何纳入以上分类
统计学模型是机器学习模型的另一种说法,这些模型高度依赖概率论和统计学的数学原理,从变量中寻找关系。
10-6 线性回归
预测因子和响应变量之间的关系,用数学公式表示:
y = \beta_{0} + \beta_{1}x_{1}+\beta_{2}x_{2}+...+\beta_{n}x_{n} \\ y是响应变量\\ x_{i}是第i个变量\\ \beta_{0}是截距\\ \beta_{i}是第i个变量的相关系数
【例子】根据每小时共享单车的使用情况,预测每小时共享单车租用的数量。
第一步 : 先做 温度 和 租用数量 的散点图 , 然后 绘制最佳拟合线。
# 画散点
df.plot(kind = 'scatter', x = '温度' , y = '租用数量' ,alpha = 0.2 )
# 作拟合
import seaborn as sns
sns.lmplot(x = '温度' , y = '租用数量',data = df ,aspect = 1.5 ,scatter_kws = {'alpha':0.2})
plt.show()
图像显示两者成正比,表面上看共享单车租用量在随气温的升高而增加。
第二步:用相关系数量化变量间的关系,验证它们是否和表面上看起来的规律一致。
bikes[['数量','温度']].corr()
# 得到 0.3944
可见,两个变量有微弱的正相关关系。回到线性回归公式:
y = \beta_{0} + \beta_{1}x_{1}+\beta_{2}x_{2}+...+\beta_{n}x_{n} \\
因为点很多,而线是直线,点全部在线上的可能性很小,所以会退而求其次,让尽可能多的点落在上面。
线性回归模型中,只要给定x和y风速概率分布,模型就可以找出最佳的贝塔系数(beta coefficients) 也叫模型系数(model coefficients).
每个点和直线上同x对应的y轴距离,叫残差(residual)
残差做平方然后加和,得到残差平方和(sum of squared residuals) ,公式:
SS_{residuals} = \sum_{i=1}^{N}{(y_{i 预测} -y_{i观测} )}^{2} \\
残差平方和最小的拟合线就是最佳拟合线。
features_cols = ['温度']
X = bikes[features_cols] # 预测因子(自变量)
y = ['租用数量'] # 响应变量(因变量)
from sklearn.linear_model import LinearRegression
linreg = LinearRegression()
linreg.fit(X,y)
print linreg.intercept_ # 得到截距 6.0462 指温度为0时的租车数量。
print linreg.coef_ # 得到变量的相关系数 这里就是温度系数 [9.1705]
注意:不是所有的数据等于0时都有意义。
温度系数 是因变量除以自变量得到的 ;表明 自变量和因变量如何同时变化;9.1705说明温度每升高1摄氏度,共享单车的租用量就会增加9量。 ; 符号很重要,如果是符号,说明租用量将随着温度升高而下降。
相关性并不等于因果关系,并不能得出“温度升高导致了自行车租用量增加”的结论。
#做预测
linreg.predict(20) # 189.4570
# 当温度为20摄氏度时,可能有190辆共享单车被租用。
10-6-1 增加更多预测因子
增加预测因子只要修改一下参数,在这之前,研究一下数据集,搞清楚每个预测因子的含义。
季节:春夏秋冬分别表示为 1 2 3 4
假期:当天是否是假期
工作日:工作日或者周末
天气:晴天、薄雾、小雪、大雨。。
温度、体感温度、空气湿度、风速、租客中未注册人数,租客中已注册人数,租用量。
feature_cols = ['温度' ,'季节','天气','空气湿度']
X = bikes[feature_cols]
y = bikes['count']
linreg = LinearRegression()
linreg.fit(X,y)
zip(feature_cols , linreg.coef_)
# 可以看到 各变量的系数
# 每一组都意味着 保持其他预测因子不变 变量增加1个单位,共享单车量增加 linreg.coef_ 辆
温度: 7.89
季节:22.53
天气:6.6
空气湿度:-3.11
可以看到 天气增加(逐渐阴天)、季节增加(逐渐冬天),单车租量也在增加,这不太符合预期。
#做散点图看一下
feature_cols = ['温度' ,'季节','天气','空气湿度']
sns.pairplot(bikes , x_vars = feature_cols , y_vars = '出租数量' ,kind = 'reg')
plt.show()
打完散点之后发现天气的趋势是下行的,和模型得到的结果相反。那就是有些预测因子可能是没用的噪声,并不能帮助预测。怎么才能知道谁有没有用,要使用一些更高级的指标。
10-6-2 回归指标
量化线性回归模型有效性,有三种主要的指标:
平均绝对误差(MAE):误差绝对值之和的平均值。 \frac{1}{n}\sum^{n}_{i=1}\left| y_{i} - y_{ihat} \right|
均方误差(MSE):误差平方之和的平均值。 \frac{1}{n}\sum^{n}_{i=1}( y_{i} - y_{ihat})^{2}
均方根误差(RMSE):误差平方之和的平均值的平方根。 \sqrt{\frac{1}{n}\sum^{n}_{i=1}( y_{i} - y_{ihat})^{2}}
n 是观测数量 \\ y_{i}是实际值\\ y_{ihat}是预测值
true = [9,6,7,6]
pred = [8,7,7,12]
from sklearn import metrics
import numpy as np
print 'MAE:' , metrics.mean_absolute_error(true,pred)
print 'MSE:' , metrics.mean_squared_error(true,pred)
print 'RMSE:' , np.sqrt(metrics.mean_absolute_error(true,pred))
# 输出:
MAE:2.0
MSE:9.5
RMSE:3.0822
MAE是绝对值,容易理解 ;MSE是平方和 更合理;RMSE是MSE开方 容易解释和说明 最常用。
对回归模型,推荐使用RMSE,这些指标都属于损失函数(loss function) , 因此 越小越好。
第三步:用均方根误差(RMSE)判断哪些预测因子有助于预测,哪些干扰了预测。
流程: 生成变量X和y ——对线性回归模型进行拟合——基于X,使用模型生成一列预测值——根据预测值和实际值,计算模型的均方根误差(RMSE)
from sklearn import metrics
feature_cols = ['温度'] # 只有温度
X = bikes[feature_cols] # 生成变量
linreg = LinearRegression()
linreg.fit(X,y) # 拟合
y_pred = linreg.predict(X)# 预测
np.sqrt(metrics.mean_squared_error(X,y_pred)) # 求RMSE
# 得到 166.45
from sklearn import metrics
feature_cols = ['温度','空气湿度'] # 温度基础上增加湿度
X = bikes[feature_cols] # 生成变量
linreg = LinearRegression()
linreg.fit(X,y) # 拟合
y_pred = linreg.predict(X)# 预测
np.sqrt(metrics.mean_squared_error(X,y_pred)) # 求RMSE
# 得到 157.79
同样的套路,多加几个变量进行预测
发现RMSE又减小了,是我们想看到的结果。
但这样做可能导致过拟合(overfitting),解决办法是把数据集拆成训练集和测试集。
步骤:分拆数据集-训练模型-测试模型-参数调优-选择最佳模型-对所有数据进行训练-使用新数据预测
1 把数据集拆成 训练集(trainging set ) 和 测试集(test set)
2 用训练集训练,然后测试集测试模型。
3 达到最佳状态(满足模型评价指标),将模型推广到整个数据集
4 模型做好准备迎接新数据
做上边这些工作的目的是让 模型预测陌生数据时,误差最小。
--暂停 2021-4-6 13:08
-- 2021-4-6 15:19
# 代码实现
from sklearn.cross_validation import train_test_split
feature_cols = ['温度']
X = bikes[feature_cols]
y = bikes['租车数量']
X_train , X_test , y_train , y_test = train_test_split(X,y)
linreg = LinearRegression()
linreg.fit(X_train,y_train)
y_pred = linreg.predict(X_test)
np.sqrt(metrics.mean_squared_error(y_test ,y_pred))
# 输出: 166.91
# 代码实现
from sklearn.cross_validation import train_test_split
feature_cols = ['温度','工作日'] # 增加
X = bikes[feature_cols]
y = bikes['租车数量']
X_train , X_test , y_train , y_test = train_test_split(X,y)
linreg = LinearRegression()
linreg.fit(X_train,y_train)
y_pred = linreg.predict(X_test)
np.sqrt(metrics.mean_squared_error(y_test ,y_pred))
#输出: 166.95
虽然得到了RMSE,但是这个数字,没有对比,也就没办法衡量它是不是好。
一个评价方法是和空模型(null model)对比 。
空模型就是什么参数都不用,直接瞎蒙,回归模型的空模型设置方法是 用每小时平均租车量作为瞎猜的结果,看它的RMSE。
average_bike_rental = bikes['租车数量'].mean()
average_bike_rental
# 输出均值 : 191.57
num_rows = bikes.shape[0]
null_model_predictions = [average_bike_rental]*num_rows
np.sqrt(metrics.mean_squared_error(y,null_model_predictions))
# 输出RMSE : 181.13613
可以看到,空模型的RMSE是大于之前构造的模型的,说明模型有效,然而击败空模型是机器学习的最低要求。
到这,线性回归就告一段落,引出Logistic 回归风速概率分布,
Logistic 和 线性回归核心思想很相似,但它用来做分类。
10-7 Logistic 回归 (logistic regression)
既然是分类算法,为什么要叫回归模型?
因为它是在用回归方法处理分类问题。
使用线性回归模型时,用以下函数拟合:
y = \beta_{0} + \beta_{1}x \\ y是响应变量\\ \beta是模型参数\\ x是输入变量 (可以有多个输入变量)\\
使用Logistic分类“类别1” 时 ,用以下函数 :
P(y=1|x) = \frac{e^{\beta_{0}+\beta_{1}x}}{1+e^{\beta_{0}+\beta_{1}x}} \\ P(y=1|x)表示给定数据x的前提下,响应变量属于“类别1”的条件概率。\\ e是欧拉常数(Euler's number) 约等于2.718
函数右边部分叫做 Logistic函数(Logistic function)
欧拉常数经常用来对自然增长和衰减进行建模。
为什么要用这个,而不用线性回归方程,一个是因为希望模型返回的范围能在0-1之间,就像真实的概率分布。 另一个是因为Logistic 更平滑,更符合分类问题的特质。
10-8 概率、几率和对数几率
3000人中有1000人购物:
任意一人购物的概率:P(购物) = 1000/3000
任意一人购物的几率:Odds(购物) = 1000/2000
概率和几率的转换 : Odds = P / 1-P
table = pd.DataFrame({'probability':[0.1,0.2,0.25,0.5,0.6,0.8,0.9]})
table['odds'] = table.probability /(1 - table.probability)
table
概率增加,几率更快增加
从转换公式也可以看出来, 概率接近1时,几率接近无穷大。但是永远大于0, 所以这样就也没法用。
接着想到自然数和对数。
table = pd.DataFrame({'probability':[0.1,0.2,0.25,0.5,0.6,0.8,0.9]})
table['odds'] = table.probability /(1 - table.probability)
table ['log_odds'] = np.log(table.odds)
table
对数几率(Log-odds ,也称为Logit )
没有上下界并且有正负,很理想,这也正是Logistic回归的起点。
Logistic 回归的数学原理
假设p表示数据点属于特定类别的概率,那么Logistic回归的公式可表示为:
log_{e}(\frac{p}{1-p}) = \beta_{0} + \beta_{1}x \\
对这个公式求解,就得到值介于0-1 之间的Logistic 曲线,转化一下:
p = \frac{e^{\beta_{0}+\beta_{1}x}}{1+e^{\beta_{0}+\beta_{1}x}}\\
要更好理解logistic,就要理解概率和几率的区别。
线性回归中, \beta_{1} 表示x发生1单位变化,响应变量的变化值。Logistic回归中, e^{\beta_{0}} 表示x发生1单位变化,几率的变化值。
假设 Logistic 得到参数 \beta_{1} = 0.693 ,可以计算几率 e^{\beta_{1}} =2 , 说明样本属于该分类的几率是其他分类的两倍 。
除了二元分类,Logistic 还可以做多分类问题,用1对多方法。
回到共享单车例子:
# 先用空模型预测 分类问题通常将最常出现的类别作为预测结果
bikes['租车数量是否大于均值'] = bikes['租车数量'] >= average_bike_rental
bikes['租车数量是否大于均值'].value_counts(normalize = True)
#输出 : False 0.599853 True 0.400147 那就是全都预测低于,这样正确率就是将近60%
# 用Logistic回归模型预测指定时间的共享单车租用量是否高于平均值
from sklearn.linear_model import LogisticRegression
feature_cols = ['温度']
X = bikes[feature_cols]
y = bikes['租车数量是否大于均值']
X_train,X_test ,Y_train,Y_test = train_test_split(X,y)
logreg = LogisticRegression()
logreg.fit(X_train, y_train)
logreg.score(X.test,y_test)
# 0.6565
到这发现做预测的时候只能使用定量特征,如果确知分类特征和响应变量间有关联,怎么处理? —— 哑变量
10.9 哑变量(dummy variables)
用于将分类特征转换成定量特征。
哑变量指定性数据被重新编码后的新数据。 把定性数据的每个值在转换成1 和 0 组成的新数据。
# 测试功能 x[11]+x[12]是在取小时,然后int 类型转换
datetime = '1234-56-78 ab:cd:ef'
datetime[11]+datetime[12]
# 得到小时 放到一列
bikes['hour'] = bikes['datetime'].apply[lambda x:int(x[11]+x[12])]
# 划分时间段
def when_is_it(hour) :
if hour >= 5 and hour < 11 :
return 'morning'
elif hour >=11 and hour <16 :
return 'afternoon'
elif hour >=16 and hour <18 :
return 'rush_hour'
else:
return 'off_hours'
# 分类放到一列
bikes['when_is_it'] = bikes['hour'].apply(when_is_it)
bikes[['when_is_it','租车数量是否大于均值']].head()
bikes.groupby('when_is_it').租车数量是否大于均值.plot(kind = 'bar')
得到条形图,发现租车量受时间段影响。
# 生成哑变量
# 第一种方法
when_dummies = pd.get_dummies(bikes['when_is_it'],prefix x='when_')
# 第二种方法
when_dummies = when_dummies.iloc[:,1:] # 不要第一列
生成哑变量之后,分类特征转换成数值,使用Logistic回归模型:
X = when_dummies
y =bikes.租车数量是否大于均值
logreg = LogisticRegression()
logreg.fit(X_train,y_train)
logreg.score(X_test,y_test)
# 0.6851 比单独用温度特征评分高
汇总特征,再做一次
new_bike = pd.concat([bikes[['温度','空气湿度']],when_dummies],axis =1 )
X = new_bike
y =bikes.租车数量是否大于均值
X_train, X_test ,Y_train,Y_test = train_test_split(X,y)
logreg = LogisticRegression()
logreg.fit(X_train,y_train
logreg.score(X_test,y_test)
# 0.7182 评分提升了
10-10 总结
更重要的是我们将学习把数据科学应用到显示世界的新方法。
发表评论