特集「がん統計データおよびその解析」［研究ノート］

統計数理（2011）
第 59 巻 第 2 号 217–237
c 2011 統計数理研究所


年齢 時代平面上における癌死亡リスクの視覚化
1

2

加茂 憲一 ・冨田 哲治 ・佐藤 健一

2（受付 2010 年 12 月 28 日；改訂 2011 年 7 月 19 日；採択 7 月 21 日）

要

旨

癌の動向に対する適切な統計的手法を用いた評価は，癌対策の根幹をなす重要な情報となる．
本論文では，癌死亡リスクを年齢と時代を座標とする空間内の曲面（リスク曲面）として視覚化することを目的とするが，実際には 2 次元の年齢 時代平面上において癌死亡リスクの高低を表現することになる．このように，ある平面上の関数曲面を推定する方法としては地理的加重回
帰，あるいは離散変数にも適用可能な地理的加重一般化線形モデルが広く利用されており，本目的に対しては，年齢と時代の組み合わせを仮想的な位置情報とみなすことで癌死亡リスク曲面の推定が可能となる．このモデルはノンパラメトリックであり自由度が極めて高い．これに対して，別のパラメトリックモデルも考察する．本論文では，年齢と時代に関する交互作用を用いたパラメトリックモデルによる癌死亡リスク曲面の推定も試みる．更に，癌の時系列解析
や予測において良く用いられる年齢・時代・コホートモデルに基づく解析結果を用いて癌死亡
リスク曲面を推定することも可能であり，その結果も併せて紹介し，各モデルより得られた結果を比較検討する．解析には統計フリーソフトウエア R を用い，解析スクリプトおよび Web
より入手可能なデンマーク肺癌死亡データを用いた解析例も併せて紹介する．本研究の目的である「視覚化」について達成されたかどうかに対する絶対的な評価は，現時点では不可能である．しかし，特性が既知であるデータに対して本手法を適用し，実際にその特性が視認できるかという点での検証は可能である．そこで，その検証の一例として「昭和一桁生まれ世代の高リスク」という出生コホート効果の示唆される日本の肝臓癌死亡データに適用し，その効果が視認できるかに着目した結果を示す．最後に，これらのモデルにおける数理的な面に着目し，その特性や問題点をまとめ，今後の改良点や応用について議論する．
キーワード： 癌，リスク曲面，視覚化，地理的加重一般化線形モデル，交互作用モデ
ル，年齢・時代・コホートモデル，コホート効果．
1.

はじめに

ある集団の癌の経年変動を解析し精確に評価することは，今後の癌対策の立案にも繋がる重要なテーマである．本論文では癌死亡や罹患に関するリスクを，年齢と時代を座標とする空間
内の曲面（以下，このような曲面を「リスク曲面」と表す）として視覚化する数理モデルの構築
およびその解釈を与えることを目的とする．しかしながら，この曲面は 3 次元空間内のものであるため表現が難しい．そこで，年齢 時代からなる 2 次元平面上に色の濃淡と等高線を用いて曲面を表現することによる視覚化を試みる．本研究においては，特に癌死亡に着目する．
1
2

札幌医科大学 医療人育成センター：〒060–8556 北海道札幌市中央区南 1 条西 17 丁目
広島大学 原爆放射線医科学研究所：〒734–8551 広島県広島市南区霞 1–2–3

218

統計数理

第 59 巻

第 2 号 2011

このような研究目的に対しては，地理的加重一般化線形モデル（Geographically weighted generalized linear model：以下「GW モデル」と表す）の適用が期待される．GW モデルは地理情報を有するデータの解析で多用されるが，本研究テーマにおいては年齢と時代の組み合わせを仮想的な位置情報とみなす．そして，カーネル平滑化の要領で，固定した年齢と時代の組み合
わせ周辺のデータに対する局所的な重み付きポアソン回帰を繰り返すことにより，癌死亡リスク曲面が推定される．
一方で，ノンパラメトリックな GW モデルに対して，年齢と時代の交互作用項を導入したパ
ラメトリックなモデル（以下「交互作用モデル」と表す）にも着目する．パラメトリックなモデ
ルは自由度が低いが，自由度の高い GW モデルによる結果をどの程度再現可能か，癌死亡リス
クの視覚化という点に注目し比較する．
更には，経時的に観測された癌データの解析や予測において広く用いられている年齢・時代・
コホートモデル（Age-period-cohort model：以下「APC モデル」と表す）についても，同様の研究目的に対して得られた結果を紹介する．
本論文は以下のように構成されている．まず第 2 章においては，本論文で取り扱う 3 種類のモデル，GW モデル，交互作用モデルおよび APC モデルの理論的な背景を紹介する．また，読者による解析結果の再現を可能にするために，統計フリーソフトウエア R のパッケージ Epiに実装されているデンマークにおける男性の肺癌死亡データ lungDK を用いた解析結果をスク
リプトと共に付加する．次に，第 3 章においては，実際に日本における癌死亡データに適用した結果を示す．本論文では，「昭和一桁生まれ世代の高リスク」という出生コホート効果の存在が示唆されている肝臓癌に着目し，その効果が上記モデルにより視覚化できるかという点に注目する．最後に，第 4 章においては，第 2，3 章での解析例の結果をもとに本論文で取り扱った 3 手法の比較検討を行い，今後の実用上の課題について議論する．
2.

年齢・時代平面上における癌死亡リスクの視覚化

本章では，癌死亡リスクを，年齢と時代を座標とする空間内における曲面として視覚化することを目的とした理論について説明する．癌死亡の発生は，全人口に対して稀な事象と考えられるので，人口（あるいは人年）をオフセットとするポアソン回帰モデルを適用する．具体的には，時代 p において年齢 a である人口（あるいは人年）を zap 人とし，そこから yap 人の癌死亡が発生したとすると，(2.1)

yap ∼ P oisson(zap λap )となる．ここで λap は時代 p において年齢 a である人の癌死亡リスクを表す未知パラメータである．従って，本研究においては，この未知パラメータ λap の年齢 時代に関する特性を視覚化する．
2.1 地理的加重一般化線形モデル
通常の回帰モデルにおいては未知パラメータを定数とするが，GW モデルでは，これを位置情報によって変化する関数として扱うのが特徴である（Brunsdon et al., 1998; Fotheringham etal., 2002; 谷村, 2010 など）．地理的に近いパラメータ間には関連性がある一方で，解析対象の全地域でパラメータが一様でないことも多いためである．
モデルに含まれる未知パラメータは，観測地点間の距離をもとに局所的な重み付けを行い推定される．この重みを定義する関数は距離低減関数と呼ばれ，推定値が周辺地域から受ける影年齢 時代平面上における癌死亡リスクの視覚化

