C++ Эффективно вычисляет текущую медиану

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

Быстрый выбор используется для вычисления k-го наименьшего элемента в несортированном списке, и эту концепцию также можно использовать для поиска медианы в несортированном списке.

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

6 ответов

Решение

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

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

Затем начинается основной алгоритм. Каждое последующее число во входной последовательности сравнивается с текущей медианой и добавляется в левую кучу, если она меньше текущей медианы, или в правую кучу, если она больше текущей медианы; если входное число равно текущей медиане, оно добавляется к той куче, которая имеет меньшее количество, или к произвольной куче, если они имеют одинаковое количество. Если это приводит к тому, что счетчики двух куч различаются более чем на 1, корень большей кучи удаляется и вставляется в меньшую кучу. Затем текущая медиана вычисляется как корень большой кучи, если они различаются по количеству, или как среднее значение корней двух кучи, если они имеют одинаковый размер.

Код в Scheme и Python доступен в моем блоге.

У Джеффа МакКлинтока средняя оценка. Требуется сохранить только два значения. В этом примере выполняется перебор массива выборочных значений (загрузка ЦП). Кажется, сходятся относительно быстро (около 100 образцов) к оценке медианы. Идея заключается в том, что на каждой итерации медианы в направлении входного сигнала с постоянной скоростью. Скорость зависит от того, какой величиной вы оцениваете медиану. Я использую среднее значение в качестве оценки величины медианы, чтобы определить размер каждого приращения медианы. Если вам нужно, чтобы медиана была точной до 1%, используйте шаг размером 0,01 * в среднем.

float median = 0.0f;
float average = 0.0f;

// for each sample
{
    average += ( sample - average ) * 0.1f; // rough running average.
    median += _copysign( average * 0.01, sample - median );
}

введите описание изображения здесь

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

Это займет время O (lg n) на вставку и время O (lg n) на медиану, в общей сложности O (n lg n) времени плюс пространство O (n).

Вот сбалансированная древовидная структура C++, которая обеспечивает возможность запроса по индексу в отсортированном списке. Поскольку он поддерживает все значения в отсортированном порядке, он не так эффективен, как подход с двумя кучами, но обеспечивает дополнительную гибкость. Например, это также может дать вам рабочий квартиль.

template <typename T>
class Node
{
public:
    T key;
    Node* left;
    Node* right;
    size_t size;

    Node(T k) : key(k)
    {
        isolate();
    }

    ~Node()
    {
        delete(left);
        delete(right);
    }

    void isolate()
    {
        left = NULL;
        right = NULL;
        size = 1;
    }

    void recount()
    {
        size = 1 + (left ? left->size : 0) + (right ? right->size : 0);
    }

    Node<T>* rotateLeft()
    {
        Node<T>* c = right;
        Node<T>* gc = right->left;
        right = gc;
        c->left = this;
        recount();
        c->recount();
        return c;
    }

    Node<T>* rotateRight()
    {
        Node<T>* c = left;
        Node<T>* gc = left->right;
        left = gc;
        c->right = this;
        recount();
        c->recount();
        return c;
    }

    Node<T>* balance()
    {
        size_t lcount = left ? left->size : 0;
        size_t rcount = right ? right->size : 0;
        if((lcount + 1) * 2 < (rcount + 1))
        {
            size_t lcount2 = right->left ? right->left->size : 0;
            size_t rcount2 = right->right ? right->right->size : 0;
            if(lcount2 > rcount2)
                right = right->rotateRight();
            return rotateLeft();
        }
        else if((rcount + 1) * 2 <= (lcount + 1))
        {
            size_t lcount2 = left->left ? left->left->size : 0;
            size_t rcount2 = left->right ? left->right->size : 0;
            if(lcount2 < rcount2)
                left = left->rotateLeft();
            return rotateRight();
        }
        else
        {
            recount();
            return this;
        }
    }

    Node<T>* insert(Node<T>* newNode)
    {
        if(newNode->key < key)
        {
            if(left)
                left = left->insert(newNode);
            else
                left = newNode;
        }
        else
        {
            if(right)
                right = right->insert(newNode);
            else
                right = newNode;
        }
        return balance();
    }

    Node<T>* get(size_t index)
    {
        size_t lcount = left ? left->size : 0;
        if(index < lcount)
            return left->get(index);
        else if(index > lcount)
            return right ? right->get(index - lcount - 1) : NULL;
        else
            return this;
    }

    Node<T>* find(T k, size_t start, size_t* outIndex)
    {
        if(k < key)
            return left ? left->find(k, start, outIndex) : NULL;
        else if(k > key)
            return right ? right->find(k, left ? start + left->size + 1 : start + 1, outIndex) : NULL;
        else
        {
            if(outIndex)
                *outIndex = start + (left ? left->size : 0);
            return this;
        }
    }

