MATLAB: Generating random numbers from 0 – 1 with limit on the sum

Question

Hi,

I am trying to generate values from 0 to 1 for several variables (x1 x2 x3).

These results must assure that x1*a1 + x2*a2 + x3*a3 <= D. a1 a2 a3 and D are other variables I define. How can I do this?

Thank you very much,

Answer 1

Ok, lets talk about the problem as defined. It is possible to solve the problem, but things are not always simple in high dimensions. In a high number of dimensions, the solution may be highly computationally intensive.

The set of numbers typically generated by rand are uniform, in 3 dimensions or in n dimensions. In n-dimensions, we can think of the points as living in a hypercube, an n-cube.

But if you require that the weighted sum

x1*a1 + x2*a2 + x3*a3 <= D

be limited to some constant D, then you actually want a UNIFORM random sample that lives in a truncated cube, but with one part of the cube sliced away, as if eaten by rats. :)

So, consider the problem in just 3 dimensions. I'll pick some arbitrary values for a1, a2, a3, D to see what happens here. For example, suppose that we chose the a's (and D) as:

a = [2;3;4];
D = 4;

I'm using a vector here, because numbered variables are just SO bad a programming solution. I won't do it. LEARN TO USE VECTORS. Were this a 10-dimensional problem, this solution would be far more readable and writable. In fact, as you see below, to change the solution to 10 dimensions from 3 dimensions, the change was trivial.

I've also used a column vector for a. So if X is a vector, of length 3 that comes from rand:

X = rand(1,3);

Well, actually, lets generate a large sample, say 1000 such samples, in 3 dimensions.

N =1000;
X = rand(N,3);

So now each row of X may be interpreted as one point in the unit 3-cube. We could plot them in 3-d, but many of them are outside of the linear constraint plane. In fact, only 353 of the 1000 samples we generated were within the constraint.

sum(X*a < D)
ans =
   353

Time for a picture.

ind = X*a < D;
plot3(X(ind,1),X(ind,2),X(ind,3),'ro')
hold on
plot3(X(~ind,1),X(~ind,2),X(~ind,3),'g+')
box on
view(42.5,19.6)

Ok, so it looks like roughly 35% of the 1000 points lie underneath the constraint plane. In fact, if D were a little smaller (in this case with these values for a, no larger then 2) then the region of interest would be described by a tetrahedron. But I chose D=4 above. In fact, it is possible to create the region of interest as a solid in 3 dimensions.

So it is a somewhat complex solid. Not just a simple tetrahedron. In fact, the code that I used to generate the picture created 14 different tetrahedrae in the simplicial complex that describes the volume. You can see the edges of them in that plot, where they touch the surface. Ok, better, using transparency...

So a bit messy.

Regardless ... the simple solution to this problem, to generate perhaps 1000 points that lie in the required region, is to just generate too many random samples, then throw away those that lie beyond the constraint.

So, lets suppose that we desire to sample N points, randomly and uniformly from the above volume, so that no point exceeds the required limit D. Simple (sort of) is a rejection sampling scheme:

a = [2;3;4];
D = 4;
N = 1000;
Xfinal = zeros(N,3);
Nfound = 0;
Nneed = N - Nfound;
while Nfound < N
  X = rand(N,3);
  % rejection step
  ind = find(X*a <= D);
  if ~isempty(ind)
    if numel(ind) > Nneed
      ind = ind(1:Nneed);
    end
    Xfinal(Nfound + (1:numel(ind)),:) = X(ind,:);
    Nfound = Nfound + numel(ind);
    Nneed = N - Nfound;
  end
end

Here, I've plotted the final set of 1000 points, with the required volume overlaid in transparent form to see that indeed, they do satisfy the constraint.

The nice thing is, the above solution is rather easily converted to work in a higher number of dimensions. The downside is that in a higher number of dimensions, the rejection step may be come slow. The curse of dimensionality nails you here. There are other, more sophisticated solutions that can be employed, but even there, they get nasty in a high number of dimensions.