# Coverage Probability

6 messages
Open this post in threaded view
|

## Coverage Probability

 Hello. I'm allready this far. I have a function which is calculating the lower (l) and upper (u) limit for a confidence interval for the odds ratio. For example for 5 simulated 2x2 tables the upper and lower limits are: > u [1] 2.496141 7.436524 8.209161 4.313587 3.318612 > l [1] -0.9718608  1.1000713  1.5715373  0.1135158 -0.2700517 With (l[1]; u[1]) being the confidence interval for the odds ratio for the first simulated table and so on. Now I want to compute the coverage probability. For that I've created a function which is return 1 if the odds ratio is in the interval and 0 if it isn't. cover <- function(theta, u, l){ if(theta >= l && theta   <= u){z=1} if(theta < l || theta   > u){z=0}; return(z) } This works but unfortunately not if I want to summarize the function and divide it with the sample size to get the coverage probability. I tried it this way for(for(x in 1:5) {a = (sum(cover(theta, u[x], l[x]))/5; return(a)} Maybe someone can help me. Thank you
Open this post in threaded view
|

## Re: Coverage Probability

 Hi hubinho, You are almost there.  Try this slightly modification of your function: # theta, u and l are vectors of the same length foo <- function(theta, u, l) mean(theta >= l & theta <= u, na.rm = TRUE) foo(theta, u, l) HTH, Jorge.- On Mon, Mar 19, 2012 at 12:55 PM, hubinho <> wrote: > Hello. > > I'm allready this far. I have a function which is calculating the lower (l) > and upper (u) limit for a confidence interval for the odds ratio. > > For example for 5 simulated 2x2 tables the upper and lower limits are: > > > u > [1] 2.496141 7.436524 8.209161 4.313587 3.318612 > > l > [1] -0.9718608  1.1000713  1.5715373  0.1135158 -0.2700517 > > With (l[1]; u[1]) being the confidence interval for the odds ratio for the > first simulated table and so on. > > Now I want to compute the coverage probability. For that I've created a > function which is return 1 if the odds ratio is in the interval and 0 if it > isn't. > > cover <- function(theta, u, l){ > if(theta >= l && theta   <= u){z=1} > if(theta < l || theta   > u){z=0}; return(z) > } > > This works but unfortunately not if I want to summarize the function and > divide it with the sample size to get the coverage probability. > > I tried it this way > > for(for(x in 1:5) {a = (sum(cover(theta, u[x], l[x]))/5; return(a)} > > Maybe someone can help me. Thank you > > > > -- > View this message in context: > http://r.789695.n4.nabble.com/Coverage-Probability-tp4485511p4485511.html> Sent from the R help mailing list archive at Nabble.com. > > ______________________________________________ > [hidden email] 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]] ______________________________________________ [hidden email] mailing list https://stat.ethz.ch/mailman/listinfo/r-helpPLEASE do read the posting guide http://www.R-project.org/posting-guide.htmland provide commented, minimal, self-contained, reproducible code.
Open this post in threaded view
|

## Re: Coverage Probability

 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) }
Open this post in threaded view
|

## Re: Coverage Probability

 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. > > ______________________________________________ > [hidden email] 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]] ______________________________________________ [hidden email] mailing list https://stat.ethz.ch/mailman/listinfo/r-helpPLEASE do read the posting guide http://www.R-project.org/posting-guide.htmland provide commented, minimal, self-contained, reproducible code.