Characterizing the harmonic oscillator creation and annihilation operators in a rotationally invariant way I am interested in a characterization of the creation and annihilation operators that is in some sense invariant under $O(n)$ rotations of $\mathbb{R}^n$:
Background
The Harmonic Oscillator on $\mathbb{R}^n$ is the differential operator 
$$ H := \sum_{k=1}^n \left[x_k^2-\frac{\partial^2}{\partial x_k^2}\right] = |x|^2 + \nabla.$$
It is not hard to see that the $L^2(\mathbb{R}^n)$ eigenvalues are exactly $\{n,n+2,n+4,\dots\}$. Furthermore, the annihilation operator is an operator on Schwartz functions on $\mathbb{R}^n$ $$C_k := \frac {1}{\sqrt {2}}\left( x_k + \frac{\partial}{\partial x_k}\right)$$ and the creation operator as its adjoint (with the $L^2$ inner product) $$ C_k^\dagger = \frac {1}{\sqrt {2}}\left( x_k -\frac{\partial}{\partial x_k}\right).$$ If we let $V_{n+2m}$ be the set of eigenfunctions with eigenvalue $n+2m$, we can show that $$C_k^\dagger V_{n+2m} \subset V_{n+2(m+1)}$$
$$C_k V_{n+2m} \subset V_{n+2(m-1)}$$ (where we define $V_r=0$ if $r$ is not an eigenvalue) and $V_n$ is spanned by $e^{-|x|^2/2}$. It turns out that $V_{n+2m}$ is isomorphic to the space of degree $m$ homogeneous polynomials in $n$ variables, which I'll denote $\mathcal{P}^n_m$, by the isomorphism $p \mapsto p(C^\dagger)$ i.e. $x_1x_2 \mapsto C_1^\dagger C_2^\dagger$, etc. All of this and more can be found here starting on page 86 (with some slightly different notation than I've used here). 
Motivation
One of the problems with this whole business is that even though $H$ and thus $V_{n+2m}$ are invariant under rotations of $\mathbb{R}^n$, the $C_k$ and $C_k^\dagger$ are not. We made an arbitrary choice of coordinates when we defined them. This leads to a non-cannonical choice of basis for the eigenspaces, and has been giving me problems in my research. My question is thus:
Question
Even though there is no cannonical choice of basis for $V_{n+2m}$, is there some characterization of the creation and anihilation operators that is invariant under $O(n)$ rotations of $\mathbb{R}^n$. That is, what is an $O(n)$-invariant characterization of the space of operators $$ \text{span}_\mathbb{R}\{C_1^\dagger,C_2^\dagger,\dots, C^\dagger_n\}$$
If this is not possible, then as an alternative answer, I am interested in insight into how rotations and the operators interact. 
 A: Let me try to answer this question.  I apologise if my notation is slightly different, since I will work in some more generality, since the equivariance properties of the creation and annihilation operators are actually more transparent, I believe, relative to the the general linear group instead of the orthogonal group.  Also the fact that this is the harmonic oscillator is a red herring.  In the case of the harmonic oscillator, we have to introduce more structure, reducing the group of symmetries.
It is also in this case, where the grading to be discussed below coincides (up to a choice of scale) with the energy of the system.