219

響を距離が遠くなるに従って減少させるように設定される．例えばガウス型の距離低減関数は，
   2d
W (d) = exp −(2.2)
hで定義される（Brunsdon et al., 1998）．ここで d は地点間の距離を表し， h は関数形を決定するパラメータでありバンド幅と呼ばれる．バンド幅は，Stone（1974, 1977）
による交差検証法を用いて最適化される．すなわちn(yi − ŷi[−i] )2i=1を目的関数として，これを最小にするように求められる．ここで yi は地点 i における観測値，
ŷi[−i] は地点 i を取り除いたときの yi の推定量
（ジャックナイフ推定量）である．
GW モデルにおける位置情報として 2 次元平面上の座標 (x, y) を考えると，未知パラメータ関数は 3 次元空間内の曲面となる．経年的な癌死亡データは，死亡年別・年齢階級別という形式で構成されることが一般的であるので，GW モデルにおける (x, y) 座標に対応する仮想的な座標として（年齢，時代）を考察することにより，癌死亡リスク曲面の推定が可能となる．いま，（2.1）における未知パラメータを，年齢 a と時代 p の関数として(2.3)

log λap = β0 (a, p)と仮定する．この β0 (a, p) はノンパラメトリックに推定される．具体的にはカーネル平滑化の要領で，固定した年齢と時代の組み合わせ近傍データに対する局所的な重み付きポアソン回帰を繰り返し適用する．その際には，年齢・時代に関する近さを（2.2）における d とし，重み関数を（2.2）で定義する．
以上は，説明変数無しのモデルであるが，GW モデルにおいては説明変数をモデルに組み込
むことも可能である．その一例として，一般的な癌死亡データが性別となっていることに着目し，性別を変数として組み込むことを考えよう． m を性別を表すダミー変数（男性なら 1，女
性なら 0）として(2.4)

log λap = β0 (a, p) + β1 (a, p)mとモデルを構築する．ここで β1 (a, p) は性差の効果（ダミー変数は男性を 1 としているので，女性をコントロールとした際の男性の相対分）を表す未知パラメータ関数である．
2.2 年齢と時代に関する交互作用モデル
GW モデルにおいて，未知パラメータは（2.3）の通り陽に与えられないため，ノンパラメト
リックに推定することとなる．ノンパラメトリックモデルは関数形に制約が少ないため適用範
囲が広い一方で，パラメトリックモデルは適用範囲が狭い分様々な出力結果が得られる．そこ
で，前節の GW モデルによる結果をパラメトリックモデルでどの程度再現可能であるかに着目する．ここでは，佐藤 他（2009）や Satoh and Yanagihara（2008）による経時データに対する変化係数の推測法および，冨田 他（2010）による空間データに対する変化係数曲面の推測法を参考に，年齢と時代を座標とする空間内での癌死亡リスク曲面の推測法を考察する． GW モデル
における未知パラメータ（2.3）に対して，年齢と時代に関する q 個の基底関数を表す説明変数
x(a, p) = (x1 (a, p), . . ., xq (a, p)) を考え，次の線形構造を仮定したモデルを考える：
β0 (a, p) = θ x(a, p).

ただし θ = (θ1 , . . . , θq ) は q 次元未知パラメータベクトルである．基底関数 x(a, p) として
3 次元データの曲面補間で用いられる多項式曲面（Ripley, 1981）を利用する場合，1 次式なら

220

統計数理

第 59 巻

第 2 号 2011ば x(a, p) = (1, a, p, ap) と，2 次式ならば (1 + a + a2 )(1 + p + p2 ) の各項を基底とするので
x(a, p) = (1, a, p, a2 , p2 , ap, a2 p, ap2 , a2 p2 ) と設定する．1 次式の基底における ap あるいは 2 次式の基底における ap, a2 p, ap2 , a2 p2 の 4 項は，年齢と時代の交互作用を表す部分となっている．
また，GW モデルと同様に交互作用モデルにおいても，年齢と時代以外の変数も組み込むことが可能である．モデル（2.4）と同様に，性別を変数として組み込む場合は，性差に関する曲面
β1 (a, p) に関して線形構造
β1 (a, p) = φ x(a, p)
を仮定し，
(2.5)

log λap = β0 (a, p) + β1 (a, p)m = θ x(a, p) + φ x(a, p)mというモデルを構築すれば良い．
交互作用モデルはパラメトリックモデルであり，様々な仮説検定が可能であるという利点がある．例えば，性別を変数とするモデル（2.5）
において実際に性差が存在するかどうかを調べるためには，帰無仮説「性差を表す曲面＝ 0」に対する検定を行えば良い．つまり，帰無仮説 H0は(2.6)

for ∀a, p

H0 : β1 (a, p) = 0

である．冨田 他（2010）における変化係数曲面の推測法を適用すると，仮説（2.6）に関する検定
統計量は

T1 = φ̂ Ω̂

(2.7)

−1

φ̂,

Ω = cov(φ)となり，T1 の漸近分布が自由度 q のカイ 2 乗分布に従うことを用いる．ここで q は基底の数（例えば，2 次式の基底では 9 つ）である．
あるいは，モデル（2.5）において性差の一様性に着目する場合は，帰無仮説「性差を表す曲面＝定数」つまり(2.8)

H0 : β1 (a, p) = c = 0

for ∀a, pに対する検定を行えば良い．仮説（2.8）に対する検定統計量は
(2.9)



−1

T2 = φ̂[−1] Ω̂[−1] φ̂[−1] ,Ω[−1] = cov(φ[−1] )となり，T2 の漸近分布が自由度 q − 1 のカイ 2 乗分布に従うことを用いる．ここで φ̂[−1] は φ̂の第 1 成分を除いた q − 1 次元ベクトルである．
2.3 年齢・時代・コホートモデル
癌統計データ，特に経時測定データを解析するにあたって多用される手法に APC モデル（例
えば Holfold, 1983）がある．これは，特に癌の動向の将来予測において多用される．この将来予測に関しては，短期予測（例えば，Tiwari et al., 2004 や Pickle et al., 2007）
と長期予測
（例えば，
Kaneko et al., 2003 や大島 他, 2004）のどちらからも異なる観点からのニーズがあるが， APC
モデルは主に長期予測において用いられる．その理由は，コホート効果を明確にすることにより，同一コホートの時系列内での動きを的確に表現できるからである．ここでは，この APC
モデルを癌死亡リスクの視覚化に適用することについて考察する．
APC モデルとは，癌の死亡率を年齢，時代およびコホートの 3 成分で記述するものであり，（2.1）における未知パラメータ λap の対数値がlog λap = β + βa + βp + βc

