« 眠り姫問題のプログラム | トップページ

2017年4月20日 (木)

過剰診断論に関するシミュレーション 中間報告

「過剰診断論」にいくつか疑問を持ち、シミュレーションを書いてその疑問を明らかにしたいと考えた。その過程で、がん患者には、基礎寿命とがん寿命があって、どちらかが尽きたとき死ぬというモデルを考えることになった。

まだシミュレーション開発の途中だが、「基礎寿命論」についてある程度開発したところで、行き詰まりを感じたため、とりあえず現在の進捗を報告することにした。
過剰診断論


NATROM 氏によると「過剰診断論」とは「治療しなくても症状を起こしたり、死亡の原因になったりしない病気を診断すること」と定義するようだ。最近では、 NATROM 氏は甲状腺がんについて話題にしていた。

がんと診断された群は、当然に何らかの治療を受けているだろう。そうした診断群における予後が、そうでない罹患群と比較してよくなってないどころか、むしろ悪くなっているというとき、診断は、がん治療の侵襲や、がんと知ったことによる生活の変化(QOL の悪化)などの「悪い結果」をもたらしうる…というのが、過剰診断の問題らしい。

治療によって延命されるならまだしも、そうではないのに余計な診断をしている……というのが過剰診断論の批判のように私はとらえた。だったら、それは「診断」が問題なのではなく。治療が有効でないことの問題ではないかと思うのだが、(英語の Wikipedia の over diagnosis の項目ではそういう批判も書いているのだが、)あくまで「診断」を問題にするようだ。

そこで NATROM 氏は「リードタイムバイアス」や「レングスバイアス」を紹介する。

>
検診で発見された時点(長井さんの場合は60歳)から、検診を受けなかったら自覚症状が生じて診断されたであろうという時点(長井さんの場合は 65歳)までの期間を「リードタイム」という(長井さんの場合はリードタイムは5年間)。検診がまったく予後を改善しなくても、リードタイムの分だけ生存期間は延び、よって生存率も改善する。これを「リードタイムバイアス」という。

検診で発見可能になってから自覚症状が出るまでの期間(遅山さんは62歳~ 68歳で6年間)が長いがんは検診でより発見されやすい。そしてこうしたがんは進行がゆっくりなので予後が良い。実際には検診が予後を改善しなくても、もともと予後のよいがんが検診でより多く発見される傾向があるがゆえに、がん検診が予後を改善するように誤認するのがレングスバイアスである。

リードタイムバイアスおよびレングスバイアスをわかっていれば、「がん死」であろうと「甲状腺全摘」であろうと、予後を改善させない無効な検診が検診で発見された患者における「有害アウトカム」が生じた患者の割合を見かけ上小さく見えさせることを理解できるだろう。


予後が改善しているように見えても、それは見かけだけかもしれないとして、過剰診断論を補強し、がん検診が予後を改善しないことが知られている甲状腺がんについて、2011年3月11日の東日本大震災の福島第一原発の事故を受けての、がん検診の徹底を NATROM 氏は批判していた。

NATROM 氏にならって、疑い深くなってみてそれを NATROM 氏に向けると、「診断」があっても「治療」をしなければ、重大な害は生じにくいので、「治療」ではなく「診断」を問題にするのは、「診断」が増えると困ることがあるからではないかと疑うことができる。

診断が増えることで、統計的に何かが有意になることはないだろうか。例えば、治療が効かないことがはっきりする。広域被害の有無が確定する。最初は治療に意味がないようにみえて、のちに医師の技能差があらわになる。…とか。

ただ、ここでこの文だけを読んでいる人に注意しておくと NATROM はニセ科学批判クラスタに属しているとみるべき人物だが、原発関連企業から何かの利益をもらってポジションを取るような人物ではないということだ。

統計的に明らかになる可能性があることには、原発震災の被害がなかったことが確定するということもありうる。それは、もしかすると被災者には酷な内容を含意するかもしれない。そういうことを NATROM 氏が案じていた可能性もあるだろう。

