Как стать автором
Обновить

Распределение для рекордов, циклических перестановок и количества дел, которые вы успеете сделать за свою жизнь

Уровень сложностиСредний
Время на прочтение12 мин
Количество просмотров4.4K

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

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

Сюжет первый. Рекорды

В одном околоматематическом паблике появилась любопытная задачка.

В затылок друг другу выстроились восемь человек. Более высокие загораживают более низких так, что тех не видно. Чему равно математическое ожидание числа людей, которых видно?

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

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

Так субботняя прогулка привела меня к компьютеру, и через часик была готова нехитрая программка, моделирующая людей в очереди, набором чисел. Она честно перебирает все возможные варианты, коих в очереди из nчеловек будет k^n, если люди имеют k различных значений роста, вычисляет количество видимого народа, и даже учитывает распределение людей по росту и ведёт статистику. С её помощью быстро удалось просчитать разные варианты, и они говорили о том, что, всё-таки, для ответа нам нужно знать то, каким именно будет распределение вероятности для высоты людей. То есть, ответ зависит от этого распределения и одним числом его не выразить.

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

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

H_8 = 1 + \frac12 + \frac13 +  ... + \frac18.

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

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

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

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

Я рассматривал варианты всех возможных очередей с людьми, рост которых округлён до некоторых дискретных значений. При этом, в очереди могли оказываться люди с одинаковым (при округлении до 5 или 10 см) ростом. Если мы примем, что рост не округлён и мы имеем дело с непрерывным распределением, то точное совпадение по росту имеет вероятность, эффективно равную нулю, а значит, в очереди все люди будут иметь различный рост. В таком случае, уже не важно, как именно он распределён: мы всегда можем ранжировать их по росту и поставить каждому в соответствие натуральное число, сохраняющее отношение порядка.

Пример ранжирования. Главное не рост, а кто кого выше! Ярко-зелёным цветом выделены рекордные значения, "загораживающие" остальные, если смотреть слева направо.
Пример ранжирования. Главное не рост, а кто кого выше! Ярко-зелёным цветом выделены рекордные значения, "загораживающие" остальные, если смотреть слева направо.

Если мы заменим настоящий рост людей на их порядковый номер в ранжировании по росту, то ответ на вопрос: "Кто кого загораживает?" или, что то же самое: "Кто является рекордсменом?" не изменится. Таким образом, мы можем смело забыть о конкретном распределении, лишь бы оно было непрерывным, и свести задачу от перебора k^nвариантов, к перебору перестановок из n различимых элементов, которых, как известно, n! штук. Итак, задача наша формализуется следующим образом:

Найти распределение для числа рекордов в перестановках длины n.

В комбинаторике полезно поиграть с исследуемой системой на малых числах. Посмотрим на то, как рекорды распределяются в перестановках длины 2, 3 и 4:

Гистограммы для числа рекордов в перестановках. Рекордные величины, возникающие при чтении чисел слева направо, выделены жирным шрифтом.
Гистограммы для числа рекордов в перестановках. Рекордные величины, возникающие при чтении чисел слева направо, выделены жирным шрифтом.

Пусть R(n,k) это количество перестановок длины n таких, чтобы количество рекордов в них было равно k. Выясним, чему будет равно R(n+1,k), то есть, сколько будет способов получить k рекордов, если мы добавим в перестановку ещё один элемент. Если его величина превышает все уже задействованные числа, то мы можем поставить его в конце последовательности и количество способов получить k рекордов будет равно R(n,k-1), k-тым рекордом будет этот последний элемент. Если же его величина оказывается где-то посередине, то существует ровно n способов его поставить так, чтобы в любом из уже рассмотренных R(n,k) вариантов он оказался загорожен каким-нибудь рекордом и не нарушил числа k. Таким образом, получается, что

R(n+1,k) = nR(n,k)+R(n,k-1).

Если в последовательности всего один элемент, то R(1,1)=1. Кроме того, есть пара тривиальных граничных значений для отсутствия рекордов и отсутствия людей в очереди:R(n,0) = R(0,k)=0.

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

\begin{eqnarray*}&{{n}\brack{k}} =  (n-1) {{n-1} \brack {k}} + {{n-1} \brack {k-1}},\\& {{0}\brack{0}} = 1,\quad {{n}\brack{0}} = {{0}\brack{k}} = 0.\end{eqnarray*}

Всех перестановок, как мы помним, n! штук, так что можно от количества перестановок с указанным числом рекордов, перейти к вероятностям.

Вот как выглядит распределение рекордов и его среднее значение для первой дюжины значений n.
Вот как выглядит распределение рекордов и его среднее значение для первой дюжины значений n.

