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, 2p+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.

Share 