《统计分析软件：使用R与Python》 课件 第5章- R 语言随机抽样和随机数

5.1随机变量分布5.2随机抽样5.3生成已知分布的随机数5.4随机数的应用第5章R语言随机抽样和随机数R语言随机抽样和随机数随机变量分布

离散变量、连续变量及其分布生成已知分布的随机数利用现有R函数、逆变换法、舍选法随机数的应用随机抽样简单随机抽样、分层抽样、整群抽样5.1随机变量分布数据的数字特征刻划了数据的主要特征,而要对数据的总体情况作全面的描述,就要研究数据的分布.常见的离散型分布有:两点分布、二项分布和Poisson分布等.常见的连续型分布均匀分布、正态分布、指数分布、χ2分布、t分布和F分布等.在R软件中,提供了计算典型分布的分布函数、分布律或概率密度函数,以及分布函数的反函数的函数.分布R名称自变量binomialB(n,p)binom(size,prob)Poissonπ(λ)pois(lambda)GeometricGeo(p)geom(prob)hypergeometrichyper(N,m,n)

H(N,m,n)negativebinomialnbinom(size,prob)NB(r,p)multinomialmultinom(size,prob)离散型随机变量及其分布uniformU(a,b)unif(min,max)Gaussian(normal)N(µ,

2)norm(mean,sd)ExponentialE(λ)

exp(rate)gammaГ(λ,β)gamma(shape,scale)BetaBe((a,b)beta(shape1,shape2,ncp)‘Student’

t(n)t(df,ncp)Fisher–SnedecorF(n1,n2)f(df1,df2,ncp)Pearsonχ2(n)chisq(df,ncp)Weibullrweibull(shape,scale)Cauchycauchy(location,scale)InversegammaInvgamma(shape,scale)Logisticlogis(location,scale)Log-normallnorm(meanlog,sdlog)分布R函数附加参数连续型随机变量及其分布上述各类表示分布的函数前分别加d——概率函数（分布列或概率密度）(distribution)p——分布函数(probability)q——分位数(quantile)r——取随机数(random)pfunc(x,...)qfunc(p,...)

0<p<1)rfunc(n,p1,p2,...)

dfunc(x,...)R中概率分布的四种状态#X~B(20,0.2),P(X=5)

>dbinom(5,20,0.2)[1]0.1745595#X~π(5),P(X=10)>dpois(10,5)[1]0.01813279随机变量取某个值的概率#X~B(100,0.2),P(46≤X≤54)sum(dbinom(46:54,100,0.2)pbinom(54,100,0.2)-pbinom(45,100,0.2)#X~B(20,0.2),P(X≤5)pbinom(5,20,0.2)#X~B(20,0.2),P(X>5)

1-pbinom(5,20,0.2)#X~π(5),P(0≤X≤16)

sum(dpois(0:16,5));ppois(16,5)#X~N(0,1),P(X≤2)

pnorm(2)#X~N(170,62),P(X≤184)pnorm(184,170,6)随机变量在某区间内取值概率#标准正态的下0.95分位数.qnorm(0.95)#标准正态的下0.025分位数.qnorm(0.025)#自由度为4的卡方分布下0.05分位数.qchisq(0.05,4)#自由度为5的t分布下0.975分位数.

qt(0.975,5)求分位数-R中给的是下分位数求分位数在收集数据过程中，绝大多数情况下，并不采取普查的方式获取总体中所有样本的数据信息，而是以各类抽样方法抽取其中若干代表性样本来进行数据获取和分析。在获得待分析数据集后，需要再次通过抽样技术选取出训练集和测试集，以便比较选择出最优的统计模型。这里主要介绍简单随机抽样、分层抽样、整群抽样三种基本抽样方法。简单随机抽样，R包自带；分层抽样，整群抽样，在sampling包中找到相应方法实现。5.2随机抽样1.简单随机抽样R语言简单随机抽样的函数是sample,句法是：

sample(x,size,replace=FALSE,prob=NULL)