Я позволю себе называть это распределение именем Стирлинга и обозначу Stirling(n). Было бы разумным называть гармоническим распределением, но это имя уже занято. Перечислим основные характеристики распределения Стирлинга.

Носитель распределения: X \in 1,2,..,n.

Функция вероятности:

\mathbb{P}[X] = \frac{1}{n!}{n \brack k}.

Среднее:

\mathbb{E}[X] = H_n ,\quad H_{n} = 1 + \frac{1}{2} + \frac{1}{3} + \dots + \frac{1}{n}.

Дисперсия:

\mathbb{D}[X] = H_n - H_{n,2},\quad H_{n,2} = 1 + \frac{1}{2^2} + \frac{1}{3^2} +\dots + \frac{1}{n^2}.

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

Сразу замечу, что то, что среднее для этого распределения выражается гармоническими числами, является его ключевым свойством, о котором мы подробно поговорим ниже. То, что это так, несложно доказать индукцией по n.Однако, есть и другие способы доказательства, один из них хорошо и подробно изложен в статье Артура Бенджамина с соавторами A Stirling Encounter with Harmonic Numbers, которую я рекомендую тем, кто хочет познакомиться с числами Стирлинга первого рода в контексте комбинаторики.

Сюжет второй. Стохастические цепочки

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

Постановка задачи была такой. Пусть в нашем распоряжении имеется n дней. В этот срок необходимо выполнить некоторую работу, состоящую из цепочки последовательных этапов. Говоря "последовательных" я имею в виду, что i-тый этап можно выполнить только по завершении i-1этапа, никакого распараллеливания процессов не предполагается, так что два дела в один день выполнить не удастся. При этом исполнитель придерживается безалаберной, но вполне реалистичной житейской стратегии: в какой-то, произвольно выбранный день, выполняется первый этап работы; второй этап выполняется в один из оставшихся дней, третий — в какой-либо день после второго, и так далее, до тех пор, пока не наступит дедлайн. Я хотел найти вероятность p_n(k) для того, что цепочка будет иметь максимальную длину k. Под максимальной длиной подразумевается, что k+1 этап выполнить уже не удастся, или, что то же самое, что выполнение k-ого этапа приходится на последний n-ный день.

Мы будем рассуждать, рассматривая размещение точек в обратном порядке. Любая цепочка завершается выполнением последнего этапа в последний день. Вероятность того, что первый и последний этап придётся именно на день с номером n, будет равна 1/n, ведь для нашего безалаберного исполнителя все дни одинаковы. Дальше можно действовать индуктивно. Для произвольного k последний этап обязательно должен случиться в последний день, это, как мы знаем, может случиться с вероятностью 1/n. Потом мы можем выполнить предпоследний этап в любой из предыдущих дней, скажем в день с номером m, сведя при этом задачу к случаю k-1 этапов и n-m дней. Выбор m ограничен сверху числом k-2, поскольку два этапа — последний и предпоследний — уже на местах. Сложив вероятности всех таких вариантов, мы получим вероятность p_n(k). Таким образом, у нас уже есть способ получить точное решение искомой задачи, но для этого нужно знать решения всех входящих в неё подзадач:

\begin{align*}p_0(k)&= p_n(0) = 0,\\p_n(1) &= \frac1n,\\ p_n(k) &= \frac1n\left[p_{n-1}(k-1)+p_{n-2}(k-1)+\dots+p_{k-2}(k-1)\right].\end{align*}

Полученную рекуррентную формулу мне, не без труда, удалось упростить, приведя к форме:

p_n(k) = \frac{n-1}{n}p_{n-1}(k)+\frac{1}{n}p_{n-1}(k-1),

которая соответствует распределению Стирлинга.

То, что ожидаемое число последовательных дел, которые можно успеть выполнить, имея ограниченное время, выражается гармоническим числом, наводит на философские мысли. Гармонические числа растут неограниченно,но чрезвычайно медленно, и имея серьёзные планы на пять, десять или на пятьдесят лет, нужно иметь в виду, что если выполнять их не регулярно, а придерживаясь стратегии "когда-нибудь обязательно сделаю", то и за всю жизнь можно успеть сделать не так уж много.

График роста гармонических чисел
График роста гармонических чисел

Сюжет третий. Циклы в перестановках

В любой справочной или вводной статье, посвящённой числам Стирлинга первого рода, говорится, что они фигурируют в комбинаторике, а именно, отражают количество перестановок, имеющих ровно k циклов. Напомню, что любую перестановку можно однозначно представить в виде произведения непересекающихся циклов. Например, перестановка (3\ 2\ 4\ 1\ 6\ 5)содержит в себе такие циклы: 1 \to 3\to 4\to 1 и 5\to 6\to 5 и одну неподвижную точку 2, которую можно считать тривиальным циклом. На иллюстрации она показана в форме подстановки, для наглядности.

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

