紫梦沁香
0%

发表于 更新于 分类于 阅读次数： Valine：
本文字数： 472 阅读时长 ≈ 1 分钟

一. 公式

1. 排列组合

（1） 排列

$P_m^n = m(m - 1)(m - 2)···(m - n + 1)$

$P_m^n = \frac{m!}{(m - n)!}$

（2） 组合

$C_m^n = \frac{P_m^n}{n!} = \frac{m(m - 1)(m - 2)···(m - n + 1)}{n!}$

$C_m^n = \frac{m!}{n!(m - n)!}$

$C_m^n = C_m^{m-n}$

$C^n_{m+1} = \frac{(m+1)m(m-1)···(m-n+2)}{n!} = \frac{[(m - n + 1) + n]m(m-1)···(m-n+2)}{n!} = \frac{m(m-1)···(m-n+1)}{n!} + \frac{m(m-1)···(m-n+2)}{(n-1)!} = C_m^n + C_m^{n-1}$

2.条件概率

（1） 条件概率

$P(B|A) = \frac{P(AB)}{P(A)}$， P(A) > 0

（2） 乘法公式

$$

（3） 完备事件组

（4） 全概率公式

（5） 贝叶斯公式

发表于 更新于 分类于 阅读次数： Valine：
本文字数： 3.5k 阅读时长 ≈ 6 分钟

算法思想

观察下图，可以发现每个特征之间的量纲不同，没有办法放在一起比较，学习出来的模型==可解释性比较差==（不会影响预测结果）。而且对于KNN算法来说，样本之间的距离极大部分会由发现时间（天）来引导，肿瘤大小（厘米）几乎可以忽略不计。

为此我们就需要采用归一化，把有量纲的表达式改成无量纲的表达式，将数据的所有特征都映射到一个统一的尺度下。（对于其它算法如线性回归则不会影响分类或预测结果）

而且归一化后可以加速训练，左侧因为进行了归一化，因此等高线几乎是正圆，等高线上的法线指向圆心，因此能快速收敛到最小值，而右图没有进行归一化则需进行多次迭代。

没有进行归一化，学习率也得按照特征取值范围最小的那个维度来取，不然很容易出现nan。

注：由于进行了归一化，我们的模型都是按照归一化的数据进行训练的，因此测试数据也要进行归一化，否则没有意义了。

（1） 均值方差归一化（standardization）
$$
X_{scale} = \frac{X - X_{mean}}{X_{std}}
$$
经过该标准化处理的数据，均值为0，方差为1，符合正态分布。

适用于数据分布没有明显的边界，有可能存在极端的数据值。

正好学了概率论，就来简单推导下吧。
$$
设经过均值方差归一化后X^的数学期望和方差分别为E^X, D^X \qquad (X^ = \frac{X - EX}{\sqrt{DX}})
\
E^X = E(\frac{X - EX}{\sqrt{DX}}) = \frac{1}{\sqrt{DX}} * (EX - EX) = 0
\
D^X = D(\frac{X - EX}{\sqrt{DX}}) = \frac{1}{DX} * D(X - EX) = \frac{1}{DX} * DX = 1
$$
（2）最值归一化（Normalization）
$$
X_{scale} = \frac{X - X_{min}}{X_{max} - X_{min}}
$$
经过该标准化处理的数据，数据分布在(0, 1)。

适用于分布有明显边界的情况，受 outliner影响较大。

例如像素范围（0 ~ 255），学生成绩（0 ~ 100）。


实例

（1） 最值归一化

1
2
3
4
5
6
7
8
9
class MinMaxScaler:
    
    def fit(self, X):
        self.X = X
        self.max = np.max(X, axis = 0)
        self.min = np.min(X, axis = 0)
        
    def transform(self, X):
        return (X - self.min) / (self.max - self.min)

原始数据均值为46.995，标准差为27.97007284581147。

1
2
3
4
5
6
X = np.random.randint(0, 100, (100, 2))
print(np.mean(X))			# 46.995
print(np.std(X))			# 27.97007284581147

plt.scatter(X[:, 0], X[:, 1])
plt.show()

