means summaryプロシジャのclassステートメント指定変数の欠損値には気を付けてって話

年末、インフルにやられてダウンしてました。
締めの記事としてはたいしたことない内容ですが、means summaryのclass変数に欠損値があるのをうっかりすると、危険だという話です。

data Q1;
CL1=1;CL2=1;X=1;output;
CL1=1;CL2=1;X=1;output;
CL1=1;CL2=.;X=1;output;
CL1=2;CL2=2;X=1;output;
CL1=2;CL2=2;X=1;output;
CL1=2;CL2=.;X=1;output;
run;











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

まずCL1のみをclass変数に指定してみます。

proc summary data=Q1;
 types CL1;
 class CL1;
 var X;
 output out=A1 n=N;
run;






CL1が1,2の場合ともにN=3と正しい結果になっています。

続いて、CL2もclass変数に加えたうえで、
typesでCL1 CL2それぞれ単品でクラスした場合の結果のみだします

proc summary data=Q1;
 types CL1 CL2;
 class CL1 CL2;
 var X;
 output out=A2 n=N;
run;









はい、CL2の結果については1,2共にN=2で問題ないのですが
CL1の方について,あれ!!N=2で結果がさっきと違う!!
これなんですよね。複数のclass変数を指定して、あるクラスのみを見た時の結果と
そのclass変数単体で指定した場合の結果が異なるという現象です。

これは、classステートメントに複数変数を指定した場合、いずれか一つ以上の変数にnullが存在する場合、そのオブザベーションが全ての集計対象から外れるという仕様によるためです。

これは知らないと結構、恐ろしい仕様です。

nullが入っていても集計から外さないためにはnull自体を値とみなしてクラス変数として成立させる
という方法があります。
具体的にはclassステートメントにmissingオプションをつけるわけです。

proc summary data=Q1;
 types CL1 CL2;
 class CL1 CL2 / missing;
 var X;
 output out=A3 n=N;
run;

すると










こうなるわけです。1obs目がCL2の欠損値がひとつのクラスとして集計されている部分です。

しかし、これはこれで微妙に危険なところがあります。何故かというと、例えば
typesステートメントをつけずに全集計をデータセットに吐く場合

proc summary data=Q1;
 class CL1 CL2 / missing;
 var X;
 output out=A4 n=N;
run;



















こうなるわけですが、このデータセットの1obs目と2obs目について、CL1 CL2の値だけ
見ると両方null nullで見分けがつきませんよね。
_TYPE_の値が違うことから別層での集計であることがわかります。
1obs目はクラスなしの全集計、2obs目はCL2のnullクラスでの集計です。
例えばこのようなデータセットか、CL1 CL2の値で抽出して臨むべき集計結果を取得しようとする
場合、別層のものも同じ条件でかかってしまう恐れがあるわけです。

結局のところ、クラス変数にnullが入る場合で
missing使う場合は、結果をきちんと想像できるようにする。
それが不安な場合は
①class変数を単体で指定して、何回も回す
②欠損値を99などのダミー値にして回す
の方がまだいいかもですね。

また、たまにN数を_FREQ_からとっているプログラムをみますが、
もしN=と_FREQ_が常に同値と思ってやっているなら、それはヤバいです。

data Q2;
CL1=1;CL2=1;X=.;output;
CL1=1;CL2=1;X=.;output;
CL1=1;CL2=.;X=.;output;
CL1=2;CL2=2;X=.;output;
CL1=2;CL2=2;X=.;output;
CL1=2;CL2=.;X=.;output;
run;

proc summary data=Q2;
 types CL1 CL2;
 class CL1 CL2 / missing;
 var X;
 output out=A5 n=N;
run;










_FREQ_は純粋に読み込まれたオブザベーション数で、N=はそのうち、
varで指定された変数の値が非欠損値であるものの数です。

以上。今年も有難うございました!!

今年は無料のSAS Universityがリリースされ、アクセス履歴をみるとacドメイン(大学)が少し増えてきていて、Rが猛威を振るう中、少しでも新規の若いSASユーザーが増えるというのは、僕にとっては幸せを感じるできごとです(Rもできた方が絶対いいと思うけど)。

