机器学习 —— 用于超参数随机化搜索的几个分布

Posted 2021-02-05 massquantity

机器学习中超参数搜索的常用方法为 Grid Search，然而如果参数一多则容易碰到维数诅咒的问题，即参数之间的组合呈指数增长。如果有 (m) 个参数，每个有 (n) 个取值，则时间复杂度为 (Theta(n^m))。 Bengio 等人在 《Random Search for Hyper-Parameter Optimization》 中提出了随机化搜索的方法。他们指出大部分参数空间存在 “低有效维度 (low effective dimensionality)” 的特点，即有些参数对目标函数影响较大，另一些则几乎没有影响。而且在不同的数据集中通常有效参数也不一样。 在这种情况下 Random Search 通常效果较好，下图是一个例子，其中只有两个参数，绿色的参数影响较大，而黄色的参数则影响很小：

Grid Search 会评估每个可能的参数组合，所以对于影响较大的绿色参数，Grid Search 只探索了3个值，同时浪费了很多计算在影响小的黄色参数上； 相比之下 Random Search 则探索了9个不同的绿色参数值，因而效率更高，在相同的时间范围内 Random Search 通常能找到更好的超参数 (当然这并不绝对)。 另外，Random Search 可以在连续的空间搜索，而 Grid Search 则只能在离散空间搜索，而对于像神经网络中的 learning rate，SVM 中的 gamma 这样的连续型参数宜使用连续分布。

在实际的应用中，Grid Search 只需为每个参数事先指定一个参数列表就可以了，而 Random Search 则通常需要为每个参数制定一个概率分布，进而从这些分布中进行抽样。然而对什么样的参数应该选择什么样的分布？这就大有讲究了，如果选的分布不恰当可能就永远找不到合适的参数值了，本文主要介绍一些超参数搜索的常用分布以及它们的特点和使用范围。这些分布都出自 scipy.stats 模块，共同特点是提供了 rvs 方法用于独立随机抽样。

Randint 分布

Randint 分布的概率质量函数 (PMF) 为：

[ f(x) = frac{1}{high - low} ]

其中 (x = low,,...,high - 1) ，下面画出随机抽样10000次后各个取值的分布图：

np.random.seed(42)
randint = sp.stats.randint(low=-10, high=11)
randint_distribution = randint.rvs(size=10000, random_state=42)

start = randint.ppf(0.01)
end = randint.ppf(0.99)
x = np.arange(start, end+1)

randint_dict = dict(zip(*np.unique(randint_distribution, return_counts=True)))   # 计算各个数的频次
randint_count = list(map(lambda x: x[1], sorted(list(randint_dict.items()), key=lambda x: x[0])))

