R语言学习笔记4_参数估计

最新推荐文章于 2022-04-20 20:08:43 发布

skye学习的小窝

最新推荐文章于 2022-04-20 20:08:43 发布

阅读量5.4k

点赞数 9

分类专栏： R

本文链接：https://blog.csdn.net/skye_s_/article/details/119650412

版权

R 专栏收录该内容

9 篇文章 23 订阅

订阅专栏

本文介绍了统计学中的参数估计方法，包括矩估计和极大似然估计。通过实例展示了如何使用R语言进行单参数和多参数的估计，并探讨了正态总体参数的区间估计，如均值和方差的置信区间。同时，讨论了两正态总体参数的区间估计，如均值差和方差比的置信区间。此外，还涉及了样本容量的确定方法。

摘要由CSDN通过智能技术生成

四、参数估计

根据样本推断总体的分布和分布的数字特征称为统计推断。
参数估计有两类，一类是点估计，以某个统计量的样本观测值作为未知参数的估计值；另一类是区间估计，用两个统计量所构造的区间来估计位置参数【给出了估计的可信度】。

4.1 矩估计和极大似然估计法

4.1.1 矩估计

若总体X的k阶矩存在，则样本的k阶矩依概率收敛到总体的k阶矩，样本矩的连续函数收敛到总体矩的连续函数----->用样本矩作为总体矩的估计量

矩估计可能是不唯一的，通常采用低阶矩给出未知参数的估计
在总体分布未知的情况下，也可以用样本均值估计总体均值，用样本方差估计总体方差
没有固定的R程序求出矩估计，可利用R的计算功能根据具体问题编写相应的R程序

例1：通常事件的成败机会比 g(θ)=θ/1-θ 是人们感兴趣的参数。对某个篮球运动员记录其一次在比赛中投篮命中与否，观测数据如下：
1 1 0 1 0 0 1 0 1 1 1 0 1 1 0 1
0 0 1 0 1 0 1 0 0 1 1 0 1 1 0 1
编写相应的R函数估计这个篮球运动员投篮的成败比

> x<-c(1,1,0,1,0,0,1,0,1,1,1,0,1,1,0,1,0,0,1,0,1,0,1,0,0,1,1,0,1,1,0,1)
> theta<-mean(x)
> t<-theta/(1-theta)
> t
[1] 1.286

例2：下面的观测值为来自指数分布的一个样本，估计参数λ
0.17834 0.33181 1.20810 0.08954 0.33990 0.68148 0.02528 0.34818 1.20790 2.62448

> x<-c(0.17834,0.33181, 1.20810, 0.08954, 0.33990, 0.68148, 0.02528, 0.34818, 1.20790, 2.62448)
> lambda<-1/mean(x)
> lambda
[1] 1.421
使用二阶矩进行矩估计：
> lambda<-1/sd(x)
> lambda
[1] 1.256

实际上，上面的数据是模拟参数为2的指数分布，可见低阶矩更为精确。

4.1.2 极大似然估计

单参数 optimize( )

optimize(f = , interval = , lower = min(interval), upper = max(interval), maximum = TRUE, tol = .Machine$double.eps^0.25, ...)

f 是似然函数，interval 是参数θ的取值范围，lower 是θ的下界，upper 是θ的上界，maximum=T是求极大值，tol表示求值的精度，省略号是对f 的附加说明

多参数 optim( ) 、nlm( )

