Create a nullvector to have a nxn-matrix to calculate MWB

[cleokatt]

Hi I try to calculate a MWB; In mathematics, particularly matrix theory, a band matrix or banded matrix is a sparse matrix whose non-zero entries are confined to a diagonal band, comprising the main diagonal and zero or more diagonals on either side.

And it requires a nxn-matrix ,but my data looks like this:  And I think the best way to calculate it when the matrix don't have the same dimensions is to add a nullvector, in my case when Mk1=scale1 

data pseudo_data;
input SCALE $ mk2 mk3 mk4 mk5 mk6 mk7;
datalines;
scale1 0.620073 0.585710 0.742480 0.654074 0.152773 0.828327
scale2 0.030940 0.081121 0.505554 0.495800 0.024913 0.874568
scale3 0.493812 0.464746 0.675302 0.900584 0.240976 0.112287
scale4 0.194194 0.483504 0.650082 0.795761 0.112740 0.145614
scale5 0.427196 0.986343 0.778883 0.631170 0.623107 0.586194
scale6 0.037071 0.375161 0.217043 0.893637 0.779019 0.161867
scale7 0.186883 0.830461 0.038186 0.000778 0.336085 0.586100
;
run;

/* ex rclassorgcount dataset */
data rclassorgcount;
input SOMETHING_THAT;
datalines;
100
200
150
250
300
350
400
;
run;

%let rcmigrdist=mk2 mk3 mk4 mk5 mk6 mk7; /* # of mk's*/
%let K=7; /*scale*/
%let K_Max=10; /*how many mk's it could possibly be depending on the dataset you running*/
%let SOMETHING_THAT=SOMETHING_THAT; /*something else in my real code; count(SOMETHING) as SOMETHING_THAT*/

proc iml;
    use pseudo_data;
        read all var {&rcmigrdist.} into rc_migr_dist; /* Load matrix with movements shares p_ij */
    close pseudo_data;

    use rclassorgcount;
        read all var {&SOMETHING_THAT.} into rc_orig_count_vec; /* Load vector with number of humans per movement i, Ni */
    close rclassorgcount;

    /* Calc Matrix_normal_u */
    M_norm_u = 0;
    do i = 1 to &K. - 1;
        sum_p = 0;
        maxval = max(abs(i - &K.), abs(i + (&K_Max. - &K.) - 1));
        N_i = rc_orig_count_vec[i];
        do j = i + 1 to &K.;
            sum_p = sum_p + (rc_migr_dist[i, j] / 100); /* Assuming percentages are given, convert to decimal */
        end;
        /*print sum_p;*/
        temp = maxval * N_i * sum_p;
        M_norm_u = M_norm_u + temp;
    end;
    create M_norm_u var {"Matrix_normal_u"}; 
    append from M_norm_u;
    close M_norm_u;
quit;

And you see that the log says: 

153 end;
ERROR: (execution) Invalid subscript or subscript out of range.

How can I add the nullvector for mk1 in this case?   there is more matrix operations  then later 

M_norm_u	Upper MWB	UMWB_skalad	M_norm_l	Lower MWB	LMWB_skalad
0	0	0	A	B	C

(A,B,C is calculated numbers that the matrix could calculated but not the other where there is 0)

(i dont post it here) beacuse I think there should be a step before all the actual calculations, like before the <actual> proc iml-step?

Thanks!

[Ksharp]

Since you are asking a question about IML, better post it at IML forum:
https://communities.sas.com/t5/SAS-IML-Software-and-Matrix/bd-p/sas_iml
And calling @Rick_SAS

[sbxkoenk]

I have just moved this topic to "SAS/IML and Matrix Computations" - board (where it belongs).

Koen

[cleokatt]

Oh, thank you very much 🙂 

[Rick_SAS]

1. What is a MWB? Reference?

2. In your program, the rc_migr_dist matrix is a 7x6 matrix. The error is because you are trying to access the 7th column of a matrix that has only 6 columns. The error is in the loop that iterates over j:
do j = i + 1 to &K.;

3. Your matrix is not banded, so I don't understand why your started your post talking about banded matrices.

4. Your matrix does not look like it contains percentages, so I don't understand the comment /* Assuming percentages are given, convert to decimal */. However, if you want to divide the matrix by 100, do it OUTSIDE the loops after you read the matrix:

    use pseudo_data;
        read all var {&rcmigrdist.} into rc_migr_dist; /* Load matrix with movements shares p_ij */
    close pseudo_data;
    rc_migr_dist = rc_migr_dist / 100; /* Assuming percentages are given, convert to decimal */

5. You can make the error go away by changing the upper limit of the inner loop to &K-1 or ncol(rc_migr_dist)

/* No need to loop:
        sum_p = 0;
        do j = i + 1 to &K.-1;
            sum_p = sum_p + rc_migr_dist[i,j];
        end;
*/
sum_p = sum( rc_migr_dist[i, i+1:&K-1] );

6. I think you would be wise to read the dimensions of the data inside IML by using nrow(rc_migr_dist) and ncol(rc_migr_dist). This is preferable to using a hard-coded macro such as K=7.

[Rick_SAS]

I did a little research, and the OP might be asking about the MWB statistic on p. 24-25 of this paper: Instructions for reporting the validation results of internal models - IRB Pillar I models for credi...

However, the paper is analyzing Markov transition matrices whereas the matrix in the OP example is not a transition matrix because the rows do not sum to 1. To learn about transition matrices in IML, see

