The proposed algorithm estimates the covariance matrix with a low-rank decomposition. The rank of the matrix is equal to the number of components \(d_z\)  in the latent space. which is a model hyper-parameter to be determined. We use the Bayesian Information Criterion (BIC) to estimate this parameter. The BIC arises from the Laplace approximation to the model posterior \(p(M | D_k)\)  [1]}, where \(M\)  is the complete model including the latent dimension \(d_z\) . This results in a Bayesian estimate of \(d_z\) : \(d_z = \text{argmax}_{d_z} ( \log p(D_{k}) - \text{BIC}/2)\) , where BIC is a function of the total number of unknown parameters, which penalizes the log likelihood with a model complexity penalty term. In the proposed model, we use ELBO lower bound to the the likelihood by following [2]}. The unknown parameters of the model are \(\lbrace \Theta _{kl}\rbrace _{l=1}^L\) , \(\mathbf {\mu }_k\) , and \(\mathbf {\Sigma _k}\) . Hence, the total number of parameters is \(\text{dof} = K \times (d_{l} + K)\) , which is used in the BIC expression as \(0.5 \times \text{dof} \times \log (I_k) \) . To illustrate the BIC model selection for the proposed model, we simulate three datasets with true latent space dimensions 4, 8, and 12, respectively, and then train multiple models while varying the dimensions \(d_z\)  over \(\lbrace 2,3,\dots ,12\rbrace \)  as the search range. We repeat the experiment 10 times to report the performance. The panel on top of Fig. REF  shows the average BIC values obtained after convergence of the variational EM algorithm. We see that maximum BIC is obtained on the vicinity of the true ranks for all the datasets. On the other hand, the panel on the bottom shows the RMSE values of the estimated covariance matrices. One can see that the lowest error is achieved by maximizing \(d_z\) , which is likely to result in over-fitting the model. This over-fitting is mitigated by using BIC estimate of \(d_z\) .
