特殊欠損値を使ってnullを区別する話。

null、欠損値とは何か?例えば、何かデータを収集していて、そのデータの一部(変数が)が不明であったり、未知であった場合、そこにはnullが入ります。
また、たとえば独身男性から収集したデータで、配偶者の名前が入る変数があった場合など、そこはnullになります(定義不能・適応不能のため)。

ひとえにNullといってもその意味は分けれたりするわけです。
これをSASで表現するのはどうするでしょうか?


data Q1;
ID='001';SEX='F';ANSWER=7;output;
ID='002';SEX='F';ANSWER=11;output;
ID='003';SEX='F';ANSWER=.;output;
ID='004';SEX='M';ANSWER=6;output;
ID='005';SEX='F';ANSWER=4;output;
run;











というデータがあるとします。
これは、ある調査のデータで、質問の対象は女性のみで、解答として
1から10までの数字がANSWERに入るとします。

ここで、

ID 002のように、解答の有効範囲を超えた値をnullに、
そしてID 004のように男性なのに解答してしまった値もnullに変えます。

data A1;
 set Q1;
 if ANSWER>10 then ANSWER=.A;
 if SEX='M'   then ANSWER=.B;
run;











このようにピリオドの後に任意の文字をいれることで
数値型のフィールドに、判別可能な欠損値を格納できます。

見た目上、AやBといった文字が入ってますが、その性質は普段の.と全く同じです。

ですが例えば

proc format;
 value ANSF .A='範囲外解答による欠損'
            .B='定義外解答による欠損'
            .='未回答による欠損';
run;

のように区別してフォーマットが設定できるので

proc freq data=A1;
 tables ANSWER/missing;
 format ANSWER ANSF.;
run;











みたいににして欠損値をカテゴライズできるわです。



同IDグループの中で、n個前のオブザベーションを参照したり、n個先のオブザベーションを参照したりする処理がハッシュオブジェクトなら割と簡単にできる話

※ちょっとコード修正しました。1obsごとにハッシュオブジェクト定義するようになっていて、効率悪い書き方してました。

SASで、1つ前のオブザベーションが変数の値が欲しい!とか、1つ先が欲しいとかってよくある話です。n個前のということであれば、retainステートメントや、或いはlag関数などで実現できます。

また、n個先であれば、例えばSAS公式の
「次オブザベーションの値を参照する方法」
http://www.sas.com/offices/asiapacific/japan/service/technical/faq/list/body/ba212.html

や、他にpoint=オプションを使った方法(いずれ紹介するかも)などが利用できます。

ただ、例えば、同じIDグループの中で次のオブザベーションといった制限がつくと、少し面倒になります。
(retainであればbyの部分と合わせて、first.で値をリセットしていけばいいし、先をとるためのmergeなら、予めマージキーとしてID内連番を作成して、それに細工をすればできます)

しかし、扱うデータセットのサイズが何万obsとかの巨大なものでないのであれば、ハッシュオブジェクトを使うとかなり楽にでき、さらにより複雑な条件であっても容易に組み込むことができるので、個人的にはお勧めです。

ただし、できれば、ハッシュオブジェクトを勉強して、理屈をある程度理解した上で利用するようにした方がいいと思います。
何をしているのか全くわからずに、コピペからカスタマイズすると、少し危険に思います。

過去のハッシュオブジェクト関連の記事
http://sas-tumesas.blogspot.jp/search/label/ハッシュオブジェクト


それではやってみましょう。

data Q1;
ID='001';X=1;output;
ID='001';X=3;output;
ID='001';X=5;output;
ID='001';X=7;output;
ID='001';X=9;output;
ID='001';X=11;output;
ID='002';X=2;output;
ID='002';X=4;output;
ID='002';X=6;output;
ID='002';X=8;output;
run;


















というデータがあって、同じIDの中で、1obs前のXの値をX_m1に、1obs先のXを_p1、2つ先をX_p2という変数にいれるとします。