年齢 時代平面上における癌死亡リスクの視覚化221で表現されているものである．ここで β は総平均効果，βa ，βp および βc はそれぞれ年齢 aの効果，時代 p の効果およびコホート c の効果を表す．3 つの変数 a, p および c には線形従属な関係 a + c = p が存在するため，前述の効果を一意に推定することができない．この問題はAPC モデルにおける識別問題と呼ばれており，様々なアプローチによる克服が試みられており（Holfold, 1983; 中村, 2000, 2005; 丹後, 2002 など），上記 3 成分が分解可能，あるいは推定可能な部分のみを用いた予測などが行われている．ここでは，3 成分の推定結果を再結合して得られる λap の推定値を用いて癌リスクの視覚化を試みる．
2.4 デンマーク肺癌死亡データを用いた解析例
前節で考察対象とした GW モデル，交互作用モデルおよび APC モデルは，統計フリーソフ
トウエア R に解析のための関数が用意されている．そこで，本章では読者の解析マニュアルとして役立つことを想定し，R のパッケージ Epi に実装されているデンマーク肺癌死亡データ解析結果を，スクリプトと併せて紹介する．このデータは lungDK というファイル名で実装されており 9 種類の項目から成るが，今回の解析で用いるのはこの中の 4 項目
（年齢，観測年，死亡数，人年）のみである．Web 上 http://staﬀ.pubhealth.ku.dk/ bxc/APC/MPIDR-2009/data/に
これら 4 項目のみのデータが公開されており，今回はそれを用いた．男性は lung5-M.txt，女
性は lung5-F.txt が対応するファイルである．以下，スクリプトの説明の簡便化のため，男性
データを変数名 Male，女性データを変数名 Female としておく．また Male と Female を合わせて，性別を表すダミー変数 male（男性＝ 1，女性＝ 0）を付加したデータセットを変数名 Whole
として準備しておく．ここで，Male および Female における項目 A，P，D および Y はそれ
ぞれ，年齢，観測年，死亡数および人年を表す．図 1 は，性別・年齢階級別の肺癌死亡数の経年変動を表すグラフであり，（a）
が男性を，
（b）
が女性を表す．これらを年齢と時代を座標とす
る平面上で表現すると（c）および（d）となる．ここで（c）は男性，（d）は女性であり，リスクの高
低を色の濃淡で表現している．色の濃い箇所が高リスクである．
2.4.1 地理的加重一般化線形モデルの適用
まず最初に GW モデルによる解析結果を紹介する． GW モデルの解析において最初に行うのは，最適なバンド幅の決定である．男性の肺癌死亡データ（Male）
を用いる場合，以下によりバンド幅が決定される．
library(spgwr); library(akima)ap <- cbind(Male$A,Male$P)
bw <- ggwr.sel(D~1+offset(log(Y)),coords=ap,Male,family=poisson)

まず，GW モデル解析において用いるパッケージ spgwr と，線形補完のためのパッケージ akimaを用意する．更に，位置情報（今回は年齢と時代）を変数名 ap として準備する．次に最適なバンド幅を関数 ggwr.sel により決定し，変数名 bw として保存する．いま，最適なバンド幅は，男性 0.883，女性 0.908 と推定された．ここで，距離低減関数としてはガウス型（2.2）が用いられている．このバンド幅の推定値（bw ）を用いて次のスクリプトを実行すると図 2 が描ける．ここで（a）は男性を，（b）は女性を表す．
gw <- ggwr(D~1+offset(log(Y)),coords=ap,Male,bandwidth=bw,family=poisson)
rr <- interp(Male$A,Male$P,exp(gw$SDF$X.Intercept.),linear=T)
image(rr)contour(rr,add=T)

関数 ggwr によって未知パラメータ λap および相対癌死亡リスク exp(λap ) を推定し，変数名

統計数理

222

図 1.

第 59 巻

第 2 号 2011

デンマークにおける肺癌の死亡数．デンマークにおける肺癌の性別・年齢階級（5 歳階級）別死亡数．（a）が男性，（b）が女性を表す．また（c），（d）はそれぞれ男性と女性について，癌死亡率を色の濃淡を用いて年齢 時代平面上に表現したものである．色の濃い箇所が高リスクである．

rr として格納する．次に，関数 image を用いて年齢 時代平面上に癌死亡リスク曲面の高低を濃淡で表現した後に，関数 contour により等高線を付加する．図 2 によると，男性に関しては（75 歳，1980 年）辺りが最大癌死亡リスクとなっており，等高線は楕円形の一部となった．一方で女性に関しては，グラフ内にピークが見当たらないが，おおよそ（65 歳，1990 年）に最大値が見受けられる．等高線の形状は男性の楕円形よりはむしろ円形に近い．
以上は男女別に解析を行った結果であるが，性別を変数とした包括的なモデル
（2.4）
を考察す
る場合，解析は以下のスクリプトで行う．
library(spgwr); library(akima)ap <- cbind(Whole$A,Whole$P)年齢 時代平面上における癌死亡リスクの視覚化

図 2.

デンマーク肺癌の GW モデル解析結果（性別）．デンマーク肺癌データに対する GW モ
デルによる推定結果．（a）と（b）はそれぞれ男性と女性の結果であり，癌死亡リスク曲
面を年齢 時代平面上に表現したものである．

図 3.

デンマーク肺癌の GW モデル解析結果（性差）．デンマーク肺癌データに対して，性別を変数として組み込んだ GW モデルによる癌死亡リスクの性差の推定結果．女性の癌死
亡リスクをコントロールとし，男性の相対分を年齢 時代平面上に表現したものである．

223

bw <- ggwr.sel(D~1+offset(log(Y))+male,coords=ap,Whole,family=poisson)
gw <- ggwr(D~1+offset(log(Y))+male,coords=ap,Whole,bandwidth=bw,family=poisson)
rr.d <- interp(Whole$A,Whole$P,exp(gw$SDF$male),duplicate="mean")
image(rr.d); contour(rr.d,add=T)

