/*---------------------------------------------------------------*/ /* Name: CANPLOT.SAS */ /* Title: Canonical discriminant structure plot. */ /* Plots class means on two canonical variables, confidence */ /* circles for those means, and variable vectors showing the */ /* correlations of variables with the canonical variates. */ /* Doc: http://www.math.yorku.ca/SCS/sasmac/canplot.html */ /*---------------------------------------------------------------*/ /* Author: Michael Friendly */ /* Created: 24 Nov 1991 13:17:10 */ /* Revised: 23 Mar 1998 09:51:55 */ /* Version: 1.2 */ /* - Added hsym, legend control */ /* - Added canx= cany= to plot other canonical dimensions */ /* - Warning message for <3 groups (only 1 dimension) - no plot */ /*---------------------------------------------------------------*/ /* 修正者 : 鈴木督久 */ /* Revised: 1998.6 */ /* - 群平均と構造ベクトルに日本語表示 */ /* - 構造ベクトルに鏃（矢尻）を表示 */ /* - 座標に平均（０）リファレンス線を表示 */ /* - 個体のラベル表示．( ID = 変数名）を指定した場合のみ */ /* - discrim の誤判別率出力を追加して分離の目安を与える */ /* - discrim で判別結果をデータセットに出力 ( _disc_ ) */ /* - candisc で１変量分散分析を追加して変数選択の目安を与える */ /* - 信頼限界のほか任意の円半径データを指定できる: circleds = */ /* */ /* ＜フィッシャーのあやめのデータをＡ４横に印刷する指定例＞ */ /* goptions vsize = 17.5 cm ; * 頁余白用の指定 ; */ /* %canplot( data=dat.irisj, */ /* class=species, var=_numeric_ , */ /* vorder = %str( order=(-5 to 5 by 1 ) ), */ /* title = f=hwdmx001 h=1 "あやめデータ", */ /* hlength=25 cm, device=winprtm ) */ /* */ /* */ /* ＜１８種類（群）のブランドがあり平均のみをプロットする例＞ */ /* legend1 label = none value = none ; */ /* %canplot( data=x, class=br, var=v1-v20, */ /* device = winprtm, */ /* vlength = 16 cm, */ /* vorder = %str( order=(-2 to 2 by 1 )), */ /* horder = %str( order=(-2 to 2 by 1 )), */ /* symbols = none, */ /* scale = 2.5, */ /* legend = legend1 , */ /* title = f=hwdmx001 h=1 "ブランドイメージ", */ /* hlabel = 0.5, hmean = 0.6 ) */ /* */ /* */ /* ＜各群の売上高シェアを面積とする円を描くためのデータ作成例＞ */ /* data share ; */ /* p = uriage / total ; */ /* size = sqrt( p / 3.14 ) ; */ /* keep &class size ; */ /* run ; */ /* このデータセットを作成しておいて，circleds = share と指定する */ /* */ /*---------------------------------------------------------------*/ %macro canplot ( data = _last_, /* data set to analyze */ class = , /* class variable */ var = , /* classification variables */ scale = 3.5, /* scale factor for variable vectors */ conf = .99, /* confidence probability for can means */ out = _dscore_, /* data set containing discrim scores */ anno = _danno_, /* data set containing annotations */ plot = YES, /* or NO to suppress plot */ id =, /* 個体を識別する変数名・ラベル表示に使う*/ add =, /* １件追加したいデータセット名(x,y変数) */ fontolab = hwdmx001,/* 個体ラベルのフォント．1はMSゴシック */ sizeolab = 0.6 ,/* 個体ラベルのサイズ */ rotate =, r=, /* 回転オプション,proc factor で回転 */ circleds =, /* 分類変数値別の半径値を含むデータセット*/ /* 分類変数名は &class と同一とする */ /* 半径変数名は size に固定する */ /* これが無指定なら信頼限界の円を描く */ /* 例えばシェア（比率p）を面積で表現する */ /* 場合は半径変数 sizeは以下のように作る */ /* size = sqrt( p / 3.14 ) */ device = win, /* グラフ出力デバイス(win, winprtm, ) */ /* win は画面，win以外はプリンタ出力 */ /* デフォルトの axis にオプションを追加 */ horigin = 0 cm, /* グラフ印刷時(device ^= WIN )の左余白 */ vorigin = 0 cm, /* グラフ印刷時(device ^= WIN )の下余白 */ hlength = 60 pct, /* グラフ印刷時(device ^= WIN )の水平軸長*/ vlength = 60 pct, /* グラフ印刷時(device ^= WIN )の垂直軸長*/ /* A4横の場合は hlength = 25 cm 程度 */ horder =, /* %str( order = ( ------ ) ) ,と指定する*/ vorder =, /* %str( order = ( ------ ) ) ,と指定する*/ hmajor =, /* 水平軸の主要目盛の定義 */ hminor =, /* 水平軸の補助目盛の定義 */ vmajor =, /* 垂直軸の主要目盛の定義 */ vminor =, /* 垂直軸の補助目盛の定義 */ /* AXIS, LEGEND をマクロ外で指定する場合 */ haxis =, /* AXIS statement for horizontal axis */ vaxis =, /* and for vertical axis- use to equate axes */ legend =, /* LEGEND statement */ hsym = 0.4, /* height of plot symbols */ hval = 0.6, /* 軸目盛の高さ */ fval = simplex, /* 軸目盛のフォント */ hmean = 0.8, /* プロットする群（平均）ラベルの高さ */ fmean = HWDMX002, /* プロットする群（平均）ラベルのフォント*/ hlabel = 0.6, /* プロットする変数ラベル高さ */ varfont = HWDMX004, /* プロットする変数ラベルフォント */ sizevtop = 0.5, /* ベクトル頭（鏃）のサイズ */ markvtop = C, /* ベクトル頭（鏃）MARKERフォントの値 */ /* 「C」は「▲」．「S」は「＋」 */ /* 「K」は細い「↑」．「G」は太い「↑」 */ can1fmt =, /* dim1 へのformt指定 */ can2fmt =, /* dim2 へのformt指定 */ vref = 0, /* 垂直線 */ href = 0, /* 水平線 */ lvref = 33, /* 垂直線の線種番号 */ lhref = 33, /* 水平線の線種番号 */ lcircle = 1, /* 信頼楕円の線種 */ canx = can1, cany = can2, name = CANPLOT, /* name for graphic catalog entry */ deltmpd = no, /* 一時データセットの削除 */ colors = RED GREEN BLUE BLACK PURPLE YELLOW BROWN ORANGE, can1lab = f = hwdmx001 h = 0.8 "群平均信頼領域（信頼係数=&conf）と構造ベクトル", can2lab = ' ', title = f = hwdmx001 h = 1 "正準判別分析による Biplot", symbols = circle plus triangle star x square y ); %let plot = %upcase(&PLOT); %if %length(&var) = 0 %then %do; %put CANPLOT: You must specify VAR= variables list ; %goto DONE ; %end; %if %length(&class) = 0 or %length(%scan(&class,2)) > 0 %then %do; %put CANPLOT: You must specify one CLASS= variable ; %goto DONE ; %end; %if &rotate ne %str() or &r ne %str() %then %do ; %if &r ne %str() %then %let rotate = &r ; %end ; /*----- データセット中の変数ラベルを取得する -------------*/ %if &data = _LAST_ %then %let data = %scan( &sysdata , 2 ) ; %put NOTE: dataset name = &data . ; data format ; if 0 then set &data ; keep start label fmtname type ; length start $8 label $40 ; fmtname = 'VLABELF' ; type = 'C' ; array vv( * ) &var ; do i = 1 to dim( vv ) ; call vname( vv(i), start ) ; call label( vv(i), label ) ; output ; end ; stop ; run ; proc format cntlin = format ; run ; title '正準判別分析' ; proc candisc data = &data outstat = _dstat_ out = &out ncan = 2 anova ; classes &class ; var &var ; run ; title '線形判別分析による判別力' ; proc discrim data = &data short crossvalidate ncan = 2 out = _disc_ ; classes &class ; var &var ; run ; title '見かけの判別率' ; proc freq data = _disc_ ; tables &class * _into_ / nocol norow nopercent ; run ; /*--- 回転オプションの選択 ---*/ %if &rotate ne %str() %then %do ; title '構造行列の回転' ; proc factor data = _dstat_ outstat = _dstat_ method = score r = &rotate all ; run ; title ; data _coeff_ ; drop _type_ ; set _dstat_ ; if _type_ = 'PATTERN' ; run ; proc score data = &data out = &out score = _dstat_ ; var &var ; run ; %end ; %else %do ; data _coeff_ ; drop _type_ ; set _dstat_ ; /* get standardized coefficients */ if _type_ = 'STRUCTUR' ; run ; %end ; proc transpose data = _coeff_ out = _coeff_ ; run ; title '構造行列' ; proc print data = _coeff_; run ; title '群平均値' ; proc sort data = &out ; by &class ; run ; proc means data = &out noprint ; var &canx &cany ; by &class ; output out = _means_ mean = &canx &cany n = n ; run ; %if &circleds ne %str() %then %do ; data &circleds ; merge &circleds _means_( keep = &class &canx &cany ) ; by &class ; run ; %end ; data &anno ; set _means_ end = eof ; length text $60 color function $8 ; retain xsys '2' ysys '2'; drop &canx &cany n ; x = &canx ; y = &cany ; color = scan("&colors", _n_ ) ; /* mark the class mean */ text = left( &class ) ; hsys = '4' ; size = &hmean ; style = "&fmean" ; function = 'LABEL ' ; output ; %if &circleds = %str() %then /* 信頼限界円を描く */ %do ; /* draw confidence region */ size = sqrt( cinv(&conf, 2, 0) / n ) ; /* radius */ line = &lcircle ; do a = 0 to 360 by 10 ; /* draw a "circle" */ ang = a * arcos(-1) / 180 ; /* convert to radians */ xp = size * cos( ang ) ; yp = size * sin( ang ) ; x = xp + &canx ; /* translate to means */ y = yp + &cany ; if a = 0 then FUNCTION = 'MOVE ' ; else FUNCTION = 'DRAW ' ; output; end; %end ; %else %do ; /*-------- 指定した半径の円を描く ---*/ if eof then do while( not endobs ) ; style = 'empty ' ; xsys = '2'; ysys = '2' ; hsys = '2' ; function = 'pie ' ; rotate = 360 ; color = 'black' ; angle = 0 ; set &circleds end = endobs ; x = &canx ; y = &cany ; output ; end ; %end ; if eof then do ; /* save number of groups */ call symput('NGP', put(_n_, best5.) ) ; end ; run ; %if &ngp < 3 %then %do; %put WARNING: CANPLOT cannot produce a plot for &NGP groups.; %goto done; %end; data _vector_ ; set _coeff_ ; where ( &canx ^= . ) ; length function $8 text $60 ; retain xsys '2' ysys '2' position '6' ; x = 0 ; y = 0 ; function = 'MOVE ' ; output ; x = &scale * &canx ; y = &scale * &cany ; function = 'DRAW ' ; output ; /*---- 変数ラベルがあれば表示する ------*/ text = _NAME_ ; text = put( text, $vlabelf. ) ; size = &hlabel ; style = "&varfont" ; select ; when( x >= 0 & y >= 0 ) position = '3' ; when( x >= 0 & y < 0 ) position = '9' ; when( x < 0 & y >= 0 ) position = '1' ; when( x < 0 & y < 0 ) position = '7' ; end ; function = 'LABEL' ; output ; /*------ 鏃（やじり）を描く ------*/ function = 'LABEL' ; style = 'MARKER ' ; text = "&markvtop" ; size = &sizevtop ; position = '5' ; ab = sqrt( x ** 2 + y ** 2 ) ; /* 斜辺の長さを求める */ sin = y / ab ; /* SINを求める */ radian = arsin( sin ) ; /* 逆サイン関数でラジアンを求める */ angle = abs( radian ) * 57.29578 ; /* ラジアンを角度(絶対値）に変換 */ select ; /* 象限の条件分けで０度〜360度に */ when( sign(x) = 1 & sign(y) = 1 ) angle = angle + 270 ; when( sign(x) = -1 & sign(y) = 1 ) angle = 90 - angle ; when( sign(x) = -1 & sign(y) = -1 ) angle = angle + 90 ; when( sign(x) = 1 & sign(y) = -1 ) angle = 270 - angle ; when( sign(x) = 0 & sign(y) = 1 ) angle = 0 ; when( sign(x) = 0 & sign(y) = -1 ) angle = 180 ; when( sign(x) = 1 & sign(y) = 0 ) angle = 270 ; when( sign(x) = -1 & sign(y) = 0 ) angle = 90 ; end ; output ; run ; %if &id ne %str() %then %do ; data &id ; /* Label the observation */ set &out ( keep = &id &canx &cany ) ; color = 'BLACK' ; xsys = '2' ; ysys = '2' ; x = &canx ; y = &cany ; text = &id ; position = '6' ; /* 左寄せ */ style = "&fontolab" ; size = &sizeolab ; function = 'LABEL ' ; if n( x, y ) = 2 ; if &id ne ' ' then output ; run ; %end ; %if &add ne %str() %then %do ; data &add ; /* １件の追加表示 */ set &add ( keep = x y ) ; color = 'BLACK' ; xsys = '2' ; ysys = '2' ; position = '6' ; /* 左寄せ */ style = "&fontolab" ; size = &sizeolab ; function = 'LABEL ' ; text = '▲' ; run ; %end ; data &anno ; set &anno _vector_ &id &add ; run ; %if %substr( %upcase( &plot ), 1, 1 ) = Y %then %do; goptions horigin = &horigin vorigin = &vorigin device = &device ; %let hlength = %str( length = &hlength ) ; %let vlength = %str( length = &vlength ) ; %if &hminor ne %str() %then %let hminor = %str( minor = &hminor ) ; %if &hmajor ne %str() %then %let hmajor = %str( major = &hmajor ) ; %if &vminor ne %str() %then %let vminor = %str( minor = &vminor ) ; %if &vmajor ne %str() %then %let vmajor = %str( major = &vmajor ) ; %gensym(n=&ngp, h=&hsym, symbols=&symbols, colors=&colors) ; %if %length( &haxis ) = 0 %then %do ; axis2 label = ( &can1lab ) &horder &hminor &hmajor value = ( f = &fval h = &hval ) &hlength ; %let haxis = axis2 ; %end; %if %length( &vaxis ) = 0 %then %do ; axis1 label = ( &can2lab ) &vorder &vminor &vmajor value = ( f = &fval h = &hval ) &vlength ; %let vaxis = axis1 ; %end; %if &legend = %str() %then %do ; legend1 label = ( h = &hmean f = &fmean ) value = ( h = &hmean f = &fmean ) frame ; %let legend = %str( legend = legend1 ) ; %end; %else %if %upcase( &legend ) = NOLEGEND %then %let legend = nolegend ; %if &can1fmt ne %str() %then %let can1fmt = %str( &canx &can1fmt ) ; %if &can2fmt ne %str() %then %let can2fmt = %str( &cany &can2fmt ) ; title &title ; proc gplot data = &out ; plot &cany * &canx = &class / anno = &anno frame href = &href vref = &vref lvref = &lvref lhref = &lhref vaxis =&vaxis haxis = &haxis &legend hm = 1 vm = 1 name = "&name" des = "canplot of &data" ; %if &can1fmt ne %str() | &can2fmt ne %str() %then %do ; format &can1fmt &can2fmt ; %end ; run ; quit; title ; %end; %if %upcase( &deltmpd ) = YES %then %do ; proc datasets lib = work memtype = data nolist ; delete _coeff_ _dstat_ _means_ _vector_ ; run ; quit ; %end ; %done: goptions device = win ; %mend CANPLOT ; /*----------------------------------------------------* | Macro to generate SYMBOL statement for each GROUP | *----------------------------------------------------*/ %macro gensym ( n = 1 , h = 1.5 , i = none, symbols = %str(- + : $ = X _ Y), colors = BLACK RED GREEN BLUE BROWN YELLOW ORANGE PURPLE ) ; /*-- note: only 8 symbols & colors are defined */ /*-- if more than 8 groups symbols and colors are recycled */ %local chr col k ; %do k = 1 %to &n ; %if %length( %scan( &symbols, &k, %str( ))) = 0 %then %let symbols = &symbols &symbols ; %if %length( %scan( &colors, &k, %str( ))) = 0 %then %let colors = &colors &colors ; %let chr = %scan ( &symbols, &k, %str( )) ; %let col = %scan ( &colors, &k, %str( )) ; symbol&k h = &h v = &chr c = &col i = &i ; %end ; %mend gensym ;