SAS

GTLでフォレストプロットを作成する (1)

フォレストプロットはリスク比やハザード比を視覚化するときによく使われます。臨床研究では頻出の作図ですので押さえておくと後々楽です。

フォレストプロットとは

理学療法士学会で非常にわかりやすい解説がありますので、こちらを参照するとよいと思います。

リスク比とその信頼区間をマーカーとエラーバーで表現し(いわゆるドットプロットです)、それらを縦積みするかたちになっています。またプロットの左右に集計結果や信頼区間のテーブルを配置することが多いようです。マーカーとエラーバーはscatterplotで表現できますが、左右のテーブルはaxistableステートメントで作成します。また複雑な表の形式になるためLayout latticeステートメントを使用することになります。layout griddedステートメントでは行と列の幅を調整できないので今回は使用しません。

メタ解析の場合複数調査を統合して得られた結果(pooled result)を表示します。この統合結果はひし形のマーカーで表現することが一般的です。ひし形の両端がリスク比あるいはハザード比の信頼区間に相当します。このため統合結果のマーカーは通常の散布図マーカーとして作図できません。

統合結果マーカー
統合結果マーカー

ひし形ではなく通常のエラーバー付きマーカーでも良いのですが、統合結果は標本数が多くなるため信頼区間が小さくなりやすく、エラーバーがマーカに隠れやすいです。ひし形で表現するのは理にかなっていると感じました。

Rのパッケージと異なりSASにはforest plotを作図するステートメントは存在しないため、この統合結果マーカーはpolygonplotステートメントで自作する必要があります。

作図例

上記のリンクに掲載されているフォレストプロットをGTLで作図してみます。


引用:http://jspt.japanpt.or.jp/ebpt_glossary/forest-plot.html

作図データセット

変数COL1-COL8はtableのデータを格納しています。すべて文字列型に変換しています。個別研究と統合解析の間に空のレコードを挿入して見やすくしています。
変数odds_ratio, upper, lower, weightはドットプロット作図用です。weightの大きさをドットプロットのマーカーサイズとして反映させます。

grpは個別の研究と統合解析の文字列とマーカーを切り替えるために用意しました。
idはデータIDで、Y軸の変数として使います。

作図データセット

id2, x,yは統合結果マーカー作図用の変数です。polygonplotでひし形の図形を作成するために頂点の座標を格納しています。ID2は図形識別IDです。

作図コード

proc template; 
define statgraph forest; 
begingraph; 

discreteattrmap name="dot"; 
   value "1" / markerattrs=GraphData1(symbol=squarefilled); 
   value "2" / markerattrs=GraphData2(symbol=diamondfilled) textattrs=(weight=bold); enddiscreteattrmap; 

discreteattrvar var=grp attrvar=_grp attrmap="dot"; 

/*フットノート*/ 
entryfootnote halign=left "Heterogenity: Chi_square = 22.20, df=7 (p=0.002), I"{sup "2"}"=68%" / ; 
entryfootnote halign=left "Test for overall effect: Z=2.00 (p=0.05)"; 

layout lattice /rows=3 rowweights=(7.5 7.5 85) 
                rowdatarange=union rowgutter=0; 

/* タイトル行1 */ 
layout lattice / columns=7 
                 columnweights=(3 2.8 2.8 1.4 3 1 6) 
                 rowdatarange=union 
                 columngutter=0 
                 pad=0 border=false ; 

   entry halign=center textattrs=(size=11 weight=bold)""; 
   entry halign=center textattrs=(size=11 weight=bold)"intervention"; 
   entry halign=center textattrs=(size=11 weight=bold)"control"; 
   entry halign=center textattrs=(size=11 weight=bold)""; 
   entry halign=center textattrs=(size=11 weight=bold)"Odds ratio"; 
   entry halign=center textattrs=(size=11 weight=bold)""; 
   entry halign=center textattrs=(size=11 weight=bold)"Odds ratio"; 

/* 罫線 */
   drawline x1=15 y1=0 x2=100 y2=0 / 
      drawspace=layoutpercent 
      lineattrs=(color=black thickness=1); 
   drawline x1=0 y1=100 x2=100 y2=100 / 
      drawspace=layoutpercent 
      lineattrs=(color=black thickness=1); 
