Entropy of an image What is the most information/physics-theoretical correct way to compute the entropy of an image? I don't care about computational efficiency right now - I want it theoretically as correct as possible.
Lets start with a gray-scale image. One intuitive approach is to consider the image as a bag of pixels and compute
$$
H = - \sum_k p_k log_2(p_k)
$$
where $K$ is the number of gray levels and $p_k$ is the probability associated with gray level $k$.
There are two problems with this definition:


*

*It works for one band (i.e. gray-scale), but how should one extend it in a statistically correct way to multiple bands? For example, for 2 bands, should one base oneself on $(X_1,X_2)$ and thus on PMF using $P(X_1=x_1,X_2=x_2)$? If one has many ($B$>>2) bands then $P(X_1=x_1, ..., X_B=x_B) \sim 1/N^B \rightarrow H_{MAX}$, which seems wrong.

*Spatial information is not taken into account. For example, the images below (custody of John Loomis) have the same $H$, although clearly they do not convey the same information. 



Anyone care to explain or give advice, or refer me to some decent reference material about the subject? I am mainly interested in a theoretically correct approach of the second problem (i.e. spatial information).
 A: There is none, it all depends on the context and your prior information. Entropy has many interpretations such as "measurement of order" or "measurement of information", but instead of looking at the interpretations you could just look at what it actually is. Entropy is just a way of expressing the number of states of a system. A system with many states has a high entropy, and a system with few states has a low entropy.
You, and the article you link to - states that the two images have the same entropy. This is not correct (for me).
The article correctly calculates the entropy is.
$$
H = - \sum_k p_k log_2(p_k)
$$
For the first image any pixel can have any gray value, $$p_k = \frac{1}{M} = 2^{-n}$$
Therefore the entropy is:
$$
H = - \sum_k p_k log_2(p_k) = - \sum_k 2^{-n} log_2(2^{-n}) = - log_2(2^{-n}) = n
$$
However, This is not the case for the second image.
The entropy can still be calculated as:
$$
H = - \sum_k p_k log_2(p_k)
$$
but you can not simple say $p_k = \frac{1}{M} = 2^{-n}$, because when you have found $p_1$ to be a value, you know that $p_2, p_3, p_4 \ldots p_{many}$ is the same value.
Therefore, the two images do not have the same entropy.
It might sound counter intuitive that entropy depends on how you look at the problem. However, you probably know it from compression. The maximum compression of a file is dictated by the Shannon's source coding theorem which sets an upper limit for how well a compression algorithm can compress a file. This limit depends on the entropy of the file. All modern compressors will compress a file close to this limit.
However, if you know the file is an audio file you can compress it using FLAC instead of some generic compressor. FLAC is lossless so all information is preserved. FLAC can not get around the Shannon's source coding theorem, that's math, but it can look at the file in a way which reduces the entropy of the file, thus do a better compression.
Identically, when I look at you second image I see that the pixels are sorted by gray value, and therefore it doesn't have the same entropy to me as the image with random noise.
A: “What is the most information/physics-theoretical correct way to compute the entropy of an image?“
An excellent and timely question.
Contrary to popular belief, it is indeed possible to define an intuitively (and theoretically) natural information-entropy for an image.
Consider the following figure:

