How to verify that an operator is a tensor operator? Suppose $\boldsymbol R,\boldsymbol P$ are the common position and momentum operators, and $\boldsymbol L=\boldsymbol R\times \boldsymbol P$ is the orbital angular momentum .
$\boldsymbol K$ is a linear combination of $\boldsymbol R,\boldsymbol P$, and Jordan product of any of their cross product, e.g. $\frac12(\boldsymbol R\times (\boldsymbol R\times \boldsymbol P)-(\boldsymbol R\times \boldsymbol P)\times \boldsymbol R)$, whose coefficients are scalar functions of $ R^2,P^2, \boldsymbol R\cdot \boldsymbol P$, independent from coordinates, i.e. terms like $\boldsymbol R\cdot \boldsymbol{\hat x}$ are not involved.
Let $\boldsymbol J$  be the collection of three angular momentum operators, then a vector operator $\boldsymbol V$ satisfies
$$
[J_i,V_j]=\text i\hbar\epsilon_{ijk}V_k
$$
Question: Can we generally proof that $\boldsymbol K$ satisfies the defnition of vector operators, or is there a counter example? In other words, is
$$
[L_i,K_j]=\text i\hbar\epsilon_{ijk}K_k
$$
true for all $\boldsymbol K$s? One example is the Runge-Lenz operator of a Hydrogen atom
$$
\boldsymbol B=\frac{\boldsymbol P\times\boldsymbol L-\boldsymbol L\times \boldsymbol P}{2me^2}-\frac{\boldsymbol R}R
$$
which is confirmed to be a vector operator.
Note: The key is to justify the fact that an operator directly constructed by a classic vector is indeed a vector operator, so simply referring to the definition can't be the proof.
 A: What defines a tensor operator is the commutation relation of its components with the angular momenta.  Thus, if $[L_a,K_b]=i\hbar\epsilon_{abc}K_c$ you have a vector operator by definition.
More generally you want something like
\begin{align}
[L_\pm,T^{\lambda}_\mu]&=\sqrt{(\lambda\mp\mu)(\lambda\pm \mu+1)}T^{\lambda}_{\mu\pm 1}\\
[L_0,T^{\lambda}_\mu]&=\mu T^{\lambda}_\mu\, ,
\end{align}
which tends to be easier to work with and goes beyond vectors ($\lambda$ is not restricted to $1$).  It is also useful to use
spherical components where, for instance,
$$
T^1_1=-\frac{V_x+iV_y}{\sqrt{2}}
$$
It is a good bet that, if you start with a classical vector with components $x^a p_y^b$ with $a$ or $b=0$ or $1$ (or such types of products), the corresponding observable will also be a vector operator.  If you have higher powers in $x$ or $p$ (i.e. both $a$ and $b\ge 2$), there may be ordering issues and you have to be quite careful.  (You want to be aware that a product of hermitian operators remains hermitian only if the operators commute.)
Polynomials in $x,y,z$ alone, or in $p_x,p_y,p_z$ alone, can be decomposed into tensors by expressing them in terms of spherical harmonics.  Thus,
\begin{align}
x+iy&= r \sin\theta e^{i\phi} = c r Y^1_{-1}(\theta,\phi)\, ,\\
(x+iy)(x-iy)&= a r^2 + b r^2 Y^2_{0}(\theta,\phi)\, ,\\
\end{align}
for some constants $a,b,c$
and are thus linear combinations of tensor operators, with the latter containing a scalar and the component of a  tensor with $\lambda=2$.
In the case of more sophisticated products like the ones you have, you need to use the definition of composite tensors to properly combine the individual vectors using Clebsch-Gordan technology.  Thus, again using spherical coordinates:
$$
Z^{j}_m =\sum_{qq'} T^{\lambda}_q W^{\lambda'}_{q'} C_{\lambda,q;\lambda'q'}^{jm}
$$
so it then becomes a matter of checking if the complicated expressions are properly combined using CGs (again note the use of spherical components).
A good source for this is the text by

Baym, Gordon. Lectures on quantum mechanics. CRC Press, 2018.

