В середине июля закончился контест по машинному обучению ML Boot Camp V от Mail.Ru. Нужно было предсказать наличие сердечно-сосудистых заболеваний по результатам классического врачебного осмотра. Метрикой являлась логарифмическая функция потерь. Полное описание задачи доступно по ссылке.

Знакомство с машинным обучением для меня началось с ML Boot Camp III где-то в феврале 2017 года и какое-то подобие представления что с такими задачами делать начинает складываться у меня только сейчас. Многое из сделанного в 5 контесте — результат в первую очередь изучения собрания статей на kaggle и обсуждений и примеров кода оттуда же. Ниже — слегка переработанный отчет о том, что пришлось сделать, чтобы занять 3 место.

Данные задачи

Датасет сформирован из 100.000 реальных клинических анализов. Даны возраст, рост, вес, пол, артериальное давление верхнее и нижнее, холестерин и глюкоза в крови.

Дополнительно есть «субъективные» данные — что пациенты сами о себе сообщили, отвечая на вопросы о курении, потреблении алкоголя и физической активности. Эта часть данных была еще и подпорчена организаторами, так что никаких особых надежд с ними у меня связано не было.

Исходные данные содержали заведомо нереальные значения — были люди с ростом 50 см в 30+ лет, с давлением вроде 16020, с отрицательными давлениями давлений. Объяснялось это ошибками при ручном вводе данных анализов.

Инструменты

Задача решалась на python с использованием стандартных для такого случая библиотек:

• pandas – чтение-запись и обработка табличных данных (на самом деле много чего еще, но в данном случае остальное не понадобилось);
• NumPy – операции над массивами чисел;
• scikit-learn – набор инструментов для машинного обучения, включая базовые алгоритмы МО, разбиение данных, валидацию;
• XGBoost – одна из наиболее ходовых реализаций градиентного бустинга;
• LightGBM – альтернатива XGBoost;
• TensorFlow + Keras — библиотека для обучения и использования нейронных сетей и обертка к ней;
• Hyperopt – библиотека для оптимизации функций в заданном пространстве аргументов;

Csv vs pickle

Для хранения данных при долгих расчетах сначала использовал csv пока не понадобилось сохранять вместе более сложные структуры чем отдельные таблицы. Очень хорошо проявил себя модуль pickle — все нужные данные сохраняются или читаются в 2 строки кода. Позже сохранять стал в сжатые файлы:

with gzip.open('../run/local/pred_1.pickle.gz', 'wb') as f:
    pickle.dump((x, y), f)

Репозиторий

Весь код, относящийся к соревнованию лежит на github. В репозитории в old/ спрятаны старые скрипты, реальной пользы от которых не было, и оставленные только потому, что результаты их работы также отправлялись на проверку. Из-за ошибок в коде промежуточные результаты их выполнения позже оказались непригодны для использования, так что на итоговые решения никакого влияния эта часть кода не оказала.

Первые 2 недели

Первые 2 недели чистил данные, подставлял их в оставшиеся от прошлых соревнований модели и все это без особого успеха. Для каждого нового сабмита весь код из одного из существующих скриптов копировался целиком в новый и там уже редактировался. Результат – на исходе второй недели сказать с ходу что делает какой-то из последних скриптов и какие функции в нем реально используются, а какие просто занимают место, уже не мог. Код был громоздкий, трудночитаемый и мог проработать несколько часов и упасть не сохранив вообще ничего полезного.

Вторые 2 недели

Когда через 2 недели после начала копировать старые скрипты и немного их менять стало слишком сложно — пришлось начать полную переделку всего кода. Он был разбит на общие части — базовые классы и конкретные их реализации.

Общая идея новой организации кода — конвейер: данные > признаки > модели 1 уровня > модели 2 уровня. Каждый этап реализует отдельный файл скрипта, который при запуске проделывает все положенные вычисления, сохраняет их, промежуточные результаты и свои данные. Скрипт для каждого следующего этапа импортирует код предыдущих этапов и получает данные для обработки от их методов. Идея за всем этим — чтобы можно было запустить скрипт для одной из финальных моделей, тот запустил скрипты для моделей уровнем ниже, те вызвали нужные им генераторы признаков, которые в свою очередь запустили нужный вариант очистки данных. Задача каждого скрипта — проверить существует ли файл, куда он должен сохранять свои результаты и если нет — выполнить необходимые расчеты и сохранить данные.

