Этот файл содержит описание и результаты всех проделанных нами работ. Работы подготовили: Орлова Татьяна и Юлия Ивлева, 401-И

• Метрические алгоритмы классификации
  • K ближайших соседей (KNN)
  • K взвешенных ближайших соседей (kwKNN)
  • Парзеновское окно (PW)
  • Потенциальные функции (PF)
  • STOLP
• Байесовские классификаторы
  • Нормальный дискриминантный анализ
  • Подстановочный алгоритм (Plug-in)
  • Линейный дискриминант Фишера (ЛДФ)
• Линейные классификаторы
  • Метод стохастического градиента
  • Адаптивный линейный элемент (adaline)
  • Персептрон Розенблатта
  • Логистическая регрессия
  • Вывод по линейным классификаторам

Гипотеза компактности утверждает, что схожим объектам соответствуют схожие ответы.

Метрические алгоритмы классификации основаны на гипотезе компактности, которая говорит о том, что схожим объектам соответствуют схожие ответы.

Функцию, определяющую "схожесть" объектов называют мерой близости. Эта функция определяется следующим образом:

(функция расстояния)

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

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

Для оценки близости объекта u к классу y алгоритм использует следующую функцию: , где i обозначает порядок соседа по расстоянию к точке u.

Программно, функция веса реализуется следующим образом:

mc.KNN.w = function(i, k) +(i <= k)

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

Программно алгоритм реализуется следующим образом:

mc.KNN = function(sortedDistances, k) {
    orderedDistances = 1:length(sortedDistances)
    names(orderedDistances) = names(sortedDistances)

    weights = mc.KNN.w(orderedDistances, k)
    weightsByClass = sapply(unique(names(weights)), mc.sumByClass, weights)

    bestClass = names(which.max(weightsByClass))
}

где sortedDistances – массив расстояний от неизвестной u до всех точек обучающей выборки. names(sortedDistances) – наименование классов для каждой точки выборки.

Реализация алгоритма доступна по ссылке

Протестируем алгоритм на выборке ирисов Фишера:

Алгоритм дает маленькую погрешность (5 ошибок) при небольшом k, однако при росте k, начиная с k > 95, ошибка стремительно растет. Это объясняется тем, что функция оценки близости никак не учитывает порядок элементов, а учитывает лишь их наличие. Из-за этого те объекты, что находятся очень далеко от точки u влияют на классификацию с такой же силой, что и объекты, находящиеся в непосредственной близости.

Плюсы:

• прост в реализации
• неплохие результаты при правильно подобраном k

Минусы:

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

Для оценки близости объекта u к классу y алгоритм использует следующую функцию: , где i обозначает порядок соседа по расстоянию к точке u, а — строго убывающая функция веса. Последней алгоритм взвешенного KNN отличается от KNN.

Мы же будем применять следующую функцию веса:

Программно функция веса реализуется следующим образом:

mc.kwKNN.w = function(i, k) +(i <= k) * (k + 1 - i) / k

Сам же алгоритм ничем, кроме функции веса, от KNN не отличается. Подставив mc.kwKNN.w вместо mc.KNN.w получим необходимый алгоритм.

Реализация алгоритма доступна по ссылке

Протестируем алгоритм на выборке ирисов Фишера:

Алгоритм демонстрирует маленькую погрешность при любом k, отклонение максимальной ошибки к минимальной всего лишь в двух объектах.

Алгоритм взвешенных k соседей качественно отличается от обычного KNN тем, что учитывает порядок объектов при классификации. Тем самым, ближние к точке u объекты будут влиять на нее гораздо сильнее, чем дальние. Однако из-за этой особенности, при большом k дальние точки быстро обесцениваются, поэтому функцию следует выбирать осторожно.

Плюсы:

• прост в реализации
• неплохие результаты при любом k

Минусы:

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

Для оценки близости объекта u к классу y алгоритм использует следующую функцию:

, где — функция ядра.

Чаще всего применяются 5 типов ядер:

• Прямоугольное
• Треугольное
• Квартическое
• Епанечниково
• Гауссовское (нормальное распределение)

Программная реализация ядер:

mc.kernel.R = function(r) 0.5 * (abs(r) <= 1) #прямоугольное
mc.kernel.T = function(r)(1 - abs(r)) * (abs(r) <= 1) #треугольное
mc.kernel.Q = function(r)(15 / 16) * (1 - r ^ 2) ^ 2 * (abs(r) <= 1) #квартическое
mc.kernel.E = function(r)(3 / 4) * (1 - r ^ 2) * (abs(r) <= 1) #епанечниково
mc.kernel.G = function(r) dnorm(r) #гауссовское

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

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