この場合，バンド幅の推定値は 0.914 であった．また exp(β1 (a, p)) は性差の相対癌死亡リスクを表す部分であり，それを年齢 時代平面上に表したものが図 3 である．ただし，ダミー変数
は男性を 1 としているので，図 3 は女性をコントロールとした際の男性の相対分を表す．この

統計数理

224

第 59 巻

第 2 号 2011

等高線についても，楕円形の一部に近い形であった．
2.4.2 年齢と時代に関する交互作用モデルの適用
次に，交互作用モデルを考える．基底関数として 3 次式をフルモデルとし，Akaike（1973）による赤池情報量規準（AIC）に基づくステップワイズ変数減少法により，最適な変数の組み合わせを探索した．それのためには，以下のスクリプトを実行すれば良い．
library(MASS)
res <- glm(D~offset(log(Y))+(A+I(A^2)+I(A^3))*(P+I(P^2)+I(P^3)),Male,family=poisson)stepAIC(res)
AIC に基づくステップワイズ関数 stepAIC はライブラリ MASS に用意されている．男性については a3 p の項を除いたモデル，女性についてはフルモデルが最適であった．この結果に基づき，癌死亡リスク曲面（女性）を描くには以下のスクリプトを実行すれば良い．
library(akima)res <- glm(D~offset(log(Y))+(A+I(A^2)+I(A^3))*(P+I(P^2)+I(P^3)),Female,family=poisson)basis <- function(a,p){c(1,a,a^2,a^3,p,p^2,p^3,a*p,a*p^2,a*p^3,a^2*p,a^2*p^2,a^2*p^3,a^3*p,a^3*p^2,a^3*p^3)}
theta <- res$coef
rr.int <- function(theta,data){aa <- seq(min(data$A),max(data$A),by=1)
表 1.

デンマーク肺癌の交互作用モデルにおける未知パラメータ推定結果．デンマーク肺癌
データに対して，年齢と時代に関する交互作用モデルを適用したときの，未知パラメー
タの推定結果．表中の 1 は定数項を，a は年齢の項を，p は時代の項を意味し，これらのクロス表で変数を表す．例えば，男性における変数 ap2 に対する未知パラメータの推定値は，2 行 3 列の −5.107 × 10−2 である．また，基底のフルモデルを 3 次関数に設定したモデル選択の結果，男性における a3 p の項は不要であったので，記号 “−” で表している．また，性差とは包括的なモデル（2.5）における男性の相対分に対する未知パラメータの推定値である．

年齢 時代平面上における癌死亡リスクの視覚化

図 4.

225

デンマーク肺癌の交互作用モデル解析結果（性別）．デンマーク肺癌データに対して，年
齢と時代に関する交互作用モデルによる推定結果．ただし基底は，男性は 3 次関数から
a3 p の項を除いたもの，女性は 3 次関数である．（a）と（b）はそれぞれ男性と女性の結果であり，癌死亡リスク曲面を年齢 時代平面上に表現したものである．

pp <- seq(min(data$P),max(data$P),by=1)
rr <- matrix(NA,length(aa),length(pp))
for(i in 1:length(aa))for(j in 1:length(pp))
rr[i,j] <- exp(theta %*% basis(aa[i],pp[j]))image(aa,pp,rr)
contour(aa,pp,rr,add=T) }rr.int(theta,data=Female)

変数名 rr.int は，データと未知パラメータの推定値を引数とする自作関数であり，癌死亡リ
スク曲面を表現するものである．未知パラメータの推定結果および図 2 に倣った年齢 時代平面上における癌死亡リスク曲面の表現をそれぞれ表 1 と図 4（（a）は男性，（b）は女性）に示す．
このモデルはパラメトリックなモデルであるので自由度は低いが，自由度の極めて高い GW モ
デルによる結果である図 2 の特性を端的に捉えている．
次に，性差を変数とした包括的なモデル（2.5）における β1 (a, p) を年齢 時代平面上に視覚化することを考える．関数 β1 (a, p) の推定結果を図 5 に示すが，これは以下のスクリプトを実行することにより描くことができる．
library(akima)res <- glm(D~offset(log(Y))+(A+I(A^2)+I(A^3))*(P+I(P^2)+I(P^3))*male,Whole,family=poisson)
basis <- function(a,p){c(1,a,a^2,a^3,p,p^2,p^3,a*p,a*p^2,a*p^3,a^2*p,a^2*p^2,a^2*p^3,a^3*p,a^3*p^2,a^3*p^3)}
theta.d <- res$coef[grep("male",names(res$coef))]
rr.int(theta.d,data=Whole)

統計数理

226

図 5.

第 59 巻

第 2 号 2011

デンマーク肺癌の交互作用モデル解析結果（性差）．デンマーク肺癌データに対して，性別を変数として組み込んだ交互作用モデルによる癌死亡リスクの性差に関する推定結果．
女性の癌死亡リスクをコントロールとし，男性の相対分を年齢 時代平面上に表現したものである．ただし，基底は 3 次関数である．

次に，2.3 節で触れた 2 種類の仮説検定（性差の存在性と一様性）を考える．これらの検定を行うには，以下のスクリプトを実行すれば良い．
A0 <- Whole$A-mean(Whole$A)
P0 <- Whole$P-mean(Whole$P)
res <- glm(D~offset(log(Y))+(A0+I(A0^2)+I(A0^3))*(P0+I(P0^2)+I(P0^3))*male,
Whole,family=poisson)
theta.d <- res$coef[ grep("male",names(res$coef))]omega.d <- vcov(res)[grep("male",names(res$coef)),grep("male",names(res$coef))]
t1 <- theta.d %*% solve(omega.d) %*% theta.d
p1 <- pchisq(t1,df=length(theta.d),lower.tail=F)
t2 <- theta.d[-1] %*% solve(omega.d[-1,-1]) %*% theta.d[-1]p2 <- pchisq(t2,df=length(theta.d)-1,lower.tail=F)

