The infinitely negative energy density that we find if we use a non-normal ordered Hamiltonian can be understood as a result of the Dirac sea. Moreover, we only find an infinitely negative energy density for fermions as a result of the anticommutation relations:
$$H = \int \frac{d^3p}{(2\pi)^3} E_p \left(a_p^\dagger a_p - \frac12\{a_p,a_p^\dagger\}\right),$$ where the second term is an infinitely large constant.
In contrast, for bosons we find if we use a non-normal ordered Hamiltonian an infinitely positive energy density:
$$H = \int \frac{d^3p}{(2\pi)^3} E_p \left(a_p^\dagger a_p + \frac12[a_p,a_p^\dagger]\right),$$ where again the second term is an infinitely large constant.
Is it possible to explain this through a similar idea as the Dirac sea does it for fermions?
PS: One obvious problem is that there is no exclusion principle for bosons. However, even for fermions the exclusion principle alone is not sufficient, to quote from Wikipedia
[the] Pauli exclusion does not definitively mean that a filled Dirac sea cannot accept more electrons, since, as Hilbert elucidated, a sea of infinite extent can accept new particles even if it is filled.
