R统计学(01): 伯努利分布、二项分布

想获取本文完整代码的下载链接,可关注微信公众号”R语言和Python学堂”,并回复发文日期”20181112″。

在接下来的一些博客将介绍各种离散概率分布(Discrete probability distribution)。

所谓”离散”是指随机变量只能取有限个或者可列举无穷多个值。常见的离散分布有二项分布、泊松分布、几何分布、负二项分布和多项分布等。

在介绍具体分布之前,先介绍离散分布中的两个重要函数:

  • 概率质量函数(probability mass function,简称PMF):是离散随机变量X在各特定取值上的概率P(x),其总和为1。与概率密度函数f(x)不同,概率质量函数是对离散随机变量定义的,本身就代表该值的概率;而概率密度函数是对连续随机变量定义的,本身不是概率,它在某区间内的积分才是概率。

  • 累积分布函数(cumulative distribution function, 简称CDF):定义为F(x)=P(X≤x),是单调递增的,且满足:F(-∞)=0和F(+∞)=1。对离散分布而言,它是所有小于等于x的值出现的概率之和。

下图是离散分布的一个例子,随机变量的可能取值为1,2和3,取值概率分别为0.25,0.5和0.25,红色虚线是其累积分布函数F(x)。

《R统计学(01): 伯努利分布、二项分布》 离散分布

好了,回到今天的主题,下面将介绍伯努利分布和二项分布。

1. 伯努利分布

在现实生活中,许多事件的结果往往只有两个。例如:抛硬币,正面朝上的结果只有两个:国徽或面值;检查某个产品的质量,其结果只有两个:合格或不合格;购买彩票,开奖后,这张彩票的结果只有两个:中奖或没中奖;拨打女朋友电话:接通或没接通。。。以上这些事件都可被称为伯努利试验

《R统计学(01): 伯努利分布、二项分布》 雅各布·伯努利

伯努利试验是单次随机试验,只有”成功(值为1)”或”失败(值为0)”这两种结果,是由瑞士科学家雅各布·伯努利(1654 – 1705)提出来的。

其概率分布称为伯努利分布(Bernoulli distribution),也称为两点分布或者0-1分布,是最简单的离散型概率分布。我们记成功概率为p(0≤p≤1),则失败概率为q=1-p,则:

  • 其概率质量函数为:

《R统计学(01): 伯努利分布、二项分布》

  • 其期望值为:

《R统计学(01): 伯努利分布、二项分布》

  • 其方差为:

《R统计学(01): 伯努利分布、二项分布》

2. 二项分布

2.1 定义

假设某个试验是伯努利试验,其成功概率用p表示,那么失败的概率为q=1-p。进行n次这样的试验,成功了x次,则失败次数为n-x,发生这种情况的概率可用下面公式来计算:

《R统计学(01): 伯努利分布、二项分布》

我们称上面的公式为二项分布(Binomial distribution)的概率质量函数。其中

《R统计学(01): 伯努利分布、二项分布》

组合公式,表示从n个不同元素中取出x个元素的所有组合的个数。感叹号!是阶乘运算符,例如:5!=5x4x3x2=120,在R中可用内置的factorial()函数来计算阶乘。此外R内置choose()函数可用来直接计算组合数。看个例子:

##20选5
> factorial(20)/(factorial(5)*factorial(20-5))
[1] 15504
> choose(20, 5)
[1] 15504

从二项分布公式可知,概率分布由试验次数n和”成功”概率p决定,因此二项分布的概率质量函数可以简写为X~B(n, p)。

2.2 性质

二项分布的均值和方差分别为np和npq,这个结论这里就不推导了,具体可参考https://en.wikipedia.org/wiki/Binomial_distribution

二项分布的另一个性质是其分布形状的变化规律。从二项分布概率质量函数P(x)可知,概率分布只与试验次数n和成功概率p有关,其分布形状的变化规律为:

  • “成功”概率p越接近0.5(也即”成功”概率与”失败”概率越接近),二项分布将越对称。保持二项分布试验的次数n不变,随着成功概率p越接近0.5,二项分布逐渐对称,且近似于均值为np、方差为npq的正态分布。(见下图的第一排3个子图)

  • 对于任意”成功”概率p,无论其距离0.5有多远,随着试验次数n的增加,二项分布与均值为np、方差为npq的正态分布越来越接近。(见下图的第二排3个子图)

以上两个二项分布形状变化规律,可明显由下图观察出来。图中的横轴代表试验”成功”的次数;纵轴代表次数对应的概率;红线是均值为np、方差为npq的正态分布曲线

《R统计学(01): 伯努利分布、二项分布》 不同n和p下的二项分布

由此可见,二项分布是一个概率分布族,随着试验次数n和成功概率p的不同而不同,且它与正态分布关系密切。

绘制上图的R代码为:

