Syleir’s note

2020.4.1より統計検定やE資格の勉強の進捗を報告しています。統計検定準1級、E資格、G検定取得しました!当ブログへのリンクはご自由にどうぞ。

MENU

【統計的因果推論】回帰分断デザイン(Regression Discontinuity Design:RDD)を解説する

1. はじめに

因果関係の分析を行うとき、理想的には無作為化比較実験(RCT)を行うのが理想的です。
syleir.hatenablog.com

しかし、実際には適用が難しいケースが多々あります。RCTができない環境で、観察データからなんとか因果関係を見出そうとします。
そのような研究の中で、比較的因果関係を示しやすい、そんな手法が回帰分断デザイン(Regression Discontinuity Design: RDDです。

この記事では、RDDの基本概念から具体的な分析方法、注意点までを解説します。

2. 回帰分断デザイン(RDD)とは?

因果関係を示すには、RCTを行うのが理想です。上の図のように、母集団をランダムに2群に分け、処置群と非処置群(対照群)に分け、介入を行って因果関係を推定します。

RDDは、一言で言えば、閾値の前後で処置群と非処置群が分けられる状況で、因果関係を分析する手法、と言えます。

皆様、中高大とテスト漬けだったあの日々を思い出してみましょう。
あるいは、資格試験を思い出してみましょう。

合格点が60点だったとして、59点の人と、60点の人、知識量にあまり差はないと思いませんか?時の運、問題セットで覆る些細な差です。
なのに59点の人だけ再試・補習・浪人をさせられ、60点の人には何もありません。
不公平ですよね?59点だった時にこの再試は意味があるのか?と考えることもあるでしょう。そういう時に考えるのが回帰分断デザイン(RDD)です。

2.1 基本的な概念

回帰分断デザイン(RDD)は、特定の閾値を境に処置群と非処置群が分けられる状況で、因果効果を推定する手法です。

2.2 具体例

試験で合格点が60点と設定されている場合、59点と60点の受験者間に知識量の大きな差は存在しない可能性があります。しかし、59点の学生には補習が課され、60点の学生には何も行われません。この「閾値付近のランダム性」を利用して、補習の効果を推定するのがRDDです。

介入群(Treatment Group):60点以下の得点の学生
対照群(Control Group):60点以上の得点の学生

上の例では、
介入群:不合格者
対照群:合格者
です。

また、
介入:補習
因果効果:次回の得点上昇率

としています。

3. RDDの仮定

  • カットオフ周辺では、ランダムに介入が発生している

この例では、60点付近の実力を持つ受験者はランダムに61点、あるいは59点を取り、介入が起きていると考えています。
つまり、閾値付近の個体では、背景要因がほぼ同じであると考えることができます。

ここで、これらのカットオフ前後の試験者を抽出し、介入による点数差が生まれれば介入による効果があると考えることができます。

  • 介入の割り付けのルール、およびカットオフ値が明確

RDD閾値を超える(または下回る)場合に、そしてその場合に限って介入が導入されることを前提にしています。
介入の基準、それに対応する閾値RDDを行う前に特定していなければいけません。

  • 介入前の変数の操作ができない

割り当てするかどうかを操作できるような環境があってはいけません。
例えば、60点以下の人が増えそうだから、配点を変更するなどの操作があってはいけません。

  • カットオフ値の前後で潜在的な介入がない

例えば、60点以下の人に補習に加えて、塾・予備校に通うことを義務付けたりした場合、補習の効果のみを判定することはできません。
背景要因に不連続性がないようにデザインすることが好ましいです。

  • カットオフ値の前後において、アウトカムに影響を与えるそのほかの因子が大きく変わらない

同様に、介入以外の被験者の要因もカットオフ前後で変わりがないことが重要です。

4 RDDの分析

通常、ITSと同様に、折れ線回帰で分析を行います。

通常、結果はこのように表示します。
折れ線に有意なずれがある場合、介入効果ありと判断します。


ITSはRDDのより特殊な例と考えることができ、カットオフ前後を時系列で横軸を取っているのにすぎません。

いつかpythonでの実装を行います。

5.参考:曖昧なRDD

この例は60点以下の群に補習という介入を行っていますが、例によっては介入群に割り当てられても全例が介入を受けない場合もあります。
例としては、大学受験合格後の辞退、留学可能かをGPAで足切りした後の辞退などです。
このような例では曖昧なRDD(fuzzy RDD)ということを行います。

6. RDDの利点

  • 因果関係の推定が可能

すでに書きましたが、RDDは観察研究ですが、ランダム化比較試験(RCT)が不可能または非倫理的な状況でも、準実験デザインとして因果関係の推定ができます。

閾値付近では交絡因子がほぼ均等に分布していると仮定されます。
これにより、未測定の交絡因子の影響のバランスをとることができます。未測定の交絡因子のバランスをとることができることについては、RCTの強いメリットであり、それが活かせるのはこの研究デザインの大きなメリットです。例えば、傾向スコアマッチングなどではこれができません。

  • 仮定を満たしているかの確認が比較的容易

RDDの仮定がいくつか必要ですが、これらは、データから直接検証できます。
例えば、ヒストグラムを用いて割り当て変数が閾値付近で不連続性を示していないか確認できます。

Yusuke Sasabuchi, Introduction to Regression Discontinuity Design, Annals of Clinical Epidemiology, 2022, 4 -1, p. 1-5,

7. RDDのlimitation

  • 閾値以外に介入を行ったときの一般化ができない

例えば、70点の群に補習を行ったらどうかはわかりません。

  • どの範囲を取るかの選択が難しい

前後5点を取るか、10点を取るか、この辺にはドメイン知識が必要です。
体感的には55点と65点は差がないと言っても良いかもしれませんが、50点と70点は同じ背景と言ってはいけないような気がします。

  • 必要なデータ数が多い

カットオフ前後の値だけを取るので、必要なデータサンプルが増えます。精度の高い結果を得るには、十分なサンプルサイズが増えます。

8. まとめ

本記事では、回帰分断デザイン(Regression Discontinuity Design: RDD) について解説しました。RDDは、無作為化比較実験(RCT)が実施できない状況でも、比較的信頼性の高い因果効果の推定を可能にする強力な手法です。

10.関連書籍

因果推論の導入におすすめです。

最近出た本ですが、因果推論関連書籍の中でもかなり出来が良く、面白いです。
ぜひ。

ちょっとだけですが、RDDの記載もあります。

modern epidemiologyの和訳です。3ページほど記載がありますがちょっと読みにくいです。


関連記事です。
普段はこんなことを書いています。

syleir.hatenablog.com
syleir.hatenablog.com
syleir.hatenablog.com
syleir.hatenablog.com

【無刻印化】Keychron Q0 MaxをHHKBとお揃いにする【Part2】

syleir.hatenablog.com

この記事の後編です。

はじめに

前編では、Keychron Q0 Maxを購入したことをお話ししました。単体で十分高機能なテンキーですが、この記事ではいろいろ調整、魔改造していきます。

前回の記事では、一般的な方にもおすすめできる記事ですが、ここからはある種の思想、信仰、信奉、崇拝が入ります。ご注意ください。

 

 

無刻印化という思想

一般的にテンキーは「数字入力専用の補助デバイス」と認識されがちですが、Keychron Q0 Maxはその枠を大きく超えた存在です。4つのレイヤーに任意のキーやマクロを割り当てることができ、キーマップを柔軟に変更することで、ユーザーの思い描く理想の入力環境を構築できます。

特にNumlockキーの存在には、多くのユーザーが疑問を抱くのではないでしょうか。数字入力モードの切り替えに過ぎないこのキーは、多くの場合、不要であり、邪魔ですらあります。Keychron Q0 Maxでは、このような無駄を排除することが容易ですが、ここで新たな課題が浮かび上がります。

キーキャップの刻印問題です。

このKeychron Q0 Maxのキーキャップには記号が刻まれています。記号がメインのデザイン性には疑問を感じます。そして、キーマップ変更が前提である場合、キーの刻印そのものが無意味に思えてきます

 

キーマップが完全にカスタマイズ可能であれば、キーキャップに余計な情報を刻む必要はありません。むしろ、無刻印であることが美しさと実用性を両立させる答えとなります。これは、無刻印キーボードの代名詞とも言えるHHKBが示すひとつの完成形です。

HHKBの思想をリスペクトしつつ、Keychron Q0 Maxにおいても無刻印化を実現する。そのためには、以下の2つの工程が必要です。

  1. キーキャップの交換 
  2. キースイッチの換装 

今回は、Keychron Q0 Maxを無刻印化し、HHKBに通ずる思想を反映させたカスタマイズを実行します。

テンキーという枠を超えたデバイス、Keychron Q0 Max。その可能性を追求し、理想のカスタマイズを目指します。

1.キーキャップの交換

キーキャップの交換自体はそこまで難しくありません。
キーキャップを外して、つける。それだけです。Keychron Q0 Maxにはキーキャップ、キースイッチの引き抜き器具がついていますが、やや使いにくいので、こちらがおすすめです。

 

1.1.キーキャップの選択肢

ここで、無刻印を実現するキーキャップの選択肢について検討します。できれば市販のもので、かつ国内で完結すると嬉しいと思ってました。

  1. HHKB Studio 無刻印キーキャップセット
    信心を発揮するならこれですが、問題はテンキーの+キー、Enterに対応するものがないことです。HHKBはテンキーレスキーボードなので、テンキーにのみ存在するキーキャップがなさそうです。また0キーが2Uサイズですが、HHKBには2Uのキーがないので、選択肢としてはイマイチです。

     

     

  2. Majestouch用 キーキャップセット
    こちらはHHKBと同じPBT素材で、2024年2月に発売された比較的新しいキーキャップです。HHKBと比較するとややつるつるした肌触りですが、色味としてはややHHKBより黒く、Keychronの黒にはより映えるキャップです。テンキー付きなので、こちらを利用して改造することにします。
  3. Razer Phantom Keycap Upgrade Set

    自分は使いませんが、Keychron Q0 Maxにはバックライトがついているので、このようなRazer Phantomシリーズのように、バックライトがついているときだけ刻印が見え、そうでない時は消失する半透明キーキャップも選択肢にはなるかなと思います。完全無刻印はやや恐ろしい方々のためのものです。私は無刻印新教なのでこちらを認めることができますが、正教の方々からはバッシングの嵐かもしれません。黒もあります。

     

これらの選択肢を吟味して、今回は、Majestouch用 キーキャップセットを注文しました。かなり質感がよく、おすすめです。

こちらが交換の途中になります。もともとのキーキャップはカーブがきつく、やや高いですが、比較してみると、ややなだらかで、低いです。

自分は長時間使用するので、キーキャップは低い方が手首が疲れず、ありがたいです。

ただ、もともとのキースイッチも25000円するテンキーのキースイッチなので、良いのですが、あと一息といったところで、静音赤軸なので静かで打ちやすく、軽いのですが、より静かでスコスコ叩けるキースイッチが欲しいと思ってしまう自分がいました。

そこで、さらにキースイッチ換装も行うことにしました。

2. キースイッチの換装

もともとデフォルトの静音赤軸もまずまずの打鍵感ですが、理想のテンキーに近づけていきます。

キースイッチの換装は、Keychron Q0 Maxを自分だけの理想形に仕上げるための重要なカスタマイズです。標準搭載の静音赤軸は十分に静かで軽い打鍵感がありますが、自分が選択したのはHHKB Studio軸です。

 

2.1 HHKB Studio軸とは?

HHKB Studio軸は、その打鍵感と静音性で高く評価されているスイッチです。柔らかく、それでいてしっかりとした押し心地があり、長時間使っても疲れにくいのが特徴です。

  • 打鍵感: 吸い付くような柔らかい打鍵感

  • 静音性: 夜中に使っても周囲に気を使わなくて済む静かさ

 

HHKB studioはHHKBが静電容量無接点方式からメカニカルスイッチを採用した最新のキーボードですが、市販のキーボードの中では最高級の部類に分類されます。メカニカルキーボードですが、その打鍵感は評判です。ただし、静電容量無接点方式の打鍵感とは似て非なるものなので、また別の良さ、というような感じです。

 

 

間違いなくよいキースイッチであることには疑念の余地がありませんが、問題は高価なことです。

キースイッチは10個単位から販売されております。通常のキースイッチと違って少数販売であり、かつお値段1個あたり3300円です。Keychron Q0 Maxはありがたいことにたくさんのキーが搭載されておりますので、全キーを交換するには30個購入が必要になります。9900円になります。高いですよね?

でもよく考えてみましょう。PFU様は、、故障用の交換を目的としてキースイッチを販売してくださっています。もともとキースイッチ交換用の販売ではないのです。これがなければHHKB Studioを購入し、分解してつけなければならないところでした。HHKB Studioは44000円です。なんと1/4。お値打ちですねPFU様、ありがとうございます。

ありがたく購入させていただきましょう。これが信心です。

なお、インターネットでは、こちらのKailh Deep-Sea SilentがかなりHHKB Studio軸の打鍵感に近いという評判です。お値段は単価で1/4位まで抑えることができます。

 

 

ただ、自分は買わずにHHKB Studio軸しか買っていないので、比較はできません。これが信心です。スイッチも黒ってかっこいいんですよね。

2.2 キースイッチの換装手順

  1. キーキャップを外す

  2. キースイッチを外す

  3. 新しいスイッチの取り付ける

  4. キーキャップを取り付ける

これだけです。実際には器具さえあれば、苦労なくできます。簡単にスイッチ交換ができるキーボードをホットスワップ対応と言いますが、本当に簡単です。ぜひやってみましょう。

2.3 換装後の感触

換装後のKeychron Q0 Maxは、打鍵感がまるで別物になりました。指先に伝わる柔らかさと反応の良さ、静音性のバランスが絶妙です。

もともとの静音赤軸も静かですが、体感半分程度の音量です。音量は対数で認識されるとされます。ってことは数字以上にすごいってことです。

換装後は自然とテンキーに手が伸び、作業が楽しくなるほどの打鍵感です。また、見た目が本当に美しくなりました。愛用機のHHKBと並べてもこの統一感で素晴らしいです。自分はデバイスは黒派ですが、HHKB雪モデルとKeychron Q0 Max白でもかなり良い構成になるのではないでしょうか。

余談ですが、自分はテンキーは左側に置く派です。トラックパッドやマウスが右側、テンキーが左側、キーボードが正中。これが人間工学です。フルサイズのキーボードはこれができなくて困りますね。これができるのが、キーボードとテンキーを分離させる最大のメリットです。

3. 総費用

さて、懺悔を込めて総費用を供養しておきます。

Keyboard:Keychron Q0 Max ¥24,970-
Switch:HHKB studio軸 ¥3,300 * 3 = ¥9,900-
Keycap:Filco Majestouch 無刻印 ¥4,027-
Keycap puller(Option):¥4,683-

合計¥43,580-

数字の入力をテンキーで行っておりますが、あまりの入力のしやすさに涙が出てきます。涙の理由はそれだけではないような気がしますが笑

向こう10年はこの装備で戦っていけそうです。

まとめ

Keychron Q0 Maxは、単なる「数字入力デバイス」ではなく、快適な作業環境を提供するパートナーへと進化しました。

  • 無刻印キーキャップ: 見た目がシンプルで、どんな配置でも違和感なし。

  • HHKB Studio軸: 理想の打鍵感、静音性。HHKBへのお布施に最適。

カスタマイズによって、Keychron Q0 Maxはもはや「テンキー」の枠を超えた存在で、万能な左手デバイスとしても使用できます。ちょっと手をかけると、デバイスにも愛着が湧きます。ぜひ皆様も最高の相棒を作ってみてください。

HHKBの布教記事はこちらです。

syleir.hatenablog.com

 

普段は統計に関する記事を書いています。

syleir.hatenablog.com

 

 

【散財】最高のテンキーを購入した【Keychron Q0 Max】

はじめに

ここ数年、HHKB Professional HYBRID Type-Sという無刻印キーボードを使用しています。3年前に購入してから、全ての記事、コードをこのキーボードでタイピングしており、大きな問題なくここまで来ていました。

 

syleir.hatenablog.com

 

 

このキーボードには文字がなく、タッチタイピング(無理やり)習得し使うことになりますが、数字が連続するときにホームポジションから遠い数字キーを連続して打つため、やや正確性が落ち、ストレスになります。

例えば、研究やスライド作成でのデータ入力、暗証番号や口座番号への入力などです。特に、文字が表示されないタイプのパスワードの入力などは、可能ではありますが、毎回若干の不安に駆られながら入力を行っています。さらに、今回、Blenderという3Dモデリングのソフトを導入し、テンキー操作が必要となったので、外付けテンキーを探していました。

条件としては1つ、HHKBに見合うくらいかっこいいテンキーであることです。
今回、Keychron Q0 Maxというテンキーを購入したので、布教して行きたいと思います。

 

こんな人におすすめ

  • HHKBユーザーでデザイン性の高いテンキーを求めている人
  • データ入力や3Dソフト、動画編集で効率を上げたい人
  • 打鍵感や静音性を自分好みにカスタマイズしたい人

競合との比較

Keychron Q0 Max は、2024年6月28日に発売されたテンキーです。2024年というのは、高級テンキー業界に革新が起きた年で、ここ数年まともな製品が出ていなかった、高級路線のテンキーが、Keychron Q0 Maxの他にも登場した年です。

その一つが、Realforce RT1 テンキーボードです。こちらのテンキーも、Realforceから、プレミア価格が付いているRealforce 23U(2009年発売)の後継品として、15年ぶりに発売されたテンキーで、テンキー界隈では話題沸騰しており、あまりの人気に品薄となり、入手困難となっています。一応1月には品薄が解消される予定のようです。

prtimes.jp

 

比較表

項目 Keychron Q0 pro Realforce RT1
キー数 26キー+ロータリーエンコーダー 23キー
端子 USB-C USB-A
無線対応

あり

(Bluetooth 5.1/2.4GHz対応)

なし
ケーブル 取り外し可能 本体一体型ケーブル(約80cm)
キースイッチ カニカルスイッチ 静電容量無接点方式
カスタマイズ性 高い(ホットスワップ対応) 低い
静音性 スイッチの選択次第 良好
キーキャップ PBTダブルショット PBT昇華印刷
バックライト あり(オフ設定可能) なし
重量 684g 300g
価格 24,970円 13,200円

どちらのキーボードを選択する?

個人的には、打鍵感、打ち心地については是が非でも求めたいところであり、HHKBに慣れてしまった自分としては、同様のキースイッチ方式である静電容量無接点方式を採用するRealforce RT1を使用したいところではありますが、どうしても延長も短縮もできない80cm固定の一体型ケーブルと、旧石器時代のUSB-A端子が令和の時代に強制されることが許せません。ケーブル断線で高級キーボードがおじゃんになるのは辛いので、今回はKeychron Q0 proを採用することにしました。値段については目を瞑ることにしました。今回は持ち運びせず、自宅での使用を想定しているので、重量については気にならず、キースイッチの選択で、可能な限り静音性がよく、打鍵感の良いキーボードに魔改造してやろうと思い至りました。

 

Keychron Q0 Max を実際に使ってみた感想

Keychron Q0 Maxを導入してから数週間が経ちましたが、期待以上の製品でした。HHKBと並べて使っても違和感なく、デザイン性の高さが最大のポイントです。さらに機能面でも、テンキーに求めるものをすべて満たしてくれました。テンキーに求めるものは0−9までの数字が楽に打てればいいだけなのでハードルはすごく低いのですが。

1. デザインが良い

Keychron Q0 Maxは、CNC加工のフルアルミニウムボディを採用しており、質感が高いです。HHKB Professional HYBRID Type-S の墨モデルの色と並べても、フレームの黒との統一感がよく、非常に美しいデスクに仕上がっています。

2. 使いやすい

  • USB-C / Bluetooth / 2.4GHz対応
    有線接続はもちろん、無線接続にも対応しているため、作業環境に合わせて柔軟に使えます。充電も50時間程度持ち、そこまで頻回の充電を要すことはなく、USB-Cの充電器は机にたくさん転がっているので、充電も楽でした。

  • 取り外し可能なUSB-Cケーブル
    Realforce RT1の固定ケーブルとは違い、Keychron Q0 Maxではケーブルが着脱可能です。断線の心配もなく、長期間の使用にも安心感があります。

  • ロータリーエンコーダー
    スライダー操作など、テンキー以上の役割を果たせるロータリーエンコーダーが思ったより便利でした。Blenderでの3D作業でも役立ちますが、一番役立っているのはYouTube視聴時の音量調整です。

3. カスタマイズ性とホットスワップ対応

Keychron Q0 Maxのもう一つの大きな強みは、ホットスワップ対応 である点です。
スイッチの交換が簡単にできるので、静音性と打鍵感を自分好みにカスタマイズできます。静音性を意識したスイッチを選択することで、静電容量無接点方式に引けを取らない打鍵感のよさと静音性にを手にいれることができました。カスタマイズについては、次回の記事で書こうと思います。

4. 静音性

HHKB Professional HYBRID Type-Sの静電容量無接点スイッチには静音性と軽いタッチが特徴ですが、Keychron Q0 Maxのデフォルトの赤軸も、非常に静音性がよく、スイッチの交換は実は全くマストではありません。遜色ない静音性が実は元から付いています。デフォルトでは、他にバナナ軸、茶軸の選択が可能ですが、静音性については赤軸が最も優れていると思います。

5. テンキーとしての実用性

Keychron Q0 Maxは、テンキーに加えて9キー追加されており、レイヤーを切り替えてアプリやソフトごとに割り当てを変更したり、マクロ導入をすることができ、以下のようなシチュエーションで便利です。一般的な左手用デバイスに準じたことが、テンキーで可能になります。

  • 数値入力作業
    エクセルやデータ入力が驚くほど快適になります。
  • パスワード登録
    セキュリティにやや懸念がありますが、パスワードをマクロ登録しておくことで高速にPCやサイトへのログインができます。
  • Blender
    ショートカットを割り当てることで操作効率が上昇します。
  • プログラミング
    マクロを導入し、よく書くコードを登録できます。

6. 価格

Keychron Q0 Maxは、24,970円と決して安価な製品ではありません。自分も当初はその値段にビビり、テンキーにその値段を払ったことに対し冷静になったことは幾度となくあります。コストパフォーマンスが心配になったことは幾度となくあります。しかし、コストパフォーマンスについては、満足度が非常に高ければ、それで改善するのです。Keychron Q0 Maxは、その機能と品質を考えると満足度は非常に高いです。値段には目を瞑って一度購入してみましょう。

まとめ

Keychron Q0 Maxは、その高級感、カスタマイズ性、機能性において、HHKBユーザーが満足する数少ないテンキーといえます。イキって無刻印キーボードを購入したけど、数字はちょっと不安というそこのあなた、生活が変わります。
無線対応ホットスワップ対応など、Realforce RT1にはない柔軟性を備えており、現代の作業環境に合ったテンキーと言えるでしょう。

Keychron Q0 Max は、デスク環境を美しく、機能的にしてくれるテンキーです。HHKBと一緒に使うことで、作業効率がさらに向上することは間違いありません。

 

 

ここまでは、比較的万人向けの記事です。次回は、このテンキーの魔改造について、思想たっぷりの記事を書こうと思います。

続きはこちら

syleir.hatenablog.com

左室収縮末期壁応力(ESWS)を解説:心臓超音波検査における理論と実践

はじめに

どうもお久しぶりです。いつもは統計やDSについて書くことが多いですが、今回は医学分野、特に心臓超音波検査について書いていきます。


本記事では、新生児医療や循環管理において重要な左室収縮末期壁応力(ESWS: End-Systolic Wall Stress)について詳しく掘り下げます。
特に新生児医療において、ESWSは、後負荷を考える指標としてよく用いられます。

左室収縮末期壁応力( \sigma)は以下の式で表されます:

\sigma = \dfrac{1.35 \cdot P_{\text{es}} \cdot D_{\text{es}}}{4 \cdot h_{\text{es}} \cdot \left( 1 + \dfrac{h_{\text{es}}}{D_{\text{es}}} \right)}

近年、新生児医療では、ESWSとmVcfc(the rate-corrected mean velocity of fiber shortening)を組み合わせて循環管理を行うことが提案されています。
しかし、ESWSは少し難しい概念です。
本記事では次の疑問に答えていきます。

  1. なぜESWSが後負荷の指標となるのか?
  2. 応力とは何か?
  3. ESWSの計算式はどのように導出されるのか?
  4. 実臨床での適用とその課題は?

応力とは?

応力とは、物体内部の単位面積あたりに働く力を指します。また、応力と抵抗の関係性を考えると、応力は血管壁や心筋に対する負荷を表す一方、抵抗は血液の流れに対する障害を示します。両者は密接に関連しており、血管系全体の機能を理解する上で重要な要素となります。例えば、血管抵抗が増加すると、心臓は同じ血流量を維持するためにより大きな応力を発生させる必要があります。

左室における後負荷とは,左室が収縮し血液を拍出しようとする際に直面する抵抗,すなわち心室壁応力のことです。これを直接見るための指標が、ESWSです。ESWSは、左室収縮末期の応力を表し、心臓の後負荷を定量化する重要な指標となります。この指標を用いることで、心臓の機能状態をより正確に評価し、適切な治療方針を立てることが可能になります。

ラプラスの法則とESWSの関係

ここで、唐突にラプラスの法則を登場させます。
ラプラスの法則とは、流体を含む容器の壁にかかる応力と、容器の半径および内圧との関係を説明する物理法則です。

\sigma = \dfrac{P \cdot r}{2h}

ここで、

  •  \sigmaは壁張力(応力)
  •  Pは内圧
  •  rは容器の半径
  • hは壁の厚さ

この法則は、心臓の壁応力を理解する上で重要な基礎となり、ESWSの概念にも密接に関連しています。

ここで、エコーでは、半径ではなく、直径で見るので、半径でなく直径で表すことにします。半径 r を直径 D の半分として表現すると、ラプラスの法則は以下のように書き換えることができます:

 \sigma = \dfrac{P \cdot D}{4h}

  • ただし、 D = 2r

ラプラスの法則とESWSの式の比較

 \sigma = \dfrac{1.35 \cdot P_{\text{es}} \cdot D_{\text{es}}}{4 \cdot h_{\text{es}} \cdot \left( 1 + \dfrac{h_{\text{es}}}{D_{\text{es}}} \right)}

この式は、ESWSの式と比較すると、基本的な構造が似ていることがわかります。主な違いは、ESWSの式には補正係数や追加の項が含まれており、より精密に心臓の形状と機能を反映していること、および、es = end-systolic = 収縮末期という特定のタイミングを考慮しているということです。                                                 

ラプラスの法則の導出の際、容器が薄壁で球形または円筒形であり、壁の厚さが半径に比べて十分に小さいということを利用していますが、新生児の心臓は小さく、壁の厚さが半径に比べて無視できません。これを考慮したのがESWSの式であり、実際にESWSの式において、壁の厚さが半径に比べて十分に小さいという仮定を利用すると大局的にESWSはラプラスの法則に収束します。

壁の厚さは直径に対しても十分に小さいということを仮定すると、直径に対する壁の厚さの比( h/D)は非常に小さな値になります。このとき、ESWSの式の分母にある(1 + h/D)の項は、ほぼ1に近似できます。つまり、この項を無視すると、ESWSの式はラプラスの法則とほぼ同じ形になり、1.35という係数を除けば、基本的な構造が一致することがわかります。1.35の意味についてはまた後ほど示します。

まとめると、ESWSは、ラプラスの法則の壁厚が半径に比べて無視できない場合を考えていることになります。

式を定性的に評価してみる

さて、ここでESWSの式を定性的に評価してみましょう。

 \sigma = \dfrac{1.35 \cdot P_{\text{es}} \cdot D_{\text{es}}}{4 \cdot h_{\text{es}} \cdot \left( 1 + \dfrac{h_{\text{es}}}{D_{\text{es}}} \right)} = 1.35 \cdot \dfrac{P_{\text{es}} \cdot D_{\text{es}}}{4 \cdot h_{\text{es}}} \cdot \dfrac{1}{\left(1 + \dfrac{h_{\text{es}}}{D_{\text{es}}} \right) }

となるのでした。ESWSは後負荷を直接的に表す指標でしたから、ベースのラプラスの法則に近似できる部分を考えると、後負荷は、圧力、直径に比例し、壁厚に反比例することがわかります。

後負荷の指標として、簡易的に収縮期血圧を利用するのも、後負荷が圧力に比例していることが理由なのでした。

また、高血圧などで、血圧が高くなっている際、心臓側で後負荷を下げるためにできることとして、直径を減らし、壁厚を増やすことです。

心臓の求心性肥大

これが、本態性高血圧や大動脈弁狭窄症などにおける求心性肥大の説明として古典的に教科書でよく扱われてきた説明でした。

ESWSの式の導出

Grossman W, Jones D, McLaurin LP. Wall stress and patterns of hypertrophy in the human left ventricle. J Clin Invest. 1975;56(1):56-64.

こちらの論文に記載されているモデルで考えます。

心臓を楕円球と仮定し、真ん中の対称線で二分割します。(一意ではない表現ですが、許して下さい)断面は球である仮定をします。

モデル図

図のように、単位面積あたりにかかる応力が均一に \sigma_m、心臓内圧がP 、中心から内壁までの距離 R_i、中心から外壁までの距離 R_o、壁厚h( = R_o - R_i)とします。

そうすると、この半楕円球にかかる垂直方向の力を考えると、応力の和=内圧の垂直方向の和が成り立たなければいけません。

応力の和

応力の和は、単位面積あたりの応力が均一に \sigma_mであること、および、その応力がかかる部分の面積が、下図の青色部分の面積となります。

(外側の円の面積)ー(内側の円の面積)であることから、この面積は
 \pi R_o ^ 2 - \pi R_i ^ 2 = \pi (R_o ^ 2 - R_i ^ 2)
となるので、応力の和は、
 \sigma_m \pi (R_o ^ 2 - R_i ^ 2)
となります。

圧力の和

一方、心臓の壁にかかる圧力による、図下向き成分の力の和は、
内面積にかかる圧力として計算して良いので、(導出は省略します)
 P  \pi R_i ^ 2となります。

式変形

よって、力の釣り合いから、
 \sigma_m \pi (R_o ^ 2 - R_i ^ 2) = P  \pi R_i ^ 2
が成立し、
 \sigma_m =  \dfrac{P  R_i ^ 2}{R_o ^ 2 - R_i ^ 2}
と計算できます。これを計算していくと、
  = \dfrac{P  R_i ^ 2}{\left( R_o - R_i \right) \left(R_o  + R_i \right) }

 h = R_o - R_iを利用して、R_oを消去して、

  = \dfrac{P  R_i ^ 2}{h \left(2 R_i + h \right)}

  = \dfrac{P  R_i ^ 2}{ 2 h R_i \left( 1 + \dfrac{h}{2 R_i} \right) }

  = \dfrac{P  R_i}{ 2 h \left( 1 + \dfrac{h}{2 R_i} \right) }

です。ここで、収縮末期(end systolic:es)における状態を考えると、
 P P_{es} h h_{es}となります。
また、超音波検査では、半径ではなく直径を考える(例:LVDd, LVDsなど)ので、収縮末期の直径をD_{es}
とすれば、
 2 R_iD_{es}となります。

よって、収縮末期の張力を \sigma_{es}とすれば、
 \sigma_{es} = \dfrac{P_{es} \cdot  \dfrac{D_{es}}{2}}{ 2 \cdot h_{es} \cdot \left( 1 + \dfrac{h_{es}}{D_{es}} \right) }
 = \dfrac{P_{es} \cdot D_{es}}{ 4 \cdot h_{es} \cdot \left( 1 + \dfrac{h_{es}}{D_{es}} \right)} となります。

実臨床に落とし込む

さて、ここまできたらあとは一息です。
実臨床でこれを扱うためには、測定可能なものを代入しないといけません。
 h_{es}D_{es}については心臓超音波検査における傍胸骨左縁長軸断層像をMモードを利用して、描出することで、
 h_{es}→LVDWs (left ventricular posterior wall end-systolic thickness:収縮末期左室後壁厚)

D_{es}→LVDs(left ventricular end-systolic diameter)

と置き換えることで計算できます。
P_{es}については、観血的動脈圧測定をAラインなどで行っていれば、波形から計算可能ですが、

Rowland DG, Gutgesell HP. Use of mean arterial pressure for noninvasive determination of left ventricular end-systolic wall stress in infants and children. Am J Cardiol. 1994;74(1):98-99.

この論文において、収縮末期血圧は新生児、小児において平均血圧(MAP)で代用して良いことが示唆されておりますので、Aラインがない状態でも、非観血的動脈圧測定における平均血圧で代用します。

平均血圧の単位は \mathrm{mmHg}ですが、ESWSの単位は \mathrm{g/cm^2}で表しますので、相互変換するために 1\ \mathrm{g/cm^2} =  1.35\ \mathrm{g/cm^2}であることを利用して、
(760 mmHg = 1013 hPaなどを利用すると簡単です)


ESWS = \dfrac{1.35 \cdot P_{\text{es}} \cdot D_{\text{es}}}{4 \cdot h_{\text{es}} \cdot \left( 1 + \dfrac{h_{\text{es}}}{D_{\text{es}}} \right)}

として計算できるようになりました。

ESWSの課題とモデリングの欠点

以上考えてきたように、ESWSは
・楕円球を考えた単純なモデリングであること
・壁厚が均一であること
・局所的な応力の不均一な分布を無視していること
・VSDなどの心室内シャントの影響を考慮していないこと
・収縮末期の応力のみを後負荷の指標と考えていること
・LVDs, LVDWsが正しく測られていない可能性があること
・収縮末期血圧を平均血圧で代用していること

などがlimitationの可能性たりうることがわかります。例えば、これはverifyされていますが、新生児期の肺高血圧が強い楕円形の左心室でも同様の推論が可能かなど、一向の余地がある前提が多いです。
一方で、心室内形態も考慮した、直接的な後負荷を考えるための指標であることには間違いなく、これらの前提を理解して使えば強力なツールとなり得るでしょう。

実際の使用法について

実際の使用法・管理については、語るに耐えませんので素晴らしい成書や論文に譲らせて下さい。

これらの本が詳しく書かれている成書かと思います。

生理学を深掘りした本では、こちらが最近読んで面白かったです。


おわりに

本記事では、ESWSの理論的背景と実臨床での応用について解説しました。
ESWSは後負荷を直接評価できる有力な指標であり、特に新生児医療や心疾患管理において有用です。しかし、その使用にはモデルの限界や測定誤差を理解し、適切に応用する必要があります。
理論、前提から見つめ直すことで、どのような脆弱性があり、どのような時に正当化されるかを考えることができます。少しでも参考になればと思います。では。


普段はもう少し疫学よりのことを書いています。
syleir.hatenablog.com

syleir.hatenablog.com

通算10万PV到達しました!【大感謝】

はじめに

いつも当ブログをお読みいただきましてありがとうございます。こんにちは、Syleirです。

 

 

先日当ブログは通算10万PVを達成することができました。これもひとえに読者の皆様ひとりひとりが当ブログに訪れてくれること、また2回3回と見にきてくださっているおかげです。本当にありがとうございます。

この記事ではこれまでにやってきたことやこれからやっていきたいことをつらつらと書いていきたいと思います。

 

目次

 

syleir.hatenablog.com

こんな記事を書いていたこともありました。この頃から3年近く経過し、コンスタントに1日100PVみていただけるようになりました。ありがたいことです。アクセスに繋がる記事は前回から大きく変わっています。

通算記事数

この記事を除いて、今までに書いた記事数は89記事になりました。どれもこれも思い入れがありますが、統計検定に関わる記事を多めに書いてきました。

syleir.hatenablog.com

記事数は多くないですが、記事の平均字数は10,000近く、これは果たしてブログと呼んで良いものかどうかわからないくらいの記事のボリュームがあります。

PVに貢献している記事Top3

1位

syleir.hatenablog.com

この記事は統計検定準1級に合格した時、どのように勉強したか、実際の内容について、まとめたものです。また、最後の現地試験であり、かつ過去最難関のセットだったものの実地レポートです。かなり力を入れて書いたので、これはたくさんアクセスしてくださってとても嬉しいです。こちらは前回から不動の1位です。統計検定の対策ってそう変わるものではありませんので時流の流れに動じずコツコツとやっていくのが大事です。

2位

syleir.hatenablog.com

これが2位になるの、このブログの読者層はどうなっているんでしょうか。褒めています。

統計的因果推論の中でもかなりコンテンツとしては奥底のコンテンツになりますが、これを解説した記事が人気です。確かに意外と難しいので、これを手軽に解説した記事ってあまりないんですよね。

でも、本当はこの次の記事、

syleir.hatenablog.com

これが書きたくてこの記事を書いたんです。伸びて欲しい記事が伸びないの、インターネットあるあるですね。

 

3位

syleir.hatenablog.com

3位は比較的因果推論の中でも初歩的な分野です。これが伸びるのは納得です。かなりまとまっているのでおすすめです。

 

2位、3位の記事の関連やさらなる学習はこの本がおすすめです。

 

到達までの期間

10万PVに達するまで、当ブログの開設日は2020年4月なので、延べ4年2ヶ月かかっていることになります。

大学の学部が終わってしまうのかなりびっくりですね。

なんで続いているのか自分でもわかりません(笑)何度も失踪しかけていますが、なんとか続けられています。読者の方々からのご質問やコメント、大変励みになっています。ありがとうございます。

 

今後の目標

今後の目標としては、やはり200000PVになりますが、なんとしても、充実して質の良いコンテンツを提供することを心がけていきたいと思います。また、日本語の質や良質な言語化、図解についても引き続き意識していきたいです。

読者のご希望があれば統計検定のコンテンツも出したいですが、今後は自己学習も兼ねて、因果推論よりの解説記事を追加していきたいとも考えております。

書籍には出せない強みを継続して出せるように頑張ります。

 

おわりに

改めて皆様のご支援・ご愛読に心から感謝しております。これからもSyleir's noteをどうぞよろしくお願いします。

 

【生存時間解析】ログランク検定の実践【Part6】【統計検定準1級•1級】

はじめに

ついに想定している最終回です。
前回まではカプランマイヤー曲線の作成を行ったのでした。
syleir.hatenablog.com

今回は

これらをベースに咀嚼してから執筆しています。ぜひ原著にも当たってください。

Kaplan-Meier法で生存曲線が書けるようになった

カプランマイヤー法で生存曲線を推定して描出したものがカプランマイヤー曲線でした。
カプランマイヤー曲線は、実際のデータから生存関数を推定したいというモチベーションで作成されます。治療群と対照群で生存曲線を作った時に、それらの生存曲線は全く同じ比でイベントや打ち切りが起きていない限り、当然違うものになります。
この生存曲線の違いが、統計的に有意なものかどうかを知りたいというのが、自然な次のモチベーションです。

ログランク検定

ログランク検定(Log-rank test)は「生存曲線に差があるか」をテーマとした検定です。
「生存曲線に差がある」とはどういうことでしょうか。

2本の生存曲線に差があるか?

明らかにこれらの生存曲線は当然違います。ここでいう、生存曲線に差があるかの「差」というのは、「2本の曲線がどれくらい違うか大域的に示すような値」のことであって視覚的な差のことではありません。
つまり、もう少し補うと、生存曲線に「統計的な」差があるかです。この統計処理にもたくさんの種類がありますが、その一つがログランク検定になります。


ログランク検定は生存曲線に差があるか?という問いに答える検定です。
検定というからには帰無仮説と対立仮説があります。
ログランク検定における帰無仮説は、「全時点で、2つの生存関数が等しい」です。
この帰無仮説
・打ち切りは群と独立に起こること
・観測された総イベント数、期待されるイベント数が大きいこと

を仮定するとログランク検定が実現されます。

ログランク統計量

ログランク検定における、「2本の曲線がどれくらい違うか大域的に示すような値」のことを「ログランク統計量」といいます。
このログランク統計量は自由度1の\chi ^2分布に従うことが知られていますので、ログランク統計量の導き方さえわかれば、あとは統計検定2級、あるいは統計検定準1級前半の知識で解くことができます。

分割表によるカイ二乗検定を復習する

ワークブック28章にありますが、変数の値の組み合わせごとの頻度を表にしてまとめたものを分割表と言います。一旦脱線しますが、今後の解説のために一度復習しておきます。
既知の方はスキップして構いません。

生存 死亡 合計
新薬 16 4 20
プラセボ 12 8 20
合計 28 12 40

この図は2×2の分割表で、

\dfrac{(観測度数-期待度数)^2}{期待度数}
をすべての群で足し合わせたものが、
が、漸近的に自由度1の \chi^2分布に従うことが知られています。

期待度数とは?

期待度数というのは、先ほどのこの表を

生存 死亡 合計
新薬 20
プラセボ 20
合計 28 12 40

このように空欄にしたとき、空欄をもっともらしく埋めた値のことです。
例えば、新薬群のうち、生存したものは、40*(20/40)*(28/40)とするのが期待度数です。
同様に、他のものも埋めていくと、このようになります。

生存 死亡 合計
新薬 40\times \dfrac{20}{40} \times \dfrac{28}{40} 40\times \dfrac{20}{40} \times \dfrac{12}{40} 20
プラセボ 40\times \dfrac{20}{40} \times \dfrac{28}{40} 40\times \dfrac{20}{40} \times \dfrac{12}{40} 20
合計 28 12 40

となり、計算すると
期待度数としては、

生存 死亡 合計
新薬 14 6 20
プラセボ 14 6 20
合計 28 12 40

となります。
これから、 \chi^2 = \dfrac{(観測度数-期待度数)^2}{期待度数}
を計算すると、
 \chi^2 =  \dfrac{(16-14)^2}{14} +\dfrac{(12-14)^2}{14} + \dfrac{(4-6)^2}{6} + \dfrac{(8 - 6)^2}{6}
 =  \dfrac{4}{14} +\dfrac{4}{14} + \dfrac{4}{6} + \dfrac{4}{6}
 =  \dfrac{4}{14} +\dfrac{4}{14} + \dfrac{4}{6} + \dfrac{4}{6}
 =  \dfrac{4}{7} +\dfrac{4}{3}
 =  \dfrac{4}{7} +\dfrac{4}{3}
 =  \dfrac{40}{21}  \fallingdotseq 1.90

となります。

観測値 B_1 B_2 合計
A_1 x_{11} x_{12} x_{1\cdot}
A_2 x_{21} x_{22} x_{2\cdot}
合計 x_{\cdot 1} x_{\cdot 2} x_{\cdot\cdot}

と抽象化した時に、
 \chi^2 = \dfrac{x_{\cdot\cdot} (x_{11} x_{22} - x_{12} x_{21})^2}{x_{1\cdot}x_{2\cdot}x_{\cdot 1}x_{\cdot 2}}としても同じになることも知られています。

これらの \chi^2をピアソンのカイ2乗検定統計量と言います。
これが帰無仮説のもとでカイ2乗分布に分布収束することが示されています。
詳細は、

7.4 適合度検定をご参照ください。証明もあります。

この例では、
 \chi^2 = \dfrac{40 \times (16 \times 8 - 4 \times 12)^ 2}{20 \times 20 \times 28 \times 12}
 \fallingdotseq 1.90
として同じ値が出てきます。

より正確には、フィッシャーの正確検定により検定を行います。

ログランク検定とは一体何をやっているのか?

ログランク検定はノンパラメトリックな検定である

まず最初に理解する必要があるのは、ログランク検定は生存曲線の関数の形を仮定する必要のないノンパラメトリックな検定ということです。実際のデータをもとに、ログランク統計量という、「2本の曲線がどれくらい違うか大域的に示すような値」を算出しますが、ここに曲線がどのような関数であるかを指定する余地はありません。

ログランク統計量の計算方法

上に書いたイベント数が多いこと、という過程を早速無視して、手計算可能なレベルで以下の例からログランク検定をやってみましょう。

t 新薬 死亡数 打ち切り リスクセット プラセボ 死亡数 打ち切り リスクセット
1 0 0 6 2 0 6
2 2 0 6 1 0 4
3 0 0 4 1 0 3
4 0 1 4 0 0 2
5 1 0 3 1 0 2
6 1 0 2 0 0 1
7 0 0 1 0 0 1

この例を使ってログランク統計量を計算していきましょう。

①どちらかの群のイベント発生数と期待イベント発生数の差をイベント発生時点ごとに計算

イベント発生ごとに分割表を作成します。

t = 1

⚪︎イベント発生数

生存 死亡 合計
新薬 6 0 6
プラセボ 4 2 6
合計 10 2 12

⚪︎期待イベント発生数

生存 死亡 合計
新薬 6
プラセボ 6
合計 10 2 12

これを埋めて、

生存 死亡 合計
新薬 5 1 6
プラセボ 5 1 6
合計 10 2 12

となります。t = 1での新薬群の実際のイベント数と期待イベントの差は0 - 1 = -1です。

上の分割表のピアソンのカイ2乗検定統計量と違って、4箇所すべてを計算する必要はありませんが、常に同じところを計算し、実際のイベント数と期待発生数の引く順番も常に揃える必要があります。
今回は(死亡, 新薬)かつ実際のイベント数ー期待発生数で統一して計算します。
揃えてさえいれば、どの群、どちらの引き算で計算しても、正負が逆転するだけで、同じ絶対値の点数になるのでどちらでも構いません。

t = 2

⚪︎イベント発生数

生存 死亡 合計
新薬 4 2 6
プラセボ 3 1 4
合計 7 3 10

⚪︎期待イベント発生数

生存 死亡 合計
新薬 5.2 1.8 6
プラセボ 2.8 1.2 4
合計 7 3 10

t = 2での死亡群の実際のイベント数と期待イベントとの差は2 - 1.8 = 0.2です。

t = 3

⚪︎イベント発生数

生存 死亡 合計
新薬 4 0 4
プラセボ 2 1 3
合計 6 1 7

⚪︎期待イベント発生数

生存 死亡 合計
新薬 3.43 0.57 4
プラセボ 2.57 0.43 3
合計 6 1 7

t = 3での死亡群の実際のイベント数と期待イベントとの差は0 - 0.57 = -0.57です。

t = 4

t = 4打ち切りはありますが、イベントはないので、計算を行いません。

t = 5

飽きた人はスクロールしてください。
⚪︎イベント発生数

生存 死亡 合計
新薬 2 1 3
プラセボ 1 1 2
合計 3 2 5

⚪︎期待イベント発生数

生存 死亡 合計
新薬 1.8 1.2 3
プラセボ 1.2 0.8 2
合計 3 2 5

t = 5での生存群の実際のイベント数と期待イベントとの差は1 - 1.2 = -0.2です。

t = 6

⚪︎イベント発生数

生存 死亡 合計
新薬 1 1 2
プラセボ 1 0 1
合計 2 1 3


⚪︎期待イベント発生数

生存 死亡 合計
新薬 1.33 0.67 2
プラセボ 0.67 0.33 1
合計 2 1 3

t = 6での生存群の実際のイベント数と期待イベントとの差は1 - 0.67 = 0.33です。

t = 7

イベントが発生していないので計算は起こりません。

②これらの値の総和をとる

次に、これらの値の総和を取ります。
この教科書で定義される通りに文字を定義していくと、
i = 1, 2(群:1 = 新薬, 2 = プラセボ)
f = 1, 2, 3, 4, 5(時間:イベントが起きたタイミングの昇順、今回はf=1,2,3,4,5t=1,2,3,5,6に対応)
m_{if}(群i, 時間fでの実観測イベント数)
e_{if}(群i, 時間fでの期待イベント数)
となり、上で計算したことをまとめていくと次のようになります。

t f m_{1f} e_{1f} m_{1f} - e_{1f}
1 1 0 1 -1
2 2 2 1.8 0.2
3 3 0 0.57 -0.57
4 0 0 0
5 4 1 1.2 -0.2
6 5 1 0.67 0.33
7 0 0 0
合計 4 5.24 -1.24

さらに文字が増えてごちゃごちゃしますが、
「群iの、観測度数ー期待度数の値をすべてのイベントが起きた時間について合計したもの(-1.24のことです)」
 O_{i} - E_{i} = \sum_{f} (m_{if} - e_{if})と定義すると、
群1(新薬群)についての  O_{1} - E_{1} = -1.24
となります。

③分散を計算する

ここまできたらログランク統計量まであと一息です。
定義としては、
Log-rank statistic  = \dfrac{(O_1 - E_1)^2}{Var(O_1 - E_1)}
です。(O_1 - E_1) は計算済みですから、Var(O_1 - E_1)を計算すればよく、これは観察度数ー期待度数の分散のことです。

この分散は、
 \hat{Var}(O_1 - E_1)
 = \sum_f \dfrac{n_{1f}n_{2f}(m_{1f} + m_{2f})(n_{1f} + n_{2f} - m_{1f} - m_{2f})}{(n_{1f} + n_{2f})^2 (n_{1f} + n_{2f} - 1)}

として推定が可能です。こちらの教科書の教科書のシグマの対象の添字にはミスがありますので持っている方はご参照ください。
各イベントが発生した時間fについて、
n_{1f}:新薬群のリスクセットの人数
n_{2f}プラセボ群のリスクセットの人数
m_{1f}:新薬群の死亡数
m_{2f}プラセボ群の死亡数
です。

当てはめるだけなので簡単に計算できますが、手計算は人間の仕事ではなく、現実的ではないので、簡単に計算でき、出題が想定される例をご紹介します。

各時間でイベントが1人にしか起こらないとき

fでイベントが1人にしか起こらないとは、すべてのfm_{1f} + m_{2f} = 1のことです。

ここではm_{1f} = 1, m_{2f} = 0としましょう。
 \dfrac{n_{1f}n_{2f}(m_{1f} + m_{2f})(n_{1f} + n_{2f} - m_{1f} - m_{2f})}{(n_{1f} + n_{2f})^2 (n_{1f} + n_{2f} - 1)}
=  \dfrac{n_{1f}n_{2f} (n_{1f} + n_{2f} - 1)}{(n_{1f} + n_{2f})^2 (n_{1f} + n_{2f} - 1)}
=  \dfrac{n_{1f}n_{2f} }{(n_{1f} + n_{2f})^2}
= \dfrac{n_{1f}}{n_{1f} + n_{2f}} \dfrac{n_{2f}}{n_{1f} + n_{2f}}

ここで、 p_f =  \dfrac{n_{1f}}{n_{1f} + n_{2f}}と置くと、 1 - p_f =  \dfrac{n_{2f}}{n_{1f} + n_{2f}}
となり、

= p_f(1 - p_f)
となり、これは二項分布 B(1, p_f)の分散に相当しますから、
各イベント発生時点でのイベント数が、1例の場合にはその群の期待イベント発生確率をパラメータとする2項分布の分散をイベント発生時点ごとに計算し、合計したものが観察度数ー期待度数の分散になる
ことがわかります。
これくらいであれば手計算でもやる気の出る範囲内であり、これは2016年の統計検定1級で既出です。

分散を期待度数で近似するとき

現実には、各fでイベントが1人にしか起こらないことはありません。
もっとも、プログラミングをできない環境にいることはさらにありえません。手計算でログランク検定をやらなくてはいけないのは、統計検定の会場にいる時か、インフラのない無人島に漂着した時くらいです。でも残念ながら、我々はこのケースに該当してしまうのでやらなければなりません。

一般にログランク統計量
Log-rank statistic  = \dfrac{(O_1 - E_1)^2}{Var(O_1 - E_1)}は、分散を使わずに、
 \fallingdotseq \sum_i \dfrac{(O_i - E_i)^2}{E_i}
 =  \dfrac{(O_1 - E_1)^2}{E_1} +  \dfrac{(O_2 - E_2)^2}{E_2}(2群の場合)

近似できることが知られています。
O_1 - E_1 = O_2 - E_2
ですし、群すべてのイベント数の平均をEとすると
 E = \dfrac{E_1 + E_2}{2} = \dfrac{O_1 + O_2}{2}
が成り立ちますから、
 E_2 = O_1 + O_2 - E_1
 = 総イベント数 - E_1
とするともう一度計算しなくても簡単に求まります。
この例では
 E_2 = 9 - 5.24 = 3.76
です。

④ログランク統計量を計算し、検定する

Log-rank statistic  \fallingdotseq \sum_i \dfrac{(O_i - E_i)^2}{E_i}
 =  \dfrac{(O_1 - E_1)^2}{E_1} +  \dfrac{(O_2 - E_2)^2}{E_2}
 =  \dfrac{(-1.24)^2}{5.24} +  \dfrac{(-1.24)^2}{3.76}
 =  \dfrac{(-1.24)^2}{5.24} +  \dfrac{(-1.24)^2}{3.76}
 \fallingdotseq 0.29 + 0.41 = 0.60

となります。近似したものにせよ、近似しないものにせよ、漸近的に
ログランク統計量は自由度1のカイ2乗分布\chi(1)に従うことが知られており、
\chi_{0.05}(1) = 3.84 < 0.60ですから、有意水準0.05で検定した時、
「本例での2本の生存曲線に差がないという帰無仮説は棄却されない」という結果になります。

ログランク統計量とピアソンのカイ2乗統計量を見比べてみる

ログランク統計量
Log-rank statistic  \fallingdotseq \sum_i \dfrac{(O_i - E_i)^2}{E_i}
と近似できるのでした。

ピアソンのカイ2乗統計量についても
 \chi^2 = \dfrac{(観測度数-期待度数)^2}{期待度数}
をすべての群で足し合わせたもの
ですから、全く同じ形をしています。
実は、ログランク統計量というのは、すべてのイベント発生時について、分割表を書き、総和を取っていった後に計算してみるとピアソンのカイ2乗統計量の形に近似できるので、自由度1のカイ2乗分布に従う、というように解釈できます。

おわりに

このパートで当初想定していた、ログランク検定までの初歩の記事執筆はおしまいになります。
余力のある方は比例ハザードモデルあたりも理解していただければと思います。
今回もなんとか失踪せずに書き上げることができました!
少しでも初歩の理解につながっていれば幸いです。

そういえば、なかなかわかりやすい和書がなかったのですが、

最近出て話題のこの本、今日届いたのですが、生存時間解析のパートわかりやすくて良かったです。
あれ、この記事いらなくないですか???

ここまで読んでくださった方で、お時間ある方、
あと、フィードバックのコメントをたくさんお待ちしています。
次はこの記事書いてほしいとかのリクエストも(可能な範囲で)対応します。
コメント、励みになるのでぜひよろしくお願いします。

syleir.hatenablog.com
統計検定対策に役立つ記事のまとめです。他の記事も気合を入れて作っているのでぜひ読んでください。

はてなブログproに課金することにしました【大赤字】

はじめに

このブログも2020年4月に開始して早いもので4年が経ちます。ブログ立ち上げ当初とは自分の生活は大きく変わりましたが、なんとか失踪せずに執筆を続けられています。これもサイトに訪問してくださって、記事をコメントをいただける皆様のおかげです。ありがとうございます。

はてなブログproに課金

今日、はてなブログproに課金することにしました。はてなブログproは月1000円の課金で、2年分前払いをすると月600円まで下げることができます。本ブログの広告での収益は小学校の遠足のおやつ代程度で損益分岐点に達していないですが、趣味としてはかなり楽しめていること、また記事を書いていくことで自分のoutputにもつながっている実感があることから課金することとしました。

はてなブログproでできること

はてなブログproになって嬉しいことは、トップページで自分の記事が一覧表示できるようになります。今までは、最新の記事が全部表示され、全部終わった後にまたその前の記事が全部表示され、というような形で、記事の視認性がかなり悪く、検索性に欠け、読者の方々に不親切な設計でした。課金することで、このような表示となり、過去記事を読みやすくなっていると思います。
ぜひトップページに飛んで読んでみてください。

記事の一覧表示ができるように!

自分の記事は一般的なブログとは違って、1記事10000字を超えることが珍しくなく、一覧表示ができないとかなりストレスフルな形になっていたと思います。これを機にぜひ色々な過去記事を読んでいただき、少しでも赤字を軽減してくださると助かります笑

そのほか独自ドメインを作ったり、10個まで他のブログを作ったりできるようになるようですが、現状は時間がないのでその予定はありません。

あとプロフィール画面の自分の名前の隣にproが付きます!かっこいいですね!

つまり一覧表示のため(と一部の広告を消すため)に課金していることになります。なんと殊勝な心がけでしょうか。

今後の目標

とりあえず、課金してしまったからにはあと2年は継続して記事を出し続けられるように頑張ります。
もう少しで累計100000PVなので、そこは目標に頑張りたいです。
また、この記事で累計88記事目なので、100記事は頑張って達成したいところです。
いつもお読みいただいてありがとうございます。今後ともよろしくお願いします。