data A1;
 if _n_ = 1 then do;
   declare hash hx();
    hx.defineKey('ID','NO');
    hx.definedata('_X');
    hx.definedone();
  do while (^FL);
   set Q1 end=FL;
   by ID;
   if first.ID then NO=0;
   NO+1;
   _X=X;
   hx.add();
  end;
 end;
   set Q1;
   by ID;
   if first.ID then RNO=0;
   RNO+1;
   
   /*同ID内の一つ前のデータを参照*/
   NO=RNO-1;
   rc=hx.find();
   X_m1=ifn(rc=0,_X,.);

   /*同ID内の一つ先のデータを参照*/
   NO=RNO+1;
   rc=hx.find();
   X_p1=ifn(rc=0,_X,.);

   /*同ID内の二つ先のデータを参照*/
   NO=RNO+2;
   rc=hx.find();
   X_p2=ifn(rc=0,_X,.);

 keep ID X X_m1 X_p1 X_p2;
run;


結果は


















このように、ハッシュオブジェクトは、データステップ中の処理と独立したルックアップテーブルを作成して、データステップ中に値のやりくりができるので、使いこなせば、プログラミングの自由度がかなりアップします。



西暦2万1年問題

さて、普段は知ってても知らなくても、どーでもいいようなことばっかり書いているこのブログですが、今回ばかりは全てのSASプログラマーにとって極めて深刻な問題について取り上げます。
姿勢を正して、真剣に読んでください。

data A1;
X=mdy(1,1,20000);
Y=mdy(1,1,20001);
format X WORDDATE32.;
run;

をやると結果はどうなるでしょう?

mdy関数を知らない方は、SAS忘備録の記事
「日付値を簡単に作るMDY関数の紹介」
をよんでください

WORDDATEフォーマットは
「 January 1 2014」みたいな形式です。

ただ、mdy関数もフォーマットも、この深刻な問題の根本にはあまり関係ないです。


結果は






です。

関数の指定は何も間違っていないのに、何故欠損なのか?

実はSASは日付として、西暦2万年12月31日にあたる範囲までしか、サポートしてないからです。
2万1年目はないんです。

今後、現在のバージョンのSASで作られたシステムを使い続けると、西暦2万1年の1月1日を
迎えた途端、全てのSASのtoday関数などが欠損になります。
それによって世界が大混乱するわけですね。

これがSASの西暦2万1年問題です。

today関数が欠損を返しても、ミサイルとか発射されないように、エラー分岐に気を付けてください。


やっぱり、どーでもよかった



retainステートメントに配列指定できる話

もしかしたら、初等技術なのかもですが、最近初めて知った話です。

例えば












があって、
Xで重複を消して、Yを使って、Zを横展開しろと言われたら
迷うことなく

proc transpose data=Q1 out=A1(drop=_NAME_) prefix=Z_;
 var Z;
 by X;
 id Y;
run;









ですね。

一番好きなプロシージャはと訊かれたらtransposeと答える僕なので
(訊かれたこと一度もありませんが)、これ以外のやり方がぱっと思いつかないんですよね。

ただ、最近知ったのが

proc sql noprint;
 select max(Y) into:MY
 from Q1;
quit;

data A2;
 set Q1;
 by X;
 array Z_{&MY};
 retain Z_;
 if first.X then call missing(of Z_{*});
  Z_{Y}=Z;
  if last.X then output;
 drop Y Z;
run;

で、この場合結果は同じです。
Yの値が例えば1と4だけなどで飛び値になっている場合、transposeだとZ_1とZ_4だだけが作られますが、こちらの方法だとZ_2とZ_3もnullで作成されるので、ケースによって使い分けれそうです。

何が面白いって、retain 配列名が通るんですね!確かに、通りそうですけど、意外と
その発想がなかったです。


あとついでに関係ないけど、これまた最近知った

