SAS

GTLでヒストグラムを作図する

  /

データの分布を可視化するときはヒストグラムをよく使います。特にデータの分布が正規分布に従うかどうかを直感的に判断したいときは便利です。

GTLでヒストグラムを作図するときはhistgramステートメントかhistgramparmステートメントを使用します。

Contents
  1. histgramステートメント
  2. mirrored histogramの作成

histgramステートメント

histgramステートメントは作図したい変数を指定するだけですぐにヒストグラムを作成してくれます。

それに対しhistgramparmステートメントの場合は、あらかじめ算出したビンの中心点および頻度またはパーセントを格納した変数を指定する必要があります。このビンの中心を格納した変数はあらかじめソートしなければなりません。当然ですが頻度あるいはパーセントを格納した変数は負の値を格納することはできません。

ビンの設定はデータの範囲とオブザベーション数から自動で設定されますので、ビンの設定を定義しなくてもOKです。

縦軸の値はデフォルトはパーセントですが、頻度や密度を指定することも可能です。

試しにsashelpのデータセットを使って新生児の体重の分布を作図してみます。

proc template;
define statgraph hist;
begingraph;
   layout overlay;
      histogram weight;
   endlayout;
endgraph;
end;
run;

proc sgrender data=sashelp.bweight template=hist;
run;

実行結果（ヒストグラム）

ビンの設定

次にビンの設定を追加してみます。histgramparmステートメントではビンの設定はできません。

ビンの設定説明
BINSTART最初のビンの中心点を定義する
BINWIDTHビンの幅を定義する
NBINSビンの個数を定義する
proc template;
define statgraph hist;
begingraph;
   layout overlay;
      histogram weight /
         binstart=200
         binwidth=100;
   endlayout;
endgraph;
end;
run;

実行結果（ビンの設定を追加）

応答変数の設定

応答変数（縦軸の変数）は特に設定しない場合パーセントとなりますが、頻度や密度に指定することができます。

変更する場合はscaleオプションをpercent, countあるいはpropotion(比率)のいずれかを指定します。まあどれを選んでもヒストグラムの形状は変わりませんが。

グループの設定

groupオプションを設定すると指定したクラス変数ごとにヒストグラムを作成します。

先ほどのヒストグラムを男女に分けて出力してみましょう。boyという変数は男子の場合1が格納されるのでこれを利用します。

この時ビン同士が重なってしまいわかりにくくなる場合は、透過度や塗りつぶしの設定を変えて見やすくしたほうがいいです。

proc template;
define statgraph hist;
begingraph;
   layout overlay;
      histogram weight / 
         binstart=200
         binwidth=100
         group=boy
         fillattrs=(transparency=0.5)
         filltype=ALPHAGRADIENT
         name="hist" ;
      discretelegend "hist";
   endlayout;
endgraph;
end;
run;

実行結果（男女に分けて出力）

若干男子のほうが重いですが、性別差はほとんどないとみてよさそうですね。

mirrored histogramの作成

前項ではクラス変数ごとにヒストグラムを作成しましたが、これだと２つの群が重なってしまいあまり見やすくありません。この場合はmirrored histogramを検討したほうが良いでしょう。

mirrored histogramは片方の群を度数（または密度）=0を境にして上下反転させることで群同士を重ならないように工夫したヒストグラムです。臨床試験分野ですと2群の傾向スコアの分布を可視化するときに用いられることがあります。
mirrored histogramはSASのオプションで作成することができません。layout latticeを併用する方法と度数をマイナスに反転させる方法の2通りがありますが、レイアウトの柔軟性を考えると後者のほうが利便性が高いと思います。

今回は2群の傾向スコアを2通りの方法でヒストグラムにしてみます。処理群とコントロール群の傾向スコアの分布をmirrored histogramとして表示します。ビン幅は0.05刻みに設定します。

layout lattice の使う方法

こちらはプロットエリアを2分割し、各群のヒストグラムを作図します。各群ごとに傾向スコアの変数を用意しなくてもeval関数をifn関数の併用することで特定の群のデータのみを表示させることが可能です。
下のヒストグラムはY軸を反転（reverse=True)とすることで上下反転させることができます。

histogramステートメントのbinaxisオプションをfalseにすると軸目盛に依存せずに任意のビンを設定できます。Trueだと目盛りに基づいてビンが設定されます。

