https://github.com/cran/simecol
Tip revision: 27972f79606ed158414d982ce17b7f8dbe5bb6e9 authored by Thomas Petzoldt on 07 October 2021, 05:40:02 UTC
version 0.8-14
version 0.8-14
Tip revision: 27972f7
pcuseries.R
`alpha2rho` <-
function(alpha) {
(alpha^2 -1 - 2* alpha * log(alpha)) / (alpha - 1) ^2
}
`rho2alpha` <-
function(rho) {
if (abs(rho)>1) stop("rho must be within interval (0, 1)")
f <- function(alpha, rho) {
(rho - alpha2rho(alpha))^2
}
if (rho==1) {
Inf ## must be handled in subsequent procedures
} else if (rho==-1) {
0
} else if (rho==0) {
1
} else {
mm <- optimize(f, interval=c(1e-12, 1e12), rho=rho)
mm$minimum
}
}
`pcu` <-
function(x, alpha = rho2alpha(rho), rho){
if (is.infinite(alpha)) alpha <- 1e99
n <- length(x)
y <- runif(n)
b <- alpha+(alpha-1)*(alpha-1)*y*(1-y)
c <- 2*y*(1-y)*(alpha*alpha*x-x+1)-2*alpha*y*(1-y)+alpha
d <- sqrt(alpha*alpha+4*alpha*(1-alpha)*(1-alpha)*x*y*(1-x)*(1-y))
y <- (c-(1-2*y)*d)/(2*b)
}
`pcuseries` <-
function(n, alpha = rho2alpha(rho), rho, min=0, max=1) {
x<-numeric(n)
x[1] <- runif(1)
for(i in 2:n){
x[i] <- pcu(x[i-1], alpha)
}
x * (max - min) + min
}