2011-12-19 1 views
2

Я ищу алгоритм инициализации k-means++. Следующие два шага алгоритма приводят к неравномерному вероятностей:Как выбрать значение из списка с неравномерными вероятностями?

для каждой точки данных х, вычислить D (х), расстояние между х и ближайший центр, который уже был выбран.

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

Как выбрать с указанным взвешенным распределением вероятности в C++?

+0

Является ли «x» скалярным количеством или вектором? –

+0

Попробуйте здесь: http://stats.stackexchange.com/ Или здесь: http://www.physicsforums.com/forumdisplay.php?f = 78 –

+0

x - 2-й пункт, извините. – zebra

ответ

3

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

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

/* 
* xset is an array of points X, 
* cdf is a preallocated array of the same size 
*/ 
function prepare_cdf(X[] xset, float[] cdf) 
{ 
    float S = 0; 
    int N = sizeof(xset); 
    for i = 0:N-1 
    { 
     float weight = /* calculate D(xset[i])^2 here */ 
     // create cumulative sums and write to the element in cdf array 
     S += weight; 
     cdf[i] = S; 
    } 

    // now normalize so the CDF runs from 0 to 1 
    for i = 0:N-1 
    { 
     cdf[i] /= S; 
    } 
} 

function select_point(X[] xset, float[] cdf, Randomizer r) 
{ 
    // generate a random floating point number from a 
    // uniform distribution from 0 to 1 
    float p = r.nextFloatUniformPDF(); 
    int i = binarySearch(cdf, p); 
    // find the lowest index i such that p < cdf[i] 

    return xset[i]; 
} 

Вы называете prepare_cdf один раз, а затем вызвать select_point столько раз, сколько вам нужно для создания случайных точек.

1

я бы следующий подход:

  • итерация над данными точками, сохраняя их D-квадрат-й в double distance_squareds[] или std::vector<double> distance_squareds или этажерке и хранение суммы их D-квадрате-х лет в double sum_distance_squareds ,
  • использовать the drand48 function для выбора случайного числа в [0.0, 1.0) и умножить его на sum_distance_squareds; сохраните результат в random_number.
  • перечислите значение distance_squareds, добавив вместе значения (снова), и как только общее количество выполнения достигнет или превысит random_number, верните точку данных, соответствующую только что добавленному квадрату D.
  • из-за ошибки округления, возможно, что вы закончите цикл без возврата; если это так, просто верните первую точку данных, или последнюю, или что-то еще. (Но не волнуйтесь, это должно быть очень редкий крайний случай.)
3

Дискретные распределения намного проще сделать в C++ с random заголовком и с помощью std::discrete_distribution. Это пример:

#include <iostream> 
#include <map> 
#include <random> 

int main() 
{ 
    std::random_device rd; 
    std::mt19937 gen(rd()); 
    std::discrete_distribution<> d({20,30,40,10}); 
    std::map<int, int> m; 
    for(int n=0; n<10000; ++n) { 
     ++m[d(gen)]; 
    } 
    for(auto p : m) { 
     std::cout << p.first << " generated " << p.second << " times\n"; 
    } 
} 

и это пример вывода:

0 generated 2003 times 
1 generated 3014 times 
2 generated 4021 times 
3 generated 962 times 
0

Здесь у вас есть что-то, что может помочь вам, используя (номер ..) массив с заданным распределением вероятностей (вероятностный ..) он будет генерировать для вас (числа) с этими вероятностями (здесь он их будет считать).

#include <iostream> 
#include <cmath> 
#include <time.h> 
#include <stdlib.h> 
#include <map> 
#include <vector> 
using namespace std; 
#define ARRAY_SIZE(array) (sizeof(array)/sizeof(array[0])) 

int checkDistribution(double random, const map<double, vector<int> > &distribution_map) 
{ 
    int index = 0; 
    map<double, vector<int> >::const_iterator it = distribution_map.begin(); 
    for (; it!=distribution_map.end(); ++it) 
    { 
     if (random < (*it).first) 
     { 
       int randomInternal = 0; 
       if ((*it).second.size() > 1) 
        randomInternal = rand() % ((*it).second.size()); 
       index = (*it).second.at(randomInternal); 
       break; 
     } 
    } 
    return index; 
} 

void nextNum(int* results, const map<double, vector<int> > &distribution_map) 
{ 
    double random = (double) rand()/RAND_MAX; 
    int index = checkDistribution(random,distribution_map); 
    results[index]+=1; 
} 

int main() { 

    srand (time(NULL)); 
    int results [] = {0,0,0,0,0}; 
    int numbers [] = {-1,0,1,2,3}; 
    double prob [] = {0.01, 0.3, 0.58, 0.1, 0.01}; 
    int size = ARRAY_SIZE(numbers); 
    // Building Distribution 
    map<double, vector<int> > distribution_map; 
    map<double, vector<int> >::iterator it; 
    for (int i = 0; i < size; i++) 
    { 
     it = distribution_map.find(prob[i]); 
     if (it!=distribution_map.end()) 
      it->second.push_back(i); 
     else 
     { 
      vector<int> vec; 
      vec.push_back(i); 
      distribution_map[prob[i]] = vec; 
     } 
    } 
    // PDF to CDF transform 
    map<double, vector<int> > cumulative_distribution_map; 
    map<double, vector<int> >::iterator iter_cumulative; 
    double cumulative_distribution = 0.0; 
    for (it=distribution_map.begin();it!=distribution_map.end();++it) 
    { 
     cumulative_distribution += ((*it).second.size() * (*it).first); 
     cumulative_distribution_map[cumulative_distribution] = (*it).second; 
    } 

    for (int i = 0; i<100; i++) 
    { 
     nextNum(results, cumulative_distribution_map); 
    } 
    for (int j = 0; j<size; j++) 
     cout<<" "<<results[j]<<" "; 
    return 0; 
}