この方法のメリットは別途ヒストグラムの集計をしなくてもプロットステートメントだけで直接作図できる点ですが、密度＝０の部分に空白ができてしまいます。これはY軸の目盛りの値を表示するスペースを確保するためのものなので上下のプロットエリアをぴったりくっつけることはできません。

/*データセットps trt=群変数(0=control, 1=treated) ps=傾向スコア*/
proc template;
define statgraph mirror1;
begingraph; /*塗りつぶしの設定*/ 

discreteattrmap name="_map" / discretelegendentrypolicy=attrmap;
   value "Control" / fillattrs=GraphData1;
   value "Treated" / fillattrs=GraphData2;
enddiscreteattrmap;

discreteattrvar var=trt attrvar=_grp attrmap="_map";
   layout lattice / rows=2 columndatarange=union rowgutter=0;
   /*列の軸設定*/ 

   columnaxes;
      columnaxis / offsetmin=0.02 offsetmax=0.02
                   label="Propencity Score"
                   linearopts=(
                      viewmin=0 viewmax=1 
                      tickvaluesequence=(start=0 end=1 increment=0.1)
                   );
   endcolumnaxes;
   /*処理群*/ 
   layout overlay / 
      yaxisopts=(display=(line ticks tickvalues)
      linearopts=(
         viewmin=0 viewmax=0.2 
         tickvaluesequence=(
            start=0 end=0.2 increment=0.05)
         )
       );
      histogram eval(ifn(trt=1, ps,.)) /
         binstart=0.025
         binwidth=0.05
         group=_grp
         scale=proportion
         binaxis=false
         name="hist";
   endlayout;

    /*コントロール群*/
   layout overlay / 
      yaxisopts=(display=(line ticks tickvalues)
      reverse=true
      linearopts=(
         viewmin=0 viewmax=0.2 
         tickvaluesequence=(
            start=0 end=0.2 increment=0.05)
         )
       );
      histogram eval(ifn(trt=0, ps,.)) /
         binstart=0.025
         binwidth=0.05
         group=_grp
         scale=proportion
         binaxis=false;
   endlayout;

   /*共通軸ラベル*/
   sidebar / align=left spacefill=false;
      entry "Density" / 
         rotate=90
         textattrs=(size=10) pad=10;
   endsidebar;

    /*凡例*/
   sidebar / align=bottom spacefill=false;
      discretelegend "hist" /
         displayclipped=true ;
   endsidebar;
endlayout;
endgraph;
end;
run;

proc sgrender data=ps template=mirror1;
run;

実行結果

barchartparamを使う方法

もう一つの方法はあらかじめヒストグラムデータセットを作成し、棒グラフとして作図する方法です。集計データセットの密度をマイナスに反転させることで
上下反転したヒストグラムを作成することができます。データセットを加工してマイナスに反転させるのも良いですが、eval関数とifn関数を併用することで
作図のために必要な変数加工をGTL上で実行することができます。

このままだとY軸の目盛りがマイナスとなってしまいますが、pictureフォーマットを利用することで先頭のマイナス表示を非表示にしています。

この方法ですと体裁の調整の自由度が高いのでお勧めです。なお傾向スコアの分布のFigureはだいたいこの形で発表されているようです。

/*軸目盛からマイナス表記を取り除く*/ 
proc format;
   picture axis(round default=10)
      low-high="9.99";
run;

/*ヒストグラムデータセットを作成*/
ods output histogram=hist;
proc univariate data=ps noprint ;
   histogram ps /
      midpoints=(0.025 to 0.975 by 0.05)
      vscale=proportion;
   by trt;
run;

proc template;
define statgraph mirror2;
begingraph;
   layout overlay / 
      xaxisopts=(
         type=linear
         label="Propencity Score"
         linearopts=(
            viewmin=0 viewmax=1 
            tickvaluesequence=(
               start=0 end=1 increment=0.1)
            )
         ) 
      yaxisopts=(
         type=linear
         label="Density"
         linearopts=(
            tickvalueformat=axis.
            viewmin=-0.2 viewmax=0.2
            tickvaluesequence=(
               start=-0.2 end=0.2 increment=0.05)
            )
         );

      barchartparm category=BinX response=eval(ifn(trt=0,-Bincount,Bincount)) / 
         barwidth=1
         group=trt
         name="hist";

      discretelegend "hist";
   endlayout;
endgraph;
end;
run;

proc sgrender data=hist template=mirror2; 
run;

実行結果

なおグラフパッケージ「SAS Plotter」でも作成できます。

https://picolabs.jp/sasplotter_mirrored_histogram/

RELATED POST