关灯

数学推导纯Python实现机器学习算法:线性回归

[复制链接]
admin 发表于 2019-1-20 13:11:44 | 显示全部楼层 |阅读模式 打印 上一主题 下一主题
 

很多同学在学习机器学习的时候,理论粗略看一遍之后就直接上手编程了,非常值得表扬。但是他不是真正的上手写算法,而是去直接调用 sklearn 这样的 package,这就不大妥当了。笔者不是说调包不好,在实际工作和研究中,封装好的简单易用的 package 给我们的工作带来了莫大的便利,大大提高了我们机器学习模型和算法的实现效率。但这仅限于使用过程中。

 

笔者相信很多有企图心的同学肯定不满足于仅仅去使用这些 package 而不知模型和算法的细节。所以,如果你是一名机器学习算法的学习者,在学习过程中最好不要一上来就使用这些封装好的包,而是根据自己对算法的理解,在手推过模型和算法的数学公式后,仅依靠 numpy 和 pandas 等基础包的情况下手写机器学习算法。如此一遍过程之后,再去学习如何调用 sklearn 等机器学习库,相信各位更能体会到调包的便利和乐趣。之后再去找数据实战和打比赛做项目,相信你一定会成为一名优秀的机器学习算法工程师。

 

本机器学习系列文章的两个主题就是数学推导+纯 numpy 实现。第一讲我们从最基础的线性回归模型开始。相信大家对回归算法一定是相当熟悉了,特别是咱们有统计背景的同学。所以,笔者直接上数学推导。

 

回归分析的数学推导

 本来想着上笔者的手推草稿的,但字迹过于张扬,在 word 里或者用 markdown 写公式又太耗费时间,这里就直接借用周志华老师的机器学习教材上的推导过程:

1538028852921638492.jpg
1538028860805880106.jpg
1538028868636135412.jpg
1538028882912308355.jpg

1538028890973371376.jpg

 推广到矩阵形式就是:

1538028904174076130.jpg

1538028910919602066.jpg

 

以上便是线性回归模型中参数估计的推导过程。

 

回归分析的 numpy 实现

按照惯例,动手写算法之前我们需要理一下编写思路。回归模型主体部分较为简单,关键在于如何在给出 mse 损失函数之后基于梯度下降的参数更新过程。首先我们需要写出模型的主体和损失函数以及基于损失函数的参数求导结果,然后对参数进行初始化,最后写出基于梯度下降法的参数更新过程。当然,我们也可以写出交叉验证来得到更加稳健的参数估计值。话不多说,直接上代码。

 

回归模型主体:

  1. import numpy as np
  2. def linear_loss(X, y, w, b):
  3. num_train = X.shape[0]
  4. num_feature = X.shape[1]
  5. # 模型公式
  6. y_hat = np.dot(X, w) + b
  7. # 损失函数
  8. loss = np.sum((y_hat-y)**2)/num_train
  9. # 参数的偏导
  10. dw = np.dot(X.T, (y_hat-y)) /num_train
  11. db = np.sum((y_hat-y)) /num_train
  12. return y_hat, loss, dw, db
复制代码
 

参数初始化:

  1. def initialize_params(dims):
  2. w = np.zeros((dims, 1))
  3. b = 0
  4. return w, b
复制代码

 

 

基于梯度下降的模型训练过程:

  1. def linar_train(X, y, learning_rate, epochs):
  2. w, b = initialize(X.shape[1])
  3. loss_list = []
  4. for i in range(1, epochs):
  5. # 计算当前预测值、损失和参数偏导
  6. y_hat, loss, dw, db = linar_loss(X, y, w, b)
  7. loss_list.append(loss)
  8. # 基于梯度下降的参数更新过程
  9. w += -learning_rate * dw
  10. b += -learning_rate * db
  11. # 打印迭代次数和损失
  12. if i % 10000 == 0:
  13. print('epoch %d loss %f' % (i, loss))
  14. # 保存参数
  15. params = {
  16. 'w': w,
  17. 'b': b
  18. }
  19. # 保存梯度
  20. grads = {
  21. 'dw': dw,
  22. 'db': db
  23. }
  24. return loss_list, loss, params, grads
复制代码

以上便是线性回归模型的基本实现过程。下面以 sklearn 中的 diabetes 数据集为例进行简单的训练。

 

 

数据准备:

  1. from sklearn.datasets import load_diabetes
  2. from sklearn.utils import shuffle
  3. diabetes = load_diabetes()
  4. data = diabetes.data
  5. target = diabetes.target
  6. # 打乱数据
  7. X, y = shuffle(data, target, random_state=13)
  8. X = X.astype(np.float32)
  9. # 训练集与测试集的简单划分
  10. offset = int(X.shape[0] * 0.9)
  11. X_train, y_train = X[:offset], y[:offset]
  12. X_test, y_test = X[offset:], y[offset:]
  13. y_train = y_train.reshape((-1,1))
  14. y_test = y_test.reshape((-1,1))
  15. print('X_train=', X_train.shape)
  16. print('X_test=', X_test.shape)
  17. print('y_train=', y_train.shape)
  18. print('y_test=', y_test.shape)
复制代码

 

1538029046027583585.jpg

 

执行训练:

  1. loss_list, loss, params, grads = linar_train(X_train, y_train, 0.001, 100000)