疑いを別方向に向けてみよう。それは、過剰診断論が、そもそも注目を浴びるためになされた議論ではないかというものだ。過剰診断論があるということは、生きているうちは非顕在に終った疾病について比較対象としてデータが得られているということであろう。本当は、統計的に余命が延びることが示したくてそういうデータをとったのだが、それが有意には示せなかった。そのため、そのデータをとったことにつき正当性を主張しようとして出てきたのが、センセーショナルな過剰診断論だったのではないか。…という疑いだ。

この疑いを示すには、余命について有意にならない状況で、過剰診断については有意性が主張できるようなことがなければならない。そのようなことはありうるのだろうか。

これらの疑いについてシミュレーションをして確かめてみようというのが、本稿の目的だった。診断が増える効果について、余命の有意性と過剰診断の有意性について、何か言えないかというのが本稿のシミュレーションをするキッカケだった。

しかし、後述のように、それらを示すまでにはまだ致っていない。


基礎寿命論


シミュレーションで、最終的には統計をとってどうこうしたいわけだが、その前に、人が、がんに罹ることとはどういうメカニズムかについての内生的なモデルが必要である。

あくまでも内生的なモデルなので、人に寿命をまず決めることができる。そして、がんに罹患することでその寿命が減り、治療されることで治療群は罹患非治療群よりは寿命が増えると「普通」は考える。「普通」というのは、過剰診断論の設定では、治療群なのに寿命が増えていない状況を考える必要があるからだ。「診断」そのものに効果はなく、あくまで「治療」にのみ効果があると考える。寿命が変化するのはこれだけで、あとは、寿命が決まったものに属性を加えていくだけだと考える。例えば、診断をいつ受けるかという属性の割り当てが必要であるが、本稿では「診断」まではモデル化していない。

そして、「がん」をどのようにモデル化するかであるが、上のレングスバイアスの議論にあったように、がんには進行の遅いもの、早いものがあって、過剰診断論が示すように、がんが末期に進行することなく死を迎えることもある。……ということで、がんには「がん期間」が短いものや長いものがあり、死亡時どれくらい、がんが進行していたかを示す「死亡時がんステージ」があるとモデル化して考える。ちなみに、がんには四段階または五段階で示す「ステージ」論があるが、それとここでのステージは関係なく、ここでは 0.0 から 1.0 の間をとるとする。ちなみに、死亡までの闘病期間は、「がん期間」×「死亡時がんステージ」で表せることになる。

ところで、この「がん期間」や「死亡時がんステージ」というものが、実際に観測しうるものかという問いには、否定的に答えざるを得ない。四段階または五段階のステージに関連させて、死亡時に 0.0 から 1.0 の間に、がんステージを特定することはできるかもしれない。しかし、そこから、がん期間を特定できるかというと難しいだろう。病変した組織の経過から、まれに、がん期間も推定できることもあるかもしれないが、普通は、生前の別の検診時の、がんステージの特定があって、その二点のデータでやっと、がん期間が決まるとすべきだろう。しかし、それでは、治療群では特定できる可能性はあるが、罹患非治療群では、ほぼ特定できる可能性がないことになる。

モデル化の議論に戻ろう。ここからが本稿の特徴である基礎寿命論である。がんによる死には、がん寿命による死と基礎寿命による死があると考える。がんが末期になって死ねば、がん寿命が尽きたと考える。末期にならずに他の原因で死んだのなら基礎寿命が尽きたと考える。基礎寿命と、がん寿命のうち短いほうが実寿命となると考える。

進行の遅い、がんの場合、がん期間が長くがん寿命が長いわけで、基礎寿命による死が多くなるはずである。それは死亡時がんステージがどちらかというと 0 に近い段階で死を迎えるということだ。

また、治療群は、侵襲度の高い手術を受けることになるかもしれない。それは基礎寿命を減らす方向に作用すると考える。やはり、これも死亡時がんステージを 0 に近づける方向で作用するだろう。治療が実寿命を増やす効果がない状況、すなわち過剰診断論の設定の状況でも、この基礎寿命減があるとすれば、そこには元は寿命増の効果があったのではないかと推測できる。