经过最值归一化后，均值为0.4721085858585859，标准差为0.2863025514861211。

所有值都被映射到（0， 1）之间

1
2
3
4
5
6
7
8
9
minMaxScaler = MinMaxScaler()
minMaxScaler.fit(X)
X1 = minMaxScaler.transform(X)

print(np.mean(X1))			# 0.4721085858585859
print(np.std(X1))			# 0.2863025514861211

plt.scatter(X1[:, 0], X1[:, 1])
plt.show()

（2） 均值方差归一化

1
2
3
4
5
6
7
8
9
class StandardScaler:
    
    def fit(self, X):
        self.X = X
        self.mean = np.mean(X, axis = 0)
        self.std = np.std(X, axis = 0)
        
    def transform(self, X):
        return (X - self.mean) / self.std

经过均值方差归一化后，均值几乎为0，标准差也几乎为1，但数据范围不一定在（0， 1）中。

1
2
3
4
5
6
7
8
9
standardScaler = StandardScaler()
standardScaler.fit(X)
X1 = standardScaler.transform(X)

print(np.mean(X1))				# -2.7755575615628914e-17
print(np.std(X1))				# 0.9999999999999998

plt.scatter(X1[:, 0], X1[:, 1])
plt.show()

测试

先介绍一下 train test split，正如它的字面意思。将一份数据集分成训练数据集和测试数据集，这样在使用 score 函数时，使用的都是训练时没有见过的数据集，得出来的准确率也就更可靠。

首先数据集需要打乱，因为有些数据集是根据类别来分类的，例如前 k 个都是类别1，后 k 个都是类别2，这样 train 和test中的类别分布差距很大，会导致泛化效果不好。

np.random.seed(seed) 使用同一个种子为了能更好的比较使用归一化和不使用归一化的区别。

1
2
3
4
5
def train_test_split(X, y, split_rate, seed):
    np.random.seed(seed)
    length = int(len(X) * split_rate)
    p = np.random.permutation(len(X))
    return X[p[length:]], X[p[:length]], y[p[length:]], y[p[:length]]

（1） 不使用归一化进行学习

1
2
3
4
5
6
7
8
9
10
boston=datasets.load_boston()
X=boston.data			# (506, 13)
y=boston.target			# (506, )

split_rate = 0.2
X_train, X_test, y_train, y_test = train_test_split(X, y, split_rate)

lr = LinearRegression()
lr.fit(X_train, y_train, eta = 1e-6, n_iters = 1e6)
lr.score(X_test, y_test)

eta = 1e-6

score = 31.348463775149956

（2） 使用最值归一化

1
2
3
4
5
6
7
8
9
10
11
12
13
14
boston=datasets.load_boston()
X=boston.data
y=boston.target

minMaxScaler = MinMaxScaler()
minMaxScaler.fit(X)
X = minMaxScaler.transform(X)

split_rate = 0.2
X_train, X_test, y_train, y_test = train_test_split(X, y, split_rate)

lr = LinearRegression()
lr.fit(X_train, y_train, eta = 1e-2, n_iters = 1e6)
lr.score(X_test, y_test)

eta = 1e-2

score = 25.209485526542817

（3） 使用均值方差归一化

1
2
3
4
5
6
7
8
9
10
11
12
13
14
boston=datasets.load_boston()
X=boston.data
y=boston.target

standardScaler = StandardScaler()
standardScaler.fit(X)
X = standardScaler.transform(X)

split_rate = 0.2
X_train, X_test, y_train, y_test = train_test_split(X, y, split_rate)

lr = LinearRegression()
lr.fit(X_train, y_train, eta = 1e-2, n_iters = 1e6)
lr.score(X_test, y_test)

eta = 1e-2

score = 18.92473923096906

经过对比，使用同样的训练次数，经过了归一化的预测率明显比没有经过归一化的预测率更高，而且eta取1e-2就行，不用多次尝试。

发表于 更新于 分类于 阅读次数： Valine：
本文字数： 4.5k 阅读时长 ≈ 8 分钟

算法思想

