A well-known way of encoding $N$ levels of a harmonic (bosonic) oscillator is as follows: \begin{equation} |n\rangle = |1\rangle^{\otimes n} \otimes |0\rangle^{\otimes N-n+1} \quad,\qquad 0\leq n \leq N \quad, \end{equation} with the ladder operators being defined as: \begin{equation} \label{bosenc} \begin{alignedat}{99} a &= |000\rangle\langle100|_{1,\,2,\,3} + \sum\limits_{j=2}^{N} \sqrt{j}|100\rangle\langle110|_{j-1,\,j,\,j+1} + \sqrt{N+1}|110\rangle\langle111|_{N-1,\,N,\,N+1} \quad&&,\\ a^\dagger &= |100\rangle\langle000|_{1,\,2,\,3} + \sum\limits_{j=2}^{N} \sqrt{j}|110\rangle\langle100|_{j-1,\,j,\,j+1} + \sqrt{N+1}|111\rangle\langle110|_{N-1,\,N,\,N+1} \quad&&. \end{alignedat} \end{equation} The subscripts denote on which qubits the term is acting (for the rest of qubits it's just the unity operator). Indeed, one can check that the commutator of the so-defined operators is \begin{equation} [a,\,a^\dagger] = \sum \limits_{n=0}^{N} |n\rangle\langle n| - (N+1)|N+1\rangle\langle N+1| \quad. \end{equation} (the last term being the unavoidable 'correction' due to the absence of finite-dimensional representations of the Heisenberg algebra)
An obvious downside of such an approach is that we are only using $N$ out of total of $2^{N+1}$ states of the $N+1$ qubit system. The main advantage, however, stems from the fact that the raising and lowering operators are $3$-local, and contain only $N+1$ terms each.
Is anybody aware of more efficient ways of storing bosonic degrees of freedom? I'm personally interested in encoding multiple bosonic degrees of freedom having various numbers of levels; so, not only I'm interested in "what's the minimum number of qubits required to efficiently work with $N$ bosonic levels?", but also "given a fixed number of qubits, how many bosonic oscillators (with, possibly, different number of levels) can one efficiently encode?" Obviously, 'efficiently' here implies that the raising and lowering operators have some reasonable scaling properties.
UPDATE
The approach opposite to the one presented above was suggested by Craig Gidney in comments: one can encode the $n$-th level of the oscillator with its binary representation $|(n)_2\rangle$, where $(\ldots)_2$ stands for the binary form of the number. The number of terms in the expression for ladder operators is the same as in the first approach (it is exponential in the number of qubits - but the number of qubits itself is logarithmic in the number of levels). Note, however, that the main downside of this approach is that each term in the sum is not $3$- or any fixed number-local. Instead, they all are $\log N$-local: $$ a = \underbrace{|\ldots0000\rangle}_{\text{$\log N$ qubits}}\langle\ldots0001| + \ldots $$
It would be interesting though to find a compromise between these two approaches.