How to derive the gradient of the reparameterized score function estimator? In the paper Evolution Strategies as a Scalable Alternative to Reinforcement Learning, the authors derive the following gradient of the score function estimator
$$
\begin{align}
\nabla_\psi\mathbb E_{\theta\sim p_\psi}[F(\theta)]&=\mathbb E_{\theta\sim p_\psi}[F(\theta)\nabla_\psi\log p_\psi(\theta)]\\
&=\nabla_\theta\mathbb E_{\epsilon\sim\mathcal N(0,I)}F(\theta+\sigma\epsilon)\\
&={1\over\sigma}\mathbb E_{\epsilon\sim\mathcal N(0,I)}\{F(\theta+\sigma\epsilon)\epsilon\}
\end{align}
$$
where $p_\psi$ is a multivariate Gaussian, $F$ is the objective function(e.g., return on reinforcement learning problems) and $\theta+\sigma\epsilon$ is the result of the reparameterization trick. The meaning of $\theta$ changes from samples from $p_\psi$ to the mean of the $p_\psi$ in the second step, following the same convention used by the paper. I'm wondering how the last step is derived? 
 A: In the list of equations you have written above, you are mixing up two different approaches. In particular, line 2 (their approach) does not follow from line 1 (the REINFORCE objective). Even if $p_\psi$ is a gaussian, this would only follow if the variance is fixed ($\psi = \theta$). 
As for the step line 2 to line 3, I think this is an approximation result based on Taylor series of $F(\cdot)$ about $a = \theta$
$$
\begin{align}
\mathbb E_{\epsilon\sim\mathcal N(0,I)} \left[\epsilon F(\theta+\sigma\epsilon)\right] &= \mathbb E_{\epsilon\sim\mathcal N(0,I)}\left[ \epsilon \left( F(\theta) + \sigma\epsilon \nabla F(\theta) + \dotso \right) \right]\\
&\approx \mathbb E_{\epsilon\sim\mathcal N(0,I)}\left[ \epsilon F(\theta) + \sigma\epsilon^2 \nabla F(\theta) \right]
\end{align}
$$
Note that since $\epsilon \sim \mathcal{N}(0,I)$, we have that $E[\epsilon] = 0$, $E[\epsilon^2] = \mathbf{I}$. So the above expression simplifies to 
$$
\begin{align}
\mathbb E_{\epsilon\sim\mathcal N(0,I)} \left[\epsilon F(\theta+\sigma\epsilon)\right] \approx \sigma \nabla F(\theta)\\
\end{align}
$$
and rearranging we get,$$\nabla F(\theta)   \approx \frac{1}{\sigma}\mathbb E_{\epsilon\sim\mathcal N(0,I)} \left[\epsilon F(\theta+\sigma\epsilon)\right]$$
I am not sure why they have the relationship as an equalilty. Perhaps the approximation is tight since for $k > 0$, $\mathbb E[\epsilon^{2k+1}] = 0$.
