I am trying to evaluate $E[f(x)f(x)^T]$ where $f(x) = Nx$ and $N \in \mathbb{R}^{N \times N}$, and also $p(x) = \mathbb{N} (x|μ, \sigma)$
So, as we know,
\begin{align*} E[f(x)f(x)^T] &= \int_{-\infty}^{\infty} f(x)f(x)^{T} p(x)dx \\ &= \int_{-\infty}^{\infty} Nx(Nx)^{T} p(x) dx \\ &= \int_{-\infty}^{\infty} Nxx^{T}N^{T} p(x) dx \end{align*}
Now, I am confused about the next step as I haven't evaluated $\mathbb{E}$ with a matrix before. Can I drag the $NN^T$ part outside as it's constant with respect to $x$? Can I break the order? What to do next?
Assuming the matrix doesn't depend on the random vector $x$, you can indeed take the matrix outside the expectation, since matrix multiplication is a linear operation. You just have to be careful to remember that matrix multiplication is "order dependent" i.e. not commutative:
$$ \mathbb{E}[(Nx)(Nx)^T] = \mathbb{E}[Nxx^T N^T] = N\mathbb{E}[xx^T]N^T = NC N^T $$ where $C = \mathbb{E}[xx^T]$.
Hope that helps!