Let $X$ be an $m$ by $n$ matrix. I would like to find the $n^2$ by $(mn)^2$ matrix $B$ such that
$$
\operatorname{vec}\left(X^{\top}X\right) = B \operatorname{vec}\left(\operatorname{vec}\left(X\right)\left(\operatorname{vec}\left(X\right)\right)^{\top}\right).$$
Here $\operatorname{vec}(\cdot)$ is the vector function that unrolls a matrix into a vector columnwise.
I would guess $B$ involves the Kronecker products, and some $I_m$ and $I_n$ matrices and some 1 vectors, but I am having a hard time with all the indices.
Matrix expression for $\operatorname{vec}(X^{\top}X)$
kronecker productmatrices
Best Answer
I believe $B$ takes the form $$ B = \left(I_n \otimes \left(\operatorname{vec}\left(I_m\right)\right)^{\top} \otimes I_n\right)\left( I_{mn} \otimes K \right), $$ where $K$ is the 'commutation matrix.' That is, $K$ is the matrix such that $K \operatorname{vec}(X) = \operatorname{vec}\left(X^{\top}\right).$
Letting $$ V = \left(I_{mn} \otimes K \right) \operatorname{vec}\left(\operatorname{vec}\left(X\right)\left(\operatorname{vec}\left(X\right)\right)^{\top}\right) = \operatorname{vec}\left(\operatorname{vec}\left(X^{\top}\right)\left(\operatorname{vec}\left(X\right)\right)^{\top}\right), $$ then the $i,j$ element of $\operatorname{vec}\left(X^{\top}X\right)$ should be: $$ \operatorname{vec}\left(X^{\top}X\right)_{i,j} = \sum_k X_{k,i} X_{k,j} = V_{i+mk+mn(k+mj)}, $$ where we index from zero, as suggested by @Chrystomath. The term $\left(I_n \otimes \left(\operatorname{vec}\left(I_m\right)\right)^{\top} \otimes I_n\right)$ then captures that indexing and the repeated $k$.
I have confirmed this relationship with some R code: