Спортивная аналитика

Прогнозируйте рейтинг чемпионата мира по фигурному катанию на основе результатов сезона

Часть 3: Многофакторные модели

Задний план

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

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

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

Эти скрытые оценки можно сложить вместе (аддитивная модель), перемножить (мультипликативная модель) или их комбинация (гибридная модель ):

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

Проблема

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

Из приведенного выше отчета по этим трем моделям за 10 лет (из 14), выбранных в обучающей выборке, все модели аппроксимируют результаты сезона лучше, чем базовая модель средних значений за сезон, о чем свидетельствует их более низкая среднеквадратичная ошибка. (RMSE) в среднем по сравнению с базовой моделью.

Тем не менее, в гораздо более важной метрике тау Кендалла эти модели, к сожалению, не обеспечивают какого-либо значительного увеличения этой метрики по сравнению с базовой моделью, поскольку 95% доверительные интервалы разницы в тау Кендалла между ними и базовой линией охватывают ноль для фигуристы, как мужчины, так и женщины.

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

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

Многофакторная модель

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

  • ŷ_event-skater: прогнозируемый результат данного фигуриста на данном мероприятии.
  • θ_baseline: исходный результат (постоянный для всех сезонов)
  • θ_event: скрытая оценка этого события (одинакова для фигуристов)
  • θ_skater: скрытая оценка этого фигуриста (постоянна в разных дисциплинах)

В приведенной выше формуле мы видим, что каждое событие представлено одной скрытой оценкой (θ_event), и каждый фигурист также представлен одной скрытой оценкой (θ_skater).

Однако что, если вместо одной скрытой оценки каждый фигурист представлен несколькими скрытыми факторами? Например, насколько он любит кататься по утрам или насколько она подготовлена ​​к соревнованиям на больших высотах. Точно так же каждое событие может быть представлено одними и теми же факторами: как рано утром происходит событие или на какой высоте.

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

  • ŷ_event-skater: прогнозируемый результат данного фигуриста на данном мероприятии.
  • f: латентный фактор, присутствующий в каждом виде спорта и в каждом фигуристе.
  • θ_event,f: оценка этого скрытого фактора для данного события
  • θ_skater,f: оценка этого скрытого фактора для данного фигуриста

Обнаружение скрытых факторов

Как и в предыдущих моделях, скрытые факторы в многофакторной модели могут быть найдены путем минимизации суммы квадратов разницы между истинными и прогнозируемыми результатами в течение сезона (которые будут называться остатками отсюда):

  • J: целевая функция модели. Продолжая пример с фигуристами-мужчинами в сезоне 2017 года, это будет функция 1 базового результата, всех n скрытых факторов для каждого из 9 видов спорта (выделены синим) и всех n скрытых факторов для каждого из 62 фигуристы (выделены желтым).
  • ŷ_e,s: прогнозируемый результат пары спортивного спортсмена и фигуриста на основе многофакторной модели.
  • y_e,s: истинный результат, зафиксированный для данной пары фигуристов и фигуристов в течение сезона.
  • θ_baseline: исходный балл
  • θ_e,f: оценка скрытого фактора f для события e
  • θ_s,f: оценка латентного фактора f для фигуриста s

Алгоритм градиентного спуска

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

Поиск градиентов

Как и в гибридной модели, градиент базовой оценки представляет собой просто сумму остатков для всех пар спортивного катания и конькобежцев в сезоне:

Затем давайте вычислим градиент некоторого фактора k данного события i (θ_ei,k):

  • При использовании этого градиента все квадраты остатков, не связанные с этим событием, исчезнут. Остались фигуристы, участвовавшие в этом мероприятии. Следовательно, будут учитываться только производные этих квадратов остатков.
  • Когда производная берется для каждого фигуриста s, который участвовал в этом событии, вступает в действие основное правило цепочки: внешняя производная возвращает остаток для этого фигуриста в этом событии. Затем мы берем внутреннюю производную от самого остатка, в результате чего исходный, истинный балл и баллы всех других скрытых факторов, которые не k, исчезают. Остается один балл латентного фактора k для этого фигуриста (θ_s,k), поскольку он был ранее умножен на соответствующий балл для коэффициента k этого самого события (θ_ei,k).
  • В результате для каждого фигуриста s, участвовавшего в этом соревновании, мы умножим внешнюю производную - остаток этого фигуриста в событии i - на внутреннюю производную - оценку латентного фактора k, который у него есть. Добавление этих продуктов для всех фигуристов, которые участвовали в этом мероприятии, наконец, даст нам градиент для фактора k события i.
  • Что интересно в этой формуле градиента, так это то, что невязка во внешней производной не зависит от фактора k, который нас интересует. Следовательно, этот остаток можно многократно использовать многократно, чтобы найти градиенты для всех других факторов этого события ( умножая его на разные внутренние производные). Это будет использовано позже, когда мы будем кодировать алгоритм градиентного спуска.

