В этом примере показано, как получить точные значения для биномиальных коэффициентов и найти вероятности в бросающих монету экспериментах с помощью Symbolic Math Toolbox.
Задайте символьную функцию, P(n,k)
, это вычисляет вероятность для голов, чтобы подойти точно k
времена из n
броски.
syms P(n,k)
P(n,k) = nchoosek(n,k)/2^n
P(n, k) =
Предположим, вы бросаете монету 2000 раз. Вероятностью, что головы подходят в половине бросков, является P(n, n/2)
, где n = 2000
. Результатом является рациональное выражение с большими количествами и в числителе и в знаменателе. Symbolic Math Toolbox возвращает точный результат.
n = 2000; central = P(n, n/2)
central =
Аппроксимируйте это рациональное число 10-разрядной точностью с помощью digits
и vpa
.
previous_digits = digits(10); vpa(central)
ans =
Вычислите вероятность, что количество "голов" отличается от ожидаемого значения не больше, чем двумя стандартными отклонениями.
sigma = sqrt(n/4); withinTwoSigma = symsum(P(n, k), k, ceil(n/2 - 2*sigma), floor(n/2 + 2*sigma))
withinTwoSigma =
Аппроксимируйте результат числом с плавающей запятой.
vpa(withinTwoSigma)
ans =
Сравните этот результат со следующими двумя оценками, выведенными из кумулятивной функции распределения (cdf) нормального распределения. Оказывается, что взятие средней точки между первым целым числом снаружи и первым целым числом в интервале 2D сигмы дает более точный результат, чем использование самого интервала 2D сигмы.
syms cdf(x)
cdf(x) = 1/2 * (1 + erf((x - n/2)/sqrt(sym(n/2))))
cdf(x) =
estimate1 = vpa(cdf(n/2 + 2*sigma) - cdf(n/2 - 2*sigma))
estimate1 =
estimate2 = vpa(cdf(floor(n/2 + 2*sigma) + 1/2) - ...
cdf(ceil(n/2 - 2*sigma) - 1/2))
estimate2 =
error1 = vpa(estimate1 - withinTwoSigma)
error1 =
error2 = vpa(estimate2 - withinTwoSigma)
error2 =
Восстановите точность по умолчанию для расчетов с плавающей точкой.
digits(previous_digits);
Постройте это распределение для k в Интервал. Графиком является Кривая Гаусса.
k = n/2 + (-2*sigma:2*sigma); plot(k, P(n,k), '-+'); title('2000 coin tosses'); xlabel('Number of heads'); ylabel('Probability');