<jupyter_start><jupyter_text>Методы подбора гиперпараметровДатасет: https://www.kaggle.com/c/bioresponse/data  Дисклеймер:1. Серебрянной пули – нет.2. Параметров – много, а времени – мало.3. Автоматический подбор параметров не всегда может спасти, см. предыдущий пункт.<jupyter_code>from __future__ import print_function, division, unicode_literals
from tqdm import tqdm_notebook
import numpy as np
import pandas as pd

import scipy

from sklearn.model_selection import cross_val_score
import numpy as np
import pandas as pd

import matplotlib.pyplot as plt

%matplotlib inline
df = pd.read_csv('data/train.csv', sep=',')
df.head()
df.shape
from sklearn.model_selection import train_test_split

X_train, X_valid, y_train, y_valid = \
    train_test_split(df.iloc[:, 1:].values, df.iloc[:, 0].values, test_size=0.3, random_state=9876)
plt.figure(figsize=(6, 3))
_ = plt.hist(y_train)
from sklearn import linear_model, ensemble
from sklearn.metrics import roc_auc_score<jupyter_output><empty_output><jupyter_text>__1__. GridSearch<jupyter_code>%%time

from sklearn.model_selection import GridSearchCV

params_grid = {
    'n_estimators': np.arange(50, 150, 25),
    'max_features': np.arange(0.5, 1.0, 0.1),
    'criterion':    ['gini', 'entropy'],
}

clf = ensemble.RandomForestClassifier()
model = GridSearchCV(estimator=clf, param_grid=params_grid, cv=3, verbose=5, n_jobs=16)
model.fit(X_train, y_train)

y_pred = model.best_estimator_.predict(X_valid)
roc_auc_score(y_valid, y_pred)
model.best_params_<jupyter_output><empty_output><jupyter_text>__2__. RandomizedSearch<jupyter_code>%%time

from sklearn.model_selection import RandomizedSearchCV

params_grid = {
    'n_estimators': np.arange(50, 150, 25),
    'max_features': scipy.stats.uniform(loc=0.5, scale=0.5),
    'criterion':    ['gini', 'entropy'],
}

clf = ensemble.RandomForestClassifier()
model = RandomizedSearchCV(estimator=clf, param_distributions=params_grid,
                           n_iter=40, cv=3, verbose=5, n_jobs=16)
model.fit(X_train, y_train)