スクリプト中の変数名 t1 と t2 はそれぞれ帰無仮説（2.6）と（2.8）に対する検定統計量（2.7）と（2.9）
であり，変数名 p1 と p2 がそれぞれに対する p 値である．また年齢と時代変数について，
年齢については 40 ∼ 85，時代については 1943 ∼ 1993 と大きな数値であるため，平均値を減じる基準化を行った．このことにより，分散共分散行列 Ω および Ω[−1] が特異になることを防いでいる．
性差の存在性に関する検定について，検定統計量（2.7）
は 31383(p < 0.001) であったので，性差が存在することが分かる．また，性差の一様性に関する検定について，検定統計量は 5591(p < 0.001)
であったので，性差の影響は年齢や時代に対して一様でないことも分かる．
2.4.3 年齢・時代・コホートモデルの適用
最後に APC モデルを考える．R には APC モデルに基づく推定関数 apc.ﬁt がパッケージ Epi
において既に用意されており（Carstensen, 2007），これを用いて，例えば次のスクリプトにより解析が可能である．

年齢 時代平面上における癌死亡リスクの視覚化

227

library(Epi)apc <- apc.fit(A=A,P=P,D=D,Y=Y,Male,model="factor")
apc.plot(apc)

パッケージ Epi の使用を宣言し，パッケージに含まれる関数 apc.ﬁt による結果を変数名 apc
として格納する．そして関数 apc.plot により図 6 の
（a）
と
（b）
が描かれる．ここで（a）は男性を，（b）は女性を表す．これらは 3 本の折れ線により構成されるが，左より年齢効果，コホート効
果および時代効果を表すものである．これら 3 成分の推定結果を再結合することにより癌死亡リスク曲面が構築される．それが図 6 の
（c）
と
（d）
であり，これらは以下のスクリプトにより作製される．ここで（c）は男性を，（d）は女性を表す．
library(akima)rr.apc <- function(apc,data){a <- apc$Age; p <- apc$Per; c <- apc$Coh
aa <- unique(data$A); pp <- unique(data$P)
rr <- matrix(NA,length(aa),length(pp))
for(i in 1:length(aa))for(j in 1:length(pp))
rr[i,j] <- a[a[,1]==aa[i],2]*p[p[,1]==pp[j],2]*c[c[,1]==pp[j]-aa[i],2]
image(aa,pp,rr)
contour(aa,pp,rr,add=T) }rr.apc(apc,data=Male)

自作関数 rr.apc によって，APC モデルによる 3 成分の推定結果から各（年齢，時代）の癌死亡リスクを算出し，癌死亡リスク曲面を年齢 時代平面上に表現する．
2.4.4 解析結果のまとめ
最後に，これらの解析結果に対する考察を行う．性別に描いた図については，癌死亡リスクが最大となる位置が存在し（図中に必ず存在するとは限らないが）そこより連続的に減少してゆく傾向にあった．これを等高線で表現すると，楕円形ないし円形（の一部）に近い形状であった．
もし，年齢効果と時代効果共に連続で単峰の凸性が存在するならば，等高線は楕円形（あるいは円形）となる．その長軸は，年齢効果と時代効果における凸性の強弱に依存して，年齢軸か時
代軸のどちらかに並行となる．しかし男性（図 2，4 および 6 の（a））および性差（図 3 および 5）に関しては，長軸はおおむね傾き 1 の方向であった．この方向は，同じ出生年に属するコホートに対する軸を表すものであるため，今後「同一コホート方向」と呼ぶこととする．もし癌死亡リスクに影響を与えるのが出生コホート効果のみであるならば，等高線は同一コホート方向に表れることになる．図 2，4 および 6 の男性および図 3，5 における等高線は楕円形（の一部）に近い形状，かつその長軸が同一コホート方向であるため，同じ出生年のコホートは共通の特徴を持つという出生コホート効果の存在性が示唆される結果であった．
3.

日本における肝臓癌データを用いた解析

本章では，上記の手法を実際に日本のデータに適用した結果を紹介する．癌死亡および人口のデータは，人口動態統計（厚生労働省大臣官房統計情報部編）
より提供されており，国立がん研究センターのホームページ http://www.ncc.go.jp/jp/ よりダウンロードできるが，これは観測年毎の性別・年齢階級別死亡数という形である．一般的に癌死亡や罹患（Matsuda et al., 2011）に関するデータはこの形式であることが多い．

統計数理

228

図 6.

第 59 巻

第 2 号 2011

デンマーク肺癌の APC モデル解析結果．デンマークにおける肺癌データに対する APC
モデルによる推定結果．
（a）と（b）はソフトウエア R におけるパッケージ Epi の関数
apc.plot による推定結果で，（a）は男性，（b）は女性に関する結果である．3 本の折れ線で構成されているが，左より，年齢効果，生まれ年（コホート）効果，時代効果を表す．
また，（c）と（d）はそれぞれ，（a）と（b）の結果を用いて，癌死亡リスク曲面を年齢 時代
平面上に表現したものである．
（c）は男性，（d）は女性の結果である．

本研究の目的である癌死亡リスクの視覚化について，現時点では目視による確認を行うこととなる．そこで本章では，年齢と時代に関連する特徴を持つことが先験的に知られている臓器に着目し，その特徴が視認できるかどうかをチェックすることとする．そのための一例として，肝臓癌に着目しよう．肝臓癌は 2009 年における臓器別死亡数 32725 人であり，肺（67583人），胃（50017 人）に次いで第 3 位である主要な癌の 1 つである
（ただし，大腸（42800 人）＝結腸（28692 人）＋直腸（14108 人）は，この順位から除いている）．まず，肝臓癌の性別・年齢階級別（40 ∼ 84 歳）の死亡率（粗率）の経年変動を図 7 に示す．ただし，（a）は男性を，（b）は女性を表す．これらを年齢と時代を座標とする平面上で表現すると（c）および（d）となる．ここで（c）

年齢 時代平面上における癌死亡リスクの視覚化

図 7.

229

日本における肝臓癌の死亡率．日本における肝臓癌の性別・年齢階級（5 歳階級）別死亡率（粗率）．
（a）が男性，（b）が女性を表す．また（c），（d）はそれぞれ男性と女性について，癌死亡率を色の濃淡を用いて年齢 時代平面上に表現したものである．色の濃い箇所が高リスクである．