x:带抽取对象size:想要抽取的样本数量replace:是否为有放回，默认为FALSE无放回prob:设置个抽取样本的抽样概率，默认为无取值，即等概率抽样以JS数据集为例#有放回抽取10个数据进行分析set.seed(1235)sampleindex=sample(nrow(JS),10,replace=T)sampleindex[1]2712337143628133512>JS[sampleindex,]以JS数据集为例#无放回抽取10个数据进行分析>sampleindex=sample(nrow(JS),10)>sampleindex##[1]271231428133853323>JS[sampleindex,]这里，sampleindex没有重复值。无放回时，size必须小于样本数sampleindex=sample(nrow(JS),50)Errorin(x,size,replace,prob):cannottakeasamplelargerthanthepopulationwhen'replace=FALSE'2.分层抽样分层抽样，也叫类型抽样。就是将总体单位按其属性特征分成若干类型或层，然后在类型或层中随机抽取样本单位。strata(data,stratanames=NULL,size,method=c(“srswor”,“srswr”,“poisson”,“systematic”),pik,description=FALSE)data:带抽样数据stratanames:进行分层所依据的变量名称size:各层中要抽出的观测样本数method:选择4中抽样方法，分别为无放回、有放回、泊松、系统抽样，默认为srsworpik:设置各层中样本的抽样概率description:选择是否输出含有各层基本信息的结果以JS数据集为例#载入分层抽样的包samplinglibrary(sampling)#按地区，每个职称抽取3人>strata(JS,stratanames=(“title"),size=c(3,3,3,3),method="srswor")

3.整群抽样整群抽样又称聚类抽样。是将总体中各单位归并成若干个互不交叉、互不重复的集合，称之为群；然后以群为抽样单位抽取样本的一种抽样方式。cluster(data,clustername,size,method=c("srswor","srswr","poisson","systematic"),pik,description=FALSE)data:带抽样数据clustername:用来划分群的变量名称size:需要抽取的群数method:选择4中抽样方法，分别为无放回、有放回、泊松、系统抽样，默认为srsworpik:设置各层中样本的抽样概率description:选择是否输出含有各群的基本信息的结果以JS数据集为例#载入分层抽样的包samplinglibrary(sampling)#以title为分群变量，#不放回简单随机抽样#方法抽2个群>cluster(JS,clustername=(“title"),size=2,method="srswor",description=T)

1.服从均匀分布的随机数R语言生成均匀分布随机数的函数是runif（）句法是：runif(n,min=0,max=1)

n表示生成的随机数数量，min表示均匀分布的下限，max表示均匀分布的上限；若省略参数min、max,则默认生成[0,1]上的均匀分布随机数。>runif(5,0,1)

#生成5个[0,1]的均匀分布的随机数[1]0.59930.73910.26170.50770.7199

>runif(5)

#默认生成5个[0,1]上的均匀分布随机数[1]0.27840.77550.41070.83920.7455

5.3生成已知分布的随机数5.3.1

R语言函数生成随机数2.服从正态分布的随机数正态分布随机数的生成函数是rnorm（）句法是：rnorm（n,mean=0,sd=1）

其中n表示生成的随机数数量，mean是正态分布的均值，默认为0，sd是正态分布的标准差，默认时为1;#生成一个长度为100的向量，向量中的每个值服从标准正态分布rnorm(100)#产生100个N(150,62)的随机数rnorm(100,150,6)3.服从二项分布的随机数二项分布是指n次独立重复贝努力试验成功的次数的分布，每次贝努力试验的结果只有两个，成功和失败，记成功的概率为p生成二项分布随机数的函数是：rbinom（）句法是：rbinom(n,size,prob)

n表示生成的随机数数量，size表示进行贝努力试验的次数，prob表示一次贝努力试验成功的概率#产生100个B(20,0.22)的随机数rbinom(100,20,0.2)R其他分布的随机数rexp指数分布rfF分布rgamma

Gamma分布rgeom几何分布rhyper超几何分布rlnorm对数正态分布rlogisLogistic分布rmultinom多项分布rnbinom负二项分布rpois泊松分布rtt分布rchisq卡方分布生成随机变量的方法有许多种，但是几乎所有产生随机变量的技术都是从生成一个或多个[0,1]区间均匀分布的伪随机数开始的，随后通过应用某种转换方法即可从[0,1]均匀分布随机数中生成非均匀分布的随机数。方法：逆变换法、卷积法、函数变换法、合成法、取舍法等5.3.2

