2013年5月17日金曜日

福島県健康調査 甲状腺調査 について 市町村レベルでの分析
(2013/10/25若干改訂→受診者数をウエイトとする必要はないので、それについての記述、結果を削除。またt検定の方が適切だが手抜きしてz検定しているので、若干*の数が減ることにも注意)
 
0.経緯
福島の甲状腺調査は、当然ながら個人別のデータであり、個人レベルでサイズや結節有無などと線量の関係を分析できるはずだが、されていない。環境省のデータも公開されたが、

 環境省 報道発表資料-平成25年3月29日-福島県外3県における甲状腺有所見率調査結果について(お知らせ)

 調査対象者層が異なる(0-2才が含まれていない、福島とマッチングさせるための変数も充分にとられていない。→そもそもそのような分析自体想定していない)。など、まあ調査設計がひどい=福島と比較しづらい=感じ。

 情報公開クリアリングハウスさんのおかげで、H23年(被爆量の多い地域)の市町村別データが公開された。
 → 福島県民健康管理調査 23年度実績(市町村別甲状腺検査の結果を含む)

 これらを使って出来る範囲のことを分析した。

  上記データが公開された2013/4月末には、だいたい分析できたので、福島県(医科大)の健康管理センター、環境省の担当に電話で調査方法などを確認。→togetter  福島県健康調査(基本、甲状腺)、環境省甲状腺調査について

 しかし、サンプル数が少なく、不安定な部分もある云々で公開はしてこなかった。が、昨日からのstudy2007さんの分析や、それを巡るやりとり→ 
 H23年度福島県甲状腺検査市町村別結節割合と土壌汚染の関係について - Togetter
 を拝見して、避難関係の変数を追加して分析した結果を公開する。

 結論としては、H23調査で市町村別データが公開された13市町村について、結節サイズ5.1mm以上の人数(B判定)を分析したところ、WHOの推定(1年間)甲状腺線量と正の相関があった。ただし、浪江を除くと有意ではなくなるなど、注意が必要である。
 あくまで市町村レベルでの、各種データ間での(因果ではなく)相関であることに注意されたい。

 そのうち論文として投稿するかも知れないが、ここでは速報として引用なども簡略化する。


1 市町村別のデータ
下表にあるように、二つのデータを中心に市町村別データを整理。

 ・福島県民健康管理調査 23年度実績(市町村別甲状腺検査の結果を含む)
 ・環境省 報道発表資料-平成25年3月29日-福島県外3県における甲状腺有所見率調査結果について(お知らせ)

 線量についてはWHO,放医研、福島県WBCなどを利用。詳細は下の方参照。この他、国勢調査のデータもサンプリングバイアス検討のために利用(表には示してないが)。

表 市町村別データ(生データ)
注)黄色の3町はWHOは線量推定していないため、これらの地域で最も低い値を記入してある。出所は文中に大体示した。

表 市町村別データ(割合)


2 線量について
 線量について市町村別で推定、測定、公表されているものは、WHO,福島県基本調査外部線量、福島県WBC、放医研(内部被曝)、などがある。それらの中で、WHO(発生後1年間の甲状腺線量10才)が最も包括的なはずなので、それを用いる(外部、食物、呼吸を考慮)。
 → Health risk assessment from the nuclear accident after the 2011 Great East Japan earthquake and tsunami, based on a preliminary dose estimation
 ただし20km圏内の3町については推定されていないので、次の2シナリオを設定。
・同地域のなかでもっとも汚染度の低い15mSvとする。
・次の手順で、他の線量を用いて内挿して推定する。
 WHOの推定値がある10市町村について、線量を被説明変数とし、放医研、福島県基本調査外部線量、福島県WBCの各線量で重回帰分析。それをステップワイズで変数絞り込みした結果、放医研の「甲状腺線量1才児」のみがp=0.14となった(下表)。
 以下、Rからの出力をそのまま掲載する。通常のソフトと有意水準の表示が異なることに留意されたい。