复制代码
 

 

查看训练得到的回归模型参数值:

  1. print(params)
复制代码
 

1538029098016755959.jpg

 

下面定义一个预测函数对测试集结果进行预测:

  1. def predict(X, params):
  2. w = params['w']
  3. b = params['b']
  4. y_pred = np.dot(X, w) + b
  5. return y_pred
  6. y_pred = predict(X_test, params)
  7. y_pred[:5]
复制代码



1538029129433634228.jpg





 

利用 matplotlib 对预测结果和真值进行展示:

  1. import matplotlib.pyplot as plt
  2. f = X_test.dot(params['w']) + params['b']
  3. plt.scatter(range(X_test.shape[0]), y_test)
  4. plt.plot(f, color = 'darkorange')
  5. plt.xlabel('X')
  6. plt.ylabel('y')
  7. plt.show()
复制代码



1538029167657377074.jpg





 

可见全变量的数据对于线性回归模型的拟合并不好,一来数据本身的分布问题,二来简单的线性模型对于该数据拟合效果差。当然,我们只是为了演示线性回归模型的基本过程,不要在意效果。

 

 

训练过程中损失的下降:

  1. plt.plot(loss_list, color = 'blue')
  2. plt.xlabel('epochs')
  3. plt.ylabel('loss')
  4. plt.show()
复制代码

1538029198536094591.jpg

 

 

 

封装一个线性回归类

笔者对上述过程进行一个简单的 class 封装,其中加入了自定义的交叉验证过程进行训练:

  1. import numpy as np
  2. from sklearn.utils import shuffle
  3. from sklearn.datasets import load_diabetes
  4. class lr_model():
  5. def __init__(self):
  6. pass
  7. def prepare_data(self):
  8. data = load_diabetes().data
  9. target = load_diabetes().target
  10. X, y = shuffle(data, target, random_state=42)
  11. X = X.astype(np.float32)
  12. y = y.reshape((-1, 1))
  13. data = np.concatenate((X, y), axis=1)
  14. return data
  15. def initialize_params(self, dims):
  16. w = np.zeros((dims, 1))
  17. b = 0
  18. return w, b
  19. def linear_loss(self, X, y, w, b):
  20. num_train = X.shape[0]
  21. num_feature = X.shape[1]
  22. y_hat = np.dot(X, w) + b
  23. loss = np.sum((y_hat-y)**2) / num_train
  24. dw = np.dot(X.T, (y_hat - y)) / num_train
  25. db = np.sum((y_hat - y)) / num_train
  26. return y_hat, loss, dw, db
  27. def linear_train(self, X, y, learning_rate, epochs):
  28. w, b = self.initialize_params(X.shape[1])
  29. for i in range(1, epochs):
  30. y_hat, loss, dw, db = self.linear_loss(X, y, w, b)
  31. w += -learning_rate * dw
  32. b += -learning_rate * db
  33. if i % 10000 == 0:
  34. print('epoch %d loss %f' % (i, loss))
  35. params = {
  36. 'w': w,
  37. 'b': b
  38. }
  39. grads = {
  40. 'dw': dw,
  41. 'db': db
  42. }
  43. return loss, params, grads
  44. def predict(self, X, params):
  45. w = params['w']
  46. b = params['b']
  47. y_pred = np.dot(X, w) + b
  48. return y_pred
  49. def linear_cross_validation(self, data, k, randomize=True):
  50. if randomize:
  51. data = list(data)
  52. shuffle(data)
  53. slices = [data[i::k] for i in range(k)]
  54. for i in range(k):
  55. validation = slices[i]
  56. train = [data
  57. for s in slices if s is not validation for data in s]
  58. train = np.array(train)
  59. validation = np.array(validation)
  60. yield train, validation
  61. if __name__ == '__main__':
  62. lr = lr_model()
  63. data = lr.prepare_data()
  64. for train, validation in lr.linear_cross_validation(data, 5):
  65. X_train = train[:, :10]
  66. y_train = train[:, -1].reshape((-1, 1))
  67. X_valid = validation[:, :10]
  68. y_valid = validation[:, -1].reshape((-1, 1))
  69. loss5 = []
  70. loss, params, grads = lr.linear_train(X_train, y_train, 0.001, 100000)
  71. loss5.append(loss)
  72. score = np.mean(loss5)
  73. print('five kold cross validation score is', score)
  74. y_pred = lr.predict(X_valid, params)
  75. valid_score = np.sum(((y_pred - y_valid) ** 2)) / len(X_valid)
  76. print('valid score is', valid_score)
复制代码

以上便是本节的内容,基于 numpy 手动实现一个简单的线性回归模型。

 

参考资料:

周志华 机器学习

回复

使用道具 举报

 
*滑块验证:
您需要登录后才可以回帖 登录 | 立即注册

本版积分规则


1关注

0粉丝

1603帖子

排行榜

关注我们:微信订阅号

官方微信

APP下载

全国服务热线:

4000-018-018

公司地址:上海市嘉定区银翔路655号B区1068室

运营中心:成都市锦江区东华正街42号广电仕百达国际大厦25楼

邮编:610066 Email:3318850993#qq.com

Copyright   ©2015-2016  比特趋势Powered by©Discuz!技术支持:迪恩网络