algorithm - MATLAB: Fast creation of random symmetric Matrix with fixed degree (sum of rows) -

i searching method create, in fast way random matrix a follwing properties:

• a = transpose(a)
• a(i,i) = 0 i
• a(i,j) >= 0 i, j
• sum(a) =~ degree; sum of rows randomly distributed distribution want specify (here =~ means approximate equality).

the distribution degree comes matrix orig, degree=sum(orig), know matrices distribution exist.

for example: orig=[0 12 7 5; 12 0 1 9; 7 1 0 3; 5 9 3 0]

orig =  0    12     7     5 12     0     1     9  7     1     0     3  5     9     3     0   sum(orig)=[24 22 11 17]; 

now 1 possible matrix a=[0 11 5 8, 11 0 4 7, 5 4 0 2, 8 7 2 0]

a =   0    11     5     8 11     0     4     7  5     4     0     2  8     7     2     0 

with sum(a)=[24 22 11 17].

i trying quite time, unfortunatly 2 ideas didn't work:

---

version 1:

i switch nswitch times 2 random elements: a(k1,k3)--; a(k1,k4)++; a(k2,k3)++; a(k2,k4)--; (the transposed elements aswell).

unfortunatly, nswitch = log(e)*e (with e=sum(sum(nn))) in order matrices uncorrelated. e > 5.000.000, not feasible (in particular, need @ least 10 of such matrices).

---

version 2:

i create matrix according distribution scratch. idea is, fill every row i degree(i) numbers, based on distribution of degree:

nn=orig; nnr=zeros(size(nn)); i=1:length(nn)     degree=sum(nn);     howmany=degree(i);     degree(i)=0;     full=rld_cumsum(degree,1:length(degree));     rr=randi(length(full),[1,howmany]);     ff=full(rr);     xx=i*ones([1,length(ff)]);     nnr = nnr + accumarray([xx(:),ff(:)],1,size(nnr)); end a=nnr; 

however, while sum(a')=degree, sum(a) systematically deviates degree, , not able find reason that.

small deviations degree fine of course, there seem systmatical deviations in particulat of matrices contain in places large numbers.

---

i happy if either show me fast method version1, or reason systematic deviation of distribution in version 2, or method create such matrices in different way. thank you!

---

edit:

this problem in matsmath's proposed solution: imagine have matrix:

orig =   0    12     3     1 12     0     1     9  3     1     0     3  1     9     3     0 

with r(i)=[16 22 7 13].

• step 1: r(1)=16, random integer partition p(i)=[0 7 3 6].
• step 2: check p(i)<=r(i), case.
• step 3:

my random matrix starts looks like

a =   0    7     3     6  7    0     .     .  3    .     0     .  6    .     .     0 

with new row sum vector rnew=[r(2)-p(2),...,r(n)-p(n)]=[15 4 7]

second iteration (here problem occures):

• step 1: rnew(1)=15, random integer partition p(i)=[0 b]: rnew(1)=15=a+b.
• step 2: check p(i)<=rnew(i), gives a<=4, b<=7. a+b<=11, a+b has 15. contradiction :-/

---

edit2:

this code representing (to best of knowledge) solution posted david eisenstat:

orig=[0 12 3 1; 12 0 1 9; 3 1 0 3; 1 9 3 0]; w=[2.2406 4.6334 0.8174 1.6902]; xfull=zeros(4);  ii=1:1000     rndmat=[poissrnd(w(1),1,4); poissrnd(w(2),1,4); poissrnd(w(3),1,4); poissrnd(w(4),1,4)];     kkk=rndmat.*(ones(4)-eye(4)); % remove diagonal     hhh=sum(sum(orig))/sum(sum(kkk))*kkk; % normalisation     xfull=xfull+hhh; end  xf=xfull/ii; disp(sum(orig)); % gives [16    22     7    13] disp(sum(xf));   % gives [14.8337    9.6171   18.0627   15.4865] (obvious systematic problem) disp(sum(xf'))   % gives [13.5230   28.8452    4.9635   10.6683] (which systematically different [16, 22, 7, 13] 

since it's enough approximately preserve degree sequence, let me propose random distribution each entry above diagonal chosen according poisson distribution. intuition want find weights w_i such i,j entry i != j has mean w_i*w_j (all of diagonal entries zero). gives nonlinear system of equations:

for i, (sum_{j != i} w_i*w_j) = d_i, 

where d_i degree of i. equivalently,

for i, w_i * (sum_j w_j) - w_i^2 = d_i. 

the latter can solved applying newton's method described below starting solution of w_i = d_i / sqrt(sum_j d_j).

once have w_is, can sample repeatedly using poissrnd generate samples of multiple poisson distributions @ once.

---

(if have time, i'll try implementing in numpy.)

the jacobian matrix of equation system 4 4 problem is

(w_2 + w_3 + w_4) w_1               w_1               w_1 w_2               (w_1 + w_3 + w_4) w_2               w_2 w_3               w_3               (w_1 + w_2 + w_4) w_3 w_4               w_4               w_4               (w_1 + w_2 + w_3). 

in general, let a diagonal matrix a_{i,i} = sum_j w_j - 2*w_i. let u = [w_1, ..., w_n]' , v = [1, ..., 1]'. jacobian can written j = + u*v'. inverse given sherman--morrison formula

                              a^-1*u*v'*a^-1 j^-1 = (a + u*v')^-1 = a^-1 - -------------- .                               1 + v'*a^-1*u 

for newton step, need compute j^-1*y given y. can done straightforwardly in time o(n) using above equation. i'll add more detail when chance.