[R] Re : Generate a random bistochastic matrix

[Duncan Murdoch]

On 10/16/2006 11:50 AM, Rolf Turner wrote:
> I don't think this idea has been suggested yet:
> 
> 	(1) Form all n! n x n permutation matrices,
> 	say M_1, ..., M_K, K = n!.
> 
> 	(2) Generate K independent uniform variates
> 	x_1, ..., x_k.
> 
> 	(3) Renormalize these to sum to 1,
> 		x <- x/sum(x)
> 
> 	(4) Form the convex combination
> 
> 		M = x_1*M_1 + ... + x_K*M_K
> 
> M is a ``random'' doubly stochastic matrix.
> 
> The point is that the set of all doubly stochastic matrices
> is a convex set in n^2-dimensional space, and the extreme
> points are the permutation matrices.  I.e. the set of all
> doubly stochastic matrices is the convex hull of the the
> permuation matrices.
> 
> The resulting M will *not* be uniformly distributed on this
> convex hull.  If you want a uniform distribution something
> more is required.  It might be possible to effect uniformity
> of the distribution, but my guess is that it would be a
> hard problem.

But using your solution a Markov chain converging to the uniform 
distribution would be easy:

Start at a bistochastic matrix X, e.g. at a permutation matrix drawn at 
random, or at constant 1/n, or whatever.

Choose a direction D at random by drawing two permutation matrices and 
taking the difference.

Calculate the range of scalar lambda such that X + lambda*D remains 
stochastic (no negatives, no values greater than 1).

Draw lambda uniformly on this range, and move there.

This is a hit-and-run sampler.  I suspect it would mix pretty rapidly. 
I'd guess a Propp-Wilson perfect version would be fairly easy to put 
together.

Duncan Murdoch