Rolling rata-rata algoritma dalam C

Saat ini saya bekerja pada sebuah algoritma untuk melaksanakan rolling median filter (analog dengan rolling mean filter) di C. Dari pencarian saya dari sastra, tampaknya ada dua cukup efisien cara untuk melakukannya. Yang pertama adalah untuk memilah-jendela awal dari nilai-nilai, kemudian melakukan pencarian biner untuk menyisipkan nilai-nilai baru dan menghapus yang sudah ada pada setiap iterasi.

Kedua (dari Hardle dan Steiger, 1995, JRSS-C, Algoritma 296) membangun double-berakhir tumpukan struktur, dengan maxheap pada salah satu ujung, minheap pada yang lain, dan median di tengah. Ini menghasilkan linear-waktu algoritma bukan salah satu yang adalah O(n log n).

Berikut ini adalah masalah saya: pelaksana mantan bisa dilakukan, tapi saya harus menjalankan ini jutaan kali seri, sehingga efisiensi banyak hal. Yang terakhir ini terbukti sangat sulit untuk diterapkan. Saya menemukan kode di Trunmed.c file kode untuk statistik paket R, tapi itu lebih terbaca.

Apakah ada yang tahu dari yang ditulis dengan baik C implementasi linear waktu bergulir rata-rata algoritma?

Edit: Link ke Trunmed.kode c http://google.com/codesearch/p?hl=en&sa=N&cd=1&ct=rc#mYw3h_Lb_e0/R-2.2.0/src/library/stats/src/Trunmed.c

Mengomentari pertanyaan (7)

Saya telah melihat R's src/perpustakaan/statistik/src/Trunmed.c beberapa kali aku ingin sesuatu yang serupa juga di mandiri kelas C++ / C subrutin. Perhatikan bahwa ini adalah benar-benar dua implementasi dalam satu, melihat src/perpustakaan/statistik/pria/runmed.Rd (sumber file help) yang mengatakan