そういう新ユーザーが、SAS面白いって思えるような何かを提供するのが先達の務めだとも少しは思うので(一部のデータステップマニア向けで、あんまりそういう人向けのブログではないのですが…)、来年もまあ、頑張りたいと思います。

SASユーザー躍進の一年になりますように!!

means(summary)プロシジャのclassステートメントに関連してtypes ways nwayの話

たまに質問されるのが、meansやsummaryプロシジャで

ERROR: 適切な処理量を超えるクラス変数の組み合せの種類が
       要求されました。 TYPES または WAYS ステートメント、
       あるいは NWAY オプションを使用して、
       組み合せの種類の生成を制限してください。

ってメッセージがでるんですが、どうしたらいいですか?って感じの質問です。
どうもclassステートメントととbyステートメントの違いがあんまりピンときていない感じで、
とりあえず層別化する変数を全部つっこんだらこうなったって感じが多いです。

実際、classステートメント、ちょっとわかりにくですよね。僕も難儀しました。

記事:MEANSプロシジャのnwayオプションとTYPESステートメントの話
http://sas-tumesas.blogspot.jp/2013/10/meansnawytypes.html

でnwayとtypesについては触れていたのですが、waysについては忘れていたので整理ついでにもう一度やりましょう。


data Q1;
call streaminit(777);
do i = 1 to 20;
 C1=rand('table',1/2);
 C2=rand('table',1/2);
 C3=rand('table',1/2);
 C4=rand('table',1/2);
 X =int(rand('uniform')*10);
 output;
end;
drop i;
run;























データはなんでもいいです。
今回はC1-C4という4つの層別化変数と、集計対象のXというデータセットです。

留意してほしいのが、層別化変数に全て値が入っている、欠損がないという点です。
何故かというとclass変数に欠損が入っていた場合について、これまた何回も質問される間違いがあるんですが、それについては次回説明します。

byでやるなら、

proc sort data=Q1;
 by C1 C2 C3 C4;
run;

として

proc summary data = Q1;
by C1 C2 C3 C4;
var X;
output out=A0 mean= / autoname;
run;

とすれば















ってな感じで、C1からC4について、元のデータに存在している全パターンで集計結果がでています。

一方 classで同じようにC1からC4を指定してまわすと

proc summary data = Q1;
class C1 C2 C3 C4;
var X;
output out=A1 mean= / autoname;
run;

結果は以下










































































はい、全部で74オブザベーションでてます。
ここで、はぁ?なんで?と思われた方はmeans summaryプロシジャでのclassステートメントの
理解ができていないということになります。

どうしてかというと

全体集計(絞りなし)-0クラスレベル
C1の値で絞って集計した結果、C2で絞って集計した結果、C3の値で絞って集計した結果、C4の値で絞って集計した結果 -1クラスレベル
そして、C1とC2の出現している値の全組み合わせで、絞った結果、C1とC3、C1とC4・・・といった2クラスレベルといったように出現しているクラスの全組み合わせで、全レベル分でるわけです。
なので、値のパターンがそれなりにある変数をいくつもclassにしているすると、集計する組み合わせパターンはあっという間に爆発的な数になって、冒頭のようにSASがギブアップメッセージをだすわけです。

そこで、まずはnway

proc summary data = Q1 nway;
class C1 C2 C3 C4;
var X;
output out=A2 mean= / autoname;
run;













class変数の全てを使ったフルクラスレベルの結果のみを使用します。

次にtypesを使うと欲しい組み合わせのみを取得できます。
例えばC1とC2の2クラスレベルでの層別結果だけが欲しいのであれば
(その場合、C3 C4をclassステートメントに指定する意味がなくなりますが…)

proc summary data = Q1;
class C1 C2 C3 C4;
types C1*C2;
var X;
output out=A3 mean= / autoname;
run;







って感じです。


或いは、C1とC2 C3 C4それぞれを組み合わせた2クラスレベルが欲しければ
proc summary data = Q1;
class C1 C2 C3 C4;
types C1*(C2 C3 C4);
var X;
output out=A4 mean= / autoname;
run;













といった感じです。


続いて、waysの場合は、欲しいクラスレベルのレベルを数字で指定します。
0クラスレベル(つまり全体集計)と、1つのクラス変数のみをそれぞれ使った1クラスレベルの結果が欲しければ

