From: <oehl_list_at_gmx.de>

Date: Tue 26 Oct 2004 - 21:38:21 EST

language R

Date: Tue 26 Oct 2004 - 21:38:21 EST

Best regards

Jens

K <- 100 J <- 60 N <- K+J p <- K/N n <- 50

nn <- 100000

urn <- rep(0:1, c(J,K))

x <- sapply(1:nn, function(i){

sum(sample(urn, n))

})

y <- rhyper(nn, K,J, n)

require(SuppDists)

z <- rghyper(nn, a=K, k=n, N=N)

# hypergeometric mean and variance

p*n

p*(1-p)*n*(N-n)/(N-1)

# check we have parametrized ghyper correctly sghyper(a=K, k=n, N=N)[c("Mean", "Variance")]

var(x) var(y) # wrong var(z)

version

*> K <- 100
**> J <- 60
**> N <- K+J
**> p <- K/N
**> n <- 50
**>
**> nn <- 100000
**>
*

> urn <- rep(0:1, c(J,K))

> x <- sapply(1:nn, function(i){

+ sum(sample(urn, n))

+ })

*> y <- rhyper(nn, K,J, n)
**>
*

> require(SuppDists)

**[1] TRUE
**

*> z <- rghyper(nn, a=K, k=n, N=N)
**>
**>
*

> # hypergeometric mean and variance

*> p*n
*

[1] 31.25

> p*(1-p)*n*(N-n)/(N-1)

[1] 8.107311

*>
*

> # check we have parametrized ghyper correctly

> sghyper(a=K, k=n, N=N)[c("Mean", "Variance")]

$Mean

[1] 31.25

$Variance

[1] 8.107311

*>
*

> var(x)

[1] 8.060095

> var(y) # wrong

[1] 8.61289

> var(z)

[1] 8.150472

*>
*

> version

_

platform i386-pc-mingw32

arch i386 os mingw32 system i386, mingw32 status major 2 minor 0.0 year 2004 month 10 day 04

language R

-- ______________________________________________ R-devel@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-develReceived on Tue Oct 26 21:46:58 2004

*
This archive was generated by hypermail 2.1.8
: Fri 18 Mar 2005 - 09:00:56 EST
*