On Tuesday, December 19, 2017 at 12:48:30 AM UTC+1, bartekltg wrote:
> On Monday, December 18, 2017 at 5:33:04 AM UTC+1, M.M. wrote:
> > Co by było, gdyby rozkład wzorcowy wygenerować innym generatorem,
>
> Powinno być w miarę dobrze, tyle, że to test na to, że oba generatory
> daj ten sam rozkład, nie na to, że generatory sa dobre!
Oczywiście, to tak naprawdę test na to, że oba dają ten sam rozkład.
Można tym też testować ten sam z innym zarodkiem.


> Do Twoich celów musiałbyś więc założyć, że jeden z nich jest dobry;-)
> a to może zostać odebrane jako mętne.
Tak.

> Pamiętaj, że te testy mają test i chi^2.
1) Najpierw jest test urodzinowy. Zliczane jest w kubełkach jak często
było zero powtórzeń, jedno powtórzenie, dwa... N-1.
2) Potem jest liczone chi-kwadrat.
3) Potem chi-kwadrat podstawiane jest do rozkładu - co w sumie nie jest
tutaj przydatne.
4) Ostatecznie jest liczona całka rozkładu (dystrybuanta) od zera do
chi-kwadrat - wychodzi p-value.

> Ten pierwszy też może być chi^2,
> ale nie musi (np kołgomorowa-kogośtam dla ciągłych).
Nie rozumiem. Pierwszy to jest test urodzinowy, jakaś całka, czy
cokolwiek co daje znany rozkład.


> Rozumiem, że mówisz o zastąpieniu właśnie tego pierwszego testu?
> Czy drugiego?
Chyba pierwszego, bo nie bardzo rozumiem jak można zastąpić drugi.
W pierwszym teście dwa generatory budują jakiś rozkład, tutaj konkretnie
rozkład dla testu urodzinowego. W drugim, jeden generator jest traktowany
tak, jakby był teoretycznie idealny i są liczone powyższe punkty od 2 do 4.


> Jeśli pierwszgo, to znany asymptotyczny rozkład (o dystrybuancie F)
> używamy tylko do rozsądnego podziału kubełków, generujemy właściwy
> test, wychodzi x, liczym y = F^-1(x), i ładujemy do odpowiedniego kubełka.
Właśnie tutaj miałem problem. Niby rozkład teoretyczny jest znany dla
testu urodzinowego. Liczyłem poissonem, wartość oczekiwaną brałem z tej
strony:

https://pl.wikipedia.org/wiki/Test_Birthday_Spacings

[
BFloat lambda = 1;
lambda *= n_births;
lambda *= n_births;
lambda *= n_births;
lambda /= 4;
lambda /= l_year;
]

Wartość rozkładu poissona poniżej. Zmienna poisson przechowuje
wartość dla i-tego kubełka.
[
BFloat poisson = boost::multiprecision::exp( lambda * -1.0 );
BFloat chi = 0;

for( unsigned int i=0 ; i<n_births ; i++ ) {
if( i > 0 ) {
poisson *= lambda;
poisson /= i;
}
BFloat E = poisson * loops;
BFloat tmp = sums[i] - E; // (test - E)
tmp = tmp * tmp / E; // (test - E)^2 / E
if( E >= 0.001 || sums[i] > 0 ) {
std::cout << std::setw( 5) << i << " " ;
std::cout << std::setw(11) << sums[i] << " " ;
std::cout << std::setw(16) << std::setprecision(12) << E << " ";
std::cout << std::setw(16) << std::setprecision(12) << tmp << std::endl;
}
chi += tmp;
}
]

Rozkład chi-kwadrat tutaj:
[
BFloat chi_dist( BFloat x , BFloat n ) {
return
boost::multiprecision::pow( x , n / 2 - 1 )
*
boost::multiprecision::exp( x / -2 )
/
boost::math::tgamma( n / 2 )
/
boost::multiprecision::pow( BFloat(2) , n / 2 );
}
]

Skumulowany (dystrybuanta) tutaj:
[
BFloat chi_dist_cum( BFloat x , BFloat n , const int steps=10000 ) {
BFloat sum = 0;
BFloat step = x / steps;
BFloat x1 = 0;
BFloat x2 = step;
BFloat y1 = chi_dist( x1 , n );
while( x2 <= x ) {
BFloat y2 = chi_dist( x2 , n );
sum += step * ( y1 + y2 ) / 2;
x1 = x2;
x2 += step;
y1 = y2;
}
return sum;
}
]

Ale... nie wiem czy warto się w to wczytywać, bo wolfram i arkusz
kalkulacyjny dawały bardzo podobne wartości dla chi-square, gamy,
poissona, itd.

Prędzej test urodzinowy robię źle. Jedna pętla testu jest tutaj:
[
template<class TRnd>
static void loopBirthdaySpacings(
TRnd &rnd, // generator
const unsigned int l_year, // length of year
const unsigned int n_births, // number of briths
QVector<unsigned int> &sums, // sums
QVector<unsigned int> &births, // work memory
QSet<unsigned int> &dupl // work memory
) {
births.clear();
dupl.clear();
for( unsigned int i=0 ; i<n_births ; i++ ) {
births.append( rnd() % l_year );
}
qSort(births);
unsigned int sum = 0;
for( unsigned int i=0 ; i<n_births ; i++ ) {
const unsigned int dist = i==0 ? births[0] : births[i] - births[i-1];
if( ! dupl.contains( dist ) ) {
dupl.insert( dist );
} else {
sum ++;
}
}
sums[ sum ] ++;

}
]

W sums są kubełki. Jeśli było N powtórzeń, to sums[N]++;
dla 0 <= N < ilość urodzeń.

Reszta kodu to wyświetlanie, albo zdublowanie.

W main definiuje się dwa (lub jeden z różnymi seedami) generatory jakie
chcemy przetestować:

ExactRnd exact_rnd(seed1);
TestRnd test_rnd(seed2);

birthdaySpacings<TestRnd,ExactRnd>( test_rnd , exact_rnd , 1<<20 , 1<<7 , 100000000
);

Podaje się długość roku i ile pierwszych urodzeń w roku chcemy analizować.
Ostatni parametr 100000000 to maksymalna ilość pętli.


> Teraz testujemy (np chi^2) czy oba generatory miały ten sam rozkład
> (w oryginale - czy jest to rozkład jenostajny).
Tak, taki miałem zamiar.


> Chyba standardowy chi^2 na porównanie wystarczy (ten, gdzie się porównuje
> do rozkładów brzegowych
> https://en.wikipedia.org/wiki/Chi-squared_test#Examp
le_chi-squared_test_for_categorical_data
> )
To chyba robię dobrze, bo testy dla dwóch różnych generatorów udają się.



>
> BTW, trochę zgaduję, bo nie opisałeś, co i w którym momencie robisz.
> Kod jest i z niego wszytko wynika, ale... przeczytanie kody zajęłoby
> mi znacznie więcej czasu niż informacji o podstawach w poście;>
Mam nadzieję że trochę wyjaśniłem.


> > albo tym samym z innym zarodkiem?
>
> Bardzo zły pomysł. Ten sam generator dla róznych seedów pewnie ma
> ten sam rozkład. Niekoniecznie dobry. Nic takim testem nie wykryjesz,
> poza bardzo złymi generatorami;-)
Hmmm, muszę spróbować i przemyśleć na spokojnie.



> >Otóż są problemy z liczeniem chi-square, bo inny generator może
> >dać w jakimś kubełku zero i by trzeba dzielić przez zero.
>
> To źle dobrane kubełki.
> Dobór kubełków powinien być dobrany do oczekiwanej liczby
> pomiarów.

Hmmm ale przy teście urodzinowym ilość kubełków jest zawsze
równa ilości pierwszych urodzin.

Mnie się zdaje, że albo źle robię test urodzinowy, albo rozkład
poissona słabo przybliża rozkład testu urodzinowego. Może trzeba
wziąć bardzo długi rok (np. miliard dni) i tylko kilka
najczęstszych kubełków... Nie wiem.


Pozdrawiam







>
>
> pzdr
> bartekltg