は男性，（d）は女性であり，リスクの高低を色の濃淡で表現している．色の濃い箇所が高リス
クである．
図 7 より，特に男性で顕著に 1935 年前後の出生コホート（昭和 1 桁生まれ世代）
における癌死亡リスクが他に比べて高いことが伺える．つまり，同一の出生コホートが共通の特徴を持つ，あるいは特定の出生コホートが他と異なる特徴を持つという「出生コホート効果」の存在が示唆される癌であることが分かる．この点については大島（2003）など多くの論文で触れられているように，日本における肝臓癌の多くは肝炎ウイルスの持続感染が原因であり，その抗体陽性率が特にこの出生コホートにおいて高いことがその理由として考えられる．このような肝臓癌
における出生コホート効果が，GW モデル，交互作用モデルおよび APC モデルにより視覚化可能かどうかは非常に興味深いテーマであろう．なお，Marugame and Sobue（2004）
や Imamura

統計数理

230

第 59 巻

第 2 号 2011

図 8.

日本肝臓癌の GW モデル解析結果．日本における肝臓癌データに対する GW モデルによる推定結果．（a）と（b）はそれぞれ男性と女性の結果であり，癌死亡リスク曲面を年
齢 時代平面上に表現したものである．全ての図は，デンマーク肺癌における図 2 に対
応する．

図 9.

日本肝臓癌の GW モデル解析結果（性差）．日本における肝臓癌データに対して，性別を変数として組み込んだ GW モデルによる癌死亡リスクの性差の推定結果．女性の癌死亡リスクをコントロールとし，男性の相対分を年齢 時代平面上に表現したものである．デンマーク肺癌における図 3 に対応する．

and Sobue（2004）
で指摘されているように，肝臓癌に限らず癌には出生コホート効果による差の存在が示唆されるものは少なくない．
解析に使用したデータは，前述の通り人口動態統計（厚生労働省大臣官房統計情報部編）により提供されている性別，年齢階級（5 歳階級）別，部位別の癌死亡および，オフセットとして人口について 1958 ∼ 2008 年のものである．これらは，国立がん研究センターホームページ
http://www.ncc.go.jp/jp/ よりダウンロードできる．ただし，本解析においては結果の安定性

年齢 時代平面上における癌死亡リスクの視覚化

図 10.

231

日本肝臓癌の交互作用モデル解析結果（性別）．日本における肝臓癌データに対して，
年齢と時代に関する交互作用モデルによる推定結果．ただし基底は，男性は 3 次関数から a3 p の項を除いたもの，女性は 3 次関数である．（a）と（b）はそれぞれ男性と女性の結果であり，癌死亡リスク曲面を年齢 時代平面上に表現したものである．全ての図は，デンマーク肺癌における図 4 に対応する．

を図るため，年齢については 40 歳 ∼ 84 歳に限定することにより，若年および老年の影響を排除した．
GW モデル（2.3）による性別の推定結果を図 8（ただし（a）は男性，（b）は女性）に示す．これは，前出の図 2 に対応する．まず，バンド幅の推定値は，男性で 0.263，女性で 1.116 であった．一方，性別を変数とした包括的なモデル（2.4）においてはバンド幅は 0.548 と推定された．
そして，図 3 に倣って性差（女性をコントロールとした際の男性の相対分）を年齢 時代平面上に表現したものが図 9 である．
次に，交互作用モデルを考察する．本データの解析においても，3 次式の基底をフルモデルとした．つまり推定する未知パラメータは，フルモデルにおいては (1 + a + a2 + a3 )(1 + p + p2 + p3 )の各項における係数であり，全部で 16 個である．AIC に基づく変数選択の結果，男性については a3 p の項を除いたモデル，女性についてはフルモデルが最適であった．未知パラメータの推定結果を表 2 に示し，前出の図 4 に倣って癌死亡リスク曲面を表現したものが図 10 である．
また，性差について図 5 に倣った結果が図 11 である．次に，性差が存在するかどうかを調べるための，帰無仮説（2.6）に対する検定を行うと，検定統計量（2.7）は 226050(p < 0.001) であったので，性差は存在することが分かる．次に性差の一様性については，帰無仮説
（2.8）
に対する
が 26527(p < 0.001) であったので，性差の影響は，年齢や時代に対して一様で
検定統計量
（2.9）
ないことも分かる．しかし，図 5 と図 11 を比較すると，性差の特徴はデンマークの肺癌と日本の肝臓癌で大きく異なった．デンマークの肺癌については，1900 年前後の出生コホートで極めて高い性差（男性の高リスク）が見られたような出生コホート効果の存在が特徴的であったが，日本の肝臓癌については，このような性差に関する出生コホート効果は見られなかった．
最後に APC モデルを適用した結果を図 12 に示す．これは，前出の図 6 に対応するものである．前の 2 つのモデルと同様な特性が表現できたが，年齢階級の影響でメッシュが粗くなって
いる．

統計数理

232

4.

第 59 巻

第 2 号 2011

図 11.

日本肝臓癌の交互作用モデル解析結果（性差）．日本における肝臓癌データに対して，性別を変数として組み込んだ年齢と時代に関する交互作用モデルによる癌死亡リスクの性差に関する推定結果．女性の癌死亡リスクをコントロールとし，男性の相対分を
年齢 時代平面上に表現したものである．ただし，基底は 3 次関数である．デンマーク
肺癌における図 5 に対応する．

表 2.

日本肝臓癌の交互作用モデルにおける未知パラメータ推定結果．日本における肝臓癌
データに対して，年齢と時代に関する交互作用モデルを適用したときの，未知パラメー
タの推定結果．表中の 1 は定数項を，a は年齢の項を，p は時代の項を意味し，これらのクロス表で変数を表す．例えば，男性における変数 ap2 に対する未知パラメータの推定値は，2 行 3 列の 1.458 × 10−1 である．また，基底のフルモデルを 3 次関数に設定したモデル選択の結果，男性における a3 p の項は不要であったので，記号 “−” で表している．また，性差とは包括的なモデル（2.5）における男性の相対分に対する未知パ
ラメータの推定値である．

おわりに

本論文では，癌死亡リスクを年齢と時代を座標とする空間内の曲面として視覚化するという目的に対する統計モデルおよび解析結果を紹介した．ここでは，GW モデル，交互作用モデル
および APC モデルの 3 つを用いて癌死亡リスク曲面を推定し，その結果を年齢 時代からなる

年齢 時代平面上における癌死亡リスクの視覚化

図 12.

233