plt.figure(figsize=(8,5))
plt.bar(x, randint_count, color=‘b‘, alpha=0.5, edgecolor=‘k‘,  label="random_samples")
plt.axhline(y=450, xmin=0.01, xmax=0.99, color=‘#FF00FF‘, linestyle="--")
plt.legend(frameon=False, fontsize=10)
plt.title("randint distribution", fontsize=17)
plt.show()

从上图可以看出 Randint 分布为离散型均匀分布，适用于必须为整数的参数 (比如神经网络的层数，决策树的深度)。

Uniform 分布

Uniform 是 Randint 分布的连续版本，概率密度函数为：
[ f(x) = frac{1}{high - low} ]
其中 (x in [low, high])

np.random.seed(42)
uniform = sp.stats.uniform(loc=10, scale=10)
uniform_distribution = uniform.rvs(size=10000, random_state=42)
start = uniform.ppf(0.01)
end = uniform.ppf(0.99)
x = np.linspace(start, end, num=10000)

plt.figure(figsize=(8,5))
plt.plot(x, uniform.pdf(x), ‘r--‘, lw=2, label="uniform distribution PDF")
plt.hist(uniform_distribution, bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘, normed=True, label="random samples")
plt.legend(frameon=False, fontsize=10)
plt.title("uniform distribution", fontsize=17)
plt.show()

Geometric 分布

Geometric 分布的概率质量函数为 :
[ f(x) = (1-p)^{x-1} p ]
其中 (x geqslant 1)

np.random.seed(42)
plt.figure(figsize=(8,5))
geom_distribution = sp.stats.geom.rvs(0.5, size=10000, random_state=42)
plt.hist(geom_distribution, bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘, label="random samples")
plt.legend(frameon=False, fontsize=10)
plt.title("Geometric Distribution, p = 0.5", fontsize=17)
plt.show()

Geometric 分布为离散型分布，表示得到一次成功所需要的试验次数，如果参数集中于少数几个值且可能性呈离散型单调递减，则适用此分布。

Geometric 分布概率质量函数中的 (p) 指定了一次试验成功的概率。如果改变此值则会增大或缩小采样范围。

np.random.seed(42)
plt.figure(figsize=(15,5))
plt.subplot(121)
geom_distribution = sp.stats.geom.rvs(0.8, size=10000, random_state=42)
plt.hist(geom_distribution, bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘, label="random samples")
plt.legend(frameon=False, fontsize=10)
plt.title("Geometric Distribution, p = 0.8", fontsize=17)

plt.subplot(122)
geom_distribution = sp.stats.geom.rvs(0.01, size=10000, random_state=42)
plt.hist(geom_distribution, bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘, label="random samples")
plt.legend(frameon=False, fontsize=10)
plt.title("Geometric Distribution, p = 0.01", fontsize=17)
plt.show()

Exponential 分布

Exponential 分布是 Geometric 分布的连续版本，其概率密度函数为 ：
[ f(x) = e^{-x} ]
可以看到上图中当Geometric 分布中的 (p) 非常小时，就会变得非常接近 exponential 分布。

plt.figure(figsize=(16,4))
expon_distribution = sp.stats.expon.rvs(loc=0, scale=1, size=10000, random_state=42)

plt.subplot(121)
start = sp.stats.expon.ppf(0.001)
end = sp.stats.expon.ppf(0.999)
x = np.linspace(start, end, num=10000)
plt.plot(x, sp.stats.expon.pdf(x), ‘r--‘, lw=2, label=" exponential 
distribution PDF")
plt.hist(expon_distribution, bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘, normed=True, label="random samples")
plt.legend(frameon=False, fontsize=13)
plt.title("exponential distribution", fontsize=17)

plt.subplot(122)
plt.hist(np.log(expon_distribution), bins=30, edgecolor = ‘k‘, color=‘b‘, alpha=0.5)
plt.title("log of exponential distribution", fontsize=17)
plt.axvline(x=-5, color=‘#FF00FF‘, linestyle="--")
plt.axvline(x=2, color=‘#FF00FF‘, linestyle="--")
plt.show()

从右边的 log 分布图来看，大部分值集中于 (e^{-5}) 到 (e^2) 之间，即 (0.0067 sim 7.389) 。如果有一些先验知识，知道参数在0附近，且值越大可能性越小 (如svm中的gamma)，则适用此分布。当然也可以调整位置 (loc) 和 比例 (scale) 参数来改变搜索范围。此时对应的概率密度函数为 (下面演示 loc=10，scale=10 的情况)：
[ f(x) = frac{e^{-frac{x - loc}{scale}}}{scale} ]

Reciprocal 分布

reciprocal 分布的概率密度函数为：
[ f(x, a, b) = frac{1}{x ext{log}(b/a)} = frac{1}{x ext{log}b - x ext{log}a} ]
其中 (a < x < b, ; b > a > 0)

np.random.seed(42)
plt.figure(figsize=(15,5))
reciprocal = sp.stats.reciprocal(a=0.1, b=100)
reciprocal_distribution = reciprocal.rvs(size=10000, random_state=42)

plt.subplot(121)
start = reciprocal.ppf(0.3)
end = reciprocal.ppf(0.99)
x = np.linspace(start, end, num=10000)
plt.plot(x, reciprocal.pdf(x), ‘r--‘, lw=2, label="  reciprocal 
distribution PDF")
plt.hist(reciprocal_distribution, bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘, normed=True, label="random samples")
plt.legend(frameon=False, fontsize=13)
plt.title("reciprocal distribution", fontsize=17)

plt.subplot(122)
plt.hist(np.log10(reciprocal_distribution), bins=30, color=‘b‘, alpha=0.5, edgecolor = ‘k‘)
plt.title("log of reciprocal distribution", fontsize=17)
plt.show()

上图中 reciprocal 分布的PDF和 exponential 分布比较相似，然而右边的 log 分布图却是比较平均的，可见 reciprocal 分布是一个典型的对数均匀分布，以10为底为例，线性空间中10倍的差距在对数空间中均为1，设$x_2 = 10,x_1 $：
[ egin{align*} ext{log},f(x_1, a, b) - ext{log}, f(x_2, a, b) &= ext{log} frac{1}{x_1 ext{log}(b/a)} - ext{log} frac{1}{x_2 ext{log}(b/a)} \[1ex] & = - ext{log} , [x_1 ext{log}(b/a)] + ext{log} , [x_2 ext{log}(b/a)] \[1ex] & = - ext{log}, x_1 + ext{log} , x_2 \[1ex] & = ext{log} frac{x_2}{x_1} = 1 end{align*} ]

下面用 np.random.uniform 可以模拟类似的分布。

np.random.seed(42)
log_uniform = 10 ** np.random.uniform(-1, 2, size=10000)

plt.figure(figsize=(15,5))
plt.subplot(121)
plt.hist(log_uniform, bins=30, color=‘b‘, alpha=0.5, normed=True, edgecolor=‘k‘)
plt.title("$10^{(-1 sim 2)}$ distribution", fontsize=17)

plt.subplot(122)
plt.hist(np.log10(log_uniform), bins=30, color=‘b‘, alpha=0.5, edgecolor=‘k‘)
plt.title("log of $10^{(-1 sim 2)}$ distribution", fontsize=17)
plt.show()

这种分布的好处是在不同的取值范围内也能均匀地抽样。如上图中参数 (x) 的取值范围是 (0.1 sim 100)， 即 (10^{-1} sim 10^2)，如果是一般的均匀分布中抽样，(10 sim 100) 这个范围被取样到的概率会远大于 (1 sim 10) 和 (0.1 sim 1) 这两个范围，因为前者的距离更大，但在对数均匀分布中三者的范围却是一样的，都是10的倍数，这样被抽样到的概率也就类似。下面的代码显示一个例子：

a = 0
b = 0
c = 0

reciprocal_distribution = sp.stats.reciprocal.rvs(a=0.1, b=100, size=10000, random_state=42)

for val in reciprocal_distribution:
    if val > 10 and val < 100:
        a += 1
    elif val > 1 and val < 10:
        b += 1
    elif val > 0.1 and val < 1:
        c +=1

print("10  到 100 之间取样 {} 次".format(a))  # 10 到 100 之间取样 3233 次
print("1   到 10  之间取样 {} 次".format(b))  # 1 到 10 之间取样 3392 次
print("0.1 到 1  之间取样 {} 次".format(c))   # 0.1 到 1 之间取样 3375 次

对于像 learning rate 这样的参数，我们希望 (0.01sim0.1) 和 (0.1sim1) 范围之间的抽样概率是类似的。举例来说，0.11和0.1的学习率可能相差不大，但0.01和0.02的学习率结果更可能大不相同，虽然这两个范围的绝对差异均为0.01。因此在这样的参数中不同值之间的比率更适合作为超参数变化范围。 另外实际上我们可以做到 “完全” 的对数均匀分布，这要用到 numpy 中的 logspace。 然而使用 np.logspace 的缺点是只能生成一个间隔均匀的固定数组进行采样，从而丧失了一定的随机性。

logspace = np.logspace(-1, 2, base=10, num=10000)
plt.figure(figsize=(15,4))
plt.subplot(121)
plt.hist(logspace, bins=30, color=‘b‘, alpha=0.5, normed=True, edgecolor=‘k‘)
plt.title("logspace", fontsize=17)
plt.subplot(122)
plt.hist(np.log10(logspace), bins=30, color=‘b‘, alpha=0.5, edgecolor=‘k‘)
plt.title("log of logspace", fontsize=17)
plt.show()

---

最后用scikit-learn中的 RandomizedSearchCV 来比较一下 Grid Search 和 Random Search 的效果，使用了 Kaggle 上的 HousePrices 比赛中的一个 Kernel 进行数据预处理，最后的特征数为410，使用的模型为超参数较多的 GBDT，评估指标为 RMSE：

1. 树的数量 (n_estimators)
2. 损失函数类型 (loss)
3. 学习率 (learning_rate)
4. 子采样率 (subsample)
5. 叶结点上的最少样本数 (min_samples_leaf)
6. 最大深度 (max_depth)
7. 分裂时考虑的特征数 (max_features)

下面先尝试 GridSearchCV

from sklearn.model_selection import GridSearchCV, RandomizedSearchCV

def neg_sqrt(val):  # 定义 RMSE
    return np.sqrt(-val)

model = GradientBoostingRegressor()
param_grid = {
    "learning_rate": [0.01, 0.05, 0.1],
    "n_estimators": [400, 800, 1500],
    "subsample": [0.8, 1.0],
    "max_depth": [2, 3, 4], 
    "max_features": [0.8, 1.0], 
    "min_samples_leaf": [1, 2],
    "loss": ["ls", "huber"],
    "random_state": [42],
}

grid_search = GridSearchCV(model, param_grid, scoring="neg_mean_squared_error", cv=5, verbose=3, n_jobs=-1)
grid_search.fit(X_scaled, y_log)
print("最优参数为: ", grid_search.best_params_, ‘
‘)
print("RMSE 为: ", neg_sqrt(grid_search.best_score_))

结果如下：

最优参数为:  {‘min_samples_leaf‘: 2, ‘learning_rate‘: 0.05, ‘max_depth‘: 2, ‘random_state‘: 42, ‘n_estimators‘: 1500, ‘loss‘: ‘huber‘, ‘subsample‘: 0.8, ‘max_features‘: 1.0} 

RMSE 为:  0.11852064590041982

上面过程中总共的参数组合为 (3 imes 3 imes 2 imes 3 imes 2 imes 2 imes 2 = 432) 个，接下来的RandomizedSearchCV 用了差不多的400个，其中 learning_rate 用了 reciprocal 这样的对数均匀分布，原因前文已经说了。叶结点上的最少样本数 (min_samples_leaf) 使用了 Geometric 分布，主要考虑到大部分值可能集中在 1和2左右。其他参数都使用均匀分布：

model = GradientBoostingRegressor()
param_distribution = {
    "learning_rate": sp.stats.reciprocal(a=0.01, b=0.1),  
    "n_estimators": sp.stats.randint(low=400, high=1500),
    "subsample": sp.stats.uniform(loc=0.8, scale=0.2),
    "max_depth": sp.stats.randint(low=2, high=4),
    "max_features": sp.stats.uniform(loc=0.8, scale=0.2),
    "min_samples_leaf": sp.stats.geom(p=0.6),
    "loss": ["ls", "huber"],
    "random_state": [42],
}

random_search = RandomizedSearchCV(model, param_distribution, n_iter=400, scoring="neg_mean_squared_error", cv=5, 
                                   verbose=3, random_state=42, n_jobs=-1) 
random_search.fit(X_scaled, y_log)
print("最优参数为: ", random_search.best_params_, ‘
‘)
print("RMSE 为: ", neg_sqrt(random_search.best_score_))

结果如下：

最优参数为:  {‘min_samples_leaf‘: 3, ‘learning_rate‘: 0.03181845026156779, ‘max_depth‘: 2, ‘random_state‘: 42, ‘n_estimators‘: 1476, ‘loss‘: ‘huber‘, ‘subsample‘: 0.8978905520555127, ‘max_features‘: 0.8557292928473224} 

RMSE 为:  0.11835604958840028

在这个数据集上 Random Search 的效果确实比 Grid Search 稍好，当然前提是为每个参数都选择合适的分布。

/