I am reading the paper Factor analysis and outliers: A Bayesian approach. The author starts with a factor analysis model given by $${\bf y}_i = {\bf \Lambda} {\bf z}_i + {\bf e}_i, \quad i = 1, \ldots, n,$$ where each ${\bf y}_i$ is a $p$-dimensional observation vector, each ${\bf z}_i$ is a $K$-dimensional latent factor vector, and ${\bf \Lambda}$ is a $p \times K$ full-rank matrix of factor loadings. The author assumes that the factors and the error term are Normal: $${\bf z}_i \sim \mathcal{N} ({\bf 0}, {\bf \Phi})$$ $${\bf e}_i \sim \mathcal{N} ({\bf 0}, {\bf \Psi})$$
The author assigns Wishart priors to ${\bf \Phi}^{-1}$ and ${\bf \Psi}^{-1}$: $${\bf \Phi}^{-1} \sim \mathcal{W}_K \left( {\bf \Phi}_{*}, \nu_{*} \right)$$ $${\bf \Psi}^{-1} \sim \mathcal{W}_p \left( {\bf \Psi}_{*}, n_{*} \right)$$
In the paper the author writes something I found to be quite interesting:
While classical factor analysis sets $\bf \Phi = I$ and uses a diagonal $\bf \Psi$ matrix, we impose these restrictions via the prior information matrices ${\bf \Psi}_{*}$ and ${\bf \Phi}_{*}$.
Question: What should the values of ${\bf \Psi}_{*}$ and ${\bf \Phi}_{*}$ be in order to do what the author is suggesting?
The author does not seem to state exactly how this can be done, but I may have missed it so I will continue reading it. My own research on this matter pointed me to these seemingly similar unanswered questions here and here.
UPDATE: I did some research on the Wishart distribution and if you specify that $\Psi_*$ and $\Phi_*$ are two diagonal matrices, then $\mathbb{E} [\Psi]$ and $\mathbb{E} [\Phi]$ will be two diagonal mean matrices. Perhaps, this is what the author is referring to. Still unsure, though.
UPDATE 2: I set $\Psi_*$ and $\Phi_*$ to diagonal matrices and ran simulations in R, but the results aren't what I expected. The simulated values I obtained are not diagonal, so I think I misinterpreted the author's statement. I thought that if you formulate the factor analysis model with the prior distributions above, that you can consider it the classical factor analysis model by choosing certain hyper-parameter value. But it seems that this formulation does not produce the classical factor analysis model.
UPDATE 3: The classical factor analysis model sets ${\bf \Phi} = {\bf I}$ (i.e. non-random), sets $\bf \Psi$ to be a diagonal matrix (i.e. random diagonal matrix) and assigns prior distributions to only the diagonal elements. What I understand the author's statement to mean, is that I can do the aforementioned things by using Wishart priors on $\bf \Phi$ and $\bf \Psi$ with special scale matrices $\bf \Phi_*$ and $\bf \Psi_*$.