Let $J$ be a random matrix, i.e. all elements are drawn randomly, with zero mean and $\mathrm{E}[J_{ij}^2] = \frac{1}{N}$ and $\mathrm{E}[J_{ij}J_{ji}] = \frac{\tau}{N}$, then it is given that its Gaussian measure obeys the formula
$P(J) \propto \exp\big[ -\frac{N}{2(1 - \tau^2)}\mathrm{Tr}(JJ^T - \tau JJ) \big]$, where $J_{ij}^T = J_{ji}$.
How could this formula be derived ? From standard formula of Gaussian distribution ?
Consider first (as a warm-up exercise) a $2\times 2$ case. The matrix $\mathbf{J}$ can be parametrized as $$ \mathbf{J}= \begin{pmatrix} J_{11} & J_{12} \\ J_{21} & J_{22} \\ \end{pmatrix}\ . $$ If you evaluate the term $\frac{N}{2(1-\tau^2)}\mathrm{Tr}(\mathbf{J}\mathbf{J}^T-\tau \mathbf{J}\mathbf{J})$ you obtain \begin{equation} \frac{N \left(-\tau \left(J_{11}^2+J_{12} J_{21}\right)+J_{11}^2+J_{12}^2-\tau \left(J_{12}J_{21}+J_{22}^2\right)+J_{21}^2+J_{22}^2\right)}{2 \left(1-\tau ^2\right)}\ . \end{equation} Therefore, the term $\exp\left[-\frac{N}{2(1-\tau^2)}\mathrm{Tr}(\mathbf{J}\mathbf{J}^T-\tau \mathbf{J}\mathbf{J})\right]$ can be cast in the form of a multivariate Gaussian joint probability density (jpd) $$ f_G(J_{11},J_{12},J_{21},J_{22})\propto \exp\left[-\frac{1}{2}\vec{J}^T {\bf\Sigma}^{-1}\ \vec{J}\right]\ , $$ where $\vec{J}$ is a column vector $\vec{J}^T = (J_{11},J_{12},J_{21},J_{22})$, and ${\bf \Sigma}$ a $4\times 4$ covariance matrix. I let you work out the detailed form of ${\bf\Sigma}$. From the diagonal elements of ${\bf\Sigma}$, it should be easy to see that (for instance) $\mathbb{E}[J_{11}^2]=(1+\tau)/N$. Of course, this 'proof' can be run backwards and for general $N$, so indeed you start from a multivariate Gaussian jpd for your $N^2$ matrix elements (arranged in a column vector), with a covariance matrix suitably constructed in order to accommodate your constraints on variances and covariances.