単位量当りの効果を記述する | 医療統計とStataプログラミングの部屋

このブログでは、統計解析ソフトStataのプログラミングのTipsや便利コマンドを紹介しています．

Facebook groupでは、ちょっとした疑問や気づいたことなどを共有して貰うフォーラムになっています．ブログと合わせて個人の学習に役立てて貰えれば幸いです．

さて、今回はStataとは直接的な関係はないのですが、何かのモデルを走らせて得られた結果をどのように記述するか、ということを考えてみたいと思います．

二値変数の場合には、単純に何かの要因があった場合に、なかった場合と比較してどの程度推定値に変化があったかを記述すればよいだけなので簡単です．しかし、連続変数の場合はどうでしょう？そんな、単位変化量当りの推定値の変化について解説してみたいと思います．

1．単位変化量が整数倍に変化したとき
2．曝露因子が2倍になるにつれての推定値の変化を記述したいとき
3．報告の仕方の例
まとめ

1．単位変化量が整数倍に変化したとき

「BMIが1上がったときのハザード比は○○です」

くらいだったらそれほど問題ないと思いますが、

「収縮期血圧が1mmHg上がったときのハザード比は○○です」

といわれても1mmHgの変化は誤差範囲に感じられてしまいます．そんなときは

「収縮期血圧が5mmHg上がったときのハザード比」

として記述するとわかりやすいかもしれません．（5でも10でもよいと思いますが）

このとき、新たな変数を作成するか、推定値自体を計算することで欲しい結果が得られます．

例えば収縮期血圧を表す変数を”sbp”としたとき、以下の様に新たな変数を作成します．

.generate sbp10 = sbp/10

とすると収縮期血圧10mmHgあたりの変化をみることができます．そのままモデルに投入することで推定値が得られます．

すでにモデルに投入した場合はどうでしょうか．推定値をそのまま使う場合にはそのまま10をかけ算するだけでよいのですが、ハザード比やオッズ比では推定値をexponentiate（指数変換）しているので、それほど単純ではありません．

推定値を10倍スケールにしたときのHRをHR₁₀としてみます．このとき、

HR₁₀ = exp(10 × β)　→　exp(β)^10

となるので、exp(β)、つまりもとのスケールのときのHRを10乗することで得られます．

disp exp(0.03)^10
    /* 1.3498588 */

disp exp(0.3)
    /* 1.3498588 */

なお、これらの操作は、他の変数の推定値に影響を与えないので後からモデルを作り直してその部分の数値だけ入れ替えればよいです．

2．曝露因子が2倍になるにつれての推定値の変化を記述したいとき

あるマーカーが2倍に変化したときに、HRがどのように変化するか、というのは、変数に対して底を2にした対数変換をすることで推定値が得られます．

ところで、Stataではもともとln(), log10()が関数として使えて、それぞれ自然対数、常用対数として活用できるのですが、log2はどうすればいいでしょうか．

これは、単にその数字の自然対数を２の自然対数で（あるいは常用対数を２の常用対数で）割ればよいです．

これによって底を自在に変えることができます．（Statalist参照）

変数の作成方法は、一度自然対数に変換した変数を、そのままln(2)でわって新たな変数とすれば完了です．これも他の変数の結果には影響を与えないのでその結果だけ後から入れ替える事も可能です．

log2(x) = ln(x)/ln(2)

generate l2var = ln(var)/ln(2)

varという変数を変換して新たな変数を作成することができました．

3．報告の仕方の例

さて、そのようにして作成した変数を使った解析をどのように論文中に記載していくといいでしょうか．留学先であったCRIC研究から出た論文（J Am Soc Nephrol. 2014 Feb;25(2):349-60.）を見ながら表現を参考にさせてもらいましょう．

… higher FGF-23 was independently associated with graded risk of CHF (hazard ratio [HR], 1.45 per doubling [95% confidence interval (CI), 1.28 to 1.65]; HR for highest versus lowest quartile, 2.98 [95% CI, 1.97 to 4.52]) and atherosclerotic events (HR per doubling, 1.24 [95% CI, 1.09 to 1.40]; HR for highest versus lowest quartile, 1.76 [95% CI, 1.20 to 2.59]).
Abstract 　J Am Soc Nephrol. 2014 Feb;25(2):349-60.

「FGF23の値が2倍になるにつれて、○○だけHRが上昇します」

というのをこういった表現で表しています．それと、段階的にリスクが上がっていくというのをgraded riskと表現しておくというのもいい感じですね．

本文中には、以下の様な表現もあり、参考になります．

… ascending FGF-23 quartiles were associated with a stepwise increase in risk.
… the highest versus lowest quartile of FGF-23 was associated with a 1.63-fold greater risk of atherosclerotic events (95% CI, 1.03 to 2.59), and a 23% higher risk per doubling of FGF-23.
Each doubling of FGF-23 was associated with a 45% higher risk in the primary model that adjusted for …
Page 350, Results section

こうしてみると、結構表現にバリエーションをつけてくれていることに気がつきます．

ちなみに他にもこの論文は参考になるところがいくつかあって、フォレストプロットの描き方は一部お手本にしました．このようなお手本となる論文をじっくり読み込むことは、内容の理解もさることながら実際に論文を書く側の意識も高まっていきますね．