Методы получения случайных комбинаторных объектов

Описание алгоритма

Задача:

Необходимо сгенерировать случайный комбинаторный объект размера [math] n [/math] с равномерным распределением вероятности, если в наличии есть функция для генерации случайного числа в заданном интервале.

Пусть [math] B = \{b_1, b_2 ..., b_k\} [/math] - множество различных элементов, которые могут находиться в данном комбинаторном объекте.

Будем получать элементы по порядку: сначала определим, какой элемент будет стоять на первом месте, потом на втором и так далее. Считаем, что мы построили префикс длинны [math] i [/math] : [math] P = \{a_1, a_2, \ldots, a_i\} [/math]. Будем выбирать элемент [math] a_{i+1} [/math] из множества всех возможных так, чтобы вероятность выбора элемнта [math] b \in B [/math], была пропорциональна числу комбинторных обьектов размера [math] n [/math] с префиксом [math] P + b [/math]. Для этого разобъем отрезок натуральных чисел [math] [1, s] [/math]. где [math] s [/math] - число различных комбинаторных объектов с текущим префиксом, на [math] k [/math] диапазонов так, чтобы размер диапазаона [math] d_j [/math] был равен числу объектов с префиксом [math] P + b_j [/math]. С помощью функция для генерации случайного числа получим число [math] r [/math] в интервале [math] [1, s] [/math] и добавим к префиксу [math] P [/math] элемент [math] b_j [/math] соответствующий диапазону отрезка в которм находится полученное число.

object randomObject(n: int, k: int):  // [math] n [/math] — размер комбинторного объекта, [math] k [/math] — число различных элемнтов.
  for i = 1 to n                               
    s = number(prefix)  // число комбинаторных объектов с текущим префиксом. 
    r = random(1, s)
    for j = 1 to k  
      if number(prefix + B[j]) < r  // [math] B [/math] — множество всех возможных элементов. 
        r = r - number(prefix + B[j])  // если [math] r [/math] не попало в текщий диапазон — перейдем к следующему.
      else 
        prefix[i] = b[j]
        break
  return prefix

Сложность алгоритма — [math]O(nk) [/math]. Количества комбинаторных объектов с заданными префиксами считаются известными, и их подсчет в сложности не учитывается. Стоит отметить, что подсчет количества комбинаторных объектов с заданным префиксом зачастую является задачей с достаточно большой вычислительной сложностью.

Доказательство корректности

Докажем по индукции, что вероятность получить любой перфикс [math] P [/math] равна [math] C(n)\over{S(P)} [/math], где [math] C(n) [/math] - число различных комбинаторных данного типа длины [math] n [/math], а [math] S(P) [/math] - число различных комбинаторных обьектов с таким префиксом.

