Version8における医薬統計関係の拡張（株）ＳＡＳインスティチュートジャパン
テクニカルサポート小野裕亮
作成：2000年6月23日
修正1：2000年6月28日

ナッシュビル・プロジェクト✸ Ｖ８での様々な機能拡張。例えば ...
✸ Base SAS：ODS✸ GRAPH：

ActiveX, Java
✸ AF： SAS Component Language✸ ACCESS： LIBNAME✸ WebAF,WebEIS ： Java✸ Integration Tech：COM/DCOM,CORBA✸ CONNECT： 非同期 , CPORT不要✸ Enterprise Guide： thin -client のＧＵＩ

1

Version8,SAS/STATプロダクトの構成✸✸✸✸✸

既存42プロシジャ＋ 新規13プロシジャ＋ 評価 3プロシジャ----------------------------合計 58プロシジャ

SAS/STATに新しく追加されたプロシジャ(1)✸ ノンパラメトリックな回帰 / 確率密度推定– LOESS, TPSPLINE(Thin-Palate Spline)– GAM( Generalized Additive Model:評価版）– KDE (=Kernel Density Estimation)✸ （２次元の）空間統計– VARIOGRAM, KRIDGE,SIM2D✸ 層化・集落抽出など

および その解析– SURVEYSELECT,SURVEYMEANS– SURVEYREG

2

新しく追加されたプロシジャ(2)✸ NLMIXED– 非線形 + 一般化 混合モデル✸ PLS– PLS回帰✸ MI+ MIANALYZE ( 評価版 )– Multiple Imputation
BOXPLOT （箱ひげ図）✸ STDIZE（変数の様々な方法による標準化）✸

本発表のトピック– × NLMIXED →石塚先生が説明。
✸ 次のプロシジャを中心に説明。

– ○LOGISTIC, GENMOD– MIXED, FREQ
– ( 評価版のプロシジャ ) MI

3

LOGISTICプロシジャ✸ 数多くの拡張が行われました。

– (1) CLASS– (2) CONTRAST– (3) EXACT→ 質的変数→ 対比→ 正確な条件付き…。

LOGISTICプロシジャ（１）✸ CLASSステートメントの使用例✸ PROCLOGISTIC DATA=DATA1;✸CLASS Ａ / PARAM=GLM ;✸MODEL Y=Ａ Ｂ Ａ＊Ｂ;
✸ RUN;– 交互作用もＯＫ。

4

LOGISTICプロシジャ(1)✸ 複数のコーディング方法をサポート– 各水準の効果✸ a)PARAM=EFFECT… PROC CATMOD流✸ b)PARAM=GLM… PROC GLM流✸ c)PARAM=REFERENCE … 最後のセルを０– 多項式… 多項式✸ e)PARAM=ORTHPOLY … 直交多項式✸ d)PARAM=POLYNOMIALPARAM=EFFECT（デフォルト）✸ Ａα1 α2 α3✸ ---------------------------✸ １100✸ ２010
✸ ３001✸ ４-1 -1 -1

CATMODプロシジャで採用されているコーディング方法
パラメータ推定値の和が０になるような制約。
α1+α2+ α3+α4=0最後の水準における効果は、
α4= －（α1+α2+ α3)のようにして求める。

5

PARAM=REFERENCE✸ Ａ

α1 α2 α3 最後の水準を０とおく✸ ---------------------------- コーディング✸ １
100 最後の水準を基準とした✸ ２010 時の推定値が算出される。
✸ ３001✸ ４000PARAM=GLM✸ Ａ

α1 α2 α3 α4PROC GLMの方式✸ -------------------------------✸ １100 0✸ ２010 0✸ ３001 0 （注）交互作用がある時、PARAM=REFERENCE✸ ４000 1 と主効果に対する検定が異なる。

6

LOGISTICプロシジャ(2)✸ CONTRASTで、対比を指定できるように。

– PARMA= オプションに応じて、対比の指定方法は異なるので、注意が必要？
– expをとった値も計算することができる。
✸ →オッズ比の計算も可能になりました！

LOGISTICプロシジャ(3)✸ EXACTステートメント✸ 正確な条件付き分布に基づくロジスティック分析✸ Version8.1から。応答２水準のみ。