表 線量内挿のための分析(従属変数:WHO線量、20km内3町を除く)

                   Estimate Std. Error t value Pr(>|t|)
(Intercept)           7.945     16.485    0.48     0.64
Thyroid.dose1yNIRS    1.588      0.969    1.64     0.14

Residual standard error: 22.1 on 8 degrees of freedom
  (3 observations deleted due to missingness)
Multiple R-squared:  0.251,Adjusted R-squared:  0.158 
F-statistic: 2.68 on 1 and 8 DF,  p-value: 0.14

 10%水準でも有意ではないし、R2も低いのだが、仕方ないのでこれを用いる。このパラメータと放医研の線量を用いて、3町の内挿値を求めた。

    放医研  内挿値
富岡町         10   23.82  
大熊町        20   39.70
双葉町         30   55.58

3 受診(回答)バイアス
表にあるように市町村によって受診率が異なる。被曝量が多い地域ほど心配して受診率が高くなる、といった傾向がある場合、それを考慮した分析を行わなければ、線量の推定にバイアスが生じる。
 これをチェックするため13市町村の受診率を従属変数、WHO甲状腺線量、国勢調査による総人口、男女比、平均年齢、15才未満人口の割合、人口密度を説明変数とした(割合なので)ロジスティック回帰分析を行った(検査数を重みとした)。さらに、ステップワイズ法で変数を絞り込んだ。
 その結果、人口密度、15才未満人口の割合などが有意となった。線量が有意とはならなかったことは、受診率の違いは、汚染度の高さ(への心配)ではなく、これらのデモグラフィクス要因によることを意味する。なお、線量については3町を15mSvとしても、内挿値を用いても有意とはならなかった
 以下では受診率バイアスを補正するために、傾向スコアを用いる。この手法では、R2が高いことが望ましい。このため、さらに、有意となった変数間の交互作用も導入した分析を行った。その結果、修正R20.83となった(下表)。
交互作用などがあるせいで、解釈は困難であるが人口密度の高い市部で受診率が高くなるようである。

   表 受診率のロジスティック回帰の結果
Coefficients:                                     
Estimate Std. Error t value Pr(>|t|)   
(Intercept)                          1.46e+02   4.00e+01    3.65   0.0082 **
人口密度                             1.13e-01   4.02e-02    2.80   0.0267 * 
人口性比                            -1.43e+00   4.10e-01   -3.50   0.0100 * 
年齢.総数.15歳未満人口割合          -1.16e+01   3.31e+00   -3.51   0.0099 **
人口密度:人口性比                   -1.16e-03   4.23e-04   -2.74   0.0288 * 
人口性比:年齢.総数.15歳未満人口割合  1.16e-01   3.40e-02    3.40   0.0115 * 
---Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 14.7 on 7 degrees of freedom

Multiple R-squared:  0.901, Adjusted R-squared:  0.83 F-statistic: 12.7 on 5 and 7 DF,  p-value: 0.00212
注)検査数を重みとした推定結果。

4 結節との関係
1)データ
 5.1mm以上の結節もしくは、20.1mm以上の嚢胞があるとB判定とされ2次検査が推奨される。後者についてはH23は1名しかおられないので、前者に注目する。
 市町村別の5.1mm以上の結節と線量の関係を示す。まず直感的にわかりやすい、受診者に占める、5.1mm以上の結節がみつかった者の割合を示す。円の大きさは検査数に比例している。

図表  WHO推定線量と5.1mm以上結節の割合
注)線量について3町は15mSvとしたもの。円の大きさは検査数に比例。

 割合を従属変数とした分析もできるが、人数というカウントデータなので、ポアソン回帰分析の方が適切である。よって、市町村別の線量と5.1mm以上の結節がみつかった「人数」をプロットしてみた。人口および、受診者数が異なるために、この図を読むときには注意が必要である。
 なお、赤マークは後述するモデルで推定した結果である。観測値とよくあてはまっていることがわかる。

