/*************************************************************/ /* ƒoƒCƒvƒƒbƒg : BIPLOT */ /* */ /* DATE: 1997/07/20 */ /* CODE: Tokuhisa SUZUKI */ /*************************************************************/ /*---------------------------------------------------- SAS DATASET ‚πμ¬ -----------------------------------------------------*/ data a ; input x1 - x10 id$ ; cards ; 4.58 4.96 6.02 5.43 5.39 4.66 4.23 4.22 4.14 4.39 •’“Έ 4.56 5.13 6.68 7.15 6.45 3.32 3.64 4.22 4.28 3.59 ”žŽπ 6.13 6.29 7.05 6.76 6.63 6.52 6.81 7.21 6.83 7.35 —Ξ’ƒ 6.08 6.73 6.56 5.80 5.48 5.77 6.00 5.96 5.68 5.75 g’ƒ 6.20 6.75 6.23 6.58 5.49 5.02 6.31 6.28 6.18 5.76 ΰΫΰθ 7.40 7.27 6.75 6.70 6.82 6.27 6.25 5.04 6.36 5.16 ‹“ϋ 8.05 6.84 6.72 6.34 5.90 7.82 7.02 6.29 6.73 6.11 ƒJƒ‹ 7.80 6.59 6.42 6.18 5.78 7.60 6.57 5.94 5.91 5.21 ƒIƒŒ 6.97 5.31 5.25 4.83 4.21 6.73 5.06 4.73 4.16 4.25 •²–– 6.87 6.86 6.23 5.56 4.96 5.60 5.77 4.49 3.87 3.34 ƒRƒ‰ 7.60 6.60 5.80 5.32 5.33 6.90 6.35 5.51 5.49 5.67 ƒTƒC 6.96 5.61 4.34 4.28 4.15 6.46 5.70 5.31 4.77 4.19 ƒlƒN ; run ; /*--------------------------------------------------*/ proc means data = a n mean var std min max maxdec = 3 ; run ; proc corr data = a noprob nosimple cov pearson ; run ; proc standard data = a out = x m = 0 ; run ; proc iml ; use x ; read all var _num_ into x[ colname = vars rowname = id ] ; n = nrow( x ) ; m = ncol( x ) ; call svd( u, q, v, x ) ; /* “ΑˆΩ’l•ͺ‰π */ q2 = q ## 2 ; /* “ΑˆΩ’lƒΙ‚Μ‚Qζ */ p = cusum( q2 / trace( x`* x ) ) ; /* —έΟŠρ—^—¦ */ /*----- type 1 -----*/ f1 = u * diag( q` ) ; g1 = v ; /*----- type 2 -----*/ f2 = sqrt( n ) # u ; g2 = ( 1 / sqrt( n ) ) # v * diag( q` ) ; reset noname; print q[ format = 10.5 colname = "“ΑˆΩ’lƒΙ" ] q2[ format = 10.5 colname = "“ΑˆΩ’l‚Qζ" ] p[ format = 5.2 colname = "—έΟŠρ—^—¦" ] ; print / ( u[ ,1:2] ) [colname = "U-matrix" format=10.5] ( v[ ,1:2] ) [colname = "V-matrix" format=10.5] ; print / " Type of biplot : 1",, ( f1[ ,1:2] ) [colname = "F-matrixiŒΒ‘́j" format=10.5] ( g1[ ,1:2] ) [colname = "G-matrixi•Ο—ʁj" format=10.5] ; print / " Type of biplot : 2",, ( f2[ ,1:2] ) [colname = "F-matrixiŒΒ‘́j" format=10.5] ( g2[ ,1:2] ) [colname = "G-matrixi•Ο—ʁj" format=10.5] ; _name_ = id // vars` ; _type_ = repeat( {"OBS "}, n, 1 ) // repeat( {"VAR "}, m, 1 ) ; id = concat(_type_, _name_ ) ; b1 = f1[ , 1:2 ] // g1[ , 1:2 ] ; b2 = f2[ , 1:2 ] // g2[ , 1:2 ] ; create biplot1 from b1[ rowname = id ]; append from b1[ rowname = id ]; close biplot1 ; create biplot2 from b2[ rowname = id ]; append from b2[ rowname = id ]; close biplot2 ; /*------------------ «ŽΏ‚P -------------------*/ fi = f1[ , 1:2 ] ; gi = g1[ , 1:2 ] ; fg = fi * gi` ; reset autoname ; print / '«ŽΏ‚PF‚†i‚Ζ‚‡j‚Μ“ΰΟ‚Ν‚wij‚Μ‹ίŽ—' ,, '‚†i‚Ζ‚‡j‚Μ“ΰΟ',, fg[ format = 5.2 ] , / '•½‹Ο•Ξ·ƒf[ƒ^s—ρ‚w' ,, ( x[, 1:10])[ format = 5.2 ] ; /*------------------ «ŽΏ‚Q -------------------*/ u2 = u[ , 1:2 ] ; v2 = v[ , 1:2 ] ; q2 = diag( q[ 1:2 , ] ) ; x2 = u2 * q2 * v2` ; /* ŠK”‚Q‚̍s—ρ‚w2 */ xx = x2 * x2` ; ff = fi * fi` ; print / '«ŽΏ‚QFŠK”‚Q‚Μƒf[ƒ^s—ρXX` = ŒΒ‘́isjƒXƒRƒAs—ρFF`' ,, 'ŠK”‚Q‚ΜXX`', xx[ format = 6.2 ] / 'ŠK”‚Q‚ΜFF`', ff[ format = 6.2 ] ; /*------------------ «ŽΏ‚R -------------------*/ g2i = g2[ , 1:2 ] ; gg = g2i * g2i` ; nxx = ( 1/n ) # x2` * x2 ; print / '«ŽΏ‚RF‹€•ͺŽUƒoƒCƒvƒƒbƒg‚ł́A(1/n)X`X = GG`' ,, '(1/n)X`X', nxx[ format = 6.2 ] / 'ŠK”‚Q‚ΜGG`', gg[ format = 6.2 ] ; /*------------------ –β‚S -------------------*/ call svd( u2, q2, v2, x2 ) ; /* ŠK”‚Q‚Μ‚w‚π“ΑˆΩ’l•ͺ‰π */ f = ( u2 * diag( q2` ))[ , 1:2 ] ; g = v2[ , 1:2 ] ; fg2 = f * g` ; print / '–β‚SFŠK”‚Q‚̍s—ρ‚w‚ΜƒoƒCƒvƒƒbƒg•ͺΝ',, 'ŠK”‚Q‚̍s—ρ‚w', x2[ format = 6.2 ] / 'FG`i“ΰΟj‚̍s—ρ' , fg2[ format = 6.2 ] ; /*------------------ –β‚T -------------------*/ j = j( n, 1, 1 ) ; ff = f * f` ; d1 = diag( ff ) * j*j` - 2 # ff + j*j` * diag( ff ) ; d2 = diag( xx ) * j*j` - 2 # xx + j*j` * diag( xx ) ; print / '–β‚TFŒΒ‘̍·ƒoƒCƒvƒƒbƒg‚Μ“_ŠΤ‹——£‚w‚̍sŠΤ‹——£',, 'ŠK”‚Q‚Μ‚w‚©‚η‚ΜƒoƒCƒvƒƒbƒg“_ŠΤ‹——£', d1[ format = 6.2 ] / '•ͺΝƒf[ƒ^‚w‚̍sŠΤ‚Μƒ†[ƒNƒŠƒh‹——£', d2[ format = 6.2 ] ; quit ; data biplot1 ; set biplot1 ; keep col1 col2 _type_ _name_ ; length _type_ $3 _name_ $8 ; _type_ = scan( id, 1 ) ; _name_ = scan( id, 2 ) ; rename col1 = dim1 col2 = dim2 ; run ; data biplot2 ; set biplot2 ; keep col1 col2 _type_ _name_ ; length _type_ $3 _name_ $8 ; _type_ = scan( id, 1 ) ; _name_ = scan( id, 2 ) ; rename col1 = dim1 col2 = dim2 ; run ; data anno1 ; set biplot1 ; length function text $8 ; style = 'hwdmx003' ; xsys = '2' ; ysys = '2' ; text = _name_ ; if _type_ = 'OBS' then do ; /* Label the observation */ color = 'BLACK' ; size = 0.5 ; x = dim1 ; y = dim2 ; position = '5' ; function = 'LABEL '; output ; end; if _type_ = 'VAR' then do ; /* Draw line from */ color = 'BLACK'; size = 0.3 ; x = 0 ; y = 0 ; /* the origin to */ function = 'MOVE' ; output ; x = dim1 ; y = dim2 ; /* the variable point */ function = 'DRAW' ; output ; if dim1 >= 0 then position = '6' ; /* left justify */ else position = '2' ; /* right justify */ function = 'LABEL'; output ; /* variable name */ end; run ; data anno2 ; set biplot2 ; length function text $8 ; style = 'hwdmx003' ; xsys = '2' ; ysys = '2' ; text = _name_ ; if _type_ = 'OBS' then do ; /* Label the observation */ color = 'BLACK' ; size = 0.5 ; x = dim1 ; y = dim2 ; position = '5' ; function = 'LABEL '; output ; end; if _type_ = 'VAR' then do ; /* Draw line from */ color = 'BLACK'; size = 0.5 ; x = 0 ; y = 0 ; /* the origin to */ function = 'MOVE' ; output ; x = dim1 ; y = dim2 ; /* the variable point */ function = 'DRAW' ; output ; if dim1 >= 0 then position = '6' ; /* left justify */ else position = '2' ; /* right justify */ function = 'LABEL'; output ; /* variable name */ end; run ; goptions hsize = 17 cm vsize = 20 cm horigin = 2 cm vorigin = 3 cm device = winprtm ; proc gplot data = biplot1 ; plot dim2 * dim1 / anno = anno1 frame href = 0 vref = 0 lvref = 3 lhref = 3 vaxis = axis2 haxis = axis1 vminor = 1 hminor = 1 ; axis1 length = 13 cm offset = ( 2 ) value = ( h = 0.5 ) label = ( h = 0.5 ) ; axis2 length = 13 cm offset = ( 2 ) value = ( h = 0.5 ) label = ( h = 0.5 ) ; symbol v = none ; footnote1 f = hwdmx003 h = 0.7 '}‚˜@ƒoƒCƒvƒƒbƒgiƒ^ƒCƒv‚Pj' ; run ; quit ; proc gplot data = biplot2 ; plot dim2 * dim1 / anno = anno2 frame href = 0 vref = 0 lvref = 3 lhref = 3 vaxis = axis2 haxis = axis1 vminor = 1 hminor = 1 ; axis1 length = 13 cm offset = ( 2 ) value = ( h = 0.5 ) label = ( h = 0.5 ) ; axis2 length = 13 cm offset = ( 2 ) value = ( h = 0.5 ) label = ( h = 0.5 ) ; symbol v = none ; footnote1 f = hwdmx003 h = 0.7 '}‚˜@ƒoƒCƒvƒƒbƒgiƒ^ƒCƒv‚Qj' ; run ; quit ; footnote1 ;