endlayout; 

/* タイトル行2 */ 
layout lattice / columns=9 
                 columnweights=(3 1.4 1.4 1.4 1.4 1.4 3 1 6) 
                 rowdatarange=union 
                 columngutter=0 
                 pad=0 border=false ; 
   entry halign=center textattrs=(size=10)"Study or Subgroup"; 
   entry halign=center textattrs=(size=10)"Events"; 
   entry halign=center textattrs=(size=10)"Total"; 
   entry halign=center textattrs=(size=10)"Events"; 
   entry halign=center textattrs=(size=10)"Total"; 
   entry halign=center textattrs=(size=10)"Weight"; 
   entry halign=center textattrs=(size=10)"M-H Fixed, 95"{unicode '0025'x}"CI"; 
   entry halign=center textattrs=(size=10)"year"; 
   entry halign=center textattrs=(size=10)"M-H Fixed, 95"{unicode '0025'x}"CI"; 

/* 罫線 */ 
   drawline x1=0 y1=0 x2=100 y2=0 / 
      drawspace=layoutpercent 
      lineattrs=(color=black thickness=1); 
endlayout; 

/*データ*/ 
layout lattice / columns=9 
                 columnweights=(3 1.4 1.4 1.4 1.4 1.4 3 1 6) 
                 rowdatarange=union 
                 columngutter=0 rowgutter=0 
                 pad=0 border=false ; 
/* study name */ 
  layout overlay / walldisplay=none 
     xaxisopts=(display=none) 
     yaxisopts=(display=none reverse=true); 
     
     axistable y=ID value=COL1 / 
        textgroup=_grp 
        display=(values) 
        position=0.4; 
  endlayout;
 
/* invention events */ 
   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true ); 

      axistable y=ID value=COL2 /
         display=(values) 
         position=0.5; 
   endlayout; 

/* invention total */ 
   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true); 

      axistable y=ID value=COL3 / 
         display=(values) 
         position=0.5; 
   endlayout; 

/* control events */ 
   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true ); 

      axistable y=ID value=COL4 / 
         display=(values) 
         position=0.5; 
   endlayout; 

/* control total */ 
   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true); 
      
      axistable y=ID value=COL5 / 
         display=(values) 
         position=0.5; 
   endlayout; 

/* weights */ 
   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true); 

      axistable y=ID value=COL6 / 
         display=(values) 
         position=0.5; 
   endlayout; 

/* odds ratio 95%CI */ 

   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true); 

      axistable y=ID value=COL7 / 
         display=(values) 
         position=0.5; 
   endlayout; 

/* year */ 
   layout overlay / walldisplay=none 
      xaxisopts=(display=none) 
      yaxisopts=(display=none reverse=true ); 

      axistable y=ID value=COL8 / 
         display=(values) 
         position=0.5; 
   endlayout; 

/* dot plot */ 
   layout overlay / walldisplay=none 
      xaxisopts=(type=log logopts=(viewmin=0.01 viewmax=100)) 
      yaxisopts=(display=none reverse=true); 

      scatterplot x=odds_ratio y=ID / 
         group=_grp 
         xerrorlower=lower 
         xerrorupper=upper 
         sizeresponse=weight 
         errorbarcapshape=none; 

      /*参照線*/ 
         referenceline x=1; 

      /*統合結果マーカー*/ 
         polygonplot id=ID2 x=x y=y / 
            display=(fill) 
            fillattrs=GraphData2; 
  endlayout; 
endlayout; 
endlayout; 
endgraph; 
end; 
run; 

*画像出力; 

ods listing gpath="出力パス" ; 
ods graphics /noborder imagename="画像名" width=20cm height=10cm outputfmt=(pngまたはemf) ; 

proc sgrender data=データセット template=forest; 
run;

出力結果

実行結果

これでほぼ論文投稿できるくらいの品質になったかと思います。列の幅はlayout latticeのcolumnweightsオプションで適宜修正しましょう。

GTLでフォレストプロットを作成する (2)コードの解説前回に引き続きフォレストプロットの続きです。フォレストプロットの作図方法の解説は日本語ではほとんどないので個々のコードを説明していきたい...