図表  WHO推定線量と5.1mm以上結節の人数(および推定結果)
注)線量について3町は内挿して推定したもの。円の大きさは検査数に比例。
赤はモデルでの推定値を用いて計算した予測値(内挿値)。

2)分析方法
 5.1mm以上の結節人数を被説明変数、下記を説明変数としてポアソン回帰分析を行う。なお、甲状腺ガンは女性に多いが、性別内訳については市町村別に公開されていないので用いることができない。

・log(1+線量)
 ポアソン回帰では、発生頻度パラメータをλ=exp(λ0+Σβx)のように定式化するので、logをとっておくと、この係数が1のときにはλが線量とともに線形に増加、1より小さければ逓減型、大きければ逓増型となる。LSSなどでもポアソン回帰の枠組みで推定しているが、ERR、EARで定式化し、かつ線形を仮定している。ここでは普通の統計パッケージでも推定でき、上述のように線形を前提としない推定を行った。
 参考→ LSS(Life Span Study) 13報 被爆者データのRによる分析

・受診者に占める6-10才の割合
  はじめの表にあるように0-5才、11-15才、16-18才の割合もあるが、被説明変数との単相関がもっとも高いこれのみを用いた。

 上述のstudy2007さんの分析を参考として、2変数も作成した。
・避難率
 3/12日23時59分ごろまでに避難した者の割合(国会事故調報告書住民アンケートの結果p.362第4部 被害状況と被害拡大の要因 (その1) | 国会事故調」よりデジタイザで読み取り)

・ヨウ素剤を配布したか?ダミー
 国会事故調報告書「4.4.2 防護策として機能しなかった安定ヨウ素剤」によると楢葉町、富岡町、双葉町、三春町で配布。今回のデータには三春町は含まれていない。
  第4部 被害状況と被害拡大の要因(その3) | 国会事故調


 さらに線量については、次の二通りを推定した。
・3町を15mSvとした場合
・内挿して推定した場合

3)結果
 下記には、上述の4つの説明変数すべてを導入し、3町については線量を内挿して推定した場合の推定結果を示す(ウエイトあり、なし)。
 線量のパラメータは正で1%水準で有意である。なお、3町の線量を15mSvとした場合でも有意水準には変化がなかった。
 
表 ポアソン回帰の結果 (ウエイトなし,
3町線量内挿)   
                        Estimate Std. Error z value Pr(>|z|)    
(Intercept)            -1.35e+01   3.10e+00   -4.36  1.3e-05 ***
log(1 + Thyroid.dose4)  5.57e-01   2.02e-01    2.76   0.0058 ** 
年齢別6.10p             2.25e+01   9.38e+00    2.40   0.0163 *  
避難率pevacuate          8.38e-04   2.24e-03    0.37   0.7079    
ヨウ素剤配布ダミーfg.iodine -3.51e-01   3.56e-01   -0.99   0.3242    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for poisson family taken to be 1)
    Null deviance: 18.1174  on 12  degrees of freedom
Residual deviance:  6.6378  on  8  degrees of freedom
AIC: 64.49


4)頑健性の確認
 サンプル数は13市町村と少ないので、結果の頑健性を確認するために、1市町村づつ取り除いて12サンプルとして推定した(1つ取り除いた13通りを推定)。
 すると、どれか一つでも市町村をとり除くと避難率、ヨウ素剤配布ダミーは有意ではなくなった。不安定なので、これは、これ以降は入れないことにする。
 線量のパラメータについては、浪江町を取り除いた場合のみ、有意ではなくなる。

 サンプル数が少ないこともあり、若干の不安定さがあることには注意が必要である。