data A3;
array AX{2,3,4}(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1);
 D1=dim(AX);
 D2=dim2(AX);
 D3=dim3(AX);
keep D1 D2 D3;
run;





へ~dim関数って dim2 dim3とかにして、次元ごとの数とれるんですね!
いつも1次元でdim使ってたけど、dim1ってことだったのか~

マクロ変数に数値をいれて戻すと誤差がでちゃう場合がある問題について考える話

SASプログラミング掲示板の方の話題で
http://tumesas.progoo.com/bbs/tumesas_topic_pr_20.html

話題の本筋ではないのですが
p値をマクロ変数に入れて、後でそこから値を取り出すと、本来の値と比べて凄い小さい誤差がでてるっぽいねという話がありました。

まあ、基本コンピューターの小数点表現には限界があるので、みんな、疑わしきはroundみたいな習慣がついてるので、問題ないね!ということです。
そういや、そもそも、僕の場合、何らかの計算結果をマクロ変数に入れるってこと自体、したことないかも。

けど、一応なんかないかな~と気になって調べてみると

hex16フォーマットを使うという記事がいくつか、ひっかかります
http://www.sascommunity.org/wiki/Unintentional_Loss_of_Precision とか)

読んでもあんまりピンとこないから、実際やってみました。

data _NULL_;
call symputx('A',1/2);
call symputx('B',1/3);
call symputx('C',1/7);
call symputx('D',sqrt(2));
call symputx('E',constant('pi'));
run;

と5つの計算結果をマクロ変数に入れてみます

%put A=&A B=&B C=&C D=&D E=&E;

とすると





お~、掲示板の話と違って、これはわかりやすく切れてますね。
そうか、なんもしないとデフォルトのフォーマットの有効桁数分しか文字に入らんわけか。
まあ一応、確認してみますか

data A1;
A=1/2;
B=1/3;
C=1/7;
D=sqrt(2);
E=constant('pi');

FL1=ifc(A=&A,'同じ','違う');
FL2=ifc(B=&B,'同じ','違う');
FL3=ifc(C=&C,'同じ','違う');
FL4=ifc(D=&D,'同じ','違う');
FL5=ifc(E=&E,'同じ','違う');

keep FL:;
run;

とすると





と当然の結果。

ではでは、hex16で文字化してからマクロ変数に入れて

data _NULL_;
call symputx('HA',put(1/2,hex16.));
call symputx('HB',put(1/3,hex16.));
call symputx('HC',put(1/7,hex16.));
call symputx('HD',put(sqrt(2),hex16.));
call symputx('HE',put(constant('pi'),hex16.));
run;

%put HA=&HA HB=&HB HC=&HC HD=&HD HE=&HE;




マクロ変数を戻す時に、hex16.で再び読み込んで比較すると

data A2;
A=1/2;
B=1/3;
C=1/7;
D=sqrt(2);
E=constant('pi');

FL1=ifc(A=input("&HA",hex16.),'同じ','違う');
FL2=ifc(B=input("&HB",hex16.),'同じ','違う');
FL3=ifc(C=input("&HC",hex16.),'同じ','違う');
FL4=ifc(D=input("&HD",hex16.),'同じ','違う');
FL5=ifc(E=input("&HE",hex16.),'同じ','違う');

keep FL:;
run;





となるんですね。

なんか、やっぱよくわかんないですね。
どんなフォーマット当てようが、どっかで切れてんのは一緒じゃないのかって
思ってしまうのですが…。
誰か詳しく解説してください。

掲示板の内容でもチェックしてみます。

proc univariate data=SASHELP.CLASS noprint;
 var AGE;
 output out=U1 probt=PT;
run;
data _NULL_;
 set U1;
 call symputx('VAL',put(PT,hex16.));
run;
data OUT1;
 set U1;
 VAL=input("&VAL",hex16.);
 if PT=VAL then FLG="同じ";
run;

とすると





