Does Difference of Gaussian need normalization?

[hrch3n]

Description

At line 340 of blob.py, why do we need to multiply the average standard deviation to provide scale invariance? The Difference of Gaussian is already scale-invariant.

    # computing difference between two successive Gaussian blurred images
    # multiplying with average standard deviation provides scale invariance
    dog_images = [(gaussian_images[i] - gaussian_images[i + 1])
                  * np.mean(sigma_list[i]) for i in range(k)]

[hrch3n]

G(x, y, kσ ) − G(x, y, σ ) ≈ (k − 1)σ ^2∇ ^2 G.
The Difference of Gaussian is a normalized approximated version of Laplacian of Gaussian

[jni]

Hmm, I honestly don't know the answer to this. Thank you for raising the issue @hrch3n.

@scikit-image/core can someone more familiar with blob_dog comment here? The "reference" we use for the function (Wikipedia) isn't very helpful.

[hrch3n]

Thanks for the reply. The reference I checked is "Distinctive Image Features from Scale-Invariant Keypoints", page 94.

[makslevental]

@hrch3n that paper doesn't have 94 pages? page 5 here

https://people.eecs.berkeley.edu/~malik/cs294/lowe-ijcv04.pdf

i also don't believe we should be normalizing the mutual differences

one can also follow the derivation here on page 12

http://www.cs.toronto.edu/~jepson/csc420/notes/imageFeaturesIIIBinder.pdf

so in fact if anything skimage should be normalizing by 1/(sigma_ratio - 1)

Edit: this must be code that's left over from when sigmas were chosen by linearly interpolating between min and max?

In that case multiplying by sigma is correct

Page 12 here

https://www.mdpi.com/2313-433X/4/6/73/pdf

[scottstanie]

Hi, I was going to open up an issue on this, but see it was already brought up.

The current implement does get the scaling wrong:

I was testing it on this fake image of a few blobs. In the blob detector, you'd expect the LoG or DoG filter to be weak at small sigmas, then get stronger, then weak again ones the sigma is bigger than all the objects (so in this case, the peaks are at like sigma = 30 and 60).

The current implementation does this:

The bigger the sigma, the stronger the filter response, so it would detect all peaks at whatever largest scale is max_sigma

Weirdly, there seems to be 2 ways to fix it:

1. Remove the normalization. @hrch3n and @makslevental are right the the Lowe SIFT paper says it's automatically scale-normalized, so we could just change
   (gaussian_images[i] - gaussian_images[i + 1]) * np.mean(sigma_list[i])
   to
   (gaussian_images[i] - gaussian_images[i + 1])

(that's probably easiest)

2. The docstring references the wikipedia article:
   https://en.wikipedia.org/wiki/Blob_detection#The_difference_of_Gaussians_approach
   lists this formula:

which is actually different than the references it points to haha. But when I tried that, it also worked (assuming all the sigmas are in `sigma_list` below)

    gaussian_images = np.stack([gaussian_filter(image, s) for s in sigma_list])
    gaussian_diffs = -1.0 * np.diff(gaussian_images, axis=0)
    t_diffs = np.diff(sigma_list ** 2)  # 't' is variance, sigma^2, in Lindeberg
    dog_images = [
        img * sigma ** 2 / dt
        for (img, sigma, dt) in zip(gaussian_diffs, sigma_list[:-1], t_diffs)
    ]

The first option seems easier, so if one of the other's don't get to it, I could put in a pull request to change that.

[grlee77]

closed by #4482