5)他の線量
 ここまでは、WHOの甲状腺線量10才(1年間)を用いた。WHOは1才、成人についても推定しているが、これらの相関はいずれも0.98と極めて高いため、結果は変わらない(はずである)。
 この他にも市町村別の線量を評価したものがあるので、比較のため、線量と6-10才の割合のみを説明変数として推定した(下表)。
 いずれも10%水準でも有意とはならないがNIRSについては10.5%であり、サンプル数が増加すると10%有意になる可能性はある。なお、WHOの推計は外部、食物、呼気からの線量をあわせて推計してあるので、これらを包含したものとなっているはずである。

表 他の線量を用いた推定結果(線量パラメータのみを抜き出したもの。
・放医研の1才児甲状腺線量
 データは→ 資料2 「甲状腺検査」の実施状況及び検査結果等について の24枚目(甲状腺線量の90パーセンタイル値:1才児)。もともと下記のWBCに基づいて算出されたもの。
                 Estimate Std. Error z value Pr(>|z|)
log(1 + Thyroid.dose1yNIRS)    0.322      0.198    1.62    0.105    

・福島県基本調査 外部線量
 データは→資料1 県民健康管理調査「基本調査」の実施状況について p.8から市町村別に平均値算出。
log(1 + ExtDose)    0.271      0.452    0.60  0.54974    

・福島県WBC測定値で預託線量1mSv以上と判定されたものの割合
 データは→平成23年6月27日~平成25年 3月31日まで の市町村別分布より算出。

log(1 + WBC1mp)   116.94      89.08    1.31     0.19    


6)環境省の3市調査
 環境省調査は0-2才を含んでいないため、福島調査とは比較が困難である。よって、あくまで参考値であるが、これら3市の線量を0として、年齢別分布も0-2才は含まれないことを無視して6-10才の割合を算出し、同様の分析を行った。
 ただし、福島甲状腺調査とは実施時期などが異なるため、これらの変数では表現できないなんらかの差異があると考え、コントロールダミー(fg.control :3市について1、福島県内は0)を導入した。
 コントロール変数をいれないモデル、入れたモデルのAICを比較すると、後者の方があてはまりは良好となった。つまり、これら3市とH23調査には体系的な何らかの差があるといえる。

 (ウエイトなしとした場合の)推定結果を示す。6~10才の割合の符号が、上の場合と逆転した。これは0-2才を含んでいない3市について無理矢理変数を導入したためと考えられる。ただし、線量パラメータは、いずれの場合も正で有意である。
 上記と同様、頑健性の検討を行った。上と同様、浪江を除いた場合のみ線量のパラメータは有意とならなくなった。6-10才の割合は抜け方によっては、有意にならない場合もあった。これは0-2才が含まれていないことを無視したためであろう。

表 環境省調査3市も含む(ウエイトなし、線量:3町は15mSv、3市は0mSv)
                      Estimate Std. Error z value Pr(>|z|)    
(Intercept)             -6.324      1.554   -4.07  4.7e-05 ***
log(1 + Thyroid.dose)    0.361      0.150    2.41  0.01599 *  
年齢別6.10p             -0.668      4.766   -0.14  0.88858    
fg.control               1.921      0.509    3.77  0.00016 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for poisson family taken to be 1)
    Null deviance: 38.688  on 15  degrees of freedom
Residual deviance: 16.324  on 12  degrees of freedom
AIC: 85.56