逆变换法生成随机数逆变换法逆变换法也称反函数法，是最常用、最简单的一种随机变量生成方法。基于概率积分变换定理

连续型随机数逆变换法的生成步骤例如，给定某指数分布的概率密度函数为:

注意到u和1-u同为均匀分布，因此产生n个参数为lambda的指数分布的R命令可以写成-log(runif(n))/lambdaset.seed(123)n=10lambda=3-log(runif(n))/lambda##[1]0.415420940.079290010.298032180.041470120.02045947

##[6]1.029600660.212819740.037939830.198410430.26130509>n<-1000>u<-runif(n)>x<-u^(1/3)>hist(x,prob=TRUE,main=expression(f(x)==3*x^2))>y<-seq(0,1,0.01)>lines(y,3*y^2)逆变换法仅适用于生成分布函数具有反函数的随机数（如均匀分布、指数分布等），而不适用于生成类似于正态分布这种分布函数不具有反函数的随机数，同样也不适用于生成离散分布函数的随机数。5.3.3舍选法由于大多数分布函数无法求反函数，因此逆变换方法的适用范围有限，为此，学者们提出了舍选法（接受-拒绝抽样，Acceptance-Rejectionsampling）。假设需要对分布π(x)进行采样，但是却很难直接进行，故借助另外一个容易采样的分布g(x)，用某种机制去除掉一些g(x)生成的样本，从而使得剩下的样本就是来自于所求分布π(x)的样本。舍选法的基本思想所需分布π(x)易抽样分布g(x)

舍选法的原理证明

逆变换法和舍选法小结设置随机种子，保证实验的可重复性

>x<-rnorm(10)#随机生成10个随机数>x[1]0.27387940-0.29359610-0.640519080.865472080.68148647-0.33057645[7]-1.08217114-1.658157360.78442397-0.01521926>x<-rnorm(10)#再次随机生成10个随机数>x[1]0.0634884060.302965807-1.851315056-0.048666815-0.261443363-1.413083258[7]-0.4164286291.0477428531.318716887-0.003966177>set.seed(10)#设定种子>x<-rnorm(10)#随机生成10个随机数>x[1]0.01874617-0.18425254-1.37133055-0.599167720.294545130.38979430[7]-1.20807618-0.36367602-1.62667268-0.25647839>set.seed(10)#设定种子>y<-rnorm(10)#再次随机生成10个随机数>y[1]0.01874617-0.18425254-1.37133055-0.599167720.294545130.38979430[7]-1.20807618-0.36367602-1.62667268-0.25647839

从上面的例子可以发现，前两次生成的10个随机数都不一样。后两次在设定了相同的种子前提下，生成的随机数是相同的。5.4随机数的应用随机模拟的方法生成指定数据，用来估计参数或验证统计学中的一些基本结论。估计模型参数验证大数定律验证中心极限定理估计参数

我们可以利用上述结论来估计圆周率。

N=2000set.seed(1234)x<-runif(N,-1,1)y<-runif(N,-1,1)p<-mean((x^2+y^2)<=1)pi_est<-4*ppi_est##[1]3.14为了验证大数定律，选择t分布进行验证#自由度为10的t分布n=20000a<-rt(n,10)llnum<-function(a){y<-rep(0,length(a))for(iin1:n){y[i]<-mean(sample(a,i,replace=TRUE))}data<-data.frame(size=1:n,value=y)data}lln_data<-llnum(a)colnames(lln_data)<-c("sample_size","sample_mean")library(ggplot2)ggplot(lln_data,aes(x=sample_size,y=sample_mean))+labs(title="样本均值变化趋势")+geom_line(color='darkblue')+geom_abline(intercept=5,slope=0,color='grey')

例：用R验证中心极限定理#First,wewillgenerate1000samplesandcomputethesamplemeanofeach.set.seed(1)sample_means

=

rep(NA,1000)

for(i

in

1:1000){

sample_means[i]

=

mean(rexp(40,0.2))

}

均值：这里用R语言对指数分布的均值进行模拟仿真。并将结果和中心极限定理的结果进行对比。#comparethesamplemeanfr