実は、この治療群における基礎寿命減が、すなわち QOL の低下であり、過剰診断論が本当は示そうとしていることではないかと、このシミュレーションをしながら、私は思うようになった。過剰診断論は、偽陽性の問題…「偽陽性が少なくないのに侵襲度の大きい治療をしてしまうことが問題」というのとは違うそうだが、ある種の「侵襲」を問題にしているのではないかと私は考えるからだ。

この基礎寿命論を考慮した内生的モデルのシミュレーションが本稿「中間報告」の到達したところである。


シミュレーションの実際


実装は、プログラミング言語 Perl で行った。Perl プログラムから確率分布の数値を出力し、数学用作図ソフトの gnuplot で頻度分布の図を描いて確かめるのが今回のだいたいの流れだった。

まず、はじめに作ったのが、厚生労働省が配布する生命表の Excel データから死亡確率の分布を Perl で利用しやすい形にして出力する make_cumulative_mortality.pl と、その分布を使い逆関数法でランダムにある年齢に対する余命(寿命)を出力する test_cm.pl (cm は cumulative mortality の略) だった。ここを、擬似的に数学的な分布ではなくリアルなデータから取ったことが、その後、数学的な解析の道をはばみ、シミュレーション一本を軸とするように舵をとらせることとなった。

次に、ランダムに発生させた寿命を適当に変数変換して罹患群と非罹患群のデータを作ろうと考えた。そこで、Y がランダムに決められた余命として、p Y (p = 0.9 ぐらい) が罹患群の余命としたい。罹患群の割合は n/N とすると、p * (n/N) * Y + q * ((N-n)/N) * Y = Y となるように非罹患群の q を求めればいいのではないかと考えた。

そして test_cm2.pl を作ったが、結果、そうはなっていなかった。平均に関する議論では、E(Y) = (n/N) * E(p Y) + ((N-n)/N) * E(q Y) は言えるようであるが、分布についてはそれは言えないようだった。p = 0.1 とか極端な数値にして分布を見ると、p Y と q Y は別の山を作っていた。

考え直して、エージェント発生のアルゴリズムを見直し、必ず元の分布になるような選び出し方ということで test_cm3.pl を作った。寿命 y に対するある関数 Q(y) について rand(1) < Q(y) であれば罹患群とし、そうでなければ非罹患軍とするようなしかたをすれば、発生させた y は罹患群か非罹患群になるわけで、これを逆にして罹患群と非罹患群を合わせれば、必ず元の分布に戻るはずである。

その数学的意味であるが、まず、Y の密度分布を f(y) とする。P(Y=y1 and U < Q(y1)) = P(Y=y1) * P(U < Q(y1) | Y=y1) = f(y1) * Q(y1) となると思われる。ただし、自信はない。このように密度分布が積の分布になるとしても、累積分布関数が 0 から 1 の範囲になるための定係数をどのように計算してらいいのかわからないからだ。(いちおう、積の分布であることを確かめるための test_cm3_1.pl を作ったり、寿命の分布以外の分布の例として指数分布について積の分布になっていることを確かめる test_cm3_2.pl を作ったりした。)

とにかく、そのような Q(y) で変形して、罹患群の寿命の平均が下がり、治療群の平均が罹患群よりも上がるというエージェントの作り方ができればよい。罹患群に関する Q(y) をどのような関数形にすべきか迷って、最初は減少関数に指数関数 1 - C * (exp(x/T) - 1) や log(exp(x) + 1) 型のものを考えたがうまくいかないのでやめた。結果、test_cm3_3.pl で試した Q(y) を正規分布型の関数としてかけるか、Q(y) を直線的関数…「 ̄\_」という型の関数としてかけるかでよいのではないかと考えるようになった。そして治療群については、Q2(y) を直線的関数…「_/ ̄」という型の関数をさらにかけることにした。

正規分布型の場合は、罹患が認められるのある程度トシを取ってからで、年齢が上の部分につても、それ以外と比較して含まれる率が少なくなると考えることになろう。直線型の場合は、罹患が認められるトシになっているかどうかは考えず、高齢の者が含まれにくいことのみ表されるということになるだろう。

