[37] | 1 | % Here is code -- I have added a few comments to help you. |
---|
| 2 | % N is dimension of square matrix, M is rank we want. alpha is number |
---|
| 3 | % between 0 and 1 that selects entries that are don't care. |
---|
| 4 | % Program generates random [0,1] NxN matrix A1. Z is generated from another |
---|
| 5 | % similar random matrix whose elements are then set to one where entries |
---|
| 6 | % exceed alpha, and zero otherwise. One entries determine the don't cares. |
---|
| 7 | % Hence, about 1-alpha entries are don't care. A is set to A1 in the care |
---|
| 8 | % states, zero elsewhere. Initial random entries are chosen for the |
---|
| 9 | % don't care states in the outer for loop in range -50,+50. |
---|
| 10 | % |
---|
| 11 | % Inner for loop is really to prevent an infinite loop (crummy programming). |
---|
| 12 | % The real test is the break whent he error change is small. That's also a |
---|
| 13 | % hack for now. |
---|
| 14 | % |
---|
| 15 | % XX accumulates sets of all final values of don't care values. |
---|
| 16 | |
---|
| 17 | N = 5; M = 2; alpha = 0.9; NUMSTEPS = 1000; |
---|
| 18 | A1 = rand(N,N); |
---|
| 19 | Z = max(0,sign(rand(N,N) - alpha)); |
---|
| 20 | A = A1 - A1.*Z; |
---|
| 21 | Zs = reshape(Z,1,N^2); |
---|
| 22 | index = find(Zs == 1) |
---|
| 23 | count = size(index,2); |
---|
| 24 | X = zeros(N,N); |
---|
| 25 | XX = []; |
---|
| 26 | STEPS = []; |
---|
| 27 | |
---|
| 28 | for i = 1:20 |
---|
| 29 | |
---|
| 30 | Xs = reshape(X,1,N^2); |
---|
| 31 | Xf = Xs(index); |
---|
| 32 | Xf = 100*(rand(1,count)-0.5); |
---|
| 33 | Xs(index) = Xf; |
---|
| 34 | X = reshape(Xs,N,N); |
---|
| 35 | % X holds random initial starting conditions for unknown values. |
---|
| 36 | XfA = []; |
---|
| 37 | |
---|
| 38 | for j = 1:NUMSTEPS |
---|
| 39 | B = A + X; |
---|
| 40 | % B is the matrix |
---|
| 41 | [U,D,V] = svd(B); |
---|
| 42 | BA = V(:,1:M)*D(1:M,1:M)*(V(:,1:M))'; |
---|
| 43 | % BA is the best fit rank M matrix to B. |
---|
| 44 | X = BA.*Z; |
---|
| 45 | % X is updated with missing values from this fit. |
---|
| 46 | |
---|
| 47 | Xs = reshape(X,1,N^2); |
---|
| 48 | Xfold =Xf; |
---|
| 49 | Xf = Xs(index); |
---|
| 50 | XfA = [XfA; Xf]; |
---|
| 51 | Error = sqrt(sum((Xf-Xfold).^2)); |
---|
| 52 | % Error is change in magnitude of X. |
---|
| 53 | |
---|
| 54 | A - (BA - X) |
---|
| 55 | if Error < 1e-10 |
---|
| 56 | break |
---|
| 57 | end |
---|
| 58 | end |
---|
| 59 | STEPS = [STEPS; j]; |
---|
| 60 | XX = [XX;Xf]; |
---|
| 61 | end |
---|
| 62 | % XX |
---|
| 63 | STEPS |
---|
| 64 | XX(find(STEPS ~= NUMSTEPS), :) |
---|