5 まとめ
 本稿では情報公開クリアリングハウスさんの情報公開請求によって"はじめて"公開された、H23年の市町村別データ等を用いて分析を行った。サンプル数は13もしくは16市町村と小さいものの、全データを用いると、甲状腺線量(10才)と5.1mm以上の結節人数には、有意な"相関"がみられた。
 なお、浪江を除くと、有意ではなくなるが、同町は避難指示があいまいであったため、何度も避難し、かつ高線量地域に避難させられた割合も高くなっている(国会事故調報告書 第4部 被害状況と被害拡大の要因 (その1) | 国会事故調)。線量が高いと推定することには一定の合理性があると考えられる。
 参考までに、WHOの推定は当地で収穫された作物を継続して摂取するなど、過大に(安全側)に見積もられているらしいが、例えば全地域の被曝量を10分の1にしても、パラメータの有意水準は変化しないことに注意したい。

 この検査はH23年10月からH24年3月までに行われたので、福島原発事故から7ー12ヶ月しか経過していない。「チェルノブイリでは事故後4年から小児の甲状腺ガンが増加した」という報告からみると、早過ぎる感もある。ただし、ここでの従属変数はガンではなく「5.1mm以上の結節」である。医学的な観点からどうなのかは、筆者にはわかりかねるが、当時からみるとスクリーニング機器の性能は向上したのではないだろうか。前駆段階としての結節の発生を捉えたのかもしれない。いずれにしてもより精密な線量評価、分析が必要である。

 なお、結節のサイズヒストグラムをみると、(より線量が低いはずの)H24の方がサイズも大きく、H23よりもピークが右にずれているようにもみえる。それが有意な差なのかも、H23,24調査の対象者を併せて分析すればある程度のことがいえる可能性はある。

図表 結節のサイズ分布
H23 調査

出所)H23  (平成24年9月11日開催) 資料2 「甲状腺検査」の実施状況及び検査結果について →p.17 (3.8万人)
                H24 調査


 はじめにも述べたが、1時点でのデータを用いた分析であり、因果ではなく相関であることにくれぐれも注意されたい。

 最後に要望を述べる。

・予断をもたない検査、調査、診断、分析を
 「100mSv以下では、影響がみられない/影響はない?」という予断をもって調査が行われている感がある。そのような予断をもたず検査、調査、診断、分析をすべきである。
 なお、被曝者を対象としたLSSで100mSv以下に限定するとパラメータが有意ではなくなる、ということについては、下記にまとめたので参照されたい。
 → LSS(Life Span Study) 13報 被爆者データのRによる分析

・素早い分析を
 いずれにしても、このような分析は、ちょっと学んだことがあれば、すぐにでもできるはずである。H24についても市町村別の集計はされているだろう(電話で要望したが公開予定はないとのことであった)。公開が望まれる(自分で公開請求すればいいだけの話ではあるが)。
 ここでは公開された市町村レベルでの分析しかできなかったが、調査側には個人レベルでのデータもあるはずである。それを用いた分析を行うべきである。
 個人の線量と甲状腺調査は未リンクだというが、甲状腺調査の実施対象者は、H23年度は4万人である。多忙であることは理解するが、それらの個人線量の推定とリンクを優先すべきである。
 
・情報の公開
 繰り返すが、この情報は情報公開請求があってはじめて公開された。被災された方々への対応を最優先にすべきであることはいうまでもないが、市町村別レベル程度の情報は即座に公開すべき内容である。
 プライバシー云々を強調しているが、匿名化したデータは医学分野でも多く公開されている。例えば最近公開されたRad Res誌掲載の長崎大によるWBC測定では、個人レベルのスペクトル、避難経路が掲載されている。この他、米国も原子力関係従業者の匿名化個人データを公開している。広島大学では被爆者について公開している。→放射線関連データ源 - Togetter 参照)。
 個人のプライバシー保護は当然であるが、過度の非公開は隠蔽というイメージを与える。

・適切な調査、分析方法の検討
 togetter  福島県健康調査(基本、甲状腺)、環境省甲状腺調査について でまとめたように、福島県健康調査調査センターには未だ疫学担当の方はおられないようである。また、環境省の調査は、福島調査とは、まったく別々に行われている。このため、必要な変数が測定されておらず、年齢、性別、社会経済変数など、個人レベルでマッチングすることも困難な状況となっているようである。さらに、詳細な比較分析の予定もないという。適切な体制構築が必要である。



謝辞
市町村データを公開請求、公開された情報クリアリングハウスさん、@study2007さん、および H23年度福島県甲状腺検査市町村別結節割合と土壌汚染の関係について - Togetter でやりとりされている方々から多くを学んだ。謝意を表する。