Точно так же градиент случайного фактора k данного фигуриста j - это просто фигуристский аналог градиента события выше:

Основные отличия этой формулы от предыдущей:

  1. Мы подводим итоги всех мероприятий, в которых e участвовал этот фигурист j
  2. Внутренние производные - это соответствующая оценка латентного фактора k для каждого из этих событий (θ_e,k)

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

Обновление скрытых оценок

После того, как градиенты для всех факторов в каждом виде спорта и для каждого фигуриста были рассчитаны, мы можем обновить скрытые оценки этих факторов как таковых:

  • Подобно градиентному спуску для гибридной модели, наличие скорости обучения α на этапах обновления будет контролировать скорость градиентного спуска и поможет лучше сходиться (если выбрана подходящая скорость обучения).
  • Обратите внимание, что шаги обновления для каждого фактора выполняются независимо от всех других факторов. Это также поможет нам реализовать код Python для градиентного спуска многофакторной модели.

Резюме

Учитывая формулы градиента и описанные выше правила обновления, градиентный спуск для многофакторной модели можно резюмировать как:

Из приведенного выше резюме видно, что:

  • После расчета остатка для каждой пары спортивного спортсмена (этап 2а) его можно использовать для расчета градиентов индивидуально для каждого латентного фактора k, независимо от каждого другого скрытого фактора (этап 2c-i). Другими словами, это как если бы существует только одна единственная скрытая оценка для каждого вида и каждого фигуриста, и их градиенты рассчитываются точно так же, как в однофакторной гибридной модели. После этого остатки повторно используются для вычисления градиентов для следующего скрытого фактора (с использованием того же метода).
  • Кроме того, шаги обновления для каждого фактора k также не зависят от всех других факторов (шаг 2c-ii). Следовательно, мы можем просто обновлять оценки для каждого скрытого фактора один за другим, пока не закончим одну итерацию шага 2.

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

Кодирование многофакторной модели

Давайте возьмем тот же пример игрушки, который мы использовали для гибридной модели, season_scores pandas DataFrame с 7 результатами за сезон у 4 фигуристов (MAJOROV, FERNANDEZ, GE, MURA) и 3 событиями в разных странах (CA, FR, RU):

Подобно гибридной модели, мы сначала преобразуем этот фрейм данных из длинного формата в формат сводной таблицы, где фигуристы представлены в виде строк, а события в виде столбцов и записей представляют собой оценки для каждой пары соревнований и фигуристов. Обратите внимание, что не все фигуристы участвуют во всех соревнованиях, поэтому в сводной таблице отсутствуют значения, которые представлены как NaN (не число).

season_pivot = pd.pivot_table(season_scores[['name', 'event', 'score']], values='score', index='name', columns='event')

Затем мы конвертируем сводную таблицу pandas в numpy-матрицу 4 × 3 true_scores, чтобы нам было проще управлять ею. Это фактически удалит все имена строк и столбцов, поэтому давайте сохраним имена фигуристов (строки) и имена событий (столбцы), чтобы мы могли идентифицировать их после запуска градиентного спуска.

skater_names = list(season_pivot.index) 
# ['Alexander, MAJOROV', 'Javier, FERNANDEZ', 'Misha, GE', 'Takahito, MURA']
event_names = list(season_pivot.columns) 
# ['CA', 'FR', 'RU']
true_scores = season_pivot.values
# array([[   nan,    nan, 192.14],
#        [   nan, 285.38, 292.98],
#        [226.07, 229.06,    nan],
#        [222.13, 248.42,    nan]])

Шаг 1: Инициализировать исходные, скрытые и скрытые результаты фигуристов.

Для базовой оценки, как и в гибридной модели, мы можем просто использовать метод random_sample объекта RandomState (с определенным начальным значением), чтобы вернуть случайное число от 0 до 1.