用于解决回归问题，做法就是找一条最能拟合测试数据的直线，求该拟合直线的斜率和截距就是最值问题，可以通过梯度下降法解决。预测数据也比较简单，将数据代入直线中求出结果即可。


实例

在举例之前，先介绍几种求误差的方法。

  1. MSE (Mean Square Error) 均方误差

    $\frac{1}{m}\sum\limits_{n=1}^{m} (y - y_n)^2$

  2. RMSE (Root Mean Square Error) 均方根误差

    $\sqrt{\frac{1}{m}\sum\limits_{n=1}^m (y - y_n)^2}$

  3. MAE (Mean Absolute Error) 平均绝对误差

    $\frac{1}{m}\sum\limits_{n=1}^{m} |y - y_n|$

    上述三种中，均方误差和均方根误差只是差了个根号，本质是一样的，都是求预测值与真实值的平方差。

    因为是平方差，会对极度偏颇的数据特别敏感，因此会着重去拟合这部分数据。

    而绝对值对所有的数据都是一样的，没有侧重点。

1
2
3
4
5
6
7
8
import numpy as np
import matplotlib.pyplot as plt

x = np.arange(-10, 10, 0.2).reshape([100,1])
y = 3.5 * x + np.random.normal(0, 20, 100).reshape([100,1])

plt.scatter(x,y)
plt.plot([-10, 10],[-55, 55], color = "red")

如上图所示，蓝色的离散点就是我们的测试数据，红色是要求的拟合曲线，方程为$y = ax + b$。

我们需要让这条直线的预测值和真实值的误差尽可能小，因此成本函数为

$cost = \frac{1}{m}\sum\limits_{n=1}^{m} (ax_n + b - y_n)^2$ （MSE）

此处$x_n$和$y_n$都是已知的数据，a和b才是我们要求的未知数。

该成本函数的偏导函数为

$\frac{dcost}{da} = \frac{2}{m}\sum\limits_{n=1}^{m} x_n(ax_n + b - y_n)$

$\frac{dcost}{db} = \frac{2}{m}\sum\limits_{n=1}^{m} (ax_n + b - y_n)$

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
def J(a, b, x, y):
    return 1/len(x) * np.sum((a * x + b - y) ** 2)

def DJ_a(a, b, x, y):
    return 2 / len(x) * np.sum(x * (a * x + b - y))

def DJ_b(a, b, x, y):
    return 2 / len(x) * np.sum(a * x + b - y)

def gradient_descent(a, b, eta, epsilon = 1e-6):
    i_iter = 0
    while True:
        last_a = a
        last_b = b
        a = a - eta * DJ_a(a, b, x, y)
        b = b - eta * DJ_b(a, b, x, y)
        if np.abs(J(last_a, last_b, x, y) - J(a, b, x, y)) < epsilon:
            break
    return a, b
1
2
3
4
5
6
7
8
9
10
a = 0
b = 0
eta = 0.02
epsilon = 1e-6

a, b = gradient_descent(a, b, eta, epsilon)

plt.plot([-10, 10],[-10 * a + b, 10 * a + b], color = "red")
plt.scatter(x,y)
plt.show()

效果也是挺不错的，不过这是一元的情况，可以写两个DJ_a和DJ_b，但如果是百元呢？

下面将会实现一种通用的梯度下降法。


先更改一下我们的成本函数，其中θ是系数矩阵，对于一元函数，它的shape就是(2, )，x需要再添一列（数值为1），这样x * a + 1 * b就正好是我们的直线方程预测值。

$cost = \frac{1}{m}\sum (x.dot(θ) - y)^2$ （MSE）

导函数

$\frac{dcost}{dθ} = \frac{2}{m} x.t.dot((x.dot(θ) - y))$

推导过程

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
def J(theta, x, y):
    return 1/len(x) * np.sum((x.dot(theta) - y) ** 2)

def DJ(theta, x, y):
    return 2 / len(x) * (x.T.dot(x.dot(theta) - y))

def gradient_descent(theta, x, y, eta, n_iters = 2e5, epsilon = 1e-6):
    i_iter = 0
    
    while i_iter < n_iters:
        last_theta = theta
        theta = theta - eta * DJ(theta, x, y)
        if np.abs(J(last_theta, x, y) - J(theta, x, y)) < epsilon:
            break
        i_iter += 1
        
    return theta

