混沌系工学特論 #3

情報科学研究科 井上純一
URL : http://chaosweb.complex.eng.hokudai.ac.jp/~j_inoue/

平成17年10月31日 第3回講義

デジタルデータの転送と復元
前回の復習※ 古文書/文献の復元等では時間経過による劣化を伝送路で表すものとする

デジタル画像劣化されたデジタル画像

確率モデルで表現

{τ }

伝送路

符号
ビット列
011000101001000001001{ξ }復元されたデジタル画像

ノイズの乗ったビット列

復号

{σ }

100111001010001001001

復号
011000110010001001001
復元されたビット列

劣化過程の確率モデル化
前回の復習
1-p+1+12元対称通信路のグラフ表現

pp-1

1-p-1

条件つき確率表現

P({τ } | {ξ }) = P(τ1 ,τ 2 ,...,τ N | ξ1 , ξ2 ,..., ξ N )劣化されたピクセル列

原画像のピクセル列

2元対称通信路の条件つき確率表現

各ピクセルの劣化過程が独立ならばP (τ i = ξi | ξi ) = e βτ 2cosh ( βτ ) = 1 − p

exp ( βττ iξi )P({τ } | {ξ }) = ∏2cosh( βτ )iP (τ i = −ξi | ξi ) = e − βτ 2cosh ( βτ ) = p

劣化過程の確率モデル

12⎛ 1− p ⎞⎟p⎝⎠βτ = log ⎜ベイズの公式と事後確率
前回の復習
尤度（ゆうど）
事後確率分布P ({σ } | {τ }) =原画像のピクセル列の推定値

事前確率分布P({τ } | {σ }) P ({σ } )∑{σ } P({τ } | {σ }) P ({σ })ベイズの公式

尤度は劣化過程の確率モデル
2元対称通信路であることが既知ならば未知な変数を含むexp ( h τ iσ i )P ({τ } | {σ } ) = ∏2cosh ( h )i
事前分布の選び方は全く任意である

画像データの特徴を生かして適切に選択する必要がある画像データの事前分布と事後確率
前回の復習

強磁性体の確率モデル

P ({σ } ) =(exp J ∑ <ij > σ iσ j

∑σ exp ( J ∑<ij >)σ iσ j)スナップショット

T=1/J=1.15※ どの温度のスナップショットの張り合わせで表現するか？
⇒ どのように J を決めるか？
（この方法は後に説明する）

P ({σ } | {τ } ) =T=1/J=2.15

T=1/J=3.15

ベイズの公式から求まる事後確率

exp(− H eff ), H eff = − J ∑ <ij > σ iσ j − h∑ i τ iσ i∑{σ } exp ( − H eff )画像復元問題のエネルギー関数

スナップショットの作り方
マルコフ連鎖モンテカルロ法（MCMC法） 前回の復習

E = − J ∑ <ij > σ iσ j エネルギーが増加

ノイズを利用したアルゴリズムの多変数版(1) 各時刻で任意にσi

この前後の状態をするプロセスもある確率で許すを1つを選び、その符号を変える。

σ ,σ ′ と書く。

（2） (1)の前後でのエネルギー差 :ΔE = E (σ ′) − E (σ ) を計算し、 ΔE

ΔE < 0ならば新しい状態を採用

ΔE > 0でも確率

e −ΔE / T で新しい状態を採用(3) (1)(2)を繰り返す

P ({σ }) の時間発展を調べてみる
P ({σ } ) =平衡状態では(exp J ∑ <ij > σ iσ j∑σ exp ( J ∑<ij >)σ iσ j) 1/J=1.15

最大事後確率推定
前回の復習
事後確率

エネルギー関数P ({σ } | {τ } ) =exp(− H eff ), H eff = − J ∑ <ij > σ iσ j − h∑ i τ iσ i∑{σ } exp ( − H eff )エネルギー関数最小⇒ シミュレーテッド・アニーリングで最適化

事後確率最大画像復元が最適化問題として定式化される
このパラメータTを温度としたアニーリング

