For any square matrix $Y$ let $\chi_x(Y) = det(xI -Y)$ denote its characteristic polynomial.
Say $A$ and $B$ are two $n-$dimensional symmetric matrices with constant row sums $a$ and $b$. Lets define the polynomials $p_1(x)$ and $p_2(x)$ as follows : $\chi_x(A)=(x-a)p_1(x)$ and $\chi_x(B)=(x-b)p_2(x)$. Then over uniform sampling from the permutation group $S_n$ one can show (quite a non-trivial proof) that ``finite free convolution" (denoted as $\boxplus$) satisfies the following identity,
$$\mathbb{E}_{P \sim S_n} [\chi_x(A + PBP^T)] = (x-(a+b))[p_1(x) \boxplus p_2(x)]$$
Given a $n-$dimensional symmetric matrix $M$ such that, $\chi_x(M) = (x-1)^{\frac {n}{2}}(x+1)^{\frac {n}{2}}$ define a polynomial $p$ such that $\chi_x(M)=(x-1)p(x)$. Now apparently the following identity holds for any positive integer $d$,
$$\underset{P_1,P_2,..,P_d \sim S_n}{\mathbb{E}}[\chi_x(P_1MP_1^T+ P_2MP_2^T+..+P_dMP_d^T)]\\ = (x-d)[p(x)\boxplus p(x) ..(d \text{ times})..\boxplus p(x)]$$
Can someone kindly help derive the second equality from the first?
I believe this is some kind of an induction but I am unable to get it work. As in even at $d=3$ I cant get this explicitly.