За всем этим был план сделать решение для организации моделей и данных в расчете на будущее использование, попутно отладив его на подходящей задаче. Собственно это решение и было важнейшим результатом контеста, которое сейчас плавно перерастает в небольшую библиотеку для облегчения жизни при участии в последующих соревнованиях такого рода.

Общий план

Изначально было запланировано 2 уровня моделей и надо было заготовить как можно больше различных моделей на 1й уровень. Способ этого добиться подготовить как можно больше по-разному обработанных данных, на которых будут обучены однотипные модели. Но готовить данные — дело долгое. Хоть работа с данными и является залогом успеха (достаточное количество добавленных осмысленных признаков позволяет обойтись простейшими моделями), но времени требует больше чем хотелось бы. Альтернатива — решение грубой силой, то есть относительно умеренная обработка данных и максимальное машинное время на вычисления.

Самое простое при таком подходе — обработать данные несколькими способами, придумать несколько наборов дополнительных признаков и использовать их комбинации. Получится куцый вариант Random subspace method, отличающийся от полноценного тем, что он совсем не random и признаки выбираются сразу группами. Так при небольшом количестве дополнительных признаков можно получить сотни вариантов обработанных данных (собственно количество способов чистки * (2 ^ количество групп признаков)). Предполагалось, что такой подход даст достаточно разные решения простых моделей, использующих разные подмножества признаков, чтобы каждое из них улучшало качество моделей 2 уровня.

Подготовка данных

То, что исходные данные были грязными надо было как-то учитывать. Основные подходы — или выкидывать все заведомо невозможные значения или пытаться как-то восстановить исходные данные. Так как источник таких искажений оставался неизвестен почти до конца — пришлось готовить данные несколькими способами, после чего обучать на них разные модели.

Каждый из вариантов обработки данных реализован классом, при обращении возвращающим полный датасет с соответствующими изменениями. Так как обработка данных на этом этапе по времени проходит достаточно быстро — промежуточные результаты сохранялись только для относительно долгого варианта (2) — восстановления субъективных признаков с помощью xgboost. Остальные данные генерировались по запросу.

Варианты обработки:

1. Исходные данные в которых испорченные значения субъективной части теста заменены на 0.0001, чтобы привести их к числовому виду, но отличать от нетронутых.
2. Испорченные субъективные признаки заменены: потребление алкоголя — 0, активность — 1. Далее по остальным колонкам данных было «восстановлено» курение.
3. В данных с восстановленными субъективными признаками почищены экстремальные значения давлений.
4. В данных с восстановленными субъективными признаками (из п.2) почищены экстремальные значения давлений, веса и роста.
5. В данных только с очищенными давлениями (из п.3) дополнительно почищены вес, рост, давления.
6. Данные с очищенными давлениями дополнительно преобразованы — любое отдельное неправдоподобное значение роста, веса или давления заменены на NaN.

Признаки

Из обработанных данных генерировались дополнительные признаки. Осмысленных из них было совсем немного — индекс массы тела, ожидаемые значения давления в зависимости от пола, веса и возраста по какой-то старой формуле и т.п. Гораздо больше дополнительных колонок с данными было получено автоматически достаточно простыми способами.

Дополнительные признаки генерировались из разных вариантов обработанных данных, но часто одним и тем же способом. Так как часть признаков могла требовать слишком много времени на повторный расчет их значений — колонки признаков сохранялись отдельно. Вычисление признаков в скриптах было реализовано аналогично очистке данных — в каждом скрипте определялся метод, возвращающий дополнительные колонки признаков.

Группы дополнительных признаков:

