Thanks, I tried someting like this, but computation takes times for large
matrices
btransf <- function(y,X=length(y)^4) {
N <- length(y)
bm <- matrix(rep(1/N,N^2),N,N)
for(j in 1:X){
coord <- sample(1:N,4,replace=T)
d <-
runif(1,0,min(bm[coord[1],coord[2]],bm[coord[3],coord[4]]))
bm[coord[1],coord[2]] <- bm[coord[1],coord[2]]-d
bm[coord[3],coord[4]] <- bm[coord[3],coord[4]]-d
bm[coord[1],coord[4]] <- bm[coord[1],coord[4]]+d
bm[coord[2],coord[3]] <- bm[coord[2],coord[3]]+d }
y.btransf <- bm%*%y
y.btransf <- y.btransf+(mean(y)-mean(y.btransf))
as.vector(y.btransf) }
the fonction is designed to perform a mean-preserving transformation of a
vector.
----- Message d'origine ----
De : Richard M. Heiberger <[EMAIL PROTECTED]>
À : Florent Bresson <[EMAIL PROTECTED]>; r-help@stat.math.ethz.ch
Envoyé le : Lundi, 16 Octobre 2006, 14h58mn 13s
Objet : Re: [R] Generate a random bistochastic matrix
bistochastic.3x3 <- function() {
B <- matrix(0, 3, 3)
## 2 df
tmp.1 <- runif(3)
B[1,] <- tmp.1/sum(tmp.1)
## 1 df
tmp.2 <- runif(2)
B[2:3, 1] <- (1-B[1,1]) * tmp.2/sum(tmp.2)
## 1 df
B[2, 2] <- runif(1, max=min(1-B[1,2], 1-B[2,1]))
## Fill in the rest
B[2,3] <- 1-sum(B[2, 1:2])
B[3,2] <- 1-sum(B[1:2, 2])
B[3,3] <- 1-sum(B[1:2, 3])
B
}
B <- bistochastic.3x3()
apply(B, 1, sum)
apply(B, 2, sum)
To extend this to larger than 3x3 requires the same kind of
conditional generation of alternating rows and columns of the
matrix. The hard part is the extension of the two-way conditioning
I illustrated in the B[2, 2] line.
Rich
______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