そのような関数形に適当な数値を与えることで、元の寿命の平均が 82.9…だったのを、罹患群は直線型のとき 79.9…に、正規分布型のとき 78.3…になった。治療群は直線:直線のとき 81.7…に、正規分布:直線のとき 80.2…になった。リアルな分布を元にし、あやふやな理論を使っているため、これらは経験的にしか求まっていない。変化をもう少し大きな数値にしたかったのだが、このような関数形での変型では、これぐらいが相場のようだ。

そして、基礎寿命論を加味する test_cm4.pl を作った。ここでは、がんで死ぬかどうかをまずランダムに決定し、非治療群ががんで死ぬ場合は、実寿命からの余命を基礎寿命とし、実寿命をがん寿命とし、死亡時がんステージは一律に 1 とした。非治療群ががん以外で死ぬ場合は、実寿命を基礎寿命とし、がん期間を一様乱数に定数をかけたものとして求め、それについてがん寿命を決定し、死亡時がんステージを 1 から 0 の間の乱数とした。治療群ががんで死ぬ場合は、がんで死ぬ場合は、実寿命からの余命をランダムに発生させて、基礎寿命論に沿うように、その分布が非治療群の余命の分布より余命が減る方向にゆがむようにして、それを基礎寿命と設定した。他は、非治療群と同様である。

しかし、test_cm4.pl のモデルはうまくいっているように見えなかった。実寿命に比べて、基礎寿命の平均が大きくなりすぎたり、がん寿命による死の割合を、治療群について少なくすると、過敏に反応して、基礎寿命が少なくなり過ぎたりした。

そこで、全面的にモデルを見直したのが、test_cm4_1.pl のモデルである。 test_cm4.pl では、がんで死ぬかどうかを死亡時がんステージとは独立に乱数を発生していたのであるが、test_cm4_1.pl では、死亡時がんステージをまず発生させ、それが 0.8 以上である場合、がん死であると判断するようにした。がん期間についても、test_cm4.pl では最大がん期間に対する一様分布にしていたが、test_cm4_1.pl では平均がん期間に対する指数分布に改めた。

そして基礎寿命減は、死亡時がんステージが 0 に近くなることで表現しようとした。がん期間が長くなると、その分、がん以外で死ぬ可能性が大きくなる、すなわち基礎寿命が減る、すなわち死亡時がんステージが 0 に近づく効果があると考える。がん期間を t、がんステージを s とすると ds/dt は負になるべきであろう。一方、がん期間とがんステージの積で表される闘病期間は、がん期間が増えるほど増すべきだろう。f(t) = t * s とすると、df/dt は正になるべきだということになる。

これらを満たすには、平均死亡時がんステージが、s(t) = a/sqrt(x + b) + c という関数形であればよいと考えた。死亡時がんステージの平均を動かすためには、一様分布から変形しないといけない。test_cm4_1.pl ではそれを右台形分布・左台形分布で表すことにし、test_cm4_2.pl では、それをベータ分布で表すことにした。test_cm4_1.pl で左右台形分布を使うには平均が 1/3 から 2/3 に含まれていないといけないので、そうなるように s(t) の a, b, c を定めた。test_cm4_2.pl では、基礎寿命論が確証できるよう a, b, c を定めようとした。

しかし、このようなモデルにしてみてもうまくいっているように見えない。過剰診断論の設定、つまり、治療群と非治療群で実寿命に差がない設定をしてみると、死亡時がんステージの違いの分だけ、平均がん期間を長くして、平均罹患年齢を治療群と非治療群でそろえる必要が出てきて、その違いはわりとはっきり出るのだが、基礎寿命を見てみると、大きな違いが出てこなかった。

過剰診断論の設定の場合、治療群でも実寿命の延びはない。しかし、基礎寿命論から来る死亡時ステージの変化はある。しかし、基礎寿命の違いはなかなか観測されない。それは、実寿命の変化が結局のところないから。がん寿命による死亡を決めて、それ以下の死亡ステージのものを集めても、基本、実寿命は死亡時ステージと独立なためその群では実寿命で決まる「基礎寿命」に変化はないからだと考えられる。がん死亡率が、治療群では小さくなり、大きい基礎寿命になりがちな、がん死亡群が少なくなるため、基礎寿命が短く見える効果はあるが、それほど大きな効果ではない。

