# probability of bivariate max (error?)

a guest
Mar 11th, 2020
147
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. r_to_theta <- function (r) { 2/pi * asin(r) }
2. round(digits=3, sapply(seq(0,0.5, by=0.05), function(r) { r_to_theta(r) }))
3. # [1] 0.000 0.032 0.064 0.096 0.128 0.161 0.194 0.228 0.262 0.297 0.333
4.
5. p_max_bivariate <- function(n,r) {
6.   ## for the special-case r=0-0.5:
7.   if (r > 0.5 || r < 0.0) { stop('Formula valid only for 0 <= r <= 0.5') }
8.     else {
9.         theta <- r_to_theta(r)
10.         t <- theta / (theta - 1)
11.         library(gsl) # for hyperg_2F1
12.
13.         t * (hyperg_2F1(1, n+1, n+2, t) / (n+1)) -
14.          2 * n * t^2 * (hyperg_2F1(1, n+2, n+3, t) / ((n+1)*(n+2))) -
15.          (2*n*t)/((n+1)^2) +
16.          (1/n)
17.         }
18. }
19. r=0.5; sapply(2:10, function(n) { p_max_bivariate(n, r) })
20. # [1] 0.540620144 0.378139568 0.292465038 0.238966583 0.202214466 0.175349622 0.154830414 0.138634576 0.125520185