Пусть C(n,k)— это число циклов длины k в перестановке длины n. Увеличим длину перестановки на единицу, и вычислим C(n+1,k). При этом возможны два случая:

  1. Дополнительный элемент может оказаться неподвижной точкой, и для того, чтобы число k осталось неизменным, мы должны прибавить этот тривиальный цикл только к тем перестановкам длины n, в которых число циклов равно k-1. Количество таких перестановок C(n,k-1).

  2. Новый элемент может войти в перестановку, имеющую k циклов, став частью одного из них. Это можно сделать вставив новый элемент между любыми двумя элементами перестановки длины n, в том числе, между последним и первым элементами, которые могут входить в один цикл. Таким образом, получается n различных вариантов добавления элемента, не нарушающих числа циклов, в каждую из C(n,k) перестановок с k циклами. Так мы получаем ещё nC(n,k) перестановок.

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

C(n+1,k) = nC(n,k)+C(n,k-1).

Осталось заметить, что перестановка с единственным элементом имеет единственный цикл C(1,1) = 1, и что не существует перестановок без циклов, а также циклов без элементов в перестановке: C(n,0) = C(0,k) = 0.

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

Цветом обозначены элементы классов эквивалентности по количеству рекордов. Перестановки здесь приведены в виде композиции циклов.
Цветом обозначены элементы классов эквивалентности по количеству рекордов. Перестановки здесь приведены в виде композиции циклов.

Что же объединяет эти сюжеты?

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

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

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

\begin{align*}&\begin{array}{rclcl}(x)^1 &=& x\\(x)^2 &=&  x(x+1) &=& x^2+x\\(x)^3 &=& x(x+1)(x+2) &=& x^3+3x^2+2x,\\(x)^4 &=& x(x+1)(x+2)(x+3) &=& x^4+6x^3+11x^2+6x\end{array}\\&(x)^n\ \ \ =\ \ \ x(x+1)(x+2)...(x+n-1)(x+n) =\\ &\quad\quad\ \ \,=\ \ (x)^{n-1}(x+n)\end{align*}

Сравните коэффициенты в этих многочленах с гистограммами, полученными нами для перестановок. В том, что это, действительно, числа Стирлинга первого рода, легко убедиться, раскрыв скобки в выражении, определяющем (x) ^n:

(x)^n = (x)^{n-1}(x+n) = x(x)^{n-1} + n(x)^{n-1}.

Если вычислить коэффициент при x^kв (x)^n,обозначив его, как s(n,k), то получится знакомое выражение:

s(n,k) = s(n-1,k-1) + ns(n-1,k).

Поделив восходящий факториал (x) ^nна (1) ^n=n! , получим многочлен в котором коэффициенты при x^k равны вероятностям для числа k в распределении Стирлинга:

\frac{(x)^n}{n!} =\sum\limits_{k=0}^{n}\frac{1}{n!}{n \brack k}x^k.

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

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

\mathbb{E}[Stirling(n)] = 1 + \frac12 + \frac13 +  ... + \frac1n.

Математическое ожидание линейно, это значит, что матожидание для суммы случайных величин равно сумме матожиданий для каждой из них. Это наводит на мысль: не является ли распределение Стирлинга распределением для суммы случайных величин, каждая из которых имеет математическое ожидание, равное 1/n?

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

Например, многочлен с коэффициентами, соответствующими распределению Бернулли с вероятностью p, имеет вид px+(1-p). Свёрткой n бернуллиевых величин будет произведение этих многочленов, или, согласно биномиальной формуле:

(px+(1-p))^n = \sum\limits_{k=0}^n\binom{n}{k}p^k(1-p)^{n-k}x^k.

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

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

\frac{(x)^n}{n!} = (1\cdot x + 0)\left(\frac12x+\frac12\right)\left(\frac{1}{3}x+\frac23\right)...\left[\frac{1}{n}x+\left(1-\frac{1}{n}\right)\right].

Теперь видно, что каждый множитель здесь — это производящая функция вероятности для распределения Бернулли с параметрами 1, 1/2, 1/3, ... 1/n . Таким образом, мы можем сделать важный вывод:

Сумма случайных величин, подчиняющихся распределениям Бернулли с гармонично убывающими вероятностями, подчиняется распределению Стирлинга:

X_k \sim Bernoulli\left(\frac1k\right),\quad \sum\limits_{k=1}^{n} X_k\sim Stirling(n).
Демонстрация получения распределении Стирлинга через свёртку с распределениями Бернулли.
Демонстрация получения распределении Стирлинга через свёртку с распределениями Бернулли.

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

