生体システムの階層適合解析

生体システムの階層適合解析
博士（工学）論文提出先
大阪大学大学院基礎工学研究科
平成 23 年 9 月
堀尾 秀之
生体システムの階層適合解析
堀尾 秀之
内容梗概
生体は呼吸や消化など様々な機能が集合して生命活動を維持している．各生体
機能システムはさらに，独自機能を持ち，相互に関連を持って機能する複数のサ
ブシステムによって構成されている．生体機能の中でも血液循環システムは，他
の生体機能に対して物質の運搬を担っていることから重要な機能システムである．
血液循環システムは，血液を媒体として，加圧し駆出する心臓，経路およびの加
圧補助を行う体幹部四肢の血管，そしてガス交換・物資交換を行う末端の毛細血
管が階層状に関連し構成されている．
血液循環システムを対象として解析を行う場合，サブシステム毎の機能の複雑
さや分布は異なっている為，共通の計測手法を適用した場合，計測時間や情報量
が冗長もしくは不足してしまうといった問題がある．本研究では，多階層の血液
循環システムの各階層の運動等に適した計測情報・解析手法の開発を目指す．
各サブシステムで開発した手法を以下に述べる
心臓は血液を全身へ送り出す駆出機能を担う．心臓は全身へ血液を送り出すた
めに単純な収縮ではなく全体がねじれながら収縮している．ねじれながら収縮し
ているため，心臓の収縮は特定断面内だけではなく，三次元的な方向を示す．MRI
の速度計測は計測断面に対して任意に設定した方向の速度が計測可能であり，複
数回計測し，位置合わせを行うことで三次元的な速度を取得することが可能であ
るといった利点がある．速度情報は運動を直接解析可能であることからも MRI の
速度計測を用いる．心臓の複雑なねじれ運動を解析するために心臓全体のねじれ
を示すねじれ角・ねじれ率の解析と心臓各部の局所ねじれを示すひずみの解析を
行う．ねじれ角・ねじれ率の解析からねじれ運動の時間変化を，ひずみ解析から
各部のねじれ分布の時間変化を獲得し心臓のねじれ運動を解析することが可能と
なった．
i
四肢の血管は，周囲の筋肉の収縮・弛緩による血管圧迫によって，血液への加
圧補助機能を担う．加圧補助機能を解析するために血管の圧迫動態の観察と筋肉
の収縮状態の判定を同時に行う．非侵襲で高コントラストな形態情報を取得でき
るなどの利点から，MRI 計測を用いる．計測画像の信号強度ヒストグラムは筋
肉・血管・脂肪起因の信号が重合しており，各体組織由来の信号強度のばらつき
はガウス分布となると考えられる．断面像の信号強度ヒストグラムを混合ガウス
分布へ近似する．近似した混合ガウス分布の標準偏差等の差異から筋肉組織に相
当するガウス分布を判定する．筋肉の収縮・弛緩で筋肉組織の分布の比較から，
血管圧迫の観察を可能とした．また，筋肉組織に相当するガウス分布の重み成分
から計測画像の収縮弛緩状態の判定が可能となった．
毛細血管はガス交換・物質交換機能を担う．ガス交換・物質交換には，血液が
供給・回収される血液循環の有無が重要であり，血液循環の有無を判定するため
に，血液の圧変化による血管の脈動を測定する．血管の脈動は単純な膨張・収縮
運動であるため，血管の直径を簡易に計測できる透過画像を用いる．血液に含ま
れるヘモグロビンは近赤外線に対する吸光率が高いことから，光源として近赤外
線 LED を用いて近赤外線透過画像を計測する．体組織毎の近赤外線吸光率の差
による透過画像の輝度値の差異から血管部分を抽出する．抽出した血管部分に対
して雑音除去と位置合わせを行い，血管太さの時間変化を検出する．検出した血
管の時間変化から，脈動周期を取得することが可能となった．
本研究では，機能・分布の異なるサブシステムである心臓・血管・毛細血管に
対して，異なる計測・解析手法を適用することで，多階層システムである血液循
環システム解析を行った．本手法のように各サブシステムの機能等に適した解析
手法を開発することで，機能と構造を考慮した生体システム全体の最適な機能解
明手法に寄与すると考える．
キーワード
医用画像計測, MRI, 近赤外線, 画像処理, 収縮解析, 捻転運動
ii
Adaptive Multi-layered Analyses for Biological
Systems
Hideyuki Horio
Abstract
A human has a lot of biological systems, for examples, respiration，digestion，
metabolism，circulation, and so on. Circulation is one of biological functional
systems and consisting of some subsystems organized with heart, blood vessels,
peripheral vessels and so on. The motion of the organizations is different, i.e., a
heart deforms once in a second and vessels change the diameter slightly. On the
other hand, a human has only one heart and there are innumerable vessels in a
human body. We have no methodology to evaluate each organization, but we need
to prepare the analysis methods for the multi-layered organizations adaptively.
The thesis presents the development of the adaptive multi-layered analyses for
blood circulation system featuring a heart, blood vessels and peripheral vessels.
A heart beats once in a second and deforms with shrinking and twisting in
order to pump blood to the whole human body. First, the method was proposed
for a heart to analyze the twist performance of a left ventricle. MR–PC (magnetic
resonance phase contrast was selected to measure the three dimensional velocity
of the cardiac muscles on the short axis section of the heart. For the detailed
analyses, two methods were developed. For the former method, time change
of the global twist was analyzed by using the axis of rotational movement of the
heart. As the result, the twist rate of the endocardium shows the maximum value
in the end systole. For the latter, time change of the local twist was evaluated
by calculate in the principal strain and maximum shear strain. The experimental
results demonstrate that the twist between sections was large in the early systole
iii
and diastole and that the twist on a section was large in the end systole. From
these results, the method has the possibility to evaluate the pumping function of
a heart.
A blood vessel in limbs helps to pressure blood for the circulation to periphery.
The shrinking with the surrounding muscles enables to pressure blood. Second,
the method was prepared for a blood vessel T1-weighted MR image was obtained
with fist held and relaxed, respectively. The intensity histogram of the MR image
was approximated with four Gauss distribution profiles using the EM (expectation
maximization) algorithm. The held and relaxed condition could be distinguished
with comparing between statistical values of two profiles corresponding to blood
vessels and muscles. The experimental results indicate that the method is thought
to monitor the muscle performance to pressure blood.
The peripheral vessel has the function to transfer and collection O2, CO2
gas. The amount of gas exchange is proportional to heart beat. Third, the
method was developed to count beat rate with transmitted near infrared ray
image sequence. The temporal characteristics of the capillary thickness were
investigated with the image processing for the image dataset. Then, applying
the Fourier transform to the characteristics, the beat of peripheral vessel was
estimated. The estimated value shows the good agreement with the value with
the conventional pulse monitor. This confirms that the method has the potential
to check the performance of peripheral vessel function.
As mentioned above, by adapted analyses for each subsystem of the circulation, we acquired information on the biological functional in necessary minimum
time. We think that it becomes possible to analyze the entire blood circulation
efficiently by using the analyses of this thesis. In addition, our methodology will
be expected to lead to other functional systems and analyses in the future.
Keywords:
Medical Imaging, MRI,Near Infrared, Image Processing, Contraction Analysis,
Twist Movement
iv
目次
1
序論
2
血液循環システム
3
4
13
2.1
血液循環システムの概要 . . . . . . . . . . . . . . . . . . . . . . .
14
2.2
心臓の特徴
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
15
2.3
血管の特徴
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
19
波動計測
23
3.1
波動 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
24
3.2
波動の物理現象 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
28
3.3
生体計測技術 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
32
3.4
生体画像計測 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
36
MR 速度画像を用いた捻転解析
41
はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
42
4.2 MR–PC 法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
43
4.3
計測画像 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
46
4.4
三次元速度場の獲得
. . . . . . . . . . . . . . . . . . . . . . . . .
48
4.5
回転の中心軸を用いたねじれ解析 . . . . . . . . . . . . . . . . . .
50
4.6
ひずみの方向によるねじれ解析 . . . . . . . . . . . . . . . . . . .
63
4.7
おわりに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
75
4.1
5
3
MR 形態画像を用いた収縮解析
77
5.1
はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
78
5.2
筋収縮による信号強度分布変化 . . . . . . . . . . . . . . . . . . .
79
5.3
計測対象および計測条件 . . . . . . . . . . . . . . . . . . . . . . .
80
v
目次
6
7
5.4
解析の流れ
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
83
5.5
信号分布を用いた血管圧迫解析 . . . . . . . . . . . . . . . . . . .
84
5.6
信号強度を用いた筋収縮解析
. . . . . . . . . . . . . . . . . . . .
89
5.7
考察 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
93
5.8
おわりに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
93
近赤外線透過像を用いた拍動解析
95
6.1
はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
96
6.2
近赤外線透過画像計測装置 . . . . . . . . . . . . . . . . . . . . . .
97
6.3
透過画像の動的変化解析 . . . . . . . . . . . . . . . . . . . . . . . 102
6.4
脈拍情報の定量的評価 . . . . . . . . . . . . . . . . . . . . . . . . 114
6.5
おわりに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116
結論
117
謝辞
123
参考文献
125
研究業績
133
vi
図目次
1.1
多階層システム例 . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
1.2
生体多階層システム例 . . . . . . . . . . . . . . . . . . . . . . . .
6
1.3
血液循環システム概要図 . . . . . . . . . . . . . . . . . . . . . . .
8
1.4
各サブシステムの運動と分布の差異 . . . . . . . . . . . . . . . . .
9
1.5
計測時間と獲得情報量の差と最適化による変化 . . . . . . . . . . .
10
2.1
血液循環システムの概要 . . . . . . . . . . . . . . . . . . . . . . .
14
2.2
心臓の外形および各部位 . . . . . . . . . . . . . . . . . . . . . . .
15
2.3
心室および心房と弁位置 . . . . . . . . . . . . . . . . . . . . . . .
16
2.4
短軸断面の各部と心筋線維の走向方向
. . . . . . . . . . . . . . .
17
2.5
血管の模式図 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
19
2.6
毛細血管でのガス交換模式図
. . . . . . . . . . . . . . . . . . . .
20
3.1
疎密波の周波数による分類 . . . . . . . . . . . . . . . . . . . . . .
26
3.2
電磁波の周波数による分類 . . . . . . . . . . . . . . . . . . . . . .
27
3.3
放射の模式図 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
28
3.4
吸収の模式図 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
29
3.5
反射の模式図 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
30
3.6
核磁気共鳴の模式図
. . . . . . . . . . . . . . . . . . . . . . . . .
31
3.7
可視光を用いて計測された眼底画像 . . . . . . . . . . . . . . . . .
32
3.8
超音波ドプラを用いて計測された心臓画像 . . . . . . . . . . . . .
33
3.9 X 線 CT により得られた腹部計測画像 . . . . . . . . . . . . . . . .
34
3.10 MRI により得られた脳計測画像 . . . . . . . . . . . . . . . . . . .
35
3.11 X 線により得られた頭部透過画像 . . . . . . . . . . . . . . . . . .
36
3.12 MRI により得られた再構成断層像 . . . . . . . . . . . . . . . . . .
37
3.13 MRI により得られた三次元再構成像 . . . . . . . . . . . . . . . .
38
vii
図目次
viii
4.1
MR–PC 法のパルスシーケンス . . . . . . . . . . . . . . . . . . .
43
4.2
移動しない水素原子に対して速度エンコード磁場を印加した場合 .
44
4.3
移動している水素原子に対して速度エンコード磁場を印加した場合 45
4.4
3 方向の速度分布と対応する形態画像 . . . . . . . . . . . . . . . .
47
4.5
注目領域の抜きだし処理 . . . . . . . . . . . . . . . . . . . . . . .
49
4.6
回転の中心軸と断面
. . . . . . . . . . . . . . . . . . . . . . . . .
50
4.7
断面内ねじれ . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
51
4.8
断面間ねじれ . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
51
4.9
計測画像と左室心筋抽出結果
. . . . . . . . . . . . . . . . . . . .
52
4.10 心内膜側と心外膜側への分割 . . . . . . . . . . . . . . . . . . . .
⃗ と各点の関係 . . . . . . . . . . . . . . . . . . . . .
4.11 回転中心軸 C
4.12 速度 V⃗ と円筒座標系速度 . . . . . . . . . . . . . . . . . . . . . . .
53
4.13 左心室内腔に仮定した中心軸
. . . . . . . . . . . . . . . . . . . .
56
4.14 収縮初期，心尖部の回転中心
. . . . . . . . . . . . . . . . . . . .
57
4.15 収縮初期，心基部の回転中心
. . . . . . . . . . . . . . . . . . . .
57
4.16 専門医による心臓の中心 . . . . . . . . . . . . . . . . . . . . . . .
58
4.17 本手法で得られた回転の中心
. . . . . . . . . . . . . . . . . . . .
58
4.18 専門医による中心位置と本手法による中心位置の比較 . . . . . . .
58
4.19 心内膜と心外膜の回転角の時間変化 . . . . . . . . . . . . . . . . .
59
4.20 心内膜側に対する心外膜側のねじれ角
. . . . . . . . . . . . . . .
59
4.21 心基部に対する心尖部のねじれ角 . . . . . . . . . . . . . . . . . .
61
4.22 心基部に対する心尖部のねじれ率 . . . . . . . . . . . . . . . . . .
61
4.23 垂直ひずみとせん断ひずみ . . . . . . . . . . . . . . . . . . . . . .
63
4.24 捻転運動の矢印表示 . . . . . . . . . . . . . . . . . . . . . . . . .
66
4.25 断面内ねじれ数値モデル . . . . . . . . . . . . . . . . . . . . . . .
68
4.26 断面内ねじれ数値モデルにおける第 1 主ひずみの方向 . . . . . . .
69
4.27 断面内ねじれ数値モデルにおける最大せん断面の法線方向 . . . .
69
4.28 断面間ねじれ数値モデル . . . . . . . . . . . . . . . . . . . . . . .
70
4.29 断面間ねじれ数値モデルにおける第 1 主ひずみの方向 . . . . . . .
71
4.30 断面間ねじれ数値モデルにおける最大せん断面の法線方向 . . . .
71
4.31 左心室における第 1 主ひずみ方向の時間変化 . . . . . . . . . . . .
72
54
54
図目次
4.32 左心室における最大せん断面法線の時間変化 . . . . . . . . . . . .
73
4.33 ねじれ方向が z 軸方向を示した領域の割合時間変化 . . . . . . . .
74
5.1
筋収縮の模式図 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
79
5.2
筋収縮による血管圧迫の模式図 . . . . . . . . . . . . . . . . . . .
79
5.3 Gradient-Echo 法のパルスシーケンス図 . . . . . . . . . . . . . . .
80
5.4
計測画像 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
82
5.5
解析の流れ
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
83
5.6
計測画像とヒストグラム . . . . . . . . . . . . . . . . . . . . . . .
85
5.7
混合ガウス分布近似
. . . . . . . . . . . . . . . . . . . . . . . . .
86
5.8
筋線維と周辺毛細血管領域の輝度値範囲 . . . . . . . . . . . . . .
87
5.9
計測画像の領域分割結果 . . . . . . . . . . . . . . . . . . . . . . .
88
5.10 領域抽出結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
90
5.11 分割領域のヒストグラムと近似結果 . . . . . . . . . . . . . . . . .
90
5.12 全計測結果に対する解析結果
. . . . . . . . . . . . . . . . . . . .
92
6.1
脈拍計測場所 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
96
6.2
近赤外線透過画像計測装置の概観 . . . . . . . . . . . . . . . . . .
97
6.3
画像計測台
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
98
6.4 CCD カメラ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
98
6.5
近赤外線照明機 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
99
6.6
波長特性 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
6.7
中指計測画像例 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
6.8
解析処理の流れ . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102
6.9
注目領域の抜き出し
. . . . . . . . . . . . . . . . . . . . . . . . . 103
6.10 閾値処理の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104
6.11 計測画像の注目領域に対する閾値処理画像 . . . . . . . . . . . . . 105
6.12 膨脹・縮退処理の概要 . . . . . . . . . . . . . . . . . . . . . . . . 106
6.13 膨脹・縮退処理前後の画像 . . . . . . . . . . . . . . . . . . . . . . 107
6.14 細線化処理前後の画像 . . . . . . . . . . . . . . . . . . . . . . . . 109
6.15 血管太さの時間変化 . . . . . . . . . . . . . . . . . . . . . . . . . 112
6.16 太さの時間変化をフーリエ変換した強度グラフ . . . . . . . . . . . 113
ix
図目次
6.17 本手法と脈拍計の測定結果比較 . . . . . . . . . . . . . . . . . . . 114
x
7.1
本論文の計測と解析手法 . . . . . . . . . . . . . . . . . . . . . . . 117
7.2
血液循環システムに対する獲得情報 . . . . . . . . . . . . . . . . . 119
7.3
呼吸システムの機能と分布 . . . . . . . . . . . . . . . . . . . . . . 120
7.4
呼吸システムに対して獲得が予想される情報 . . . . . . . . . . . . 120
7.5
階層適合解析の展開例 . . . . . . . . . . . . . . . . . . . . . . . . 121
表目次
2.1
世代による心拍数変化 . . . . . . . . . . . . . . . . . . . . . . . .
18
2.2
血液循環のサブシステムの運動と分布量の違い . . . . . . . . . . .
21
3.1
画像種による計測時間と情報量の違い
. . . . . . . . . . . . . . .
38
4.1
速度画像仕様 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
46
4.2
形態画像仕様 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
46
5.1
計測条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
81
5.2
屈筋領域の混合ガウス分布のパラメータ . . . . . . . . . . . . . .
91
5.3
伸筋領域の混合ガウス分布のパラメータ . . . . . . . . . . . . . .
91
5.4
全解析結果における重み成分の平均と分散 . . . . . . . . . . . . .
92
6.1
使用機器仕様 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
6.2
計測画像仕様 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
1
第1章
序論
自身は複数の製品を製造・販売している企業で開発員として所属している．製
品の開発・販売を行う企業の価値は開発した各製品の機能や利益等による評価に
よって決定されることが多い．製品の機能や利益率は多数の開発員や販売員の努
力や能力の結晶であり，各人の能力の高さではなく，開発集団としての能力の高
さが重要である．製品の価値を高める為には集団の力が必要であるが、集団の力
を向上させるためには自身のような開発員個人の能力や連携力の向上が重要であ
る．末端開発員である自身の能力の評価や向上は一見して企業全体の価値向上に
直接関係ないようであるが，開発集団，製品を通して企業全体へ影響を及ぼして
いる．
天体の軌道といった自然現象から，電力供給などの社会インフラ，人体の生命
活動に至るまで世界中の様々な現象・活動は，サブシステムと呼ばれる複数の事
象が集まり，相互に関係することで一つの集合であるシステムとして機能してい
る．さらにシステム毎にサブシステムの構成・関連性は異なっている．
さまざまなシステムの中でも，サブシステム毎に機能だけでなく，分布が異な
り，より末端のサブシステムではより分布密度が高くなり，特定のサブシステム
を集積・分解しても別のサブシステムとして機能できないシステムが存在する．
上記のようなサブシステムの階層毎に分布構造や機能が異なり，かつ特定のサブ
システムが別のサブシステム機能を担うことが出来ないシステムを，構造上の階
層を持つ多階層システムと定義する．
3
第 1 章 序論
図 1.1 に多階層システムである電力システム，天体軌道システム，血液循環シ
ステムを示す．
Electric Power
System
Solar System
Electric
Generation Plant
Sun
Transformer Station
Transformer Station
Pole Pole
Transformer
Pole
Transformer
Pole
Transformer
Transformer
Pole Pole
Transformer
Pole
Transformer
Pole
Transformer
Transformer
Saturn
Earth
Titan
Moon
Rhea
Consumer Consumer
(a)
Consumer Consumer
(b)
Blood Circulation
1
Heart
Vessel
Vessel
BloodBlood
capillary
Blood
capillary
Capillary
BloodBlood
capillary
Blood
capillary
Capillary
1
(c)
図 1.1 多階層システム例
(a) 太陽軌道システム (b) 電力システム
(c) 血液循環システム
1
図 1.1(a) に示したように電力システムにおいて，発電所で発電された電力は，
変圧所等で変圧され，最終的に一般家庭等で消費される．発電所から変圧所，一
4
般家庭へと末端になるほど分布量は増大し，一般家庭を集積しても発電所の機能
を担うことは出来ないことから，電力システムは構造上の階層をもつ多階層シス
テムである．
図 1.1(b) に示したように太陽系の軌道は恒星である太陽を中心として，惑星で
ある地球や土星は固有の円周軌道上を，固有の円周周期で公転している．さらに
各惑星を中心としてそれぞれの衛星である月やレア・タイタンも固有の円周起動
上を固有の円周周期で公転している．上記のように天体軌道も構造上の階層をも
つ多階層システムとすることが可能である．
図 1.1(c) に示したように生体機能の一つである血液循環システムでは，心臓か
ら駆出された血液は，四肢の血管を経由することで全身に供給され，毛細血管で
代謝物交換・ガス交換を行う．心臓から四肢の血管，毛細血管へと分布量は増大
し，毛細血管では心臓の機能を担うことは不可能である．これらの特徴から血液
循環システムも多階層システムと見なすことができる．
構造上の階層をもつ多階層システムを解析する場合，各階層の分布構造を考慮
した計測・解析が必要となる．例えば電力システムにおいては，発電所と一般家
庭に対して同じ計測・解析を行うのではなく，発電所に対しては最大発電量や安
定性等の各所の詳細情報が重要である．しかし，発電所に対して一般家庭におい
ては個別の詳細情報ではなく，正常の電気が供給されているかなどの情報が重要
である．
上記のように多階層システムにおいては分布密度が低く，重要な機能を担うサ
ブシステムの計測・解析は詳細に行い，分布密度が高く，末端機能を担うサブシ
ステムに対しては一定範囲を計測し，必要最小限の解析を行うことが重要である．
5
第 1 章 序論
人体は，運動や呼吸，消化・吸収など複数のシステムが機能することによって，
生命活動を維持している．人体の生命機能においても血液循環だけでなく多階層
システムと見なすことが可能な機能は多岐にわたる．図 1.2 に生体における多階
層システムの例を示す．
Blood Circulation
Respiration
Heart
Costa and Diaphragm
Vessel
Vessel
Bronchi
Bronchi
BloodBlood
capillary
Blood
capillary
Capillary
BloodBlood
capillary
Blood
capillary
Capillary
Alveoli
Alveoli
Alveoli
Alveoli
Alveoli
Alveoli
(a)
(b)
Digestion Absorption
Postural Control
Palate
1
1
Central Nerve
Belly
Pancreas
Liver
Semicircular Canal
Intestinal
Intestinal
villus
Intestinal
villus villus
Oculus
Feet
Intestinal
Intestinal
villus
Intestinal
villus villus
(c)
(d)
図 1.2 生体多階層システム例
(a) 血液循環システム (b) 呼吸システム
(c) 消化吸収システム (d) 姿勢制御システム
1
6
1
図 1.2(a) の血液循環システムは先に示した通り，多階層システムである．
図 1.2(b) に示したように呼吸システムにおいて，横隔膜等で肺全体を減圧し，
気管支を経由し空気が肺胞へ供給される．肺胞ではガス交換が行われることで血
液に酸素が供給される．横隔膜，気管支，肺胞の分布量や機能の独立性から，呼
吸システムも多階層システムと見なす．
また，図 1.2(c) の消化吸収システムは，咀嚼を行う口蓋や，消化酵素を生成を
担うすい臓などの臓器，吸収を行う腸壁などから構成された生体機能である．咀
嚼は一箇所で行われるのに対して，消化酵素を生成する臓器は複数あり，さらに，
吸収を行う腸壁の腸柔毛は無数に分布している．各機能や分布の違いから消化吸
収システムも多階層システムである．
さらに，図 1.2(d) の姿勢制御システムでは，周囲の状態を知覚する足や眼，平
行状態を認識する三半規管，各種情報を統合する中枢神経系から構成される生体
機能である．中枢神経系と各感覚器の機能や分布の違いから立位制御システムも
多階層システムと見なすことができる．
上記のように生体機能においても構造上の階層をもつ多階層システムは多岐に
渡っている．
人体の機能の一つである血液循環は，呼吸機能のガス交換によって赤血球に保
持された酸素の運搬や，消化・吸収機能によって体内に吸収された物質の運搬な
どを行う重要な機能であり血液循環システムに着目する．
7
第 1 章 序論
血液循環システムを図 1.3 に示す．血液循環は血液を供給・回収する一つの循
環システムであり，血液は心臓から体幹部・四肢の血管を経由して末梢の毛細血
管へと循環し，各サブシステムが階層的に関係して構成される．
Blood Circulation
Heart
Vessel
Vessel
Blood capillary
Blood capillary
Blood Capillary
Blood capillary
Blood capillary
Blood Capillary
図 1.3 血液循環システム概要図
血液循環における各サブシステムの運動と分布の例を図 1.4 に示す．心臓は血
液循環システムの中枢であり，駆出機能を担う唯一のサブシステムである．心臓
は収縮によって血液へ加圧し，全身へ血液を供給している．心臓の収縮運動は，
全体がねじれながら収縮するといった複雑な運動を行っている．
四肢の血管は血液を全身へ供給する経路機能を担うサブシステムであり，全身
の複数の部分に存在する．また四肢の血管は血液への加圧補助機能も担っている，
血液への加圧補助機能は心臓の収縮だけでは全身からの血液を再び心臓まで輸送
8
1
することが困難であることから重要な機能である．血液への加圧補助機能は血管
周囲の筋肉の収縮・弛緩によって血管が圧迫されることによって実現している．
末梢の毛細血管は筋肉・臓器へのガス交換・物質交換を担うサブシステムであ
り，全身に無数に分布している．毛細血管において，ガス交換・物質交換を行う
ためには，酸素や栄養素を含んだ血液が供給され，老廃物等を含んだ血液が回収
される血液循環が行われていることが重要である．血液循環の有無は，血液の圧
力変化による血管の脈動で識別することが可能である．
以上のように血液循環システムの各サブシステムは，分布量や機能も異なる．
Many
Spot
Blood capillary
Vessel
Heart
Few
Movement
Simple
Complex
図 1.4 各サブシステムの運動と分布の差異
1
9
第 1 章 序論
多階層システムの機能を明らかとするため，全体を一つとして解析・評価する
手法や，各サブシステムの機能と相互関係から全体を解析する手法が提案されて
いる．各サブシステムの機能解析では対象に適した解析を選択することが重要で
ある．
生体システムを対象とした機能解析を行うために様々な計測技術が開発され，
計測時間や獲得できる情報量などの特徴が異なっている．超音波計測における計
測時間と獲得可能な情報量を図 1.5 に示す．さらに，解析を行い，計測手順を最
適化することで，より短時間で計測を行ったり [2]，もしくは三次元といったより
詳細な情報を獲得することが可能である [3]．
Fast
Measurement
Speed
Modify protocol and analysis
Slow
US
Information
Little
Much
図 1.5 計測時間と獲得情報量の差と最適化による変化 [2][3]
超音波計測は短時間で形態情報を獲得可能であるが，高コントラスト画像を獲
得することは不可能である．さらに，計測手順の変更や，画像処理を行うことで，
三次元形態画像やより短時間で形態情報を獲得することが可能となる．しかしな
がら，三次元速度分布や生体代謝物質の分布を直接計測することは不可能である．
上記のように計測手順・解析処理を最適化しても必要な情報が取得することが
不可能な場合がある．したがって，解析対象に応じて計測技術の選択・解析処理
1
10
の最適化を行う必要がある．
血液循環システムの各サブシステムに対して一意の解析手法を適用した場合，
ねじれながら収縮する心臓の解析が不足，もしくは単純な拍動である末梢毛細血
管の解析が冗長になってしまうといった問題がある．
本研究では，血液循環システムの中枢から末端にかけて，各サブシステムの構
造上の分布を考慮した計測解析を適用する手法を提案する．具体的には唯一の臓
器である心臓に対しては，三次元再構成画像を用い，収縮におけるねじれ・ひず
み解析を行う．つぎに，複数が同時に分布している四肢の血管に対しては，再構
成断層像を用い，周囲筋肉の収縮・弛緩解析を行う．さらに広範囲に無数に分布
している末梢血管に対しては，透過画像を用い，拍動解析を行う．各サブシステ
ムに対応した解析手法を適用することで多階層システムである血液循環システム
の機能解析を目指す．第 2 章では血液循環システム概要を述べ，第 3 章では解析
に用いる波動計測について述べる．さらに，第 4∼6 章で提案の解析手法と解析
結果について述べる．第 7 章で本論文をまとめる。
11
第2章
血液循環システム
生体は生命を維持するために呼吸や消化・吸収，循環など様々な機能を有して
いる．各機能はそれぞれが一つのシステムとして機能している．更に各システム
は複数のサブシステムから構成されている．生体が有する機能のなかでも血液循
環は，各臓器・細胞の活動を維持する上で不可欠なシステムである．血液循環シ
ステムは血液に加圧する心臓，運搬経路となる動脈・静脈，酸素や栄養・老廃物
などを細胞と交換する毛細血管などのサブシステムから構成され，一つのシステ
ムとして機能している．本章では，2.1 節で血液循環システムの概要を述べ，2.2
節と 2.3 節で血液循環システムを構成する各サブシステムである心臓と血管の特
徴について述べる．
13
第 2 章 血液循環システム
2.1
血液循環システムの概要
図 2.1 に血液循環システムを示す．循環の媒体は血液であり，血液は心臓によっ
て加圧される．心臓によって加圧された血液は血管を経路として，全身の臓器へ
酸素や栄養素を運搬する．さらに，各臓器から二酸化炭素や老廃物を回収し，ガ
ス交換や老廃物のろ過を行う臓器に運搬される．ガス交換を行う肺への循環は肺
循環と呼ばれ，全身への酸素・栄養運搬を行う循環は体循環と呼ばれる．
Muscle
Blood Vessel
Lung
(Gas Exchange)
Heart
Liver
(Nutrient Exchange)
Blood Capillary
図 2.1 血液循環システムの概要
心臓の拍動による血液の圧の変化は様々な部位の血管網での同期した圧変化を
引き起こし，血管網での圧変化を知ることにより心臓の拍動リズムを知ることが
で出来る．拍動リズムを計る方法には血管音で計測する聴診法と血管の拍動を計
測する触診法等がある [4]．脈拍の測定には首の頚動脈，上腕や手首，指等の動脈
が計測対象とされている．
循環の媒体である血液の成分は赤血球，血小板，白血球の細胞成分が
46 ％，血
1
漿が 54 ％を占める．赤血球は酸素運搬を行なう直径約 7µm，厚み 3µm 程度の細
胞成分である．赤血球は酸素が透過できる半透過性を持ち，乾燥重量換算で 90
％程度のヘモグロビンと水分で構成されている．ヘモグロビンは鉄を含むタンパ
ク質であり，酸化することにより酸素の運搬が行なわれる [5]．酸素と結合したヘ
モグロビンを酸化ヘモグロビン，酸素を分離したヘモグロビンを還元ヘモグロビ
ンという．
14
2.2. 心臓の特徴
2.2
心臓の特徴
心臓は血液循環システムにおいて，血液に加圧を行うサブシステムであり，全
身へ血液を供給するために，非常に効率的かつ複雑な収縮運動を行っている．心
臓の収縮拡張運動は体液循環の源であり，心臓の運動不全による心疾患は主要な
死亡原因の一つでもある [6]．
図 2.2 に心臓の外形および各部位の名称を示す．心臓の形状は丸みを帯びた円
錐状をしており，円錐の先端部分は人体の前方左下方向を向いている．また心臓
各部の名称として，円錐の先端部分は心尖部，円錐の底の部分は心基部と呼ばれ
る．さらに，心基部と心尖部をつなぐ直線を心長軸，心長軸に垂直な面を短軸断
面と呼ぶ．心臓の大きさは個人によって異なり，人間のこぶし程度の大きさであ
る [4][7]．
Base
Long Axis
Apex
Short Axis Plane
図 2.2 心臓の外形および各部位
15
第 2 章 血液循環システム
図 2.3 に心臓の各心室と心房の位置および弁の位置と血液の流れを示す．心臓
は左心と右心の二つに分けられ，各部分は心房と心室で構成される．心房は心臓
に灌流してきた血液を受け入れる部屋として働き，心室は心臓から血液を送り出
すポンプとして働く．左心室と左心房の境界には僧帽弁があり，右心室と右心房
の境界には三尖弁がある．この二つの弁によって心室から心房への血液の逆流を
防いでいる．左心室と左心房で構成される左心は，全身に酸素と代謝基質を供給
し，二酸化炭素と代謝産物を除去する大循環を駆動し，右心室と右心房で構成さ
れる右心は，ガス交換のための小循環を駆動している [4][7]．
Right Atrium
Left Atrium
Arterial Blood
Venous Blood
Mitral Valve
Tricuspid Valve
Left Ventricle
Right Ventricle
図 2.3 心室および心房と弁位置
16
2.2. 心臓の特徴
図 2.4 に心臓に対する短軸断面と，心内膜と心外膜，そして心筋線維の走向方
向を示す．心筋壁は心内膜，心筋層，心外膜の三層で構成される．心内膜と心外
膜は心筋層を内と外から滑らかな膜として覆っている．心筋層は心筋線維と結合
組織である心臓骨格で構成され，心筋線維の走向方向は心臓の外膜側と内膜側で
異なっている [4][7]．外膜側において，心筋線維は心尖部側からみて時計回りに螺
旋状に走向している．一方，内膜側における心筋線維は反時計回りに走向してい
る．さらに，内膜側と外膜側の心筋線維の走向方向は約 120deg 異なる [8][9]．
Epicardium
Midcardium
Endocardium
Short Axis Plane
Base
Endocardium
Epicardium
Apex
Cardiac Fiber
Fiber Direction
図 2.4 短軸断面の各部と心筋線維の走向方向
17
1
第 2 章 血液循環システム
心臓の運動周期には収縮期と拡張期がある．収縮期のはじめに心筋の緊張が心
室内に圧上昇を引き起こし，血液が動脈へ駆出される．収縮期に続く拡張期では
心筋の弛緩が起こり，心房と心室間の弁は開かれることで心房から心室に血液が
流入する [4][7]．
心臓周期における心臓の収縮・弛緩は，心長軸方向，心臓全体の重心方向，心
長軸を中心とした円周方向に生じ，さらに，心長軸方向に対して心筋のねじれが
生じている．ねじれは心臓の拡張や収縮に大きな影響を与えていると考えられる
[10][11]．
年齢による心拍数の変化を表 2.1 に示す．心臓の収縮・弛緩周期である心拍数
は，世代による基礎代謝の違いから変化する．また拍動数は日常生活でも睡眠時
や入浴時等で変化することが知られている．
表 2.1 世代による心拍数変化 [12]
Generation Pulse Rate(Times/min)
18
Newborn
120∼140
Infant
90∼130
Adult
70∼80
2.3. 血管の特徴
2.3
血管の特徴
血管は血液循環システムにおいて，主に運搬経路を担うサブシステムである．
血管には各臓器や組織で血管網があり，部位によって血管網や太さにも違いがあ
る．たとえば，腎臓では尿を漉すという目的に最適な球状血管網を構成し，手な
どの筋肉の血管は筋線維の並びに平行に走っている [13]．
図 2.5 に血管の概要を示す．血管は血液の通る粘弾性管血管壁に被われており，
ある程度の伸縮性を持つ．さらに，太い大動脈は直径 30mm，細い毛細血管は 1µm
と太さも多様である [4]．血管は内膜・中膜・外膜から構成され，動脈血管は拍動
性の圧力変化に対しても形状を保てる構造である．また，静脈血管は場所によっ
て静脈弁を持ち，血液の逆流を防ぐ構造となっている．
Adventitia
Media
(External Elastic Lamina)
(Smooth Muscle Cell)
Intima
(Internal Elastic Lamina)
図 2.5 血管の模式図
さらに，心臓の拍動だけでは全身の血液を再び心臓へ運搬することは不可能で
あることから，血液への加圧補助も血管の重要な役割である．しかし，血管自体
に収縮機能はないため，血管周囲の筋肉の収縮に伴う血管の圧迫によって血液の
加圧を行っている．筋肉による加圧機能は筋ポンプ機能と呼ばれる．
19
1
第 2 章 血液循環システム
図 2.6 に毛細血管におけるガス交換と物質交換の概念を示す．毛細血管は全身
に分布し，内皮細胞だけで構成され，細胞とのガス交換や物資交換を行っている．
ガス交換では，血中の酸化ヘモグロビンから還元された酸素は血漿内に溶け出し，
各細胞へ供給される．また，細胞内の二酸化炭素は炭酸となって血漿内に溶け出
すことで回収を行っている．また，物質交換では，血中の栄養素と細胞内の老廃
物を行っている．
Metabolic
Product
CO2
Nutrient
CO2
O2
Metabolic
Product
図 2.6 毛細血管でのガス交換模式図
末端の血管では，ガス交換・物質交換を行うために，酸素と栄養分を含む血液
が供給され，二酸化炭素や老廃物を含む血液が回収される血液循環が重要である．
末端血管における血流の有無は血流の圧変化による血管の脈動から判断すること
が可能である．
上記のように血液循環システムを構成する各サブシステムは運動の複雑さや分
布量が異なっている．各サブシステムの運動の複雑さと分布量の違いを表 2.2 に
示す．
20
2.3. 血管の特徴
表 2.2 血液循環のサブシステムの運動と分布量の違い
Blood Flow Subsystem
Movement Number of Spot
Heart
Complex
Few
Large Vessel
Moderate
Medium
Peripheral Vessel
Simple
Many
したがって，血液循環を解析する場合，計測時間を短縮して必要な情報を取得
するには各サブシステムの運動分布に適した手法を選択することが重要である．
21
第3章
波動計測
波動とは，エネルギが伝播する現象を示しており，計測対象による反射や吸収，
または対象自体から放射された波動から対象の形状や状態，特性を計測すること
が可能である．生体を対象とした波動計測は切開などの破壊なしに内部構造・状
態の計測が可能であることから，運動機能の解析，疾病の診断，運動効果の評価
など多岐に渡って適用されている．本章では，波動の概要と物理現象，および波
動を用いた生体計測技術について述べる．3.1 節では波動の概要について述べ，3.2
節では波動に関する物理現象を述べる．さらに 3.3 節で生体計測技術について,3.4
節で画像計測について述べる．
23
第 3 章 波動計測
3.1
波動
波動はエネルギが空間を伝播する現象を示していることから，位置 x と時間 t
の関数によって記述される．位置と時間の関数を Ψ とし，時間微分と空間微分と
の関係を式 (3.1) に示す．
1 ∂ 2Ψ
(3.1)
s2 ∂t2
式 (3.1) において s は波動の速度であり，真空中の電磁波では光速 c となる．ま
△Ψ =
た，△ はラプラシアンを示す．式 (3.1) は波動方程式と呼ばれ，波の記述に用い
られる．
波動方程式を簡略化するため x 方向の位置のみに依存するとし、波動方程式の
一般解を求める．Ψ は位置と時間の関数であることから，位置の関数 Ψx と Ψt を
用いて Ψ = Ψx Ψt と表されるものとすると，式 (3.1) は
△(Ψx Ψt ) =
1 ∂ 2 (Ψx Ψt )
s2
∂t2
(3.2)
となり，△ は位置の 2 階微分であることなどから
1
∂ 2 Ψt
Ψ
x
s2
∂t2
1
1 1 ∂ 2 Ψt
△ Ψx = 2
Ψx
s Ψt ∂t2
Ψt △ Ψx =
(3.3)
(3.4)
と式変形することが可能である．式 (3.4) の両辺で変数が分離されていることか
ら，式 (3.4) は定数 Υ に等しいと置くことができる．
Υ=
1 1 ∂ 2 Ψt
1
△ Ψx = 2
Ψx
s Ψt ∂t2
(3.5)
したがって，式 (3.1) は，
となる．
24
Υ =
1
Ψx
△ Ψx
(3.6)
Υ =
1 1 ∂ 2 Ψt
s2 Ψt ∂t2
(3.7)
3.1. 波動
Υ < 0 において両式の一般解は任意定数 Ak (k = 0, 1, 2, 3) を用いて
Ψx (x) =
A0 cos(Υx) + A1 sin(Υx)
(3.8)
Ψt (t) = A2 cos(Υst) + A3 sin(Υst)
(3.9)
となり，さらに波長 l の定常波とすると固定端 x0 = 0 xl = l での境界条件
Ψx (x0 ) = 0, Ψx (xl ) = 0 から
A0 = 0, A1 sin(Υx) = 0, Υ =
kπ
(k = 0, 1, ..)
l
(3.10)
が求まり，
Ψx (x) =
A1 sin( kπl x)
Ψt (t) = A2 cos( kπl st) + A3 sin( kπl st)
となり，波動方程式は角周波数 ω =
(3.11)
(3.12)
kπ
s の関数となる．したがって，波動は周波
l
数によって記述することが可能である．波動の中でも運動エネルギが伝播する波
動は疎密波，電磁波エネルギが伝播する波動は電磁波と呼ばれ，それぞれ特性が
異なる．
25
第 3 章 波動計測
3.1.1
疎密波の特徴
疎密波は音波とも呼ばれ，伝播する媒体によって，波長が異なる．音波の周波
数による分類を図 3.1 に示す．したがって，周波数によって分類され，音波の中で
も人間が認識できる周波数は可聴音と呼ばれる．また，工学的には人間がコミュ
ニケーション以外で使用する音を超音波と定義する [14]．可聴音は発声計測に用
いられ，超音波はソナーや医療用超音波計測などに用いられている．
20 Hz
20 kHz
Audible Sound
×
1 10 1
×
1 10 2
×
1 10 3
×
1 10 4
Frequency (Hz)
図 3.1 疎密波の周波数による分類
26
3.1. 波動
3.1.2
電磁波の特徴
音波に対して電磁波の速度はほぼ一定である．電磁波の波長による分類を図 3.2
に示す．電磁波の波長は 10pm 以下から 1mm 以上と広く，ガンマ線から電波と
呼ばれ，適用範囲は広い [15]．X 線は透過率が高く，内部構造の計測などに用い
られ，赤外線は対象の温度測定などに用いられている．
×10
12
3
×
3 1016 Hz
Hz
×
3 1018 Hz
X-Ray
Radio
Wave
Ultraviolet
Ray
Infrared
Ray
×10
1
×
12
1
1 10 14
4
×10
Gamma
Ray
14
7.5
×10
16
1
×10
18
×10 Frequency (Hz)
14
Visible Ray
図 3.2 電磁波の周波数による分類
1
27
第 3 章 波動計測
3.2
波動の物理現象
波動はエネルギの伝播であり，物質と様々な相互作用を起こす．代表的な物理
現象を以下に示す．
放射
放射の概念と式を図 3.3 と式 (3.13) に示す．EHigh と ELow はエネルギ準位を示
し，h はプランク定数，ω は角周波数を示す．物質はエネルギ準位 EHigh からよ
り低いエネルギ状態 ELow に遷移した時に電磁波を放出する．エネルギ準位間の
遷移による電磁波の放出を放射と呼ぶ．
Energetic State
EHigh
h
2
ω
π
Energetic State
ELow
図 3.3 放射の模式図
EHigh − ELow =
28
hω
2π
(3.13)
3.2. 波動の物理現象
吸収
吸収の概要と式を図 3.4 と式 (3.14) に示す．Ein ，Eout は電磁波のエネルギとす
る．また，l は透過した距離，ϑ は物質と波長によって決まる吸収係数を示す．電
磁波は物質を透過する際，物質にエネルギを奪われ強度が減衰する．透過による
エネルギ減衰を吸収と呼ぶ．エネルギの減衰量は透過した物質と距離に依存する．
Transmitted Wave Energy
Eout
Incident Wave Energy
Ein
l
図 3.4 吸収の模式図
Eout = Ein exp(−ϑl)
(3.14)
29
第 3 章 波動計測
反射
図 3.5 と式 (3.15) に反射の概念と屈折の式を示す．θin と θout はそれぞれ入射角
と反射角を示し，ϵin は ϵout それぞれ入射側の屈折率と透過側の屈折率を示す．電
磁波は屈折率の異なる物資の境界面に入射した時に一部は透過して，一部は反射
する．反射波において境界面に対する入射角と反射角は等しいが，透過波は屈折
し，入射角とは異なる角度で進む．
Incident Wave
Refractive Index
ϵin
Reflected Wave
θ
in
Refractive Index
ϵout
θ
out
Reflected Wave
図 3.5 反射の模式図
ϵin sin θin = ϵout sin θout
30
(3.15)
3.2. 波動の物理現象
共鳴
図 3.6 と式 (3.16) に核磁気共鳴の概要とラーモア周波数の式を示す．ω は共鳴
角周波数，B は磁場強度，γ は磁気回転比を示す．物質はそれぞれ状態や形状に
よって固有の周波数を持ち，同じ周波数の電磁波を照射することによって，対象
のエネルギ状態が遷移する現象を共鳴と呼ぶ．共鳴のなかでも核磁気共鳴現象に
おける周波数は磁場強度と対象の原子核に依存し，ラーモア周波数と呼ばれる．
ω = γB
Magnetic Field
Intensity
B
図 3.6 核磁気共鳴の模式図
ω = γB
(3.16)
以上のように波動に関係する物理現象は多岐に渡る．波動を利用した計測にお
いては取得したい情報に応じた物理現象が選択されている．
31
第 3 章 波動計測
3.3
生体計測技術
波動を用いた生体計測技術は多数開発され応用されている．各計測技術は用い
る波動や物理現象が異なっている．代表的な生体計測技術の例を以下に示す．
光学計測
図 3.7 に光学計測によって得られた眼底血管を示す [16]．光学計測は紫外線や
可視光，赤外線を用いた計測を示す．光を計測対象に照射し，透過光・反射光を
計測することで対象の形状や反射率や吸光率等の情報を獲得する．光学計測は実
時間で計測が可能であり，計測装置の準備が比較的簡便などの利点があるが，表
面や浅部のみで深部の計測が困難などの欠点もある．
図 3.7 可視光を用いて計測された眼底画像 [16]
32
3.3. 生体計測技術
超音波計測
図 3.8 に超音波ドプラを用いた計測により得られた心臓画像例を示す [17]．超
音波計測は疎密波を用いた計測を示す．音波照射から計測対象内の組織境界面で
の反射波受信までの時間差から照射位置からの距離を獲得する．また，照射対象
の移動によって反射波の周波数が変化するドプラ効果を応用することで、対象の
移動速度を獲得することも可能である．超音波計測は，非侵襲的に体深部の情報
を短時間で計測可能などの利点があるが，体組織間の信号コントラストが低いこ
とや骨を超えた領域の計測が不可能などの欠点もある．
図 3.8 超音波ドプラを用いて計測された心臓画像 [17]
33
第 3 章 波動計測
X 線計測
図 3.9 に X 線 CT 計測により得られた腹部画像例を示す [18]．X 線計測は紫外
線からガンマ線の間の周波数の電磁波である X 線を用いた計測を示す．X 線を計
測対象に照射し，透過した X 線量を計測することで内部構造を獲得する．また，
複数方向から計測し，再構成することで断層像の獲得や，断層像を複数計測する
ことで三次元再構成像を獲得することも可能である．X 線は透過率が高く、体組
織間の吸収率も異なるため，体深部の情報を高コントラストで計測可能などの利
点がある．しかし，X 線は生体の DNA を傷つけるなどの侵襲性があり，高頻度
の計測や胎児への計測が行えないなどの欠点もある．
図 3.9 X 線 CT により得られた腹部計測画像 [18]
34
3.3. 生体計測技術
MRI 計測
図 3.10 に MRI 計測により得られた脳画像例を示す [19]．MRI 計測は赤外線よ
りも低い周波数の電磁波であるラジオ波を用いた計測を示す．MRI 計測では，対
象を高磁場内に配置し、磁場強度と対象原子核に固有の周波数でラジオ波を照射
することで核磁気共鳴現象を引き起こす．ラジオ波の照射を停止した後，複数の
傾斜磁場を印加し，対象からの緩和信号を計測し、再構成することで内部構造等
を獲得する．MRI 計測は形態情報だけでなく，緩和信号の減衰時定数の違いや，
周波数分布から機能情報を獲得することも可能であり，計測が非侵襲的であるな
どの利点がある．しかし，高磁場を用いることから骨固定用ボルトなど体内の金
属が画像アーチファクトになり、計測不可能であることや，画像再構成に必要な
計測回数が多く，時間を要するなどの欠点もある．
図 3.10 MRI により得られた脳計測画像 [19]
35
第 3 章 波動計測
3.4
生体画像計測
生体計測の中でも画像計測は計測対象の分布が取得可能であり，画像種によっ
て，計測時間や情報量も多様である．得られる画像の種類を以下に示す．
透過画像
図 3.11 に X 線により得られた頭部透過画像例に示す [20]．計測対象に波動を照
射し，透過した波を収集し，配置することで対象の分布を画像化する．一度の照
射で情報の取得が完了するため，短時間での計測が可能である．しかし，波動の
進行方向の情報がすべて合算されることから三次元情報の提示が不可能である．
図 3.11 X 線により得られた頭部透過画像 [20]
36
3.4. 生体画像計測
再構成断層像
図 3.12 に MRI により得られた再構成断層像例を示す [21]．計測物に対して複数
の角度・位置から波動を照射して，透過，または輻射した電磁波を収集する．収
集した信号を配置，もしくはフーリエ変換を行うことによって計測対象の再構成
断層像を取得する．計測画像を取得する為には，複数回の照射が必要であり，計
測時間が長くなるが，計測断面における内部構造を取得することが可能である．
図 3.12 MRI により得られた再構成断層像 [21]
37
第 3 章 波動計測
三次元再構成像
図 3.13 に MRI により得られた三次元再構成像例を示す [22]．断層像を計測断
面の奥行き方向に複数回計測して，計測画像間を関連付けを行うことで三次元情
報を獲得する．組織間のつながりなどや立体構造などの詳細な情報を取得するこ
とが可能である．しかし，複数回の断面像計測が必要であることから，断面像計
測よりも長い計測時間が必要となる．
図 3.13 MRI により得られた三次元再構成像 [22]
上記のように透過画像，再構成断面像，三次元再構成像の順に計測時間は延長
するが，取得できる情報量も増加する．各画像種における計測時間と情報量の違
いを表 3.1 に示す．
表 3.1 画像種による計測時間と情報量の違い
38
Image Assortment
Measurement Time
Amount of Information
Transfer Imaging
Short
Little
Cross-Sectional Imaging
Medium
Moderate
3D Reconstruction Imaging
Long
Much
3.4. 生体画像計測
本研究は，血液循環システムに対して，各サブシステムの運動に対応した画像
種を取得して，解析に用いることで，計測時間や情報量の冗長性を回避して，多
階層に渡るシステムの解析を目指す．
まず心臓の収縮は特定面上だけではないため，三次元の情報が必要となるため，
三次元再構成画像を用いる．また，ねじれ収縮の詳細な解析には収縮速度が重要
であることから MR–PC 法を用いる．
次に四肢の血管解析では血管だけなく，筋組織の分布も同時に必要であるため，
再構成断層像計測を用いる．また，MRI は血液と柔組織を高コントラストで獲得
可能であることから MRI 計測を用いる．
最後に抹消血管は全身に無数に分布していることから，簡易に計測可能である
透過像計測を用いる．また，血液中のヘモグロビンは他の体組織よりも近赤外線
に対する吸光率が高いことから近赤外線を用いる．
39
第4章
MR 速度画像を用いた捻転解析
本章では，血液循環システムにおいてポンプ機能を担うサブシステムである心
臓に関する解析方法を示す．心臓は血液へ効率的な加圧を行うため、全体がねじ
れながら収縮している．複雑な心臓のねじれ収縮を詳細に解析することを目指し
た．まず、4.2 節では計測に用いた MR Phase–Contrast 法の原理について述べ、
4.3 で計測した形態画像と速度画像について述べる．次に 4.4 節で計測した速度画
像からの三次元速度場の獲得について述べる．4.5 節では心臓全体のねじれ解析
について述べ、4.6 節では心臓局所のねじれ解析について述べる．最後に、4.7 節
で本手法についてまとめる．
41
第 4 章 MR 速度画像を用いた捻転解析
4.1
はじめに
心臓は収縮・拡張によって血液への加圧を行うポンプ機能を担うサブシステム
である．また，人体の各臓器のなかでも，心臓は収縮拡張運動を絶えず繰り返す
特異な器官である [4][5]．心臓の収縮拡張運動では，様々な運動が同時に起きてお
り，非常に複雑である．
心臓の収縮拡張は，心臓を構成する心筋線維の収縮によって行われるが，心筋
線維の収縮率が 20%程度であるのに対し，全身に血液を送り出す左心室内腔の
容積変化率（駆出率）は 70%にも達することが解っている．この様な心臓全体の
効率的な収縮様式は，心臓の長軸周りに左心室を「ねじる」動きに深く関係があ
ると考えられている [11]．心臓のねじれる運動に対して、数値モデルを用いて解
析する手法 [23]–[27] や，計測画像における特徴点の経時的変化解析 [28]–[31] や、
MRI-tagging 法によるタグの経時的変化解析から心収縮を解析する手法 [32]–[35]
が行われている．特徴点やタグを追跡する形状変化解析に対して、速度情報を用
いた収縮解析は，運動を直接解析できるといった点で優れていると考えられる．
速度情報を用いた解析として、超音波ドプラ計測を用いて収縮速度を獲得し、収
縮を解析する手法 [36]–[40] などが行われている．しかし、心臓の収縮方向は，特
定の断面内だけではないため，一断面上の情報だけでは心臓の収縮の詳細な解析
には不十分である．したがって、心臓の収縮解析には三次元再構成された情報を
用いることが重要である．
MRI における Phase–Contrast 法 (以下 MR–PC 法) は計測面に関わらず任意方
向の速度が計測出来るといった利点があり，心臓の新たな解析手法として注目さ
れている [41]．本手法では，計測された三次元速度場から心臓の収縮様式にとっ
て重要な左心室のねじれ運動の解析手法を開発する．
42
4.2. MR–PC 法
4.2
MR–PC 法
本節では MR–PC 法の原理について説明する．図 4.1 に MR–PC 法におけるパ
ルスシーケンスを示す．MR–PC 法では，スライス選択後に逆勾配の二つの傾斜
磁場を用いて速度エンコードを行う．
TR
RF
Slice
Encode
Speed Encode
Gs
Phase
Encode
Frequency
Encode
FID Signal
TE
図 4.1 MR–PC 法のパルスシーケンス
1
43
第 4 章 MR 速度画像を用いた捻転解析
図 4.2，4.3 に傾斜磁場印加による水素原子の位相変化を示す．求めたい速度方
向を S とし，S 方向の傾斜磁場を速度エンコード磁場を Gs とする．各図におい
て円は各位置にある水素原子を示し、円内の矢印は位相方向を示す．
図 4.2 に S 方向に沿って移動していない水素原子のスピン位相変化を示す．
手順 1: 図 4.2(a) に一つ目の Gs の印加によるスピン位相の変化を示す．S 方向の
傾斜磁場 Gs を短時間印加する．Gs を印加後の核スピンの位相は，印加さ
れた磁場強度に依存して変化する．
手順 2: 図 4.2(b) に二つの Gs 間のスピン位相を示す．水素原子位置は変化しない．
手順 3: 図 4.2(c) に二つ目の Gs の印加によるスピン位相の変化を示す．逆勾配の
傾斜磁場 Gs の短時間印加によって，スピン位相は変化する．
手順 4: 図 4.2(d) に二つ目の Gs を印加した後のスピン位相を示す．移動していな
い水素原子のスピン位相は，二つの傾斜磁場 Gs によって，Gs を印加され
る前の状態に戻る．
Magnetic Field
Intensity
Gs
S
(a)
Phase Change
(b)
Magnetic Field
Intensity
Gs
S
(c)
Phase Change
(d)
図 4.2 移動しない水素原子に対して速度エンコード磁場を印加した場合
(a)Gs の印加による位相変化 (b) 水素原子の位置変化
(c) 逆勾配 Gs の印加による位相変化 (d) 位置不変のため、Gs による位相変化なし
44
4.2. MR–PC 法
次に，図 4.3 に速度エンコード方向に沿って移動している水素原子のスピン位
相変化を示す．
手順 1: 図 4.3(a) に一つ目の Gs の印加によるスピン位相の変化を示す．S 方向の
傾斜磁場 Gs を短時間印加する．Gs を印加後の核スピンの位相は，印加さ
れた磁場強度に依存して変化する．
手順 2: 図 4.3(b) に二つ目の Gs が印加される前のスピン位相を示す．水素原子位
置は変化する．
手順 4: 図 4.3(d) に二つ目の Gs を印加した後のスピン位相を示す．移動によって，
水素原子の核スピンの位相は Gs を印加される前の状態に戻らない．
移動量によって移動前後に受ける磁場強度差は異なるため、位相変化は速度と対
応する．Gs1 の磁場勾配の方向を速度エンコード方向と呼ぶ [42]．
Magnetic Field
Intensity
Gs1
S
(a)
Phase Change
(b)
Moving
Magnetic Field
Intensity
Gs2
S
(c)
Phase Change
(d)
図 4.3 移動している水素原子に対して速度エンコード磁場を印加した場合
(a)Gs の印加による位相変化 (b) 水素原子の位置変化
(c) 逆勾配 Gs の印加による位相変化 (d) 位置変化と Gs による位相変化
45
第 4 章 MR 速度画像を用いた捻転解析
4.3
計測画像
MRI–PC 法では、信号の強度を保存した形態画像と，信号の位相を保存した速
度画像が同時に得られる．表 4.1 と 4.2 に速度画像と形態画像の仕様を示す．計
測は Siemens 社製 1.5T MRI 装置 MAGNETIM Avanto を用いて行った．心臓短
軸断面に対して，画像の左右方向，上下方向，垂直方向の三方向に速度エンコー
ドを施し，時系列画像を取得した．心臓の位置を同じとする為に息止めを行い、
心電同期計測を行った．
表 4.1 速度画像仕様
Image Size
208 × 256 pixel
Image Resolution
1.25 mm × 1.25 mm/pixel
Temporal Resolution
72 msec
Number of Sequence
9
Number of Slices
11
Slice Thickness
5 mm
Slice Interval
5 mm
Image Depth
13 bit
Velocity Range
± 250 mm/sec
表 4.2 形態画像仕様
46
Image Size
208 × 256 pixel
Image Resolution
1.25 mm × 1.25 mm/pixel
Temporal Resolution
72 msec
Number of Sequence
9
Number of Slices
11
Slice Thickness
5 mm
Slice Interval
5 mm
Image Depth
8 bit
4.3. 計測画像
図 4.4 に計測画像と心臓の位置を示す．左心室中央部の短軸断面の収縮初期に
おける計測画像の左右方向，上下方向，垂直方向の 3 方向に速度エンコードを施
し，得られた各速度画像と各形態画像の例を示す．形態画像では輝度の高さがプ
ロトン密度の高さを示す．速度画像では、垂直方向では手前方向を、上下方向で
は上方向、左右方向では左方向を正方向とし、正方向の速度を高輝度、負方向の
速度を低輝度で示す．
Outside the Body
Left Ventricle
Right Ventricle
(a)
(d)
(b)
(e)
(f)
(c)
(g)
図 4.4 3 方向の速度分布と対応する形態画像
(a) 計測画像における心臓の位置
1
(b) 垂直方向の速度画像 (c) 垂直速度画像に対応する形態画像
(d) 上下方向の速度画像 (e) 上下速度画像に対応する形態画像
(f) 左右方向の速度画像 (g) 左右速度画像に対応する形態画像
47
第 4 章 MR 速度画像を用いた捻転解析
4.4
三次元速度場の獲得
心臓の位置は呼吸によって変化するため、計測時に息止めを行ったが，吸気量
によって心臓の位置は変化する．吸気量は同程度で行ったことから位置変化は小
さいと考えられ，移動量は 10mm 以下であると考えらる [7]．従って，吸気量の差
による影響を取り除く為に形態画像間で位置合わせを行い，速度画像を対応付け
ることで，三次元速度場を生成する．
計測で得られた画像において，位置 (o, p) における形態画像の画素値を Ik (o, p)(k =
x, y, z) とする．I の添字 k は速度エンコードの方向を示し，x は画像の左右方向，
y は画像の上下方向，z は画像の垂直方向を示す．x 方向の画像に対して、y 方向
画像と z 方向画像の位置合わせを行う．以下では，z 方向の x 方向の画像に対す
る位置合わせ手法の手順について述べる．
手順 1
z 方向の形態画像から左心室が含まれる様に画像サイズ N × N の注目領域 Hz を
抜き出す．領域 Hz の位置は用手的に決め，左上の点を (o0 , p0 ) とする．
手順 2
x 方向の形態画像から N × N の注目領域 Hxab を抜き出す．注目領域の左上の点
を (o0 + a, p0 + b)(a = −10, −9, ..., 9, 10 b = −10, −9, ..., 9, 10) とする．
手順 3
Hz と Hxab の全画素について，画素値差の絶対値の総和を計算する．
∑
| (Hz (o, p) − Hxab (o, p)) |
(4.1)
o,p
手順 4
式 (4.1) が最小となる点 (o0 + a′ , p0 + b′ ) を左上とした領域と Hz を対応付ける．
48
4.4. 三次元速度場の獲得
(o0, p0)
N
N
Magnitude Image in z-Direction
(o0 -10, p0 +10)
(o0 +10, p0 -10)
Registration
Magnitude Image in x-Direction
図 4.5 注目領域の抜きだし処理
位置合わせによる移動量
位置合わせによる b′ と a′ の最大量はそれぞれ，5pixel と 3pixel である．表 4.2
に示した形態画像の画像解像度から心臓の縦横の移動量はそれぞれ 6.25mm と
3.75mm であり、仮定通り 10mm 以下であった．
1
49
第 4 章 MR 速度画像を用いた捻転解析
4.5
回転の中心軸を用いたねじれ解析
本節では心臓全体のねじれの解析手法を示す．三次元速度場から心臓各部の回
転速度を解析し、ねじれの時間変化を求める．
4.5.1
ねじれ運動の定義
図 4.6 に回転物体の中心軸と物体の断面の例を示す．回転運動の中心軸に対し
て直交する面を物体の断面とする．
Center Axis
Cross Section
図 4.6 回転の中心軸と断面
本論文では同一断面上の 2 領域における中心軸周りの角速度が異なる現象を断
面内ねじれ，別断面上の 2 領域における中心軸周りの角速度が異なる現象を断面
内ねじれと定義する．断面内ねじれと断面間ねじれの詳細定義を以下に示す．
断面内ねじれ
断面内の 2 領域 Hin と Hout のねじれをねじれ角 θ で定義する．図 4.7 に断面内
ねじれを示す．Hin と Hout における平均角速度をそれぞれ ϕin ，ϕout とする．
1
50
4.5. 回転の中心軸を用いたねじれ解析
Center Axis
Center Axis
Hin
Hin
φ
Hout
in
φ
Hout
Angle of Twist :
θ
out
図 4.7 断面内ねじれ
このとき，単位時間後のねじれ角 θ は ϕin − ϕout となる．
断面間ねじれ
図 4.8 に断面間ねじれを示す．異なる断面の 2 領域 H0 と H1 のねじれをねじれ
率 τ を求める．H0 と H1 間の距離を l，H0 と H1 における平均角速度をそれぞれ
ϕ0 ，ϕ1 とする．また，H0 と H1 における回転軸に沿った並進速度を L0 とする．
H0
H0
L0
φ
l + L0 – L1
0
l
L1
H1
φ
1
H1
φ -φ
1
0
図 4.8 断面間ねじれ
平均角速度差 ϕ0 − ϕ1 が単位時間当たりのねじれ角となる．また，断面間の距
離 l と相対並進速度 L0 − L1 が単位時間後の断面間距離となる．従って，ねじれ
率 τ は式 (4.2) となる．
τ=
(ϕ0 − ϕ1 )
l + (L0 − L1 )
(4.2)
51
第 4 章 MR 速度画像を用いた捻転解析
4.5.2
左心室部分の抽出
図 4.9–(a) に計測によって得られた形態画像を示す．形態画像から心臓を含む
ように注目領域を抜き出した．図 4.9–(b) に抜き出した注目領域を示す．注目領
域に対して，心筋部分の境界を保ち、雑音を低減するため 3×3 の中央値フィル
タを用いた [43]．図 4.9–(c) に注目領域に中央値フィルタを用いた結果を示し，図
4.9–(d) に用手的に心臓領域を抜き出した結果を示す．図 4.9–(d) において，白色
部分が左室心筋部分，灰色部分が左室内腔，黒色部分が背景を示す．
(a)
(b)
(c)
(d)
図 4.9 計測画像と左室心筋抽出結果
(a) 計測画像 (b) 計測画像の注目領域
(c) 雑音低減後の注目領域 (d) 心筋部分と内腔部分抽出結果
52
4.5. 回転の中心軸を用いたねじれ解析
4.5.3
各断面の領域分割
図 4.10 に心筋領域の分割の概要を示す．心内膜側と心外膜側では心筋線維の走
向方向が 120deg 回転している．従って心内側と心外膜側では運動に差異がある
と考えられる．本研究では抽出した心筋領域を心内膜側と心外膜側に分割して処
理を行う．
抽出した心筋領域では心筋の厚みが 10mm であった．内腔との境界から 5mm
の範囲を心内膜側とし，外側の境界から 5mm の範囲を心外膜側とした．
Epicardium
Endocardium
Left Ventricle
図 4.10 心内膜側と心外膜側への分割
また，断面内の解像度に対して断面間の距離が大きく，断面毎の心臓形状が異
なるため，断面間の位置合わせが行えない．従って，左心室全体に対して一つの
中心軸を求めるのではなく，各断面の領域毎に回転の中心軸を決定した．
53
第 4 章 MR 速度画像を用いた捻転解析
回転速度の抽出
4.5.4
⃗
図 4.11 と図 4.12 に直交座標系と円筒座標系の関係を示す．直交座標系速度 V
⃗ と軸に沿った並進速度 L
⃗ に対する角速度 ϕ
⃗ を抽出する．
から中心軸 C
⃗ (Vx , Vy , Vz ) とし，中心
直交座標系の点 Q(x, y, z) における直交座標系速度を V
⃗ に対する回転速度を Ω
⃗ ，半径速度を R
⃗ とする．また，C
⃗ の XY 平面
軸ベクトル C
⃗ に下ろした垂線と C
⃗ の交点を Q0 とし，QXY
上の点を QXY 、点 Q(x, y, z) から C
から Q0 までの距離を l とする．
⃗ を，点 QXY を通り単位ベクトル ⃗e に平行なベクトルとする．
中心軸ベクトル C
⃗ を示す．また，図 4.12 に V⃗ と Ω
⃗ ，R
⃗ とL
⃗
図 4.11 に Q,QXY ,Q0 と l，中心軸 C
を示す．
Z
Q
Q0
Ω
L
R
X
Y
V
C
Q
QXY
⃗ と各点の関係 図 4.12 速度 V⃗ と円筒座標系速度
図 4.11 回転中心軸 C
このとき，l は式 (4.3) で示される．
−−−−→
l = (QXY Q) · ⃗e
(4.3)
−−→
また，Q0 Q は式 (4.4) で示される．
−−→
−−−−→
Q0 Q = QXY Q − l⃗e
(4.4)
⃗ に沿った並進速度成分 L
⃗ は，V⃗ と ⃗e を用いて式 (4.5) で示される．
中心軸 C
⃗ = (V⃗ · ⃗e)⃗e
L
54
(4.5)
4.5. 回転の中心軸を用いたねじれ解析
−→
−−→
⃗ は，V⃗ と −
また半径方向速度成分 R
Q0 Q の内積を | Q0 Q | で割った値で示される．
−−→
· Q0 Q
⃗ = V⃗ −
R
−→
| Q0 Q |
(4.6)
⃗ のノルム | ϕ
⃗ | は，
さらに，中心軸に対する角速度成分 ϕ
⃗ L
⃗ −R
⃗ |
⃗ |= | V −−−
|ϕ
→
| Q0 Q |
(4.7)
となる．
⃗ の回転方向は，C
⃗ における正の方向から見たとき，
また，角速度 ϕ

 Clockwise
⃗=
ϕ
 Counterclockwise
⃗ × R)
⃗ · ⃗e < 0
if (ϕ
⃗ × R)
⃗ · ⃗e > 0
if (ϕ
(4.8)
となる．
55
第 4 章 MR 速度画像を用いた捻転解析
4.5.5
中心軸の獲得
以下の手順で各領域における回転運動の中心軸を算出する．中心軸の算出は、
中心軸を様々に変更しながら領域内各点の角速度の分散を算出し、分散が最小と
なる軸を求めることで決定した．短軸断面上の左室内腔に位置する画素を Hin と
し，内腔にある全画素数を Nin とする．また，分割した領域内にある画素を Harea
とし，分割した領域内の全画素数を Narea とする．
手順 1 : Hin の 1 画素を通り，断面に対する軸の傾きが α，断面の横軸に対す
る角度が β の直線を中心軸と設定する．図 4.13 に設定した中心軸を示す．
α
Center Axis
β
Harea
Lumen of Left Ventricle : Hin
図 4.13 左心室内腔に仮定した中心軸
⃗ の平均 ϕ と分散 σ 2 を計算する．
手順 2 :手順 1 の軸に対する Harea の角速度 ϕ
ϕ=
2
σ =
1
Narea
∑−1
Narea
k=0
ϕ⃗k
(4.9)
(ϕ − ϕ⃗k )2
(4.10)
1
Narea
∑−1
Narea
k=0
手順 3 :Hin の画素において，α と β を，0 ≤ α < 90deg, 0 ≤ β < 360deg の範
囲で 5deg ずつ変化させながら中心軸を設定し，分散の計算を行う．σ 2 が最小で
ある点 Hin0 を同領域の回転運動の中心とする．
手順 4 : 手順 3 で得られた中心軸に対する ϕ を同領域の平均角速度とする．ま
⃗ の平均を同領域の軸に沿った並進速度 L とする．
た，心筋上各点の L
56
4.5. 回転の中心軸を用いたねじれ解析
4.5.6
本手法で得られた回転運動の中心軸の検証
4.5.6 では機能情報から算出された回転運動の中心と専門医によって形態情報か
ら判断された左心室の中心との比較を行う．図 4.14，4.15 に収縮初期における心
基部と心尖部における回転の中心を示す．図 4.14，4.15 において赤色の十字が心
内膜における運動の中心を示し，青色の十字は心外膜における運動の中心を示す．
図 4.14 収縮初期，心尖部の回転中心 図 4.15 収縮初期，心基部の回転中心
57
第 4 章 MR 速度画像を用いた捻転解析
専門医によって決定された中心との比較
本手法で得られた心臓の回転運動の中心軸の妥当性を検証を行う．図 4.16 に専
門医が経験を元に用手的に決定した左心室の中心を示す．図中の黒色の十字が左
心室の中心を示す．また，図 4.17 に同断面に本手法を用いて得られた中心を示す．
図 4.17 おいて赤色の十字が心内膜における運動の中心を示し，青色の十字は心外
膜における運動の中心を示す．
図 4.16 専門医による心臓の中心 図 4.17 本手法で得られた回転の中心
専門医によって決定された中心と心内膜側と心外膜で得られた中心との比較を
行う。図 4.18 に専門医による中心を原点とし，本手法で得られた中心を示す．x，
y 軸は人体における左右・前後方向を示す．赤色と青色の十字は心内膜と心外膜に
おける平均的な運動の中心位置を示す．専門医による中心と心内膜と心外膜で決
定した中心の差はそれぞれ平均 4.23mm，5.52mm であった．画像の分解能から誤
差は 5pixel 未満であり、本手法で専門医による決定と同程度の結果が得られた．
図 4.18 専門医による中心位置と本手法による中心位置の比較
中心比較から医師は心臓全体から中心を決定していると考えられる．
58
4.5. 回転の中心軸を用いたねじれ解析
断面内ねじれの時間変化と考察
4.5.7
本手法で得られた断面内ねじれを示す．まず，心内膜側と心外膜側の回転角度
を示す．次に心内膜側と心外膜側のねじれ角を示す．
図 4.19 に心内膜側と心外膜側の回転角度を示す．図において横軸は収縮初期を
時刻 0 とした時間を示す．また，縦軸は収縮初期を回転角度 0 とし，心尖部から
見て反時計回りを正とした回転角度を示す．
30
24
18
12
6
0
End Systole
-6
-12
0
100 200 300 400 500
Base Epicardium
24
Apex Endocardium
Angle (deg)
Angle (deg)
30
Apex Epicardium
Base Endocardium
18
12
6
0
End Systole
-6
600
700 800
-12
0
100 200
Time (msec)
300 400 500
600
700 800
Time (msec)
(a)
(b)
図 4.19 心内膜と心外膜の回転角の時間変化
(a) 心先部 (b) 心基部
図 4.20 に心内膜側と心外膜側のねじれ角を示す．図 4.20 において横軸は収縮
初期を時刻 0 とした時間を示す．また，縦軸は収縮初期をねじれ角 0 とし，心基
部に対する心尖部のねじれ角を示す．
1
24
Twist Angle (deg)
1
Base
Apex
18
End Systole
12
6
0
0
100 200
300 400 500
600
700 800
Time (msec)
図 4.20 心内膜側に対する心外膜側のねじれ角
59
1
第 4 章 MR 速度画像を用いた捻転解析
心臓は収縮期において心尖部から見て反時計回りに回転し，拡張期には時計回
りに回転する．また．収縮初期を基準とした収縮末期の回転角は+15deg∼–5deg
と報告されている [44]．本手法で得られた心臓の回転運動は収縮期に反時計回り
に回転し，拡張期に時計回りに回転していることから先行研究とも一致する．ま
た，収縮末期の回転角度は+25deg∼–5deg であり，先行研究と同程度の結果が得
られた．本手法は，スペックルや tagging 法による標識等を用いる手法とは異な
り，目印の有無や速度の計測方向に制限されないことから，1 心周期に渡る計測
や，心内膜側と心外膜側への分割が可能であるといった利点がある．
さらに，心内膜側と心外膜側に領域を分割することで，心内膜側と心外膜側の
ねじれ角を獲得することができた．心内膜側と心外膜側のねじれ角においては心
基部が心尖部よりも大きなねじれ角を示した．本結果は心筋線維方向が心尖部を
中心に回転し始めていることから，線維方向の差異の大きさによる影響であると
考えられる．
60
4.5. 回転の中心軸を用いたねじれ解析
4.5.8
断面間ねじれの時間変化と考察
本手法で得られた断面間ねじれを示す．まず，心尖部と心基部のねじれ角を示
す．次に心尖部と心基部のねじれ率を示す．
図 4.21 に心尖部と心基部のねじれ角を示す．図 4.21 において横軸は収縮初期
を時刻 0 とした時間を示す．また，縦軸は収縮初期をねじれ角 0 とし，心基部か
らみて反時計回りを正としたねじれ角を示す．
Twist Angle (deg)
18
Epicardium
12
Endocardium
6
0
-6
-12
-18
End Systole
0
100 200 300 400 500
600
700 800
Time (msec)
図 4.21 心基部に対する心尖部のねじれ角
図 4.22 に心尖部と心基部のねじれ率を示す．図 4.22 において横軸は収縮初期
を時刻 0 とした時間を示す．また，縦軸は収縮初期をねじれ率 0 とし，心基部に
Twist Rate (deg/mm)
対する心尖部のねじれ率を示す．
90
70
60
50
40
30
20
10
0
Epicardium
Endocardium
1
End Systole
0
100 200
300 400 500
600
700 800
Time (msec)
図 4.22 心基部に対する心尖部のねじれ率
61
第 4 章 MR 速度画像を用いた捻転解析
先行研究において心尖部と心基部のねじれ角は収縮末期で最大となり，5∼18deg
と報告されている [44]．本手法で得られた心外膜側における心尖部と心基部のね
じれ角は収縮末期で最大の 18deg となった．従って先行研究と同程度の結果が得
られた．
また，本手法によって，心尖部と心基部のねじれ率は収縮末期で最大，かつ拡張
末期で 0 となり，心外膜側が心内膜側よりも大きいねじれ率を示すことが解った．
心内膜側と心外膜側のねじれ角が拡張初期で最大値に達するのに対して，心尖
部と心基部のねじれ率は収縮末期で最大となる．さらに，心尖部と心基部のねじ
れ率において心内膜側のねじれ率が心外膜側よりも高い．これらの結果から心臓
の効率的な収縮様式にとって心内膜側における心尖部と心基部のねじれ率が関係
していると考えられる．
4.5.9
ねじれ解析に対する考察
本手法は，解析において標識等の目印を用いないことから心内膜側と心外膜側
への領域分割が可能であり，三次元速度場を用いることから運動を直接解析する
ことができる．さらに，本手法は各領域間のねじれ角やねじれ率を定量的に評価
できるといった利点を持つ．しかしながら，解析には中心軸の決定が不可欠であ
るといった問題点も残る．この問題点を解決する為に次章では中心軸の決定を必
要としないねじれ解析手法について示す．
62
4.6. ひずみの方向によるねじれ解析
4.6
ひずみの方向によるねじれ解析
本節では心臓各部に働くひずみからねじれ運動を解析する手法を示す．三次元
速度場から算出したひずみを解析することで局所ねじれの時間変化を解析する．
4.6.1
ひずみの定義
図 4.23 に面に働く垂直ひずみとせん断ひずみを示す [45]．ひずみとは対象物体
の形状変形を指し，対象物体に力が作用することで起こる．対象物体内に微小な
立方体を仮定する．微小立方体の各面に力が作用すると，図 4.23 のように作用し
た力によって各面に変形が生じる．微小立方体の各面の変形は 3 軸方向に分解で
きる．従って，微小立方体に起こるひずみは 3 平面の 3 軸方向分解によって 9 成
分で示される．各面に対して垂直に起こるひずみを垂直ひずみ，水平方向に起こ
るひずみをせん断ひずみと呼ぶ．
XY Plane
z
y
Stress
x
YZ Plane
XZ Plane
Normal Strain
Shear Strain
図 4.23 垂直ひずみとせん断ひずみ
微小直方体の 3 軸を x,y,z とし，各軸に直交する面を Y Z,XZ,XY 面とする．ま
た，Y Z 平面に起こり，y 軸方向に分解したひずみを εxy としたとき，ひずみの 9
成分は式 (4.11) に示すような 2 階のテンソルになる．この ε をひずみテンソルと
呼ぶ [46]．


 εxx εxy εxz

ε=
 εyx εyy εyz

εzx εzy εzz





(4.11)
63
1
第 4 章 MR 速度画像を用いた捻転解析
4.6.2
速度場を用いたひずみの獲得
三次元速度を単位時間当たりの各領域の変位 Uk (k = x, y, z) とみなし，変位解
析によって各領域に生じるひずみを求める．各領域の変形は微小であると考え，
ひずみテンソルには Cauchy の微小ひずみテンソルを用いる．Cauchy の微小ひず
みテンソルを式 (4.12) で示す [47]．
εkj =
1
(uk,j + uj,k ) (k = x, y, z j = x, y, z)
2
(4.12)
uk,j (k = x, y, z j = x, y, z) は式 (4.13) で示される．
uk,j =
∂Uk
(k = x, y, z j = x, y, z)
∂j
(4.13)
式 (4.12) を式 (4.11) に代入することで，ひずみテンソル ε は 3 × 3 の正方行列

ux,x


 (uy,x + ux,y )
ε=


2
 (u + u )
z,x
x,z
2
(ux,z + uz,x )
2
(uy,z + uz,y )
2
uz,z
(ux,y + uy,x )
2
uy,y
(ux,y + uy,x )
2








(4.14)
となる．
作成した正方行列のうち対角成分が垂直ひずみとなり，対角成分以外がせん断
ひずみとなる．
4.6.3
主ひずみの獲得
式 (4.14) で示したひずみテンソル ε を対角化することでせん断ひずみが 0 とな
る座標変換を行う．

 λ0 0
⃗ 0, W
⃗ 1, W
⃗ 2 ]−1 ε[W
⃗ 0, W
⃗ 1, W
⃗ 2] = 
 0 λ
