How to create a function that returns a p-value min(1, 2*min(P{N ≤ i}, P{N ≥ i})) in R?

user_13543233 0
Viewed 46

I am trying to create a function that performs the sign test. The function should take in vector x and median m and returns a p-value defined by min(1, 2*min(P{N ≤ i}, P{N ≥ i})), where N~Bin(n, 0.5), n is the amount of values which are not equal to m, and i is the amount of values in n which are smaller than m.

Here's what I've gotten so far:

test<-function(x,m) {
  n<-length(x[x!=m])
  i<-length(n[n<m])
  min(1,c(pbinom(q=i,size=n,p=0.5),1-pbinom(q=i,size=n,p=0.5)))
}

However, when I test it out with given values, it gives the wrong answer:

test(x=1:3, m=2)
[1] 0.25

test(x = 1:5, m = 2) 
[1] 0.0625

The right answers should be 1 and 0.625. I don't know which step I have done wrong.

p-value min r function
1 Answers

You're on the right track: Aside from two minor errors you made (forgot the factor 2, problem in your calculation of i: Note that in your script, n is an integer, not a vector!), there's a more subtle but important error in your calculation of the third term inside the min:

We have 1 - P(N<=i) = P(N<i) which is different from P(N<=i) for discrete-valued probability distributions (like the binomial distribution).

You can fix this by subtracting a tiny value eps from i in the last term and using the approximation P(N <= i-eps) ≈ P(N < i):

test <- function(x , m){
  # amount of entries that are not m
  n <- sum(x!=m)
  # amoutn of entries that are smaller than m
  i <- sum(x<m)
  
  eps = 1e-10

  min(1,
      2*pbinom(q=i,size=n,p=0.5),
      # to get P(N >= i), we need 1 - P(N < I), not 1 - P(N <= i) 
      2*(1-pbinom(q=i-eps,size=n,p=0.5))
  )
}
user_47559700
Related