optim(par, fn, gr = NULL, method = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN", "Brent"), lower = -Inf, upper = Inf,  control = list( ), hessian = FALSE,...)

函数nlm( )仅使用牛顿-拉夫逊算法求函数的最小值点；函数optim( )提供method选项中的六种方法中的一种进行优化。

nlm(f, p, ..., hessian = FALSE, typsize = rep(1, length(p)), fscale = 1, print.level = 0, ndigit = 12, gradtol = 1e-6,  stepmax = max(1000 * sqrt(sum((p/typsize)^2)), 1000),  steptol = 1e-6, iterlim = 100, check.analyticals = TRUE)

例：一地质学家为研究密歇根湖的湖滩地区的岩石成分，随机地自该地区取出100个样品，每个样品有十块石子，他记录了每个样品中属石灰石的石子数，得到的数据如下：

样品中的石子数 0 1 2 3 4 5 6 7 8 9 10
样品个数 0 1 6 7 23 26 21 12 3 1 2

假设这100次观测相互独立，求这地区石子中的石灰石比例p的极大似然估计

样品中的石子数	0	1	2	3	4	5	6	7	8	9	10
样品个数	0	1	6	7	23	26	21	12	3	1	2

显然，每个样品中的石子数服从二项分布binom(10,p),下面根据100次观测估计参数p:
> f <- function(p)(p^517)*(1-p)^483   #似然函数
> optimize(f,c(0,1),maximum = T)
$maximum
[1] 0.517
$objective
[1] 1.664e-301
因此，该地区石子中石灰石的比例p的最大似然估计为0.517

4.2 单正态总体参数的区间估计

4.2.1 均值μ的区间估计

1）方差σ²已知时μ的置信区间

> # 求方差已知时均值的置信区间
> z.test<- function(x,n,sigma,alpha,u0=0,alternative="two.sided"){
+   options(digits = 4)
+   result<- list( )
+   mean<- mean(x)
+   z<- (mean-u0)/(sigma/sqrt(n))
+   p<- pnorm(z,lower.tail = F)
+   #把计算结果放到resul里
+   result$mean<-mean
+   result$z<-z
+   result$p.value<-p
+   #假设检验
+   if(alternative=="two.sided")
+     result$p.value<- 2*pnorm(abs(z),lower.tail = F)
+   else if (alternative=="greater")
+     result$p.value<- pnorm(z)
+   #求置信区间
+   result$conf.int<-c(
+     mean-sigma*qnorm(1-alpha/2,mean=0,sd=1,lower.tail = T)/sqrt(n),
+     mean+sigma*qnorm(1-alpha/2,mean=0,sd=1,lower.tail = T)/sqrt(n)
+   )
+   result
+ }

例：一个人10次称自己的体重（单位500g）：175,176,173,175,174,173,173,176,173,179，假设此人体重服从正态分布，标准差为1.5，求体重的置信水平为95%的置信区间。

> x<-c(175,176,173,175,174,173,173,176,173,179)
> result<-z.test(x,10,1.5,0.05)
> result$conf.int
[1] 173.8 175.6

2）方差σ²未知时μ的置信区间

t.test(x, y = NULL, alternative = c("two.sided","less","greater"), mu=0, paired = F, var.equal = F, conf.level = 0,95, ...)

若仅出现数据x，进行单样本t检验;否则进行二样本t检验。
alternative=“two.sided”是缺省值（默认），表示求置信区间；alternative="less"表示求置信上限；alternative="greater"表示求置信下限。
mu表示均值，仅在假设检验中起作用。

在上例中如果不知道方差，就需要用函数t.test()来求置信区间
> x<-c(175,176,173,175,174,173,173,176,173,179)
> t.test(x)$conf.int
[1] 173.3 176.1

4.2.2 方差σ²的区间估计

#卡方检验:方差的置信区间
chisq.var.test<- function(x,var,alpha,alternative="two.sided"){
  options(digitis=4)
  results <- list( )
  n<- length(x)
  v<- var(x)
  result$var<- v
  chi2<-(n-1)*v/var
  result$chi2<-chi2
  p<- pchisq(chi2,n-1)
  result$p.value <-p
  if(alternative=="less")
    result$p.value<-pchaisq(chi2,n-1,lower.tail=F)
  else if(alternative=="two.sided")
    result$p.value<- 2*min(pchisq(chi2,n-1),pchisq(chi2,n-1,lower.tail = F))
  result$conf.int<- c(
    (n-1)*v/qchisq(alpha/2,df=n-1,lower.tail = F),
    (n-1)*v/qchisq(alpha/2,df=n-1,lower.tail = T))
  result
}

4.3 两正态总体参数的区间估计

4.3.1 均值差μ1-μ2的置信区间

1）两方差都已知时两均值差的置信区间

> #两正态总体均值差的区间估计（方差已知）
> two.sample.ci<- function(x,y,conf.level=0.95,sigma1,sigma2){
+   options(digits = 4)
+   m= length(x);n=length(y)
+   xbar=mean(x)-mean(y)
+   alpha=1-conf.level
+   zstar=qnorm(1-alpha/2)*(sigma1/m+sigma2/n)^(1/2)
+   xbar+c(-zstar,+zstar)
+ }

2）两方差都未知但相等时两均值差的置信区间
利用t.test(x,y,var.equal=TRUE)可求

4.3.2 两方差比的置信区间

var.test(x, y, ratio=1, alternative=c("two.sided","less","greater"), conf.level=0.95, ...)

4.4 单总体比率p的区间估计

prop.test(x, n, p=NULL, alternative=c("two.sided","less","greater"), conf.level=0.95, correct=TRUE)

correct=TRUE是否做连续型矫正，不矫正的区间长度＜矫正后的区间长度

binom.test(x, n, p=NULL, alternative=c("two.sided","less","greater"), conf.level=0.95)

例：从一份共有3042人的人名录中随机抽200人，发现38人的地址已变动，试以95%的置信水平，估计这份名录中需要修改地址的比例。

用正态分布来近似
> prop.test(38,200,correct = TRUE)

	1-sample proportions test with continuity correction

data:  38 out of 200, null probability 0.5
X-squared = 76, df = 1, p-value <2e-16
alternative hypothesis: true p is not equal to 0.5
95 percent confidence interval:
 0.1395 0.2527
sample estimates:
   p 
0.19 
以95%的置信水平认为这份名录中需要修改地址的比例p落在（0.1395，0.2527）中，点估计为0.19

用二项分布来近似
> binom.test(38,200)

	Exact binomial test

data:  38 and 200
number of successes = 38, number of trials = 200, p-value <2e-16
alternative hypothesis: true probability of success is not equal to 0.5
95 percent confidence interval:
 0.1381 0.2513
sample estimates:
probability of success 
                  0.19

4.5 两总体比率差p1-p2的区间估计

例：据一项市场调查，在A地区被调查的1000人中有478人喜欢品牌K，在B地区被调查的750人中有246人喜欢品牌K，试估计两地区人们喜欢品牌K比例差的95%置信区间。

> like<-c(478,246)
> people<-c(1000,750)
> prop.test(like,people)

	2-sample test for equality of proportions with continuity correction

data:  like out of people
X-squared = 39, df = 1, p-value = 4e-10
alternative hypothesis: two.sided
95 percent confidence interval:
 0.1031 0.1969
sample estimates:
prop 1 prop 2 
 0.478  0.328

可以看出，A地区喜欢品牌K的人更多，且A、B两地区喜欢品牌K的比例之差的95%的置信区间为（0.1031，0.1969）

4.6 样本容量的确定

4.6.1 估计正态总体均值时样本容量的确定

1）总体方差σ²已知

size.norm1<- function(d,var,conf.level){
	alpha = 1-conf.level
	((qnorm(1-alpha/2)*var^(1/2))/d)^2
}

d是允许的最大绝对误差
2）总体方差σ²未知

size.norm2<- function(s,alpha,d,m){
	t0<- qt(alpha/2,m,lower.tail=FALSE)
	n0<- (t0*s/d)^2
	t1<- qt(alpha/2,n0,lower.tail=FALSE)
	n1<- (t1*s/d)^2
	while(abs(n1-n0)>0.5){
	n0<- (qt(alpha/2,n1,lower.tail=FALSE)*s/d)^2
	n1<- (qt(alpha/2,n0,lower.tail=FALSE)*s/d)^2
	}
	n1
}

m是事先给定的一个很大的数

4.6.2 估计比例 p时样本容量的确定

size.bin<- function(d,p,conf.level=0.95){
	alpha=1-conf.level
	((qnorm(1-alpha/2))/d)^2*p*(1-p)
}

例：某市一所重点大学历届毕业生就业率为90%，试估计应届毕业生就业率，要求估计误差不超过3%，试问在α=0.05下要抽取应届毕业生多少人？

> size.bin(0.03,0.9,0.95)
[1] 384.1

skye学习的小窝

关注

9
点赞
踩
52

收藏

觉得还不错? 一键收藏
0
评论
R语言学习笔记4_参数估计

目录四、参数估计4.1 矩估计和极大似然估计法4.1.1 矩估计4.1.2 极大似然估计单参数 optimize( )多参数 optim( ) 、nlm( )4.2 单正态总体参数的区间估计4.2.1 均值μ的区间估计4.2.2 方差σ^2^的区间估计4.3 两正态总体参数的区间估计4.3.1 均值差μ1-μ2的置信区间4.3.2 两方差比的置信区间4.4 单总体比率p的区间估计4.5 两总体比率差p1-p2的区间估计4.6 样本容量的确定4.6.1 估计正态总体均值时样本容量的确定4.6.2 估计比例 p
复制链接

扫一扫