random matrix with constraints:

Hello, I want to generate a random matrix M * N, with two constraints: 1-the sum of the lines is fixed 2-the sum of the columns is fixed. for exmaple: matrix 6*4, sum of lines = 1 and the sum of the 4 columns respectively: 2.43 1.68 1.53 1.18 thank you in advance

[Torsten]

Not possible since 2.43, 1.68, 1.53 and 1.18 must add to 6.

Best wishes

Torsten.

[zelmat mohamed]

I want to fill the matrix in blue, why it's not possible?

[zelmat mohamed]

I correct the sums!

[Torsten]

Any constraints on the matrix elements, e.g. between 0 and 1 ?

[zelmat mohamed]

yes the elements are between 0-1.

Answer 1

[Torsten]

I did not test it, but this code at least promises to do the job:

https://de.mathworks.com/matlabcentral/fileexchange/36070-generate-random-points-in-multi-dimensional-space-subject-to-linear-constraints?focused=5229329&tab=function

If you have access to MATHEMATICA, you can use

http://reference.wolfram.com/language/ref/RandomPoint.html

where the "region" is given by your constraints.

Best wishes

Torsten.

[zelmat mohamed]

Thank you

Answer 2

[John D'Errico]

Suppose I gave you ANY 6x4 matrix. Compute the sum of all of the rows. Then compute the sum of those sums. Clearly the result must be the sum of all elements in the matrix.

So if the row sums are [1 1 1 1 1 1], then the sum of all elements in the matrix MUST be 6.

However, you want the column sums of the matrix to be [2.43 1.68 1.53 1.18].

Just as with the row sums, the sum of the column sums of ANY matrix is also the sum of all of the elements.

sum([2.43 1.68 1.53 1.18])
ans =
         6.82

Since this is not 6, there is no possible way in any universe you will ever contrive, such that your matrix can ever have both row and column sums as you want to see. It is mathematically impossible.

In fact, you can prove mathematically that it is impossible, by knowing nothing more than the simple fact that addition is commutative.

[zelmat mohamed]

thank you, yes I took a bad example, but if we assume that the sum on the line is equal to that on the columns,

 it is possible to generate a random matrix with matlab?

[zelmat mohamed]

I changed the sums of the columns (sum (2,43 1,28 1,53 0,76) = 6) , look at the new matrix, my question is this: is there a code to solve this problem (fill the matrix with random values while respecting the two constraints)

[John D'Errico]

You actually have 24 random variables here, and 10 linear constraints, not TWO constraints. I've actually posted a code on the file exchange, randFixedLinearCombination, that solves a similar problem, subject to ONE general linear constraint. But not 10 linear constraints. So the code I posted will not solve this problem.

You also don't specify what the distribution of points should be. I'd probably assume when someone does not say, but I'd assume most people are thinking about a somewhat uniform set. In fact, it is fairly difficult to generate TRULY uniformly distributed points on such a domain.

A simple iterative rescaling scheme might work. Although it would provably NOT result in uniformly distributed points at all. (Have I said that often enough?)

But here is a little function that does what you want.

% ***************************************************

function [A,iter] = randrowcolconstr(rowsum,colsum,contol)
% generate random quasi-uniformly distributed arrays, subject to both row and column sum constraints
%
% rowsum: row sum constraint vector
% colsum: column sum constraint vector
% contol: constraint tolerance error
%
% Requirement: sum(rowsum) must be the same as sum(colsum)

if abs(sum(rowsum) - sum(colsum)) > contol
  error('Row and column sums are incompatible')
end

nrows = numel(rowsum);
ncols = numel(colsum);

% ensure that rowsum and colsum are vectors in the proper orientation
rowsum = rowsum(:);
colsum = colsum(:).';

% initial random guess
A = rand(nrows,ncols);

iter = 0;
while any(abs(rowsum - sum(A,2)) > contol) || any(abs(colsum - sum(A,1)) > contol)

    % count of iterations, just to see how long it took
    iter = iter + 1;

    % simple iterative rescaling of the rows and columns
    A = A.*(rowsum./sum(A,2));
    A = A.*(colsum./sum(A,1));

end

% ***************************************************

Now, to test it out:

rowsum = [ones(6,1)];
colsum = [2.43,1.28,1.53,0.76];

% tolerance supplied on the constraint error
contol = 1.e-14;

[A,iter] = randrowcolconstr(rowsum,colsum,contol);

A
A =
    0.42407      0.38333      0.11902     0.073579
    0.41896      0.31257      0.21557     0.052901
    0.32853      0.17579     0.076306      0.41938
    0.47763      0.13715      0.34272     0.042504
    0.44073     0.089205      0.45387     0.016198
    0.34009      0.18195      0.32252      0.15544

[rowsum,sum(A,2)]
ans =
          1            1
          1            1
          1            1
          1            1
          1            1
          1            1

[colsum;sum(A,1)]
ans =
       2.43         1.28         1.53         0.76
       2.43         1.28         1.53         0.76

But, are the points uniformly sampled under some meaning of the word uniform? Of course not. This sampling scheme will work. But it will not produce uniformly distributed samples.

Checking the iteration count, it took only 13 iterations to achieve the desired tolerance. That seems to be fairly consistant, taking roughly 13-15 iterations in a few tests.

So the above code is pretty fast, and as long as the rowsum vector sums to the same amount as the colsum vector, it should converge reliably. I did not test that, but here...

sum(rowsum)
ans =
     6

sum(colsum)
ans =
     6

As I said, A is NOT going to be truly uniformly distributed over the domain. It will look moderately uniform at a glance though. The uniformity failure will be that you will tend to see fewer points than expected in the extreme corners of the hyper-dimensional cube. But doing that in a way that is provably uniform will take considerably more effort.

Personally, I'd call the code I posted here a hack, but it will produce something that is probably acceptable, and do so quite efficiently. I do tend to be a perfectionist in my code though.

[zelmat mohamed]

thank you very much, I implemented the function, but during the test it showed me a dimension problem, I did not understand why

[John D'Errico]

You have an old MATLAB release. I used a capability of MATLAB that was introduced R2017. It might work in R2016. It can be repaired using BSXFUN.

[zelmat mohamed]

yes I'm working with the 2016 version, what should I do to solve this problem?

[John D'Errico]

You already know how to solve the problem, you did accept Torsten's answer. Fixing the code I wrote merely requires you learn what bsxfun is, and how to use it.

A = bsxfun(@times,A,rowsum./sum(A,2));

Do something similar for the next line in that code.

[zelmat mohamed]

I admit that your solution is much better :), thank you very much, it worked