– 十分統計量で条件付けた分布のもとで、、、✸ ２つの検定のｐ値および “mid-p”値– exact probability test
– exact conditional scores test
✸ 推定( ESTIMATE オプション )7

最も単純な例（２×２表の解析）
データ-------------YX
Yes A群Yes A群No A群
Yes B群
No B群
No B群
No B群No B群

デザイン
計画行列Ｘの各列 x_i------------------------ に対応した
yx0x1 十分統計量t_i１11１11t_i = y` x_i
０11
切片x1の十分統計量
１10t0=y’ x0 = ３０
10０
10
効果x1の十分統計量
０10t1=y’ x1 = ２０10

最も単純な例（２×２表の解析）

A

Yes No 計
21 ３B

14５計

35８✸ 行和（赤色）だけでなく、列和（水色）も固定した時の(1,1)セルの分布

8ＥＸＡＣＴステートメントを手計算すると、、、
T1
何通り？
✸ ------------------------✸010
3C 0×5C3
130✸
3C 1×5C2215✸
3C 2×5C131✸
3C 3×5C0✸

2つの考え方＝2つの検定方法✸ (1)“フィッシャーの直接確率検定” 風– Prob(T1=u | T0=3)を計算していって、– Prob ( T1=2 | T0=3) 以下となる確率を足していく。
✸ (2)“コクラン-マンテル-ヘンツェル”風[ t1 － E(T1|T0=3) ]2✸ -------------------------------✸Var(T1| T0=3)（＊単変量の場合）✸– よりも大きくなるものが生じる確率を足していく。

9ObsxA
1234Count
0
10
130
215
3
1Score2.52000
0.031111.524447.00000

Prob0.17857
0.535710.26786
0.017860.26786以下の確率を足すと、、、
0.17857+ 0.26786 + 0.01786
1.52444以上になる確率を求めると、、
0.17857+ 0.26786 + 0.01786（＊２つの検定は、少なくとも２×２表の時は同じ結果）

Mid-p value✸ 条件付き正確検定✸ →ｐ値の分布は離散✸ → 保守的になる傾向がある。
✸ Mid -p … 「保守的な傾向」を“薄める”ために考え出された方法。
– 0.17857+ 0.26786 / 2 + 0.01786

10

Test
Score
ProbExact
Mid
Statistic PValue PValue
1.5244 0.4643
0.33040.2679 0.4643
0.3304

ESTIMATEオプション✸ 条件付き尤度に基づく推定値・信頼区間– ESTIMATE= BOTH, ODDS, PARM
✸ (1-α)%の信頼区間– 帰無仮説をαで棄却するようなβ0の値を。
–H0: β=β0- H1: β>β0– および H0: β=β0+ H1: β<β0+– （両側の場合は、それぞれα/2で）。

11Exact Odds Ratios (2-tail)
95% Confidence
Est
Limitsp-Valuex A 5.784
0.158 586.994
0.5714

Exact Odds Ratios (1-tail)One-sided 95%ConfidenceOne-sided
Est
Limitsp-Valuex A 5.784
0.245 Infinity
0.2857
Cf. 推定値… PHREG , 信頼区間 …FREQ /EXACT OR
One-sided p-Value … Fisherの正確検定

GENMODプロシジャ✸ “相関がある”二値応答データの解析✸ ２次のモーメントまでを考慮したＧＥＥ推定の他に、Alternating Logistic Reg(交互ロジスティック回帰）が行えるように。
– Type3検定が行えるように。
– 「相関構造」ではなく、「（対数）オッズ比」の構造を指定する。
– 基本は全被験者に対して共通のオッズ比構造

12

LOGOR= オプション✸ PROC GENMOD….;✸CLASS …; MODEL …..
✸/ DIST=BIN LINK=LOGIT;✸REPEATED✸SUBJECT=id(drug) /✸WITHIN = time✸
LOGOR = EXCH ;✸ RUN;LOGOR= オプション✸ LOGOR =EXCH– … すべての時点間で等しいオッズ比。
– （パラメータ数 : １個）✸ LOGOR = FULLCLUS– … すべての時点間でばらばらのオッズ比。
– （パラメータ数： t(t-1)/2個)✸ LOGOR= …ユーザ指定など...

13

対 数 オ ッズ の 構 造 を E xch an ga ble に し た 結 果標 準 誤 差 は Em p irica l
S tan d ar
E s tima te d Erro rPara met er
In ter cept

9 5%
C on fid en ce
Lim its

Z Pr > |Z |1.3 863

0 .790 6 -0. 1632

2.9 358

1. 75

0. 0795

1 .250 0 -2. 4500

2.4 500 -0. 00

1. 0000time

1-0. 0000

….

.

…

…

省略…

…

…

…

.

…

…

…

…

…

…

d rug

1

0.8 109

1 .317 6 -1. 7715

3.3 934

0. 62

0. 5383d rug

2

0.0 000

0. 0000…

0.0 000

0 .000 0.

.

d rug *tim e 1 1-1. 7918

1 .860 3 -5. 4378

1.8 543 -0. 96

0. 3355

d rug *tim e 1 2-2. 1609

1 .691 3 -5. 4757

1.1 540 -1. 28

0. 2014d rug *tim e 1 3-0. 2719

1 .696 2 -3. 5964

3.0 525 -0. 16

0. 8726….

.. . .

…

…

…

…

…

…

.. . .

…

…

省略…

…

.

.

0 .361 7 -1. 8296 -0. 4116 -3. 10

0. 0019

A lp h a 1-1. 1206…

S core S tat is tics For T yp e 3 GE E
A na lys is
C h iS q ua re Pr > C h iS q

S o urce

DFt ime4

1. 76

0. 7800

d ru g10. 01

0. 9327d ru g* tim e4

2. 28

0. 6842

FREQプロシジャ✸

EXACTステートメント✸ MCオプションによる乱数シミュレーション✸ 処理時間の上限を決定するMAXTIME=オプション✸ NPAR1WAYプロシジャも同様。

✸ 二項分布 Bin(n,p) の確率ｐに対する信頼区間が一元表の時にも計算できるように。
– （方法自体は、RISKDIFFオプションと同じ）

14

FREQプロシジャの変更点✸ 共通オッズ比と共通相対リスクのマンテル-ヘンツェル流の推定値– Version6までは、「帰無仮説のもとでの分散推定値」によって信頼区間を計算していた。
– → Version8で、ようやく改良。

MIXEDプロシジャ✸ 小標本時の標準誤差・検定– DDFM=KENWARDROGER✸ Kenward & Roger (1997)の論文に基づく方法によって、標準誤差・検定統計量（Ｆ値, ｔ値）・ ｐ値を計算。
✸ 注）自由度の計算だけでなく、固定パラメータ間の標準誤差なども調整される。

15

Kenward & Roger(1997)✸ REML推定における標準誤差の求め方– 「σの推定値」をプラグインして、– 応答変数Ｙの分散共分散構造 Σ(σ) の推定値–＝ Σ（σのL推定値）– （ＧＬＳ推定の枠組みに添って、）Σの推定値を既知のものとして、βの分散共分散行列Vを推定する。
– 欠点A) σの推定に伴う変動を考慮していない。
– 欠点B) Vにおける小標本に伴うバイアスを考慮していな
い。

Kenward & Roger(1997)
これらの欠点を補うために、近似を高める方法を提案。
✸ そして、Ｆ分布の近似も良くなるような方法も。既存の検定と対応させて、、、✸ 自由度１の時：サタースウェイト近似✸ 正確な検定の時：Hotteling T^2統計量✸ になるような計算を提案。
✸✸

シミュレーションではカイ２乗近似と比較。

16

MI & MIANALYZE(V8.1～ ただし、評価版)✸ ３種のMultiple Imputation をサポート– PROC MI✸ 乱数によって欠損値を補ったデータセット(complete dataset)を複数作成。Imputeする方法
としては、３種類。

– PROC MIANALYZE✸ 様々なプロシジャによる推定結果から標準誤差などを計算。

MIにおける３種類の方法✸ monotoneな欠測パタン– 1) 通常の回帰に基づく方法✸ regression method– 2) ロジスティック回帰に基づく方法✸propensity score method✸ 多変量正規に基づく– 3)MCMC法

17

Monotoneな欠測値パタン✸ ＩＤ

Y0 Y1 Y2 Y3 Y4✸ ---------------------------------✸ 001 ○○ ○ ○ ○
✸ 002 ○
○ ○ × ×✸ 003 ○
○ ○ ○ ×✸ 004 ○
○ × × ×– 「×の後ろは、必ず×」MAR✸ どの手法も、Missing At Randomでないと駄目。
– 欠測するかどうかが、その欠測部分の値によっては左右されない。
– 例えば、トービットモデルに従っているようなデー
タは駄目。→ある以上のものが欠測するようなデータなどは駄目。

18

方法１：回帰法✸ 次のモデルでY(1), Y(2), Y(3)… の順番に代入していく。
✸ Y(t) の欠測値部分=b0 + b１ Y(1) ＋ … ＋ b(t-1) Y(t-1)✸+ s×乱数✸– b0,b1,…,b(t-1)およびsは、Y(t)の非欠測部分から推定

方法２：propensity score
method「欠測」 or 「非欠測」 を応答変数として、ロジスティックモデルを推定。
✸ ロジスティックモデルの予測値が近いもの同士で、幾つかのグループに分ける。
✸ グループ化されたもののなかで、非欠測値を無作為抽出して、それを欠測値に代入。
✸19

方法３：多変量正規に基づく
ＭＣＭＣ法✸

何らかの方法で、Σの初期値を与える。
– デフォルトはＥＭアルゴリズムによる推定値✸ 「非欠測部分」および「Σの推定値」に基づき、欠測部分に値を代入。
✸ Σを、代入された部分を含めて推定する。
Ｂ）に戻る。
– B) ,C)の処理は事前に何回か行われる。

MIANALYZE✸ complete datasetsからの推定– 各データ毎に推定値Q(i)およびその分散推定値U(i)を求める。
✸ 点推定値= Q(i)の平均– within-imputation variance = U(i)の平均– between-imputation variance = Q(i)の分散✸ 点推定値の分散（total variance)✸= within ＋ (１＋1/m) between

20

Multiple-Imputation Variance Information
Variance
Variable

BetweenTotal DF

Relative
Increasein Variancex1

0.031937 0.065363 0.098897 165

0.513030

0.346931

x20.030456 0.031908 0.063887

1.002241

0.513231

Within76

Fraction
Missing
Information

GLM,REG,TTEST,UNIVARIATEでも...
✸ 信頼区間を出力する機能が追加。

✸ GLM: 線形従属の関係を表示できるように。
✸ TTEST: １標本,対応のあるｔ検定など。

21

Base SASでも...
✸ UNIVARIATE– トリム平均、ウィンザー化平均– MADなどによる標準偏差の推定– 信頼区間（平均・分散・パーセント点…)– ヒストグラム・ＱＱプロット（SAS/QCから移植）✸ MEANS– CLASSDATA= オプション– パーセント点も計算できるようになった。

その他の統計関連における拡張 (GUI)✸ Enterprise Guide (ｸﾗｲｱﾝﾄ-ｻｰﾊﾞｰ環境)✸ SAS/STATプロダクト✸ Analyst Application (英語）– 記述統計・グラフ– １標本・２標本の平均– REG･PRINCOMP･ＧLM・MIXED・LIFETEST・PHREG– １標本や２標本問題における検定・信頼区間・同等性などの検出力計算

22補足：ＯＤＳに伴う変更点（１）✸ ODSの導入により、outputウィンドウへの出力形式が一部変更されています。

–“PRINTTOプロシジャ”を利用したV6時代のプログラムが、うまく動かなくなる可能性があります。
– 対処方法： ODSを使えば、直接データセットに出力できるようになります。

補足：ＯＤＳに伴う変更点（２）✸ GENMOD,MIXEDにおける変更– Version8との互換性があるステートメント✸ ○MAKEステートメント– → しかし、変数名が変更されています。

– V8で廃止されたステートメントやプロシジャ✸ × %LET PRINT_=OFF;%LET_DISK_=ON;✸ × PROC OUTPUT → 将来のバージョンでは、これに相当した機能が正式サポートされる予定。

23

情報✸ 米国SUGIでのV8に関する論文– http://www.sas.com/service/techsup/news/– からダウンロード可能(PDF形式)✸ “What’s New in Version8”...
✸ “What’s New in Version7” …

✸ SＵＧＩ-Japan…

***********

以上です。

24