SAS/IML programmers know that the VECDIAG matrix can be used to extract the diagonal elements of a matrix. For example, the following statements extract the diagonal of a 3 x 3 matrix:

proc iml; m = {1 2 3, 4 5 6, 7 8 9}; v = vecdiag(m); /* v = {1,5,9} */ |

But did you know that you can also *assign* the diagonal elements without using a loop? Because SAS/IML matrices are stored in row-major order, the elements on the diagonal of an *n* x *p* matrix have the indices 1, *p*+1, 2*p*+2, ... *np*. In other words, the following statements assign the diagonal elements of a matrix:

start SetDiag(A, v); diagIdx = do(1,nrow(A)*ncol(A), ncol(A)+1); A[diagIdx] = v; /* set diagonal elements */ finish; run SetDiag(m, 3:1); /* set diagonal elements to {3,2,1} */ run SetDiag(m, .); /* set diagonal elements to missing */ |

I've used this trick in several past blog posts (for example, in my post about how to compute the log-determinant of a matrix), but recently I saw some SAS/IML code that used loops to set the diagonal elements, as follows:

do i = 1 to ncol(m); m[i,i] = v[i]; /* set diagonal elements (not efficient) */ end; |

As I have said before, to increase efficiency, replace loops over matrix elements with an equivalent vector operation. Assigning the diagonal elements of a matrix by using this trick is another example of how to avoid looping in the SAS/IML language.

## 3 Comments

Pingback: An efficient way to increment a matrix diagonal - The DO Loop

Pingback: Compute nearest neighbors in SAS - The DO Loop

Pingback: Difference operators as matrices - The DO Loop