• Markov transition matrices in SAS/IML - The DO Loop
• The probability of reaching a terminal state in a Markov chain - The DO Loop (sas.com) 

[cleokatt]

Hi, 

Yes I know. But how can I solve the issue with nx(n-1) matrix? Add a null vector? That was my thought 

[cleokatt]

@Rick_SAS

[Rick_SAS]

If there are 7 states, then each cell indicates the proportion of accounts that transitioned from State=i into State=j during the time period. If no account transitioned to State=7, then the last column of your 7x7 matrix should contain all zeros.

If "add a nullvector" means "add a column of 0s," then I agree. Otherwise, explain a "nullvector" in this context. In matrix theory, a null vector is a vector v such that A*v=0 for the matrix A. But I don't think you are using that definition(?).

[cleokatt]

A={0,0,0,0,0}

sorry 😞

in my case : mk1={0,0,0,0,0,0,0]

[Rick_SAS]

I suggest you read in the names of the variables dynamically (instead of using a macro) and then add the first column if necessary:

proc iml;
    use pseudo_data;
        read all var _NUM_ into rc_migr_dist[colname=varNames]; /* Load matrix with movements shares p_ij */
    close pseudo_data;
    nr = nrow(rc_migr_dist);
    print varNames;
    if upcase(varNames)[1] ^= "MK1" then do;
       rc_migr_dist = j(nr, 1, 0) || rc_migr_dist;
       varNames = "MK1" || varNames;
    end;
    print rc_migr_dist[r=varNames c=varNames];

Be aware that this empirical estimate of the transition matrix will be singular. You are specifying a system that will never allow a transition to the first state. If this is not what you want, then use a small positive probability (such as 1E-3) for the first column and adjust the remaining probabilities so that the rows sum to 1.

[cleokatt]

the query you posted. Just to be sure, hehe. That should be before my code;

proc iml;
    use pseudo_data;
        read all var {&rcmigrdist.} into rc_migr_dist; /* Load matrix with movements shares p_ij */
    close pseudo_data;

    use rclassorgcount;
        read all var {&SOMETHING_THAT.} into rc_orig_count_vec; /* Load vector with number of humans per movement i, Ni */
    close rclassorgcount;

or after? instead? 🙂 

[Rick_SAS]

The changes that we've discussed would make your program look like this:

proc iml;
    use pseudo_data;
        read all var _NUM_ into rc_migr_dist[colname=varNames]; /* Load matrix with movements shares p_ij */
    close pseudo_data;
    nr = nrow(rc_migr_dist);
    if upcase(varNames)[1] ^= "MK1" then do;
       rc_migr_dist = j(nr, 1, 0) || rc_migr_dist;
       varNames = "MK1" || varNames;
    end;
    if nr^=ncol(rc_migr_dist) then 
      print "ERROR: Matrix is not square";

    rc_migr_dist = rc_migr_dist / 100; /* ??? Is this necessary? */
    print rc_migr_dist[r=varNames c=varNames];

    use rclassorgcount;
        read all var {&SOMETHING_THAT.} into rc_orig_count_vec; /* Load vector with number of humans per movement i, Ni */
    close rclassorgcount;

    /* Calc Matrix_normal_u */
    M_norm_u = 0;
    do i = 1 to &K. - 1;
        maxval = max(abs(i - &K.), abs(i + (&K_Max. - &K.) - 1));
        N_i = rc_orig_count_vec[i];
        sum_p = sum( rc_migr_dist[i, i+1:nr] );
        /*print sum_p;*/
        temp = maxval * N_i * sum_p;
        M_norm_u = M_norm_u + temp;
    end;
    create M_norm_u var {"Matrix_normal_u"}; 
    append from M_norm_u;
    close M_norm_u;

[cleokatt]

thanks! I will have a look, but is 

[colname=varNames]; 

where 

      read all var _NUM_ into rc_migr_dist[colname=varNames]; /* Load matrix with movements shares p_ij */

is it necessary? do I have to have it? 

---

@Rick_SAS wrote:

The changes that we've discussed would make your program look like this:

proc iml;
    use pseudo_data;
        read all var _NUM_ into rc_migr_dist[colname=varNames]; /* Load matrix with movements shares p_ij */
    close pseudo_data;
    nr = nrow(rc_migr_dist);
    if upcase(varNames)[1] ^= "MK1" then do;
       rc_migr_dist = j(nr, 1, 0) || rc_migr_dist;
       varNames = "MK1" || varNames;
    end;
    if nr^=ncol(rc_migr_dist) then 
      print "ERROR: Matrix is not square";

    rc_migr_dist = rc_migr_dist / 100; /* ??? Is this necessary? */
    print rc_migr_dist[r=varNames c=varNames];

    use rclassorgcount;
        read all var {&SOMETHING_THAT.} into rc_orig_count_vec; /* Load vector with number of humans per movement i, Ni */
    close rclassorgcount;

    /* Calc Matrix_normal_u */
    M_norm_u = 0;
    do i = 1 to &K. - 1;
        maxval = max(abs(i - &K.), abs(i + (&K_Max. - &K.) - 1));
        N_i = rc_orig_count_vec[i];
        sum_p = sum( rc_migr_dist[i, i+1:nr] );
        /*print sum_p;*/
        temp = maxval * N_i * sum_p;
        M_norm_u = M_norm_u + temp;
    end;
    create M_norm_u var {"Matrix_normal_u"}; 
    append from M_norm_u;
    close M_norm_u;

---