Программная реализация алгоритма:

mc.PW = function(distances, u, h) {
    weights = mc.PW.kernel(distances / h)
    classes = unique(names(distances))

    weightsByClass = sapply(classes, mc.sumByClass, weights)

    if (max(weightsByClass) == 0) return("") #ни одна точка не попала в окно

    return(names(which.max(weightsByClass)))
}

где distances – расстояние от точки u до всех точек выборки, names(distances) – наименование классов точек выборки.

Реализация алгоритма доступна по ссылке

Протестируем алгоритм на выборке ирисов Фишера:

Прямоугольное ядро:

Треугольное ядро:

Квартическое ядро:

Епанечково ядро:

Алгоритм хорошо себя показывает при h лишь в небольшом диапазоне. Этот диапазон зависит от плотности классифицируемых точек. Если h сделать слишком маленьким, количество точек, способных к классификации заметно уменьшится. Если h сделать слишком большим, он начнет учитывать при классификации слишком дальние точки. Однако стоит отметить, что в некоторых ядрах веса дальних точек гораздо меньше, и качество алгоритма при больших h падает не так сильно. Неудачным ядром (по мнению авторов программной реализации) является прямоугольное, так как веса всех точек внутри окна одинаковые.

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

Гауссовское ядро:

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

Плюсы:

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

Минусы:

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

Для оценки близости объекта u к классу y алгоритм использует следующую функцию:

, где — функция ядра.

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

Простыми словами: алгоритм для каждого обучающего объекта x строит окружность, радиуса h и силы воздействия (потенциала) .

Программная реализация функции классификации:

mc.PF = function(distances, potentials, h) {
    weights = potentials * mc.PF.kernel(distances / h)
    classes = unique(names(distances))

    weightsByClass = sapply(classes, mc.sumByClass, weights)

    if (max(weightsByClass) == 0) return("") #ни одна точка не попала в окно

    return(names(which.max(weightsByClass)))
}

где potentials и h – массивы потенциалов и радиусов для каждой точки выборки соответственно.

Однако, прежде, чем использовать потенциалы, их необходимо подобрать. Этот процесс также входит в обязанности PF.

Описание: Изачально потенциалы заполняются нулями. Далее, пока количество ошибок классификации не достигнет нужного предела, выбираем случайно точку x из выборки. Если для нее классификация выполняется неверно, увеличиваем потенциал на 1 и пересчитываем общее количество ошибок.

Функция вызывается следующим образом:

mc.PF.potentials(points, classes, h, maxMistakes)

где points и classes – массив точек и их классов соответственно, h – массив предпосчитанных радиусов каждой точки, maxMistakes – максимальное количество допустимых ошибок.

Реализация алгоритма доступна по ссылке

Протестируем алгоритм на выборке ирисов Фишера. .Так как целью работы является подбор , а не h, то выбирем их любым образом (пускай для красных точек h = 1, так как они отделены от остальной обучающей выборки; остальным присвоим h = 0.25).

Алгоритм подобрал оптимальные параметры (4 ошибки) за 16 проходов:

Алгоритм заметно сложнее в понимании и реализации, чем предшествующие алгоритмы. К тому же подбор силы потенциалов занимает долгое время. К сожалению, это время оценить невозможно, так как алгоритм является случайным (следующий x в итерации берется случайным образом). При запуске алгоритма на одной и той же выборке, он показывает разное время работы, разную итоговую ошибку, а также разное количество итераций.

Плюсы:

• результат зависит от 2n параметров

Минусы:

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

Выделяют несколько выдов объектов обучения:

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

Алгорим STOLP исключает из выборки выбросы и неинформативные объекты, оставляя лишь нужное количество эталонных. Таким образом улучшается качество классификации, сокращается объем данных и уменьшается время классификации объектов. Другими словами STOLP — алгоритм сжатия данных.

Он использует функцию отступа:

является весовой функцией и зависит от выбранного алгоритма классификации.

Программно функция веса реализуется следующим образом:

mc.STOLP.w = function(sortedDistances, indecies, k) {
    if (length(indecies) == 0) return(0)

    orderedDistances = 1:length(sortedDistances)
    names(orderedDistances) = names(sortedDistances)
    weights = mc.kwKNN.w(orderedDistances, k)

    weights = weights[indecies]
    if (length(weights) == 0) return(0)

    weightsByClass = sapply(unique(names(weights)), mc.sumByClass, weights)
    max(weightsByClass)[1]
}

где indecies – индексы точек, гдя которых ее необходимо подсчитать, если indecies < 0 – индексы точек, которые необходимо выкинуть перед подсчетом. Функция основана на функции веса KNN, поэтому передается дополнительный параметр k и используется mc.kwKNN.w. Реализация STOLP для других алгоритмов классификации будет немного отличаться.