劣化画像データは与えられている{})(P ({σ } ) =∑ exp (T { J ∑ σ σ + h∑ τ σ })exp T −1 J ∑ <ij > σ iσ j + h∑ i τ iσ i−1σ<ij >iji ii※ J および h の決め方は次回。今回は既知とする

アニーリング法による画像復元の実例
前回の復習
エネルギー誤り率

シミュレーテッド・アニーリング
によるエネルギー関数の時間変化

シミュレーテッド・アニーリング
によるビット誤り率の時間変化

アニーリング・スケジュール

T (t ) ∼ 3
モンテカルロ・ステップ

原画像

J,h は真の値に固定

t
モンテカルロ・ステップ

劣化画像

復元画像

確率モデルのパラメータ推定
事後確率は次のように書けるが・・・{})(P ({σ } ) =∑ exp (T { J ∑ σ σ + h∑ τ σ })exp T −1 J ∑ <ij > σ iσ j + h∑ i τ iσ i−1σ<ij >+1-1

どの J で指定されるスナップショットの張り合わせで画像を表現するのか

未知

この決め方を今回説明する

iji ii1-pp pP p+1-11-p

劣化率 h がどのくらいなのか

用いることのできるのは劣化画像についての情報のみ

未知ハイパーパラメータの推定
ビット誤り率1⎡ 1⎤
ˆPb ( J , h ) = ⎢1 − ∑ i ξiξi ⎥2⎣ N⎦
ビット誤り率を最小にするパラメータを選べばよい？

この「エネルギー関数」を使うには、原画像についての情報が必要

観測データのみを含むコスト関数を選ぶべきである

周辺尤度最大化周辺尤度

L ( J , h;{τ } ) = log tr{σ } Ph ({τ } | {σ } ) PJ ({σ } )観測データ= log tr{σ }eJ∑ <ij>σ iσ j +h ∑ iτ iσ i+ log tr{σ }eJ∑ <ij>σ iσ j− log 2cosh(h)一般的にKL情報量の非負性より⎡⎣ L ( J * , h* ;{τ } ) ⎤⎦ ≥ ⎡⎣ L ( J , h;{τ } ) ⎤⎦{τ }{τ }真のパラメータ値で最大周辺尤度をコスト関数とみなしてパラメータに関して最大化すればよい

EM法による周辺尤度の最大化
周辺尤度を最大化する代わりに、次の尤度関数の事後確率での期待値を最大化

Q ( J , h | J t , ht ) = tr{σ } PJt ,ht ({σ } | {τ } ) log Ph ({τ } | {σ } ) PJ ({σ } )
尤度関数

事後確率
EMアルゴリズムJ t +1 = arg max Q ( J , h | J t , ht )Eステップ（期待値計算）とMステップ（最大化）の繰り返しで構成される

J

ht +1＝arg max Q ( J , h | J t , ht )h一方の変数hが既知の場合

um ( J t +1 ) =tr{σ }ε Be
tr{σ }eJ t +1J t +1

∑ <ij>σ iσ j

∑ <ij>σ iσ j

厳密解あり

∑ <ij>σ iσ j +h* ∑ iτ iσ i厳密解ないのでMCMC法で計算
tr{σ }ε Be== u p ( J t , h* )J t ∑ σ iσ j + h* ∑ τ iσ ii
tr{σ }e <ij>Jt

ε B = −(1/ 2 N )∑ <ij > σ iσ j

EM法における期待値と最大化
um ( J t +1 ) =tr{σ }ε Be
tr{σ }eJ t +1J t +1∑ <ij>σ iσ j

∑ <ij>σ iσ j=tr{σ }ε Be
tr{σ }eJt

Jt∑ <ij>σ iσ j +h* ∑ iτ iσ i∑ <ij>σ iσ j +h* ∑ iτ iσ i

= u p ( J t , h* )厳密解がないのでMCMC法で計算交点が次ステップのパラメータの値⇒モンテカルロステップに注意（緩和の問題）
具体的に原画像 ：☆ 磁石モデルのスナップショット☆ 「井上」の印鑑画像

厳密解に関する適用例を見て行く実行結果（原画像：イジングモデル）

固定点

h の時間変化
初期値

J の時間変化

ハイパーパラメータの時間発展相空間での軌道実行結果（原画像：「井上」の印鑑画像）

固定点

h の時間変化J の時間変化
初期値

ハイパーパラメータの時間発展相空間での軌道熱揺らぎを用いた画像の復元
周辺分布

P (σ i | τ )による各画素ごとの期待値

MPM推定値^

ξi = sgn(< σ i >) = sgn[∑ σ =±1σ i P(σ i | τ )]画素 jの推定値i

統計力学的な観点からは、エネルギー関数Hで記述される磁性体の局所磁化
原画像と復元画像の重なり

σ exp[− H (σ | τ )]∑σ< σ >=
=miiZiの符号を推定値とする

各画素の期待値

MAP⇒ ノイズレベルが大きい場合MAP推定より高精度で復元できる

新聞など画像の印刷

微視的にみると新聞の印刷は白と黒の2値化画像となっている。
白と黒だけでどのように「濃淡」を実現しているのであろうか？

画像印刷技術 : ハーフトーン処理
ハーフトーン処理 : 1画素をさらに細かなセルに分割し、このセルを白と黒で塗り分けたときに、黒で塗りつぶされた数で元の濃淡画素の階調値を表現する

どのように黒を分布させれば良いの？
（http://www.nime.ac.jp/~inaba/LR/IP/3-5.html より）

このハーフトーン処理を統計力学を用いて定式化、性能評価できないものか？
⇒ ハーフトーン処理のある側面に着目して問題設定を行う

閾値を用いたセル分割
Q 階調の濃淡画像を次のように、(Q-1) 枚の白黒画像に分割する⎧1(ξi ≥ k )Θ(ξi − k ) = ξi ,k ⎨
ξi = ξi,1 + ξi,2 +⎩0(ξi < k )Q−1+ ξi,Q−1 = ∑ξi,kk =1

3値濃淡画像の分解例

2階調
画像2枚3階調画像

標準画像（レナ）での2値分解例

8階調濃淡画像
各ピクセル毎の1の値の数で濃淡値を表現⇒ ハーフトーン処理になっている

7枚の2値化画像

性能評価のための問題設定
画像A

画像A

印刷したい Q値濃淡画像 A ⇒ ハーフトーン処理 ⇒ 劣化 ⇒ 画像修復
& 逆ハーフトーン処理 ⇒ Q値濃淡画像 B
画像Aと画像Bを比較することで、ハーフトーン処理を施した画像に対するインクのにじみ（ノイズ）等による劣化の耐性を評価する

Q状態イジングモデル
濃淡画像の確率モデル
P(ξ ) = exp[−( β s / 2 N )∑ ij (ξi − ξ j ) 2 ] / Z s
磁化の温度依存性（MCMCシミュレーションによる）絶対ゼロ度では4状態がエネルギー的に縮退

Q=4、T=0.35 におけるMCMCのスナッショット

劣化過程の確率モデル
ハーフトーン処理済み画像の劣化（各ハーフトーン・スライスk毎に）1Pk ({τ i ,k }|{ξi ,k }) =exp[−2πτ2(τ−τξ)∑i i,k 0 i,k

ハーフトーン未処理画像の劣化

1
P({τ i }|{ξi }) =exp[−'2πτ2τ22(τ−ξ)∑i i i2τ

'2]]

熱揺らぎを用いた画像復元を行い、得られる画像と原画像の2乗誤差を基準としてその性能を評価する

MCMCシミュレーション(Q=4)
システムサイズ１00ｘ100 約15000モンテカルロステップ

ハミング距離

原画像

劣化画像

復元画像

エネルギー関数βH eff = − h∑ i ,k (σ i ,k − τ i ,k ) 2 − ( m2N)∑ ij (σ i − σ j ) 2

真のモデルパラメータの値で最良

平均場アニーリング法（詳細は講義ノート#2 pp.26-27を参照）

局所磁化に関する反復式を解く

mijt +1 =ε(t )tσexp(H∑σ ij
MF / T )Z MF

1 N= ∑ | mij(t +1) − mij( t ) |< 10−8N ( ij )画素数オーダの計算量
t
H MF
=βm

(mit, j +1 + mit, j −1 + mit+1, j + mit−1, j )2−βmσ ij2 − h∑k (σ ij ,( k ) )2 + 2h∑k τ ij ,( k )σ ij ,( k )← 各ステップでこの条件を満たしたら状態更新

T = T − ΔT , Tint = 2, ΔT = 0.01
温度をゼロまで落とせばMAP解が得られる← 温度の制御平均場アニーリング法
標準画像 (chair)

BDD 劣化BDD 復元「しみ」がクラスタとして残る

原画像

Q-Ising 劣化

Q-Ising 復元

平均場アニーリング法
標準画像 (house)

BDD 劣化BDD 復元原画像

Q-Ising 劣化

Q-Ising 復元

画像の種類による性能比較※ より現実的な状況設定での解析は和歌山工業高等専門学校 雑賀洋平先生と共同研究を実施中

NNP1とは最隣接する画素が同じ画素値をとるものの全画像に占める割合でありO(原画像)、D(劣化画像)、R(修復画像)である

解析的性能評価（この手の解析評価の詳細は次回）両者でフェアな比較を行うため、原画像と劣化画像間の2乗誤差がQ-Ising, BDD（逆ハーフトーン）で同じになるようにノイズをかけてある

ハーフトーン処理画像の劣化からの濃淡画像修復 (BDD)
最良値どうしを比較するとBDDの方が高い性能を示す
劣化への耐性という観点からはハーフトーン処理は好ましい