For the *inward* pointing unit normal, it seems we arrive at the ($L^2(\Omega)$) energy identity, $$\frac{d}{dt}\|u\|^2 + \|\nabla u\|^2 + \int_\Gamma u^2 u dx =0.$$
Assuming the "solution" $u$ is nonnegative and continuous in time, we rewrite the boundary integral using the MVT for Integrals (a suitable version of this) as follows: for each $t>0$, there is $\xi(t)\in\Gamma$ in which $$\int_\Gamma u^2(t,x) u(t,x) d\sigma=u(t,\xi(t))\int_\Gamma u^2(t,x)d\sigma.$$
Hence, $$u(t,\xi(t))=\frac{\langle u^3(t),1 \rangle_\Gamma}{\|u(t)\|^2_\Gamma}.$$
The identity now reads $$\frac{d}{dt}\|u(t)\|^2 + \|\nabla u(t)\|^2 + u(t,\xi(t))\|u(t)\|^2_\Gamma =0.$$
Let $\lambda>0$ be the best constant so that the Poincare/Sobolev-inequality holds,
$$\lambda\int_\Omega u^2 dx \leq \int_\Omega |\nabla u|^2dx + \int_\Gamma u^2 d\sigma.$$
Define $m(t):=\min\{1,u(t,\xi(t))\}\geq0$.
Considering the above, we arrive at the estimate $$\frac{d}{dt}\|u(t)\|^2+\lambda m(t)\|u(t)\|^2\leq 0.$$