【数理统计实验(一)】统计量近似分布的随机模拟

简介: 【数理统计实验(一)】统计量近似分布的随机模拟

该篇文章以实例的形式演示了利用R语言完成关于统计量近似分布的随机模拟:不同分布的正态性检验(正态分布、泊松分布);似然估计量的渐进正态过程;方差估计量的优劣性比较。

1 模拟实验一:不同分布的正态性检验

   背景:总体分别为均匀分布,指数分布,两点分布(二阶矩存在)等情况下,利用随机模拟的方法画出样本均值的直方图,并进行正态性检验,给出p值

   要求:(1)至少选两个总体,一个离散型,一个连续型。

   (2)报告体现需出样本容量不同时,样本均值的趋近于正态分布的过程。

   函数定义:

fenbuhanshu<- function (r,distpar,e,s,n,N)#r-随机数函数,distpar-参数取值范围,e-均值,s-标准差,n-样本数,N-抽取样本重复次数               #定义函数
{  for (i in n) {
  if (length(distpar)==2)
  {x <- matrix(r(i*N, distpar[1],distpar[2]),nc=i)}
  else 
  {x <- matrix(r(i*N, distpar), nc=i)}
  x <- (apply(x, 1, sum) - i*e )/(sqrt(i)*s)
  p=shapiro.test(x)$p.value                       #正态性检验
  hist(x,col='green',probability=T,main=paste("n=",i,",p=",round(p,4)),ylim=c(0,max(.4, density(x)$y)))
  lines(density(x), col='red', lwd=3)}}

   (1)连续型$N \sim U(0,1) $。

   运行程序:

op=par(mfrow=c(2,2))
fenbuhanshu(runif,c(0,1),0.5,1/sqrt(12),c(3,30,300,500),1000)  #分别做3,30,300,500次实验
par(op)

   运行结果:

   (2)离散型$N \sim P(1) $

   运行程序:

op=par(mfrow=c(2,2))
limite.central(rpois,1,1,1,c(3,30,300,500),1000)
par(op)

   运行结果:

2 模拟实验二:似然估计量的渐进正态过程

   背景:总体为:

image.png

  运行程序:

t=par(mfrow=c(2,2))                               #两行两列画布
f=function(theta){
  logL=m*log(theta)+(theta-1)*sum(log(x))
  return(logL)}                                      #返回联合密度函数的ln值
n=c(5,50,500,5000)                                   #分别模拟5,50,500,1000次                                 
for (m in n){
  a=c()
  for (i in 1:1000){
    y=runif(m,min=0,max=1)
    x=sqrt(y)
    a[i]=1/(optimize(f,c(0,100),maximum = TRUE)$maximum)}  #求最大似然函数值
  p=shapiro.test(a)$p.value                        #Shapiro-Wilk方法进行正态检验
  hist(a,col='green',probability=T,main=paste("n=",m,",p=",round(p,4)),ylim=c(0,max(.4, density(a)$y)))
  lines(density(a), col='red', lwd=3)}
par(t)

   运行结果:

   由图可知,似然估计量随着样本越来越大有着渐进正态的过程。

3 模拟实验三:方差估计量的优劣性比较

   image.png

   利用R随机模拟方法,比较三个估计量分别在样本容量不同的情况下的优劣。给出表格和图形。

   运行程序:

library(ggplot2)
s1=c()
s2=c()
s3=c()
op=par(mfrow=c(3,1))                               #绘制3行1列画布
for (i in c(5,50,100)){                            #样本容量分别设置为5,50,100
  for (j in 1:100){
    x=rnorm(i, mean=0,sd=1)                     #随机选取均值为0,方差为1的样本                  
    y1=sum((x-mean(x))^2)/(i-1)
    s1[j]=y1
    y2=sum(x^2)/i
    s2[j]=y2
    y3=sum(x^2)/(i+2)
    s3[j]=y3}
  print(mean(s1))
  print(mean(s1)-1)
  print(mean(s2))
  print(mean(s2)-1)
  print(mean(s3))
  print(mean(s3)-1)
  plot(density(s1),ylim=c(0,3),main=paste("n=",i))  
  lines(density(s2),col='red')
  lines(density(s3),col='green')}                     #分别绘制3根密度曲线
par(op)

   结果:

   由表格及结果图可知, image.png 最接近方差,即为最优估计量。


相关文章
|
4月前
|
算法 数据挖掘
WinBUGS对多元随机波动率SV模型:贝叶斯估计与模型比较
WinBUGS对多元随机波动率SV模型:贝叶斯估计与模型比较
|
2月前
|
资源调度 数据可视化 算法
贝叶斯统计是一种基于贝叶斯定理的统计学方法,它不同于传统的频率派统计(或称为经典统计)。
贝叶斯统计是一种基于贝叶斯定理的统计学方法,它不同于传统的频率派统计(或称为经典统计)。
|
4月前
R语言进行数值模拟:模拟泊松回归模型的数据
R语言进行数值模拟:模拟泊松回归模型的数据
R语言进行数值模拟:模拟泊松回归模型的数据
|
4月前
|
数据可视化
R语言极值理论:希尔HILL统计量尾部指数参数估计可视化
R语言极值理论:希尔HILL统计量尾部指数参数估计可视化
|
4月前
|
算法 vr&ar Python
R语言隐马尔可夫模型HMM连续序列重要性重抽样CSIR估计随机波动率模型SV分析股票收益率时间序列
R语言隐马尔可夫模型HMM连续序列重要性重抽样CSIR估计随机波动率模型SV分析股票收益率时间序列
|
4月前
R语言中固定与随机效应Meta分析 - 效率和置信区间覆盖
R语言中固定与随机效应Meta分析 - 效率和置信区间覆盖
|
4月前
|
数据可视化 C语言
使用R语言随机波动模型SV处理时间序列中的随机波动率
使用R语言随机波动模型SV处理时间序列中的随机波动率
|
4月前
|
Python 数据可视化 索引
PYTHON用GARCH、离散随机波动率模型DSV模拟估计股票收益时间序列与蒙特卡洛可视化
PYTHON用GARCH、离散随机波动率模型DSV模拟估计股票收益时间序列与蒙特卡洛可视化
PYTHON用GARCH、离散随机波动率模型DSV模拟估计股票收益时间序列与蒙特卡洛可视化
|
4月前
|
数据可视化 算法 区块链
R语言泊松过程及在随机模拟应用可视化
R语言泊松过程及在随机模拟应用可视化
|
4月前
R语言用GAM广义相加模型研究公交专用道对行程时间变异度数据的影响
R语言用GAM广义相加模型研究公交专用道对行程时间变异度数据的影响