I have the following derivation here
Where M is a square matrix, 
is also a square matrix (diagonal). But I do not know why the derivative of 
is 
.
I only know before that the derivative of tr$(X^2)$ with respect to $X$ is 2$X$ from the matrix cookbook, but I find it no way to deal with 
.
It can be solved by taking matrix derivative, Let $Y = M \Psi M^\mathrm{T}$ and $J = \text{tr}[( M \Psi M^\mathrm{T})^2]=\text{tr}(Y^2)$, we take derivative w.r.t $Y$ and $M$, respectively,
$$ \begin{aligned} \mathrm{d} J &= \mathrm{tr}[\mathrm{d}Y^2 ] \\ &=\mathrm{tr}[(\mathrm{d}Y) Y + Y (\mathrm{d}Y)] \\ &=\mathrm{tr}[(2Y (\mathrm{d}Y)] \\ \end{aligned} $$ so,
$$ \frac{\partial J}{\partial Y} = 2 Y^\mathrm{\mathrm{T}} $$ and from matrix derivaty, we get, $$ \begin{aligned} \mathrm{d} J &= \mathrm{tr}[\frac{\partial^\mathrm{\mathrm{T}} J}{\partial Y} \mathrm{d}Y ] \\ &=\mathrm{tr}[\frac{\partial^\mathrm{\mathrm{T}} J}{\partial Y} \mathrm{d}(M \Psi M^\mathrm{\mathrm{T}}) ] \\ &=\mathrm{tr}[\frac{\partial^\mathrm{\mathrm{T}} J}{\partial Y} (\mathrm{d}M) \Psi M^\mathrm{\mathrm{T}} + \frac{\partial^\mathrm{T} J}{\partial Y} M \Psi (\mathrm{d}M^\mathrm{T})] \\ &=\mathrm{tr}[\Psi M^\mathrm{T}\frac{\partial^\mathrm{T} J}{\partial Y} (\mathrm{d}M) + \frac{\partial^\mathrm{T} J}{\partial Y} M \Psi (\mathrm{d}M^\mathrm{T})] \\ &=\mathrm{tr}[\Psi M^\mathrm{T}\frac{\partial^\mathrm{T} J}{\partial Y} (\mathrm{d}M)] + \mathrm{tr}[\frac{\partial^\mathrm{T} J}{\partial Y} M \Psi (\mathrm{d}M^\mathrm{T})] \\ &=\mathrm{tr}[\Psi M^\mathrm{T}\frac{\partial^\mathrm{T} J}{\partial Y} (\mathrm{d}M)] + \mathrm{tr}[(\mathrm{d}M) \Psi M^\mathrm{T} \frac{\partial J}{\partial Y}] \\ &=\mathrm{tr}[\Psi M^\mathrm{T}\frac{\partial^\mathrm{T} J}{\partial Y} (\mathrm{d}M)] + \mathrm{tr}[ \Psi M^\mathrm{T} \frac{\partial J}{\partial Y} (\mathrm{d}M)] \\ &=\mathrm{tr}[ (\Psi M^\mathrm{T}\frac{\partial^\mathrm{T} J}{\partial Y} + \Psi M^\mathrm{T} \frac{\partial J}{\partial Y} )\mathrm{d}M] \\ \end{aligned} $$ so, $$ \begin{aligned} \frac{\partial J}{\partial M} &= \frac{\partial J}{\partial Y} M \Psi + \frac{\partial^\mathrm{T} J}{\partial Y} M \Psi \\ &= 4(M \Psi M^\mathrm{T}) M \Psi \end{aligned} $$