R语言上机实验三

理学院实验报告

班级：学号：姓名：实验编号：

实验三：概率和分布的R实现

一、实验目的与要求：

１、会用R给出常见分布的概率密度、概率、分位数和随机数。

２、会利用sample命令进行随机抽样，prod，choose命令计算概率。

3、会利用R绘制各类分布的图形。

4、会利用choose，prod命令计算古典概率。

二、实验内容：

1.从一副扑克牌(52张)中随机抽5张，求下列概率

(1) 抽到的是10，J，Q，K，A；

> 4/choose(52,5)

[1] 1.539077e-06

(2) 抽到的是同花顺。

> 9*choose(4,1)/choose(52,5)

[1] 1.385169e-05

注：同花顺是指5张同一色牌能按从小到大连续排序，如2<3<4<5<6，3<4<5<6<7，…，10

2.模拟随机游动:

(1)从两点分布中产生1000个随机数;

> x<-rbinom(1000,1,0.5)

> x

(2)用函数ifelse( )将上面随机数中的0替换成-1;

> ifelse(x==0,-1,1 )

(3)用函数cumsum( )作出累积和; > y<-ifelse(x==0,-1,1 )

> cumsum(y)

(4)使用命令plot( ) 作出随机游动的示意图. > plot(cumsum(y))

3.在同一个图形中画出统计的四大分布密度曲线(dnorm, dchisq, dt, df)，注意不同分布有不

同的线型、颜色和宽度，还有图形都要在同一方框中，最后用图例说明(legend)。

> curve(dnorm(x,0,1),xlim=c(-1,5),ylim=c(0,0.5),col=1,lwd=1,lty=1)

> curve(dchisq(x,1),xlim=c(-1,5),ylim=c(0,0.5),lwd=2,lty=2,col=2,add=T)

> curve(dt(x,1),xlim=c(0,8),ylim=c(0,0.5),lwd=3,lty=3,col=3,add=T)

> curve(dt(x,1,1),xlim=c(0,8),ylim=c(0,0.5),lwd=4,lty=4,col=4,add=T)

> legend('topright',c("dnorm","dchisp","dt","df"),lty=c(1,2,3,4),col=c(1,2,3,4),lwd=c(1,2,3,4))

> curve(dnorm(x,0,1),xlim=c(-1,5),ylim=c(0,0.5),col=1,lwd=1,lty=1)

> curve(dchisq(x,1),xlim=c(-1,5),ylim=c(0,0.5),lwd=2,lty=8,col=2,add=T)

> curve(dt(x,1),xlim=c(0,8),ylim=c(0,0.5),lwd=5,lty=3,col=7,add=T)

> curve(dt(x,1,1),xlim=c(0,8),ylim=c(0,0.5),lwd=4,lty=4,col=4,add=T)

> legend('topright',c("dnorm","dchisp","dt","df"),lty=c(1,8,3,4),col=c(1,2,7,4),lwd=c(1,2,5,4))

>

4. 除本章给出的标准分布外, 非标准的随机变量X的抽样可通过格式点离散化方法实现.

设p (x )为X 的密度函数, 其抽样步骤如下

(1) 在X 的取值范围内等间隔地选取N 个点x 1, x 2,…, x N , 例如取N =1000; (2) 计算p (x i ); i = 1, 2, …, N ;

(3) 正则化p (x i ); i =1, 2,…,N , 使其成为离散的分布律, 即每一项除以

∑=N

i i

x p 1

)(;

(4) 按离散分布抽样方法使用命令sample( )从x i , i = 1, 2, … ,N 有放回地抽取n 个数, 例如 n =1000.

注：前面４小步是用来编一个函数，功能是对给定的概率密度产生随机数，形式应与ｒｎｏｒｍ差不多。

试以标准正态分布为例来说明. 为与R 中的正态抽样函数rnorm( )进行比较, 将作图区域分为左右两部分,

(1) 使用rnorm( )抽取n =1000个标准正态随机数, 并在左侧区域画出相应的直方图和核

密度估计曲线;

(2) 用格子点离散化抽样方法完成抽样, 并在右侧区域画出相应的直方图和核密度估计

曲线, 离散化所用的N =1000, n =1000, 取点范围为[-4, 4]. > N=seq(-4,4,length=1000) > f<-function(x)

+ dnorm(x)/sum(dnorm(x)) > f1=f(N)

> result=sample(N,replace=T,size=1000,prob=f1) > rn=rnorm(1000)

> op<-par(mfrow=c(1,2)) > hist(rn,probability=T)

> lines(density(rn),col="green",lwd=4) > hist(result,probability=T)

> lines(density(result),col="red",lwd=4) > par(op)