 Solution

For one Fourier mode, put $\mu =\hat{\mathbf{k}}\cdot \mathbf{e}$. The collisionless equation becomesUse the unweighted Legendre polynomial expansion specified in the paper: $\Theta ={\sum }_{j\ge 0}(-i{)}^j{\Theta }_j{P}_j(\mu )$. The Legendre polynomial recurrence relation saysThe coefficient of $P_{\ell }$ in $ik\mu \Theta$ receives contributions from $j=\ell -1$ and $j=\ell +1$. Dividing by $(-i{)}^{\ell }$ gives respectively $-k\ell {\Theta }_{\ell -1}/(2\ell -1)$ and $k(\ell +1){\Theta }_{\ell +1}/(2\ell +3)$. The metric sources are $\dot{\varphi }{P}_0$ and $-ik\psi P_1$, the latter becoming $k\psi$ after division by $(-i)$. Hence the neutrino Boltzmann hierarchy isThe lower-neighbour term is absent for $\ell =0$. The low moments areThese coefficients depend on the expansion convention. In the convention with a $(2\ell +1)$ factor, the temperature multipoles are ${\Theta }_{\ell }/(2\ell +1)$; using those multipoles without converting them would give incorrect factors in the scalar neutrino anisotropic stress.