y_pred = model.best_estimator_.predict(X_valid)
roc_auc_score(y_valid, y_pred)
model.best_params_<jupyter_output><empty_output><jupyter_text>Вопрос: чем `RandomizedSearchCV` лучше `GridSearchCV`?  __3__. Sequential Model Optimization Понятия:* objective function – оптмизируемая функция, $f$;* surrogate model – модель, которая приближает objective function;* acquisition function – функция выгоды, по которой выбираем следующий набор гиперпараметров. $$ \mathbf{P}(y|\theta, D) \rightarrow \max_\theta $$ Зададим acquisition function, как $a_f = \mathbf{E}[u(x)|D]$.* $u(x) = \mathbf{I}(f(x) \leq f')$ – probability of improvement;* $u(x) = max(f' - f(x), 0)$ – expected improvement. Принцип работы:1. Построим surrogate model для objective function.2. Найдем новую точку (значение гиперпаметров), которые максимизируют acquisition function.3. Посчитаем значение objective function в новой точке.4. Обновим surrogate model.5. Продолжаем выполнять шаги 2 – 4.  __3.1__. Гауссовские процессы* [`bayesian-optimization`](https://github.com/fmfn/BayesianOptimization);* [`scikit-optimize`](https://scikit-optimize.github.io/);* [`gpyopt`](http://sheffieldml.github.io/GPyOpt/).[Ссылка](https://papers.nips.cc/paper/4522-practical-bayesian-optimization-of-machine-learning-algorithms.pdf) на статью.<jupyter_code>%%time

from bayes_opt import BayesianOptimization

def model_auc_score(**params):
    params['n_estimators'] = int(round(params['n_estimators']))
    params['max_depth'] = int(round(params['max_depth']))
    
    clf = ensemble.GradientBoostingClassifier(**params)
    score = cross_val_score(clf, X_train, y_train, scoring='roc_auc', cv=3, n_jobs=3).mean()
    return score

bounds = {
    'learning_rate': (0.001, 0.1),
    'n_estimators':  (50, 200),
    'subsample':     (0.5, 1.0),
    'max_depth':     (3, 8),
}

optimizer = BayesianOptimization(
    f=model_auc_score,
    pbounds=bounds,
    random_state=9872,
)

optimizer.maximize(init_points=10, n_iter=30)
optimizer.max
%%time

from skopt import gp_minimize
from skopt.utils import use_named_args
from skopt.space import Real, Integer, Categorical

space = [
    Real(0.001, 0.1, prior='log-uniform', name=str('learning_rate')),
    Integer(50, 200, name=str('n_estimators')),
    Real(0.5, 1.0, prior='uniform', name=str('subsample')),
    Integer(3, 8, name=str('max_depth')),
]

@use_named_args(space)
def objective(**params):
    clf = ensemble.GradientBoostingClassifier(**params)
    score = cross_val_score(clf, X_train, y_train, scoring='roc_auc', cv=3, n_jobs=3).mean()
    return -score

res = gp_minimize(objective, space, n_random_starts=10, n_calls=30, verbose=1, random_state=9872)
res.x
res.fun
from skopt.plots import plot_evaluations

_ = plot_evaluations(res, dimensions=[s.name for s in space])
from skopt.plots import plot_objective

_ = plot_objective(res, dimensions=[s.name for s in space])<jupyter_output><empty_output><jupyter_text>Недостатки:* проблема холодного старта;* плохо работает в случае категориальных и дискретных параметров;* сильно подвержен проклятию размерности. Для решения проблемы холодного старта неплохо подойдут RandomizedSearch и GridSearch.  __3.2__. Деревья Чисто технически ничто не мешает заменить расчет параметров, среднее и матрицу ковариации, гауссовского процесса на что-то другое, например, на деревья.* [`scikit-optimize`](https://scikit-optimize.github.io/).<jupyter_code>from skopt import forest_minimize, gbrt_minimize<jupyter_output><empty_output><jupyter_text>**Tree-Structured Parzen Estimators**Идея: перейдем от рассмотрения $ \mathbf{P}(y|\theta, D) $ к $ \mathbf{P}(\theta|y, D) $.$$ \mathbf{P}(\theta|y, D) = \left\{\begin{matrix}l(x); & y < y^*  \\ g(x); & y \geq y^*\end{matrix}\right. $$где* $l(x)$ – генеративное распределение, построенное на основе наблюдений $\{x_i | f(x_i) < y^*\}$;* $g(x)$ – генеративное распределение, построенное на основе остальных наблюдений;* $y^*$ – выборочный квантиль порядка $\gamma$.Распределения $l(x)$ и $g(x)$ моделируются с помощью ядерной оценки плотности (KDE или Parzen density estimation).Оказывается, что в этом случае для expected improvement справедливо$$ a_f \propto \left(\gamma + \frac{g(x)}{l(x)} (1 - \gamma) \right)^{-1} $$Вопрос: что означает такое выражение для expected improvement?Библиотеки:* [`hyperopt`](http://hyperopt.github.io/hyperopt/);* [`optuna`](https://optuna.org/).<jupyter_code>%%time

from hyperopt import fmin, tpe, hp, Trials

space = {
    'learning_rate': hp.loguniform('learning_rate', np.log(0.001), np.log(0.1)),
    'n_estimators':  hp.uniformint('n_estimators', 50, 200),
    'subsample':     hp.uniform('subsample', 0.5, 1.0),
    'max_depth':     hp.choice('max_depth', range(3, 9)),
}

def objective(params):
    clf = ensemble.GradientBoostingClassifier(**params)
    score = cross_val_score(clf, X_train, y_train, scoring='roc_auc', cv=3, n_jobs=3).mean()
    return -score

trials = Trials()
res = fmin(objective, space, algo=tpe.suggest, max_evals=40, trials=trials,
           rstate=np.random.RandomState(9872), verbose=1)
res
f, axes = plt.subplots(nrows=1, ncols=len(space), figsize=(18, 5))

for ax, val in zip(axes, space):
    xs = [t['misc']['vals'][val] for t in trials.trials]
    ys = [-t['result']['loss'] for t in trials.trials]
    xs, ys = zip(*sorted(zip(xs, ys)))
    ys = np.array(ys)
    ax.scatter(xs, ys, s=20, linewidth=0.01, alpha=0.75)
    ax.set_title(val)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14283800080396622, 'fraction_numerical': 0.021573093930054937, 'mean_word_length': 5.564643799472296, 'pattern_counts': {'":': 0, '<': 27, '<?xml version=': 0, '>': 25, 'https://': 6, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '535073', 'license': 'MIT', 'n_tokens_mistral': 2963, 'n_tokens_neox': 2983, 'n_words': 594, 'path': '1_term/made_2019_ml/seminars/week03_hyperopt/hyperopt.ipynb', 'repo_name': 'Ilyabasharov/made_mail.ru', 'size': '7463'}