平均がん期間を延ばすと、平均死亡時ステージを少なめにする効果があり、さらに平均がん期間を延ばす必要があるため、数字が大きくなりやすく、平均がん期間推定は数字として拾いやすいようだ。

シミュレーションにおける、がん期間延長の原因は、基礎寿命減 + 治療による寿命増が考えられる。過剰診断論の設定のもとでも、治療効果はあると考えるべきで、だからこそ推定がん期間が伸びていると考えられる。もしこれが現実で観測されているなら、がん期間増の分(のある関数分)だけ、実寿命が延びるはずなのに、実寿命が延びていないということだから、その分基礎寿命減があるはずとなるだろう。

治療の侵襲による基礎寿命減がないとした場合の寿命を、平均死亡時がんステージの治療群と非治療群の差に、がん期間をかけたものを基礎寿命に足すことで、推定してみたりもした。が、本当にこの推定で正しいかは自信がない。


シミュレーションで使った名前付けについて


シミュレーションのプログラムでは、変数名や関数名に「英語っぽい」名前を付けた。私は疫学の知識がないため、実際の英語表記をどうすればよいかわからなかったのを適当に付けたため、滑稽に感じる人もいるかもしれない。謝罪しておく。ただ、プログラムでは、短い名前のほうがわかりやすいということもあるので、ある程度、許してもらえればありがたい。

cumulative mortality
死亡確率の累積分布。

life rest
余命。本当は life expectancy というらしい。

base life
基礎寿命。basic lifespan のほうがよかったかも。

cancer life
がん寿命。cancer lifespan のほうがよかったかも。

cancer term
がん期間。cancer span のほうがよかったかも。

of cure
治療群のデータ。あまりこういう言い方は見ないけど。

of disease
罹患非治療群のデータ。

cancer stage at death
死亡時がんステージ。

pre base life
治療の侵襲による基礎寿命がないとした場合の寿命。

over diagnosis
過剰診断。これは合ってるはず。



結果と今後の課題


過剰診断論の設定、つまり、治療群と非治療群で実寿命に差がない設定においても、治療の侵襲による基礎寿命減があるために、実際には、がん期間(または闘病期間)増に相当する寿命増が隠れている場合もありうることをシミュレーションを通じて示した。これがこの中間報告までの成果といえる唯一のことだろう。

しかし、がん期間などは現実に観測するのは難しく、がん期間や死亡時がんステージは、実寿命の関数になっている可能性もあるが、それとは独立に乱数で発生させるモデルを使っていて、現実的意味を見出しがたくなっている。このあやふやな根拠をもとに、シミュレーションをこれから続けていっても、どれほど意味があるのだろうかと私は悩むようになった。医学・疫学の専門家でも、シミュレーションの専門家でもない私には、知識がどうしようもなく不足している。

今回は、寿命の分布をもとにしただけだったが、この先は、人口推計をもとに、各年齢のエージェントを発生させて、そこで「統計」をとってみることが考えられる。そこで、総務省統計局の人口推計年齢男女別人口の Excel データから人口の分布を(100歳以上のデータを補完しながら) Perl に使いやすい形式で取り出す make_population_table.pl と、その分布にしたがって、ランダムな年齢のエージェントを発生させる test_pop.pl を書いたが、それ以上のことはしていない。

「統計」をシミュレートするためには、現時点の年齢を問題にするだけでなく、そこから何十年かたって寿命(余命)が尽きるエージェントが出てくる必要がある。ある年齢以上の統計という形にもできるが、通常は、その何十年かの間に出生があると考えるべきで、その出生をどう扱うかが問題となる。最初にエージェントを発生させたあと、一定数毎年エージェントを新たに足していくか、マイナスの年齢のエージェントを作っておくか、いろいろ方法はある。

そうやって発生させたエージェントに属性を割り当てて、統計をとるわけであるが、先ほど書いたように、あやふやな根拠でそれにどれだけ意味があるかわからない上に、統計学の多変量解析は私が不得意だった学問で、どうにも行き詰まりを感じている。最初に書いた疑問についてクリアにするとはどういうことなのか、今のところ判然としない。現代の統計においては計画が大事で、とにかく診断数が増えることに意義を見出すようなものではないことも、このシミュレーションを作る際の学習で思い出したということもある。

