As @NinadMunshi commented, the trick here is through a rotation of coordinates, so that the original coordinate vector $r$ and the new coordinate vector $r'$ are related by
$ r = R r' $
were $R = [u_1, u_2, u_3] $. Taking the third column vector of $R$ to be
$u_3 = [1, 1, 1]/\sqrt{3} $
Then $ x + y + z = [1, 1, 1] r = [1, 1, 1] R r' = \sqrt{3} z' $
And using spherical coordiantes for $r'$ we have
$ x' = \rho \sin(t) \cos(s) , y' = \rho \sin(t) \sin(s) , z' = \rho \cos(t) $
Then the integral becomes
$I = \displaystyle \int_{\rho = 0}^{1} \int_{s=0}^{2\pi} \int_{t=0}^{\pi} (\sqrt{3} \rho \cos(t) )^{2n} \rho^2 \sin(t) dt ds d\rho $
Integrating with respect to $\rho$ and $\phi$, reduce the integral to
$I =(2 \pi) \dfrac{(3)^n}{2n+3} \displaystyle \int_0^\pi \cos^{2n}(t) \sin(t) dt $
And this evaluates to
$I = \dfrac{(4 \pi) (3)^n}{(2n+3)(2n + 1)} $