random_state = np.random.RandomState(seed=42)
baseline = random_state.random_sample() # 0.37

В этом примере игрушки предположим, что для каждого вида спорта и каждого фигуриста у нас есть 2 скрытых фактора. Следовательно, мы можем инициализировать все скрытые оценки для 4 фигуристов в виде матрицы 4 × 2. Точно так же мы можем инициализировать все скрытые оценки для 3 событий в виде матрицы 2 × 3.

skater_scores = random_state.random_sample((4, 2))
# array([[0.95071431, 0.73199394],
#        [0.59865848, 0.15601864],
#        [0.15599452, 0.05808361],
#        [0.86617615, 0.60111501]])
event_scores = random_state.random_sample((2, 3))
# array([[0.70807258, 0.02058449, 0.96990985],
#        [0.83244264, 0.21233911, 0.18182497]])

Шаг 2а: Рассчитайте остаток для каждой пары спортивного спортсмена и фигуриста.

Причина инициализации этих скрытых оценок в матричных формах заключается в том, что прогнозируемые оценки для каждой пары соревнований и фигуристов могут быть рассчитаны сразу путем умножения этих векторов (с использованием оператора умножения матриц @) и добавления базового показателя поверх этого:

predicted_scores = skater_scores @ event_scores + baseline
# array([[1.65705782, 0.54954103, 1.42974207],
#        [0.92831034, 0.41999206, 0.98355296],
#        [0.53334684, 0.39008461, 0.53640179],
#        [1.48824946, 0.52001014, 1.32395061]])

Чтобы проиллюстрировать это более наглядно, на диаграмме ниже показаны две скрытые матрицы (размером 4 × 2 и 2 × 3 соответственно), умноженные вместе с базовой линией, добавленной сверху, чтобы сформировать матрицу 4 × 3 прогнозируемых баллов для каждой пары спортивного спортсмена и фигуриста. .

Красным выделено то, как с помощью этой операции был вычислен прогнозируемый результат для пары спортивного спортсмена и фигуриста (CA-MAJOROV). Обратите внимание на то, как умножение двух матриц скрытых оценок автоматически выполняет суммирование по всем факторам пары спортивное соревнование-фигурист, взяв скалярное произведение соответствующей строки этого фигуриста и соответствующего столбца этого события в соответствующих матрицах скрытых оценок.

После того, как прогнозируемый балл для каждой пары спортивного спортсмена-фигуриста рассчитан, найти остатки так же просто, как вычесть сводную сводную таблицу истинных баллов из прогнозируемой матрицы баллов:

residuals = predicted_scores - true_scores

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

Шаг 2b: Рассчитайте градиент для базовой оценки и обновите ее.

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

Следовательно, для вычисления базового градиента (baseline_gradient) мы используем ту же функцию np.nansum для суммирования остатков существующих пар событий-фигуристов в причине, игнорируя несуществующие пары, которые являются NaN в матрице residuals. Наконец, базовая оценка обновляется путем вычитания ее градиента, умноженного на скорость обучения alpha.

alpha = 0.0005
baseline_gradient = np.nansum(residuals)
baseline = baseline - alpha * baseline_gradient

Шаг 2c: Для каждого фактора k вычислите его градиент в каждом виде и для каждого фигуриста и обновите каждый результат с помощью соответствующего градиента.

Ниже приведен код всего шага:

# For each factor k
for k in range(2):
    # Extract latent scores for that factor
    skater_scores_k = skater_scores[:, [k]]
    event_scores_k = event_scores[[k], :]
    # Step 2c-i: Find gradients of factor k for all skaters & events
    event_gradients_k = np.nansum(residuals * skater_scores_k, axis=0, keepdims=True)
    skater_gradients_k = np.nansum(residuals * event_scores_k, axis=1, keepdims=True)
    # Step 2c-ii: Update scores of factor k for all skaters & events
    event_scores[[k], :] = event_scores_k - alpha * event_gradients_k
    skater_scores[:, [k]] = skater_scores_k - alpha * skater_gradients_k

Как упоминалось ранее при выводе градиентного спуска, после вычисления остатков (шаг 2а) мы обрабатываем каждый фактор как полностью отдельный от всех других факторов. Следовательно, мы можем вычислять градиенты и обновлять оценки для каждого скрытого фактора один за другим с помощью цикла for k in range(2) (с 2 в качестве количества факторов в нашем примере игрушки): первый фактор будет иметь k=0, а второй k=1.