1. Простейшие осмысленные признаки — ИМТ, пульсовое давление, усредненные значения давлений вида $\$inline\$\backslash frac\{ap\backslash \_hi+x*ap\backslash \_lo\}\{x+1\}\$inline\$$ для разных значений x. Также взяты приближенные формулы вычисления давлений по возрасту/весу и для каждого пациента вычислены ожидаемые давления (формулы вида $\$inline\$ap\backslash \_X\; =\; a\; +\; b\; *\; age\; +\; c\; *\; weight\$inline\$$). Вычислено на основе необработанных значений.
   
   
2. То же, что в п.1, но дополнительно на основе имеющихся давлений сделана попытка восстановить вес пациента. Для каждого из предсказанных таким образом признаков добавлена разница с «реальным» значением. Вычислено на основе необработанных значение.
3. Текстовое представление колонок необработанных данных, разбитое посимвольно — сначала с выравниванием слева, потом справа. Символы заменены их числовыми значениями (ord()). Где строка была слишком короткая и на все колонки не хватило — ставилось -1.
   
   
4. То же, что в п.3, но полученные колонки бинарно закодированы (one-hot encoding).
   
   
5. Данные из п.4, но пропущенные через PCA – тяжкое наследие недавно прошедшего на kaggle соревнования от mercedes, где модели с этим шаманством неплохо выглядели на паблике и печально на привате.
   
   
6. Для всех колонок необработанных исходных данных кроме возраста рассчитаны средние значения целевой колонки. Для этого сначала значения давлений, роста и веса поделил на 10 и округлил, получив из них таким образом категориальные признаки. Далее разбил данные на 10 фолдов и для каждого из них по 9 другим фолдам посчитал для каждой категории взвешенное среднее значение целевой колонки (болен-не болен). Там, где среднее вычислить было не из чего, проставил просто глобальное среднее.
   
   
7. То же что в п.6, но средние вычислялись еще и для признаков из п.2.
   
   
8. То же что в п.7, но в качестве исходных данных взяты очищенные по варианту №5.
   
   
9. То же что в п.7, но в качестве исходных данных взяты очищенные по варианту №3.
   
   
10. Необработанные данные кластеризованы методом k-means, количества кластеров выбраны произвольно — 2, 5, 10, 15, 25. Номер кластера для каждого из таких случаев — бинарно закодирован.
    
    
11. То же что в п.10, но использовались данные, очищенные по варианту №3.

Модели

Так как модели могут работать очень долго (десятки часов), падать с ошибками или прерываться специально — необходимо сохранять не только окончательные результаты, но и промежуточные данные. Для этого каждой модели задается базовое имя. Далее из имени модели и имени назначенного данным получается имя файла, где эти данные будут лежать. Все сохранение-загрузка идут через базовые методы модели, что и обеспечивает единообразное хранение промежуточных данных. В планах на будущее — держать данные не в файлах, а в какой-нибудь базе. Недостаток использованной реализации — можно забыть обновить имя при копировании модели и получить неопределенное состояние для данных исходной модели и ее копии.

Если у модели сохранены результаты ее вычислений — при обращении остается только их прочитать и вернуть вызывающей стороне. Если есть только промежуточные результаты — их тоже не приходится повторно считать. Это сильно экономит время, особенно когда речь идет о вычислениях на несколько часов.

Основное разделение сохраняемых моделями данных — по времени жизни этих данных. Для каждой такой группы данных предусмотрен свой базовый путь для сохранения. Всего таких 3 группы:

1. Временные, которые не будут использованы при следующих запусках, например лучшие веса нейросетей для отдельных фолдов;
2. Данные модели, которые при следующих запусках понадобятся — почти все остальное;
3. Глобально полезные данные, которые могут пригодиться нескольким моделям, например дополнительные признаки;

Интерфейс у всех моделей общий и позволяет не только запускать их отдельно как обычные скрипты, но и грузить как модули python. Если некоторой модели нужны результаты работы других моделей — она их грузит и выполняет. Как результат — описание каждой из моделей 2 уровня удалось свести к списку имен моделей, чьи результаты надо объединить, и признаку необходимости отбора признаков жадным алгоритмом.

Среди моделей были основанные на нейросетях, которые на выходе могли дать очень уверенные 0 или 1 или очень близкие к крайним значения. Так как в случае ошибки подобная самоуверенность штрафуется logloss-ом очень сильно — значения всех моделей при сохранении обрезались так, чтобы до 0 или 1 оставалось не менее 1e-5. Проще всего оказалось добавить np.clip(z, 1e-5, 1-1e-5) и забыть об этом. В результате обрезались данные всех моделей, но большинство и так давали результаты в диапазоне примерно 0.1-0.93.

hyperopt

Для подгонки параметров моделей пришлось задействовать hyperopt (подробности). Результаты улучшались, но при этом долгое время для особенно медленных моделей я выставлял количество попыток порядка 20. И дня за 2 до конца в статье нашел упоминание hyperopt bootstrapping — по умолчанию первые 20 запусков моделей производится со случайными параметрами, что можно увидеть в исходниках. Пришлось срочно пересчитывать часть моделей.