Программная реализация функции отступа:

mc.STOLP.M = function(points, classes, u, class) {
    k = 4
    dist = mc.distances(points, u)
    names(dist) = classes
    sortedDistances = sort(dist)

    #Отсутуп относительно данного класса
    uIndecies = which(sapply(names(sortedDistances), function(v) any(v == class)))
    m1 = mc.STOLP.w(sortedDistances, uIndecies, k)

    #Отступ относительно других классов
    m2 = mc.STOLP.w(sortedDistances, - uIndecies, k)

    #Отступ
    m1 - m2
}

где class – класс Xi-й точки выборки.

Эти функции используются как вспомогательные для основного алгоритма STOLP.

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

Алгоритм вызывается следующей функцией:

mc.STOLP = function(points, classes, mistakes)

где mistakes – допустимое количество ошибок.

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

Сам алгоритм громоздкий, его реализация доступна по ссылке

Протестируем алгоритм на выборке ирисов Фишера. В качестве алгоритма классификации рассмотрим взвешанный KNN и его весовую функцию: .

Между соседними типами объектов не существуют четких разграничений, поэтому границу выбросов установим сами. Будем называть выбросом объект, отступ которого < -1.

Как мы знаем из предыдущих работ, алгоритм kwKNN дает лучший результат (6 ошибок) при практически любом k (возьмем, k = 30).

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

Попытаемся улучшить алгоритм с помощью STOLP, сократив количество ошибок в 2 раза.

С поставленной задачей алгоритм справился за 10 проходов, допустив всего 2 ошибки на эталонной выборке:

Как можно заметить, выборка сократилась со 150-ти до 8-ти объектов (в 18 раз), тем самым ускорив время классификации.

Плюсы:

• в разы уменьшает размер выборки
• способен улучшить качество классификации

Минусы:

• сложный в реализации

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

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

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

Разделяющая поверхность между классами t и s – это геометрическое место точек таких, что максимум апостериорной вероятности достигается одновременно при y = s и y = t.

.

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

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

, в которой

– объект, состоящий из n признаков,

– математическое ожидание,

– ковариационная матрица (положительно определенная, симметричная, невырожденная).

1. Если признаки некореллированы, то есть , то плотности распределения имеют форму эллипсоидов, параллельных осям координат:

2. Если признаки имеют одинаковые дисперсии , линии уровня имеют форму эллипсоидов:

3. Если матрица не диагональна, то линии уровня – эллипсоиды, повернутые относительно оси координат:

Саму программу можно посмотреть здесь: shinyapps.io

Подстановочный алгоритм относится к нормальному дискриминантному анализу.

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

estimateMu = function(points) {
    rows = dim(points)[1]
    cols = dim(points)[2]
    mu = matrix(NA, 1, cols)
    for (col in 1:cols) {
        mu[1, col] = mean(points[, col])
    }
    return(mu)
}

.

estimateCovarianceMatrix = function(points, mu) {
    rows = dim(points)[1]
    cols = dim(points)[2]
    covar = matrix(0, cols, cols)
    for (i in 1:rows) {
        covar = covar + (t(points[i,] - mu) %*% (points[i,] - mu)) / (rows - 1)
    }
    return(covar)
}

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

Программа доступна по ссылке: shinyapps.io.

Алгоритм довольно точно восстанавливает ковариационную матрицу и мат. ожидание объектов выборки.

Однако при малом количестве количестве объектов, точность падает:

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

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

Так как матрицы равны, можем оценить их все вместе:

estimateCovarianceMatrix = function(xy1, mu1, xy2, mu2) {
    rows1 = dim(xy1)[1]
    rows2 = dim(xy2)[1]
    rows = rows1 + rows2
    cols = dim(xy1)[2]
    sigma = matrix(0, cols, cols)

    for (i in 1:rows1)
        sigma = sigma + (t(xy1[i,] - mu1) %*% (xy1[i,] - mu1))

    for (i in 1:rows2)
        sigma = sigma + (t(xy2[i,] - mu2) %*% (xy2[i,] - mu2))

    return(sigma / (rows + 2))
}

Разделяющая поверхность имеет вид:

, где

,

.

getLDFCoeffs = function(covar, mu1, mu2) {
    invCovar = solve(covar)
    alpha = invCovar %*% t(mu1 - mu2)
    beta = (mu1 %*% invCovar %*% t(mu1) - mu2 %*% invCovar %*% t(mu2)) / 2
    list("alpha" = alpha, "beta" = beta)
}

Программа доступна по ссылке: shinyapps.io

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

