Bagaimana saya dapat dengan mudah menghasilkan nomor acak mengikuti distribusi normal dalam C atau C ++?

Saya tidak ingin menggunakan Boost.

Saya tahu bahwa Knuth membicarakan hal ini panjang lebar tetapi saya tidak memiliki buku-bukunya sekarang.

Ada banyak metode untuk menghasilkan bilangan terdistribusi Gaussian dari RNG biasa .

The Box-Muller transform umumnya digunakan. Ini menghasilkan nilai dengan benar dengan distribusi normal. Perhitungannya mudah. Anda menghasilkan dua nomor acak (seragam), dan dengan menerapkan rumus padanya, Anda mendapatkan dua nomor acak terdistribusi normal. Kembalikan satu, dan simpan yang lain untuk permintaan berikutnya untuk nomor acak.

C ++ 11

Penawaran C ++ 11 std::normal_distribution, yang akan saya lakukan hari ini.

C atau lebih lama C ++

Berikut beberapa solusi dalam urutan kompleksitasnya:

1. Tambahkan 12 angka acak seragam dari 0 ke 1 dan kurangi 6. Ini akan cocok dengan mean dan deviasi standar variabel normal. Kelemahan yang jelas adalah rentangnya dibatasi hingga ± 6 - tidak seperti distribusi normal sebenarnya.

2. Transformasi Box-Muller. Ini tercantum di atas, dan relatif mudah diterapkan. Namun, jika Anda memerlukan sampel yang sangat tepat, ketahuilah bahwa transformasi Box-Muller yang dikombinasikan dengan beberapa generator seragam mengalami anomali yang disebut Neave Effect ¹ .

3. Untuk presisi terbaik, saya sarankan untuk menggambar seragam dan menerapkan distribusi normal kumulatif terbalik untuk mendapatkan variasi yang terdistribusi normal. Berikut adalah algoritma yang sangat baik untuk distribusi normal kumulatif terbalik.

^{1. HR Neave, "Tentang menggunakan transformasi Box-Muller dengan generator bilangan pseudorandom kongruensial perkalian," Statistik Terapan, 22, 92-97, 1973}

Metode yang cepat dan mudah adalah dengan menjumlahkan sejumlah angka acak yang terdistribusi merata dan mengambil rata-ratanya. Lihat Teorema Batas Pusat untuk penjelasan lengkap mengapa ini berhasil.

Saya membuat proyek C ++ open source untuk benchmark pembuatan nomor acak yang didistribusikan secara normal .

Ini membandingkan beberapa algoritma, termasuk

• Metode teorema batas pusat
• Transformasi Box-Muller
• Metode kutub Marsaglia
• Algoritma ziggurat
• Metode pengambilan sampel transformasi terbalik.
• cpp11randommenggunakan C ++ 11 std::normal_distributiondengan std::minstd_rand(sebenarnya Transformasi Box-Muller dalam dentang).

Hasil floatversi presisi tunggal ( ) pada iMac [email protected], clang 6.1, 64-bit:

Untuk kebenarannya, program memverifikasi mean, deviasi standar, kemiringan dan kurtosis sampel. Diketahui bahwa metode CLT dengan menjumlahkan 4, 8 atau 16 bilangan yang seragam tidak memiliki kurtosis yang baik dibandingkan dengan metode lainnya.

Algoritma Ziggurat memiliki performa yang lebih baik dari yang lain. Namun, ini tidak cocok untuk paralelisme SIMD karena memerlukan pencarian tabel dan cabang. Box-Muller dengan set instruksi SSE2 / AVX jauh lebih cepat (x1.79, x2.99) dibandingkan versi non-SIMD dari algoritma ziggurat.

Oleh karena itu, saya akan menyarankan menggunakan Box-Muller untuk arsitektur dengan set instruksi SIMD, dan mungkin ziggurat sebaliknya.

PS benchmark menggunakan LCG PRNG paling sederhana untuk menghasilkan nomor acak terdistribusi seragam. Jadi mungkin tidak cukup untuk beberapa aplikasi. Tetapi perbandingan kinerja harus adil karena semua implementasi menggunakan PRNG yang sama, jadi benchmark terutama menguji kinerja transformasi.

Berikut contoh C ++, berdasarkan beberapa referensi. Ini cepat dan kotor, lebih baik Anda tidak menemukan kembali dan menggunakan perpustakaan pendorong.

#include "math.h" // for RAND, and rand
double sampleNormal() {
    double u = ((double) rand() / (RAND_MAX)) * 2 - 1;
    double v = ((double) rand() / (RAND_MAX)) * 2 - 1;
    double r = u * u + v * v;
    if (r == 0 || r > 1) return sampleNormal();
    double c = sqrt(-2 * log(r) / r);
    return u * c;
}

Anda dapat menggunakan plot QQ untuk memeriksa hasil dan melihat seberapa baik plot tersebut mendekati distribusi normal nyata (beri peringkat sampel Anda 1..x, ubah peringkat menjadi proporsi dari jumlah total x yaitu. Berapa banyak sampel, dapatkan nilai z dan plot mereka Garis lurus ke atas adalah hasil yang diinginkan).

Gunakan std::tr1::normal_distribution.

Namespace std :: tr1 bukan bagian dari pemacu. Ini adalah namespace yang berisi tambahan pustaka dari C ++ Technical Report 1 dan tersedia di kompiler Microsoft terbaru dan gcc, terlepas dari boost.

Ini adalah cara Anda membuat sampel pada compiler C ++ modern.

#include <random>
...
std::mt19937 generator;
double mean = 0.0;
double stddev  = 1.0;
std::normal_distribution<double> normal(mean, stddev);
cerr << "Normal: " << normal(generator) << endl;

Anda dapat menggunakan GSL . Beberapa contoh lengkap diberikan untuk mendemonstrasikan bagaimana menggunakannya.

Silakan lihat di: http://www.cplusplus.com/reference/random/normal_distribution/ . Ini cara paling sederhana untuk menghasilkan distribusi normal.

Jika Anda menggunakan C ++ 11, Anda dapat menggunakan std::normal_distribution:

#include <random>

std::default_random_engine generator;
std::normal_distribution<double> distribution(/*mean=*/0.0, /*stddev=*/1.0);

double randomNumber = distribution(generator);

Ada banyak distribusi lain yang dapat Anda gunakan untuk mengubah keluaran mesin bilangan acak.

Saya telah mengikuti definisi PDF yang diberikan di http://www.mathworks.com/help/stats/normal-distribution.html dan menghasilkan ini:

const double DBL_EPS_COMP = 1 - DBL_EPSILON; // DBL_EPSILON is defined in <limits.h>.
inline double RandU() {
    return DBL_EPSILON + ((double) rand()/RAND_MAX);
}
inline double RandN2(double mu, double sigma) {
    return mu + (rand()%2 ? -1.0 : 1.0)*sigma*pow(-log(DBL_EPS_COMP*RandU()), 0.5);
}
inline double RandN() {
    return RandN2(0, 1.0);
}

Ini mungkin bukan pendekatan terbaik, tetapi cukup sederhana.

Daftar FAQ comp.lang.c berbagi tiga cara berbeda untuk dengan mudah menghasilkan angka acak dengan distribusi Gaussian.

Anda dapat melihatnya di: http://c-faq.com/lib/gaussian.html

Implementasi Box-Muller:

#include <cstdlib>
#include <cmath>
#include <ctime>
#include <iostream>
using namespace std;
 // return a uniformly distributed random number
double RandomGenerator()
{
  return ( (double)(rand()) + 1. )/( (double)(RAND_MAX) + 1. );
}
 // return a normally distributed random number
double normalRandom()
{
  double y1=RandomGenerator();
  double y2=RandomGenerator();
  return cos(2*3.14*y2)*sqrt(-2.*log(y1));
}

int main(){
double sigma = 82.;
double Mi = 40.;
  for(int i=0;i<100;i++){
double x = normalRandom()*sigma+Mi;
    cout << " x = " << x << endl;
  }
  return 0;
}

Terdapat berbagai algoritma untuk distribusi normal kumulatif terbalik. Yang paling populer dalam keuangan kuantitatif diuji di http://chasethedevil.github.io/post/monte-carlo--inverse-cumulative-normal-distribution/

Menurut pendapat saya, tidak ada banyak insentif untuk menggunakan sesuatu yang lain selain algoritma AS241 dari Wichura : ini adalah presisi mesin, dapat diandalkan, dan cepat. Hambatan jarang terjadi dalam pembuatan bilangan acak Gaussian.

Selain itu, ini menunjukkan kelemahan pendekatan seperti Ziggurat.

Jawaban teratas di sini mendukung Box-Müller, Anda harus menyadari bahwa ia memiliki kekurangan yang diketahui. Saya mengutip https://www.sciencedirect.com/science/article/pii/S0895717710005935 :

dalam literatur, Box-Muller kadang-kadang dianggap agak inferior, terutama karena dua alasan. Pertama, jika seseorang menerapkan metode Box-Muller ke bilangan dari generator kongruensial linier yang buruk, bilangan yang diubah memberikan cakupan ruang yang sangat buruk. Plot angka yang berubah dengan ekor spiral dapat ditemukan di banyak buku, terutama dalam buku klasik Ripley, yang mungkin adalah orang pertama yang melakukan pengamatan ini "

1) Cara intuitif secara grafis untuk menghasilkan bilangan acak Gaussian adalah dengan menggunakan sesuatu yang mirip dengan metode Monte Carlo. Anda akan menghasilkan titik acak dalam kotak di sekitar kurva Gaussian menggunakan pembuat bilangan acak-semu di C. Anda dapat menghitung apakah titik itu di dalam atau di bawah distribusi Gaussian menggunakan persamaan distribusi. Jika titik itu ada di dalam distribusi Gaussian, maka Anda mendapatkan bilangan acak Gaussian sebagai nilai x dari titik tersebut.

Metode ini tidak sempurna karena secara teknis kurva Gaussian berlanjut menuju tak terhingga, dan Anda tidak dapat membuat sebuah kotak yang mendekati tak terhingga dalam dimensi x. Tetapi kurva Guassian mendekati 0 dalam dimensi y dengan cukup cepat jadi saya tidak perlu khawatir tentang itu. Batasan ukuran variabel Anda di C mungkin lebih menjadi faktor pembatas keakuratan Anda.

2) Cara lain adalah dengan menggunakan Teorema Batas Pusat yang menyatakan bahwa ketika variabel acak independen ditambahkan, mereka membentuk distribusi normal. Dengan mengingat teorema ini, Anda dapat memperkirakan bilangan acak Gaussian dengan menambahkan sejumlah besar variabel acak independen.

Metode ini bukan yang paling praktis, tetapi itu diharapkan saat Anda tidak ingin menggunakan pustaka yang sudah ada sebelumnya. Ingatlah bahwa jawaban ini datang dari seseorang dengan sedikit atau tanpa pengalaman kalkulus atau statistik.

Metode Monte Carlo Cara paling intuitif untuk melakukannya adalah dengan menggunakan metode monte carlo. Ambil kisaran yang sesuai -X, + X. Nilai X yang lebih besar akan menghasilkan distribusi normal yang lebih akurat, tetapi membutuhkan waktu lebih lama untuk menyatu. Sebuah. Pilih nomor acak z antara -X sampai X. b. Pertahankan dengan probabilitas di N(z, mean, variance)mana N adalah distribusi gaussian. Jatuhkan sebaliknya dan kembali ke langkah (a).

Lihatlah apa yang saya temukan.

Perpustakaan ini menggunakan algoritma Ziggurat.

Komputer adalah perangkat deterministik. Tidak ada keacakan dalam perhitungan. Selain itu, perangkat aritmatika di CPU dapat mengevaluasi sum atas beberapa set bilangan bulat terbatas (melakukan evaluasi dalam bidang hingga) dan himpunan bilangan rasional nyata yang terbatas. Dan juga melakukan operasi bitwise. Matematika mengambil kesepakatan dengan lebih banyak set hebat seperti [0,0, 1,0] dengan jumlah poin tak terhingga.

Anda dapat mendengarkan beberapa kabel di dalam komputer dengan beberapa pengontrol, tetapi apakah itu memiliki distribusi yang seragam? Saya tidak tahu. Namun jika diasumsikan sinyalnya merupakan hasil akumulasi nilai variabel random independen dalam jumlah besar maka Anda akan menerima variabel random terdistribusi mendekati normal (Terbukti dalam Teori Probabilitas)

Ada algoritma yang disebut - generator acak semu. Seperti yang saya rasakan, tujuan generator acak semu adalah untuk meniru keacakan. Dan kriteria kebaikan adalah: - distribusi empiris terkonvergensi (dalam arti tertentu - pointwise, uniform, L2) ke teoritis - nilai yang Anda terima dari generator acak tampaknya tidak bergantung. Tentu saja itu tidak benar dari 'sudut pandang yang sebenarnya', tetapi kami menganggapnya benar.

Salah satu metode yang populer - Anda dapat merangkum 12 irv dengan distribusi seragam .... Tapi jujur ​​saja selama penurunan Teorema Central Limit dengan bantuan Fourier Transform, Taylor Series, diperlukan asumsi n -> + inf beberapa kali. Jadi misalnya teoritis - Secara pribadi saya tidak mengerti bagaimana orang melakukan penjumlahan 12 irv dengan distribusi seragam.

Saya memiliki teori probilitas di universitas. Dan khususnya bagi saya itu hanya soal matematika. Di universitas saya melihat model berikut:

double generateUniform(double a, double b)
{
  return uniformGen.generateReal(a, b);
}

double generateRelei(double sigma)
{
  return sigma * sqrt(-2 * log(1.0 - uniformGen.generateReal(0.0, 1.0 -kEps)));
}
double generateNorm(double m, double sigma)
{
  double y2 = generateUniform(0.0, 2 * kPi);
  double y1 = generateRelei(1.0);
  double x1 = y1 * cos(y2);
  return sigma*x1 + m;
}

Begitulah cara melakukannya itu hanya contoh, saya kira ada cara lain untuk mengimplementasikannya.

Bukti bahwa itu benar dapat ditemukan dalam buku ini "Moscow, BMSTU, 2004: XVI Probability Theory, Example 6.12, p.246-247" dari Krishchenko Alexander Petrovich ISBN 5-7038-2485-0

Sayangnya saya tidak tahu tentang adanya terjemahan buku ini ke dalam bahasa Inggris.