逻辑回归
Logistic Regression 虽然被称为回归，但其实际上是分类模型，并常用于二分类。Logistic Regression 因其简单、可并行化、可解释强深受工业界喜爱。Logistic 回归的本质是：假设数据服从这个分布，然后使用极大似然估计做参数的估计。逻辑回归的原理是用逻辑函数把线性回归的结果(-∞,∞)映射到(0,1)。 Logistic 分布是一种连续型的概率分布，其分布函数和密度函数分别为：
F
(
x
)
=
P
(
X
≤
x
)
=
1
1
+
e
−
(
x
−
μ
)
/
γ
f
(
x
)
=
F
′
(
x
)
=
e
−
(
x
−
μ
)
/
γ
γ
(
1
+
e
−
(
x
−
μ
)
/
γ
)
2
F(x)=P(X\leq x)=\frac{1}{1+e^{-(x-\mu)/\gamma}}\\ f(x)=F'(x)=\frac{e^{-(x-\mu)/\gamma}}{\gamma(1+e^{-(x-\mu)/\gamma})^2}
F(x)=P(X≤x)=1+e−(x−μ)/γ1f(x)=F′(x)=γ(1+e−(x−μ)/γ)2e−(x−μ)/γ 其中，
μ
\mu
μ表示位置参数，
γ
\gamma
γ 为形状参数。我们可以看下其图像特征： %pip install sympy
import sympy
import numpy as np
import matplotlib.pyplot as plt
# -----------构造数据--------------
# u=0 σ=1
x_1, m_1 = [], []
y_1, n_1 = [], []
for i in np.arange(-10, 10, 0.1):
x_1.append(i)
y_1.append(np.exp(-(i * i / 2.0))/np.sqrt(2*np.pi))
x = sympy.symbols('x')
y = sympy.exp(-(x * x / 2.0))/sympy.sqrt(2*np.pi)
Y = sympy.integrate(y)
res = Y.subs(x, i) - Y.subs(x, -float('inf'))
m_1.append(i)
n_1.append(res)
# u=0 σ=2
x_2, m_2 = [], []
y_2, n_2 = [], []
for i in np.arange(-10, 10, 0.1):
x_2.append(i)
y_2.append(np.exp(-(i * i / (2.0*2*2)))/(np.sqrt(2*np.pi)*2))
x = sympy.symbols('x')
y = sympy.exp(-(x * x / (2.0*2*2))) / (sympy.sqrt(2 * np.pi)*2)
Y = sympy.integrate(y)
res = Y.subs(x, i) - Y.subs(x, -float('inf'))
m_2.append(i)
n_2.append(res)
# u=1 σ=1
x_3, m_3 = [], []
y_3, n_3 = [], []
for i in np.arange(-10, 10, 0.1):
x_3.append(i)
y_3.append(np.exp(-((i-1) ** 2 / 2.0)) / np.sqrt(2 * np.pi))
x = sympy.symbols('x')
y = sympy.exp(-((x-1) ** 2 / 2.0)) / sympy.sqrt(2 * np.pi)
Y = sympy.integrate(y)
res = Y.subs(x, i) - Y.subs(x, -float('inf'))
m_3.append(i)
n_3.append(res)
# --------------可视化概率密度函数---------------------------
plt.title("正态分布概率密度") # 设置标题
plt.rcParams['font.sans-serif'] = 'KaiTi' # 设置字体
plt.rcParams['axes.unicode_minus'] = False # 正常显示负号
plt.plot(x_1, y_1, c='red', label='μ=0 σ=1') # 添加数据
plt.plot(x_2, y_2, c='green', label='μ=0 σ=2') # 添加数据
plt.plot(x_3, y_3, c='blue', label='μ=1 σ=1') # 添加数据
plt.xlabel('x', loc='right') # x轴标签
plt.ylabel('y', loc='top') # y轴标签
plt.legend(framealpha=1, frameon=True) # 添加图标
plt.show() # 展示
# --------------可视化分布函数---------------------------------
plt.title("正态分布分布函数")
plt.rcParams['font.sans-serif'] = 'KaiTi'
plt.rcParams['axes.unicode_minus'] = False
plt.plot(m_1, n_1, c='red', label='μ=0 σ=1')
plt.plot(m_2, n_2, c='green', label='μ=0 σ=2')
plt.plot(m_3, n_3, c='blue', label='μ=1 σ=1')
plt.xlabel('x', loc='right')
plt.ylabel('y', loc='top')
plt.legend(framealpha=1, frameon=True)
plt.show()
sigmoid函数是一个s形曲线，就像是阶跃函数的温和版，阶跃函数在0和1之间是突然的起跳，而sigmoid有个平滑的过渡。从图形上看，sigmoid曲线就像是被掰弯捋平后的线性回归直线，将取值范围(−∞,+∞)映射到(0,1) 之间，更适宜表示预测的概率，即事件发生的“可能性” 。 Logistic 分布是由其位置和尺度参数定义的连续分布。Logistic 分布的形状与正态分布的形状相似，但是 Logistic 分布的尾部更长，所以我们可以使用 Logistic 分布来建模比正态分布具有更长尾部和更高波峰的数据分布。在深度学习中常用到的 Sigmoid 函数就是 Logistic 的分布函数在
μ
=
1
,
γ
=
0
\mu=1,\gamma=0
μ=1,γ=0的特殊形式。 数据说明：多元变量（二元），数据共有3列。前2列表示2个exam score，由2个考试分数共同决定该条样本是“通过考试”还是“没通过”考试，第3列0表示没通过，1表示通过 import pandas as pd # 加载数据用
import numpy as np # 数值计算，矩阵处理用
from matplotlib import pyplot as plt # 绘图用
from scipy.optimize import minimize # 优化器用
from sklearn.preprocessing import PolynomialFeatures # 多项式计算用
def sigmoid(z):
return 1./(1+np.exp(-z)) #
# 计算损失函数
def compute_loss(theta, X, y):
if np.ndim(theta) == 1:
theta = theta[:, np.newaxis]
z_x = X.dot(theta) # h=θ^T dot X=θ0*x0+θ1*x1
# print(z_x.shape) # (100,3)x(3,1)=>(100,1)
h_x = sigmoid(z_x) # =>(100,1)
m = y.size # 100，用于对loss求平均
J_loss = -1./m * (np.log(h_x).T.dot(y) + np.log(1-h_x).T.dot(1-y)) # =>(1,1) # 损失函数核心点
if np.isnan(J_loss):
return np.inf
return J_loss[0][0]
# 绘制数据点，便于查看数据分布
def plot_data(data, label_pos, label_neg, axes=None):
# 获取正负样本的index
pos = data[:,2] == 1
neg = data[:,2] == 0
if axes is None:
axes = plt.gca()
axes.scatter(data[pos][:,0], data[pos][:,1], marker='+', c='k', s=60, linewidth=2, label=label_pos)
axes.scatter(data[neg][:,0], data[neg][:,1], marker='x',c='y', s=60, linewidth=2, label=label_neg)
axes.legend()
data_dir = './data/'
data = np.loadtxt(data_dir+'data1.txt', delimiter=',')
print('data shape:', data.shape) # (100,3)
print(data[:5])
plot_data(data, 'y=1', 'y=0')
plt.gca()matplotlib.pyplot.gca — Matplotlib 3.7.1 文档，matplotlib库的pyplot模块中的gca()函数用于获取与给定关键字args匹配的当前图形上的当前Axes实例,可以实现移动，自定义坐标轴的效果。函数gcf(）与gca()分别得到当前的figure与axes。（get current figure, get current axes). 利用sca()函数实现把一个axes对象作为当前的axes，axes对象为figure对象的一个属性，当切换了axes,figure肯定也切换了。 # Implementation of matplotlib function
import matplotlib.pyplot as plt
import numpy as np
import matplotlib.gridspec as gridspec
from mpl_toolkits.axes_grid1 import make_axes_locatable
plt.close('all')
arr = np.arange(100).reshape((10, 10))
fig = plt.figure(figsize =(4, 4))
im = plt.imshow(arr, interpolation ="none", cmap ="plasma")
divider = make_axes_locatable(plt.gca())
cax = divider.append_axes("left", "15 %", pad ="30 %")
plt.colorbar(im, cax = cax)
fig.suptitle('matplotlib.pyplot.gca() function Example', fontweight ="bold")
plt.show()
# 把θ0与x0=[1,1,...,1] 与特征x1,x2,..一起组成一个大X
X = np.c_[np.ones(data.shape[0]), data[:,0:2]] # 只要了特征
print(X[:5], X.shape) # (100, 3)
y = np.c_[data[:,2]] # 只要了标签
print(y[:5], y.shape) # (100,1)
loss = compute_loss(np.array([[0],[0],[0]]), X, y) # 100个3维度的特征值传到损失函数中会进行计算获得预测值，然后与y做计算出loss
print('\nloss:', loss)
[[ 1. 34.62365962 78.02469282]
[ 1. 30.28671077 43.89499752]
[ 1. 35.84740877 72.90219803]
[ 1. 60.18259939 86.3085521 ]
[ 1. 79.03273605 75.34437644]] (100, 3)
[[0.]
[0.]
[0.]
[1.]
[1.]] (100, 1)
loss: 0.6931471805599453
手动实现梯度计算
# 计算梯度
def compute_gradient(theta, X, y):
if np.ndim(theta) == 1:
theta = theta[:, np.newaxis]
# 计算model输出
z_x = X.dot(theta) # θ^T dot X=θ0*x0+θ1*x1
# print(z_x.shape) # (100,3)x(3,1)=>(100,1)
h_x = sigmoid(z_x) # =>(100,1)
m = y.size # 100 对应loss中求平均时的m
# 计算梯度并更新参数 X.T.dot(h_x-y):(3,100)x(100,1)=>(3,1)
grad = 1./m * X.T.dot(h_x-y) # (3,1) # 复习一下矩阵求导
return grad.flatten() # (3,)
initial_theta = np.zeros(X.shape[1]) # (3,) # 初始化参数，也可随机初始化
initial_theta = initial_theta[:, np.newaxis] # (3,1)
print(initial_theta.shape) # (3,1)
loss = compute_loss(initial_theta, X, y)
grad = compute_gradient(initial_theta, X, y)
print('loss:', loss)
print('grad:', grad)
(3, 1)
loss: 0.6931471805599453
grad: [ -0.1 -12.00921659 -11.26284221]
使用梯度下降求解参数，利用scipy.optimize.minimize()函数，指定最小化的目标即损失函数的计算， 然后是梯度的计算，优化参数theta，输入输出(X, y)，即可实现最小化损失函数。
from scipy.optimize import minimize
res = minimize(compute_loss, initial_theta, args=(X,y),jac=compute_gradient, options={'maxiter':400})
print(res)
# 查看最后优化得到的参数theta的值
final_updated_theta = res.x
print('theta:', final_updated_theta, final_updated_theta.shape) # (3,)
message: Optimization terminated successfully.
success: True
status: 0
fun: 0.20349770158950983
x: [-2.516e+01 2.062e-01 2.015e-01]
nit: 25
jac: [-2.686e-09 4.364e-07 -1.397e-06]
hess_inv: [[ 2.853e+03 -2.329e+01 -2.274e+01]
[-2.329e+01 2.045e-01 1.730e-01]
[-2.274e+01 1.730e-01 1.962e-01]]
nfev: 33
njev: 29
theta: [-25.16131634 0.2062316 0.20147143] (3,)
根据求解出的参数휃绘制出决策边界：
# 绘制决策边界
# 边界：
x1_min, x1_max = X[:,1].min(), X[:,1].max() # 注意X[0]是常数列[1,1,...,1]
x2_min, x2_max = X[:,2].min(), X[:,2].max()
xx1, xx2 = np.meshgrid(np.linspace(x1_min, x1_max), np.linspace(x2_min, x2_max)) # 生成网格坐标
# print(xx1.shape, xx2.shape) # (50, 50) (50, 50)
XX = np.c_[np.ones(xx1.ravel().shape[0]), xx1.ravel(), xx2.ravel()]
# print(XX.shape) # (2500, 3)
z_x = XX.dot(final_updated_theta) # (2500, 3)x(3,)=>(2500,)
h_x = sigmoid(z_x)
h_x = h_x.reshape(xx1.shape) # (50,50)
# print(h_x.shape)
plt.contour(xx1, xx2, h_x, [0.5], linewidths=1, colors='b') # 决策边界
plot_data(data, 'y=1', 'y=0') # # 样本数据点
plt.legend(loc=1)
X,Y=numpy.meshgrid(x,y)——生成网格点坐标矩阵。下图中，每个交叉点都是网格点，描述这些网格点的坐标的矩阵，就是坐标矩阵。输入的x，y，就是网格点的横纵坐标列向量（非矩阵） 输出的X，Y，就是坐标矩阵。 import numpy as np
import matplotlib.pyplot as plt
x = np.linspace(0,1000,20)
y = np.linspace(0,500,20)
X,Y = np.meshgrid(x, y)
plt.plot(X, Y,
color='limegreen', # 设置颜色为limegreen
marker='.', # 设置点类型为圆点
linestyle='') # 设置线型为空，也即没有线连接点
plt.grid(True)
plt.show()
plt.contour是python中用于画等高线的函数. import numpy as np
import matplotlib.pyplot as plt
x = np.linspace(-3, 3, 50) # 生成连续数据
y = np.linspace(-3, 3, 50) # 生成连续数据
X, Y = np.meshgrid(x, y)
# 生成能够在坐标系中形成点阵的数组，这个可以去参考一下别的文章
# https://lixiaoqian.blog.csdn.net/article/details/81532855 这里讲的比较详细
Z = X**2 + Y**2 # 这里将高度设置为x^2+y^2，就能画一个圆形的等高线
C=plt.contour(x, y,Z,[2,5,8,10]) # 画等高线 # 使用plt.contour(X, Y,Z,[2,5,8,10])也是没问题的
plt.clabel(C, inline=True, fontsize=10)
fig=plt.figure()
fig = plt.figure(figsize=(10,10))
ax1 = plt.axes(projection='3d')
ax1.scatter3D(X,Y,Z, cmap='Blues')
x = np.linspace(-3, 3, 50)
y = np.linspace(-3, 3, 50)
X, Y = np.meshgrid(x, y)
z = (np.exp(-X**2 - Y**2) - np.exp(-(X - 1)**2 - (Y - 1)**2))*2
fig=plt.figure()
fig = plt.figure(figsize=(10,10))
ax1 = plt.axes(projection='3d')
ax1.scatter3D(X,Y,z, cmap='Blues')
加正则项的逻辑斯蒂回归的损失函数手动实现和加正则项的梯度计算的手动实现
# 计算损失函数
def compute_loss_reg(theta, reg, *args):
if np.ndim(theta) == 1:
theta = theta[:, np.newaxis] # (28,)=>(28,1)
XX, y = args
z_x = XX.dot(theta) # h=θ^T dot X=θ0*x0+θ1*x1
# print(z_x.shape) # (118,28)x(28,1)=>(118,1)
h_x = sigmoid(z_x) # =>(118,1)
m = y.size # 118 用于求平均损失
J_loss = -1./m * (np.log(h_x).T.dot(y) + np.log(1-h_x).T.dot(1-y)) # =>(1,1)
# 注意：正则项的求和中j是从1开始的，而不是0，因为约束的是特征，而\theta_0对应的是x_0=1常数项
reg_item = reg/(2.*m) * np.sum(np.square(theta[:,0][1:]))
J_loss = J_loss + reg_item
if np.isnan(J_loss):
return np.inf
return J_loss[0][0]
# 计算梯度
def compute_gradient_reg(theta, reg, *args):
if np.ndim(theta) == 1:
theta = theta[:, np.newaxis] # (28,)=>(28,1)
XX, y = args
# 计算model输出
z_x = XX.dot(theta) # θ^T dot X=θ0*x0+θ1*x1
# print(z_x.shape) # (118,28)x(28,1)=>(118,1)
h_x = sigmoid(z_x) # =>(118,1)
m = y.size # 118 对应于平均损失中的m
# 计算梯度并更新参数 X.T.dot(h_x-y):(28,118)x(118,1)=>(28,1)
grad = 1./m * XX.T.dot(h_x-y) # (28,1)
reg_item = reg/m * np.r_[ [[0]], theta[1:] ] # (1,1) 行拼接 (27,1)=>(28,1)
grad = grad + reg_item # (28,1)
return grad.flatten() # (28,)
数据集简要说明:多元变量（二元），数据共有3列。前2列表示2个Microchip Test的值，由2个测试结果共同决定该芯片是“通过”还是“没通过”测试,第3列0表示没通过，1表示通过;数据分布是非线性的 data2 = np.loadtxt('./data/data2.txt', delimiter=',')
print('data shape:', data2.shape) # (100,3)
print(data2[:5])
plot_data(data2, 'y=1', 'y=0')
可见该数据不是直接线性可分的，不能直接使用逻辑斯蒂回归。考虑在“线性回归原理小结”中提到“多项式回归”， 这里使用“多项式特征”（sklearn.preprocessing.PolynomialFeatures）对原始数据X XX进行特征转化，对转化后的数据X ′ X’X ′再使用逻辑斯蒂回归。由于数据比之前的复杂，这次便可在学习时加入正则项以避免过拟合。 X = data2[:, 0:2]
y = np.c_[data2[:,2]]
print(X.shape, y.shape) # (118, 2) (118, 1)
poly = PolynomialFeatures(6) # 最高次项为6次
XX = poly.fit_transform(X) # X是有2个特征，XX有28个特征（含组合特征）
XX.shape # (118, 28)
# 0次项：1个，1次项：2个，2次项：3个（x1^2,x2^2,x1x2），3次项：4个
# 4次项：5个，5次项：6个，6次项：7个，一共28个特征
initial_theta = np.zeros(XX.shape[1]) # (28,)
initial_theta = initial_theta[:, np.newaxis]
print(initial_theta.shape) # (28,1)
loss = compute_loss_reg(initial_theta, 1, XX, y)
grad = compute_gradient_reg(initial_theta, 1, XX, y)
print('loss:', loss)
print('grad:', grad)
(118, 2) (118, 1)
(28, 1)
loss: 0.6931471805599454
grad: [8.47457627e-03 1.87880932e-02 7.77711864e-05 5.03446395e-02
1.15013308e-02 3.76648474e-02 1.83559872e-02 7.32393391e-03
8.19244468e-03 2.34764889e-02 3.93486234e-02 2.23923907e-03
1.28600503e-02 3.09593720e-03 3.93028171e-02 1.99707467e-02
4.32983232e-03 3.38643902e-03 5.83822078e-03 4.47629067e-03
3.10079849e-02 3.10312442e-02 1.09740238e-03 6.31570797e-03
4.08503006e-04 7.26504316e-03 1.37646175e-03 3.87936363e-02]
绘制决策边界，并查看不同的正则化项系数，对于决策边界的影响
lambda = 0 : 就是没有正则化，此时会容易过拟合lambda = 1 : 合适的正则化项系数lambda = 100 : 正则化项太激进，导致基本没拟合出决策边界 # 预测
def predict(theta, X, threshold=0.5):
z_x = X.dot(theta) # (b,3)x(3,1)=>(b,1)
h_x = sigmoid(z_x) # (b,1)
pred = h_x >= threshold
return pred.astype('int') # (b,1) b个样本的预测结果
fig, axes = plt.subplots(1, 3, sharey=True, figsize=(17,5))
print(axes.shape) # (3,)
lambda_list = [0., 1., 100.]
for i,C in enumerate(lambda_list):
# 最小化损失
res2 = minimize(compute_loss_reg, initial_theta, args=(C,XX,y),
jac=compute_gradient_reg, options={'maxiter':3000})
final_updated_theta = res2.x # (28,)
# print('theta:', final_updated_theta.shape) # (28,)
# 计算准确率 XX:(118, 28) y:(118, 1) theta:=>(28,1)
pred = predict(final_updated_theta.reshape(-1,1), XX) # (118,1)
# print(pred.shape) # (118,1)
accuracy = 100. * sum(pred.ravel() == y.ravel()) / y.size
# 绘制原始数据分布
plot_data(data2, 'y=1', 'y=0', axes=axes.flatten()[i])
# 绘制决策边界 X:(118,2) y:(118,1)
# 边界：
x1_min, x1_max = X[:,0].min(), X[:,0].max() # 注意这里的X不含常数列[1,1,...,1]
x2_min, x2_max = X[:,1].min(), X[:,1].max()
xx1, xx2 = np.meshgrid(np.linspace(x1_min, x1_max), np.linspace(x2_min, x2_max)) # 生成网格坐标
# print(xx1.shape, xx2.shape) # (50, 50) (50, 50)
X_ = np.c_[xx1.ravel(), xx2.ravel()]
# print(X_.shape) # (2500,2) # 没算常数列
XX_ = poly.fit_transform(X_) # (2500,28)
# print(XX_.shape) # (2500,28)
z_x = XX_.dot(final_updated_theta) # (2500, 28)x(28,)=>(2500,)
h_x = sigmoid(z_x) # (2500,)
h_x = h_x.reshape(xx1.shape) # (50,50)
#print(h_x.shape) # (50,50)
axes.flatten()[i].contour(xx1, xx2, h_x, [0.5], linewidths=1, colors='g') # 决策边界
# 设置标题
axes.flatten()[i].set_title('Train accuracy {:.2f}% with Lambda = {}'.format(accuracy, C))
线性回归pytorch-把线性回归实现一下。原理到实现，python到pytorch_羞儿的博客-CSDN博客
# 导包
import pandas as pd # 读取文件
import numpy as np # 做数据运算依赖包
from matplotlib import pyplot as plt # 绘图
from sklearn.linear_model import LinearRegression # 导入机器学习算法
def compute_loss(X, y, theta):
if np.ndim(theta) == 1:
theta = theta[:, np.newaxis]
h_x = X.dot(theta) # h=θ^T dot X=θ0*x0+θ1*x1
# print(h_x.shape) # (97,2)x(2,1)=>(97,1)，97个样本，每个对应两个属性，需要学习求解这两个属性的权重，
m = y.size # 97
J_loss = 1./(2*m) * np.sum(np.square(h_x-y))
# 预测值 h_x减去真实值 y的差值平方，求和除以样本2倍数，其实就是方差的二分之一
return J_loss
data_dir = './data/'
data = np.loadtxt(data_dir + 'linear_regression_data1.txt', delimiter=',')
print(data.shape) # (97, 2)
# X0 即theta0 对应的那一列
X = np.c_[np.ones(data.shape[0]), data[:,0]]# (97,)=>(97, 2)
print(X[:5], X.shape) # (97, 1)
y = np.c_[data[:,1]] # (97,)=>(97, 1)
print(y[:5], y.shape) # (97, 1)
plt.scatter(X[:,1], y, s=30, c='r', marker='x', linewidths=1)
plt.xlim(4, 24)
plt.xlabel('Population of City in 10,000s')
plt.ylabel('Profit in $10,000s')
np.ndim()，numpy数组中这四个方法的区别ndim、shape、dtype、astype、newaxis
ndim返回的是数组的维度，返回的只有一个数，该数即表示数组的维度。 shape：表示各位维度大小的元组。返回的是一个元组。 dtype：一个用于说明数组数据类型的对象。返回的是该数组的数据类型。 astype：转换数组的数据类型。 newaxis：是numpy中的一个函数，顾名思义，就是插入新维度的意思，比如将一维数组变成二维数组，二维变成三维等。 numpy.dot() 和 x.dot(y) 为矩阵乘法计算。x.dot(y) 等价于 np.dot(x,y) ———x是m × n 矩阵 ，y是n×m矩阵，则x.dot(y) 得到m×m矩阵。 np.r_是按列连接两个矩阵，就是把两矩阵上下相加，要求列数相等。np.c_是按行连接两个矩阵，就是把两矩阵左右相加，要求行数相等。
import numpy as np
x = np.arange(12).reshape(3,4)
print('x:',x, x.shape)
y = np.arange(10,22).reshape(3,4)
print('y:',y, y.shape)
z = np.c_[x,y]
print('z:',z, z.shape)
x: [[ 0 1 2 3]
[ 4 5 6 7]
[ 8 9 10 11]] (3, 4)
y: [[10 11 12 13]
[14 15 16 17]
[18 19 20 21]] (3, 4)
z: [[ 0 1 2 3 10 11 12 13]
[ 4 5 6 7 14 15 16 17]
[ 8 9 10 11 18 19 20 21]] (3, 8)
import numpy as np
x = np.arange(12).reshape(3,4)
print('x:',x, x.shape)
y = np.arange(10,22).reshape(3,4)
print('y:',y, y.shape)
z = np.r_[x,y]
print('z:',z, z.shape)
x: [[ 0 1 2 3]
[ 4 5 6 7]
[ 8 9 10 11]] (3, 4)
y: [[10 11 12 13]
[14 15 16 17]
[18 19 20 21]] (3, 4)
z: [[ 0 1 2 3]
[ 4 5 6 7]
[ 8 9 10 11]
[10 11 12 13]
[14 15 16 17]
[18 19 20 21]] (6, 4)
np.r_，np.c_和numpy.vstack，numpy.hstack的功能不完全相同。np.c_或np.r_用在两个一维数组上时，作用机制应该是先把这两个一维数组转置成列向量，然后拼接。 设置参数并进行训练 # 参数theta默认初始化为0
def gradient_descent(X, y, theta, alpha=0.01, num_iters=1000): # 梯度下降
if np.ndim(theta) == 1:
theta = theta[:, np.newaxis]
m = y.size
J_losses = [] # 存放每个step过程中loss，便于绘图查看loss随着优化更新的变化
for _ in range(num_iters):
# 计算model输出
h_x = X.dot(theta) # (50,1)
# print(h_x.shape)
# 计算梯度并更新参数 X.T.dot(h_x-y):(2,50)x(50,1)=>(2,1)
grad = 1./m * X.T.dot(h_x-y) # (2,1)
theta = theta - alpha*grad # (2,1) 学习率为0.01
#print(theta)
# 计算损失
loss = compute_loss(X, y, theta)
J_losses.append(loss)
return (theta, J_losses)
theta, J_losses = gradient_descent(X, y, [[0],[0]], num_iters=1500) # 参数theta默认初始化为0
print('theta:', theta.ravel()) # 拉平 (2,1)=>(2,)
plt.plot(J_losses)
plt.xlim(0)
plt.ylabel('loss_mse(mean_square_error)')
plt.xlabel('iter_steps')
plt.legend(labels=["loss"])
plt.show()
加载sklrean中的算法对比自己设计的梯度下降算法和损失函数 xx = np.arange(5, 25)
xx = np.c_[np.ones(len(xx)), xx]
print(xx[:5], xx.shape) # (20,2)
yy = xx.dot(theta) # (20,2)x(2,1)=>(20,1)
print(yy[:5], yy.shape)
# 拟合曲线
plt.plot(xx[:,1], yy, label='linear regression (gradient descent)')
# 原始数据点
plt.scatter(X[:,1], y, s=30, c='r', marker='x', linewidths=1)
# sklearn 库中的线性回归模型拟合的曲线（与自己实现的对比）
regr = LinearRegression()
regr.fit(X[:,1].reshape(-1,1), y.ravel()) # X[:,1]的shape为(97,),reshape=>(97,1), y (97,1)=>拉平(97,)
yy_regr = regr.intercept_ + regr.coef_*xx[:,1]
plt.plot(xx[:,1], yy_regr, label='linear regression (sklearn)')
plt.xlim(4,24)
plt.legend(loc=4)
plt.show()
分析
逻辑回归是在线性回归的基础上加了一个 Sigmoid 函数（非线形）映射，使得逻辑回归称为了一个优秀的分类算法。本质上来说，两者都属于广义线性模型，但他们两个要解决的问题不一样，逻辑回归解决的是分类问题，输出的是离散值，线性回归解决的是回归问题，输出的连续值。 线性回归是在实数域范围内进行预测，而分类范围则需要在 [0,1]，逻辑回归减少了预测范围；线性回归在实数域上敏感度一致，而逻辑回归在 0 附近敏感，在远离 0 点位置不敏感，这个的好处就是模型更加关注分类边界，可以增加模型的鲁棒性。 逻辑回归和最大熵模型本质上没有区别，最大熵在解决二分类问题时就是逻辑回归，在解决多分类问题时就是多项逻辑回归。 与 SVM相比，都是分类算法，本质上都是在找最佳分类超平面；都是监督学习算法；都是判别式模型，判别模型不关心数据是怎么生成的，它只关心数据之间的差别，然后用差别来简单对给定的一个数据进行分类；都可以增加不同的正则项。
LR 是一个统计的方法，SVM 是一个几何的方法；SVM 的处理方法是只考虑 Support Vectors，也就是和分类最相关的少数点去学习分类器。而逻辑回归通过非线性映射减小了离分类平面较远的点的权重，相对提升了与分类最相关的数据点的权重； 损失函数不同：LR 的损失函数是交叉熵，SVM 的损失函数是 HingeLoss，这两个损失函数的目的都是增加对分类影响较大的数据点的权重，减少与分类关系较小的数据点的权重。对 HingeLoss 来说，其零区域对应的正是非支持向量的普通样本，从而所有的普通样本都不参与最终超平面的决定，这是支持向量机最大的优势所在，对训练样本数目的依赖大减少，而且提高了训练效率； LR 是参数模型，SVM 是非参数模型，参数模型的前提是假设数据服从某一分布，该分布由一些参数确定（比如正太分布由均值和方差确定），在此基础上构建的模型称为参数模型；非参数模型对于总体的分布不做任何假设，只是知道总体是一个随机变量，其分布是存在的（分布中也可能存在参数），但是无法知道其分布的形式，更不知道分布的相关参数，只有在给定一些样本的条件下，能够依据非参数统计的方法进行推断。所以 LR 受数据分布影响，尤其是样本不均衡时影响很大，需要先做平衡，而 SVM 不直接依赖于分布； LR 可以产生概率，SVM 不能；LR 不依赖样本之间的距离，SVM 是基于距离的；LR 相对来说模型更简单好理解，特别是大规模线性分类时并行计算比较方便。而 SVM 的理解和优化相对来说复杂一些，SVM 转化为对偶问题后，分类只需要计算与少数几个支持向量的距离，这个在进行复杂核函数计算时优势很明显，能够大大简化模型和计算。
好文链接
发表评论