При одинаковых выборках подстановочный алгоритм визуально близок к ЛДФ и разницу можно увидеть лишь при малом количестве объектов выборки.

Причем, чем меньше объектов в выборке, тем очевидней превосходство ЛДФ над подстановочным алгоритмом при одинаковых ковариационных матрицах.

Пусть и . Алгоритм

является линейным алгоритмом классификации. Если f>0, то алгоритм a относит x к классу +1. Иначе к классу -1.

Отсюда называется уравнением разделяющей поверхности.

Параметр иногда опускают. Однако в таком случае разделяющая поверхность (в нашем случае с 2мя признаками – прямая), соответствующая уравнению , будет всегда проходить через начало координат. Чтобы избежать такого обобщения, будем полагать, что среди признаков x есть константа , тогда роль свобоного коэффициента играет параметр . Тогда разделяющая поверхность имеет вид .

Величина называется отступом объекта относительно алгоритма классификации. Если , алгоритм совершает на объекте ошибку.

– монотонно невозрастающая функция потерь, мажорирует пороговую функцию . Тогда минимизацю суммарных потерь можно рассматривать как функцию вида

Для минимизации применяется метод градиентного спуска.

В начале выбирается некоторое начальное приближение вектора весов w. Не существует единого способа инициализации весов. Хорошей практикой считается инициализировать веса случайными малыми значениями: , где n – количество признаков x.

Далее высчитывается текущая оценка функционала

Затем запускается итерационный процесс, на каждом шаге которого вектор w изменяется в сторону наиболее быстрого убывания Q. Это направление противоположно вектору градиента . Соответственно веса меняются по правилу:

или

,

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

Критерий останова основан на приблизительной оценке Q методом экспоненциальной скользящей средней:

,

где – ошибка алгоритма на случайном элементе ,

– параметр сглаживания, полагаем его равным .

Алгоритм может остановиться в двух случаях:

1. Он не допускает ошибки ни на каком элементе;
2. Значение Q стабилизировано, то есть .

Примечание. Градиентный метод чувствителен к масштабу измерения признаков. Если норма x большая, итерационный процесс может оказаться парализованным. Чтобы этого не произошло, рекомендуется нормализовать признаки:

, где – j-й признак.

В обобщенном виде алгоритм стохастического градиента выглядит следующим образом:

stoh = function(xl, classes, L, updateRule) {
    #изначальная настройка алгоритма
    rows = dim(xl)[1]
    cols = dim(xl)[2]
    w = runif(cols, -1 / (2 * cols), 1 / (2 * cols))
    lambda = 1 / rows

    # начальное Q
    Q = 0
    for (i in 1:rows) {
        margin = sum(w * xl[i,]) * classes[i]
        Q = Q + L(margin)
    }
    Q.prev = Q

    iter = 0
    repeat {
        iter = iter + 1

        # выбрать объекты с ошибкой
        margins = rep(0, rows)
        for (i in 1:rows) {
            xi = xl[i,]
            yi = classes[i]
            margins[i] = sum(w * xi) * yi
        }
        errorIndecies = which(margins <= 0)

        #выходим, если выборки полностью разделены
        if (length(errorIndecies) == 0) break

        # выбираем случайный ошибочный объект          
        i = sample(errorIndecies, 1)
        xi = xl[i,]
        yi = classes[i]

        # высчитываем ошибку
        margin = sum(w * xi) * yi
        error = L(margin)

        # обновляем веса
        eta = 1 / iter
        w = updateRule(w, eta, xi, yi)

        # новое Q
        Q = (1 - lambda) * Q + lambda * error

        # выходим, если Q стабилизировалось
        if (abs(Q.prev - Q) / abs(max(Q.prev, Q)) < 1e-5) break

        # выходим, если слишком много итераций (алгоритм парализован)
        if (iter == 20000) break

        Q.prev = Q #запоминаем Q на предыдущем шаге
    }

    return(w)
}

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

Имеет квадратичную функцию потерь и дельта-правило правило обновления весов .

Имеет кусочно-линейную функцию потерь и правило Хебба для обновления весов .

Также является оптимальный байесовским классификатором из-за своих довольно сильных вероятностных предположений.

Имеет логистическую функцию потерь и логистическое правило обновления весов , где – сигмоидная функция.

Для сравнения классификаторов была разработана программа.

Рассмотрим примеры работы:

1. Разделенные между собой параллельные друг другу выборки.

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

2. Рядом стоящие выборки

3. Смешанные выборки (точного ответа нет)

4. Иногда алгоритмы параллизуются и не способны завершиться, как, например, персептрон на рисунке ниже

В данном примере существует точный ответ, однако, ни один алгоритм до него не дошел.

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