Модели 1 уровня

В общий код моделей 1 уровня попал выбор входных данных для каждой модели — обязательно один вариант очистки исходных данных и 0 или более групп признаков. Сборка данных и признаков в общий датасет реализована в общем для моделей коде. Это сократило код отдельных моделей до задания конкретных исходных данных и дополнительных признаков.

Сделать общий код для оптимизации не хватило времени, так что отдельные базовые модели 1 уровня по-прежнему сильно копируют друг друга. Всего их получилось 2 разновидности:

• нейросети (keras)
• деревья (XGBoost, LightGBM, rf, et)

Главное отличие использовавшихся моделей на основе нейросетей — отсутствие подгонки гиперпараметров. Для остальных моделей использовался hyperopt.

Нейросети

Сколько-нибудь серьезного подбора параметров нейросетей не делал, так что и результаты их были хуже чем у бустинга. Позже в чате видел упоминание устройства сети наподобие 64-64 с активацией leaky relu и дропаутом на 1-5 нейронов в каждом слое, которая давала относительно приличный результат.

У себя нейронные сети использовал примерно следующего вида:

• вход;
• несколько сотен нейронов (обычно 256);
• какая-нибудь нелинейность, дропаут (где был — брал значения порядка 0.7, потому что считал, что слишком много параметров и сеть переобучается); если при обучении модель расходилась в nan-ы — добавлял batch normalization — подробности здесь или здесь;
• сотня-другая нейронов (64-128);
• нелинейность;
• десяток-другой нейронов (16);
• нелинейность;
• 1 выходной нейрон с классической сигмоидой на выходе.

Подобное устройство перекочевало из предыдущих соревнований почти без изменений. По отдельности нейросети проявляли себя не очень, но были оставлены, чтобы их результаты использовать при расчете моделей 2 уровня.

Выбор функций активации для внутренних слоев очень простой — из доступного набора исключил все варианты сигмоиды (из-за близких к 0 градиентов у границ значений), «чистый» ReLU (из-за того, что нейрон, начавший выдавать на выходе 0, из обучения выпадает) и брал что-то из оставшихся. Изначально это было Parametric Relu, в последних моделях начал брать Scaled Exponential Linear Units. Каких-либо значимых различий от такой замены заметить не удалось.

Как и для других моделей, данные для нейросетей бились на фолды с помощью KFold из sklearn. Для обучения на каждом разбиении приходилось строить модель заново, так как я слишком поздно узнал как можно переинициализировать веса в слоях сети без ее пересоздания.

Сети обучались до тех пор, пока улучшалось качество предсказаний на части данных, выделенных для валидации. При этом сохранялись веса для сети каждый раз, когда результат на валидации улучшался. Для этого стандартные callback-и keras – на сохранение состояния сети с лучшими результатами на валидации, на раннее прекращение обучения если за заданное количество проходов по обучающей выборке результаты валидации не улучшались и на уменьшение learning rate если результаты не улучшались несколько проходов.

Получалось, что если обучение сети заходило в тупик (локальный минимум) и результаты не улучшались за несколько проходов по данным — уменьшался learning rate и если это не помогало — еще через сколько-то ходов обучение останавливалось. После обучения загружалось наилучшее состояние весов сети за все время обучения.

При этом я достаточно поздно заметил проблему, возникающую при попытке несколько раз использовать один и тот же набор экземпляров callback-ов при обучении нескольких сетей. В этом случае состояние callback-ов при начале обучения новой сети не сбрасывается автоматически в начальное. Как результат learning rate для каждой новой сети снижался все больше и больше, а лучшие результаты не сохранялись, если они не были лучшими чем вообще все ранее полученные на всех сетях, где использовались те же callback-и.

Модели, основанные на деревьях

Использовалось 2 варианта «древесных» моделей — основанные на bagging random forest и extra trees и 2 реализации градиентного бустинга — XGBoost и LightGBM. Оба варианта случайного леса плохо проявили себя и на кросс-валидации и на паблике, так что остались просто из-за того, что на них пришлось потратить много машинного времени и была надежда, что они пригодятся при объединении результатов моделей. LightGBM и XGBoost проявили себя гораздо лучше и большая часть предсказаний на первом уровне была получена от них.

