I'm following the third recipe of [this answer](https://stats.stackexchange.com/questions/313052/calculate-entropy-of-sample) to estimate the Shannon entropy of my samples using histograms.  My expectation was, increasing the sample size should lead to a better estimation of the *true* entropy. To test that, I sampled from a Gaussian distribution $N(0, \sigma^2)$ with [known](https://en.wikipedia.org/wiki/Mutual_information#Linear_correlation) entropy $H(X) = 0.5[1 + \log(2\pi \sigma^2)]$:

```python
import numpy as np
from scipy.stats import entropy # by default in natural log
import matplotlib.pyplot as plt

sigma = 2
H0 = 0.5*(1 + np.log(2*np.pi*sigma**2)) # theoretical value in nats

Hs = []
ns = np.logspace(2,6,5, dtype=int)
for n in ns:
    X = np.random.normal(0, sigma, size=n)
    nbins = int(n/20)  # <--- This is causes divergence, read UPDATE 2
    hist, bin_edges = np.histogram(X, nbins, density=True)
    Hs.append(entropy(hist)) # will be automatically normalized to sum to 1

plt.plot(ns, Hs, '-x', label='estimation')
plt.hlines(H0, ns[0], ns[-1], linestyle='--', color='k', label='exact H')
plt.legend()
plt.yticks(np.arange(11))
plt.xscale('log')
plt.xlabel('num samples')
plt.ylabel('sample entropy')
```

To my surprise, they diverge:

[![Enter image description here][1]][1]

Am I missing something very simple?

I found a [relevant question](https://math.stackexchange.com/questions/395121/how-entropy-scales-with-sample-size) that is concerned about the case where the probability distribution is either "peaked" or "flat". The upshot of the most voted answer is that one has to normalize the entropy by $\log n$ where $n$ is the number of samples. Although I find this normalization sensible, I still observe an offset between the theoretical value and the estimation:

```python
plt.plot(ns, np.array(Hs)/np.log(ns), '-x', label='normalized estimation') # normalized
```
[![enter image description here][2]][2]

What causes this discrepancy?

As @sextus-empiricus noted below, scaling the number bins with the sample size leads to a coarse density estimation. In fact, **fixing the number of bins, resolves the divergence**. However, no matter what number of bins one chooses, there's still some offset between the estimated (normalized or non-normalized) entropy and the theoretical value, as you can see below:

[![Enter image description here][3]][3]

What is the reason behind this gap? The curve displacement by `nbins` clearly suggest that there's a relationship between the two, but I found the reason elusive.

  [1]: https://i.sstatic.net/UZWDQ.png
  [2]: https://i.sstatic.net/61vMp.png
  [3]: https://i.sstatic.net/uKsVI.png