Давайте сначала рассмотрим первый фактор (k=0):

  • Во-первых, мы извлекаем скрытые оценки этого фактора. Это будет первая строка в скрытой матрице событий (event_scores[[0], :]) и первый столбец в скрытой матрице фигуристов (skater_scores[:, [0]]). Причина появления этих забавных квадратных скобок 0 в том, что мы хотим сохранить скрытые оценки для первого фактора как вектор-строку размера (1, 3) по соревнованиям и как вектор-столбец размера (4, 1) для фигуристов. Без этих скобок numpy просто свернет эти фрагменты в одномерные массивы.

  • С этими извлеченными баллами для первого фактора по всем дисциплинам (event_scores_k) и для всех фигуристов (skater_scores_k) проблема теперь сводится к однофакторной проблеме гибридной модели. Следовательно, расчет градиентов для всех оценок этого фактора точно такой же, как и в гибридной модели:
event_gradients_k = np.nansum(residuals * skater_scores_k, axis=0, keepdims=True)
skater_gradients_k = np.nansum(residuals * event_scores_k, axis=1, keepdims=True)

Я не буду вдаваться в подробности этого шага, поскольку он был полностью объяснен во второй части проекта, включая аргументы np.nansum: axis, который управляет направлением суммирования градиента (по фигуристам или по соревнованиям) и keepdims, который не позволяет numpy сворачивать градиенты до 1-D и сохранять их векторные формы строк / столбцов для более позднего шага обновления. Вкратце, я просто резюмирую эти операции на диаграммах ниже, сначала для расчета градиентов (для первого скрытого фактора) во всех событиях:

Затем для расчета градиентов (для первого латентного фактора) у всех фигуристов:

  • После того, как градиенты для первого скрытого фактора рассчитаны, их можно использовать для обновления оценок этого фактора для всех фигуристов и соревнований.
event_scores[[0], :] = event_scores_k - alpha * event_gradients_k
skater_scores[:, [0]] = skater_scores_k - alpha * skater_gradients_k

После обновления скрытых оценок для первого фактора цикл затем переходит ко второму фактору (k=1) и повторяет то же самое: вычислить градиенты для этого фактора для всех видов спорта и фигуристов (используя ту же матрицу residuals из шага 2а). , и обновите скрытые оценки, используя эти градиенты. Ниже приведены извлеченные однофакторные векторы, над которыми мы будем работать, когда k=1:

Наконец, весь шаг 2 повторяется много раз, переходя между вычислением остатков, поиском градиентов и обновлением скрытых оценок, пока RMSE после каждого цикла - np.sqrt(np.nanmean(residuals**2)) - не стабилизируется. После 1000 итераций RMSE многофакторной модели для этого игрушечного примера составляет поразительные 0,004 (по сравнению с 4,30 в однофакторной гибридной модели в части 2) с разницей от итерации к итерации 1e-5.

Фактически, при умножении матриц скрытых фигуристов и событий, чтобы получить окончательные прогнозируемые результаты в течение сезона, совсем не удивительно, что RMSE модели настолько низок, поскольку прогнозируемые оценки (слева) практически идентичны истинные оценки (справа):

Наконец, скрытые оценки для двух факторов для всех фигуристов могут быть получены в конце в pandas DataFrame, с ранее сохраненными именами фигуристов, добавленными обратно в качестве индекса строки, а коэффициенты (0 и 1) в виде столбцов:

pd.DataFrame(skater_scores, index=skater_names)
#                             0         1
# Alexander, MAJOROV  10.226857  3.122723
# Javier, FERNANDEZ   16.009246  4.594796
# Misha, GE           11.919684  7.280747
# Takahito, MURA      14.059611  3.020423

Алгоритм работает! Однако, когда я применил его к реальным данным, потребовалось довольно много времени для его запуска, особенно когда задействовано много факторов. Конечно, я также записывал кучу промежуточных значений по пути, но давайте посмотрим, как мы можем ускорить градиентный спуск, чтобы итерация модели стала менее болезненной.

Оптимизация градиентного спуска с помощью numpy broadcasting

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