Let $E$ be an $n$-dimensional real vector space and let $E^*$ denote its dual.  Then on $H = E \oplus E^* \oplus \mathbb{R}K$ one defines a Lie algebra by the following relations:
$$ [x,y]= 0 = [\alpha,\beta] \qquad [x,\alpha] = \alpha(x) K = - [\alpha,x] \qquad [K,*]=0$$
for all $x,y \in E$ and $\alpha,\beta \in E^*$.  This is called the Heisenberg Lie algebra of $E$, denoted $\mathfrak{h}$.
The automorphism group of $\mathfrak{h}$ is the group $\operatorname{Sp}(E\oplus E^*)$ of linear transformations of $E\oplus E^*$ which preserve the symplectic inner product defined by the dual pairing:
$$\omega\left( (x,\alpha), (y,\beta) \right) = -\alpha(y) + \beta(x).$$
Let $\mathfrak{a} < \mathfrak{h}$ denote the abelian subalgebra with underlying vector space $E \oplus \mathbb{R}K$.  One can induce a $\mathfrak{h}$-module from an irreducible (one-dimensional) $\mathfrak{a}$-module as follows.  Let $W_k$ denote the one-dimensional vector space on which $E$ acts trivially and $K$ acts by multiplication with a constant $k$.  Then letting $U$ be the universal enveloping algebra functor, we have that
$$ V_k = U\mathfrak{h} \otimes_{U\mathfrak{a}} W_k$$
is an $\mathfrak{h}$-module.  The Poincaré-Birkhoff-Witt theorem implies that $V_k$ is isomorphic as a vector space to the symmetric algebra of $E^*$, which we may (as we are over $\mathbb{R}$) identify with polynomial functions on $E$.
The subgroup of $\operatorname{Sp}(E\oplus E^*)$ which acts on $V_k$ is the general linear group $\operatorname{GL}(E)$ and hence $V_k$ becomes a $\operatorname{GL}(E)$-module.  In fact, $V_k$ is graded (by the grading in the symmetric algebra of $E^*$ or equivalently the degree of the polynomial):
$$V_k = \bigoplus_{p\geq 0} V_k^{(p)}$$
and each $V_k^{(p)}$ is a finite-dimensional $\operatorname{GL}(E)$-module isomorphic to $\operatorname{Sym}^p E^*$.
Every vector $x \in E$ defines an annhilation operator: $A(x): V_k^{(p)} \to V_k^{(p-1)}$ via the contraction map
$$E \otimes \operatorname{Sym}^p E^* \to \operatorname{Sym}^{p-1} E^*$$
whereas every $\alpha \in E^*$ defines a creation operator: $C(\alpha): V_k^{(p)} \to V_k^{(p+1)}$ by the natural symmetrization map
$$E^* \otimes \operatorname{Sym}^p E^* \to \operatorname{Sym}^{p+1} E^*.$$
Both of these maps are $\operatorname{GL}(E)$-equivariant, and this is perhaps the most invariant statement I can think of concerning the creation and annihilation operators.
A: The SO(n) invariant subspaces of polynomials of the creation operators are completely characterized by the theory of spherical harmonics, see for example, the following  exposition. 
Summary of the construction, Let:
$ d(z) = z_1^2+ . . . + z_n^2$
$ \Delta = (\frac{\partial}{\partial z_1})^2 + . . . + \frac{\partial}{\partial z_n})^2$
Let $\mathcal{P}_m$ be the space of homogeneous polynomials of z of degree m, and 
$\mathcal{H}_m = \mathcal{P}_m  \cap ker(\Delta)$
Then the SO(n) invariant subspaces are given by:
$\mathcal{P}^l_m = d^l \mathcal{H}_m $
Reply to the comment:
$SO(n)$ acts on the linear space spanned by the creation operators according to the fundamental vector representation. The creation operators commute among themselves,
thus, one can identify their span (from the point of view of invariant theory) with $C^n$.
The components of the vector $z$ in the answer can be identified with the creation perators. The only invariants of the $SO(n)$ action are functions of $d(z)$ (and its complex conjugate). In addition, one can construct invariant subspaces of polynomials of the creation operators, according to the method given above. The action of these polynomial subspaces on the cyclic vector generate the irreducible $SO(n)$ subspaces of eigenfunctions of the  isotropic harmonic oscillator Hamiltonian.
Remark: 
The homogeneous subspaces $V_{n+2m}$ are not $SO(n)$ irreducible and they decompose according to corollary 1 in the reference.
A: As a physicist let me try to give my opinion on this question. First, the polynomial representation you mentioned is called the Fock-Bargmann representation (in physics literiture) which can be found in Bargmann's paper in 60's (Bargmann was Einstein's assistant). This representation is useful for its convinience in dealing with Gaussian states. Second, while the statement that $H$ is invariant under ratation is probabaly ture,  $V_{n+2m}$ does change under such rotation. What happens here is that $H$ has some degeneracy (the energy quantum for each mode $k$ is the same), rotations in the degenerate subspace don't change $H$. Third, may be the most convinient way to study the changing of states under rotation is to use the so called coherent state (eigenstates of the annihilation operators) instead of number states (eigenstates of $H$). Such coherent states which form an over complete non-orthogonal basis can be parametrized by a vector $\vec{\alpha}$, and such vectors rotate just like $C$. 
