Hi,
Roughly reading the code, I find this statement "phat <- x / m" is
probably incorrect since this will give you the set of 1000000 observed x
values /1000000. I redefine the function cover with three inputs: lambda
for the parameter of the poisson distribution, sample.size and
significance.level. The output is 1 or 0, depending on whether lambda is
inside the confidence interval or not. With 5% level of significance I
expect to get 95% of the time the parameter will be included in the
confidence interval. The two runs below shows 96% and 94.8%, pretty close.
> cover <- function(lambda, sample.size, significance.level)
+ {
+ x <- rpois(sample.size,lambda)
+ estimate <- mean(x)
+ lower <- estimate - qnorm(1 - significance.level/2) *
sqrt(estimate/sample.size)
+ upper <- estimate + qnorm(1 - significance.level/2) *
sqrt(estimate/sample.size)
+ if (lambda > lower & lambda < upper){1}else{0}
+ }
> mean(sapply(1:100, function(x)cover(2.5,1000000,0.05)))
[1] 0.96
> mean(sapply(1:1000, function(x)cover(2.5,100000,0.05)))
[1] 0.948
--
View this message in context:
http://r.789695.n4.nabble.com/Coverage-probability-for-a-Poisson-parameter-tp4702535p4702537.html
Sent from the R help mailing list archive at Nabble.com.
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