proc summary data = Q1;
class C1 C2 C3 C4;
ways 0 1;
var X;
output out=A5 mean= / autoname;
run;











となります。

あまり、いい説明ではありませんでしたが、これが冒頭のTYPES WAYS NMAYいずれかを使って結果を絞ってくれというSASのお願いに繋がるのでした。

以上

プロシジャのアウトプットをカタログとして保存し、参照する方法

実は結構、下書きのまま、公開せずに寝かせてしまっている記事が結構あります。
後でもっと内容を良くしたり、自分で理解を深めたり、書き足してから公開しようと思って保存するのですが結局自分がそのネタに飽きたり、何を書き足すつもりだったのかを忘れてしまいます。

ただ最近、SASYAMAさんのブログは質より量なのがいいところですね!とまっすぐな目で言ってくれた人がいました。

そうですね!確かに!


今回の話は、プロシジャの出力をカタログとして保存することで、参照しやすくする話です。

探索的に、プロシジャ何回も回しながらプログラムいじってると、
あれ、さっき見てた出力どれだっけ?とか、出力が多かったから全部消したけど、あの結果みたかったのにもう一回実行するの時間かかるから嫌だな~とか、外部ファイルにいちいちだして確認するの面倒とか、そういうことがあると思います。

そういう時に、後でちらちら参照する出力をカタログとして保存しとくと、SAS上でぱっと確認できるので楽という話です。


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







というデータがあって

まず、プロシジャ出力を保存するカタログを以下のように指定します。

proc printto print=work.oricatalog.Univariate_1;
run;

今回はworkにoricatalogというカタログを作って
univariate_1という名前で保存します。
proc printtoはログやアウトプットを外部ファイル保存する時によく使うプロシジャですね。

proc univariate data=Q1;
 var X Y;
run;

proc printto;
run;

ここで閉じます(printtoの空指定実行は出力先を元に戻す効果を持つ)。

ついでにもう1つcorrプロシジャのアウトプットも保存してみましょう。

proc printto print=work.oricatalog.Corr_1;
run;
proc corr data=Q1;
 var X Y;
run;
proc printto;
run;

そしてWorkを見てみると

カタログができて、開けると

これらは出力結果のlstファイルで、ダブルクリックするとSAS内のウインドウで
テキストエディタ(例はNOTEPAD)で開かれる。




パラメータを変えながら、解析結果を比較したい場合などに、前の結果を
このように別ウインドウで開いて見れるのは地味に便利に感じてます。






call sortc(sortn sorth)ルーチンで降順ソートする小技

call sortc(sortn sorth)ルーチンについては何度か取り上げています。

・横方向(列方向)にcall sortn,call sortcルーチンを使ってソート処理をかける
http://sas-tumesas.blogspot.jp/2013/10/call-sortncall-sortc.html

・CHOOSEN(C)関数とCALL SORTN(C)ルーチンで、横に可変的に増える変数をソートして指定した順番にくる変数を取得する
http://sas-tumesas.blogspot.jp/2014/02/choosenccall-sort.html

call sort系ルーチンには降順指定がないため、基本的には昇順しかできません。
それをどうしたら降順にできるかという話です。
※もしかしたら、既にどこかの記事で話したか、或いはamatsuさんからコメント等で教えて貰った内容かもしれません。最近、自分でも何を投稿したか忘れてて、もし重複してたらすみません。

今、

data Q1;
X1='S';
X2='Y';
X3='O';
X4='G';
X5='I';
run;



があって

data A1;
set Q1;
call sortc(of X1-X5);
run;

とすると






このように昇順になります。

これを降順にするには

data A2;
set Q1;
call sortc(of X5-X1);
run;

というように指定を逆向きにすれば






と、結果的に降順になります。

amatsuさんが
記事:変数指定は「V100-V1」のように逆にもできる。
http://sas-boubi.blogspot.jp/2014/06/v100-v1.html

で紹介されていたことと理屈は同じです。

引数に配列を使う場合は、配列に要素を指定する場合の、指定順番を逆にすれば
OKです。

meansプロシジャのweightステートメントを使用した時の分散とか標準偏差の話

最近まで知らなくて、凄いびっくりしたのですが