Stirling(n-1) \circ Bernoulli\left(\frac1n\right) \sim Stirling(n),

где кружком обозначен оператор свёртки. Именно этот итерационный процесс показан на анимации.

Вот как выглядит комбинаторная схема суммирования четырёх случайных величин, подчинённых распределениям Бернулли с вероятностями 1,1/2? 1/3 и 1/4.
Вот как выглядит комбинаторная схема суммирования четырёх случайных величин, подчинённых распределениям Бернулли с вероятностями 1,1/2? 1/3 и 1/4.

Эта связь с распределением Бернулли и есть основное свойство распределения Стирлинга, которое кроется во всех трёх задачах и "уши" которого "торчат" в выражении для среднего значения. Похоже, что в каждой из них мы имеем дело с суммами случайных величин, подчиняющиеся распределению Бернулли с параметром 1/n. Давайте выясним, что это за величины для каждого из трёх сюжетов.

Циклы

Применительно к циклам в перестановках, нетрудно показать, что вероятность в случайной перестановке длины n встретить цикл длины m равна 1/m.

Действительно, среди всех возможных n! перестановок, (n-m)! содержат некоторый конкретный цикл длины m. При этом существует \binom{n}{m}способов наполнить этот цикл элементами и m!/m способов переставить в нём элементы, получая уникальные циклы (деление на m исключает одинаковые циклы, отличающиеся друг от друга лишь циклической перестановкой). Итого, вероятность обнаружить цикл длины m равна:

\frac{(n-m)!}{n!} \binom{n}{m}\frac{m!}{m} = \frac{(n-m)!}{n!} \frac{n!}{m!(n-m)!} \frac{m!}{m} = \frac{1}{m}.

Рекорды

Вероятность того, что в перестановке длины n рекордом будет конкретное число m \in [1,n], равна 1/(n-m+1).Убедиться в этом можно пересчитав сколько раз рекордными оказывались числа 1, 2, 3 и 4 в примерах, приведённых на иллюстрации, показывающей гистограммы для числа рекордов в перестановках .

Число m может стать рекордом, только если все числа превышающие его, оказываются по правую сторону от его позиции в перестановке, при этом, нас не интересует расположение чисел не превышающих m. Перестановку, в которой число m является рекордом, можно схематично показать так: (\dots m \dots M_1\dots  M_2\dots  M_{n-m+1}\dots). Здесь буквами M_i условно обозначены числа превышающие m, количество которых равно n-m+1.

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

(m-1)!\binom{n}{n-m+1},

а всего таких классов будет (n-m)! Перемножив эти два выражения, мы получим общее количество перестановок, в которых число m является рекордом:

(m-1)!\binom{n}{n-m+1}(n-m)! = \frac{n!}{n-m+1}.

Поделив его на n!, получаем искомую вероятность 1/(n-m+1).

Например, в перестановке длиной 4, вероятность того, что число 4 будет рекордом равна 1, для 3 она равна 1/2 , для 2 — 1/3, и для 1 — 1/4.

Стохастические цепочки

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

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

Для каждой цепочки эти вероятности перемножаются, в свою очередь, вероятности для цепочек с одинаковым числом этапов, складываются. Так получается распределение вероятности для длин цепочек в случае четырёх дней. Умножив эти вероятности на 4! = 24, мы получим последовательность чисел Стирлинга для n = 4: 6, 11, 6 и 1.

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

\begin{align*} 1\ \text{день}&: \frac{1}{12}+\frac{1}{24}+\frac{1}{12}+\frac{1}{24} &= \frac{1}{4}\\ 2\ \text{день}&: \frac{1}{8}+\frac{1}{24}+\frac{1}{8}+\frac{1}{24} &= \frac{1}{3}\\ 3\ \text{день}&: \frac{1}{4}+\frac{1}{12}+\frac{1}{8}+\frac{1}{24} &= \frac{1}{2}\\ 4\ \text{день}&: \frac{1}{4}+\frac{1}{12}+\frac{1}{8}+\frac{1}{4}+\frac{1}{24}+\frac{1}{12}+\frac{1}{8}+\frac{1}{24} &= 1 \end{align*}

Переводя все наши задачи на классический язык теории вероятностей, можно все их свести к такой канонической форме:

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

Вероятность того, что конкретный номер mокажется в полученной последовательности равна 1/m, в силу линейности среднего, ожидаемая длина последовательности выражается гармоническим числом H_n, а длина последовательности подчиняется распределению Стирлинга.

Теги:
Хабы:
Всего голосов 26: ↑26 и ↓0+26
Комментарии7

Публикации

Истории

Ближайшие события