\details{
  Apart from the end values, the result \code{y = runmed(x, k)} simply has
  \code{y[j] = median(x[(j-k2):(j+k2)])} (k = 2*k2+1), computed very
  efficiently.

  The two algorithms are internally entirely different:
  \describe{
    \item{"Turlach"}{is the Härdle-Steiger
      algorithm (see Ref.) as implemented by Berwin Turlach.
      A tree algorithm is used, ensuring performance \eqn{O(n \log
        k)}{O(n * log(k))} where \code{n 
Komentar (2)

Saya tidak't menemukan yang modern pelaksanaan c++ struktur data dengan order-statistic sehingga akhirnya menerapkan ide-ide di atas coders link yang disarankan oleh MAK ( Pertandingan Editorial: gulir ke bawah untuk FloatingMedian).

Dua multisets

Ide pertama partisi data ke dalam dua struktur data (tumpukan, multisets dll) dengan O(ln N) per menyisipkan/menghapus tidak memungkinkan kuantil untuk diubah secara dinamis tanpa biaya besar. I. e. kita dapat memiliki rolling rata-rata, atau rolling 75% tetapi tidak keduanya pada saat yang sama.

Segmen pohon

Kedua ide menggunakan segmen pohon adalah O(ln N) untuk menyisipkan/penghapusan/pertanyaan tapi lebih fleksibel. Terbaik dari semua "N" adalah ukuran dari data range. Jadi jika anda bergulir rata-rata memiliki jendela dari satu juta item, tapi data anda bervariasi dari 1..65536, maka hanya 16 operasi yang diperlukan per gerakan rolling jendela dari 1 juta!!

C++ kode ini mirip dengan apa Denis diposting di atas ("di Sini's algoritma sederhana untuk terkuantisasi data")

GNU Statistik Orde Pohon

Sebelum menyerah, saya menemukan bahwa stdlibc++ berisi statistik orde pohon!!!

Ini memiliki dua operasi penting:

iter = tree.find_by_order(value)
order = tree.order_of_key(value)

Lihat sebuah berkas++ manual policy_based_data_structures_test (cari "membagi dan bergabung").

Saya telah dibungkus pohon untuk digunakan dalam fasilitas header untuk compiler yang mendukung c++0x/c++11 gaya parsial typedefs:

#if !defined(GNU_ORDER_STATISTIC_SET_H)
#define GNU_ORDER_STATISTIC_SET_H
#include 
#include 

// A red-black tree table storing ints and their order
// statistics. Note that since the tree uses
// tree_order_statistics_node_update as its update policy, then it
// includes its methods by_order and order_of_key.
template 
using t_order_statistic_set = __gnu_pbds::tree<
                                  T,
                                  __gnu_pbds::null_type,
                                  std::less,
                                  __gnu_pbds::rb_tree_tag,
                                  // This policy updates nodes'  metadata for order statistics.
                                  __gnu_pbds::tree_order_statistics_node_update>;

#endif //GNU_ORDER_STATISTIC_SET_H
Komentar (2)

I've dilakukan [implementasi C][1][ di sini][2]. Beberapa rincian lebih lanjut dalam pertanyaan ini: https://stackoverflow.com/questions/5527437/rolling-median-in-c-turlach-implementation.

Contoh penggunaan:


int main(int argc, char* argv[])
{
   int i,v;
   Mediator* m = MediatorNew(15);

   for (i=0;i
Komentar (2)

Saya menggunakan tambahan ini rata-rata estimator:

median += eta * sgn(sample - median)

yang memiliki bentuk yang sama seperti yang lebih umum berarti estimator:

mean += eta * (sample - mean)

Di sini eta kecil belajar menilai parameter (misalnya 0.001), dan sgn() adalah signum fungsi yang mengembalikan salah satu dari{-1, 0, 1}. (Penggunaan konstan eta seperti ini jika data non-stasioner dan anda ingin melacak perubahan dari waktu ke waktu; jika tidak, untuk sumber stasioner menggunakan sesuatu seperti eta = 1 / n untuk berkumpul, di mana n adalah jumlah sampel yang dilihat sejauh ini.)

Juga, aku dimodifikasi rata-rata estimator untuk membuatnya bekerja untuk sewenang-wenang quantiles. Secara umum, kuantil fungsi memberitahu anda nilai yang membagi data menjadi dua fraksi: p dan 1 - p. Berikut perkiraan nilai ini secara bertahap:

quantile += eta * (sgn(sample - quantile) + 2.0 * p - 1.0)

Nilai p harus berada [0, 1]. Ini pada dasarnya menggeser sgn() fungsi's simetris output {-1, 0, 1} untuk bersandar ke satu sisi, partisi data sampel ke dalam dua merata berukuran sampah (fraksi p dan 1 - p dari data yang kurang dari/lebih besar dari kuantil perkiraan, masing-masing). Perhatikan bahwa untuk p = 0.5, hal ini akan mengurangi rata-rata estimator.

Komentar (2)

Berikut ini's algoritma sederhana untuk terkuantisasi data (bulan):

""" median1.py: moving median 1d for quantized, e.g. 8-bit data

Method: cache the median, so that wider windows are faster.
    The code is simple -- no heaps, no trees.

Keywords: median filter, moving median, running median, numpy, scipy

See Perreault + Hebert, Median Filtering in Constant Time, 2007,
    http://nomis80.org/ctmf.html: nice 6-page paper and C code,
    mainly for 2d images

Example:
    y = medians( x, window=window, nlevel=nlevel )
    uses:
    med = Median1( nlevel, window, counts=np.bincount( x[0:window] ))
    med.addsub( +, - )  -- see the picture in Perreault
    m = med.median()  -- using cached m, summ

How it works:
    picture nlevel=8, window=3 -- 3 1s in an array of 8 counters:
        counts: . 1 . . 1 . 1 .
        sums:   0 1 1 1 2 2 3 3
                        ^ sums[3] < 2  a shorter list, len(y) = len(x) - window + 1
    """
    assert len(x) >= window, (len(x), window)
    # np.clip( x, 0, nlevel-1, out=x )
        # cf http://scipy.org/Cookbook/Rebinning
    cnt = np.bincount( x[0:window] )
    med = Median1( nlevel=nlevel, window=window, counts=cnt )
    y = (len(x) - window + 1) * [0]
    y[0] = med.median()
    for j in xrange( len(x) - window ):
        med.addsub( x[j+window], x[j] )
        y[j+1] = med.median()
    return y  # list
    # return np.array( y )

def pad0__( x, tolen ):
    """ pad x with 0 s, numpy array or list """
    n = tolen - len(x)
    if n > 0:
        try:
            x = np.r_[ x, np.zeros( n, dtype=x[0].dtype )]
        except NameError:
            x += n * [0]
    return x

#...............................................................................
if __name__ == "__main__":
    Len = 10000
    window = 3
    nlevel = 256
    period = 100

    np.set_printoptions( 2, threshold=100, edgeitems=10 )
    # print medians( np.arange(3), 3 )

    sinwave = (np.sin( 2 * np.pi * np.arange(Len) / period )
        + 1) * (nlevel-1) / 2
    x = np.asarray( sinwave, int )
    print "x:", x
    for window in ( 3, 31, 63, 127, 255 ):
        if window > Len:  continue
        print "medianfilter: Len=%d window=%d nlevel=%d:" % (Len, window, nlevel)
            y = medianfilter( x, window=window, nlevel=nlevel )
        print np.array( y )

# end median1.py
Komentar (0)

Rolling rata-rata dapat ditemukan dengan mempertahankan dua partisi angka.

Untuk menjaga partisi menggunakan Min Heap dan Max Heap.

Max Heap akan berisi angka-angka yang lebih kecil dari sama dengan median.

Min Heap akan berisi angka lebih besar dari sama dengan median.

Menyeimbangkan Kendala: jika total jumlah elemen yang bahkan kemudian kedua tumpukan harus memiliki elemen yang sama.

jika jumlah dari unsur-unsur yang aneh kemudian Max Heap akan memiliki satu lagi elemen dari Min Heap.

Rata-rata Elemen: Jika Kedua partisi memiliki jumlah yang sama dari unsur-unsur maka rata-rata akan menjadi setengah dari jumlah maksimal elemen dari partisi pertama dan menit elemen dari partisi kedua.

Jika rata-rata akan menjadi max elemen dari partisi pertama.

Algoritma-
1 - Mengambil dua Tumpukan(1 Min Heap dan 1 Max Heap)
Max Heap akan berisi lebih setengah jumlah elemen
Min Heap akan berisi kedua setengah jumlah elemen

2 - Membandingkan nomor baru dari aliran dengan atas Max Heap, jika lebih kecil atau sama menambahkan bahwa jumlah max heap. Jika tidak, tambahkan nomor di Min Heap.

3 - jika min Heap memiliki elemen lebih dari Max Heap kemudian hapus atas unsur Min Heap dan tambahkan di Max Heap. jika max Heap memiliki lebih dari satu elemen dari Tumpukan Min kemudian hapus atas unsur Max Heap dan tambahkan di Min Heap.

4 - Jika Kedua tumpukan memiliki jumlah yang sama dari unsur-unsur yang lalu rata-rata akan menjadi setengah dari jumlah maksimal elemen dari Max Heap dan min elemen dari Min Heap. Jika rata-rata akan menjadi max elemen dari partisi pertama.


public class Solution {

    public static void main(String[] args) {
        Scanner in = new Scanner(System.in);
        RunningMedianHeaps s = new RunningMedianHeaps();
        int n = in.nextInt();
        for(int a_i=0; a_i < n; a_i++){
            printMedian(s,in.nextInt());
        }
        in.close();       
    }

    public static void printMedian(RunningMedianHeaps s, int nextNum){
            s.addNumberInHeap(nextNum);
            System.out.printf("%.1f\n",s.getMedian());
    }
}

class RunningMedianHeaps{
    PriorityQueue minHeap = new PriorityQueue();
    PriorityQueue maxHeap = new PriorityQueue(Comparator.reverseOrder());

    public double getMedian() {

        int size = minHeap.size() + maxHeap.size();     
        if(size % 2 == 0)
            return (maxHeap.peek()+minHeap.peek())/2.0;
        return maxHeap.peek()*1.0;
    }

    private void balanceHeaps() {
        if(maxHeap.size() < minHeap.size())
        {
            maxHeap.add(minHeap.poll());
        }   
        else if(maxHeap.size() > 1+minHeap.size())
        {
            minHeap.add(maxHeap.poll());
        }
    }

    public void addNumberInHeap(int num) {
        if(maxHeap.size()==0 || num 
Komentar (3)

Hal ini mungkin menunjukkan bahwa ada kasus khusus yang memiliki sederhana solusi yang tepat: ketika semua nilai-nilai dalam aliran adalah bilangan bulat dalam sejumlah (relatif) kecil kisaran yang ditetapkan. Misalnya, menganggap mereka semua harus berada di antara 0 dan 1023. Dalam hal ini hanya mendefinisikan sebuah array 1024 elemen dan menghitung, dan yang jelas semua nilai-nilai ini. Untuk masing-masing nilai dalam aliran kenaikan yang sesuai bin dan menghitung. Setelah aliran berakhir menemukan tempat sampah yang berisi hitungan/2 nilai tertinggi - dengan mudah dicapai dengan menambahkan berturut-turut sampah mulai dari 0. Menggunakan metode yang sama nilai sewenang-wenang urutan peringkat dapat ditemukan. (Ada minor komplikasi jika mendeteksi bin saturasi dan "upgrade" ukuran dari tempat penyimpanan untuk jenis yang lebih besar saat berlari akan dibutuhkan.)

Kasus khusus ini mungkin tampak buatan, tetapi dalam praktek hal ini sangat umum. Hal ini juga dapat diterapkan sebagai pendekatan untuk bilangan real jika mereka berbohong dalam jarak dan "cukup baik" tingkat presisi yang diketahui. Ini akan tahan selama cukup banyak setiap set pengukuran pada kelompok "dunia nyata" benda-benda. Misalnya, ketinggian atau bobot dari sekelompok orang. Tidak cukup besar set? Itu akan bekerja sama dengan baik untuk panjang atau bobot dari semua (individu) bakteri di planet - asumsi seseorang bisa memasok data!

Sepertinya aku salah membaca asli - yang sepertinya ingin jendela geser rata-rata bukan hanya rata-rata dari panjang sungai. Pendekatan ini masih bekerja untuk itu. Beban pertama N aliran nilai-nilai untuk jendela awal, maka untuk N+1th streaming nilai kenaikan yang sesuai bin sambil mengurangi sampah yang sesuai dengan 0th aliran nilai. Hal ini diperlukan dalam hal ini untuk mempertahankan N terakhir nilai untuk memungkinkan pengurangan, yang dapat dilakukan secara efisien oleh siklis pengalamatan array dari ukuran N. Karena posisi rata-rata hanya bisa berubah dengan -2,-1,0,1,2 pada setiap langkah dari jendela geser, isn't perlu untuk jumlah semua sampah hingga rata-rata pada setiap langkah, hanya menyesuaikan "rata-rata pointer" tergantung pada sisi mana(s) sampah yang dimodifikasi. Misalnya, jika kedua nilai baru dan salah satu yang dihapus jatuh di bawah rata-rata saat ini maka itu doesn't perubahan (offset = 0). Metode yang rusak ketika N menjadi terlalu besar untuk menahan terletak dalam memori.

Komentar (0)

Untuk mereka yang ingin berjalan rata-rata di Jawa...PriorityQueue adalah teman anda. O(log N) masukkan, O(1) rata-rata saat ini, dan O(N) hapus. Jika anda mengetahui distribusi data anda, anda dapat melakukan banyak hal lebih baik dari ini.


public class RunningMedian {
  // Two priority queues, one of reversed order.
  PriorityQueue lower = new PriorityQueue(10,
          new Comparator() {
              public int compare(Integer arg0, Integer arg1) {
                  return (arg0 < arg1) ? 1 : arg0 == arg1 ? 0 : -1;
              }
          }), higher = new PriorityQueue();

  public void insert(Integer n) {
      if (lower.isEmpty() && higher.isEmpty())
          lower.add(n);
      else {
          if (n 
Komentar (4)

Jika anda memiliki kemampuan untuk referensi nilai-nilai sebagai fungsi dari titik-titik dalam waktu, anda bisa nilai sampel dengan penggantian, menerapkan bootstrapping untuk menghasilkan dinyalakan rata-rata nilai dalam interval kepercayaan. Ini akan membiarkan anda menghitung diperkirakan rata-rata dengan efisiensi yang lebih besar daripada terus-menerus pemilahan nilai-nilai yang masuk ke dalam struktur data.

Komentar (0)

Berikut ini adalah salah satu yang dapat digunakan ketika output yang tepat adalah tidak penting (untuk keperluan display dll.) Anda perlu totalcount dan lastmedian, ditambah newvalue.

{
totalcount++;
newmedian=lastmedian+(newvalue>lastmedian?1:-1)*(lastmedian==0?newvalue: lastmedian/totalcount*2);
}

Menghasilkan cukup hasil yang tepat untuk hal-hal seperti page_display_time.

Aturan: input stream harus halus pada urutan halaman yang menampilkan waktu dalam hitungan (>30 dll), dan memiliki non-nol rata-rata.

Contoh: page load time, 800 item, 10 ms...3000ms, rata-rata 90ms, median nyata:11ms

Setelah 30 input, median kesalahan ini umumnya <=20% (9ms..12ms), dan akan kurang dan kurang. Setelah 800 input, kesalahan adalah +-2%.

Lain pemikir dengan solusi yang sama di sini: https://stackoverflow.com/questions/11482529/median-filter-super-efficient-implementation/15150968#15150968

Komentar (0)

Berikut ini adalah implementasi java

package MedianOfIntegerStream;

import java.util.Comparator;
import java.util.HashSet;
import java.util.Iterator;
import java.util.Set;
import java.util.TreeSet;

public class MedianOfIntegerStream {

    public Set rightMinSet;
    public Set leftMaxSet;
    public int numOfElements;

    public MedianOfIntegerStream() {
        rightMinSet = new TreeSet();
        leftMaxSet = new TreeSet(new DescendingComparator());
        numOfElements = 0;
    }

    public void addNumberToStream(Integer num) {
        leftMaxSet.add(num);

        Iterator iterMax = leftMaxSet.iterator();
        Iterator iterMin = rightMinSet.iterator();
        int maxEl = iterMax.next();
        int minEl = 0;
        if (iterMin.hasNext()) {
            minEl = iterMin.next();
        }

        if (numOfElements % 2 == 0) {
            if (numOfElements == 0) {
                numOfElements++;
                return;
            } else if (maxEl > minEl) {
                iterMax.remove();

                if (minEl != 0) {
                    iterMin.remove();
                }
                leftMaxSet.add(minEl);
                rightMinSet.add(maxEl);
            }
        } else {

            if (maxEl != 0) {
                iterMax.remove();
            }

            rightMinSet.add(maxEl);
        }
        numOfElements++;
    }

    public Double getMedian() {
        if (numOfElements % 2 != 0)
            return new Double(leftMaxSet.iterator().next());
        else
            return (leftMaxSet.iterator().next() + rightMinSet.iterator().next()) / 2.0;
    }

    private class DescendingComparator implements Comparator {
        @Override
        public int compare(Integer o1, Integer o2) {
            return o2 - o1;
        }
    }

    public static void main(String[] args) {
        MedianOfIntegerStream streamMedian = new MedianOfIntegerStream();

        streamMedian.addNumberToStream(1);
        System.out.println(streamMedian.getMedian()); // should be 1

        streamMedian.addNumberToStream(5);
        streamMedian.addNumberToStream(10);
        streamMedian.addNumberToStream(12);
        streamMedian.addNumberToStream(2);
        System.out.println(streamMedian.getMedian()); // should be 5

        streamMedian.addNumberToStream(3);
        streamMedian.addNumberToStream(8);
        streamMedian.addNumberToStream(9);
        System.out.println(streamMedian.getMedian()); // should be 6.5
    }
}
Komentar (1)

Jika anda hanya memerlukan merapikan rata-rata cepat/cara mudah adalah dengan mengalikan terbaru dengan nilai x dan nilai rata-rata sebesar (1-x) maka tambah mereka. Hal ini kemudian menjadi rata-rata baru.

edit: Bukan apa yang pengguna meminta dan tidak valid secara statistik, tetapi cukup baik untuk banyak kegunaan. I'll meninggalkannya di sini (terlepas dari downvotes) untuk pencarian!

Komentar (3)