data Q1;
X=1;output;
X=1;output;
X=1;output;
X=2;output;
X=2;output;
X=2;output;
X=2;output;
X=3;output;
X=3;output;
X=3;output;
run;



















proc means data=Q1 sum  mean var std;
 var X;
run;










となります。

これと同じことを

data Q2;
X=1;WGT=3;output;
X=2;WGT=4;output;
X=3;WGT=3;output; 
run;









というようにウェイト値にして

proc means data=Q2 sum  mean var std;
 var X;
 weight WGT;
run;

weightステートメントに指定すれば、当然同じ結果になると思い込んでいたのですが









うわっ、分散と標準偏差がえらいでかいな!ってびっくりです。

で、教えて貰ったのが

proc means data=Q2 vardef=wdf sum  mean var std;
 var X;
 weight WGT;
run;

です。

これで実行すると









おお、よかった。

どういうことかというと

SASのデフォルトだと、ウェイトが反映された偏差平方和をウェイトが反映されてないN-1で割って分散だしちゃんですよね。





え~、なにそれ。weightステートメント使ってたらデフォルトをwgtに切り替えるぐらいしてよ。メッセージも無しに、そんな使い道のない謎値だされてもって思いました。

医薬品開発の業界にいたときは、freqで何か検定とかするときにweightステートメント使いましたが、あんまり要約統計量出す時に使った記憶がないです。
最近ウェイトバック集計とかする機会が増えて、そこで教わった次第です。


thenやwhenの後に何も処理しなくてもOKな話

かなりちょっとした話を。

data Q1;
do X=1 to 5;
 output;
end;
run;

データは何でもいいんですが、例えば

data A1;
 set Q1;
 if X =1 then ;
 else if X>=2 then Y=1;
run;



data A2;
 set Q1;
 select;
 when(X=1) ;
 when(X>=2) Y=1;
 end;
run;

のように、thenやwhenの後に何も書かなくてもエラーにならないです。
特にwhenの場合、明示的にスルーしてます感がだせます。

あと分岐後の処理をマクロで展開したりしなかったりする場合、空白でもエラーにならないことを知っておくと便利かも。



糸谷新竜王爆誕の話 ついでに#byvalと put (=)とか

糸谷新竜王が誕生しました!本当に嬉しいです!大学で哲学の授業も受けたくて文学部入って心理学専攻した僕からしたら、院で哲学専攻しながら(現在休学中)竜王を獲得した糸谷竜王はもう憧れです。
嬉しいな~。色々と癖強いところも、もう好きすぎです。

あ、一応、SASの話もしときます。
SAS University Edition(もう面倒くさいから、これからはSASウニもしくは雲丹と表記します)

雲丹のアウトプットは通常のSASのと若干見栄えが違いますが、まあだいたい似たようなものなので製品版をお使いの方は、実際動かして出力してみてください。


data Q1;
input JCLASS $ DAN NAME:$20.;
cards;
A 9 久保利明
A 9 三浦弘行
A 8 広瀬章人
A 8 阿久津主税
B1 8 木村一基
B1 9 藤井猛
B1 8 山崎隆之
run;











といったデータがあって(今期の順位戦在籍クラスと、現時点の段位です)

proc sort data = Q1;
 by JCLASS DAN;
run;

ods html file="/folders/myfolders/output_1.html";
proc print data=Q1 noobs;
 by JCLASS DAN;
run;
ods html close;

とすると



























となります。

このby値の表示の部分を好きにカスタマイズできます。
具体的にはtitleステートメント内に #byval(by変数) とすることでそこの部分に
by値が入ります。

例えば

ods html file="/folders/myfolders/output_2.html";
options nobyline;
proc print data=Q1 noobs;
 title "順位戦:#byval(JCLASS) #byval(DAN)段";
 by JCLASS DAN;
run;
ods html close;



って感じです。

もう一つ小技、putステートメントで

data _null_;
set Q1;
file print;
put JCLASS= DAN= NAME=;
run;

とすると













といった出力がでます。

これを

data _null_;
set Q1;
file print;
put (JCLASS DAN NAME) (=);
run;


ともかけます。
最近知りました。

もちろん

data _null_;
set Q1;
file print;
put (_all_) (=);
run;

でもOK