I need to generate a list of random variables $\bf{x}$ subject to constraints that can be expressed in the form $\bf{E}x=b$ where $\bf{E}$ is an $m \times n $ matrix if $\bf{x}$ has $n$ entries. In all the cases I'm dealing with, $n >> m$, for example $n$ will be around 14,000 and $m$ will be 50. I'm not sure what method I will use for random sampling, either normal or uniform, it is not clear which is best for the problem I'm trying to solve, but I need all the variables to be sampled from distributions with the same mean and range/variance.

What I've been doing to solve this is reducing $\bf{E}$ to row-echelon form, setting all the variables corresponding to columns to the right of the last pivot to random values, and then solving the remaining square matrix equality.

There is a problem however, to solve the square matrix equality, I subtract the already set values from the right hand side. Unfortunately, the variances add as well, so my last 50 values tend to vary hugely, which is unfortunately unacceptable in this problem.

Is there a better way to do this? I cannot think of a way to fix the current method I am using. I use R.

  • 2
    $\begingroup$ Unfortunately, you won't be able to do this, unless you get really lucky with your constraint matrix. For example, consider a constraint matrix with two rows, one nonzero entry in the first row that constrains $x_1 = 0$ and two nonzero entries in the second row that constrain $x_2+x_3 = 1$. Obviously, $x_1$ will have a different mean than at least one of $x_2$ and $x_3$, and unless you set the variances for $x_2$ and $x_3 = 0$, a different variance as well. $\endgroup$
    – jbowman
    Jun 16, 2012 at 13:50

2 Answers 2


This paper and R package completely solved my problem. It uses the Markov Chain Monte Carlo method, which relies on the fact that if you can find an initial solution of the constraint, through linear programming, you can find an arbitrary number of them by using a matrix that when multiplied by $E$, the constraints, gives zero. Read about it here:


and here is the package:



Might seem trivial (and not terribly machine efficient), but consider repeating the process until you get a suitable answer? Preferably only modifying a smaller subset each time.

Can you create a "distance" measure for how far you are away from your ideal answer? It might help you "optimize"?

  • $\begingroup$ I might try that out. One problem is that I would still need to make sure the shape of the distribution would be the same. Also time constraints will be prohibitive. $\endgroup$
    – Mike Flynn
    Jun 17, 2012 at 7:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.