• Любой комбинаторный объеут имеет пустой перфикс, следовательно [math] S(\varnothing)=C(n) [/math]. Вероятность получить любой префикс [math] P [/math] длины [math] 1 [/math] равна [math] S(\varnothing)\over{S(P)} [/math], что равно [math] C(n)\over{S(P)} [/math].
• Пусть вероятность получить префикс [math] P [/math] длины [math] l [/math] равна [math] C(n)\over{S(P)} [/math]
• Тогда вероятность получить из [math] P [/math] любой префикс [math] P' [/math] длины [math] l+1 [/math] равна [math] C(n)\over{S(P)} [/math][math]\cdot[/math][math] S(P)\over{S(P')} [/math] , что равно [math] C(n)\over{S(P')} [/math]

Битовые вектора

Рассмотрим алгоритм получения случайного битового вектора. В битовом векторе может находиться только два типа элементов: [math] 0 [/math] и [math] 1 [/math], следовательно [math] k = 2 [/math]. Заметим что для любого префикса длины [math] l [/math] число возможных комбинаторных объектов одинаково и равно, следовательно на каждом шаге алгоритма небходмо выбирать с равной вероятностью [math] 0 [/math] или [math] 1 [/math]

vector<int> randomBitVector(n: int):  // [math] n [/math] — размер битового вектора.
  for i = 1 to n                               
    r = random(0, 1)
    v[i] = r
  return prefix

Сложность алгоритма — [math] O(n) [/math], так как в случае двоичных векторов [math] k [/math] постоянно и равно [math] 2 [/math].

Правильные скобочные последовательности

Рассмотрим алгоритм получения случайной правильной скобочной последовательности. Правильная скобочная пследовательность состоит из двух типов элементов: открывающей и закрывающей скобок, следовательно [math] k = 2 [/math].

Рассмотрим "полуправильные" скобочные последовательности т.е. такие что всякой закрывающей скобке соответствует парная открывающая, но не все открытые скобки закрыты. Такую последовтеьность можно охарактеризовать двумя числами: [math] l [/math] — длина скобочной последовательности и [math] b [/math] — баланс (т.е. разность между количеством открывающих и закрывающих скобок). Заметим что любой префикс правильной скобочной последовательности является "полупраильной" скобочной последовательностью, и что для любого префикса [math] P [/math] длины [math] l [/math] число различных ПСП длины [math] n [/math] равно числу "полуправильных" скобочных последовательностей длины [math] n-l [/math] с таким же балансом как у [math] P [/math].

Научимся считать [math] number(l, b) [/math] — число последовательностей длины [math] l [/math] и баланса [math] b [/math]). Если [math] l = 0 [/math], то ответ понятен сразу: [math] number[0][0] = 1 [/math], все остальные [math] number[0][b] = 0 [/math]. Пусть теперь [math] i \gt 0 [/math], тогда переберём, чему мог быть равен последний символ этой последовательности. Если он был равен [math] '(' [/math], то до этого символа мы находились в состоянии [math] (l-1,b-1) [/math]. Если он был равен [math]')'[/math], то предыдущим было состояние [math](l-1,b+1)[/math]. Таким образом, получаем формулу:

[math]number[l][b] = number[l-1][b-1] + number[l-1][b+1][/math]

(считается, что все значения [math] number[l][b] [/math] при отрицательном [math]j[/math] равны нулю). Этот преподсчет можно выполнить за [math]O(n^2)[/math].

Будем строить префикс следующим образом: на каждом шаге интервал случайных чисел [math] [0, s] [/math] (где [math] s = number[n - l][b] [/math]) , будет разбиваться на два диапазона размерами [math] number[n - l - 1][b + 1] [/math] и [math] number[n - l - 1][b - 1] [/math] , и к префиксу будет прибавляться [math]'('[/math] или [math]')'[/math] если случайное число попало в первый или второй диапазон соответственно.

string randomBracketSequence(n: int):  // [math] n [/math] — длина скобочной последовательности. 
  b = 0
  l = 0
  for i = 1 to n                               
    s = number[n - l][b] 
    r = random(1, s)
     if number[n - l - 1][b + 1] >= r
       l = l + 1
       b = b + 1
       result = result + '('
     else
       l = l + 1
       b = b - 1
       result = result + ')'
  return result

Итоговая сложность алгоритма — [math] O(n) + O(n^2) [/math] на преподсчет.

Разбиения на множества

Разбиение на [math] k [/math] подмножеств

Рассмотрим множество первых [math] n [/math] натуральных чисел: [math] N_n = \{1, 2, ..., n\} [/math]. Необходимо разбить его на [math] k [/math] непустых подмножеств [math] \{B_1, B_2, ..., B_k\} [/math] с равным распределением вероятности.

Будем строить разбиение таким образом, чтобы в результате подмножества [math] \{B_1, B_2, ..., B_k\} [/math] оказались отсортированы в лексикографическом порядке (т.е. чтобы для любых [math]i, j| 1 \leqslant i \lt j \leqslant k [/math] наименьший элемент [math] B_i [/math] был меньше наименьшего элемента [math] B_i [/math]). Для этого будем по очереди добавлять каждое число от [math] n [/math] до [math] 1 [/math] в одно из подмножеств и для каждого из подмножеств начиная с [math] B_n [/math] и заканчивая [math] B_1 [/math] будем выбирать какой элемент будет добавлен в него последним(т.е. будет минимальным).

На каждом шаге префиксом считаем текущее разбиение. Оно характеризуется двумя значениями: [math] l [/math] — число добавленных элементови и [math] m [/math] — число подмножеств для которых определен последний элемент. Заметим, что количество разбиений на подмножества с заднным префиксом равно числу способов разбить еще не добавленные элементы на еще не законченные подмножества так, чтобы они оказались лексикографически упорядочены, то есть равно числу разбиений [math] n-l [/math] элементов на [math] k-m [/math] непустых подмножеств, что равно [math]\lbrace{n-l\atop k-m}\rbrace[/math] (т.е числу Стирлинга второго рода). Таким образом из префикса [math] P [/math] мы можем получить следующий префикс [math] P' [/math] двумя способами:

• Добавить текущий элемент ([math] n-l [/math]) в одно из [math] k-m [/math] незаполненых подмножеств. В таком случае число обьектов с префиксом [math] P' [/math] будет равно [math]\lbrace{n-l-1\atop k-m}\rbrace[/math] .
• Сделать текущий элемент последним в подмножестве [math] B_{k-m} [/math] . В таком случае это подмножество станет законченым, следовательно число обьектов с префиксом [math] P' [/math] будет равно [math]\lbrace{n-l-1\atop k-m-1}\rbrace[/math].

Таким образом на каждом шаге интервал случайных чисел [math] [0, s] [/math] (где [math] s = [/math][math]\lbrace{n-l\atop k-m}\rbrace[/math]) , будет разбиваться на два диапазона размерами [math]\lbrace{n-l-1\atop k-m-1}\rbrace[/math] и [math] (k-m)\cdot [/math][math]\lbrace{n-l-1\atop k-m}\rbrace[/math]. Если случайно сгенерированное число попадет в первый диапазон, то сделаем [math] n-l [/math] последним элементом в подмножестве [math] B_{k-m} [/math] . Иначе добавим [math] n-l [/math] в случайно выбранное из незаконченных подмножеств ([math] \{B_1, B_2, ..., B_{k-m}\} [/math]).

int[] randomSetPartition(n: int k: int):  // [math] n [/math] — количество элементов в множестве, [math] n [/math] — число подмножеств на которые нужно разбить исходное множество. 
  l = 0
  m = 0
  for i = n to 1                               
    s = stirling(n - l, k - m)  // [math] stirling(a, b) [/math] - функция возвращующая число Стирлинга второго рода для заданных аргументов. 
    r = random(1, s)
     if stirling(n - l - 1, k - m - 1) >= r
       result[i] = k - m  // [math] result[i] [/math] - номер подмножества в котором находится элемент [math] i [/math]. 
       l = l + 1
       m = m + 1
     else
       p = random(1, k - m)  // Случайным образом выбираем номер одного из незаконченных подмножеств. 
       result[i] = p
       l = l + 1
  return result

Так как на каждом шаге интервал случайных чисел разделяется только на на два диапазона а всего шагов — [math] n [/math] то итоговая сложность алгоритма — [math] O(n) + O(n^2) [/math] на преподсчет чисел Стирлинга второго рода(если предпосчитать их динамически).