Однако вместо обработки факторов по одному - с помощью дорогостоящего цикла for на шаге 2c - мы можем вычислять градиенты и обновлять оценки всех скрытых факторов одновременно. Это возможно, если снова использовать трансляции numpy, на этот раз с дополнительным измерением факторов (а не только фигуристами и событиями). Код для измененной версии шага 2c ниже:

# 2c-i: Calculate gradients for all factors
# Reshape matrices
reshaped_residuals = residuals[np.newaxis, :, :]
reshaped_event_scores = event_scores[:, np.newaxis, :]
reshaped_skater_scores = skater_scores.T[:, :, np.newaxis]
# Calculate gradients
event_gradients = np.nansum(residuals * reshaped_skater_scores, axis=1)
skater_gradients = np.nansum(residuals * reshaped_event_scores, axis=2).T
# 2c-ii: Update latent scores for all factors
event_scores = event_scores - alpha * event_gradients
skater_scores = skater_scores - alpha * skater_gradients

Что за странная индексация матриц с np.newaxis? Почему в наивной реализации keepdims=False вместо True? Прочтите ниже, чтобы узнать, как они работают.

Изменение формы матриц

Напомним, что у нас есть инициализированная матрица скрытых оценок для всех соревнований и фигуристов, которые представляют собой массивы event_scores размера (2, 3) и skater_scores размера (4, 3) соответственно. Из шага 2а у нас также есть матрица остатков (residuals) размера (4, 3), в которой каждая строка представляет фигуриста, а каждый столбец - событие. Все матрицы показаны на схеме ниже:

В скрытых матрицах выше эти два фактора склеены. Однако алгоритм градиентного спуска требует, чтобы они рассматривались отдельно. Следовательно, нам нужно найти способ разделить два фактора в каждой из скрытых матриц. В частности, для событий нам нужны 2 вектора-строки размером 1 × 3, а для фигуристов нам нужны 2 вектора-столбца размером 4 × 1, причем каждый вектор представляет один фактор.

Вот как видоизменяются эти скрытые матрицы:

reshaped_event_scores = event_scores[:, np.newaxis, :]

Напомним, скрытая матрица для событий (event_scores) представляет собой двумерный массив чисел размера (2, 3). Следовательно, event_scores[:, np.newaxis, :] по-прежнему сохранит исходные размеры - через : - при добавлении новой оси между исходными осями - через np.newaxis. В результате получается трехмерный массив размером (2, 1, 3). Вы можете думать об этом как о двух разных уровнях массива (1, 3), каждый из которых представляет один фактор для всех событий (выделен красным ниже):

reshaped_skater_scores = skater_scores.T[:, :, np.newaxis]

Мы можем сделать то же самое для скрытой матрицы skater, которая изначально является двумерным массивом numpy размера (4, 2). Однако мы хотим разделить это на 2 разных слоя массива (4, 1), каждый из которых представляет один фактор для всех фигуристов. Для этого мы используем skater_scores.T[:, :, np.newaxis]: поскольку skater_scores.T имеет размер (2, 4), трехмерный массив после добавления новой оси в конце будет иметь размер (2, 4, 1). Два уровня скрытых оценок фигуристов представлены ниже (красным):

reshaped_residuals = residuals[np.newaxis, :, :]

Наконец, матрица остатков, изначально представлявшая собой двумерный массив размера (4, 3), преобразована в residuals[np.newaxis, :, :], который представляет собой трехмерный массив размера (1, 4, 3) из-за новой оси, добавленной впереди. Вы можете думать о матрице остатков как о той же самой, что и раньше, но теперь существующей в одном плоском слое в 3-D.

Шаг 2c-i: рассчитайте градиенты для всех факторов

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

event_gradients = np.nansum(reshaped_residuals * reshaped_skater_scores, axis=1)
  • reshaped_residuals * reshaped_skater_scores: с измененными трехмерными массивами для остатков (1, 4, 3) и скрытых оценок фигуристов (2, 4, 1), когда мы умножаем их вместе, оба транслируются с помощью numpy. Остаточная матрица, ранее лежавшая в трехмерной плоскости, будет дублирована на другой слой, чтобы сформировать матрицу (2, 4, 3). Между тем, оценки скрытых фигуристов в каждом слое измененной матрицы скрытых фигуристов будут дублироваться три раза по горизонтали, чтобы также сформировать матрицу (2, 4, 3):

  • Как только трансляция завершена, обратите внимание, что каждый слой теперь снова становится просто поиском градиентов для одного фактора, независимо от других факторов. Следовательно, градиенты для каждого события для всех факторов будут просто знакомыми np.nansum(reshaped_residuals * reshaped_skater_scores, axis=1): аргумент axis=1 предназначен для суммирования фигуристов, которые участвовали в каждом мероприятии, в то время как аргумент keepdims=False сворачивает результат суммирования обратно в 2 -D массив размера (2, 3) для всех градиентов событий:

skater_gradients = np.nansum(reshaped_residuals * reshaped_event_scores, axis=2).T
  • То же самое можно сделать, чтобы найти градиенты скрытых оценок фигуристов сразу по всем факторам. Однако на этот раз преобразованные трехмерные массивы остаточной (1, 4, 3) и скрытой матрицы событий (2, 1, 3) умножаются вместе, и обе будут транслироваться: матрица reshaped_residuals снова будет дублирована на другую. слой, чтобы сформировать (2, 4, 3) матрицу. Между тем, оценки скрытых событий в каждом слое матрицы reshaped_event_scores будут дублироваться четыре раза по горизонтали, чтобы также сформировать матрицу (2, 4, 3):

  • После того, как эти две матрицы умножаются вместе, сумма берется по всем событиям, в которых участвовал каждый фигурист (axis=2), чтобы вернуть градиенты фигуриста, в то время как аргумент keepdims=False сворачивает градиенты фигуриста обратно в двумерный массив размера (2, 4 ). Однако мы знаем, что скрытая матрица фигуриста для матрицы имеет размер (4, 2), поэтому мы должны транспонировать матрицу градиента фигуриста (.T), прежде чем использовать ее для обновления оценок:

Шаг 2c-ii: обновите скрытые баллы по всем факторам

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

alpha = 0.0005
event_scores = event_scores - alpha * event_gradients
skater_scores = skater_scores - alpha * skater_gradients

Оптимизированный алгоритм градиентного спуска, использующий трансляцию numpy для вычисления градиента и обновления скрытых оценок для всех факторов одновременно, полностью представлен ниже (для 4 фигуристов, 3 видов соревнований и 2 скрытых факторов):

# Step 1: Initialize baseline score, and scores of all factors
random_state = np.random.RandomState(seed=42)
baseline = random_state.random_sample()
skater_scores = random_state.random_sample((4, 2))
event_scores = random_state.random_sample((2, 3))
# Step 2: Repeat until convergence
for i in range(1000):
    # a. Calculate residual for every event-skater pair
    predicted_scores = skater_scores @ event_scores + baseline
    residuals = predicted_scores - true_scores
    
    # b. Calculate baseline gradient and update baseline score
    baseline_gradient = np.nansum(residuals)
    baseline = baseline - alpha * baseline_gradient
    
    # c. Calculate gradient and update latent scores for all factors
    # Reshape matrices
    reshaped_residuals = residuals[np.newaxis, :, :]
    reshaped_event_scores = event_scores[:, np.newaxis, :]
    reshaped_skater_scores = skater_scores.T[:, :, np.newaxis]
    
    # c-i: Calculate gradients for all factors
    event_gradients = np.nansum(residuals * reshaped_skater_scores, axis=1)
    skater_gradients = np.nansum(residuals * reshaped_event_scores, axis=2).T
    
    # c-ii: Update latent scores for all factors
    event_scores = event_scores - alpha * event_gradients
    skater_scores = skater_scores - alpha * skater_gradients

Оптимизация градиентного спуска с умножением матриц

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

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

  • Чтобы найти все градиенты событий по всем факторам, нам просто нужно умножить (2, 4) транспонирование латентной матрицы фигуриста на (4, 3) остаточную матрицу. Обратите внимание, как ориентация матриц выравнивается, чтобы дать идеальную форму (2, 3) для матрицы градиента событий.

  • Точно так же, чтобы найти все градиенты фигуристов по всем факторам, мы умножаем остаток формы (4, 3) на (3, 2) транспонирование скрытой матрицы событий. Это дает правильную форму (4, 2) для матрицы градиента фигуриста.

Однако обратите внимание, что для того, чтобы это работало, операция умножения матриц numpy должна игнорировать значения NaN в остаточной матрице. К сожалению, это не так - np.nanmatmul не существует. Тем не менее, мы всегда можем заменить значения NaN нулями с помощью np.nan_to_num, а затем продолжить умножение матриц как обычно:

residuals = np.nan_to_num(residuals)
event_gradients = skater_scores.T @ residuals
skater_gradients = residuals @ event_scores.T

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

Сравнение трех реализаций

Давайте посмотрим, превосходит ли эти оптимизированные версии алгоритма градиентного спуска наивную поэтапную реализацию с использованием цикла for.

Сложность времени

Используя волшебную функцию %timeit записных книжек Jupyter, я измеряю время, необходимое каждой реализации для выполнения одной итерации градиентного спуска в этом примере игрушки (4 скейтера, 3 события). Я также варьирую количество факторов и измеряю соответствующее время для обеих реализаций, чтобы увидеть, насколько они эффективны по времени при увеличении количества факторов.

Результат можно увидеть на прилагаемом графике:

  • При наличии только одного скрытого фактора все 3 реализации занимают примерно одинаковое количество времени (100 мкс или около того) для выполнения одной итерации градиентного спуска.
  • По мере увеличения числа факторов время, необходимое наивной реализации для завершения одной итерации, также линейно увеличивается. Это имеет смысл, поскольку цикл for на шаге 2c по существу просто переходит от одного фактора к другому.
  • Напротив, транслируемая реализация намного быстрее, чем наивная, когда количество факторов увеличивается: всего с 10 факторами она уже почти в 10 раз быстрее, в то время как при большом количестве факторов она постоянно в 100 раз быстрее - обратите внимание на запуск время нанесено в десятичном логарифмическом масштабе на графике выше.
  • Еще лучше реализация с умножением матриц, которая примерно в два раза быстрее транслируемой. Причина, по которой эти оптимизированные версии работают намного быстрее, вероятно, проистекает из их эффективных np.nansum и np.matmul операций по сравнению с более медленным for циклом в наивных реализациях.

Космическая сложность

Однако одним из потенциальных недостатков этих оптимизированных методов является увеличенная пространственная сложность: например, от широковещательной передачи к более крупным трехмерным матрицам. Мы можем проверить это, построив график объема памяти, потребляемой при запуске каждой реализации одну за другой. Мы по-прежнему используем игрушечный пример с 4 фигуристами и 3 событиями, но количество факторов увеличено до более худшего сценария из 100 000 факторов.

Из приведенного выше графика использования памяти во времени (с использованием пакета memory_profiler) мы видим, что:

  • Из почти 6,7 секунд, в течение которых были запущены обе реализации, простая версия выполнялась в течение первых 6,5 секунд, транслируемая версия занимала следующие 0,15 секунды, а версия с матричным умножением занимала только последние 0,02 секунды. Это согласуется с эталонным временем, которое мы нашли ранее.
  • Наивная реализация потребляла всего 10 MiB памяти (1 MiB = 1.024 мегабайт) во время своего запуска. Интересно, что профилировщик строк показывает, что цикл for по 100 000 факторов практически не потребляет памяти! Вероятно, это связано с тем, что для каждого фактора мы работаем только с крошечной матрицей numpy 4 × 3 для вычисления его градиентов.
  • Напротив, транслируемая реализация потребляла около 30 МБ памяти, что не так плохо, как ожидалось, с учетом ее трехмерных матриц. Однако версия матричного умножения потребляла только половину этого - около 15 МБ.

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

Результат многофакторной модели

Когда мы применяем выбранный выше алгоритм градиентного спуска к знакомому примеру фигуристов-мужчин в сезоне 2017 года, мы снова можем отслеживать невязки и RMSE модели по итерациям алгоритма. Они отображаются на анимированной панели инструментов ниже для первых 150 итераций градиентного спуска:

  • По мере того, как выполняется больше итераций градиентного спуска, остаток для каждой пары соревнований и фигуристов значительно уменьшается, о чем свидетельствует море серого в конце на остаточной тепловой карте (по сравнению с той однофакторной модели частично 2 проекта).
  • Это подтверждается очень низким значением RMSE модели 2,41 после 150 итераций. Это более чем в три раза меньше, чем у предыдущих моделей (8,84 для среднего, 8,88 для мультипликативного, 8,86 для гибридного), и более чем в четыре раза меньше, чем у базовой модели средних за сезон (10,3). Другими словами, многофакторная модель делает феноменальную работу по приближению результатов сезона.
  • Фактически, если мы позволим алгоритму действительно сходиться, он достигнет фактически несуществующей RMSE 0,03 после 1000 итераций с разницей от итерации к итерации 1e-4. Однако это неудивительно, поскольку это гораздо более сложная модель, чем предыдущие, с гораздо большим количеством параметров, которые необходимо настроить, чтобы наилучшим образом приблизиться к результатам сезона.