日本肝臓癌の APC モデル解析結果．日本における肝臓癌データに対する APC モデルによる推定結果．
（a）と（b）はソフトウエア R におけるパッケージ Epi の関数 apc.plot
による推定結果で，（a）は男性，（b）は女性に関する結果である．3 本の折れ線で構成されているが，左より，年齢効果，生まれ年（コホート）効果，時代効果を表す．また，（c）と（d）はそれぞれ，（a）と（b）の結果を用いて，癌死亡リスク曲面を年齢 時代平面上に表現したものである．
（c）は男性，（d）は女性の結果である．全ての図は，デンマー
ク肺癌における図 6 に対応する．

2 次元平面上に濃淡と等高線を用いて視覚化した．
まず始めに GW モデルについて考察する．このモデルの概念は，カーネル平滑化に基づいて
おり，距離低減関数を用いて近隣情報に重みをつけた推定を行うため，自由度は非常に高い．
従って，リスクの集積性など，局所的な特性を表現することが可能となる．今回，年齢 時代を
仮想的な位置情報とみなすことにより癌死亡データへの適用が可能となった．また，年齢と時代以外の変数を組み込むことも可能であり，例えば性別を変数として組み込むと図 3 や 9 のよ
うに性別の効果（性差）についても，年齢 時代平面上に視覚化が可能である．

234

統計数理

第 59 巻

第 2 号 2011

次に，交互作用モデルについて考察する．このモデルは，汎用的な一般化線形モデルがベースとなっているのが利点である．実際 R においては，GW モデルについて推定関数 ggwr，APC
モデルについては推定関数 apc.ﬁt が用意されているが，現時点では他の多くのソフトウエアには実装されていない．しかし，一般化線形モデルについては殆どの統計ソフトウエアに実装されているので，交互作用モデルについては，多くの読者は手持ちのソフトウエアで実行可能であろう．また，今回用いた多項式基底では自由度が極めて低いものの，信頼区間の構築や検定が実行可能であるのも利点である．実際の解析結果（デンマークの肺癌と，日本の肝臓癌）を見ても，3 次式をフルモデルとする基底でほぼ特徴は再現できていたので，癌死亡リスクのおおまかな特徴を知るにはこのようなパラメトリックなモデルで十分であると考えられる．
最後に APC モデルについて考察する． APC モデルの最大の利点は，年齢・時代・コホート
効果の定量的な評価が可能な点にある．例えば，今回解析例として取り上げた日本の肝臓癌における昭和一桁生まれの高リスク出生コホート効果については，図 12（a）のように明確なピークとして捉えられる．しかし，本論文の目的である，年齢 時代平面上における癌死亡リスクの視覚化という点に関しては，説明変数が飽和状態となるほどの柔軟なモデルであるが，結果に関しては GW モデルほどの柔軟さは見られない．一方で，モデルのシンプルさという点に関しては交互作用モデルほどではないという位置付けになる．
次に，これらの手法を出生コホート効果の存在が示唆されている日本の肝臓癌に適用した結果について考察を行う．大島（2003）など多くの先行研究により，1935 年前後（昭和一桁世代）の出生コホートにおいて，他に比べて癌死亡リスクが高いことが指摘されており，この点を表現できるかが 1 つの鍵であった． GW モデルによる結果である図 8 においては，男性における1935 年前後の出生コホートにおける特異的な高リスク効果を視認することができた．それは，この出生コホートが辿る同一コホート方向において，あたかも山の尾根のように高リスクなラインが走っているからである．性差について年齢 時代平面上に表現した図 9 は，45 歳–1995 年
周りにピークを持ちそこから均等に低くなる傾向が見られた．また，交互作用モデルにおいても（3 次関数をフルモデルとして）1935 年出生コホート効果および，その他全体的な特性の再現が可能であった．ただし自由度が低いため，1935 年出生コホート効果は若干丸められて表現されている．
最後に APC モデルに関しては，同様の特性が観測されたが，年齢階級設定の影響によりメッシュが粗くなる．
以上，GW モデル，交互作用モデルおよび APC モデルについて，R のスクリプトを交えながら説明し，日本の肝臓癌データに適用した結果についても議論を行った．本研究の目的は，統計モデルによる癌死亡リスクの視覚化であり，癌死亡リスク曲面を年齢 時代平面上に表現することを試みた．実解析として，日本の肝臓癌における「昭和一桁世代の高リスク」というコホート効果の存在も，具体的な年の確認までには至らなかったが視覚化した．これらの結果を踏まえて，3 つのモデルを 2 つの視点
（結果の精密さ，モデルの汎用性）から比較する．最も精密かつ柔軟な結果（リスク曲面の推定結果）が得られたのは GW モデルであった．一方で，最も柔軟性を欠く結果は交互作用モデルによるものであったが，ソフトウエアの利用や検定の可能性等を鑑みると，最も汎用性の高いモデルであるといえる．
これらの推定結果および解釈はあくまで目視に依存するものであり，解析者の主観や先入観がある程度混入することが考えられる．そこで，主観に頼らず自動的に癌死亡リスクの特性を検出する手法の開発が今後のテーマの 1 つとなる．
また，今後の研究テーマとして癌の時系列解析で多用される「予測」に着目する．この点に関して，GW モデルにおいては外挿の問題が発生する．本来 GW モデルはホットスポットと呼ばれる癌死亡リスクの集積する地点の検出など，考察対象地域内における特性を検出するこ年齢 時代平面上における癌死亡リスクの視覚化

235とに優れており，予測を対象地域外に適用（外挿）することは基本的に想定されていない．また交互作用モデルにおいても，自己回帰モデルのように同一コホート内の経時的な相関を記述す
るなどの工夫が必要である．従って，これら 2 つのモデルについては，現時点では予測に適したモデルとは言い難い．一方で APC モデルについては，本来予測において多用される手法でもあり，予測方法はある程度確立されている．従って GW モデルと交互作用モデルに関して
は，予測の可能性についての議論およびモデルの改良が今後のテーマである．
本研究においては，癌死亡リスクの視覚化，特にそれを年齢と時代を座標とする空間内における曲面として表現する手法に関する考察を行った．地図を描く際に有効な GW モデル，汎用性の高い交互作用モデル，癌の時系列解析で頻出の APC モデルを用いた解析を 2 種類
（デン
マークの肺癌と日本の肝臓癌）行った．本論文では，モデルの特性や結果の解釈および今後の発展の可能性等を総合的に判断して，これらの優劣を完全に決する段階には至らなかった．そこで，現時点では各モデル固有の長所や短所を勘案しつつ，実際の解析においては，手法を限定せず複数を併用して評価を行うのが有効であると考えられる．
謝

