Hi hubinho,
You need to initialize the for() loop and then store the u and l values
properly:
# parameters
n1 <- 10
n2 <- 10
y <- 100
alpha <- 1
z<-1.96
# creating B 2x2 tables
B <- 50
u <- l <- vector('numeric', B)
for (i in 1:B){
theta <- i
x1 <- exp(alpha +theta)/ (1+ exp(alpha +theta))
x2 <- exp(alpha)/ (1+ exp(alpha))
n11 <- rbinom(y, 10, x1)
n12 <- n1 - n11
n21 <- rbinom(y, 10, x2)
n22 <- n2 - n21
# upper and lower limit gart interval
gartu <-function(z,d,e, f, g){log(((d+.5)*(g+.5))/((e+.5)*(f+.5)))+
z*sqrt(1/(d+.5)+1/(e+.5)+1/(f+.5)+1/(g+.5))}
gartl <-function(z,d,e, f, g){log(((d+.5)*(g+.5))/((e+.5)*(f+.5)))-
z*sqrt(1/(d+.5)+1/(e+.5)+1/(f+.5)+1/(g+.5))}
# store results
u[i] <- gartu(z, n11[i],n22[i],n12[i],n21[i])
l[i] <- gartl(z, n11[i],n22[i],n12[i],n21[i])
}
# coverage
theta <- 1:B
foo <- function(theta, u, l) mean(theta >= l & theta <= u, na.rm = TRUE)
foo(theta, u, l)
# [1] 0.14
HTH,
Jorge.-
On Mon, Mar 19, 2012 at 2:25 PM, hubinho <> wrote:
> Thank you very much. This was, was i needed. Unfortunately I have one
> futher
> problem with this Code. I don't only need the coverage probability for one
> but for a range of different odds ratios. (for example [1;30]). I tried it
> with a loop but I get an error. I think again, that I'm almost there but
> having a little mistake. The complete code is:
>
> #setting values
>
> n1 <- 10
> n2 <- 10
> y <- 100
> alpha <- 1
> z<-1.96
>
> # creating 2x2 table
>
> for (i in 1:30)
>
> {
>
> theta <- i
> x1 <- exp(alpha +theta)/ (1+ exp(alpha +theta))
> x2 <- exp(alpha)/ (1+ exp(alpha))
>
>
> n11 <- rbinom(y, 10, x1)
> n12 <- n1 - n11
> n21 <- rbinom(y, 10, x2)
> n22 <- n2 - n21
>
> # upper and lower limit gart interval
>
> gartu <-function(z,d,e, f, g){log(((d+.5)*(g+.5))/((e+.5)*(f+.5)))+
> z*sqrt(1/(d+.5)+1/(e+.5)+1/(f+.5)+1/(g+.5))}
> gartl <-function(z,d,e, f, g){log(((d+.5)*(g+.5))/((e+.5)*(f+.5)))-
> z*sqrt(1/(d+.5)+1/(e+.5)+1/(f+.5)+1/(g+.5))}
>
>
> u <- gartu(z, n11[i],n22[i],n12[i],n21[i])
> l <- gartl(z, n11[i],n22[i],n12[i],n21[i])
>
> foo <- function(theta, u, l) mean(theta >= l & theta <= u, na.rm = TRUE)
> foo(theta, u, l)
> }
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/Coverage-Probability-tp4485511p4485865.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help@r-project.org mailing list
> 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.
>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
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.