    Node<T>* remove_by_index(size_t index, Node<T>** outNode)
    {
        size_t lcount = left ? left->size : 0;
        if(index < lcount)
            left = left->remove_by_index(index, outNode);
        else if(index > lcount)
            right = right->remove_by_index(index - lcount - 1, outNode);
        else
        {
            *outNode = this;
            size_t rcount = right ? right->size : 0;
            if(lcount < rcount)
                return left ? right->insert(left) : right;
            else
                return right ? left->insert(right) : left;
        }
        return balance();
    }

    Node<T>* remove_by_value(T k, Node<T>** outNode)
    {
        if(k < key)
        {
            if(!left)
                throw "not found";
            left = left->remove_by_value(k, outNode);
        }
        else if(k > key)
        {
            if(!right)
                throw "not found";
            right = right->remove_by_value(k, outNode);
        }
        else
        {
            *outNode = this;
            size_t lcount = left ? left->size : 0;
            size_t rcount = right ? right->size : 0;
            if(lcount < rcount)
                return left ? right->insert(left) : right;
            else
                return right ? left->insert(right) : left;
        }
        return balance();
    }
};

template <typename T>
class MyReasonablyEfficientRunningSortedIndexedCollection
{
private:
    Node<T>* root;
    Node<T>* spare;

public:
    MyReasonablyEfficientRunningSortedIndexedCollection() : root(NULL), spare(NULL)
    {
    }

    ~MyReasonablyEfficientRunningSortedIndexedCollection()
    {
        delete(root);
        delete(spare);
    }

    void insert(T key)
    {
        if(spare)
            spare->key = key;
        else
            spare = new Node<T>(key);
        if(root)
            root = root->insert(spare);
        else
            root = spare;
        spare = NULL;
    }

    void drop_by_index(size_t index)
    {
        if(!root || index >= root->size)
            throw "out of range";
        delete(spare);
        root = root->remove_by_index(index, &spare);
        spare->isolate();
    }

    void drop_by_value(T key)
    {
        if(!root)
            throw "out of range";
        delete(spare);
        root = root->remove_by_value(key, &spare);
        spare->isolate();
    }

    T get(size_t index)
    {
        if(!root || index >= root->size)
            throw "out of range";
        return root->get(index)->key;
    }

    size_t find(T key)
    {
        size_t outIndex;
        Node<T>* node = root ? root->find(key, 0, &outIndex) : NULL;
        if(node)
            return outIndex;
        else
            throw "not found";
    }

    size_t size()
    {
        return root ? root->size : 0;
    }
};

Алгоритм скользящего окна медиана:

Медиана - это отсортированный массив, где вы берете из него среднее значение.

простая скользящая реализация с очередью (dqueue) и sorted_array (любая реализация, двоичное дерево, skiparray).

d_queue - это массив, в котором вы можете нажать на хвост и сдвинуть (поп) из передней части массива.

sorted_array - это массив, в который вы вставляете по порядку в позиции, найденной с помощью бинарного поиска.

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

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

  1. Вставка:

    • Вставить в указанном месте в sorted_array,
    • и поместите значение в очередь.
  2. Удаление:

    • Если первый элемент d_queue находится за пределами окна или если в другой очереди вы можете иметь индексы, индекс слишком старый, то:
      • удалить первый элемент из d_queue (s),
      • и двоичный поиск его в отсортированном массиве и удалить его.
  3. Чтобы иметь медиану:

    • Используйте значение (я) в середине sorted_array.
    • Если длина sorted_array равна, используйте элемент посередине.
    • Если длина sorted_array нечетная, используйте в среднем два элемента в середине.
#include<cstdio>
#include<iostream>
#include<queue>
#include <vector>         
#include <functional>

typedef priority_queue<unsigned int> type_H_low;
typedef priority_queue<unsigned int, std::vector<unsigned int>, std::greater<unsigned int> > type_H_high;

size_t signum(int left, int right) {
    if (left == right){
        return 0;
    }
    return (left < right)?-1:1;
}

void get_median( unsigned int x_i, unsigned int &m, type_H_low *l, type_H_high *r) {

    switch (signum( l->size(), r->size() )) {
    case 1: // There are more elements in left (max) heap
        if (x_i < m) {
            r->push(l->top());
            l->pop();
            l->push(x_i);
        } else {
            r->push(x_i);
        }
        break;

    case 0: // The left and right heaps contain same number of elements
        if (x_i < m){
            l->push(x_i);
        } else {
            r->push(x_i);
        }
        break;

    case -1: // There are more elements in right (min) heap
        if (x_i < m){
            l->push(x_i);
        } else {
            l->push(r->top());
            r->pop();
            r->push(x_i);
        }
        break;
    }

    if (l->size() == r->size()){
        m = l->top();
    } else if (l->size() > r->size()){
        m = l->top();
    } else {
        m = r->top();
    }

    return;
}

void print_median(vector<unsigned int> v) {
    unsigned int median = 0;
    long int sum = 0;
    type_H_low  H_low;
    type_H_high H_high;

    for (vector<unsigned int>::iterator x_i = v.begin(); x_i != v.end(); x_i++) {
        get_median(*x_i, median, &H_low, &H_high);
        std::cout << median << std::endl;
    }
}
Другие вопросы по тегам