We can see that the differential image has a more compact histogram, therefore its Shannon information-entropy is lower.  So we can get lower redundancy by using second order Shannon entropy (i.e. entropy derived from differential data).  If we can extend this idea isotropically into 2D, then we might expect good estimates for image information-entropy.
A two dimensional histogram of gradients allows the 2D extension.  
We can formalise the arguments and, indeed, this has been completed recently.  Recapping briefly:
The observation that the simple definition (see for example MATLAB’s definition of image entropy) ignores spatial structure is crucial. To understand what is going on it is worth returning to the 1D case briefly. It has been long known that using the histogram of a signal to compute its Shannon information/entropy ignores the temporal or spatial structure and gives a poor estimate of the signal’s inherent compressibility or redundancy. The solution was already available in Shannon’s classic text; use the second order properties of the signal, i.e. transition probabilities. The observation in 1971 (Rice & Plaunt) that the best predictor of a pixel value in a raster scan is the value of the preceding pixel immediately leads to a differential predictor and a second order Shannon entropy that aligns with simple compression ideas such as run length encoding. These ideas were refined in the late 80s resulting in some classic lossless image (differential) coding techniques that are still in use (PNG, lossless JPG, GIF, lossless JPG2000) whilst wavelets and DCTs are only used for lossy encoding.
Moving now to 2D; researchers found it very hard to extend Shannon’s ideas to higher dimensions without introducing an orientation dependence. Intuitively we might expect the Shannon information-entropy of an image to be independent of its orientation. We also expect images with complicated spatial structure (like the questioner’s random noise example) to have higher information-entropy than images with simple spatial structure (like the questioner’s smooth gray-scale example). It turns out that the reason it was so hard to extend Shannon’s ideas from 1D to 2D is that there is a (one-sided) asymmetry in Shannon’s original formulation that prevents a symmetric (isotropic) formulation in 2D. Once the 1D asymmetry is corrected the 2D extension can proceed easily and naturally.
Cutting to the chase (interested readers can check out the detailed exposition in the arXiv preprint at https://arxiv.org/abs/1609.01117 ) where the image entropy is computed from a 2D histogram of gradients (gradient probability density function).
First the 2D pdf is computed by binning estimates of the images x and y derivatives. This resembles the binning operation used to generate the more common intensity histogram in 1D. The derivatives can be estimated by 2-pixel finite differences computed in the horizontal and vertical directions. For an NxN square image f(x,y) we compute NxN values of partial derivative fx and NxN values of fy.  We scan through the differential image and for every pixel we use (fx,fy) to locate a discrete bin in the destination (2D pdf) array that is then incremented by one. We repeat for all NxN pixels. The resulting 2D pdf must be normalised to have overall unit probability (simply dividing by NxN achieves this). The 2D pdf is now ready for the next stage.
The computation of the 2D Shannon information entropy from the 2D gradient pdf is simple. Shannon’s classic logarithmic summation formula applies directly except for a crucial factor of one half which originates from special bandlimited sampling considerations for a gradient image (see arXiv paper for details). The half factor makes the computed 2D entropy even lower compared to other (more redundant) methods for estimating 2D entropy or lossless compression.
I’m sorry I haven’t written the necessary equations down here but everything is available in the preprint text. The computations are direct (non-iterative) and the computational complexity is of order (the number of pixels) NxN . The final computed Shannon information-entropy is rotation independent and corresponds precisely with the number of bits required to encode the image in a non-redundant gradient representation. 
By the way, the new 2D entropy measure predicts an (intuitively pleasing) entropy of 8 bits per pixel for the random image and 0.000 bits per pixel for the smooth gradient image in the original question.
A: Essentially the idea of entropy is something like "number of micro-states consistent with the macrostate".
I think the comment by sean507 and the answer by bottiger both point to a common framework. If you represent the image space by a generative model, $p[\,I,h\,]$, then for a given image $I$ you can (in principle) compute a posterior over the hidden states $p[\,h\mid I\,]$ (see also here). Then you can (in principle) compute the entropy of the posterior.
So I would agree that any "entropy", even in the "most theoretically correct sense", would seem to depend on both the representation used, and the generative model linking "microstates" ($h$) to "macrostates" ($I$). 
A: $$
H = - \sum_k p_k log_2(p_k)
$$
does NOT work in practice, for the simple reason that it's almost impossible to determine $P_k$.  You think that you can do it, as you've done by considering the number of grey levels.  Pk is not that.  Pk is all possible combinations of grey levels.  So you have to create a multi dimensional probability tree considering 1, 2, 3... combinations of pixels.  If you read Shannon's work you see him do this calculation for plain English considering a tree depth of 3 letters.  It then gets unwieldy without a computer.
You proved this yourself with statement 2.  That's why your entropy calculation returns the same level of entropy for the two images, even though one is clearly less ordered than the other.
There is also no such concept of spatial distribution within entropy calculation.  If there was, you'd also have to calculate entropy differently for temporally distributed samples.  And what would you do for an 11 dimensional data array?  For informational entropy; it is measured in bytes.
Just simply compress the images using a compression algorithm.  It will output an estimate of the entropy in bytes.  It will do this for any image or literally anything else that can be digitised, such as music or Shakespearean plays.
So.  Your random image contains approximately 114 KBytes, and your ordered image contains approximately 2.2 KBytes.  This is what you'd expect, but you already knew this kinda because you saw the image file sizes were of this size.  I have reduced the compressed size by 33% to allow for future improvements in compression algorithms.  I can't see them improving beyond this as the improvement curve is becoming asymptotic to a true underlying value.
P.S.  For interest, Shakespeare only produced 1 MByte of entropy in his entire life's work, calculated by this technique.  Most of it's quite good though.
A: The entropy of an image can be computed generalizing classical 1D
methods to the 2D case.
Obviously it must take into account “spatial information”, not just gray levels.
Classical approaches used in 1D systems take into account, at least potentially, configurations at all the scales (all the n-grams, also for very large n).
The fact that the correct Shannon entropy (of a sequence, or of an image) can be estimated only considering all the possible "n-grams" is a key point that did not sufficiently permeate the literature of image processing. This fundamental idea can be found in the seminal paper "Prediction and entropy of printed English" of Shannon (1952), and it is at the base of the classical methods used for 1D systems, that can be found in:
T. Schürmann and P. Grassberger, "Entropy estimation of symbol sequences", Chaos 6, 414 (1996).
A recent analysis of 2D generalisations of these methods can be found in: F.N.M. de Sousa Filho, V. G. Pereira de Sá, and E. Brigatti, "Entropy estimation in bidimensional sequences", Phys. Rev E, 105, 054116 (2022):
https://www.researchgate.net/publication/361742112_Entropy_estimation_in_bidimensional_sequences
where the block-entropies method and lossless compression methods are taken into account.
