title 'Basic matrix operations in Proc IML'; options nodate ls=79; /*----------------------------------------------------* | Basic operations in Proc IML | *----------------------------------------------------*/ proc iml; reset print log; *-- Scalars; x = 12.3; y = {57}; name = 'Bob'; msg = "Hi there,"; print msg name x y; *-- Creating matrices and vectors; a = { 2 4, 3 1}; b = { 4 5, 0 1}; c = { 1 2 3 }; d = { 1, 2, 3}; e={'a' 'b' 'c', 'd' 'e' 'f'}; *-- Simple matrix operations; sum = a + b; diff = a - b; *-- Elementwise product; times = a # b; *-- Matrix product; prod = a * b; sqr = a##2; asq = a * a; *-- Comparison operations; max = a<>b; min = a> b; page; *--------------------------------; * Matrix-generating functions ; *--------------------------------; ident = I(3); ones = j(1,5); zero = j(3,2,0); a = j(2,2,2); b = j(2,2,1); ab = block(a,b); a = shape({5 12},3,3); b = shape({5 12},3,3,-1); c = repeat({5 12},3,3); index = 1:5; col = t(4:6); rindex= 5:1; series= do(12,72,12); vars = 'XX1': 'XX7'; page; *-- Diagonal matrices; d = diag( {1 2 4} ); s = diag( {1 2, 3 4} ); v = vecdiag( a ); page; *-----------------; * Subscripts ; *-----------------; *-- Selecting rows or columns; col1 = d[,1]; row2 = d[2,]; row21 = d[{2 1},]; *-- Selecting a submatrix; print ab; row12 = ab[{1 2},]; a = ab[1:2, 1:2]; *-- Assigning a row, column or element; a[1,] = {6 7}; a[,2] = { 0, 10}; a[2,2]= 0; page; *----------------------------------; * Subscript reduction operators ; *----------------------------------; a = { 0 1 2, 5 4 3, 7 6 8 }; colsum = a[+,]; rowsum = a[,+]; colmax = a[<>,]; rowmin = a[,><]; colmean= a[:,]; colprod= a[#,]; colssq = a[##,]; page; *-----------------; * Missing data ; *-----------------; x = { 1 2 . , . 5 6 , 7 . 9 }; y = { 4 . 2 , 2 1 3 , 6 . 5 }; sum = x + y; times = x # y; max = x <> y; *-- NB: missing treated as zero; colsum= x[+,]; page; *--------------------------------------; * Working with vectors and matrices ; *--------------------------------------; coffee = { 4 2 2 3 2, 3 3 1 2 2, 2 1 0 4 5 }; days = { Mon Tue Wed Thu Fri }; names = { 'Lenny', 'Linda', 'Sue'}; print coffee; print coffee[r=names c=days]; *-- Calculate daily and weekly cost at $.50/cup; daycost = .50 # coffee; ones = j(5,1); weektot = daycost * ones; weektot = daycost[,+]; daytot = daycost[+,]; total = daycost[+,+]; print coffee[r=names c=days] weektot[format=dollar7.2] , daytot[c=days f=dollar8.2] ' ' total[f=dollar7.2]; rowmean = coffee[,:]; colmean = coffee[:,]; most = weektot[<>]; index= weektot[<:>]; who = names[index]; *-- Sum of squares and cross-products matrix; n = nrow(coffee); sscp = t(coffee) * coffee; *-- Deviations from column means; dev = coffee - J(n,1) * colmean; *-- Corrected SSCP matrix; sscp= t(dev) * dev; *-- Variance covariance matrix; cov = sscp / (n-1); *-- Standard deviations; std = sqrt( diag(cov) ); *-- Correlation matrix; corr= inv(std) * cov * inv(std); print corr[r=days c=days f=7.2]; quit;