2015-10-30 5 views
1

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

vector<vector<double>> matrix; 

Каждая ячейка может иметь значение, равное или большее 0, и это значение представляет возможность выбранной ячейки. В частности, например, ячейка со значением равна 3 имеет в три раза больше вероятность того, чтобы быть выбраны по сравнению с ячейкой со значением 1.

нужно выбрать N клетки матрицы (0 = < < = N общее количество ячеек) случайным образом, но в зависимости от их вероятности выбора.

Как это сделать?

Алгоритм должен быть как можно быстрее.

+2

std :: discrete_distribution предназначен для этой ситуации. – user515430

+0

@ user515430 В каком направлении? – Nick

+0

@ Ник, если вы считаете свою матрицу единым массивом измерений, вы можете напрямую использовать [discrete_distribution] (http://en.cppreference.com/w/cpp/numeric/random/discrete_distribution) (+1 to @ user515430) – fjardon

ответ

2

Я опишу два метода, А и В.

работает во время приблизительно N * number of cells, и использует пространство O(log number of cells). Хорошо, когда N невелик.

B работает со временем примерно (number of cells + N) * O(log number of cells) и использует пробел O(number of cells). Таким образом, хорошо, когда N является большим (или даже «средним»), но использует намного больше памяти, на практике в некоторых режимах он может быть медленнее по некоторым причинам.


Метод A:

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

Затем образец, как это:

  1. Выберите случайный [i,j] запись матрицы (по выбору i,j каждый равномерно случайным образом из диапазона целых чисел 0 - n-1).

  2. Выберите одноразовый случайный номер p в диапазоне [0, 1].

  3. Проверьте, есть ли matrix[i][j] > p. Если да, верните пару [i][j]. Если нет, вернитесь к шагу 1.

Почему это работает? Вероятность того, что мы закончим на шаге 3 с каким-либо конкретным выходом, равна вероятности того, что была выбрана вероятность [i][j] (это одинаково для каждой записи), умножая вероятность того, что число p было достаточно маленьким. Это пропорционально значению matrix[i][j], поэтому выборка выбирает каждую запись с правильными пропорциями. Возможно также, что на шаге 3 мы вернемся к началу - это что-то уклонение? В принципе, нет. Причина в том, предположим, что мы произвольно выбираем номер k, а затем рассмотрим распределение алгоритма, обусловленного остановкой точно после k раундов.Исходя из предположения, что мы останавливаемся на k'-м раунде, независимо от того, какую величину мы выбрали, распределение, которое мы пробоем, должно быть в точности в соответствии с приведенным выше аргументом. Поскольку, если мы устраним случай, когда p слишком мал, другие возможности имеют правильные пропорции. Так как распределение идеально подходит для каждого значения k, которое мы можем установить, и общее распределение (не обусловленное на k) является средним распределением для каждого значения k, и полное распределение также отлично.

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

Вероятность того, что мы остановимся на шаге 3, может быть определена с учетом условной вероятности того, что мы остановимся на шаге 3, учитывая, что мы выбрали какую-либо конкретную запись [i][j]. По формулам для условного математического ожидания, вы получите, что

Pr[ stop at step 3 ] = sum_{i,j} (1/(n^2) * Matrix[i,j]) 

Поскольку мы предположили, матрица нормализуется, эта сумма сводится только 1/n^2. Таким образом, ожидаемое количество раундов составляет около n^2 (то есть n^2 с точностью до постоянного множителя) независимо от того, какие элементы в матрице. Вы не можете надеяться сделать намного лучше, чем я думаю - это примерно столько же времени, сколько нужно, чтобы просто прочитать все записи в матрице, и трудно отбирать из дистрибутива, что вы даже не можете прочитать все ,

Примечание. То, что я описал, является способом правильного выбора одного элемента - для получения элементов N из одной матрицы, вы можете просто повторить его N раз.


Метод B:

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

В C++ это может выглядеть следующим образом:

// Make histogram 
typedef unsigned int uint; 
typedef std::pair<uint, uint> upair; 
typedef std::map<double, upair> histogram_type; 
histogram_type histogram; 
double cumulative = 0.0f; 
for (uint i = 0; i < Matrix.size(); ++i) { 
    for (uint j = 0; j < Matrix[i].size(); ++j) { 
    cumulative += Matrix[i][j]; 
    histogram[cumulative] = std::make_pair(i,j); 
    } 
} 

std::vector<upair> result; 
for (uint k = 0; k < N; ++k) { 
    // Do a sample (this should never repeat... if it does not find a lower bound you could also assert false quite reasonably since it means something is wrong with rand() implementation) 
    while(1) { 
    double p = cumulative * rand(); // Or, for best results use std::mt19937 or boost::mt19937 and sample a real in the range [0,1] here. 
    histogram_type::iterator it = histogram::lower_bound(p); 
    if (it != histogram.end()) { 
     result.push_back(it->second); 
     break; 
    } 
    } 
} 
return result; 

Здесь время, чтобы сделать гистограмма что-то вроде number of cells * O(log number of cells) так вставляя в карте занимает много времени O(log n). Вам нужна упорядоченная структура данных, чтобы получить дешевый поиск N * O(log number of cells) позже, когда вы делаете повторную выборку. Возможно, вы могли бы выбрать более специализированную структуру данных, чтобы идти быстрее, но я думаю, что есть только ограниченное пространство для улучшения.

Редактировать: Поскольку @Bob__ указывает в комментариях, в методе (B) написано, что потенциально будет некоторая ошибка из-за округления с плавающей запятой, если матрицы достаточно велики, даже используя тип double, при этом линия:

cumulative += Matrix[i][j]; 

проблема заключается в том, что, если cumulative гораздо больше, чем Matrix[i][j] за то, что точности с плавающей точкой могут обрабатывать то это каждый раз, когда это утверждение выполняется вы можете наблюдать существенные ошибки, которые аккумулируют ввести значительную неточность.

Как он предлагает, если это произойдет, самый простой способ исправить это - сначала отсортировать значения Matrix[i][j]. Вы даже можете сделать это в общей реализации, чтобы быть в безопасности - сортировка этих ребят не займет больше времени асимптотически, чем у вас уже есть.

+0

Я думаю, что вам не нужно нормализовать матрицу. Как только вы вычислили сумму, в вашей точке 2 просто выберите случайное число в [0, sum/n^2]. –

+0

Ops, возможно [0, sum]. –

+0

@Bob__: Да, я имею в виду, я думаю, это то же самое, что и нормализация, вы просто умножаете 'p' вместо деления на сумму. Это может быть немного быстрее, поскольку вам не нужно делать 'n^2' дивизии, но это будет не асимптотически быстрее, так как вам все равно придется получить сумму. –