imgaussfilt asymmetry as a linear operator

12 views (last 30 days)
Stephen Vavasis
Stephen Vavasis on 30 Nov 2021
Edited: Matt J on 30 Nov 2021
The Gaussian image blurring operation imgaussfilt is not a symmetric linear operator, as the test below shows. This is surprising to me based on my understanding of Gaussian blurring. I am wondering (a) why it is unsymmetric, and (b) given that it is unsymmetric, how I can obtain its adjoint operator. In other words, I would like a function imgaussfilt_adjoint such that sum(sum(X.*imgaussfilt(Y,k)))==sum(sum(imgaussfilt_adjoint(X,k).*Y)). Note that I am not referring to the 'symmetry' boundary option.
>> X = randn(200,200);
>> Y = randn(200,200);
>> Xblur = imgaussfilt(X,100);
>> Yblur = imgaussfilt(Y,100);
>> sum(sum(X.*Yblur))
ans =
3.6114
>> sum(sum(Xblur.*Y))
ans =
4.5137

Accepted Answer

Matt J
Matt J on 30 Nov 2021
Edited: Matt J on 30 Nov 2021
It is because of edge effects. If you add sufficient zero padding you will see symmetric behavior.
  3 Comments
Matt J
Matt J on 30 Nov 2021
It also appears that symmetry is present with any of the padding modes except for 'replicate'.
X = randn(200,200);
Y = randn(200,200);
for p=["symmetric","circular","replicate"]
Xblur = imgaussfilt(X,200,'Padding',p);
Yblur = imgaussfilt(Y,200,'Padding',p);
asymmetry=sum(X.*Yblur, 'all')-sum(Xblur.*Y , 'all')
disp ' '
end
asymmetry = -1.9429e-15
asymmetry = -3.3654e-15
asymmetry = 62.3654

Sign in to comment.

More Answers (2)

Image Analyst
Image Analyst on 30 Nov 2021
It's because they're random numbers. X and Y are not equal, and neither are the blurred versions. So why would you expect two different random matrices multiplied by each other element by element to have the identical integrated gray value? Look, I don't even get the same numbers as you
X = randn(200,200);
Y = randn(200,200);
Xblur = imgaussfilt(X,100);
Yblur = imgaussfilt(Y,100);
sum(sum(X.*Yblur))
ans = 24.7835
sum(sum(Xblur.*Y))
ans = -12.0049
  1 Comment
Stephen Vavasis
Stephen Vavasis on 30 Nov 2021
My original question, written in terms of linear algebra, is the following: if A is a symmetric matrix, then the equation x'*(A*y) == y'*(A*x) holds for any choices (including random) of x and y. The above test showed that this equality did not hold when A is the linear operator corresponding to imgaussfilt. However, another poster in this thread pointed out that I need to add the padding=0 option in order to attain symmetry of the operator.

Sign in to comment.


Matt J
Matt J on 30 Nov 2021
Edited: Matt J on 30 Nov 2021
(b) given that it is unsymmetric, how I can obtain its adjoint operator.
In the asymmetric case, a brute force solution would be to use func2mat to get the matrix form of the operator,
A=func2mat( @(X) imgaussfilt(X,10), ones(200)).';

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!