Рейтинг фигуристов с использованием индивидуальных факторов

Несмотря на то, что многофакторные модели очень хорошо аппроксимируют результаты сезона, у них есть серьезный недостаток: вместо единой оценки, которую можно использовать для ранжирования фигуристов (от самого высокого до самого низкого), у каждого фигуриста теперь есть несколько скрытых оценок, по одному для каждого. фактор. Возникает вопрос: по какому счету мы должны оценивать фигуристов?

С 5 скрытыми факторами, когда мы ранжируем фигуристов, используя оценки каждого фактора в отдельности, мы видим, что:

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

  • Это также можно увидеть, построив тау Кендалла, которые измеряют попарные корреляции между любыми двумя рейтингами от разных факторов (см. Прилагаемую тепловую карту): наивысшая корреляция тау Кендалла между прогнозируемыми рейтингами любых двух факторов составляет ничтожные 0,37 между факторами 1 и 5. .
  • Хуже того, ни один из рейтингов, созданных отдельно по каждому фактору, не подходит для прогнозирования окончательного рейтинга чемпионата мира, с самым высоким тау Кендалла, равным всего 0,45 от 1-го фактора. Это намного ниже, чем у базовой модели средних значений за сезон (0,695), не говоря уже о предыдущих моделях, которые мы построили (0,732 для аддитивной и мультипликативной модели и 0,724 для гибридной модели).

Таким образом, ни один из факторов, взятых по отдельности, не может дать фигуристам достойного звания на чемпионат мира! Вместо этого мы должны найти способ объединить скрытые оценки по различным факторам, чтобы более точно ранжировать фигуристов. В следующей части проекта я расскажу, как это сделать.

Ресурс

В многофакторной модели скрытые оценки фигуриста хранятся в матрице skater_scores размера количество фигуристов × количество факторов, в то время как скрытые оценки событий хранятся в матрице event_skaters размера количество факторов × количество событий. В результате многофакторная модель может рассматриваться как факторизация матрицы сезонных оценок, состоящей из количества фигуристов × количество соревнований, в эти две скрытые матрицы, так что результат их работы над базовым баллом - skater_scores @ event_scores + baseline - хорошо аппроксимируется исходным. матрица результатов сезона.

Фактически, алгоритм градиентного спуска, используемый в этой части проекта, почти идентичен алгоритму известной матричной факторизации FunkSVD, используемой в задаче Netflix, чтобы рекомендовать фильмы пользователям: путем факторизации матрицы рейтингов на пользовательские и кинофильмы. конкретных скрытых матриц, мы можем использовать произведение этих скрытых оценок для прогнозирования рейтинга фильмов, которые пользователь еще не смотрел. Таким образом, многофакторная модель - это, по сути, FunkSVD с пользователем = фигурист, фильм = событие и рейтинг = результат сезона.

Однако основные различия между этими двумя случаями заключаются в следующем:

  • Для этой задачи нас не интересует прогнозирование того, сколько фигурист потенциально может забить в событии прошлого сезона, а скорее, насколько высоко он будет занимать место в чемпионате мира в будущем.
  • Для простоты и небольших обучающих данных градиентный спуск, используемый в этом проекте, представляет собой пакетный градиентный спуск, в котором все обучающие данные использовались для обновления всех скрытых оценок на каждой итерации. Напротив, FunkSVD часто выполняется методом стохастического градиентного спуска: каждая итерация выбирает только одну оценку за раз, чтобы обновить скрытые оценки пользователя и фильма, связанные с этим рейтингом. Это может быть даже быстрее, чем оптимизированная реализация, представленная здесь.
  • В отличие от многофакторной модели, исходная формула FunkSVD не включает базовый член. Напротив, другие варианты формулы включают не только универсальный базовый термин, но также и базовые термины, которые относятся к конкретным видам спорта и фигуристам. Вы можете ознакомиться с этой статьей для очень удобочитаемой презентации FunkSVD и других методов матричной факторизации, распространенных в рекомендательных системах.