Введение
Метод логистической регрессии один из простейших методов классификации в машинном обучении. Тем не менее он ряд сильных сторон и поэтому его не стоит сбрасывать со счетов. Во-первых логистическая регрессия выдаёт не только предполагаемы класс объекта, но и оценку вероятности принадлежности объекта к тому или иному классу. Причём оценка эта выдаётся в виде гладкой функции. Во-вторых этот метод очень просто реализовать и понять: вот, например, реализация из документации Apache Spark:
В третьих, натренированная модель представляет собой вектор весов, которые просто интерпретировать.
Когда объясняют логистическую регрессию, чаще всего ограничиваются бинарной классификацией (с двумя классами). Например, эта тема прекрасно раскрыта в курсе Andrew Ng на Coursera. Однако нередко нужно различать более чем 2 класса. В этом случае, например, в курсе Andrew Ng предлагается стратегия one-vs-all. То есть для различения n классов надо натренировать n бинарных классификаторов по одному на каждый класс объектов. Поэтому мне кажется интересным с самого начала рассмотреть именно вариант с произвольным числом классов.
Теория (нестрогое изложение)
Пусть наши наблюдения описываются некоторым количеством признков ${\bf f} = \{f_1...f_m\}$. Обозначим буквой $Y$ класс объекта ($Y \in 1...k$). Будем моделировать вероятность с помощью экспоненциальной функции от взвешенной суммы признаков
$$P(Y=l) \propto \exp(w^{(l)}_1 f_1 + w^{(l)}_2 f_2 + ... + w^{(l)}_m f_m ).$$
То есть для каждого класса $l$ имеется свой набор весов $\bf w^{(l)}$. Вероятность должна быть нормирована на 1: $$P(Y=l) = \frac{\exp(w^{(l)}_1 f_1 + w^{(l)}_2 f_2 + ... + w^{(l)}_m f_m )}{\sum\limits_{j=1}^k \exp(w^{(j)}_1 f_1 + w^{(j)}_2 f_2 + ... + w^{(j)}_m f_m )}.$$ Задача обучения сводится к тому, чтобы подобрать такие веса $\bf w^{(l)}$, при которых будет достигаться оптимальное разделение классов. Для этого берётся "обучающее множество": некоторое количество экземпляров данных, на которых известны $Y$ и $\bf f$. Затем изменяя значения $\bf w^{(l)}$ максимизируют функцию правдоподобия на этих данных (на практике максимизируют логарифм правдоподобия, так как в этом случае вычисления проще).
Обозначим ${\bf W} = \left( {\bf w}^{(1)}, {\bf w}^{(2)}, ... , {\bf w}^{(k)}\right)$, ${\bf F}=\left({\bf f}, {\bf f} ... {\bf f} \right)$ (всего $k$ повторений ${\bf f}$), ${\bf F}\left(Y=j\right)=\left({\bf 0}, {\bf 0} ... {\bf f} ...{\bf 0} \right)$ (здесь положены равными нулю все ${\bf f}$, кроме той, что стоит на $j$-й позиции). Такая запись позволяет представить веса модели и признаки в виде векторов рамерностью $mk$, что на мой взглят облегчает объяснения.
$$ P(Y=l|{\bf F; W}) = \frac{\exp\left({\bf W \cdot F}\left(Y=l\right)\right)}{\sum\limits_{j=1}^k\exp\left({\bf W \cdot F}\left(Y=j\right)\right)} $$ Итак, запишем логарифм функции правдоподобия: $$ L\left({\bf W}\right) = \sum\limits_{i=1}^n \ln P\left(Y=Y_i|{\bf F}_i; {\bf W}\right). $$ Здесь $i$ обозначает номер примера в тренировочном множестве, $Y_i$ -- его класс, ${\bf F}_i$ -- его набор признаков. Продиференцировав вражение для нормированной верятности $P\left(Y|{\bf F; W}\right)$ получим: $$ \frac{\partial L\left({\bf W}\right)}{\partial {\bf W}} = \sum\limits_{i=1}^n{\bf F}_i\left(Y=Y_i\right) - \sum\limits_{i=1}^n\sum\limits_{j=1}^kP\left(Y=j|{\bf F}_i; {\bf W}\right) $$ Зная градиент функции можно искать её максимум, например методом градиентного подъёма. Ещё одно замечательное свойство многоклассовой логистической регрессии -- выпуклость возникающей задачи оптимизации (см Collins, стр 16). Значит мы сразу найдём глобальный максимум.
Для того чтобы контролировать переобучение к функции правдоподобия добаволяется дополнительный член: штраф за слишком большие веса модели ${\bf W}$. Мы будем рассматривать только квадратичный штраф: $$ L\left({\bf W}\right) = \sum\limits_{i=1}^n \ln P\left(Y=Y_i|{\bf F}_i; {\bf W}\right) - \frac{\lambda}{2}{\bf W^T}\cdot{\bf W}, $$ соответственно, в производной появится дополнительный член $-\lambda {\bf W}$.
$$P(Y=l) \propto \exp(w^{(l)}_1 f_1 + w^{(l)}_2 f_2 + ... + w^{(l)}_m f_m ).$$
То есть для каждого класса $l$ имеется свой набор весов $\bf w^{(l)}$. Вероятность должна быть нормирована на 1: $$P(Y=l) = \frac{\exp(w^{(l)}_1 f_1 + w^{(l)}_2 f_2 + ... + w^{(l)}_m f_m )}{\sum\limits_{j=1}^k \exp(w^{(j)}_1 f_1 + w^{(j)}_2 f_2 + ... + w^{(j)}_m f_m )}.$$ Задача обучения сводится к тому, чтобы подобрать такие веса $\bf w^{(l)}$, при которых будет достигаться оптимальное разделение классов. Для этого берётся "обучающее множество": некоторое количество экземпляров данных, на которых известны $Y$ и $\bf f$. Затем изменяя значения $\bf w^{(l)}$ максимизируют функцию правдоподобия на этих данных (на практике максимизируют логарифм правдоподобия, так как в этом случае вычисления проще).
Обозначим ${\bf W} = \left( {\bf w}^{(1)}, {\bf w}^{(2)}, ... , {\bf w}^{(k)}\right)$, ${\bf F}=\left({\bf f}, {\bf f} ... {\bf f} \right)$ (всего $k$ повторений ${\bf f}$), ${\bf F}\left(Y=j\right)=\left({\bf 0}, {\bf 0} ... {\bf f} ...{\bf 0} \right)$ (здесь положены равными нулю все ${\bf f}$, кроме той, что стоит на $j$-й позиции). Такая запись позволяет представить веса модели и признаки в виде векторов рамерностью $mk$, что на мой взглят облегчает объяснения.
$$ P(Y=l|{\bf F; W}) = \frac{\exp\left({\bf W \cdot F}\left(Y=l\right)\right)}{\sum\limits_{j=1}^k\exp\left({\bf W \cdot F}\left(Y=j\right)\right)} $$ Итак, запишем логарифм функции правдоподобия: $$ L\left({\bf W}\right) = \sum\limits_{i=1}^n \ln P\left(Y=Y_i|{\bf F}_i; {\bf W}\right). $$ Здесь $i$ обозначает номер примера в тренировочном множестве, $Y_i$ -- его класс, ${\bf F}_i$ -- его набор признаков. Продиференцировав вражение для нормированной верятности $P\left(Y|{\bf F; W}\right)$ получим: $$ \frac{\partial L\left({\bf W}\right)}{\partial {\bf W}} = \sum\limits_{i=1}^n{\bf F}_i\left(Y=Y_i\right) - \sum\limits_{i=1}^n\sum\limits_{j=1}^kP\left(Y=j|{\bf F}_i; {\bf W}\right) $$ Зная градиент функции можно искать её максимум, например методом градиентного подъёма. Ещё одно замечательное свойство многоклассовой логистической регрессии -- выпуклость возникающей задачи оптимизации (см Collins, стр 16). Значит мы сразу найдём глобальный максимум.
Для того чтобы контролировать переобучение к функции правдоподобия добаволяется дополнительный член: штраф за слишком большие веса модели ${\bf W}$. Мы будем рассматривать только квадратичный штраф: $$ L\left({\bf W}\right) = \sum\limits_{i=1}^n \ln P\left(Y=Y_i|{\bf F}_i; {\bf W}\right) - \frac{\lambda}{2}{\bf W^T}\cdot{\bf W}, $$ соответственно, в производной появится дополнительный член $-\lambda {\bf W}$.
Реализация
Для экспериментов я реализовал многоклассовый классификатор на C++ в виде консольной программы. В приоритете была простота, поэтому для нахождения минимума целевой функции при обучении используется обычный градиентный спуск с постоянным шагом. Хотя на практике применяют методы получше. Из этих же соображений для представления векторов я использовал класс valarray из стандартной библиотеки C++. Мне больше нравится библиотека Eigen, но пока что я не хочу связываться с подключением внешних библиотек.
Эксперименты с обучением на реальных данных
Чтобы попробовать метод и его реализацию на реальных данных, я взял два датасета из репозитория UCI. Один из них - известнейший датасет Iris (Ирисы Фишера), второй - датасет Wine . В первом предлагается определить тип ириса, используя в качестве параметров геометрические размеры элементов цветка (всего 4 параметра). Во втором наборе днных предлагается определить происхождение вина (категориальная переменная, которая принимает значения 1, 2, 3) по 13 параметрам, полученным в результате химического анализа вина. Надо сказать, что оба набора данных не из тех, что "бросают вызов". Практически любой классификатор может показать на них неплохие результаты без какой-то особой обработки данных напильником.
Итак, сначала убедмся, что классификатор может обучиться разделять объекты в тренировочном наборе данных. Для датасета Iris будем использовать данные "как есть", слегка поменяв формат файла (заменим запятные в качестве разделителя полей на символы табуляции) и заменив названия классов на числовые значения: "Iris setosa" => 0, "Iris virginica" => 1, "Iris versicolor" => 2. В *nix или cygwin окружении можно использовавть такую команду (в одну строчку)
>sed -e 's/,/\t/g;s/Iris-setosa/0/;s/Iris-versicolor/1/;s/Iris-virginica/2/' ./iris.data > ./iris.tsv
или произвести замену в своём любимом текстовом редакторе.
Тренируем модель
>./mlr train --verbose --modelfile irismodel.txt --datafile iris.tsv --featureincolumns 0,1,2,3 --labelincolumns 4 --numclasses 3 --learnrate 0.2
Теперь посмотрим на матрицу неточностей (confusion matrix) на тех же самых тренировочных данных:
>./mlr eval --modelfile irismodel.txt --datafile iris.tsv --featureincolumns 0,1,2,3 --labelincolumns 4
Получим такой вывод
50 0 0
0 49 1
0 1 49
В этой матрице 3x3 число в i-й строке и j-й колонке соответствует количеству объектов класса i которым классификатор назначил класс j. "Идеальный" классификатор должен назначать всем объектам класса i только лишь класс i и никакой другой, поэтому матрица неточностей должны быть диагональной. Недиагональные элементы соответствуют ошибкам классификации.
Мы видим, что классификатор ошибся лишь два раза . Можно сделать вывод что в данном случае классификатор хорошо настроился на разделение тренировочных данных. Несмотря на линейность мы можем хорошо разделять классы ирисов в пространстве признаков. Правда, пока что рано говорить, что получился хороший классификатор: возможно мы просто переобучились, и на новых данных ошибка будет велика. Поэтому проверим теперь наш классификатор с помощью процедуры кросс-валидации. Разделим данные в пропорции 70/30. Большую часть будем использовать для тренировки, а меньшую - для проверки качества.
>sort -R ./iris.tsv > iris.shuffled.tsv
>head -n 105 iris.shuffled.tsv > ./iris.train.tsv
>tail -n 45 iris.shuffled.tsv > ./iris.test.tsv
>./mlr train --verbose --modelfile irismodel.txt --datafile iris.train.tsv --featureincolumns 0,1,2,3 --labelincolumns 4 --numclasses 3 --learnrate 0.2
>
Теперь матрица неточностей выглядит так:
16 0 0
0 15 0
0 0 14
Можно сказать, что мы научили классификатор различать ирисы.
Теперь проделаем то же самое с датасетом Wine. Сначала как и прежде заменим запятые на символы табуляции в файле данных. Кроме того класс объекта представлен числами 1,2,3. У меня всё zero-based, поэтому надо замнить в первой колонке (у неё тоже индекс 0 а не 1): 1=>0, 2=>1, 3=>2.
>sed -e 's/,/\t/g;s/^1/0/;s/^2/1/;s/^3/2/' wine.dat > wine.tsv
Взглянув на колонки признаков, можно заметить что там встречаются значения разные по порядку величины. Например ... Поэтому неплохо было бы провети сначала нормализацию данных: с помощью сдвига и масштабирования привести их значения к 0, а разброс к некоторой одинаковой величине, например, 1. Под "разбросом" можно понимать либо разницу максимального и минимального значения, либо среднеквадратичное отклонение. Для нормализации данных я написал скрипт на Perl. Скрипт можно использовать так:
>perl normalize.pl stddev --in wine.tsv --fields 1,2,3,4,5,6,7,8,9,10,11,12,13 > wine.normalized.tsv
Тренируем:
>./mlr train --verbose --modelfile wine.model.txt --datafile wine.normalized.tsv --featureincolumns 1,2,3,4,5,6,7,8,9,10,11,12,13 --labelincolumns 0 --numclasses 3 --learnrate 0.2
Итак, сначала убедмся, что классификатор может обучиться разделять объекты в тренировочном наборе данных. Для датасета Iris будем использовать данные "как есть", слегка поменяв формат файла (заменим запятные в качестве разделителя полей на символы табуляции) и заменив названия классов на числовые значения: "Iris setosa" => 0, "Iris virginica" => 1, "Iris versicolor" => 2. В *nix или cygwin окружении можно использовавть такую команду (в одну строчку)
>sed -e 's/,/\t/g;s/Iris-setosa/0/;s/Iris-versicolor/1/;s/Iris-virginica/2/' ./iris.data > ./iris.tsv
или произвести замену в своём любимом текстовом редакторе.
Тренируем модель
>./mlr train --verbose --modelfile irismodel.txt --datafile iris.tsv --featureincolumns 0,1,2,3 --labelincolumns 4 --numclasses 3 --learnrate 0.2
Теперь посмотрим на матрицу неточностей (confusion matrix) на тех же самых тренировочных данных:
>./mlr eval --modelfile irismodel.txt --datafile iris.tsv --featureincolumns 0,1,2,3 --labelincolumns 4
Получим такой вывод
50 0 0
0 49 1
0 1 49
В этой матрице 3x3 число в i-й строке и j-й колонке соответствует количеству объектов класса i которым классификатор назначил класс j. "Идеальный" классификатор должен назначать всем объектам класса i только лишь класс i и никакой другой, поэтому матрица неточностей должны быть диагональной. Недиагональные элементы соответствуют ошибкам классификации.
Мы видим, что классификатор ошибся лишь два раза . Можно сделать вывод что в данном случае классификатор хорошо настроился на разделение тренировочных данных. Несмотря на линейность мы можем хорошо разделять классы ирисов в пространстве признаков. Правда, пока что рано говорить, что получился хороший классификатор: возможно мы просто переобучились, и на новых данных ошибка будет велика. Поэтому проверим теперь наш классификатор с помощью процедуры кросс-валидации. Разделим данные в пропорции 70/30. Большую часть будем использовать для тренировки, а меньшую - для проверки качества.
>sort -R ./iris.tsv > iris.shuffled.tsv
>head -n 105 iris.shuffled.tsv > ./iris.train.tsv
>tail -n 45 iris.shuffled.tsv > ./iris.test.tsv
>./mlr train --verbose --modelfile irismodel.txt --datafile iris.train.tsv --featureincolumns 0,1,2,3 --labelincolumns 4 --numclasses 3 --learnrate 0.2
>
Теперь матрица неточностей выглядит так:
16 0 0
0 15 0
0 0 14
Можно сказать, что мы научили классификатор различать ирисы.
Теперь проделаем то же самое с датасетом Wine. Сначала как и прежде заменим запятые на символы табуляции в файле данных. Кроме того класс объекта представлен числами 1,2,3. У меня всё zero-based, поэтому надо замнить в первой колонке (у неё тоже индекс 0 а не 1): 1=>0, 2=>1, 3=>2.
>sed -e 's/,/\t/g;s/^1/0/;s/^2/1/;s/^3/2/' wine.dat > wine.tsv
Взглянув на колонки признаков, можно заметить что там встречаются значения разные по порядку величины. Например ... Поэтому неплохо было бы провети сначала нормализацию данных: с помощью сдвига и масштабирования привести их значения к 0, а разброс к некоторой одинаковой величине, например, 1. Под "разбросом" можно понимать либо разницу максимального и минимального значения, либо среднеквадратичное отклонение. Для нормализации данных я написал скрипт на Perl. Скрипт можно использовать так:
>perl normalize.pl stddev --in wine.tsv --fields 1,2,3,4,5,6,7,8,9,10,11,12,13 > wine.normalized.tsv
Тренируем:
>./mlr train --verbose --modelfile wine.model.txt --datafile wine.normalized.tsv --featureincolumns 1,2,3,4,5,6,7,8,9,10,11,12,13 --labelincolumns 0 --numclasses 3 --learnrate 0.2
Смотрим матрицу неточностей:
59 0 0
0 71 0
0 0 48
20 0 0
0 18 0
0 0 15
Мы рассмотрели метод многоклассовой логистической регрессии. Применили его к реальным наборам данных и оценили способность линейного классификатора, натренированного по методу многоклассовой логистической регрессии разделять классы в этих наборах данных. Можно сделать вывод, что простые алгоритмы, такие, как линейный классификатор иногда работают очень неплохо. В продолжении я собираюсь обсудить оптимизацию приведённого кода, рассмотреть обобщения метода и попробовать применить его к наборам данных посложнее.
>./mlr eval --modelfile wine.model.txt --datafile wine.normalized.tsv --featureincolumns 1,2,3,4,5,6,7,8,9,10,11,12,13 --labelincolumns 0
59 0 0
0 71 0
0 0 48
Видим, что классификатор научился идеально разделять объекты в тренировочных данных. Можно сказать, что чем больше рамерность пространства признаков тем больше шансов, что данные можно будет разделить линейным классификатором.
Для чистоты эксперимента опять разделим данные на тестовую и тренировочную подвыборки.
>sort -R ./wine.normalized.tsv > wine.shuffled.tsv
>head -n 124 ./wine.shuffled.tsv > ./wine.train.tsv
>tail -n 53 ./wine.shuffled.tsv > ./wine.test.tsv
Посмотрим, что получилось.
>./mlr train --verbose --modelfile wine.model.txt --datafile wine.train.tsv --featureincolumns 1,2,3,4,5,6,7,8,9,10,11,12,13 --labelincolumns 0 --numclasses 3 --learnrate 0.2
>./mlr eval --modelfile wine.model.txt --datafile wine.test.tsv --featureincolumns 1,2,3,4,5,6,7,8,9,10,11,12,13 --labelincolumns 0
0 18 0
0 0 15
Внедиагональных элементов вообще не оказалось. Таким образом можно сделать вывод что классификатор научился уверенно определять происхождение вина по результатам химического анализа.
Заключение
Мы рассмотрели метод многоклассовой логистической регрессии. Применили его к реальным наборам данных и оценили способность линейного классификатора, натренированного по методу многоклассовой логистической регрессии разделять классы в этих наборах данных. Можно сделать вывод, что простые алгоритмы, такие, как линейный классификатор иногда работают очень неплохо. В продолжении я собираюсь обсудить оптимизацию приведённого кода, рассмотреть обобщения метода и попробовать применить его к наборам данных посложнее.
Ссылки
Простая реализация логистической, написанная мной специально для этого поста
Статья "Multinomial logistic regression" в англоязычной википедии
Статья "Multinomial logistic regression" в англоязычной википедии