然后我们来拟合刚刚随机的数据

1
2
3
4
5
6
7
8
9
10
11
12
x = np.hstack([x, np.ones([100, 1])])

theta = np.zeros([2])
eta = 0.02
y = y.reshape([100])
theta = gradient_descent(theta, x, y, eta)

plt.plot([-10, 10],[-10 * theta[0] + theta[1], 10 * theta[0] + theta[1]], color = "red")
plt.scatter(x[:, 0],y.reshape([100, 1]))
plt.show()

theta 		# array([3.4686148, -1.06100372])

显然，这条直线还是能较好拟合这个测试数据的。

theta[0] = 3.4686148 ≈ 3.5，也和我们预设的基本一致。

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
class LinearRegression:
    
    def __init__(self):
        self.X = None
        self.y = None
        
    def fit(self, X, y):
        self.X = X
        self.y = y
        self.len = len(X)
        self.X_b = np.hstack([self.X, np.ones([self.len, 1])])
        return self
    
    def J(self, theta):
        return np.sum((self.y - self.X_b.dot(theta)) ** 2) / self.len
    
    def DJ(self, theta):
        return self.X_b.T.dot(self.X_b.dot(theta) - y) * 2 / self.len
    
    def transform(self, eta = 0.01, n_iters = 1e5, epsilon = 1e-8):
        assert self.X_b is not None and self.y is not None, \
        "must fit before transform"
        
        theta = np.zeros(self.X.shape[1] + 1)
        i_iter = 0
        while i_iter < n_iters:
            last_theta = theta
            theta = theta - eta * self.DJ(theta)
            if np.abs(self.J(theta) - self.J(last_theta)) < epsilon:
                break
            i_iter += 1
        self.theta = theta
        self.intercept = self.theta[-1]
        self.coefficient = self.theta[:-1]
        return theta
    
    def predict(self, X):
        return X.dot(self.coefficient) + self.intercept
    
    def score(self, X, y):
        y_pred = X.dot(self.coefficient) + self.intercept
        return np.sum((y_pred - y) ** 2) / len(X)

需要注意的是，第28行中theta = theta - eta * self.DJ(theta)，不能写成theta -= eta * self.DJ(theta)，这是因为现在last_theta和theta都指向同一块内存单元，-=是原地操作，在更改了theta的同时，last_theta也被更改了，而=则是重新开辟一块新的内存空间，然后再让theta指向它。


测试

1
2
3
4
5
6
7
8
9
10
import sklearn
from sklearn import datasets

boston=datasets.load_boston()
X=boston.data
y=boston.target

lr = LinearRegression()
lr.fit(X, y)
lr.transform(eta = 1e-6, n_iters = 1e6)
1
2
3
4
5
6
7
8
"""lr.coefficient
array([-0.09218218,  0.06825598, -0.00529991,  0.46137796,  0.17983213,
        4.30941577,  0.02543425, -0.70440603,  0.17038895, -0.00977886,
       -0.09173025,  0.01804173, -0.55937148])
       
   lr.intercept 	
   0.4361364266724813
"""

因为数据是13维的，没有办法可视化，我们只能使用MSE指标来衡量，其中lr.score(X, y) = 26.11494。

为了能看得清楚点，我们只使用其中一维来进行预测。

1
2
3
X = X[:, 5]
plt.scatter(X, y)
plt.show()

这里上面有一排数值都为50，这是因为在做调查的时候，最高一个档就是50（大于50W都只能填50W），因此这部分数据我们当做无效数据，直接丢弃掉。

1
2
3
4
5
X = X[y < 50]
y = y[y < 50]

plt.scatter(X, y)
plt.show()

做一下筛选，只选y<50的数据就行了。

然后开始整体，进行预测。

1
2
3
4
5
6
7
8
lr = LinearRegression()
X = X.reshape([-1, 1])
lr.fit(X, y)
theta = lr.transform(eta = 1e-6, n_iters = 1e6)