Каждая из «древесных» моделей после подгонки параметров рассчитывалась для нескольких (обычно — 3) исходных состояний генератора случайных чисел. Все такие результаты сохранялись отдельно для использования моделями 2 уровня. Предсказания же моделей 1 уровня получалось из результата для последнего использованного состояния ГСЧ.

В LightGBM и XGBoost предусмотрена возможность останавливать обучение, если качество обучения не улучшается на валидации за заданное число итераций. За счет этого можно было просто разрешить им обучаться 10000 шагов и останавливаться когда результаты на валидации перестанут улучшаться. В результате при подборе параметров таких моделей количество деревьев подбирать не было необходимости. Для random forest и extra trees от sklearn такой возможности нет, так что пришлось переложить подбор количества деревьев на hyperopt, который, похоже, за недостаточное количество попыток с задачей полностью не справился. Можно было обучать их по одному шагу, каждый раз проверяя качество на валидации самому, но помешала лень.

Несколько сидов

Результат работы отдельных моделей сильно зависит от состояния генератора случайных чисел. Для избавления от такой зависимости для моделей 1 уровня рассчитывались результаты обучения с несколькими сидами. При этом для каждого сида результаты сохранялись отдельно. При этом после окончания конкурса выяснилось, что в качестве результата, сохраняемого отдельной моделью 1 уровня, использовался результат для последнего из сидов. Результаты для остальных все же сохранялись и использовались моделями 2 уровня.

Модели 2 уровня

За счет того, что каждая из моделей 1 уровня давала от 1 до 4 предсказаний, на 2 уровне данные содержали до 190 колонок. Исходные данные и признаки туда не попали — только уже предсказанные вероятности. Каждая из моделей 2 уровня объединяла некоторое подмножество из моделей 1 уровня (некоторые использовали еще результаты ранних моделей 2 уровня).

Все модели 2 уровня устроены примерно одинаково — загрузка по имени модулей моделей, получение результатов их работы, опциональный выбор какие колонки из загруженных использовать и подгонка параметров регрессии для объединения предсказаний.

Что не сработало — попытка добавить еще уровень и объединить предсказания всех моделей 2 уровня. Результат такого объединения на паблике оказался настолько плохим, что делать вторую попытку я даже не думал.

При выборе части предсказаний, комбинация которых давала бы наилучший результат, использовался «жадный» алгоритм — выбиралась одна лучшая колонка из доступных, после чего в цикле к ранее выбранным колонкам по одной добавлялись колонки из числа оставшихся. Добавление колонок продолжалось пока оно улучшало результаты модели 2 уровня на валидации. При этом для экономии времени в качестве модели для отбора использовался BayesianRidge, результаты которого уступали только Ridge с хорошо подогнанными параметрами. В результате такого отбора обычно оставалось порядка 20 колонок с данными.

Для итогового подсчета изначально прогонялся выбор через hyperopt всех доступных регрессоров, но более-менее прилично себя проявили только BayesianRidge и Ridge из sklearn, так что в итоге код выродился до определения насколько хорошо комбинирует модели BayesianRidge и попыток улучшить этот результат подгонкой параметров Ridge.

Валидация

Изначально валидировался на разбиении на 10 фолдов. При этом некоторые модели начали показывать cv под 0.534-0.535 при результатах на паблике за 0.543-0.544, а то и хуже. Ближе к концу конкурса, чтобы сблизить результаты при валидации и на паблике я пробовал увеличить разбиение до 30 фолдов. Выбор числа 30 основан на возможностях процессора — выбирал максимальное значение при котором вычисления одной модели требовали менее 10 часов.

При этом все равно часть моделей валидировались на уровне 0.535-0.536, что на фоне результатов на паблике порядка 0.543 на паблике вызывало сомнения в адекватности схемы валидации. Дня за 3 до окончания соревнования к разбиению на 30 фолдов пришлось добавить 30 случайных разбиений на 0.7 и 0.3 от тренировочных данных. Выбор именно 30 — та же причина что и для cv — возможности процессора. Все разбиения были зафиксированы по random_state. После этого лучший результат на валидации был порядка 0.537.

Это тоже было далеко от желаемого, но до конца оставалось несколько дней, которых хватало только чтобы дождаться пока досчитаются последние модели, так что пришлось на этом остановиться. В итоге выбирал 2 сабмита с результатами лучше 0.543 на паблике и 0.538 на валидации. Как выяснилось позже, из 12 таких сабмитов 7 давали 3 место, так что оставалось просто не промахнуться.