> par(mar=c(2.5, 4, 2, 0.1), mfrow=c(2,3))
> n <- c(8, 8, 8, 5, 10, 30)
> p <- c(0.1, 0.3, 0.5, 0.2, 0.2, 0.2)
> for(i in 1:length(p))
+ {
+     n1 <- n[i]
+     p1 <- p[i]
+     x <- 0:n1
+     px <- choose(n1, x)*p1^x*(1-p1)^(n1-x)
+     
+     ## 绘制二项分布的概率分布
+     plot(x, 
+          type='n',
+          xlim=c(0, n1),
+          ylim=c(0, max(px)+0.025),
+          ylab='P(x)',
+          main=paste0('n=',n1,', p=',p1),
+          cex.main=2)
+     rect(x-0.1, 0, x+0.1, px, col='green')
+     
+     ## 添加正态分布曲线,均值为np和方差为npq
+     xv <- seq(-1, n1, length.out=1000)
+     yv <- dnorm(xv, mean=n1*p1, sd=sqrt(n1*p1*(1-p1)))
+     lines(xv, yv, col='red', lwd=1.5)
+ }
2.3 R中的相关函数

对于二项分布,R中有四个相关函数可用,分别是:

  • dbinom(x, size, prob):返回成功x次的概率

  • pbinom(q, size, prob):返回至多成功x次的概率,即累积概率

  • qbinom(p, size, prob):返回相应分位点x,详情见下面的例子

  • rbinom(n, size, prob):返回每组试验的成功次数

这四个函数都有sizeprob参数,分别对应于二项分布的试验次数n和成功概率p。下面通过一个例子来了解如何使用它们:

假设我们玩捉泥鳅游戏,每次捉住的概率为0.1(p=0.1),总共玩6次(size=6)。

第一个问题:6次中有2次捉住泥鳅的概率有多大?这时就要用到dbinom(x, size, prob)函数,其中x参数指定成功的次数,函数返回相应概率,比如:

> size <- 6   ##试验次数
> p <- 0.1  ##成功概率

> dbinom(2, size, p)  ##成功捉住2次的概率
[1] 0.098415
> dbinom(0:size, size, p)  ##整个概率分布
[1] 0.531441 0.354294 0.098415 0.014580 0.001215 0.000054 0.000001
> sum(dbinom(0:size, size, p))  ##所有概率之和为1
[1] 1

通过dbinom()函数可以很容易绘制出概率分布图:

《R统计学(01): 伯努利分布、二项分布》

从上图可知,一次都捉不住的概率很大(大于50%)

上图的R代码如下:

> size <- 6
> p <- 0.1
> barplot(dbinom(0:size, size, p),
+         names=0:size,
+         xlab='成功次数',
+         ylab='P(x)',
+         col="red",
+         cex.lab=1.4)

第二个问题:6次中至多有3次捉住泥鳅的概率有多大?这时就要用到pbinom(q, size, prob)函数,其中q参数指定至多次数(这里是3),函数返回相应累积概率,比如:

> pbinom(3, size, p)  ##至多成功捉住3次的概率
[1] 0.99873

## 绘制整个累积分布
> barplot(pbinom(0:size, size, p),
+         names=0:size,
+         xlab='成功次数',
+         ylab='F(x)',
+         col="red",
+         cex.lab=1.4)

《R统计学(01): 伯努利分布、二项分布》

从上图和计算可知,6次中至多有3次捉住泥鳅的概率非常大(约99.873%)

第三个问题:有时想知道,90%概率下我们至多能捉住多少次泥鳅?这时就要用到qbinom(p, size, prob)函数,其中p参数指定概率(这里是0.9),函数返回相应分位点x(即F(x)≥0.9对应的最小x值),比如:

> size <- 6
> p <- 0.1
> qbinom(0.9, size, p)
[1] 2

计算结果显示,6次中90%概率下我们至多能捉住2次泥鳅

最后一个问题:重复10000组,每组6次试验,每组捉住泥鳅的次数是多少?这时就要用到rbinom(n, size, prob)函数,其中n参数指定试验组数(这里为10000),函数返回每组的成功次数,比如:

> set.seed(12)   ##设置随机数种子,使下面随机结果可重复
> size <- 6   ##每组试验次数
> p <- 0.1    ##成功概率
> n <- 10000  ##试验组数

> ns <- rbinom(n, size, p)   ##每组成功的次数
> table(ns)   ##统计成功次数
ns
   0    1    2    3    4    5 
5286 3512 1040  145   16    1 
 
> mean(ns)  ##成功次数的平均值
[1] 0.6096
> size*p  ##理论值
[1] 0.6

> var(ns)   ##成功次数的方差
[1] 0.5542433
> size*p*(1-p)   ##理论值
[1] 0.54

计算结果显示,当组数足够大,每组成功次数的平均值和方差与理论值很接近。

伯努利分布和二项分布的介绍就到此结束,希望对大家的学习有所帮助,也希望大家多多支持本公众号。

感谢您的阅读!想了解更多有关技巧,请关注我的微信公众号“R语言和Python学堂”,我将定期更新相关文章。

《R统计学(01): 伯努利分布、二项分布》

    原文作者:R语言和Python学堂
    原文地址: https://www.jianshu.com/p/59335680cc29
    本文转自网络文章,转载此文章仅为分享知识,如有侵权,请联系博主进行删除。
点赞