[W
1


(| λ0 |>| λ1 |>| λ2 |)

0 

0 


(4.15)
0 0 λ2
⃗ 0 |=| W
⃗ 1 |=| W
⃗ 2 |= 1
|W
⃗0
このとき，Wk (k = 0, 1, 2) が変換後の三つの主ひずみの方向となる．また λ0 W
⃗ 2 を第 3 主ひずみと呼ぶ．
⃗ 1 を第 2 主ひずみ，λ2 W
を第 1 主ひずみ，λ1 W
64
4.6. ひずみの方向によるねじれ解析
4.6.4
せん断ひずみの獲得
⃗ k, W
⃗ j (k = 0, 1, 2 j = 0, 1, 2 k ̸=
得られた主ひずみからせん断ひずみを求める．W
j) による平面上のせん断ひずみの大きさ ξk,j は
ξk,j =| λk − λj |
(4.16)
⃗k と W
⃗ j のそれぞれと 45deg の角度をなす方向と，W
⃗ k と −W
⃗j
で示され，方向は W
のそれぞれと 45deg の角度をなす方向の二つとなる [46]．ξk,j のうち，最大の値
を示すせん断ひずみを最大せん断ひずみと呼び，最大せん断ひずみの 2 方向によ
る平面を最大せん断面と呼ぶ．
65
第 4 章 MR 速度画像を用いた捻転解析
4.6.5
捻転運動の表示
速度場から獲得したベクトルから捻転運動の表示を行う．図 4.24(a) に第 1 主
ひずみの方向を，図 4.24(b) に最大せん断面の法線を矢印で表現した結果を示す．
両図は収縮初期の心尖部と心基部の中間に位置する短軸断面上の解析結果を示す．
第 1 主ひずみに沿った変形量が最も大きいことから，第 1 主ひずみの方向が各領
域の主な変形方向を示していると考えられる．また，最大せん断面に対する法線
がねじれの面を示すと考えられる．
(a)
(b)
(c)
図 4.24 捻転運動の矢印表示
(a) 心臓短軸断面における表示範囲
(b) 第 1 主ひずみ方向 (c) 最大せん断面の法線
66
4.6. ひずみの方向によるねじれ解析
図 4.24 に示した矢印を用いた表現手法では捻転運動の観察が困難である．本手
法では方向を簡略化し表示を行った．得られた方向を簡略化する手法を以下に示
⃗ とする．
す．以下では簡略化するベクトルを D
⃗ と，直交座標系の各 3 軸にそれぞれ平行な三つの単位ベクトル e⃗k (k = x, y, z)
D
との内積 Ḋk (k = x, y, z) を求める．
⃗ · e⃗k (k = x, y, z)
Ḋk = D
(4.17)
⃗ の簡略化した方向とする．簡略化した方向が
Ḋk の絶対値が最大となる k を D
x 軸方向である場合は赤，y 軸方向である場合は緑，z 軸方向である場合は青で表
示する．
67
第 4 章 MR 速度画像を用いた捻転解析
4.6.6
数値モデルを用いた解析手法評価
本手法によるねじれ運動解析の評価のためにねじれを持つ数値モデルを作成し
た．複数のねじれ運動が重なった数値モデルでは計測データとの比較が複雑とな
るため、断面内ねじれと断面間ねじれのみを持つ数値モデルをそれぞれ作成した．
作成した数値モデルに対して提案手法で解析を行い、比較を行う．
断面内ねじれ数値モデル
計測データと比較を行う為に断面間ねじれを持たず，断面内ねじれを持つ数値
モデルを作成した．断面内ねじれを持つ数値モデルを図 4.25 に示す．
Z
γ
ｒ
0
ｒ
0
γ
(r, ,z)
r
図 4.25 断面内ねじれ数値モデル
二つの断面が共に r = r0 を境界に角速度が変化する数値モデルを作成した．数
⃗ の条件を以下に示す．
値モデルにおける円筒座標系の点 (r, γ, z) における角速度 ϕ

⃗  = 0(r =
̸ r0 )
∂ϕ

∂r
̸= 0(r = r0 )
⃗
∂ϕ
∂γ
⃗
∂ϕ
∂z
68
=0
=0
(4.18)
4.6. ひずみの方向によるねじれ解析
⃗ i の簡略化した方向の分布を示す．
図 4.26 に W
図 4.26 断面内ねじれ数値モデルにおける第 1 主ひずみの方向
図 4.27 に最大せん断面の法線を簡略化した方向の分布を示す．
図 4.27 断面内ねじれ数値モデルにおける最大せん断面の法線方向
角速度が変化する境界を含む領域でひずみが確認できた．さらに，第 1 主ひず
みの方向では x 軸方向と y 軸方向の領域のみが確認され，最大せん断面の法線方
向では z 軸方向の領域のみが確認された．従って，最大せん断面の法線方向を用
いることで断面内ねじれを確認できることが解った．
69
第 4 章 MR 速度画像を用いた捻転解析
断面間ねじれ数値モデル
計測データと比較を行う為に断面間ねじれを持ち，断面内ねじれを持たない数
値モデルを用いる．図 4.28 に断面間ねじれを持つ数値モデルの概念を示す．
Z
γ
γ
(r, ,z)
r
図 4.28 断面間ねじれ数値モデル
数値モデルは同じ断面内では同一角速度を持ち，断面毎に違う角速度で回転す
る．数値モデルにおける円筒座標系の点 (r, γ, z) における角速度 ϕ の条件を以下
に示す．
⃗
∂ϕ
=0
∂r
⃗
∂ϕ
=0
∂γ
⃗
∂ϕ
̸= 0
∂z
70
(4.19)
4.6. ひずみの方向によるねじれ解析
⃗ i の簡略化した方向の分布を示す．
図 4.29 に W
図 4.29 断面間ねじれ数値モデルにおける第 1 主ひずみの方向
図 4.30 に最大せん断面の法線を簡略化した方向の分布を示す．
図 4.30 断面間ねじれ数値モデルにおける最大せん断面の法線方向
数値モデルの全領域でひずみが確認できた．さらに，第 1 主ひずみの方向では
z 軸方向の領域のみが確認され，最大せん断面の法線方向では x 軸方向と y 軸方
向の領域のみが確認された．従って，第 1 主ひずみの方向を用いることで断面間
ねじれを確認できることが解った．
71
第 4 章 MR 速度画像を用いた捻転解析
4.6.7
計測データにおける捻転表示
図 4.31 に計測データにおける第 1 主ひずみの方向分布を，図 4.32 に最大せん
断面の法線分布を示し，分布に関する考察を行う．計測データを一辺 5mm の立
方体領域に分割し，解析を行う．
(a)
72
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
図 4.31 左心室における第 1 主ひずみ方向の時間変化
0msec:収縮初期 (b) 72msec
(c)144msec
(d)216msec
(e)288msec:収縮末期
(f)360msec:拡張初期
(g)432msec
(h)502msec
(i)648msec:拡張末期
4.6. ひずみの方向によるねじれ解析
(a)
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
図 4.32 左心室における最大せん断面法線の時間変化
0msec:収縮初期 (b) 72msec
(c)144msec
(d)216msec
(e)288msec:収縮末期
(f)360msec:拡張初期
(g)432msec
(h)502msec
(i)648msec:拡張末期
73
第 4 章 MR 速度画像を用いた捻転解析
図 4.33 に第 1 主ひずみの方向と最大せん断面の法線において，z 軸方向を示す
領域の時間変化を示す．図 4.33 において横軸は収縮初期を 0 とした時刻を示し，
Rate of z direction area (%)
縦軸は各時間の全領域に対する z 軸方向を示した領域の割合を示す．
70
principle strain
shaer strain
60
50
40
end systole
30
20
10
00
100
200
300
400
500
600
700
Time (msec)
図 4.33 ねじれ方向が z 軸方向を示した領域の割合時間変化
第 1 主ひずみの方向での z 軸方向の領域割合は，収縮初期と拡張初期において
高く，収縮末期に低くなった．また，最大せん断面の法線方向での z 軸方向の領
域割合は，収縮末期に高く，収縮初期と拡張末期に低かった．
これらの結果から心臓は収縮拡張運動において，初期に心尖部と心基部のねじ
れが大きくなり，その後，心内膜側と心外膜側のねじれが大きくなると考えられ
る．従って，心臓のねじれ運動の時間変化を明らかに出来たと考えられる．
4.6.8
ひずみ方向解析に対する考察
本解析手法は，速度場解析を行っていることから標識等の目印を必要としない．
また，目印等の追跡も必要ではないことから複数の時系列データを必要としない
といった利点がある．さらに，三次元速度場を用いて解析を行っていることから
運動を直接解析することができる．本手法は回転の中心軸の決定を必要としない
ことから結果が中心軸の正確さに影響されないといった利点がある．加えて，本
手法では断面内ねじれと断面間ねじれの傾向を確認することができる．
74
4.7. おわりに
4.7
おわりに
本手法の目的は血液循環システムの中枢ポンプ機能を担うサブシステムである
心臓の効率的かつ複雑な捻転収縮拡張運動の詳細を明らかにすることである．捻
転運動の詳細を明らかにする為に，MR–PC 法で計測された三次元速度場に対す
る二つの解析手法を提案した．
具体的には，MR–PC 法を用いて互いに直交する速度を計測し，得られた速度
情報の位置合わせを行うことで三次元速度場を得た．得られた三次元速度場から
左心室部分を用手的に抜き出し，解析を行った．一つ目の解析手法として，心内
膜側と心外膜側の各領域における回転運動の中心軸を用いた解析手法を提案した．
解析によって各領域間のねじれ角とねじれ率の時間変化を獲得した．解析結果か
ら心臓の収縮様式には，心内膜側における心尖部と心基部のねじれ率が関係して
いると考えられる．二つ目の解析手法として，心筋上の直方体領域に生じるひず
みを用いた解析手法を提案した．解析によって各直方体領域に起こるねじれ運動
の時間変化を獲得した．解析結果から心臓のねじれ運動は，収縮初期と拡張初期
に心尖部と心基部のねじれが心内膜側と心外膜側のねじれよりも大きく，収縮末
期には心内膜側と心外膜側のねじれが心尖部と心基部のねじれよりも大きいと考
えられる．中心軸とひずみを用いる 2 種類の解析手法による結果と考察から，心
臓の捻転運動の時間変化を明らかにした．
75
第5章
MR 形態画像を用いた収縮解析
本章では，血液循環システムにおいて血液循環の経路および加圧補助機能を担
うサブシステムである四肢血管に関する解析方法を示す．血液への加圧補助機能
に重要な血管圧迫観察と筋肉収縮解析を同時に行う手法を目指した．まず 5.1 節
では本手法の目的について述べる．つぎに 5.3 節では計測条件と計測画像につい
て述べ，5.4 節で解析の流れについて述べる．さらに，5.5 節で血管圧迫観察手法
について述べ，5.6 節で筋収縮の解析方法について述べる．最後に 5.7 節で解析結
果についての考察を行い，5.8 節において本手法をまとめる．
77
第 5 章 MR 形態画像を用いた収縮解析
5.1
はじめに
四肢の血管は主に血流の経路としての役割を担うサブシステムであるが，経路
機能のほかに，加圧補助機能も担っている．四肢の血管の加圧補助機能は，心収
縮による加圧だけでは末端の血液を心臓に戻すことが困難であることから，血液
循環にとって重要な機能である．血管自体は収縮機能を持たないため，血管周囲
の筋肉収縮による圧迫・弛緩によって，加圧補助機能を実現している．
筋収縮による血圧変化解析手法として，血管にカテーテルを挿入し，圧力変化
を直接計測する手法 [48] などがあり，加圧補助機能は広く解析されている．加圧
補助機能にとって重要な筋収縮の計測として，筋電位計測や超音波計測なども行
われているが，深部筋肉の計測には電極針を刺す必要や，計測可能な方向に制限
があるなどの問題がある．筋電位計測や超音波計測に対して，MRI 計測は，非侵
襲的かつ，計測可能な方向に対する制限もなく，非常に有用な手法である．
MRI 計測では，筋収縮による拡散値変化を利用し，拡散強調画像から解析す
る手法 [49]–[52] や筋線維の収縮速度を Phase-Contrast 法で計測・解析する手法
[53]–[55] などが行われているが，加圧補助機能に重要な血管圧迫動態の観察を同
時に行うことは困難である．
本手法では，血液への加圧補助機能にとって重要な筋収縮判定と血管圧迫動態
の観察を同時に行う解析方法を開発する．
具体的には計測画像から筋肉領域と血管領域を特定し，収縮・弛緩時の分布変
化から血管圧迫動態の観察を行う．さらに，5.2 節に示す通り，筋収縮によって信
号強度分布は変化すると考えられることから，筋肉由来の信号割合変化から筋肉
の収縮状態を解析する．
78
5.2. 筋収縮による信号強度分布変化
5.2
筋収縮による信号強度分布変化
図 5.1 に筋原線維の収縮の概念を示す．筋肉は複数の筋束から構成され，筋束
は筋原線維の集まりである．筋肉の収縮は筋原線維のアクチンフィラメントとミ
オシンフィラメントの間で滑り込みが発生することで実現している [56]．
Muscular Fibril
Contraction
Myosin Filament
Actin Filament
図 5.1 筋収縮の模式図
図 5.2 に筋線維の体積増加に伴う血管圧迫の概念を示す．筋原線維の滑り込み
によって筋線維の体積は増加し，周囲の毛細血管が圧迫される [57]．筋肉の断面
上では，筋線維の占める割合が高くなる一方で，血管は圧迫され，血管の占める
割合が低下する．
Muscular Fibril
1
Vessel
Contraction
図 5.2 筋収縮による血管圧迫の模式図
したがって，筋収縮による筋線維の体積増加によって，計測画像全体において
筋線維の画素の割合は増加し，血管の画素割合は低下すると考えられる．
79
1
第 5 章 MR 形態画像を用いた収縮解析
5.3
計測対象および計測条件
図 5.3 に Gradient-Echo 法のパルスシーケンスを示す． 本手法で用いた計測条
件，計測対象を示す．弛緩状態は手を広げて力を抜いた状態とし，収縮状態は手
を握って力を入れた状態とした．また，収縮状態を長時間維持することは困難な
ため，短時間で計測が可能な Gradient-Echo 法を用いた．計測は Hitachi Medical
社製 1.5T MRI 装置 Echelon Vega を用いて行った．
TR
RF
Slice
Encode
Phase
Encode
Frequency
Encode
FID Signal
TE
図 5.3 Gradient-Echo 法のパルスシーケンス図
計測条件を表 5.1 に示す．また，図 5.4 に弛緩・弛緩状態での計測画像，および
画像内の組織分布を示す．
1
80
5.3. 計測対象および計測条件
表 5.1 計測条件
Slice plane
Axial plane of forearm
Condition
Relax condition and Strain condition
Sequence
Gradient-Echo Sequence
Flip Angle
70
TR
30msec
TE
8msec
Image resolution
1.0mm × 1.0mm/pixel
Slice thickness
5.0mm
Image size
100 × 100 pixel
81
第 5 章 MR 形態画像を用いた収縮解析
(a)
(b)
(c)
図 5.4 計測画像
(a) 組織分布例 (b) 弛緩状態の計測画像 (c) 収縮状態の計測画像
82
5.4. 解析の流れ
5.4
解析の流れ
図 5.5 に本手法の解析の流れを示す．MRI 計測では核磁気共鳴現象を用いて計
測対象のプロトン密度分布を獲得することが可能である．空気や筋肉・血液・脂
肪ではプロトン密度が異なるため，MRI 計測画像において各体組織の信号強度に
差が生じる．本手法では，計測画像の信号強度に差に着目する．各手順の目的と
処理内容を以降の節に示す．
Calculation of Histogram
Approximation to
Mixture Gaussian Distribution
Association
Gaussian Distribution with
Tissues of the Body
Segmentation to
Flexor Muscle and
Extensor Muscle
Flexor Muscle
Extensor Muscle
Calculation of Histogram
Calculation of Histogram
Approximation to
Mixture Gaussian Distribution
Approximation to
Mixture Gaussian Distribution
Comparison of
Weighting Coefficient
Comparison of
Weighting Coefficient
図 5.5 解析の流れ
1
83
第 5 章 MR 形態画像を用いた収縮解析
5.5
信号分布を用いた血管圧迫解析
背景および各体組織由来の信号強度のばらつきはガウス分布になっており，計
測画像のヒストグラムは複数のガウス分布が重合していると考えられる．まず，
計測画像のヒストグラムを混合ガウス分布へ近似し，各ガウス分布と背景・体組
織を対応づけ，計測画像を領域分割することで，血管圧迫動態の観察を行う．
5.5.1
EM アルゴリズムを用いた混合ガウス分布近似
混合ガウス分布近似は，EM アルゴリズム [58] を用いて漸化的に求める．M 個
のガウス分布から構成される混合ガウス分布 η(I) を式 (5.1) で表す．I は輝度値，
µk ，σk は平均，標準偏差であり，重み wk は式 (5.2) の条件を満たす．
η(I) =
M
−1
∑
k=0
(
(
1 I − µk
√
exp −
2
σk
2πσk2
wk
M
−1
∑
)2 )
wk = 1, 0 ≤ wk ≤ 1
(5.1)
(5.2)
k=0
式 (5.1) で示した混合ガウス分布において，対数尤度の差が最小となるように
wk , µk , σk を更新する．更新前を hk , wk , µk , σk 更新後を h′k , wk′ , µ′k , σk′ とした場合，
EM アルゴリズムでは式 (5.3)-(5.5) が成立する．Iop は座標 (o, p) の画素の輝度値
を示す．o 軸の画素数を O，p 軸の画素数を P とする．
O−1
−1
∑ P∑
µ′k =
ηk (Iop )Iop
o=0 p=0
−1
O−1
∑ P∑
(5.3)
ηk (Iop )
o=0 p=0
O−1
−1
∑ P∑
(σk′ )2 =
ηk′ (Iop )(Iop − µ′k )2
o=0 p=0
O−1
P
∑∑
(5.4)
−1p=0 ηk′ (Iop )
o=0
wk′ =
84
−1
∑ P∑
1 O−1
η ′ (Iop )
OP o=0 p=0 k
(5.5)
5.5. 信号分布を用いた血管圧迫解析
本手法では，各ガウス分布の平均値の初期値をヒストグラムのピーク，ピーク
半値幅を標準偏差の初期値とし，重みは M の逆数を初期値とした．また，全て
の変数の更新による変更量が，0.01 未満になるまで繰り返し計算を行った．
5.5.2
ヒストグラム算出と特徴確認
計測画像から画素値のヒストグラムを求めた．図 5.6 に弛緩状態における計測
画像とヒストグラムを示す．
Frequency
0.16
0.12
0.08
0.04
0
0
100
200
300
400
Intensity
(a)
(b)
図 5.6 計測画像とヒストグラム
(a) 弛緩状態計測画像 (b) 弛緩状態ヒストグラム
取得したヒストグラムにおいて，低輝度と中輝度に明瞭なピークが見られ，高
輝度になだらかなピークが見られる．各々，低輝度のピークは画像背景に起因す
ると考えられ，高輝度のピークは骨および血管からの信号に起因すると考えられ
る．したがって中輝度のピークは，筋線維および線維周辺の毛細血管に起因する
と考えられる．
85
第 5 章 MR 形態画像を用いた収縮解析
5.5.3
獲得ヒストグラムの解析
取得したヒストグラムを混合ガウス分布で近似することで，背景および各体組
織に相当する信号分布を解析する．
図 5.6 において，低輝度，高輝度のピークは各々単一のガウス分布で近似し，中
輝度のピークは単峰性であるものの裾野が広いため二つのガウス分布へ近似した．
すなわち，式 (5.6) のような 4 個の混合ガウス分布に近似した．図 5.7 に近似結果
を示す．
(
(
1 I − µk
√
η(I) = η0 (I) + η1 (I) + η2 (I) + η3 (I) =
exp −
2
σk
2πσk2
k=0
3
∑
wk
)2 )
(5.6)
Frequency
0.16
0.12
0.08
0.04
0
0
100
200
300
400
Intensity
図 5.7 混合ガウス分布近似
ここで
µ0 ≤ µ1 ≤ µ2 ≤ µ3
(5.7)
とする．分布 η1 ，η2 における平均値の初期値は同じとし，標準偏差の初期値
は半値幅とその 3 倍とした．
86
5.5. 信号分布を用いた血管圧迫解析
5.5.4
領域分割による分布変化の観察
図 5.8 に筋組織と毛細血管に相当する輝度値範囲を示す．近収縮による血管圧
迫を観察するために，計測画像の各画素を背景および体組織へ割り当てを行う．
各画素に対して座標 (o, p) の輝度値 Iop に対して，η0 ∼η3 を求め，最大となる ηk
に座標 (o, p) は属するとした．最も低信号のガウス分布である η1 を背景とし，最
も高信号のガウス分布である η3 を骨および血管とした．ここで，筋線維および周
辺の毛細血管領域の輝度値 Iop は式 (5.8) を満たす．すなわち，図 5.8 の矢印の範
囲に存在するものとする．
η0 (Iop ) ≤ η1 (Iop ), η2 (Iop ) ≥ η3 (Iop )
(5.8)
Frequency
0.16
0.12
Muscle
and
Blood Capillary
0.08
0.04
0
0
100
200
300
400
Intensity
図 5.8 筋線維と周辺毛細血管領域の輝度値範囲
抽出した領域に，再度，混合ガウス分布近似を施す．ここでは，M = 2 とし，
標準偏差が小さいガウス分布を A，大きいのをガウス分布 B とする．混合ガウス
87
第 5 章 MR 形態画像を用いた収縮解析
分布において，筋線維起因の信号は均一度が高いと考えられるため，偏差の低い
ガウス分布 A を筋線維起因の信号，分布 B を毛細血管とした．図 5.9 に弛緩状態
と収縮状態の画像に対して，背景を黒色，骨および血管を白色，ガウス分布 A を
赤色，ガウス分布 B を灰色とした結果を示す．
(a)
(b)
図 5.9 計測画像の領域分割結果
(a) 弛緩状態, (b) 収縮状態
本実験の手を握り込む運動では，屈筋領域では筋収縮が起こり，伸筋領域では
筋線維が引き伸ばされる．したがって，屈筋領域では血管・毛細血管が圧迫され
減少し，伸筋領域では筋肉領域が減少する．また，血管は細血管，毛細血管へと
分岐することから，血管付近には比較的太い細血管・毛細血管が分布していると
考えられる．灰色破線円範囲は屈筋領域にあり，弛緩状態に比べ収縮状態では血
管領域とガウス分布 B が減少し，ガウス分布 A が増加している．破線円範囲に
対して破線長方形範囲は伸筋領域にあり，ガウス分布 A が減少し，血管領域とガ
ウス分布 B が増加している．ガウス分布 B が血管領域付近に分布していること
や，屈筋・伸筋領域での増加・減少傾向からガウス分布 B を毛細血管とし，ガウ
ス分布 A を筋肉領域とした．血管領域およびガウス分布 B の変化から血管圧迫
変化を確認することができた．
88
5.6. 信号強度を用いた筋収縮解析
5.6
信号強度を用いた筋収縮解析
手を握りこむ運動では，屈筋と伸筋で収縮・弛緩状態が異なるため，屈筋領域
と伸筋領域のそれぞれで解析を行う．各領域のヒストグラムには筋肉と毛細血管
由来のガウス分布が重合していると考えれることから二つのガウス分布からなる
混合ガウス分布へ近似する．筋肉に相当するガウス分布の割合から筋収縮の解析
を行う．
5.6.1
屈筋と伸筋への領域分割
前腕には尺骨と橈骨の二つの骨があり，二つの骨を領域分割の基準として用い
る．分布 h4 に相当する，骨，血管領域に対して縮退処理，ラベリング処理を行
い，画像中心に近い二つの領域を検出し，尺骨と橈骨とした．二つの骨領域に対
して重心計算処理を施し，二つの重心を結ぶ直線の下部を屈筋領域，上部を伸筋
領域とした．
89
第 5 章 MR 形態画像を用いた収縮解析
図 5.10 に屈筋・伸筋への分割結果を示す．
(a)
(b)
図 5.10 領域抽出結果
(a) 屈筋領域, (b) 伸筋領域
また，分割した各領域においてヒストグラムを算出し，二つの混合ガウス分布
へ近似した．図 5.11 にヒストグラムと近似結果を示す．
0.020
0.030
Frequency
Frequency
0.015
0.010
0.020
0.010
0.005
0
0
0
100
200
300
01
Intensity
(a)
100
200
300
Intensity
(b)
図 5.11 分割領域のヒストグラムと近似結果
(a) 屈筋領域, (b) 伸筋領域
90
1
5.6. 信号強度を用いた筋収縮解析
5.6.2
各筋肉領域のパラメタ比較
屈筋と伸筋の各領域の混合ガウス分布パラメータを表 5.2,5.3 に記す．なお，小
数点以下第 3 桁を四捨五入した．
表 5.2 屈筋領域の混合ガウス分布のパラメータ
Weighting coefficient wA
Gaussian Distribution A
Gaussian Distribution B
0.77
Mean value µA
323.76
Variance σA
30.42
Weighting coefficient wB
0.23
Mean value µB
300.55
Variance σB
110.42
表 5.3 伸筋領域の混合ガウス分布のパラメータ
Weighting coefficient wA
Gaussian Distribution A
Gaussian Distribution B
0.71
Mean value µA
297.81
Variance σA
33.40
Weighting coefficient wB
0.29
Mean value µB
280.11
Variance σB
104.26
筋線維に対して、毛細血管分布のばらつき大きいと考えれることから，は筋線
維起因信号のガウス分布は毛細血管起因信号のガウス分布より分散が低いとし，
ガウス分布 A を筋線維起因の信号とする．筋線維起因の信号の割合を求めるため
にガウス分布 A の重みに着目する．収縮・弛緩状態における屈筋領域と伸筋領域
の重み成分の比較を行う．重み成分を屈筋領域の値で正規化した結果，伸筋領域
の重み成分は，弛緩状態で 1.01 であり，収縮状態では 0.92 であった．したがって，
弛緩状態では各領域で，筋肉領域ガウス分布の重み成分に差はなく，収縮状態で
は屈筋領域の重み成分が伸筋領域よりも高いことが解った．
91
第 5 章 MR 形態画像を用いた収縮解析
5.6.3
解析手法の評価
図 5.12 に評価結果を示す．収縮状態 4 例，弛緩状態 4 例の計測を行い，同様の
Weight
解析を行った．図 5.12 の各印は計測結果における重み成分を示す．
0.85
Relax - Flexor
Relax - Extensor
0.8
Strain - Flexor
Strain - Extensor
0.75
0.7
0.65
0.6
1
2
3
4
Trial
図 5.12 全計測結果に対する解析結果
また，解析結果における重み成分の平均と分散を表 5.4 に示す．
表 5.4 全解析結果における重み成分の平均と分散
Relax condition
Strain condition
Flexor Muscle
0.73 ± 0.05
Extensor Muscle
0.71 ± 0.02
Flexor Muscle
0.78 ± 0.05
Extensor Muscle
0.65 ± 0.02
屈筋領域と伸筋領域の重みについて，屈筋の平均値から伸筋の平均値を引いた
値を，屈筋と伸筋の標準偏差の和で割ることによって有意差を調べた．
弛緩状態では，0.29 となり，標準偏差の 0.29 倍の区間面積は正規分布表から
0.23 であることから，有意な差はなかった．次に収縮状態にでは，1.86 となり，
標準偏差の 1.86 倍の区間面積は正規分布表から 0.94 であることから，有意差 (信
頼度 90 ％以上) があることが明らかとなった．
92
5.7. 考察
5.7
考察
本手法で用いたガウス分布の分散は空間分解能により変化する．空間分解能を
低下させた場合，分散の大小関係は変化しないと想定されるが，差が小さくなり，
二つのガウス分布の分離が困難になると考えられる．したがって，解析に必要な
空間分解能の検討が今後の課題となる．本手法で用いた短い TR では拡散係数変
化が信号強度に及ぼす影響は軽微である．したがって信号強度変化は図 5.9 に示
すように，屈筋が収縮することによって，筋線維起因信号の割合が増加したこと
に起因すると考えられる．
5.8
おわりに
計測画像のヒストグラムを混合ガウス分布へ近似することで筋肉の収縮状態の
解析を目指した．弛緩状態の筋肉と緊張状態の筋肉に対して解析を行った結果，
収縮した筋肉では混合ガウス分布のうち，分散の小さいガウス分布の重みが大き
く，引っ張られた筋肉では重みが小さいことが解った．また各ガウス分布に属す
る画素の分布から血管圧迫動態も観測可能であることが解った．従って，本手法
によって，対象の筋肉が収縮・弛緩の識別と血管圧迫観察が可能であることが明
らかとなった．筋収縮による末端血流への効果検討に繋がると考えられる．また，
本手法は信号強度の分布を用いていることから高速撮影シーケンスを用いた画像
に対しても適用可能である．したがって，本手法を用いることで，運動中の経時
的な収縮変化の解析を短時間の計測で行うことが可能になると考えられる．
93
第6章
近赤外線透過像を用いた拍動解析
本章では，血液循環システムにおいて末端のガス交換，物質交換を担うサブシ
ステムである毛細血管に関する解析方法を示す．毛細血管は全身に無数に分布し
ていることから，より簡便にガス交換に重要な血液循環の有無を解析することを
目指した．まず 6.1 節では本手法の背景および目的について述べる．つぎに，6.2
節では対象の近赤外線透過画像を計測する為に作成した計測装置と計測画像につ
いて述べる．さらに 6.3 節では計測画像からの血管抽出と時間変化解析手法と結
果について述べ，6.4 で本手法と他装置の計測値の評価結果について述べる．最
後に 6.5 節で本手法をまとめる．
95
第 6 章 近赤外線透過像を用いた拍動解析
6.1
はじめに
毛細血管は血液循環システムにおいて，ガス交換・物質交換の機能を担うサブ
システムである．毛細血管においてガス交換・物質交換を行うためには，酸素・
栄養素が含まれる血液が供給され，老廃物等を含む血液が回収されることが必要
である．したがって，血液の供給・回収を行う血液循環の有無が重要となる．血
液循環の有無を判定するには，血管に直接圧力センサを当てる手法 [59] や血流速
度の計測や，血管網の圧変化を計測する手法が適用されている．血管網の圧変化
は体表で簡易に計測可能であり，血液循環の有無のほかに，脈拍と呼ばれる心臓
の拍動周期も獲得可能である．さらに，動脈血管網の圧変化から心臓による加圧
効果の評価が可能となる．
拍動周期を計る方法には血管音で計測する聴診法と血管の拍動を計測する触診
法等 [4] があり，測定には首の頚動脈，上腕や手首，指等の動脈が計測対象とされ
ているが，首や腕等は衣服で覆われており，通常，簡易には計測出来ない．首や
腕に対し指は日常生活でも露出度の高く，気楽に計測出来る部位である．図 6.1
に脈拍計測部位である上半身と手の概要を示す．
Carotid Artery
Neck
Phalanx
Brachial Artery
Radial Artery
Ulnar Artery
Digital Artery
Elbow
Wrist
Blood Vessel
Radial Artery
Ulnar Artery
(a)
(b)
図 6.1 脈拍計測場所
(a) 上半身の模式図 (b) 手の模式図
触診法を用いた計測方法では，計測部位を拘束し，加圧する必要があり，簡便
•Phalanx
に計測を行うことが出来ない．また，計測時に過度な加圧を行うと血液循環を阻
害するといった問題がある．本手法では，近赤外線を用いて動脈血管の透過画像
を計測し，血管の太さの動的変化を解析し，脈動周期を求める．
96
1
6.2. 近赤外線透過画像計測装置
6.2
近赤外線透過画像計測装置
計測対象を撮影するためにカメラの固定が可能な計測台，近赤外線照明器を作
製した．また画像はビデオキャプチャボードを用いて計算機に保存する．図 6.2
に計測装置の概観を示す．
PC
LED
Infrared
Beam
Capture Board
CCD
Camera
Optical
Filter
Measurement
Object
図 6.2 近赤外線透過画像計測装置の概観
1
97
第 6 章 近赤外線透過像を用いた拍動解析
計測台
中指第 1 関節と第 2 関節の間を計測対象としたが，計測中の対象のブレを比較
的小さくするために手のひら全体置くことの出来る計測台を製作した．
Holding Fixture
CCD Camera
LED
(a)
(b)
図 6.3 画像計測台
(a) 作成イメージ (b) 作製した計測台
1
カメラ及び計算機接続
計測対象の撮影には CMOS よりも受光感度の高い CCD 白黒ビデオカメラ (東
京電子工業株式会社製 CS8630H) を用い，フィルタは可視光フィルタ (富士フィ
ルム製 IR76) を用いる．計算機との接続にはビデオキャプチャボード (メルコ製
CBP-AV PCI バス) を使用し，画像を計算機に保存する．
図 6.4 CCD カメラ
98
6.2. 近赤外線透過画像計測装置
近赤外線照明機
光源には，動脈血管を計測するために他の体組織に比べて血液中のヘモグロビ
ンの吸光率が高い波長 870nm の近赤外線 LED を 8 個使用した．LED は計測対象
に光が可能なかぎり均一に照射されるように，基板上に格子状ではく，近接する
三つの LED の配置を正三角形とした．また計測対象である指は人によって太さ
が異なり，必要となる光量も変化するため，可変抵抗を組込み，光の強度を計測
毎に変更出来るようにした．図 6.5 に照明機の回路図を示す．
100
Ω
LED
LED
100
Ω
100
Ω
LED
LED
100
Ω
100
Ω
LED
LED
100
Ω
100
Ω
LED
LED
100
Ω
50
Ω
5V
(a)
(b)
図 6.5 近赤外線照明機
(a) 回路図 (b) 作製した近赤外線照明機
1
99
第 6 章 近赤外線透過像を用いた拍動解析
使用機器仕様
使用した機器の仕様を表 6.1 に示す．
表 6.1 使用機器仕様
CCD Camera
Sensor Resolution
768 × 494
Visible Light Filter
Cutoff Wavelength
760nm
Near Infrared LED
Wavelength
870nm
Light Intensity
100nr/sr
Resolution
320 × 240
Frame Rate
30fps
Clock Fresuency
2.6GHz
Main Memory
248 MB
Video Capture Board
CPU
Relative Spectral Distribution
Absorptive Power (%)
100
80
60
50
40
20
Wavelength (nm)
Wavelength (nm)
(a)
(b)
図 6.6 波長特性
(a) 可視光フィルタの分光特性 (b) 近赤外線 LED の相対放射強度特性
1
100
1
6.2. 近赤外線透過画像計測装置
計測画像
作成した装置を用いて，成人男性の人差し指の第一関節と第二関節の間を対象
とし，10sec 計測を行った．計測の結果，300 の計測画像を取得した．図 6.7 に可
視光で計測した反射画像と近赤外線透過画像を，計測画像仕様を表 6.2 に示す．血
液中のヘモグロビンは近赤外線の吸光率が他の体組織よりも高いため，低輝度の
部分が動脈血管部位を示す．
(a)
(b)
図 6.7 中指計測画像例
(a) 可視光反射画像 (b) 近赤外線透過画像
表 6.2 計測画像仕様
Image Size
320 × 240 pixel
Image Resolution
0.125 mm × 0.125 mm/pixel
Temporal Resolution
33 msec
Number of Sequence
300
Images Depth
8 bit
101
第 6 章 近赤外線透過像を用いた拍動解析
6.3
透過画像の動的変化解析
図 6.8 に本手法の解析の流れを示し，各手順の目的と処理内容を以降の節に示す．
Select
First Image
Decide
ROI and
Threshold
Select m th Image
Threshold
Process
Increment
m
Expansion
Degeneracy
Process
Thinning
Process
Image
Registration
Landmark
Decision
Time Variation
Analysis
Width
Width
Measurement
Time
図 6.8 解析処理の流れ
102
6.3. 透過画像の動的変化解析
6.3.1
注目領域の抜き出し
計測対象である指は第一関節付近よりも第二関節付近のほうが太くなっている
ため，共通の処理を画像全体に施すことは困難であるため，画像全体から注目領
域 (ROI) を抜きだす．さらに注目領域の抜き出しは計算回数を減少させ，処理速
度の向上につながる．本手法では血管部分の分岐点を含み，かつ画像処理速度の
向上の為に計測画像の４分の１を注目領域として抜き出した．
(ａ)
(ｂ)
図 6.9 注目領域の抜き出し
(ａ) 計測画像 (ｂ) 注目領域を設定した画像
6.3.2
閾値処理による血管抽出
計測した画像には血管部分以外に細胞部分等も存在する．血管部分の太さ変化
を評価する為に計測画像から血管部分を抽出を行う．ヘモグロビン以外の体組織
は近赤外線の吸光率が低いため，近赤外線透過画像では，血管部分が細胞部分よ
りも輝度値が低くなる．細胞部分と血管部分の輝度値が異なることを利用し，閾
値処理によって血管部分を抽出する．
閾値処理では入力画像部分の各画素に注目し．注目画素の輝度値が閾値より低
い輝度値をもつ画素の輝度値を黒色とし，閾値より高い輝度値をもつ画素の輝度
値を白色とした．人によって指の太さが違うことから，血管部分と血管以外の部
分を別ける閾値を，一定の値を設定することは出来ない．計測対象の代表例とし
103
第 6 章 近赤外線透過像を用いた拍動解析
て 1 枚目の画像を選び，注目領域にある各画素の輝度値の平均＋ 5 を閾値とした．
E = 5+
I ′ [o, p] =
−1
O−1
∑ P∑
I[o, p]/(O · P )
(6.1)
o=0 p=0


0 I[o, p] ≤ E
(6.2)
 255 I[o, p] > E
式 (6.1) において E は閾値，I[o, p] は入力画像の座標 (o, p) における輝度値，O，
P はそれぞれ画像の o 軸方向，p 軸方向の画素数を示す．式 (6.2) において I[o, p]
と I ′ [o, p] は各々入力画像と閾値処理後の座標 (o, p) における輝度値を示す．本研
究では I ′ [o, p] = 0 の領域を血管部分，I ′ [o, p] = 255 の領域をそれ以外の部分と
した．
図 6.10 に入力画像への閾値処理の概要を示す．図 6.10(a) に対して閾値処理を
施すことにより図 6.10(b) に示すように白色部分を血管部分以外，血管部分を黒
色の部分として抽出した．
o
o
35
70
199
198
199
190
219
0
0
255
255
255
255
255
45
135
58
57
204
190
165
0
255
0
0
255
255
255
55
148
57
65
70
145
65
0
255
0
0
0
255
00
151
182
57
160
65
75
70
255
255
0
255
0
0
0
154
45
57
164
205
205
217
255
0
0
255
255
255
255
155
57
70
193
200
210
218
255
0
0
255
255
255
255
p
p
(a)
(b)
図 6.10 閾値処理の概要
(a) 入力画像 (b) 注目領域への閾値処理結果
104
6.3. 透過画像の動的変化解析
図 6.11(a) に計測画像の注目領域に対して，閾値処理を施した結果を示す．血
管の形状から，血管辺縁部の凹みは血管上部の皮膚にあるシワの影響と考えられ
る．図 6.11(b) に (a) の灰色枠内の拡大を示す．
(a)
(b)
図 6.11 計測画像の注目領域に対する閾値処理画像
血管の太さ評価を行うためには，シワの影響による凹みの影響を除外する必要
がある．皮膚のシワによって欠落した部分を補い，シワを取り除く処理を以下に
示す．
105
第 6 章 近赤外線透過像を用いた拍動解析
膨張・縮退処理による雑音除去
6.3.3
上部にシワのある血管部分は輝度値が比較的高くなるため閾値処理だけでは，
検出することが困難である．指のシワによる血管部分の未検出問題の解決のため，
膨脹・縮退処理を行い未検出部分を取り除いた．
図 6.12 に x 軸方向への膨張・縮退処理の概要を示す．
o
o
o
0
0
255
255
255
255
255
0
0
0
255
255
255
255
0
0
255
255
255
255
255
0
255
0
0
255
255
255
0
0
0
0
0
255
255
0
0
0
0
255
255
255
0
255
0
0
0
255
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
255
255
0
255
0
0
0
255
0
0
0
0
0
0
255
255
0
0
0
0
0
255
0
0
255
255
255
255
0
0
0
0
255
255
255
0
0
0
255
255
255
255
255
0
0
255
255
255
255
0
0
0
0
255
255
255
0
0
0
255
255
255
255
p
p
(a)
p
(b)
(c)
図 6.12 膨脹・縮退処理の概要
(a) 入力画像 (b) 膨脹処理結果画像 (c) 縮退処理結果画像
まず閾値処によって，対象物と背景に区別した各画素に注目する．膨張処理で
は，注目画素の x 軸方向に一定幅以内の画素が対象物である時，注目画素を対象
物とする．次に膨張処理後の各画素に注目する．縮退処理では，注目画素の o 軸
方向に一定幅以内の画素が背景である時，注目画素を背景とする．

 0
if I[o + ∆o , p] =
I ′ [o, p] =
 255 otherwise
0
(−a ≤ ∆o ≤ a)

 255 if I[o + ∆o , p] = 255
I ′ [o, p] = 
0
otherwise
(−a ≤ ∆o ≤ a)
(6.3)
(6.4)
式 (6.3) は膨脹処理の式を示し，式 (6.4) は縮退処理の式を示す．I[o, p] は入力画
像の座標 (o, p) における輝度値を示す．a は任意の画素幅，対象物の輝度値は 0，
背景の輝度値は 255 である．I ′ [o, p] は処理後の座標 (o, p) の輝度値を示す．本研
究では o 軸方向に 2 画素幅，p 軸方向に 2 画素幅の膨脹・縮退処理を行なった．
106
1
6.3. 透過画像の動的変化解析
図 6.10 の注目領域に膨脹・縮退処理を行なった画像を図 6.13 に示す．膨張・縮
退処理により，皮膚のシワによって欠落していた血管部分を補うことができた．
(a)
(b)
(c)
(d)
図 6.13 膨脹・縮退処理前後の画像
(ａ) 入力画像 (ｂ) 注目領域への膨脹・縮退処理結果画像
(c) 入力画像の拡大図，(d) 膨張・縮退結果画像の拡大図
107
第 6 章 近赤外線透過像を用いた拍動解析
6.3.4
ランドマークによる位置合わせ
血管の太さの時間変化を求めるために同じ位置の血管の太さを計測する必要が
あるが，計測時に手が微動した場合，血管の位置は変化する．無拘束計測におい
て，計測対象を静止させることは困難であるため，計測画像の位置合わせによる
補正を行った．位置合わせの方法として，まず計測画像内のランドマークを決定
し，各々の画像のランドマークを一致させることによって実施する．
ランドマークの設定方法には対象物に目印となる物体などを貼り付け，ランド
マークにする方法等があるが，本研究では対象物に直接目印を付けるのではなく，
抽出した血管部分の構造的な特長である分岐点をランドマークとした．血管部分
に細線化によって獲得した血管部分の中心線から分岐点を検出した．
また分岐点は血流の合流または分岐点であることから血管の太さの変化も分岐
点以外とくらべ大きいと考え，血管の太さ計測位置として用いる．
細線化処理
血管部分の中心線を得るために血管部分と背景の境界部分から血管部分の連結
が途切れないように血管部分を１画素幅になるまで削った．本処理は細線化と呼
ばれ，いくつかの方法があるが本研究では他の細線化アルゴリズムと比較して単
純な Hilditch の方法を用いた．Hilditch の方法による細線化とは，注目画素が以
下の６条件を全てを満たしたとき削除する [60]．細線化処理の手順を以下に示す．
処理 1： 全画素の削除マークを削除する．
処理 2： 全画素について，削除条件を調べる．全ての条件を満足する場合，削除
マークをつける．
条件 1： 注目画素は対象物の画素である
条件 2： 注目画素は境界画素である
注目画素の四近傍点のうち，少なくとも１画素が背景である場合，注
目画素は境界画素である．
条件 3： 注目画素は端点ではない
注目画素の八近傍点のうち，2 画素以上が対象物である場合，注目画
素は端点ではない．本条件により，細線化の過程で線は短縮されない．
108
6.3. 透過画像の動的変化解析
条件 4： 注目画素は孤立点ではない
8 近傍が全て背景である場合，注目画素は孤立点である．
条件 5： 削除しても連結性が維持される
条件 6： 削除時に線幅が 0 とならない
近傍画素の削除マークを確認し，注目画素の削除しても，近傍画素の
連結性が維持される場合，注目画素の削除時に線幅は 0 とならない．
処理 3： 削除マークの付いた画素を背景画素とする．同時に削除マークの付いて
いる画素数を数える．
処理 4： 削除マークが付いている画素の数が 0 になるまで以上の１∼３の処理を
繰り返し，0 になった場合処理を終了する．
以上のような細線化を行なうことで血管部分の中心線をえることが出来た．以
下に図 6.14(b) に細線化を行なった結果を示す．
(ａ)
(ｂ)
図 6.14 細線化処理前後の画像
(ａ) 入力画像 (ｂ) 注目領域に細線化処理を行なった画像
109
第 6 章 近赤外線透過像を用いた拍動解析
分岐点の決定
血管部分の細線化により，血管部分の中心線を得ることが出来た．本研究では
血管の中心線全てをランドマークにするのではなく，中心線のなかの分岐点に注
目した．以下に細線化された画像から分岐点を捜し，ランドマークを決定する方
法を示す．
細線化された画像に対し，先ず全ての分岐点をランドマーク候補点とし，候補
点が複数あった場合は各点の座標を軸毎に大きい順番に並び換え，中央の値を持
つ点をランドマークとした．並び換えには候補点の数が少ないこととアルゴリズ
ムが理解しやすいことからバブルソートを使用した．バブルソートはデータを前
から順に 2 個づつ大小を比較し，より大きいデータが前になるように並び換える
ソーティングアルゴリズムである [61]．以下に処理方法を示す．
処理１:以下の二つの条件を満たす点を候補点とする．同時に候補点の数を数える．
条件１:注目画素が線の画素である．
条件２:注目画素の八近傍に三つ以上線画素が存在する．
処理２:各候補点を軸毎に座標を大きい順に並び換えを行ない，中間の値を持つ
点をランドマークとする．
110
6.3. 透過画像の動的変化解析
位置合わせ
各画像で決定したランドマークの座標が１枚目の画像のランドマークと一致す
るように画像の黒色部分と白色部分を平行移動し位置合わせを行い，移動前の座
標が存在しない場合は背景を示す 255 を移動後の輝度値とした．位置合わせ処理
の式を以下に示す．



 0 ≤ o − ∆om < O



 Im [o − ∆om , p − ∆pm ] if
′
 0≤p−∆
Im
[o, p] = 
pm < P




255
otherwise
(6.5)
式 (6.5) において，(∆om , ∆pm ) は m 枚目の画像と 1 枚目の画像のズレを示す．
′
Im [o, p] と Im
[o, p] はそれぞれ入力画像と出力画像の座標 (o, p) における輝度値
を示し，O，P は各々入力画像の o 軸方向，p 軸方向の幅を示す．
111
第 6 章 近赤外線透過像を用いた拍動解析
6.3.5
血管太さの測定
太さの評価は位置合わせを行なった画像の y 軸に平行な一列における血管部分
の画素の個数を計測し，血管の太さとした．位置合わせを行なった各け計測画像
に対して太さ計測を行なった．図 6.15 に 10sec 計測した画像における血管の太さ
Width (pixel)
変化を示す．
Time (sec)
図 6.15 血管太さの時間変化
1
112
6.3. 透過画像の動的変化解析
6.3.6
フーリエ解析
フーリエ解析は連続したデータの周波数成分をえることができる解析方法であ
る．以下に一次元でのフーリエ変換の式を示す [62]．
∫
X(ω) =
∞
−∞
x(t) exp(−jωt)dt
(6.6)
式 (6.6) において各々，j は虚数単位，ω は角周波数，x(t) は時間 t におけるデー
タを示す．更に x(t) が離散データ y[n = 0, 1, 2∼…N − 1] であるとき式 (6.6) は式
(6.7) に書き変えられる．
(
N
−1
∑
2πkn
Y [k] =
y[n] exp −j
N
n＝0
)
(6.7)
式 (6.7) において N は離散データの数を示す．本研究では 6.2 で示した通り，10sec
で 300 の画像を計測していることから N は 300，フーリエ変換後の周波数分解能
は 0.1Hz となる．図 6.16 にフーリエ変換後の強度グラフを示す．
120
Intensity
100
1.4Hz
80
60
40
20
0
0
2
4
6
8
10
12
14
Frequency (Hz)
図 6.16 太さの時間変化をフーリエ変換した強度グラフ
図 6.16 より，直流成分以外に 1.4Hz 刻みでピークが現れている．1.4Hz にピー
クを得たことから 84 times/min の拍動を計測することが出来た．2.2 節に示した
通り，被験者である成人の標準脈拍数は 70∼80 times/min であることから，標準
脈拍数の範囲内に入っていると考えられる．
113
1
第 6 章 近赤外線透過像を用いた拍動解析
6.4
脈拍情報の定量的評価
6.3 節で示したように，連続的に計測した画像から得たデータで脈拍数を計測
することが出来た．本計測手法の定量的評価を行うために以下の実験を行なった．
被験者 12 名に対して，左手中指の第１関節と第２関節の間で本手法を用い，右
手で OMRON 製血圧計 (オシロメトリック法) を用いて計測を行う．各被験者で
5 回の試行を行い，測定結果の比較を行った．図 6.17 に評価実験の結果を示す．
図 6.17 本手法と脈拍計の測定結果比較
図 6.17 の縦軸は提案手法で計測した脈拍数であり，横軸は血圧計で計測した脈
拍数である．各印は異なる被験者を示す．両手法による脈拍数の相関係数は 0.903
となり，強い相関があった．したがって，本手法で脈拍計と同程度の脈拍数を得
ることができた．
血圧計の計測では脈拍数を 1 times/min 刻みで計測している．血圧計による計
測に対して，本手法では計測時間が 10sec であることからフーリエ解析では 0.1Hz
114
6.4. 脈拍情報の定量的評価
単位での解析しか出来ず，脈拍数は 6 times/min 刻みでしか計測出来なかった．
また閾値処理では計測の途中から血管部分の抽出が出来なくなる場合もあった．
これは計測時のブレにより LED と指の距離が変化したため全体の輝度値分布が
変化したことが原因と考えられ，閾値の決定方法を改善する必要があると考えら
れる．
115
第 6 章 近赤外線透過像を用いた拍動解析
6.5
おわりに
末梢血管の動態計測では血管透過画像の解析を行なった．まず，光源に近赤外
線を利用し，血管部分が強調された画像を連続的に計測した．つぎに，計測した
画像に対し様々な画像処理を行なうことで血管の太さ変化のデータを得ることが
出来た．太さ変化のデータにフーリエ解析を施すことで脈拍情報を得ることが可
能となり，血圧計と比較して妥当性を確認した．本手法によって，単純な末梢毛
細血管の拍動を簡便に計測することが可能となった．
116
第7章
結論
本論文の目的は機能・分布の異なる複数のサブシステムから構成されている血
液循環システムの解析を行うことである．本論文で各サブシステムに対して適用
した計測手法と解析結果を図 7.1 に示す．
Many
Blood
Capillary
1.4Hz
Vessel Width
Analysis
Near Infrared
Pulse Rate
Vessel
tati
Segmen
on
Spot
Vessel Pressure
MRI
Signa
l Inte
Analy nsity
sis
Heart
End
Systole
ngle
Twist A
s
Analysi
Contraction State
MR-PC
Global Twist
Strain
Analy
sis
Local Twist
Few
Simple
Movement
Complex
図 7.1 本論文の計測と解析手法
血液循環システムを構成するサブシステムである心臓や四肢の血管，末梢血管
は各々運動動態や分布量が異なっていることから，各サブシステムに適した計測
手法を選択した．
117
1
第 7 章 結論
具体的には，ポンプ機能を担う唯一の臓器である心臓に対しては，効率的，か
つ複雑なポンプ能力の基であるねじれ収縮運動を解析する手法を提案した．まず，
任意断面の設定方向の速度情報が取得可能な MR-PC 法を用いた複数回の計測を
心臓短軸断面に対して行い，計測間の位置合わせを行うことで心臓の三次元速度
は取得した．取得した三次元速度場から心臓各部の回転速度や回転軸沿い速度を
解析し，全体のねじれの時間変化を獲得した．さらに，三次元速度場から心臓局
所の変形を算出し，変形によって生じるひずみの方向から，局所ひずみの時間変
化を獲得した．解析で得られた心臓全体・局所のねじれの時間変化から，心臓の
捻転運動の時間変化の詳細を明らかにした．
つぎに周囲の筋肉による血液への加圧補助機能を担う四肢の血管に対しては，
加圧補助に重要な筋収縮解析と血管圧迫観察を同時に行う手法を提案した．まず，
体組織分布を高コントラストで取得可能な MRI 計測で前腕輪切り断面像を計測
した．計測画像の信号ヒストグラムから筋肉と血管等の体組織由来の信号を獲得
した．つぎに，各体組織に対応する信号から計測画像の領域分割を行い，信号の
割合を解析した．解析の結果，計測画像の領域分割による血管圧迫の観察と，割
合変化を用いた筋肉の収縮解析を同時に可能とした．
さらに，ガス交換等を行い，全身に無数に分布する末梢血管に対しては，ガス
交換に重要な血液循環の有無をより簡易に解析する手法を提案した．まず，ヘモ
グロビンは近赤外線の吸効率が高いことから，近赤外線の透過画像を計測した．
計測画像に対して信号強度を基にした閾値処理・雑音除去等を行うことで末梢血
管の太さの時間変化を獲得した．太さの時間変化をフーリエ解析し，太さ変化の
周波数成分を獲得した．太さ変化の周波数成分から脈動周期を取得することが可
能となった．
上記のように血液循環のサブシステムに対して各々最適な解析手法を提案する
ことで必要な情報を最低限の時間で解析することが可能となった．
118
本論文の手法によって各サブシステムで得られた情報を図 7.2 にまとめる．
Global Twist
Local Twist
Vessel Pressure
Contraction State
Pulse Rate
図 7.2 血液循環システムに対する獲得情報
血液循環システムにおいて心臓のねじれ収縮の詳細情報は血流システム全体に
おける中枢機能を示す．また，四肢の血管における血管圧迫と筋収縮情報は血流
システム全体における補助機能を示す．さらに抹消血管における脈動情報は血流
システム全体における末端情報を示す．各階層で機能と分布を考慮して得られた
情報を統合することで血液循環システム全体の評価に繋がることが期待される．
さらに他の生体システムに対しても本論文で示した様な各階層の機能と分布を
考慮した個別の解析手法を適用することでシステム全体の評価に繋がると考えら
れる．呼吸システムに対して本論文の提案手法を適用した例を以下に示す．
呼吸システムの各サブシステムの機能と分布，解析を図 7.3 に示す．
119
第 7 章 結論
Many
Simple
Analysis
Spot
Alveoli
Normal
Analysis
Further
Analysis
Bronchi
Costa and Diaphragm
Few
Movement
Simple
Complex
図 7.3 呼吸システムの機能と分布
呼吸システムにおいて肺内部の減圧と加圧を行う肋骨と横隔膜は体内に一つで
あり，機能は複雑であることから，詳細な解析が必要である．次に気管支は外気
を肺の各部へ誘導する経路の役割を担っており，体内に複数あり，機能は比較的
単純である．最後に肺胞は血液とのガス交換を行う役割を担っており，体内に無
1
数にあり，機能は単純であることから単純な解析で十分である．呼吸システムに
おいても各サブシステムで担う機能の複雑さや分布量は異なっていることからそ
れぞれに最適な解析手法の適用が必要である．呼吸システムの各階層に適した解
析手法によって獲得が予想される情報を図 7.4 にまとめる．
Bronchi
Information
Costa and Diaphragm
Detail Information
Alveoli
Simple Information
図 7.4 呼吸システムに対して獲得が予想される情報
120
呼吸システムにおける中枢サブシステムである肋骨・横隔膜における加圧・減
圧機能の詳細情報や，末端のサブシステムである肺胞における外気循環の簡易情
報が取得される．血液循環システムと同様に呼吸システムにおいても機能と分布
を考慮した計測・解析手法によって得られた情報は呼吸システム全体の評価に繋
がると考えれる．
また，呼吸や血液循環だけでなく，消化・吸収や姿勢制御など他の生体システ
ムも機能と分布の異なる複数のサブシステムから構成されていることから他の生
体システムに対しても本論文で提案した手法を適用するでシステム全体の評価に
繋がると考えれる．以上のように本論文で提案した階層の機能と分布に対応した
解析手法の適用が様々な生体システムの解析へ繋がることを期待する．
次に，本論で述べたような階層適合解析を臓器等のマクロレベルだけでなく，
体組織や DNA 等のミクロレベルへの適用を仮定する．血液循環システムの階層
適合解析にミクロレベルの塩基配列情報を追加した例を図 7.5 に示す．
Many
DNA
Blood
Capillary
DNA Sequence
Vessel
Spot
1.4Hz
Pulse Rate
Heart
Vessel Pressure
Contraction State
End
Systole
Global Twist
Few
Simple
Movement
Local Twist
Complex
図 7.5 階層適合解析の展開例
心収縮などのマクロレベルの情報と塩基配列といったミクロレベルの情報を統
合することで生体医工学と生命科学を融合した診断・治療が導かれることを期待
する．
121
1
第 7 章 結論
さらに，血液循環や呼吸といった生体機能の評価だけではなく，生体機能を制
御する場合も各階層の機能と分布を考慮する「階層適合解析」の概念が制御機能
の向上に貢献することを願う．
122
謝辞
本研究の機会を与えて頂き、研究の方針や内容について多大なる有益な御指導
を賜りました、大阪大学大学院基礎工学研究科大城理教授に心からの深謝の意を
表すると共に、篤く御礼申し上げます．また、多くの発表や研究室外における研
究の機会を持たせて頂きましたことは、本研究を行なう上で非常に有益でした．
重ねて御礼申し上げます．
本研究をまとめるにあたり、多大な御指導、御教示を頂きました、大阪大学大
学院基礎工学研究科和田成生教授、野村泰伸教授に深く感謝いたします．特に、
和田成生教授のひずみに関しての御助言や野村泰伸教授のシステム解析に関する
御指摘は理解を深める上で非常に有益でした．篤く御礼申し上げます．
本研究の可視化を進めるにあたり、数々の有益な御助言を頂きました、大阪大
学大学院基礎工学研究科井村誠孝准教授に深い感謝の意を表します．
本研究の解析等を進めるにあたり、多大なる工学的な御指導を頂きました、大
阪大学大学院基礎工学研究科黒田嘉宏助手に深い謝辞を申し上げます．
筋収縮の研究を進めるにあたり、研究方針に対して多大なる御意見を頂きまし
た、京都大学医学部附属病院黒田知宏准教授に篤く御礼申し上げます．
末梢血管の研究を進めるにあたり、計測・解析方法に対して数々の有益な御教
示を頂きました、大阪大学大学院基礎工学研究科増田泰助手 (現 株式会社アクセ
ンス・テクノロジー) に深く感謝します．
心臓に対する研究を進めていくうえで、心臓の収縮に関して日々数々の丁寧な
御指導を頂きました、国立循環器病センタ研究機器管理室中沢一雄氏、原口亮氏
に深い謝恩を表します．
MRI の撮影や、心臓の構造・運動の理解を深めるうえで様々な御意見を頂きま
した、国立循環器病センタ病院内藤博昭医師, 東将浩医師, 佐久間利治技師に深い
謝意を表します．
日々の研究において、常日頃から御助言を頂いた研究室先輩である勝原慎介氏、
123
謝辞
互いの研究について意見しあった森下裕太氏を始めとする研究室同期の皆様、ま
た、研究や実験に協力して頂いた亀井俊智氏を始めとした研究室後輩の皆様に感
謝いたします．最後になりましたが、著者を暖かく見守り続けてくれた家族に深
く感謝します．
124
参考文献
[1] ほくでん HP「電気が発電所からご家庭に届くまで から一部抜粋
http : //www.hepco.co.jp/atoe nve ne/energy/nattoku/powerr oute.html
[2] Z. Sparchez, P. Radu, T. Zaharia, G. Kacso, B. Diaconu, I. Grigorescu and
R. Badea: B-mode and Contrast Enhanced Ultrasound Guided Biopsy of
Portal Vein Thrombosis Value in the Diagnosis of Occult Hepatocellular
Carcinoma in Liver Cirrhosis: Medical Ultrasonography, Vol.12(4):pp.286294, 2010
[3] E. Q. Barreto, H. F. Milani, E. Araujo, K. K. Haratz, L. C. Rolo, L. M.
Nardozza, H. A. Filho and A. F. Moron: Reproducibility of Fetal Heart
Volume by 3D-Sonography using the XI VOCAL Method: Cardiovascular
Ultrasound, Vol.8(17): pp.1–5, 2010
[4] 村尾誠 監訳: BEST ＆ TAYLOR の病態生理 2 –循環器–: 廣川書店, 1988
[5] 村尾誠 監訳: BEST ＆ TAYLOR の病態生理 3 –血液とリンパ 体液と尿の排
泄–: 廣川書店, 1988
[6] (財) 厚生統計協会 選択死因別死亡状況 平成 18 年 4 月分
http://www.hws-kyokai.or.jp/html-table/table05-j.html
[7] 越智淳三 訳: 解剖学アトラス–第 3 版–: 文光堂, 1990
[8] 高山 康夫: 心臓の筋層の配置配列: Vol.7(10):pp.762–771, 2006
[9] P. M. Nielsen, I. J. Grice, B. H. Smaill and P. J. Hunter: Mathematical
Model of Geometry and Fibrous Structure of the Heart: Heart and Circulatory Physiology: Vol.260:pp.1365–1378, 1991
125
参考文献
[10] 竹 内 正 明,
中 井 博 美:
心 筋 の ね じ れ–病 気 病 態 特 性–:
心 エ コ ー,
Vol.7(10):pp.831–840, 2006
[11] G. D. Buckberg: Basic Science Review the Helix and the Heart: The Journal of Thoracic and Cardiovascular Surgery, Vol.124(5):pp.863–883, 2002
[12] ぷちナスの広場 看護学
http : //www9.plala.or.jp/sophief /
[13] NHK サイエンススペシャル 驚異の小宇宙人体 しなやかなポンプ心臓・血
管: 日本放送出版協会, 1989
[14] 千原 國宏: ME 教科書シリーズ 超音波: コロナ社. 2001
[15] 小川 力，若木守明: 光工学入門–光の基礎知識の全て–: 実教出版, 1998
[16] Y. Aizuy, T. Asakuraz, K. Oginox, T. Sugitax, Y. Suzukik and K. Masudak:
Bio-speckle Flowmetry for retinal Blood Flow Diagnostics: Bioimaging,
Vol.4:pp.254–267,1996
[17] としか内科循環器科 より抜粋
http : //www.toshikanaika.com/index.html
[18] H. Kawasoe, Y. Eguchi1, T. Mizuta, T. Yasutake, I. Ozaki, T. Shimonishi,
K. Miyazaki, T. Tamai, A. Kato, S. Kudo and K. Fujimoto: Radiofrequency
Ablation with the Real-Time Virtual Sonography System for Treating Hepatocellular Carcinoma Difficult to Detect by Ultrasonography: Journal of
Clinical Biochemistry and Nutrition, Vol.40(1):pp.66–72, 2007
[19] M. Takizawa, T. Shimoda, M. Nonaka, H. Mochizuki, S. Kawasaki, H.
Takeuchi, M. Tachibana, T. Takahashi: Parallel Imaging of Head with a
Dedicated Multi-Coil on a 0.4T Open MRI: Magnetic Resonance in Medical
Sciences, Vol.4(2):pp.95–101, 2005
[20] 財団法人大阪神経外科病院 より抜粋
http : //www.oni.or.jp/
126
参考文献
[21] 医療法人社団 金内メディカルクリニック より抜粋
http : //www.kmc.or.jp/
[22] 医療法人 SOPHIA 健康増進センター より抜粋
http : //www.is − dock.com/
[23] M. Sermesant, C. Forest, X.Pennec, H. Delingette and N. Ayache: Deformable Biomechanical Models: Application to 4D Cardiac Image Analysis: Medical Image Analysis, Vol.7:pp.475–488, 2003
[24] M. Vendelim, P. H. M. Bovendeerd, J. R. Englebrecht AND T.Arts: Optimizing Ventricular Fibers: Uniform Strain or Stress, but not ATP Consumption, Leads to High Efficiency: American Journal of Physiology - Heart
and Circulatory Physiology, Vol.283:pp.1072–1081, 2002.
[25] H. U. Klemm, R. Ventura, O. Franzen, S. Baldus, K. Mortensen, T. Risius and S. Willems: Simultaneous Mapping of Activation and Motion
Timing in the Healthy and Chronically Ischemic Heart: Heart Rhythm,
Vol.3(7):pp.781–788, 2006
[26] J. Rijcken, P. H. Bovendeerd, A. J. G. Schoofs, D. H. V. Campen, and T.
Arts: Optimization of Cardiac Fiber Orientation for Homogeneous Fiber
Strain During Ejection: Annals of Biomedical Engineering, Vol.27:pp. 289–
297, 1999
[27] L. Zhuang, H. Liu, X, Liang, H, Bao, H, Hu, and P, Shi: A Simultaneous
Framework for Recovering Three Dimensional Shape and Nonrigid Motion
from Cardiac Image Sequences: Engineering in Medicine and Biology Society, Vol.6:pp.5731-5734, 2005
[28] Y. Notomi, P. Lysyansky, R. M. Setser, T. Shiota, Z. B. Popovic, M. G.
Martin, J. A. Weaver, S.J. Oryszak, N. L. Greenberg, R. D. White and
J. D. Thomas: Measurement of Ventricular Torsion by Two-Dimensional
Ultrasound Speckle:
Journal of the American College of Cardiology,
Vol.45(12):pp.2034–2041, 2005
127
参考文献
[29] E. B. Lo, P. Shi, and H. Liu: Robust Recoverv of Volumetric Cardiac Motion with Physically Constrained HOC’ Filtering: Engineering in Medicine
and Biology Society, Vol.1:pp.17–21, 2003
[30] T. H. Valle, J. Crosby, T. Edvardsen, E. Lyseggen, B. H. Amundsen,
H. J. Smith, B. D. Rosen, J. A. C. Lima, H. Torp, H. Ihlen and O. A.
Smiseth: New Noninvasive Method for Assessment of Left Ventricular Rotation Speckle Tracking Echocardiography: Circulation, Vol.112:pp.3149–
3156, 2005
[31] M. Sahling, C. Jansen, M. Arigovindan, P. Buser, S. Marsch, M. Unser and
P. Hunziker: Multiscale Motion Mapping A Novel Computer Vision Technique for Quantitative, Objective Echocardiographic Motion Measurement
Independent of Doppler First Clinical Description and Validation: Circulation, Vol.110:pp.3093–3099, 2004
[32] I. Haber, D. N. Metaxas, T. Geva and L. Axel: Three-dimensional Systolic
Kinematics of the Right Ventricle: American journal of Physiology Heart
and Circulatory Physiology, Vol.289:pp.1826–1833, 2005
[33] J. C. Baron, K. Chow, N. S. Khoo, B. T. Esch, J. M. Scott, M. J.
Haykowsky, J. V. Tyberg, and R. B. Thompson: Measurements of Changes
in Left Ventricular Volume, Strain, and Twist During Isovolumic Relaxation using MRI: American Journal of Physiology – Heart and Circulatory
Physiology, Vol.298(6):pp.1908–18, 2010
[34] E. W. Remme, K. F. Augenstein, A. A. Young, and P. J. Hunter: Parameter
Distribution Models for Estimation of Population Based Left Ventricular
Deformation Using Sparse Fiducial Markers: IEEE Transaction on Medical
Imaging, Vol.24(3):pp.381–388, 2005
[35] E. W. B. Lo, P. Shi, and H. Liu: Robust Recoverv of Volumetric Cardiac Motion with Physically Constrained H Filtering: Proceedings of the
25’ Annual Intemational Conference of the IEEE EMBS, Vol.1:pp.814–817,
2003
128
参考文献
[36] Y. Notomi, R. M. Setser, T. Shiota, M. G. Miklovic, J. A. Weaver, Z. B.
Popovic, H. Yamada, N. L. Greenberg, R. D. White and J. D. Thomas:
Assessment of Left Ventricular Torsional Deformation by Doppler Tissue
Imaging: Validation Study With Tagged Magnetic Resonance Imaging: Circulation, Vol.111:pp.1141–1147, 2005
[37] D. Mihu, D. Diculescu, N. Costin, C. M. Mihu, L. Blaga, R. Ciortea and
A. Malu?an: Applications of Doppler Ultrasound During Labor: Medical
Ultrasonography, Vol.13(2):pp.141–149, 2011
[38] H. Maruyoshi, S. Kojima, S. Kojima, Y. Nagayoshi, Y. Horibata, K.
Kaikita, S. Sugiyama and H. Ogawa: Waveform of Ophthalmic Artery
Doppler Flow Predicts the Severity of Systemic Atherosclerosis: Circulation
Journal, Vol.74:pp.1251–1258, 2010
[39] M. Hosaka, A. Takagi, T. Takagi, K. Ashihara and N. Hagiwara: Strain
Measurements During Adenosine Triphosphate Infusion Before and After
Percutaneous Coronary Intervention: Circulation Journal, Vol.74:pp.1600–
1608, 2010
[40] M. Luckie and R. Khattar: Paradoxical Systolic and Diastolic Flow Abnormalities in Hypertrophic Cardiomyopathy with Mid-cavity Systolic Obstruction: Cardiology, Vol.18(3):pp.314–317, 2011
[41] 内藤博昭: 心臓 MRI の最前線: 循環制御, Vol.24(4):pp.338–342, 2003
[42] M. F. Walker, S. P. Souza, and C. L. Dumoulin: Quantitative Flow Measurement in Phase Contrast MR Angiography: Journal of Computer Assisted Tomography, Vol.12:pp.304–313,1988
[43] 井上誠喜, 八木伸行, 林正樹, 中須英輔, 三谷公二, 奥井誠人: C 言語で学ぶ
実践画像処理: オーム社, 2002
[44] 納冨雄一: 心筋のねじれ–評価の方法–: 心エコー, Vol.7(10):pp.812–819,
2006
129
参考文献
[45] 村上敬宣: 弾性力学: 養賢堂, 1985
[46] 小林繁夫, 近藤恭平: 工学基礎講座 7 弾性力学: 培風館, 1987
[47] 井上達雄: 弾性力学の基礎: 日刊工業新聞社, 1979
[48] P. Darren. and J. J. Michael: Skeletal Muscle Blood Flow Responses to
Hypoperfusion at Rest and During Rhythmic Exercise in Humans: Journal
of Apply: Physiology, Vol.107(2):pp.429–437, 2009
[49] O. Yanagisawa, D. Shimao, K. Maruyama and M. Nielsen: Evaluation of
Exercised or Cooled Skeletal Muscle on the Basis of Diffusion-Weighted
Magnetic Resonance Imaging: European Journal of Apply Physiology,
Vol.105:pp.723-729, 2009
[50] T. Shiraishi, T. Chikui, K. Yoshiura and Yuasa: Evaluation of T2 Values
and Apparent Diffusion Coefficient of the Masseter Muscle by Clenching:
Dento Maxillo Facial Radiology, Vol.40(1):pp.35–41, 2011
[51] O. Yanagisawa, T. Kurihara, K. Okumura and T. Fukubayashi Effects of
Strenuous Exercise with Eccentric Muscle Contraction: Physiological and
Functional Aspects of Human Skeletal Muscle: Magnetic Resonance in
Medical Sciences, Vol.9(4):pp.179–86, 2010
[52] J. F. Deux, P. Malzy, N. Paragios, G. Bassez, A, Luciani, P. Zerbib,
F. T. Roudot, A. Vignaud, H. Kobeiter and A. Rahmouni: Assessment of Calf Muscle Contraction by Diffusion Tensor Imaging: European
Radiology,Vol.18:pp.2303–2310, 2008
[53] E. Castillo, J. Lima, and D.A. Bluemke: Regional Myocardial Function.
Advances in MR Imaging and Analysis: Radiographics,Vol.23:pp.127–140,
2003
[54] R. Kinugasa, D. Shin, J. Yamauchi, C. Mishra, J. Hodgson, V. Edgerton and
S. Sinha: Phase–Contrast MRI Reveals Mechanical Behavior of Superficial
130
参考文献
and Deep Aponeuroses in Human Medial Gastrocnemius During Isometric
Contraction: Journal of Applied Physiology, Vol.105(4):pp.1312–1320, 2008
[55] N. A. Curtin, F. Lou and R. C. Woledge: Sustained Performance by Red
and White Muscle Fibres from the dogfish Scyliorhinus Canicula: The Journal of Experimental Biology, Vol.213(11):pp.1921–1929. 2010
[56] 齋藤 昭彦: 骨格筋の構造: 理学療法科学 Vol.18(1):49-53, 2003
[57] M. Nakao and S. S. Segal: Muscle Length Alters Geometry of Arterioles
and Venules in Hamster Retractor: The American Jurnal of Pysiology,
Vol.268(1-2):pp.336–344,1995
[58] C. M. Bishop: Neural Networks for Pattern Recognition. §2.6–Mixture
models–: pp.59–73, Clarendon, 1995
[59] C. M. Sauer, D. H. Tomlin, H. M. Naeini, O. Gerovichev and N. V. Thakor:
Real-Time Measurement of Blood Vessel Occlusion During Microsurgery:
Computer Aided Surgery, Vol.7:pp.364–370, 2002
[60] 井上誠喜，八木伸行，林正樹，中須英輔，三谷公二，奥井誠人: C 言語で学
ぶ実践画像処理: オーム社, 2002
[61] 津田孝夫: 岩波講座ソフトウェア科学 9 数値処理プログラミング: 岩波書
店, 1988
[62] 新井仁之: フーリエ解析学: 朝倉書店, 2003
131
研究業績
論文
1. 大城 理, 堀尾 秀之, 増田 泰:
”生体血流情報の近赤外線センシングシステム”,
電気学会論文誌 E, センサ・マイクロマシン準部門, 126 (3), (114–118), 2005
2. 黒田 嘉宏, 堀尾 秀之, 増田 泰, 黒田 知宏,
大城 理, 和田成生, 原口 亮, 中沢 一雄:
”MR Phase–contrast 画像からのひずみ算出に基づく心臓の捻転解析”,
生体医工学, 46 (2), (246–253), 2008
3. 原口 亮, 堀尾 秀之, 黒田 嘉宏, 増田 泰, 黒田 知宏,
大城 理, 内藤 博昭, 東 将浩, 中沢 一雄:
”MR 位相コントラスト法による左室心筋の局所ひずみ速度解析”,
電子情報通信学会論文誌 D, 情報・システムソサイエティ部門, 91 (7), (1818–
1828), 2008
4. 堀尾 秀之, 黒田 嘉宏, 井村 誠孝, 大城 理;
”MR 画像を用いた筋収縮統計解析”,
生体医工学, 49 (1), (207–211)，2011
技術研究報告
1. 堀尾 秀之, 原口 亮, 中沢 一雄, 内藤 博昭,
東 将浩, 佐久間 利治, 増田 泰, 大城 理:
”MR Phase–contrast 法による心室壁運動の解析 –心室壁速度場の可視化と
その評価–”,
電子情報通信学会 信学技報 MBE 106(303), (1–4), 2006
2. 堀尾 秀之, 黒田 嘉宏, 黒田 知宏, 大城 理, 和田 成生, 原口 亮, 中沢 一雄:
”MR Phase–contrast 画像からの心臓の捻転運動解析”,
電子情報通信学会 信学技報 MI 134, (337–342), 2008
133
研究業績
国際学会発表
1. H. Horio, Y. Kuroda, T. Kuroda, O. Oshiro,
S. Wada, R. Haraguchi and K. Nakazawa:
”Analysis of Cardiac Torsion from MR Phase-Contrast Image ”,
Technical Committee on Medical Imaging, (337–342), Taipei, Taiwan, 2009
2. H. Horio, Y. Kuroda, T. Kuroda, O. Oshiro,
S. Wada, R. Haraguchi and K. Nakazawa:
”Analysis of Cardiac Torsion with Strain Tensor”,
World Congress 2009 –Medical Physics and Biomedical Enginering–, (0401), Munich, German, 2010
国内学会発表
1. 堀尾 秀之, 増田 泰, 大城 理:
”近赤外線を利用した血流情報の画像計測”,
生体医工学シンポジウム 2005, (3–5), (102–106), 吹田, 日本, 2005
2. 堀尾 秀之, 増田 泰, 大城 理:
”赤外線を用いた血流情報の能動的計測”,
第 25 回日本医用画像工学会大会, (OP4–1), 京都, 日本, 2006
3. 堀尾 秀之, 黒田 嘉宏, 黒田 知宏, 大城 理,
”透磁率分布を考慮した静磁場補正”,
第 28 回日本医用画像工学会大会, (OP5–06) 名古屋, 日本, 2009
4. 堀尾 秀之, 黒田 嘉宏, 井村 誠孝, 黒田 知宏, 大城 理,
”位相画像を用いた透磁率分布計測”,
生体医工学シンポジウム 2009, (3–1–04), 千葉, 日本, 2009
5. 堀尾 秀之, 黒田 嘉宏, 井村 誠孝, 大城 理,
”MRI 画像を用いた筋収縮解析”,
第 29 回日本医用画像工学会大会, 伊勢原, 日本, 2010
6. 堀尾 秀之, 黒田 嘉宏, 井村 誠孝, 大城 理,
”MR 画像を用いた筋収縮統計解析”,
生体医工学シンポジウム 2010, 札幌, 日本, 2010
その他
1. 堀尾 秀之, 増田 泰, 大城 理:
”画像計測による血流情報の獲得”,
JOINT ミーティング, 大阪, 日本, 2005
134