今回「中間報告」と書いたが、この先、続けるかは未決定である。少し頭を冷やして、何かいいアイデアが浮かばないか探ってみたい。


ターニングポイント


2017-03-07
NATROM 氏の記事に刺激され「過剰診断論」に興味を持つ。

2017-03-19
には make_cumulative_mortality.pl を作っていた。

2017-03-29
ドーキンス『利己的な遺伝子』を読了。

2017-04-12
make_population_table.pl と test_pop.pl を作る。

2017-04-18
test_cm4_2.pl で、だいたいの発生モデルはできたと考える。

2017-04-20
中間報告の公表。



参考文献


主なところのみを示す。

検診で発見されたがんの予後が良くても、検診が有効だとは言えない - NATROMの日記》。今回、過剰診断論を私が考えるキッカケとなった記事。引用もここから。原発震災から 6 年目の日の近くに公表されたため、私には刺激的に映った。

「過剰診断」とは何か - NATROMの日記》。以前に、NATROM さんが、「過剰診断」について詳しく説明した記事。

利己的な遺伝子 - 増補改題「生物=生存機械論」』(リチャード・ドーキンス 著, 日高 敏隆 & 岸 由二 & 羽田 節子 & 垂水 雄二 訳, 紀伊國屋書店, 1991年)。NATROM さんはドーキンスを尊敬しているらしいが、私は、ドーキンスはセンセーショナリズムの罪があると思う。それについては [cocolog:87123943] に書いた。

計算機シミュレーションのための確率分布乱数生成法』(四辻 哲章 著, プレアデス出版, 2010年 2013年第2版)。シミュレーションを書くとき頼った本。左右台形分布やベータ分布のアルゴリズムはこの本から引用した。ちなみに左右台形について [cocolog:87218726] でひとことしている。

gnuplot でヒストグラム(頻度分布図)を描画する - Qiita》。 gnuplot で分布を描くとき参照したが、ここの方法は、0点で +0 と -0 が集計されて突出してしまう「バグ」があるようだ。他を読むと (filter でなく bin という名前だが) filter(x,y)=floor(x/y)*y + y/2.0 とすれば良いようだ。



配布物


上で挙げたもののソースは以下にある。上で挙げた以外のソースも少し入っている。

更新: 2017-04-20
初公開: 2017年04月20日 16:22:41
最新版: 2017年04月20日 16:29:05

2017-04-20 16:22:40 (JST) in Perl, シミュレーション | | コメント (1) | トラックバック (0)

Perl」カテゴリ内の最近の記事

シミュレーション」カテゴリ内の最近の記事

批評や挨拶のためのネットコミュニティ

  • はてなブックマーク(って何?) このエントリーをはてなブックマークに追加 このエントリーを含むはてなブックマーク このエントリーを含むはてなブックマーク
  • Twitter (って何?)

トラックバック

この記事のトラックバックURL:
http://app.cocolog-nifty.com/t/trackback/93568/65174802

トラックバックのポリシー

他サイトなどからこの記事に自薦された関連記事(トラックバック)はまだありません。
» JRF のソフトウェア Tips:過剰診断論に関するシミュレーション 中間報告 (この記事)

コメント

tennis 初公開: over_diagnosis_simulation-20170420.zip。バージョン 0.01。

投稿: JRF | 2017-04-20 16:32:06 (JST)

コメントを書く



(メールアドレス形式)


※匿名投稿を許可しています。ゆるめのコメント管理のポリシーを持っています。この記事にまったく関係のないコメントはこのリンク先で受け付けています。
※暗号化パスワードを設定すれば、後に「削除」、すなわち JavaScript で非表示に設定できます。暗号解読者を気にしないならメールアドレスでもかまいません。この設定は平文のメールで管理者に届きます。
※コメントを書くために<b>ボールド</b>、<pre>詩文やソースコード</pre>、<a href="">リンク</a>、その他のHTMLタグが使えます。また、漢字[かんじ]でルビが、[google: キーワード] で検索指定が使えます。


ランダムことわざ: 温故知新。