SAS

GTLでbeeswarm plotを作図する

beeswarm plotという作図手法を最近耳にしたのでSASでの作図方法を調査したのでメモしておきました。
beeswarm plotだけでなく、箱ひげ図を重ね書きする方法もご紹介します。

beeswarm plotとは

beeswarm plot(蜂群図)とは個別のデータを散布図のマーカーとして表現し、各マーカーが重ならないように表示位置を微調整することで、全体のデータの分布を俯瞰できるようにしたグラフです。データが比較的少ない場合に有効な手法です。箱ひげ図やバイオリンプロットと異なり、各群のデータの個数の違いも可視化できる点もメリットといえるでしょう。

データ数が多い場合はマーカーが見づらくなるため、箱ひげ図やバイオリンプロットの利用を推奨します。

臨床データはデータ取得コストが高いためデータ数は少ないことが多いため、beeswarm plotは相性はよさそうです。

個々のデータをそのまま表示するため、要約統計量は表示されません。そのため箱ひげ図やバイオリンプロットと併用することもあるようです。

Rではbeeswarmパッケージを、pythonではseabornのswarmplotメゾッドを使用すれば作図できます。

beeswarm plotの作成方法

SASにはbeeswarm plotを作図するステートメントは用意されていません。scatterplotステートメントのjitterオプションで一応マーカーが重ならないように自動調整することもできますが、データによってはマーカー同士が重なるため、jitterオプションでの作図は汎用性に欠けます。

beeswarm plotのデータ表示位置調整アルゴリズムをpythonからSASへ移植するには大変です。そのためすでに公開されているSASマクロを利用してマーカーの表示位置を計算した後、GTLを利用して作図する方法を採用します。

SASマクロ(beeswarm.sas) は以下のサイトから入手可能です。

マクロを使わない場合

マクロを使わない場合、マーカー同士が重なってしまい個々のデータが見にくくなってしまいます。

scatterplot x=x y=total_bill / group=day markerattrs=(symbol=circlefilled transparency=0.5);
マクロを使わずに散布図として作図した場合

scatterplotステートメントにはjitterオプションというものがあります。これはマーカー同士が重なってしまう場合はマーカーの位置を微調整するオプションですが、データ次第ではマーカーが重なってしまい可読性が改善するとは限りません。現に今回のテストデータではマーカーは重なってしまいます。jitterの設定はjitteroptsオプションで設定できますが、少なくともR やpythonのbeeswarmと同じ体裁にすることは難しいでしょう。

scatterplot x=x y=total_bill / 
   jitter=auto 
   group=day 
   markerattrs=(symbol=circlefilled transparency=0.5);
jitterオプションを使用した場合

作図データセットの作成

上記のSASマクロの引数に使用データセット、カテゴリ変数、応答変数を指定するだけでbeeswarm plot用の変数を追加してくれます。画像サイズなどに応じて多少設定をいじる必要がありますが、
基本的に上記の3つを指定するだけで作成できるでしょう。

このときカテゴリ変数は必ず数値変数に変換する必要があります。文字変数の場合はあらかじめ数値変数に読み替えましょう。また以降の作図のためにフォーマットを作成してきます。

今回は木曜日から日曜日までの請求額の総額をbeeswarm plotで可視化してみます。%beeswarmを実行するとデータセットbeeswarmが作成されます。変数名に[_bee]がついている変数がbeeswarm plot用の変数です。

%include "beeswarm.sas"; 
*カテゴリ軸用のフォーマット; 

proc format; 
value dat 1="Thur" 
          2="Fri" 
          3="Sat" 
          4="Sun"; 
run; 

*カテゴリ変数を数値変数に読み替え; 

data import; 
set import; 
   select (day); 
     when("Thur") x=1; 
     when("Fri") x=2; 
     when("Sat") x=3; 
     when("Sun") x=4; 
   end; 
run; 

*カテゴリ変数=x 応答変数=total_bill; 

%beeswarm(data=import ,grpvar=x ,respvar=total_bill );

作図テンプレート

データセットをつくってしまえばあとは散布図として作図するだけです。マクロで作成した変数x_beeをカテゴリ変数として指定するだけです。
カテゴリ変数は数値変数のため、軸目盛に余計な値が表示されないようにtickvaluelistオプションとtickvalueformatオプションで表示する目盛を制御しています。

今回は曜日毎にマーカー色を変えたいので、groupオプションを指定しています。

*軸目盛リスト; 
%let lst=(1 2 3 4); 

proc template; 
define statgraph swarm; 
begingraph; 
   layout overlay / xaxisopts=(offsetmin=0.05 offsetmax=0.05 
      linearopts=(tickvaluelist=&lst. tickvalueformat=dat.)
   ); 

      scatterplot x=x_bee y=total_bill / 
         group=day 
         markerattrs=(symbol=circlefilled); 
    endlayout; 
  endgraph; 
end; 
run; 

proc sgrender data=beeswarm template=swarm; 
label x_bee="曜日" total_bill="請求額"; 
run;
実行結果

水平方向にプロットする場合はscatterplotステートメントの引数を入れ替えてY軸の設定を変更します。

*軸目盛リスト; 
%let lst=(1 2 3 4); 

proc template; 
define statgraph swarm; 
begingraph; 
layout overlay / 
   yaxisopts=(offsetmin=0.05 offsetmax=0.05 
      linearopts=(tickvaluelist=&lst. tickvalueformat=dat.) 
      reverse=true); 

   scatterplot y=x_bee x=total_bill / 
      group=day 
      markerattrs=(symbol=circlefilled); 
endlayout; 
endgraph; 
end; 
run; 

proc sgrender data=beeswarm template=swarm; 
label x_bee="曜日" total_bill="請求額"; 
run;
実行結果

箱ひげ図と重ね書きする

beeswarm plotは要約統計量を可視化できないため、箱ひげ図と併用するケースがあります。上記のSASマクロを使えば作図は簡単です。
boxplotステートメントを追加するだけです。またboxplotのカテゴリ変数(今回のケースはX軸変数)はマクロで作成した変数ではなく元の変数を指定するのがポイントです。

上記のテンプレートに箱ひげ図を追加すると以下のようになります。デフォルトだと見ずらいのので箱ひげ図の塗りを非表示にして、透明度を調整しています。

%let lst=(1 2 3 4); 
/* 箱ひげ図と重ね書き */ 

proc template; 
define statgraph swarm; 
begingraph; 
layout overlay / 
   xaxisopts=(type=linear linearopts=(tickvaluelist=&lst. tickvalueformat=dat.)); 

   scatterplot x=x_bee y=total_bill / 
      group=day 
      markerattrs=(symbol=circlefilled transparency=0.5); 

   boxplot x=x y=total_bill / 
      display=(median mean) 
      boxwidth=0.7 
      meanattrs=(symbol=diamondfilled color=black) 
      outlineattrs=(thickness=2 color=black) 
      medianattrs=(thickness=3 color=black) 
      whiskerattrs=(thickness=2 color=black); 
endlayout; 
endgraph; 
end; 
run; 

proc sgrender data=beeswarm template=swarm; 
label x_bee="曜日" total_bill="請求額"; 
run;
実行結果