\begin{align} I &=\int_0^1dr\int_0^{2\pi}\left(cos\theta\cdot\frac{\partial f}{\partial x}+sin\theta\cdot\frac{\partial f}{\partial y}\right)rd\theta=\int_0^1dr\int_{\color{red}{x^2+y^2=r^2}}\left({\frac{\partial f}{\partial x}}dy-{\frac{\partial f}{\partial y}}dx\right) \\ & \quad =\color{blue}{\int_0^1dr\iint_{x^2+y^2\le r^2}\left(\frac{{\partial}^2 f}{\partial x^2}+\frac{{\partial}^2 f}{\partial y^2}\right)dxdy}=\color{green}{\int_0^1dr\iint_{x^2+y^2\le r^2}\left(x^2y^2\right)dxdy} \\ & \quad \quad =\int_0^1dr\int_0^r{\rho}^5d\rho\int_0^{2\pi}cos^2\theta sin^2\theta d\theta \end{align}
Green part is given by another equation, so just forget it. I don't understand why the range(red part) become a equation like this; and for the blue part, I have totally no idea.
The original integral is a double integral on the unit disk. The second equality (the one with the red part) comes from changing the polar coordinate double integral to line integral on the circle $x^2+y^2=r^2$ where $r$ varies from $0$ to $1$. This covers the whole double integral on the disk. The next equality comes from Green's Theorem, which says that $$\int_C Pdx+Qdy=\iint_D \frac{\partial Q}{\partial x}-\frac{\partial P}{\partial y} dxdy$$
if $C$ is a positively oriented, piecewise smooth, simple, closed curve and $D$ is the region enclosed by $C$.