辞

査読者より大変有益なコメントを頂きましたこと，感謝申し上げます．

参 考 文 献Akaike, H.（1973）. Information theory and an extension of the maximum likelihood principle, 2ndInternational Symposium on Information Theory（eds. B. N. Petrov and F. Csáki）, 267–281,
Akadémiai Kiadó, Budapest.
. Geographically weighted regression ModBrunsdon, C., Fotheringham, S. and Charlton M.（1998）（3）, 431–
elling spatial non-stationarity, Journal of the Royal Statistical Society: Series D, 47
443.
Carstensen, B.（2007）. Age-period-cohort models for the Lexis diagram, Statistics in Medicine, 26,
3018–3045.
. Geographically Weighted Regression:
Fotheringham, A. S., Brunsdon, C. and Charlton, M. E.（2002）
The Analysis of Spatially Varying Relationships, Wiley, Chichester.
. The estimation of age, period and cohort eﬀect for vital rates, Biometrics, 39,
Holfold, T. R.（1983）
311–324.
Imamura, Y. and Sobue, T.（2004）. Mortality trend of colon, rectal, liver, “gallbladder and biliary
tract” and pancreas cancer in Japan by birth cohort, Japanese Journal of Clinical Oncology,
34（8），491–493.
Kaneko, S., Ishikawa, K. B., Yoshimi, I., Marugame, T., Hamashima, C., Kamo, K., Mizuno, S. and. Projection of lung cancer mortality in Japan, Cancer Science, 94（10），
Sobue, T.（2003）
919–923.
国立がん研究センター ホームページ（2011）. がん情報サービス，http://www.ncc.go.jp/jp/.
Marugame, T. and Sobue, T.（2004）. Mortality trend of mouth and pharynx, esophagus, stomach,（7）,
larynx and lung cancer in Japan by cohort, Japanese Journal of Clinical Oncology, 34
432–438.
Matsuda, T., Marugame, T., Kamo, K., Katanoda, K., Ajiki, W., Sobue, T. and Japan Cancer Surveillance Research Group（2011）. Cancer incidence and incidence rates in Japan in 2005: Basedon data from 12 population-based cancer registries in the Monitoring of Cancer Incidence inJapan（MCIJ）
project, Japanese Journal of Clinical Oncology, 41（1），139–147.

236

統計数理

第 59 巻

第 2 号 2011中村 隆
（2000）. 質問項目のコウホート分析 多項ロジット・コウホートモデル ，統計数理，48（1）,
93–119.
中村 隆（2005）. コウホート分析における交互作用効果モデル再考，統計数理，53（1）, 103–132.
大島 明（2003）. 日本での慢性肝疾患，特に肝癌の疫学的特徴，第 123 回日本医学会シンポジウム記
録集 ウイルス肝炎，13–19.
大島 明，黒石哲生，田島和雄（2004）. 『がん・統計白書』，篠原出版新社，東京．
Pickle, L. W., Hao, Y., Jemal, A., Zou, Z., Tiwari, R. C., Ward, E., Hachey, M., Howe, H. L. and Feuel,. A new method of estimating United States and state-level cancer incidence
E. J.（2007）
counts for the current calendar year, CA A Cancer Journal for Clinicians, 57（1），30–42.
Ripley, B. D.（1981）. Spatial Statistics, Wiley, New York.
佐藤健一，柳原宏和，加茂憲一（2009）. 離散分布の経時測定データにおける線形な変化係数の推測について，応用統計学, 38, 1–11.
. Estimation of varying coeﬃcients for a growth curve model,
Satoh, K. and Yanagihara, H.（2008）
American Journal of Mathematical and Management Sciences（in press）．
Stone, M.（1974）. Cross-validatory choice and assessment of statistical predictions, Journal of theRoyal Statistical Society: Series B, 36, 111–147.
. An asymptotic equivalence of choice of model by cross-validation and Akaike’s
Stone, M.（1977）
criterion, Journal of the Royal Statistical Society: Series B, 39, 44–47.
丹後俊郎（2002）. 『医学データ―デザインから統計モデルまで』，共立出版，東京．
谷村 晋（2010）. 『地理空間データ分析』，共立出版，東京．
Tiwari, R. C., Ghosh, K., Jemal, A., Hachey, M., Ward, E., Thun, M. J. and Feuer, E. J.（2004）. A
new method of predicting US and state-level cancer mortality counts for the current calendar
year, CA A Cancer Journal for Clinicians, 54（1），30–40.
冨田哲治，佐藤健一，柳原宏和（2010）. 空間データに対する交互作用モデルを用いた変化係数曲面の推測について，応用統計学, 39（2&3）, 1–12.

Proceedings of the Institute of Statistical Mathematics Vol. 59, No. 2, 217–237 (2011)

237

Cancer Mortality Risk Visualization on Age-period Plane
Ken-ichi Kamo1 , Tetsuji Tonda2 and Kenichi Satoh2
1 Center for Medical Education, Sapporo Medical University
2 Research Institute for Radiation Biology and Medicine, Hiroshima University

Cancer information obtained by a suitable method becomes a basis for planning aneﬀective cancer control program. Our purpose in this paper is to visualize cancer mortality
risk as the curved surface on the age-period plane. To achieve this objective, a geographically weighted generalized linear model (GW model) may be suitable by regarding theinformation for age and period as the geographical one on a two-dimensional plane. On
the other hand, we consider a parametric interaction model by age and period as a rivalof the non-parametric GW model. Moreover, we consider one more model, that is, an ageperiod-cohort model, which is widely used in the analysis and prediction of longitudinalcancer data. These models are described with the script which can be executed in R and
results for lung cancer mortality in Denmark for reappearance by the reader.
To verify the achievement of our objective, that is, cancer mortality risk visualization,
we also apply these methods to data of liver cancer mortality in Japan. Liver cancer issaid to have a cohort eﬀect such that the birth cohort around 1935 suﬀers a high risk. We
check whether this eﬀect can be visualized.
Finally, we discuss not only the statistical aspect, but also future problems for our
research objective. Each model has its merit, so it may be eﬀective to use plural methods
jointly to estimate the properties for cancer mortality risk.

Key words: Cancer, visualization, risk curved surface, geographically weighted generalized linear
model, interaction model, age-period-cohort model, cohort eﬀect.