Марковские модели. Примеры задач и приложений

Пример задачи. Болен или здоров

Здравствуйте и вновь добро пожаловать на занятия по теме «Машинное обучение без учителя: скрытые марковские модели на языке Python».

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

Предположим, у нас может быть два состояния: болен и здоров (в том смысле, что отсутствует грипп или простуда). Как известно, большую часть времени мы вполне здоровы, а значит, вероятность перехода от состояния «здоров» к болезни очень и очень невелика. Допустим, она равна 0,005: p(болен|здоров) = 0,005. Следовательно, вероятность перехода от здорового состояния к здоровому равна 0,995: p(здоров|здоров) = 0,995.

Далее, вероятность перехода от болезни к болезни тоже довольно высока – ведь если вы вчера заболели, то и сегодня, скорее всего, будете болеть. Однако вероятность перехода от болезни к здоровому состоянию должна быть значительно выше, чем вероятность обратного перехода, поскольку вы, надо полагать, не более столько же времени, сколько являетесь здоровым. Поэтому допустим, что p(болен|болен) = 0,8, а p(здоров|болен) = 0,2.

Теперь мы полностью определили нашу матрицу переходов состояний и можем провести некоторые вычисления. Например, какова вероятность оставаться здоровым в течение 10 дней подряд при условии, что в начальном состоянии мы здоровы? Она равна 0,9959, то есть 95,6%.

Но почему мы возводим в девятую степень, а не в десятую? Ну, рассмотрим вероятность пребывания в здоровом состоянии двух дней подряд, учитывая, что в первый день мы и так здоровы. Она равна 0,9951. Потому и показатель степени всегда на единицу меньше.

А какова вероятность быть здоровым 100 дней подряд при условии, что мы начинаем отсчёт, будучи здоровыми? Это 0,99599 или 60,9%.

А теперь – упражнение. Что вероятнее – проболеть 10 дней подряд, или проболеть 5 дней подряд, а потом 5 дней быть здоровым при условии, что начинаем мы, пребывая в состоянии болезни? Поставьте видео на паузу и попробуйте посчитать самостоятельно, перед тем как вернуться к видео.

Вероятность проболеть 10 дней подряд равна 0,89 = 13,4%. Вероятность проболеть 5 дней подряд, а потом 5 дней быть здоровым равна 0,84 * 0,2 * 0,9954 = 8%. Следовательно, первый случай более вероятен.

Ещё один вопрос. Какова вероятность того, что здоровое и болезненное состояния будут чередоваться друг с другом на протяжении 10 дней подряд при условии, что мы начинаем с состояния болезни? Интуитивно мы знаем, что это крайне невероятно, но давайте всё же рассчитаем. Поставьте видео на паузу и попробуйте посчитать самостоятельно.

Ответ: 2*10-9 – это действительно совершенно невероятно.

Пример задачи. Ожидаемое количество дней болезни

Задача: каково ожидаемое количество дней пребывания в одном и том же состоянии, например, какое ожидаемое количество дней мы будем болеть? Это сложная математическая задача, но вполне разрешимая. Надо лишь владеть операциями с бесконечными суммами, о которых обычно рассказывают на первом курсе высшей математики.

Сначала определим вероятность, пребывая в состоянии i, продолжать пребывать в состоянии i и на следующий день. Это всего лишь наше A(i,i):

p(s(t) = i|s(t-1) = i) = A(j,i).

Теперь рассмотрим распределение вероятности, которое на самом деле хотим рассчитать. Как рассчитать вероятность того, что мы будем пребывать в состоянии i в течение n переходов до момента, когда перейдём в другое состояние? Вероятность перехода из состояния i в состояние не-i равна 1 – A(i,i):

p(s(t) \neq i|s(t-1) = i) = 1 - A(j,i).

Совместное распределение вероятности, которое мы хотим смоделировать равна

p(s(1) = i,s(2) = i, \;...\;, s(n) = i,s(n+1) \neq i) = A(i,j)^{n-1} (1 - A(j,i)).