plt.plot(X, lr.predict(X), color="red")
plt.scatter(X, y)
plt.show()

lr.score(X, y) = 42.28978618974813

发表于 更新于 分类于 阅读次数： Valine：
本文字数： 1.8k 阅读时长 ≈ 3 分钟

算法思想

梯度下降法是一种非常通用的优化算法，能够为大范围的问题找到最优解。梯度下降的中心思想就是迭代地调整参数从而使成本函数最小化。（PS：梯度上升法解决最大化问题）

斜率为正，代表单调递增，说明极值点在左边；

斜率为负，代表单调递减，说明极值点在右边。

因此我们可以通过减去导数来找到最值点（假如极值点只有一个），因此被称为梯度下降法。

不过由于每个点上的导数是固定的，例如接近笔直的直线斜率就几乎是无穷大，如果我们直接减去导数，显然就不大合适。

所以就引申出了eta，也被称作学习率（learning rate）。现在我们只需减去eta * gradient就可以了，通过调整eta的大小来改变拟合的速度

eta过大可能导致无法拟合，eta过小则会拟合过慢，下面会有实例展示。


实例

这是eta适中的情况，只花了26步就找到了最小值。

1
2
3
4
5
6
start = -10.0
epsilon = 1e-5
eta = 0.4

# 26 steps
plot_gradient_tracks(start, epsilon, eta)

这里eta就偏小了，虽然也能找到最小值，但运行了875次，是上面的五十多倍。

1
2
3
4
5
6
start = -10.0
epsilon = 1e-5
eta = 0.01

# 875 steps
plot_gradient_tracks(start, epsilon, eta)

这里eta也算是适中的情况，只不过图像会比较特殊，左右反复横跳，因此我也拿出来当个样例。

1
2
3
4
5
6
start = -10.0
epsilon = 1e-5
eta = 1.5

# 21 steps
plot_gradient_tracks(start, epsilon, eta)

这是最差的情况，eta过大，导致y越来越大，已经无法取得最小值了。

1
2
3
4
5
start = -10.0
epsilon = 1e-5
eta = 2.5

plot_gradient_tracks(start, epsilon, eta, max_len = 4)


当然也会有存在多个极值点的情况，下图中如果选择左边作为起始点，基本就找不到全局最优解了，只能找到局部最优解。

不过幸好，一般使用MSE(Mean Square Error)作为成本函数，而MSE恰好是个凸函数，这就意味着连接曲线上任意两点的线段永远也不会跟曲线相交。也就是说不存在局部最小值，只有一个全局最小值。它同时也是一个连续函数，所以斜率不会产生陡峭的变化。

这两点保证的结论是：即便是乱走，梯度下降也可以趋近到全局最小值（只要等待时间足够长，学习率也不是太高）。

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
import numpy as np
import matplotlib.pyplot as plt

# 目标函数
def J(x):
    return 0.5 * x ** 2 + 3.5 * x + 2

# 目标导函数
def DJ(x):
    return x + 3.5

# 绘制导数轨迹曲线
def plot_gradient_tracks(start, epsilon, eta, max_len = 1e4):
    a = start

    historys = []
    while(max_len > 0):
        last = J(a)
        historys.append([a, J(a)])
        a = a - eta * DJ(a)
        if np.abs(last - J(a)) < epsilon:
            break
        max_len -= 1

    historys = np.array(historys)
    plt.plot(x, y)
    plt.plot(historys[:, 0], historys[:, 1])
    plt.show()
    print(historys.shape[0])
    return historys[-1]

测试

1
2
3
4
5
start = -10.0
epsilon = 1e-5
eta = 1.5

plot_gradient_tracks(start, epsilon, eta)		# array([-3.5000062, -4.125])

将目标函数配方得，$\frac{1}{2}x^2 + \frac{7}{2}x + 2 = \frac{1}{2}[(x + \frac{7}{2})^2 + \frac{33}{4}]$

可以得出在$x = -\frac{7}{2}$处，取得最小值$y = \frac{33}{8}$

与测试结果基本一致。

注：计算结果不可能与理论值相等。