Nothing
### R code from vignette source 'fullpres.Rnw'
###################################################
### code chunk number 1: wh.predict
###################################################
wh.predict <- function(x)
{
M1 <- 30269
M2 <- 30307
M3 <- 30323
y <- round(M1*M2*M3*x)
s1 <- y %% M1
s2 <- y %% M2
s3 <- y %% M3
s1 <- (171*26478*s1) %% M1
s2 <- (172*26070*s2) %% M2
s3 <- (170*8037*s3) %% M3
(s1/M1 + s2/M2 + s3/M3) %% 1
}
RNGkind("Wichmann-Hill")
xnew <- runif(1)
maxerr <- 0
for (i in 1:1000) {
xold <- xnew
xnew <- runif(1)
err <- abs(wh.predict(xold) - xnew)
maxerr <- max(err, maxerr)
}
print(maxerr)
###################################################
### code chunk number 2: congrurand (eval = FALSE)
###################################################
## congruRand(10)
###################################################
### code chunk number 3: congrurand2
###################################################
library(randtoolbox)
options(width = 40)
congruRand(10)
###################################################
### code chunk number 4: congrurandseed1 (eval = FALSE)
###################################################
## setSeed(1)
## congruRand(10)
###################################################
### code chunk number 5: congrurandseed2
###################################################
options( width =40)
setSeed(1)
congruRand(10)
###################################################
### code chunk number 6: congrurandseed3 (eval = FALSE)
###################################################
## setSeed(1)
## congruRand(10, echo=TRUE)
###################################################
### code chunk number 7: congrurandseed4
###################################################
options( width =40)
setSeed(1)
congruRand(10, echo=TRUE)
###################################################
### code chunk number 8: congrurandseed5 (eval = FALSE)
###################################################
## setSeed(1614852353)
## congruRand(5, echo=TRUE)
###################################################
### code chunk number 9: congrurandseed6
###################################################
options( width =40)
setSeed(1614852353)
congruRand(5, echo=TRUE)
###################################################
### code chunk number 10: congrurandseed5 (eval = FALSE)
###################################################
## setSeed(12)
## congruRand(5, mod = 2^8, mult = 25, incr = 16, echo= TRUE)
###################################################
### code chunk number 11: congrurandseed6
###################################################
options( width =30)
setSeed(12)
congruRand(5, mod = 2^8, mult = 25, incr = 16, echo= TRUE)
###################################################
### code chunk number 12: sfmt (eval = FALSE)
###################################################
## SFMT(10)
## SFMT(5, 2) #bi dimensional variates
###################################################
### code chunk number 13: sfmt2
###################################################
options( width =40)
SFMT(10)
SFMT(5, 2)
###################################################
### code chunk number 14: sfmt3 (eval = FALSE)
###################################################
## SFMT(10, mexp = 607)
###################################################
### code chunk number 15: sfmt4
###################################################
options( width =40)
SFMT(10, mexp = 607)
###################################################
### code chunk number 16: halton (eval = FALSE)
###################################################
## halton(10)
## halton(10, 2)
###################################################
### code chunk number 17: halton2
###################################################
options( width =40)
halton(10)
halton(10, 2)
###################################################
### code chunk number 18: halton3 (eval = FALSE)
###################################################
## halton(5)
## halton(5, init=FALSE)
###################################################
### code chunk number 19: halton4
###################################################
options( width =40)
halton(5)
halton(5, init=FALSE)
###################################################
### code chunk number 20: sobol (eval = FALSE)
###################################################
## sobol(10)
###################################################
### code chunk number 21: sobol2
###################################################
options( width =40)
sobol(10)
###################################################
### code chunk number 22: unitsquare1 (eval = FALSE)
###################################################
## par(mfrow = c(1,2))
## plot(sobol(1000, 2))
## plot(sobol(10^3, 2, scrambling=1))
###################################################
### code chunk number 23: unitsquare2
###################################################
par(mfrow = c(1,2))
plot(sobol(1000, 2), xlab ="u", ylab="v", main="Sobol (no scrambling)")
plot(sobol(10^3, 2, scrambling=1), xlab ="u", ylab="v", main="Sobol (Owen)")
###################################################
### code chunk number 24: torus (eval = FALSE)
###################################################
## torus(10)
###################################################
### code chunk number 25: torus2
###################################################
options( width =40)
torus(10)
###################################################
### code chunk number 26: torus3 (eval = FALSE)
###################################################
## torus(5, use =TRUE)
###################################################
### code chunk number 27: torus4
###################################################
options( width =40)
torus(5, use =TRUE)
###################################################
### code chunk number 28: torus5 (eval = FALSE)
###################################################
## torus(5, p =7)
###################################################
### code chunk number 29: torus6
###################################################
options( width =40)
torus(5, p =7)
###################################################
### code chunk number 30: torus7 (eval = FALSE)
###################################################
## torus(5, mixed =TRUE)
###################################################
### code chunk number 31: torus8
###################################################
options( width =40)
torus(5, mixed =TRUE)
###################################################
### code chunk number 32: torus9 (eval = FALSE)
###################################################
## par(mfrow = c(1,2))
## acf(torus(10^5))
## acf(torus(10^5, mix=TRUE))
###################################################
### code chunk number 33: torusacf
###################################################
par(mfrow = c(1,2))
acf(torus(10^5))
acf(torus(10^5, mix=TRUE))
###################################################
### code chunk number 34: unitsquare3 (eval = FALSE)
###################################################
## par(mfrow = c(1,2))
## plot(SFMT(1000, 2))
## plot(torus(10^3, 2))
###################################################
### code chunk number 35: unitsquare4
###################################################
par(mfrow = c(1,2))
plot(SFMT(1000, 2), xlab ="u", ylab="v", main="SFMT")
plot(torus(1000, 2), xlab ="u", ylab="v", main="Torus")
###################################################
### code chunk number 36: unitsquare5 (eval = FALSE)
###################################################
## par(mfrow = c(1,2))
## plot(WELL(1000, 2))
## plot(sobol(10^3, 2, scrambling=2))
###################################################
### code chunk number 37: unitsquare6
###################################################
par(mfrow = c(1,2))
plot(WELL(1000, 2), xlab ="u", ylab="v", main="WELL 512a")
plot(sobol(10^3, 2, scrambling=2), xlab ="u", ylab="v", main="Sobol (Faure-Tezuka)")
###################################################
### code chunk number 38: integralcos (eval = FALSE)
###################################################
## I25 <- -1356914
## nb <- c(1200, 14500, 214000)
## ans <- NULL
## for(i in 1:3)
## {
## tij <- sobol(nb[i], dim=25,
## scrambling=0, normal=TRUE)
## Icos <- sqrt(rowSums(tij^2/2))
## Icos <- mean(cos(Icos)) * pi^(25/2)
## ans <- rbind(ans, c(n=nb[i],
## I25=Icos, Delta=(Icos-I25)/I25))
## }
## data.frame(ans)
###################################################
### code chunk number 39: integralcos2
###################################################
options( width =40)
I25 <- -1356914
nb <- c(1200, 14500, 214000)
ans <- NULL
for(i in 1:3)
{
tij <- sobol(nb[i], dim=25, scramb=0, norm=TRUE )
Icos <- mean(cos(sqrt( apply( tij^2/2, 1, sum ) ))) * pi^(25/2)
ans <- rbind(ans, c(n=nb[i], I25=Icos, Delta=(Icos-I25)/I25 ))
}
data.frame(ans)
###################################################
### code chunk number 40: hist (eval = FALSE)
###################################################
## par(mfrow = c(1,2))
## hist(SFMT(10^3), 100)
## hist(torus(10^3), 100)
###################################################
### code chunk number 41: hist
###################################################
par(mfrow = c(1,2))
hist(SFMT(10^3), 100)
hist(torus(10^3), 100)
###################################################
### code chunk number 42: gap1 (eval = FALSE)
###################################################
## gap.test(runif(1000))
###################################################
### code chunk number 43: gap2
###################################################
options( width =40)
res <- gap.test(runif(1000), echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Gap test\n")
cat("\nchisq stat = ", stat, ", df = ",df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n(sample size : ",1000,")\n\n", sep="")
cat("length observed freq theoretical freq\n")
for(i in 1:(df+1))
cat(i,"\t", obsnum[i],"\t", expnum[i],"\n")
###################################################
### code chunk number 44: gap3 (eval = FALSE)
###################################################
## gap.test(SFMT(1000), 1/3, 2/3)
###################################################
### code chunk number 45: order1 (eval = FALSE)
###################################################
## order.test(runif(4000), d=4)
###################################################
### code chunk number 46: order2
###################################################
options( width =40)
res <- order.test(runif(4000), d=4, echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Order test\n")
cat("\nchisq stat = ", stat, ", df = ",df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n (sample size : ", 1000,")\n\n", sep="")
cat("observed number ",obsnum[1:6],"\n",obsnum[7:18],"\n", obsnum[19:24],"\n")
cat("expected number ",expnum,"\n")
###################################################
### code chunk number 47: freq (eval = FALSE)
###################################################
## freq.test(runif(1000), 1:4)
###################################################
### code chunk number 48: freq2
###################################################
options( width =40)
res <- freq.test(runif(1000), 1:4, echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Frequency test\n")
cat("\nchisq stat = ", stat, ", df = ",df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n (sample size : ",1000,")\n\n", sep="")
cat("observed number ",obsnum,"\n")
cat("expected number ",expnum,"\n")
###################################################
### code chunk number 49: serial (eval = FALSE)
###################################################
## serial.test(runif(3000), 3)
###################################################
### code chunk number 50: serial2
###################################################
options( width =40)
res <- serial.test(runif(3000), 3, echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Serial test\n")
cat("\nchisq stat = ", stat, ", df = ",df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n (sample size : ",3000,")\n\n", sep="")
cat("observed number ",obsnum[1:4],"\n", obsnum[5:9])
cat("expected number ",expnum,"\n")
###################################################
### code chunk number 51: coll1 (eval = FALSE)
###################################################
## coll.test(runif, 2^7, 2^10, 1)
###################################################
### code chunk number 52: coll2
###################################################
options( width =40)
res <- coll.test(runif, 2^7, 2^10, 1, echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Collision test\n")
cat("\nchisq stat = ", stat, ", df = ", df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n exact distribution \n(sample number : ", 1000,"/sample size : ", 128,"\n / cell number : ", 1024,")\n\n", sep="")
cat("collision observed expected\n", "number count count\n", sep="")
for(i in 1:(df + 1) )
cat(" ", i," ", obsnum[i]," ", expnum[i],"\n")
###################################################
### code chunk number 53: coll3 (eval = FALSE)
###################################################
## coll.test(congruRand, 2^8, 2^14, 1)
###################################################
### code chunk number 54: coll4
###################################################
options( width =40)
res <- coll.test(congruRand, 2^8, 2^14, 1, echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Collision test\n")
cat("\nchisq stat = ", stat, ", df = ", df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n Poisson approximation \n(sample number : ", 1000,"/sample size : ", 256,"\n / cell number : ", 16384,")\n\n", sep="")
cat("collision observed expected\n", "number count count\n", sep="")
for(i in 1:(df + 1) )
cat(" ", i-1," ", obsnum[i]," ", expnum[i],"\n")
###################################################
### code chunk number 55: poker (eval = FALSE)
###################################################
## poker.test(SFMT(10000))
###################################################
### code chunk number 56: poker2
###################################################
options( width =40)
res <- poker.test(SFMT(10000), echo=FALSE)
stat <- res$statistic
pvalue <- res$p.value
df <- res$parameter
obsnum <- res$observed
expnum <- res$expected
cat("\n\t Poker test\n")
cat("\nchisq stat = ", stat, ", df = ", df, "\n, p-value = ", pvalue, "\n", sep="")
cat("\n (sample size : ", 10000,")\n\n", sep="")
cat("observed number ", obsnum,"\n")
cat("expected number ", expnum,"\n")
###################################################
### code chunk number 57: wh1
###################################################
wh.predict <- function(x)
{
M1 <- 30269
M2 <- 30307
M3 <- 30323
y <- round(M1*M2*M3*x)
s1 <- y %% M1
s2 <- y %% M2
s3 <- y %% M3
s1 <- (171*26478*s1) %% M1
s2 <- (172*26070*s2) %% M2
s3 <- (170*8037*s3) %% M3
(s1/M1 + s2/M2 + s3/M3) %% 1
}
RNGkind("Wichmann-Hill")
xnew <- runif(1)
err <- 0
for (i in 1:1000)
{
xold <- xnew
xnew <- runif(1)
err <- max(err, abs(wh.predict(xold)-xnew))
}
print(err)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.