proc iml; /* matrix "multiplication" where missing values are propagated */ start MVMult(A, B); C = j(nrow(A), ncol(B), .); rows = loc(countmiss(A, "ROW")=0); cols = loc(countmiss(B, "COL")=0); if ncol(rows)>0 & ncol(cols)>0 then C[rows, cols] = A[rows,] * B[,cols]; return(C); finish; A = {1 2 3, . 4 1, -1 0 1}; B = {1 2 -1, 3 4 0, 0 1 .}; C = MVMult(A,B); print C; A1 = MVMult(A, I(3)); /* A*I */ A2 = MVMult(I(3), A); /* I*A */ H = MVMult(A2, I(3)); /* I*A*I */ print A1, A2, H; /* Skip missing values during multiplication. Equivalent to substituting 0 for missing */ start SkipMVMult(A, B); C = j(nrow(A), ncol(B)); do i = 1 to nrow(A); w = A[i,]`; do j = 1 to ncol(B); v = w#B[,j]; if all(v=.) then C[i,j] = .; else C[i,j] = sum(v); end; end; return(C); finish; S = SkipMVMult(A,B); print S; A0 = A; B0 = B; A0[ loc(A=.) ] = 0; /* replace missing values with 0 */ B0[ loc(B=.) ] = 0; C0 = A0*B0; print C0; A1 = SkipMVMult(A, I(3)); /* A*I */ A2 = SkipMVMult(I(3), A); /* I*A */ H = SkipMVMult(A2, I(3)); /* I*A*I */ print A1, A2, H; A = {-1 . ., 1 3 4}; M1 = SkipMVMult(A,B); M2 = MVMult(A,B); print M1 M2; B = { . 1, -1 3, -3 1}; M1 = SkipMVMult(A,B); M2 = MVMult(A,B); print M1 M2;