チコノフ(Tikhonov)の正則化

以下の本に載っていた事例
プログラミングのための線形代数
https://www.amazon.co.jp/dp/4274065782/

ある画像をぼかし行列で変換したを、でもとに復元することを考える

clear; close all;

%N=225; %Nは奇数
X = im2double(rgb2gray(imread('lena.jpg')));
N=size(X,1);
A=zeros(N,N);
high=0.4;
middle=0.24;
low=0.05;

b=[low middle high middle low];
A(1,1:3)=[high middle low];
A(2,1:4)=[middle high middle low];
for i=3:(N-2)
  A(i,i-2:(i-2)+4)=b;
end
A(N-1,N-3:N)=[low middle high middle];
A(N,N-2:N)=[low middle high];

subplot(221);
imshow(X);
title('original');

subplot(222);
imshow(A);
title('Blur map');

subplot(223);
B=A*X;
imshow(B);
title('Blurred image');

subplot(224);
C=inv(A)*B;
imshow(C);
title('Recovered image');

つまりをにかけると、もとの画像Xが復元される。

ここで=にノイズを与え、としてみる。
として、ガウシアンノイズ N(0,0.0002)を与えてみる。

subplot(131);
B=A*X;
imshow(B);
title('Blurred image (no noise)');

subplot(132);
%Y = imnoise(A*X,'salt & pepper',0.003);
Y = imnoise(A*X,'gaussian',0,0.0002);
imshow(Y);
title('Add noise to blurred image');

subplot(133);
Y2=inv(A)*Y;
imshow(Y2);
title('Recovered image');

復元後はノイズが増幅されていることがわかる。

チコノフの正則化を行う。
普通の逆行列ではなく、でもとに戻す
（※α=0のときムーア＝ペンローズの擬似逆行列になっている）

subplot(221);
B=A*X;
imshow(B);
title('Blurred image (no noise)');

subplot(222);
%Y = imnoise(A*X,'salt & pepper',0.003);
Y = imnoise(A*X,'gaussian',0,0.0002);
imshow(Y);
title('Add noise to blurred image');

subplot(223);
Y2=inv(A)*Y;
imshow(Y2);
title('Recovered image');

subplot(224);
alpha=0.03;
A2=inv(A'*A+alpha*eye(N,N))*A';
Y3=A2*Y;
imshow(Y3);
title("Recovered image (Tikonov's regularization)");

ノイズが抑制されていることがわかる。

パラメータαの値を増やすとノイズが減るが、画像がどんどん暗くなっていく。

subplot(131);
alpha=0.003;
A2=inv(A'*A+alpha*eye(N,N))*A';
Y3=A2*Y;
imshow(Y3);
title('alpha=0.003');

subplot(132);
alpha=0.03;
A2=inv(A'*A+alpha*eye(N,N))*A';
Y3=A2*Y;
imshow(Y3);
title('alpha=0.03');

subplot(133);
alpha=0.3;
A2=inv(A'*A+alpha*eye(N,N))*A';
Y3=A2*Y;
imshow(Y3);
title('alpha=0.3');