やっぱり、できましたね。

へ~

ワードクラウドを作ってみる話

ワードクラウドっていう表現方法があります

ざっくりいうと、なんかしらの頻度を文字の大きさで表現する方法で、だいたいWebページのキーワードの出現頻度なんかを表現する方法ですね。

SAS Visual Analytics だと簡単にできますみたいな情報をよく聞くんですが、そんなもん持ってないんで、SAS Baseの機能でできないかなと思って調べると面白いPaperがありました。

http://analytics.ncsu.edu/sesug/2008/SIB-096.pdf

おおっ、SAS忘備録で取り上げらているRWI (http://sas-boubi.blogspot.jp/2014/10/rwi1.html)を使ってるのか!RWIやってみたいと思っていたのでちょうどいいや!

データは、本来、テキストマイニング的なデータをよく使うんですが、頻度であればなんでも表現できるはずなので

data Q1;
input brand $ count;
cards;
Apple 754
Sharp 53
Fujitsu 23
DoCoMo 76
Samsung 54
Sony 45
HTC 22
Google 20
SonyEricsson 23
LG 13
Microsoft 12
Kyocera 9
Motorola 6
KDDI 5
Lenovo 3
NEC 3
Panasonic 4
Asus 2
Acer 2
Oppo 2
BlackBerry 1
Huawei 1
Pantech 1
;
run;




























と、これは1年間で僕のブログに携帯端末でアクセスしたユーザーが、どこのブランドの端末を使っているかの人数データです。

あ、でも明らかにApple多すぎて、これだけでかい字になっちゃいそうだなぁ。
もう面倒なので、Apple製品以外の携帯端末アクセスユーザーという趣旨で作ります!

これを論文に従って

proc freq data=Q1 noprint;
 where brand^='Apple';
 tables brand / out=OUT_1;
 weight count;
run;























こうして

%let lowPt=10;
%let highPt=40;

フォントの大きさの範囲を最低10、最大40に設定して

proc sql noprint;
 create table OUT_2 as
 select *
  ,(percent-min(percent))
  /(max(percent)-min(percent))
  *(&highPt-&lowPt)
  +&lowPt as size
 from OUT_1;
quit;






















出現頻度を文字の大きさ(size)に反映して


後は

ods escapechar='^';
ods pdf file='/folders/myfolders/wordcloud.pdf';

options orientation=landscape;
data _null_;
 length render $32760;
 do until (last_tag);
  set OUT_2 end=last_tag;
  link = catt ( " ^S={URL='#" , brand , "' font_size=" , size , 'pt linkcolor=RED}' , brand, '^S={}' );
  render = catt (render, link);
 end;
 declare odsout cloud();
  cloud.layout_absolute();
   cloud.region(width:"6in");
    cloud.table_start();
     cloud.row_start();
      cloud.format_cell(text:render);
    cloud.row_end();
   cloud.table_end();
 cloud.layout_end();
 run;

 ods pdf close;

とすると
(すみません、コードについて、詳細まで理解できていないので解説は論文をみてください)
なんかよくわかんないWarningでたりするけど、無視して、

結果は
















おおっ!なんかそれっぽい!

満足。


変数の値によって指定された変数の値をとる方法 VVALUEXによる裏技

なんのこっちゃいってタイトルですが、

今回は、簡単そうで意外とわからない処理の話です。

data Q1;
A='赤';B='い';C='歩';X='A';output;
A='青';B='ろ';C='金';X='C';output;
A='黄';B='は';C='飛';X='B';output;
run;


というデータセットがあるとします。

ここから











というデータセットが欲しいとします。

つまり、Xの中に入っている値は、Yにどの変数の値を持ってくるかを指示しているわけです。

さて、1つ1つ

if X='A' then Y=A; 
else if ・・・・と書いていいってもいいですが、仮に変数が100個とかあったら面倒ですね。

この処理を1行で書くことができます。
それは

data A1;
 set Q1;
  Y=vvaluex(X);
run; 

です。

これだけで、上記の結果がでます。


vvaluex?、vvalueじゃないんですね。この関数って本来どういうのでしたでしょうか?

本来は

data Q2;
X=1;
format X Z2.;
run;





のように、フォーマットが設定された変数に対して、その変数に設定されたフォーマット
を使って、putして文字変数を作ってくれる便利関数です(デフォのlengthはいつもの$200.)

data A2;
set Q2;
XV=vvalue(X);
XVX=vvaluex('X');
format X best.;
run;






わかりにくいですが、XVとXVXは文字型変数です。

でvvalue関数がそのまま変数を引数にとるのに対して、vvaluexは変数名の文字値を引数に
とるんですね!括弧の中のXによくみると、シングルコーテーションがついてます。

そして変数に特にフォーマットが設定されていなければそのままデフォルトのフォーマットをつかいますので、つまり値そのまんまということです。

その性質をうまく使ってるんですね!!

vvaluexというと、フォーマットでputするだけの関数と思わせておいてのこの離れ業!
面白いですね!


means(summary)プロシジャのautonameオプションの話

今、以下のデータセットがあって

data Q1;
X=1;Y=2;Z=3;output;
X=2;Y=4;Z=3;output;
X=3;Y=2;Z=6;output;
run;









で、

proc means data=Q1 noprint;
 var X Y Z;
 output out=OUT(drop=_TYPE_ _FREQ_) mean=;
run;

とすると、データセットOUTの中身は





となります。

要約統計量=の後、何も指定しなければ、元の集計対象の変数名を上書きして
結果が格納されます。

出したい要約統計量が1つならそれでいいのですが

proc means data=Q1 noprint;
 var X Y Z;
 output out=OUT(drop=_TYPE_ _FREQ_) mean= sum=;
run;

のように2つ指定しても、格納する変数がないので、




上記のようなWarningがでて(University Editionのメッセージは全部英語なんです、、。)
sum=の方は無視されます。


なので

proc means data=Q1 noprint;
 var X Y Z;
 output out=OUT_1(drop=_TYPE_ _FREQ_) 
            mean=MEAN_X MEAN_Y MEAN_Z
            sum=SUM_X SUM_Y SUM_Z;
run;

のように集計結果を格納したい新規変数名を、順番に記述してやれば



となります。


が、はっきり言って、1つ1つ指定するのは、変数が多いと大変面倒なので

proc means data=Q1 noprint;
 var X Y Z;
 output out=OUT_2(drop=_TYPE_ _FREQ_)
        mean= sum=/autoname;
run;

こうしてやります。

すると







「変数名_要約統計量名」のルールに従って、自動的に変数が作成されるのです。
まあ便利!


ちなみに、例えば、XとYの平均、Zの合計だけを出したいとかって感じで、選択したい場合は

proc means data=Q1 noprint;
 var X Y Z;
 output out=OUT_3(drop=_TYPE_ _FREQ_)
            mean(X Y)=
            sum(Z)= /autoname;
run;






OKです。

森下卓九段が今期のNHK杯で敗退してしまったので、最近テンション低めです

小技紹介:①sumabsで絶対値合計 ②リンク付きタイトル

全く関連性ないですが、最近知った小技を2つ。

data Q1;
A=1;B=.;C=-3;D=4;
X=sumabs(A,B,C,D);
run;

sumabsという関数があって、そのまんまで、絶対値の合計を出してくれます。

結果は





9.2から登場の関数ですが、そんなに需要あったんかなぁ。何の処理に使うんでしょうか。


で次は

ods html;
title1 '普通のタイトル';
title2 link="http://tumesas.progoo.com/bbs/" "リンクがついてるタイトル";
proc print data=Q1;
run;
ods html close;















のtitle2ように link=を指定すると、タイトルをクリックするとリンク先に飛ぶ仕様にすることができます。