2

我有恒定数量的样本，每个样本都有一个概率。现在我想从这个数据源中重新采样以获得相同数量的新样本，每个样本具有相同的概率。

例如：

                                          random | 0.03 | 0.78 | 0.45 | 0.70
                                          -------+------+------+------+------
sample | 0000 | 0001 | 0002 | 0003   RNG  sample | 0000 | 0003 | 0002 | 0003
-------+------+------+------+------ ====> -------+------+------+------+------
 prob. | 0.10 | 0.20 | 0.30 | 0.40         prob. | 0.25 | 0.25 | 0.25 | 0.25

就我而言，概率不会直接给出，而是作为权重给出。然而，概率可以直接从权重中导出，因为所有权重的总和是已知的（但不是恒定的）。

在一个 MATLAB 实现中，我使用了Statistics Toolbox 的randsample函数来实现这个重采样过程：

y = randsample(n,k,true,w)或y = randsample(population,k,true,w)返回一个带替换的加权样本，使用一个正权重向量w，其长度为n。i为 的条目选择整数的概率y是w(i)/sum(w)。通常，w是一个概率向量。randsample不支持无替换加权抽样。

function [samples probabilities] = resample(samples, probabilities)
    sampleCount = size(samples, 1);
    indices = randsample(1 : samplecount, samplecount, 
                         true, probabilities);
    samples = samples(indices, :);
    probabilities = repmat(1 / sample count, samplecount, 1);
end

我现在想将这部分算法移植到 iPad 2 上，用于更新重新采样512 个样本的实时 (~25fps) 数据。因此，时间效率至关重要，因为还将执行其他计算。内存不必最小化。

我研究了 Alias 方法，但似乎结构构建过程非常繁琐，可能不是最有效的解决方案。

是否有任何其他有效的方法可以满足实时要求，或者 Alias 方法是可行的方法？

4

1 回答 1

1

这是一个如何resample在 C中实现你的示例。

typedef int SampleType;
typedef double ProbabilityType;

static ProbabilityType MyRandomFunction(ProbabilityType total)
{
    static boolean_t isRandomReady = 0;
    if ( ! isRandomReady ) {
        srandomdev();
        isRandomReady = 1;
    }

    long randomMax = INT_MAX;
    return (random() % (randomMax + 1)) * (total / randomMax);
}

static void MyResampleFunction(SampleType *samples, ProbabilityType *probabilities, size_t length)
{
    ProbabilityType total = 0;

    // first, replace probabilities with sums
    for ( size_t i = 0; i < length; i++ )
        probabilities[i] = total += probabilities[i];

    // create a copy of samples as samples will be modified
    SampleType *sampleCopies = malloc(sizeof(SampleType) * length);
    memcpy(sampleCopies, samples, sizeof(SampleType) * length);

    for ( size_t i = 0; i < length; i++ )
    {
        ProbabilityType probability = MyRandomFunction(total);

        // We could iterate through the probablities array but binary search is more efficient

        // This is a block declaration
        int (^comparator)(const void *, const void *);

        // Blocks are the same a function pointers
        // execept they capture their enclosing scope
        comparator = ^(const void *leftPtr, const void *rightPtr) {

            // leftPtr points to probability
            // rightPtr to an element in probabilities

            ProbabilityType curr, prev;
            size_t idx = ((const ProbabilityType *) rightPtr) - probabilities;
            curr = probabilities[idx];                   // current probablity
            prev = idx > 0 ? probabilities[idx - 1] : 0;   // previous probablity

            if ( curr < probability )
                return 1;
            if ( prev > probability )
                return -1;

            return 0;
        };

        void *found = bsearch_b(&probability,            // the searched value
                                probabilities,           // the searched array
                                length,                  // the length of array
                                sizeof(ProbabilityType), // the size of values
                                comparator);             // the comparator

        size_t idx = ((const ProbabilityType *) found) - probabilities;
        samples[i] = sampleCopies[idx];
    }

    // now, probabilities are all the same
    for ( size_t i = 0; i < length; i++ )
        probabilities[i] = 1.0 / length;

    // Now the can dispose of the copies
    free(sampleCopies);
}

static void MyTestFunction()
{
    SampleType samples[4] = {0, 1, 2, 3};
    ProbabilityType probabilities[10] = {0.1, 0.2, 0.3, 0.4};
    MyResampleFunction(samples, probabilities, 4);
}
于 2011-12-22T16:21:47.673 回答