模型之母:简单线性回归的代码实现
关于作者:饼干同学,某人工智能公司交付开发工程师/建模科学家。专注于AI工程化及场景落地,希望和大家分享成长中的专业知识与思考感悟。
0x00 前言
在《模型之母:简单线性回归&最小二乘法》中,我们从数学的角度理解了简单线性回归,并且推导了最小二乘法。
本文内容完全承接于上一篇,我们来以代码的方式,实现简单线性回归。话不多说,码起来
0x01 简单线性回归算法的实现
首先我们自己构造一组数据,然后画图
import numpy as np
import matplotlib.pyplot as plt
x = np.array([1.,2.,3.,4.,5.])
y = np.array([1.,3.,2.,3.,5,])
plt.scatter(x,y)
plt.axis([0,6,0,6])
plt.show()
下面我们就可以根据样本真实值,来进行预测。
实际上,我们是假设线性关系为: 这根直线,然后再根据最小二乘法算a、b的值。我们还可以假设为二次函数:。可以通过最小二乘法算出a、b、c
实际上,同一组数据,选择不同的f(x),即模型,通过最小二乘法可以得到不一样的拟合曲线。
不同的数据,更可以选择不同的函数,通过最小二乘法可以得到不一样的拟合曲线。
下面让我们回到简单线性回归。我们直接假设是一条直线,模型是:
根据最小二乘法推导求出a、b的表达式:
下面我们用代码计算a、b:
# 首先要计算x和y的均值
x_mean = np.mean(x)
y_mean = np.mean(y)
# a的分子num、分母d
num = 0.0
d = 0.0
for x_i,y_i in zip(x,y): # zip函数打包成[(x_i,y_i)...]的形式
num = num + (x_i - x_mean) * (y_i - y_mean)
d = d + (x_i - x_mean) ** 2
a = num / d
b = y_mean - a * x_mean
在求出a、b之后,可以计算出y的预测值,首先绘制模型直线:
y_hat = a * x + b
plt.scatter(x,y) # 绘制散点图
plt.plot(x,y_hat,color='r') # 绘制直线
plt.axis([0,6,0,6])
plt.show()
然后进行预测:
x_predict = 6
y_predict = a * x_predict + b
print(y_predict)
5.2
0x02 向量化运算
我们注意到,在计算参数a时:
# a的分子num、分母d
num = 0.0
d = 0.0
for x_i,y_i in zip(x,y): # zip函数打包成[(x_i,y_i)...]的形式
num = num + (x_i - x_mean) * (y_i - y_mean)
d = d + (x_i - x_mean) ** 2
a = num / d
我们发现有这样一个步骤:向量w和向量v,每个向量的对应项,相乘再相加。其实这就是两个向量“点乘”
这样我们就可以使用numpy中的dot运算,非常快速地进行向量化运算。
总的来说:
向量化是非常常用的加速计算的方式,特别适合深度学习等需要训练大数据的领域。
对于 y = wx + b, 若 w, x都是向量,那么,可以用两种方式来计算,第一是for循环:
y = 0
for i in range(n):
y += w[i]*x[i]
y += b
另一种方法就是用向量化的方式实现:
y = np.dot(w,x) + b
二者计算速度相差几百倍,测试结果如下:
import numpy as np
import time
a = np.random.rand(1000000)
b = np.random.rand(1000000)
tic = time.time()
c = np.dot(a, b)
toc = time.time()
print("c: %f" % c)
print("vectorized version:" + str(1000*(toc-tic)) + "ms")
c = 0
tic = time.time()
for i in range(1000000):
c += a[i] * b[i]
toc = time.time()
print("c: %f" % c)
print("for loop:" + str(1000*(toc-tic)) + "ms")
运行结果:
c: 249981.256724
vectorized version:0.998973846436ms
c: 249981.256724
for loop:276.798963547ms
对于独立的样本,用for循环串行计算的效率远远低于向量化后,用矩阵方式并行计算的效率。因此:
只要有其他可能,就不要使用显示for循环。
0x03 自实现的工程文件
3.1 代码
还记得我们之前的工程文件吗?创建一个SimpleLinearRegression.py,实现自己的工程文件并调用
import numpy as np
class SimpleLinearRegression:
def __init__(self):
"""模型初始化函数"""
self.a_ = None
self.b_ = None
def fit(self, x_train, y_train):
"""根据训练数据集x_train,y_train训练模型"""
assert x_train.ndim ==1, \
"简单线性回归模型仅能够处理一维特征向量"
assert len(x_train) == len(y_train), \
"特征向量的长度和标签的长度相同"
x_mean = np.mean(x_train)
y_mean = np.mean(y_train)
num = (x_train - x_mean).dot(y_train - y_mean) # 分子
d = (x_train - x_mean).dot(x_train - x_mean) # 分母
self.a_ = num / d
self.b_ = y_mean - self.a_ * x_mean
return self
def predict(self, x_predict):
"""给定待预测数据集x_predict,返回表示x_predict的结果向量"""
assert x_predict.ndim == 1, \
"简单线性回归模型仅能够处理一维特征向量"
assert self.a_ is not None and self.b_ is not None, \
"先训练之后才能预测"
return np.array([self._predict(x) for x in x_predict])
def _predict(self, x_single):
"""给定单个待预测数据x_single,返回x_single的预测结果值"""
return self.a_ * x_single + self.b_
def __repr__(self):
"""返回一个可以用来表示对象的可打印字符串"""
return "SimpleLinearRegression()"
3.2 调用
下面我们在jupyter中调用我们自己写的程序:
首先创建一组数据,然后生成SimpleLinearRegression()
的对象reg1,然后调用一下
from myAlgorithm.SimpleLinearRegression import SimpleLinearRegression
x = np.array([1.,2.,3.,4.,5.])
y = np.array([1.,3.,2.,3.,5,])
x_predict = np.array([6])
reg = SimpleLinearRegression()
reg.fit(x,y)
输出:SimpleLinearRegression()
reg.predict(x_predict)
reg.a_
reg.a_
输出:array([5.2]) 0.8 0.39999999999999947
y_hat = reg.predict(x)
plt.scatter(x,y)
plt.plot(x,y_hat,color='r')
plt.axis([0,6,0,6])
plt.show()
0xFF 总结
在本篇文章中,我们实现了简单线性回归算法的代码,并且使用了向量化运算,事实证明,向量化运算能够提高运算效率。
同时我们发现,只要数学公式推导清楚了,实际写代码时没有太多难度的。
那么我们思考一个问题,在之前的kNN算法(分类问题)中,使用分类准确度来评价算法的好坏,那么回归问题中如何评价好坏呢?