オッズ比と信頼区間からp値を出す。

論文投稿準備していたら、投稿規定に「p値を書け絶対だ」（意訳）というp値縛りがありました。

書けと言われたのだから、p値を計算して書けば良いのですが、投稿準備中の論文は「データ提供元の規約によって生データを削除済み」ということで、今からp値を計算できないという状況です。

ただし、既に作っている表には、オッズ比とその95%信頼区間のみがかかれています。

そこで、オッズ比と95%信頼区間からp値を算出することにしました。

使用する例

下記がStataコードと出力です。lbwデータセットを用いて、「妊娠中の喫煙」を曝露として、「出産児の低体重」をアウトカムとしたロジスティック回帰分析を行いました。

> webuse lbw
(Hosmer & Lemeshow data)

> logistic low i.smoke

Logistic regression                                     Number of obs =    189
                                                        LR chi2(1)    =   4.87
                                                        Prob > chi2   = 0.0274
Log likelihood = -114.9023                              Pseudo R2     = 0.0207

------------------------------------------------------------------------------
         low | Odds ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
       smoke |
     Smoker  |   2.021944   .6462989     2.20   0.028      1.08066    3.783112
       _cons |   .3372093   .0724103    -5.06   0.000     .2213694    .5136667
------------------------------------------------------------------------------
Note: _cons estimates baseline odds.

Odds Ratio=2.021944（95%CI: 1.08066, 3.783112）となっています。

表示されているp値は0.028ですが、より正確な値はp=0.0276197です。この値は、matlist r(table) を実行すれば確認できます。

計算方法

手順1．オッズ比を回帰係数にする

オッズ比と信頼区間の値を回帰係数の値に戻します。これは、logを取るだけです。

> di ln(2.021944)
.70405942

> di ln(1.08066)
.07757197

> di ln(3.783112)
1.330547

回帰係数=0.70405942 （95%CI: 0.7757197, 1.330547）と計算されました。

手順2．標準偏差を計算する

信頼区間の上限値（または下限値）と点推定値の差は、標準誤差のおよそ1.96倍（下記では、invnormal(0.975)で算出）になっているはずです。これを元に標準誤差を算出可能です。

> di invnormal(0.975)
1.959964

> di abs(1.330547-0.70405942)/invnormal(0.975)
.31964239

これで、標準誤差=0.31964239と算出が出来ました。

手順3．Z値を計算する

Z値は、回帰係数/標準誤差で計算できます。

> di 0.704059/0.31964239
2.2026459

シンプルに計算ですが、z=2.2026459となりました。
このz値は、logitコマンドで算出された値と一致しています（当たり前ですが）。

手順4．z値からp値を計算する

あとは、z値をp値に変換します。

> di (1-normal(2.2026459))*2
0.02761972

p=0.02761972と算出されました。これは、元々のp=0.0276197と一致しています。無事に推定できたようです。

adoファイルにまとめました。

この手順を踏めば、関数電卓を用いて算出することが可能ですが、手間がかかるので、Stata用のadoファイルを作りました。

エラー処理とかユーザビリティとか無視しているので、やや使い辛いですが、そのうち（今世紀中）には修正したいと思います。

odds_to_pvalというコマンドです。
引数を3つ必要で、点推定値・95%信頼区間下限・95%信頼区間上限です。

> odds_to_pval 2.021944 1.08066 3.783112
p-value = .0276196026000386

adoファイルをダウンロードして、sysdirコマンドで表示されるフォルダのどこか（仮置きならPERSONALがお勧め）に入れるとStata実行時に使える様になります。

. sysdir
   STATA:  C:\Program Files\Stata17\
    BASE:  C:\Program Files\Stata17\ado\base\
    SITE:  C:\Program Files\Stata17\ado\site\
    PLUS:  C:\Users\username\ado\plus\
PERSONAL:  C:\Users\username\ado\personal\
OLDPLACE:  c:\ado\

odds_to_pval.ado

606 Bytes

ファイルダウンロードについて

ダウンロード