При этом  A(i,j)^n-1 значит, как вытекает из предыдущего примера, что мы n раз оставались в состоянии i, а (1-A(i,j) – что мы вышли из этого состояния.

Итак, как же нам решить данное уравнение? Используем для этого арифметику. Как известно, ожидаемое значение n должно быть равно

Преобразуем данное выражение:

и вычислим обе суммы по отдельности.

Первую сумму умножим на a и вычтем:

X = 1+2a+3a^2+...

aX = 0+1a+2a^2+...

—————————-

(1-a) X = 1+a+a^2+...

Это даёт нам другую бесконечную сумму, но, к счастью, куда более лёгкую для решения:

X = 1+a+a^2+...

aX = 0+a+a^2+...


(1-a)X=1.

Следовательно, сумма равна

X = \frac{1}{1-a}.

Подставляя полученное значение, получим:

(1-a)X =1+a+a^2+... = \frac{1}{1-a}.

X = \frac{1}{(1-a)^2}.

Рассмотрим теперь вторую сумму. Умножая на a и вычитая, мы получим похожую сумму, но без единицы. Это значит, что данная сумма равна

X = 1a+2a^2+3a^3+...

aX = 0+1a^2+2a^3+...


(1-a) X = a+a^2+a^3+ ... = \frac{1}{1-a} - 1 = \frac{a}{1-a}.

И окончательно решая, получим:

X = \frac{a}{(1-a)^2}.

Следующий этап – объединение полученных результатов. Получится

E(n) = \frac{1}{(1-a)^2} - \frac{a}{(1-a)^2},

то есть

E(n) = \frac{1}{1-a}.

Теперь мы можем рассчитать количество дней болезни. Подставив числа в данное выражение, получим результат – 5 дней:

E(n) = \frac{1}{1-0.8} = \frac{1}{0.2}=5.

Пример приложения. Оптимизация SEO и показателя отказов

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

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

Во-первых, как люди попадают на вашу страницу? Это ваша домашняя страница? Страница перехода по ссылке? Это самая первая страница, которая, будем надеяться, открывает собой последовательность других страниц. Аналогия с марковской моделью тут следующая. Страница – это начальное распределение, то есть π. Таким образом, в нашей марковской модели вектор π показывает, с какой из страниц пользователи предпочитают начинать.

Теперь что касается последовательности страниц. Если вы полагаете, что люди заходят по ссылке на вашу страницу, нажимают кнопку «Купить», расплачиваются и закрывают окно браузера, то вы можете проверить правильность этого предположения, вычислив вероятность такой последовательности.

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

Последняя аналогия – показатель отказов. С технической точки зрения измерить это довольно сложно, если только вы не являетесь компанией Google с собственным браузером и аналитикой по почти всем сайтам в Интернете. Это связано с тем, что как только пользователь покинул ваш сайт, вы больше не можете запускать код на его компьютере или отслеживать его действия. Но предположим, что у вас в команде есть чудесный программист, который создал идеальную систему логов, благодаря которой вы знаете, когда пользователь просто покинул страницу, а когда – получил то, что ему нужно, и лишь потом покинул сайт. Другими словами, мы знаем, счёл ли пользователь нашу страницу полезной или нет. В связи с этим мы можем подсчитать, какая страница имеет наибольший показатель отказов, вручную проанализировать этот момент и спросить у маркетологов, чем же так особенна эта страница, что люди не считают её полезной, или что на этой странице такого, что люди хотят уйти? Далее мы можем решить эту проблему с надеждой на то, что последующий анализ покажет, как исправленная страница больше не выдаёт высоких показателей отказа. В марковской модели конец последовательности может быть представлен как переход в нулевое состояние.

Далее, поскольку мы ещё напишем код для данного примера, давайте обсудим вид необходимых для этого данных. Мы будем ориентироваться на файл site_data.csv. В нём два столбца last_page_id и next_page_id. У нашего сайта есть 10 страниц с ID от 0 до 9. Начальную страницу можно представить в виде last_page_id = – 1, а next_page_id – это будет актуальная страница. Кроме того, конец последовательности можно представить двумя разными кодами – B (от bounce – отказ) и C (от close – закрыть). Поэтому если мы видим ID страницы, а затем B – это значит, что пользователь зашёл на страницу и сразу ушёл. Если же видим ID страницы, а затем C – значит, пользователь увидел страницу, оставался там некоторое время и, вероятно, нашёл для себя какую-то полезную информацию, а лишь затем закрыл окно. Это можно представить себе так, что наш программист определил время как фактор, указывающий на то, был ли отказ или страницу покинули после некоторой работы.

Перейдём теперь к коду.

Прежде всего создадим словарь transitions для хранения всех переходов. В нём будут храниться в качестве ключа стартовая страница и переходы на другие страницы, после чего подсчитываться, сколько раз случалось открытие той или иной страницы. Это даст возможность преобразовать эти данные в распределение вероятностей. Для этого создадим ещё один словарь row_sums, в котором в качестве ключа будут находиться стартовые состояния, благодаря чему мы сможем делить все переходы, начинающиеся с данных стартовых состояний, на этот row_sums.

import numpy as np

transitions = {}

row_sums = {}

Итак, сначала соберём все данные, проходя по каждой строке.

for line in open(‘site_data.csv’):

    s, e = line.rstrip().split(‘,’)

    transitions[(s, e)] = transitions.get((s, e), 0.) + 1

    row_sums[s] = row_sums.get(s, 0.) + 1

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

for k, v in transitions.iteritems():

    s, e = k

    transitions[k] = v / row_sums[s]

Далее вычислим начальное распределение состояний.

print “initial state distribution:”

for k, v in transitions.iteritems():

    s, e = k

    if s == ‘-1’:

        print e, v

Теперь подсчитаем, какая страница имеет высший показатель отказов. Для этого нужно найти в словаре всё, где конечное состояние имеет B.

for k, v in transitions.iteritems():

    s, e = k

    if e == ‘B’:

print “bounce rate for %s: %s” % (s, v)

Запустим программу.

Мы видим, что страница с ID 9 имеет наибольшее значение – то есть начинать предпочитают именно с неё. Если же мы посмотрим на показатель отказов, то обнаружим, что наивысший показатель также имеет страница с ID 9.

Пример приложения. Построение языковой модели 2-го порядка и генерация фраз

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

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

В качестве данных мы используем сборник стихов Роберта Фроста. Это обычный текстовый файл со всеми стихами, соединёнными вместе.

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

def remove_punctuation(s):

    return s.translate(None, string.punctuation)

tokens = [t for t in remove_punctuation(line.rstrip().lower()).split()]

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

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

Проделав все эти подсчёты, мы должны будем создать массив для всех возможностей. Например, если у нас в данных есть два предложения «Я люблю собак» и «Я люблю кошек», то ключом в словаре будет «Я люблю», а значениями будут [собак, кошек]. Если же в данных будет предложение «Я люблю», то значениями будут [собак, кошек, КОНЕЦ]. Для этого я создал функцию add2dict. Сначала она проверяет, если ли какое-либо значение для ключа, и если нет – создаёт массив, в противном случае добавляет значение в массив:

def add2dict(d, k, v):

    if k not in d:

        d[k] = []

    d[k].append(v)

Собрав все эти массивы возможных следующих слов, мы должны преобразовать их в распределения вероятностей. К примеру, массив [кошка, кошка, собака] станет словарём, где ключ «кошка» будет иметь значение 2/3, а «собака» – 1/3. Вот функция, которая за это отвечает:

def list2pdict(ts):

    d = {}

    n = len(ts)

    for t in ts:

        d[t] = d.get(t, 0.) + 1

    for t, c in d.iteritems():

        d[t] = c / n

    return d

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

def sample_word(d):

    p0 = np.random.random()

    cumulative = 0

    for t, p in d.iteritems(d):

        cumulative += p

        if p0 < cumulative:

            return t

    assert(False) #сюда никогда не должны попадать

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

Если вы не желаете писать код сами, а хотите лишь запустить его, то соответствующий файл называется frost.py. Он находится в папке hmm_class по адресу https://github.com/lazyprogrammer/machine_learning_examples.

Итак, начинаем с импорта библиотек numpy и string.

import numpy as np

import string

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

initial = {}

second_word = {}

transitions = {}

Теперь функция для удаления знаков пунктуации:

def remove_punctuation(s):

    return s.translate(None, string.punctuation)

Теперь функция add2dict. Если ключа k нет в словаре, то создаётся пустой массив, в противном случае – добавляется соответствующее значение к этому массиву:

def add2dict(d, k, v):

    if k not in d:

        d[k] = []

    d[k].append(v)

Теперь можем перебрать файл стихов:

for line in open(‘robert_frost.txt’):

    tokens = remove_punctuation(line.rstrip().lower()).split()

Найдём длину последовательности и запустим цикл по каждому токену в последовательности.

    T = len(tokens)

    for i in xrange(T):

        t = tokens[i]

        if i == 0:

            initial[t] = initial.get(t, 0.) + 1

        else:

            t_1 = tokens[i-1]

            if i == T – 1:

                add2dict(transitions, (t_1, t), ‘END’)

            if i == 1:

                add2dict(second_word, t_1, t)

            else:

                t_2 = tokens[i-2]

                add2dict(transitions, (t_2, t_1), t)

Справившись с этим, нормализуем распределения:

# нормализуем распределения

initial_total = sum(initial.values())

for t, c in initial.iteritems():

    initial[t] = c / initial_total

Теперь функция list2pdict, преобразующая списки возможностей в словари вероятностей:

def list2pdict(ts):

    d = {}

    n = len(ts)

    for t in ts:

        d[t] = d.get(t, 0.) + 1

    for t, c in d.iteritems():

        d[t] = c / n

    return d

Эту функцию можно использовать и для словаря вторых слов, и для словаря переходов:

for t_1, ts in second_word.iteritems():

    second_word[t_1] = list2pdict(ts)

for k, ts in transitions.iteritems():

    transitions[k] = list2pdict(ts)

Теперь напишем код для выбора слова из словаря вероятностей. Функция берёт словарь и генерирует случайное число от 0 до 1.

def sample_word(d):

    p0 = np.random.random()

    cumulative = 0

    for t, p in d.iteritems():

        cumulative += p

        if p0 < cumulative:

            return t

    assert(False)

А теперь создадим функцию для генерирования стихотворения с 4 строками.

def generate():

    for i in xrange(4):

        sentence =[]

Генерируем начальное слово и присоединяем его к предложению:

        w0 = sample_word(initial)

        sentence.append(w0)

Генерируем второе слово и тоже присоединяем к предложению:

        w1 = sample_word(second_word[w0])

        sentence.append(w1)

Теперь переходим в бесконечный цикл.

        while True:

            w2 = sample_word(transitions[(w0, w1)])

            if w2 == ‘END’:

                break

            sentence.append(w2)

            w0 = w1

            w1 = w2

        print(‘ ‘.join(sentence))

generate()

Посмотрим, что у нас получилось… Ну, довольно неплохо.

Пример приложения. Алгоритм ранжирования страниц Google

Вы можете удивиться, но марковские цепи используются в такой современном деле, как ведущая поисковая система в мире.

Основная задача в следующем. У нас есть M веб-страниц, ссылающихся друг на друга, и нам нужно пометить их важность от x(1) до x(M) – все больше 0. То есть, нам нужно дать рейтинг для каждой страницы. Как же этого добиться?

Посещения страниц можно рассматривать как последовательность, а каждую страницу – как состояние. Очень просто. Но откуда берётся рейтинг? Он берётся из предельного распределения. В конечном итоге это доля посещений, которая приходится на данную страницу.

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

Разгадка в том, чтобы использовать структуру интернет-ссылок для определения рейтинга. Основная идея состоит в том, что переход на страницу подобен голосу в пользу её важности. В первом приближении мы можем использовать подсчёт частоты переходов для измерения результатов «голосования». Разумеется, это будет некорректной матрицей вероятностей, поэтому можно разделить каждую строку на их сумму, чтобы та стала равной единице:

где n(i) – общее количество переходов с i-й страницы. Легко убедиться, что в таком случае их сумма равна n(i)/ n(i) = 1, а значит, это корректная марковская матрица.

Правда, у нас всё ещё нет уверенности, что данное предельное распределение уникально. Но обратите внимание, что вышеуказанное решение даёт ответ на сразу несколько проблем! Во-первых, предположим, что вы спамер и размещаете ссылки на свою страницу – скажем, 1000 ссылок. Да, у вас есть 1000 ссылок на свою страницу, но, к несчастью, поскольку матрица переходов должна оставаться корректной матрицей вероятностей, в которой сумма строк равняется 1, то получится, что значимость каждой из ваших ссылок будет равняться всего 1/1000.

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

Но мы все ещё не удостоверились, что предельное распределение существует и что оно уникально. Как же убедиться, что наша модель имеет уникальное стационарное распределение? Дело в том, что ещё в 1910 году кое-кто этим уже занялся и открыл так называемую теорему Перрона-Фробениуса. Она гласит, что если наша матрица переходов является марковской матрицей, то есть сумма всех строк равна 1 и все значения её элементов строго положительны (то есть нет элементов, равных нулю), то стационарное распределение существует и является уникальным. Фактически мы можем стартовать с любого начального состояния и со временем, стремящимся к бесконечности, мы всегда закончим в одном и том же стационарном, а следовательно – и предельном распределении.

Как же нам удовлетворить критерии теоремы Перрона-Фробениуса? Вернёмся к идее сглаживания, с которой мы впервые встретились при обсуждении обучения марковской модели. Основная идея в том, что элементы, равные нулю, стали не равными нулю, то есть добавить маленькую возможность достичь данных состояний, что должно порадовать спамеров.

Итак, создадим равномерное распределение вероятностей

U = \frac{1}{M}.

Это матрица размерностью MxM, где M – количество состояний.

Решением для задачи ранжирования страниц будет взять предыдущую матрицу переходов, умножить её на 0,85 и добавить матрицу равномерного распределения, умноженную на 0,15:

G = 0.85A + 0.15U.

Теперь все элементы строго положительны и удовлетворяют критериям теоремы Перрона-Фробениуса.

А в качестве упражнения предлагаю доказать, что G по-прежнему является корректной марковской матрицей.

Понравилась статья? Поделить с друзьями:
Добавить комментарий

;-) :| :x :twisted: :smile: :shock: :sad: :roll: :razz: :oops: :o :mrgreen: :lol: :idea: :grin: :evil: :cry: :cool: :arrow: :???: :?: :!: