Файл binomialDF.h
#ifndef __BINOMIAL_H__
#include "betaDF.h"
double binomialDF(double trials, double successes, double p); /* * Пусть имеется 'trials' независимых наблюдений * с вероятностью 'p' успеха в каждом. * Вычисляется вероятность B(successes|trials,p) того, что число * успехов заключено между 0 и 'successes' (включительно). */
double rev_binomialDF(double trials, double successes, double y); /* * Пусть известна вероятность y наступления не менее m успехов * в trials испытаниях схемы Бернулли. Функция находит вероятность p * успеха в отдельном испытании. * * В вычислениях используется следующее соотношение * * 1 - p = rev_Beta(trials-successes| successes+1, y). */
double binom_leftCI(double trials, double successes, double level); /* Пусть имеется 'trials' независимых наблюдений * с вероятностью 'p' успеха в каждом * и количество успехов равно 'successes'. * Вычисляется левая граница двустороннего доверительного интервала * с уровнем значимости level. */ double binom_rightCI(double n, double successes, double level); /* Пусть имеется 'trials' независимых наблюдений * с вероятностью 'p' успеха в каждом * и количество успехов равно 'successes'. * Вычисляется правая граница двустороннего доверительного интервала * с уровнем значимости level. */
#endif /* Ends #ifndef __BINOMIAL_H__ */
Файл binomialDF.cpp
/***********************************************************/ /* Биномиальное распределение */ /***********************************************************/
#include <math.h> #include <assert.h>
#include "betaDF.h"
ENTRY double binomialDF(double n, double m, double p) /* * Пусть имеется 'n' независимых наблюдений * с вероятностью 'p' успеха в каждом. * Вычисляется вероятность B(m|n,p) того, что число успехов заключено * между 0 и 'm' (включительно), т.е. * сумму биномиальных вероятностей от 0 до m: * * m * -- ( n ) j n-j * > ( ) p (1-p) * -- ( j ) * j=0 * * Вычисления не подразумевают тупое суммирование - используется * следующая связь с центральным бета-распределением: * * B(m|n,p) = Beta(1-p|n-m,m+1). * * Аргументы должны быть положительными, причем 0 <= p <= 1. */ { assert((n > 0) && (p >= 0) && (p <= 1)); if (m < 0) return 0; else if (m == 0) return pow(1-p, n); else if (m >= n) return 1; else return BetaDF(n-m, m+1).value(1-p); }/* binomialDF */
ENTRY double rev_binomialDF(double n, double m, double y) /* * Пусть известна вероятность y наступления не менее m успехов * в n испытаниях схемы Бернулли. Функция находит вероятность p * успеха в отдельном испытании. * * В вычислениях используется следующее соотношение * * 1 - p = rev_Beta(y|n-m,m+1). */ { assert( (n > 0) && (m >= 0) && (m <= n) && (y >= 0) && (y <= 1) ); return 1-BetaDF(n-m, m+1).inv(y); }/*rev_binomialDF*/
ENTRY double binom_leftCI(double n, double m, double y) /* Пусть имеется 'n' независимых наблюдений * с вероятностью 'p' успеха в каждом * и количество успехов равно 'm'. * Вычисляется левая граница двухстороннего доверительного интервала * с уровнем значимости y. */ { assert( (n > 0) && (m >= 0) && (m <= n) && (y >= 0.5) && (y < 1) ); return BetaDF(m, n-m+1).inv((1-y)/2); }/*binom_leftCI*/
ENTRY double binom_rightCI(double n, double m, double y) /* Пусть имеется 'n' независимых наблюдений * с вероятностью 'p' успеха в каждом * и количество успехов равно 'm'. * Вычисляется правая граница доверительного интервала * с уровнем значимости y. */ { assert( (n > 0) && (m >= 0) && (m <= n) && (y >= 0.5) && (y < 1) ); return BetaDF(m+1, n-m).inv((1+y)/2); }/*binom_rightCI*/
Пример
#include <stdio.h> #include "betaDF.h"
void main(void) { char str[10]; double trials, successes, p, x1, x2, y, level = 0.95;
printf("\nChecking Binomial:\n"); while(1) { printf("\nEnter number of trials (0 to exit): "); gets(str); sscanf(str, "%lf", &trials); if (trials == 0) break; printf("Enter number of successes: "); gets(str); sscanf(str, "%lf", &successes); printf("Enter probability of success: "); gets(str); sscanf(str, "%lf", &p); x1 = binomialDF(trials, successes, p); x2 = binomialDF(trials, successes-1, p); printf("\tBinomial(<=%.0lf|%.0lf,%.4f)=%lg\n", successes, trials, p , x1); printf("\tBinomial(==%.0lf|%.0lf,%.4f)=%lg\n", successes, trials, p , x1-x2); y = rev_binom(trials, successes, x1); printf("\tEstimated p=%lg\n", y); y = binom_leftCI(trials, successes, level); printf("\tCI(%.2lg): %lg", level, y); y = binom_rightCI(trials, successes, level); printf(" <= %.2lg <= %lg\n", successes/trials, y); } }
1 2
8 8 8
|