Обозначение | |
Область значений | , где m – целое |
Параметры | n – целое положительное число (испытаний), – параметр схемы Бернулли (вероятность "успеха"). Помните? величину 1-p принято обозначать буквой q. |
Плотность (функция вероятности) | Плотность дискретна: . Здесь – число сочетаний из n элементов по m, причем . |
Математическое ожидание | np |
Дисперсия | npq |
Функция распределения |
В классических учебниках по статистике для получения значений биномиального распределения часто рекомендуют использовать формулы, основанные на предельных теоремах (типа формулы Муавра-Лапласа). Необходимо отметить, что с чисто вычислительной точки зрения ценность этих теорем близка к нулю, особенно сейчас, когда практически на каждом столе стоит мощный компьютер. Основной недостаток приведенных аппроксимаций – их совершенно недостаточная точность при значениях n, характерных для большинства приложений. Не меньшим недостатком является и отсутствие сколько-нибудь четких рекомендаций о применимости той или иной аппроксимации (в стандартных текстах приводятся лишь асимптотические формулировки, они не сопровождаются оценками точности и, следовательно, мало полезны). Я бы сказал, что обе формулы пригодны лишь при n < 200 и для совсем грубых, ориентировочных расчетов, причем делаемых “вручную” с помощью статистических таблиц. А вот связь между биномиальным распределением и бета-распределением позволяет вычислять биномиальное распределение достаточно экономно.
Я не рассматриваю здесь задачу поиска квантилей: для дискретных распределений она тривиальна, а в тех задачах, где такие распределения возникают, она, как правило, и не актуальна. Если же квантили все-таки понадобятся, рекомендую так переформулировать задачу, чтобы работать с p-значениями (наблюденными значимостями). Вот пример: при реализации некоторых переборных алгоритмов на каждом шаге требуется проверять статистическую гипотезу о биномиальной случайной величине. Согласно классическому подходу на каждом шаге нужно вычислить статистику критерия и сравнить ее значение с границей критического множества. Поскольку, однако, алгоритм переборный, приходится определять границу критического множества каждый раз заново (ведь от шага к шагу объем выборки меняется), что непроизводительно увеличивает временные затраты. Современный подход рекомендует вычислять наблюденную значимость и сравнивать ее с доверительной вероятностью, экономя на поиске квантилей.
Поэтому в приводимых ниже кодах отсутствует вычисление обратной функции, взамен приведена функция rev_binomialDF, которая вычисляет вероятность p успеха в отдельном испытании по заданному количеству n испытаний, числу m успехов в них и значению y вероятности получить эти m успехов. При этом используется вышеупомянутая связь между биномиальным и бета распределениями.
Фактически, эта функция позволяет получать границы доверительных интервалов. В самом деле, предположим, что в n биномиальных испытаниях мы получили m успехов. Как известно, левая граница двухстороннего доверительного интервала для параметра p с доверительным уровнем равна 0, если m = 0, а для является решением уравнения . Аналогично, правая граница равна 1, если m = n, а для является решением уравнения . Отсюда вытекает, что для поиска левой границы мы должны решать относительно уравнение
Хочу заметить, что если не нужна совсем уж неимоверная точность, то при достаточно больших n можно воспользоваться следующей аппроксимацией [Б.Л. ван дер Варден, Математическая статистика. М: ИЛ, 1960, гл. 2, разд. 7]:
Для вычислений с помощью нижеследующих кодов потребуются файлы betaDF.h, betaDF.cpp (см. раздел о бета-распределении), а также logGamma.h, logGamma.cpp (см. приложение А). Вы можете посмотреть также пример использования функций.
#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__ */ |
/***********************************************************/ /* Биномиальное распределение */ /***********************************************************/ #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*/ |
Дата последней модификации: 25 октября 2000 г.