# Deep Cervix Model Development from Heterogeneous and Partially Labeled Image
Datasets

Anabik Pal  National Library of Medicine
National Institutes of Health (NIH)
Bethesda, Maryland, USA
<EMAIL_ADDRESS>Zhiyun Xue  National Library of Medicine
National Institutes of Health (NIH)
Bethesda, Maryland, USA
<EMAIL_ADDRESS>Sameer Antani  National Library of Medicine
National Institutes of Health (NIH)
Bethesda, Maryland, USA
<EMAIL_ADDRESS>
###### Abstract

Cervical cancer is the fourth most common cancer in women worldwide. The
availability of a robust automated cervical image classification system can
augment the clinical care provider’s limitation in traditional visual
inspection with acetic acid (VIA). However, there are a wide variety of
cervical inspection objectives which impact the labeling criteria for
criteria-specific prediction model development. Moreover, due to the lack of
confirmatory test results and inter-rater labeling variation, many images are
left unlabeled.

Motivated by these challenges, we propose a self-supervised learning (SSL)
based approach to produce a pre-trained cervix model from unlabeled cervical
images. The developed model is further fine-tuned to produce criteria-specific
classification models with the available labeled images. We demonstrate the
effectiveness of the proposed approach using two cervical image datasets. Both
datasets are partially labeled and labeling criteria are different. The
experimental results show that the SSL-based initialization improves
classification performance (Accuracy: 2.5%↑) and the inclusion of images from
both datasets during SSL further improves the performance (Accuracy: 1.5%↑).
Further, considering data-sharing restrictions, we experimented with the
effectiveness of Federated SSL and find that it can improve performance over
the SSL model developed with just its images. This justifies the importance of
SSL-based cervix model development. We believe that the present research shows
a novel direction in developing criteria-specific custom deep models for
cervical image classification by combining images from different sources
unlabeled and/or labeled with varying criteria, and addressing image access
restrictions.

###### Index Terms:

Self-supervised Learning, Federated Learning, Deep Learning, Cervical Image
classification.

##  I Introduction

Cervical cancer is the fourth common cancer in women worldwide. Regular
cervical screening can help in the early detection of the pre-cancerous
lesions and reduce premature death. Visual inspection with acetic acid (VIA)
is a commonly used low-cost cervical screening approach, but it suffers from
significant inter- and intra-reader variability. In this regard, the
availability of an automated cervical image classification system could
address this limitation [1]. It is challenging to develop a diagnostic system
using traditional classification approaches as it needs enormous effort to
develop effective hand-crafted features and are known to underperform [2, 3,
4]. Hence, the deep learning approaches [5] can be considered to build a
robust classification model. However, the development of a deep model for
robust classification needs a huge number of images and their class labels
[6]. Labeling cervical images is costly, needs multiple experts’ agreement and
requires multiple diagnostic information. Transfer learning [7], i.e.
transferring knowledge from natural images, is a commonly used approach to
overcome these data limitations. However, transferring knowledge from the same
domain may be more effective than transferring knowledge from a different
domain.

Collaboration among multiple organizations to perform centralized learning
(CL) by uniting all labeled images from all sources is another effective
solution to overcome the data scarcity. However, supervised learning using
combined images has limitations. First, different imaging devices may be used
for image acquisition resulting in variations in the visual quality across
different sources. It is challenging to estimate the proper ratio of images to
be selected from various sources for the CL. Moreover, there may be data-
sharing restrictions. To deal with data-sharing restrictions, federated
learning (FL) can be employed [8]. However, federated supervised learning
needs research efforts to address variability in class distribution among the
sources [9]. Furthermore, inclusion of a dataset having noisy labels impede
the training to produce robust models. Finally, and most importantly, the
cervical image labeling criteria vary and depends on: availability of other
diagnostic results, population under study, treatment planning, severity
grading strategy, etc. The variety in the image labeling criteria across
datasets makes the task more challenging as it restricts researchers in
performing any kind of supervised collaborative (CL or FL) learning.

In this paper, we propose a self-supervised learning (SSL) based approach to
develop a pre-trained cervix model (or cervix model). We use two cervical
image datasets which are (i) labeled in a heterogeneous manner: labeling
criteria vary across datasets and (ii) partially labeled: not all images in
the datasets are labeled. As the SSL does not require any label, it allows us
to include all available images in our datasets for cervix model development.
Both centralized SSL and federated SSL are experimented. To evaluate the
effectiveness of the developed cervix model, criteria (or dataset) specific
classification models are trained with the available labeled images. The
classification networks are initialized with the developed cervix models. Note
that, according to our survey, no image dataset is publicly available for
machine learning research towards supporting experts’ effort in visual
assessment of acetic acid applied cervix and no well-accepted classification
network is available for the present task. Hence, we confined our experiments
only with the present datasets for the chosen competing algorithms.

In summary, our work is motivated from [10] and [11]. In [10], for medical
image representation, an encoder-decoder based architecture was used to
reconstruct original images from synthetically distorted images. In [11], a
contrastive learning-based framework is proposed for natural image
representation. In contrast to these works, we utilize the framework presented
in [11] and demonstrate the power of this framework for cervical image
representation tasks in both centralized and federated learning schemes. We
believe that our work has two key novelties: (a) it is the first work
attempting to develop a cervix model from unlabelled images; (b) the first
work where Federated Self-Supervised Learning (FSSL) is demonstrated for any
medical image representation.

The remainder of the paper is organized as follows: Section II discusses the
proposed methodology. The experimental protocol and analysis of experimental
results are presented in Section III and Section IV respectively. Finally,
Section V concludes the paper.

##  II Methods

###  II-A Self-supervised learning

Self-supervised learning (SSL) is a discriminative approach for visual
representation learning. Developers define a pretext task and develop a deep
model which captures the image semantics (i.e. good initialization weights for
related domain’s downstream tasks) with zero labeling cost. In this paper, we
employ a contrastive feature learning algorithm as the pretext task- i.e. all
images will be semantically well separated, and an image and its augmented
version will be semantically closer [11]. During SSL training, every mini-
batch of size $2N$ is constructed with $N$ random images and an augmented
version of them. The training loss ($L_{i,j}$) between an image $i$ and its
augmented version $j$ is given as:

$L_{i,j}=-log\frac{exp(f_{i}^{T}f_{j}/\tau)}{\sum_{i=1}^{2N}1_{[k\neq
i]}exp(f_{i}^{T}f_{k}/\tau)}$ (1)

where $f_{i}$ is the feature vector of $i^{th}$ image; $1_{[k\neq
i]}\in\\{0,1\\}$ is an indicator function evaluating to 1 iff $k\neq i$; T
representation transpose operation and $\tau$ is a constant. The loss in a
mini-batch is computed across all pairs constructed with an image and its
augmented version.

###  II-B Federated Self-supervised Learning (FSSL)

Recently, several different data protection laws and regulations (GDPR 2018 by
EU, CCPA 2020 in the US, etc) have been created to protect information leakage
due to sensitive data sharing (especially medical/banking domain). This
motivates the researchers to develop federated learning (FL) algorithms for
robust inter-institutional collaborated deep development from the data
distributed in multiple institutions (clients) without sharing the raw data
[8].

In this paper, we experiment with Federated Self-supervised Learning (FSSL)
for cervix model development. We perform FSSL in the following two different
ways:- (a) Client-Server FSSL (CSFSSL) and (b) Peer-to-Peer FSSL (PPFSSL). The
CSFSSL contains multiple clients and a server. In this framework, firstly, the
server sends a deep model to all clients. Then independently, at every client,
the model is fine-tuned for $E$-epochs with local data. After that, the
updated models from all clients are aggregated at the server and used as the
initialized model for the next iteration. This procedure continues until the
deep model is trained. Note that, we use weight averaging for model
aggregation. On the other hand, in the PPFSSL, SSL is circularly performed-
firstly the starting client runs SSL with its own data for $E$-epochs and then
sends the trained model to the next client at which SSL is performed for
$E$-epochs using the received model as starting point and then the updated
model is sent to next client and so on. This circular computation is performed
for few iterations and the final model is shared among the clients.

###  II-C Proposed Approach: Cervical Model Development

In this paper, we use centralized SSL (CSSL) and two different variants of
FSSL presented in Sec II-B for cervix model development. The block diagram of
the CSSL, CSFSSL, and PPFSSL are shown in Fig 1 (a), Fig 1 (b) and Fig 1 (c)
respectively. All methods use contrastive feature learning-based SSL presented
in Sec II-A for cervical model development. We perform flip (horizontal and
vertical) axis, rotation, random shift, and random zoom, gamma changing, and
brightness changing to get the augmented version of an image.

Figure 1: Block diagram of the SSL training system- (a) Centralized SSL (CSSL)
(b) Client-Server FSSL (CSFSSL) and (c) Peer-to-Peer FSSL (PPFSSL).

##  III Experimental Protocol

###  III-A Dataset description

National Cancer Institute (NCI) at the US National Institutes of Health (NIH)
conducted two different cohort studies (NHS [12] and ALTS [13]) for cervical
examination. During these studies, in every visit, two images of the acetic
acid applied cervix are captured. For the present research, NCI provides us a
subset of images collected during these studies. The dataset containing NHS
images is referred to as NHS dataset and the dataset containing ALTS images is
referred to as ALTS dataset. NCI scientists labeled a subset of images
available in NHS and ALTS datasets as case (disease) or control (non-disease)
based on the availability of several screening and diagnostic information
(like visual assessment, HPV, cytology, histopathology, colposcopy, etc).
Different criteria were used to label these two datasets as the study
objectives differs111NHS was from general population study, and ALTS was from
triage study in colposcopic clinics.. For our research, we split both NHS and
ALTS datasets at the woman level into three disjoint subsets- train,
validation, and test. The training images are used for model training, the
validation images are used for training hyper-parameter selection, and the
test images are used for classification performance evaluation. The split-wise
number of patients, number of labeled images, and total available images
(including labeled images) for the datasets are given in Table I.

TABLE I: Data set splits. Split | Class | Patients | Labeled Images | Total Images
---|---|---|---|---
NHS | ALTS | NHS | ALTS | NHS | ALTS
Train | Case | 91 | 124 | 182 | 248 | 2029 | 3145
Control | 181 | 242 | 361 | 481
Valid | Case | 22 | 31 | 44 | 62 | 520 | 791
Control | 45 | 60 | 90 | 120
Test | Case | 25 | 34 | 49 | 68 |
Control | 50 | 65 | 99 | 130 |

###  III-B Network architecture

We use ResNet-50 as a backbone network architecture. For SSL approaches,
firstly, the top 1000-way classification layer is removed and a dense layer
with a ReLU activation is added. The dense layer serves as the image
representation vector. We vary the number of neurons among [64, 128, 256] for
this layer and empirically decide to set it $64$ as we obtain very close
performances. For preparing the classification model, we remove the dense
layer and put a single output neuron with sigmoid activation. The output
obtained from the sigmoid layer is the case probability. The classification
network is fully fine-tuned with the prepared cervix model.

###  III-C Competing methods

The paper aims to develop a cervix model, i.e., weight initialization method
for cervical image classification network. We compare among six different
approaches (i) Random: network weights are randomly initialized. (ii)
ImageNet: network weights are taken from the pre-trained ImageNet
classification model- knowledge is transferred from natural images. (iii)
Self-supervised Learning (SSL): SSL is employed with the available images in a
dataset. (iv) Centralized Self-supervised Learning (CSSL): All images from
both datasets are combined to train the SSL as shown in Fig 1(a). (v) Client-
Server Federated Self-supervised Learning (CSFSSL): Images are not shared,
client SSL models are aggregated in the server (see Fig 1(b)). (vi) Peer-to-
Peer Federated Self-supervised Learning (PPFSSL): Image are not shared,
circularly SSL based fine-tuning is performed among the clients (see Fig
1(c)).

###  III-D Parameter Settings

We vary the network hyper-parameters for both SSL and classification model
training and choose the best hyper-parameters based on the validation loss.
The SSL networks are trained with following parameters: learning rate $0.01$,
weight decay = $1e-5$, momentum = $0.9$, $\tau=0.1$ (See Eq 1), epochs = $50$
and the classification networks are trained with following hyper-parameters:
learning rate $0.001$, weight decay = $1e-6$, momentum = $0.09$, epochs =
$50$, batch size =4. In both SSL and classification model training, we
randomly shuffle images during batch construction. The learning stops when the
validation loss is not decreasing for $5$ epochs. We use reverse class
weighting to address the class imbalance issue in classification model
training.

###  III-E Implementation

The Keras [14] deep learning tool-kit is used for implementing the networks.
The networks are trained with 2 GeForce RTX 2080 Ti GPUs installed with an
Intel(R) Xeon(R) Gold 5218 CPU (@ 2.30GHz). We implement federated learning in
the same computing resources as sequential processes.

###  III-F Evaluation metrics

In this paper, we evaluate the performance of the dataset-specific
classification models produced from the supervised learning with the varying
initialization approaches mentioned in Sec III-C. The following four
quantitative evaluation metrics are computed for performance comparison: (1)
Accuracy (ACC), (2) Recall, (3) Precision, (4) F1-Score.

##  IV Experimental results and discussion

(a) (b)

Figure 2: Receiver Operating Curve (ROC): (a) NHS and (b) ALTS. The numeric
values represent the AUC values for the classifiers with considered
initialization.

The Receiver Operating Curves (ROC) for all supervised classifiers for both
NHS and ALTS datasets developed with varying initialization approaches are
shown in Fig 2. The quantitative classification performance (performed with
the metrics mentioned in Sec III-C) for the same is shown in Table II.
According to Table II, we find that ImageNet initialization improves the
accuracy for both datasets but recalls did not improve. The cervical image
classification model initialization with SSL based approach is more effective
than the ImageNet weights. In the SSL approach, for the NHS dataset, best
accuracy, recall, and F1$\\_$Score are obtained when $N=16$; and best
precision is received when $N=8$. In the SSL approach for ALTS dataset, best
accuracy and precision are obtained when $N=16$; best recall and F1$\\_$Score
are obtained when $N=8$. According to our experimental results, we observe a
noticeable performance improvement in using SSL. Thus the experimental results
justify the importance of SSL-based cervical model development. The CSSL based
initialization approach further improves the performance of SSL for both
datasets. This observation supports the development of the cervix model
uniting images from both datasets. We find that for CSSL, in general, $N=8$
produces the best performance in both NHS and ALTS. Hence, we evaluate the
performance of both Federated Self-Supervised Learning (FSSL) models for
$N=8$. For CSFSSL, we vary the value of $E$ i.e. local model updating for $1$,
$5$, $10$ and find that $E=1$ provides us the best result which is reported in
the ninth row of Table II. For PPFSSL, we opt two approaches. The first
approach starts SSL training with NHS (called PPFSSL$\\_$NHS) images and the
other approach starts SSL training with ALTS (called PPFSSL$\\_$ALTS) images.
We set $E=1$ for experimental similarity, and the classification performance
with this model initialization is listed in the tenth and eleventh rows of
Table II. We find that both PPFSSL and CSFSSL produce comparative performance
and in general, federated SSL produces better results than SSL with images in
its own dataset. This justifies the effectiveness of FSSL in addressing data
sharing constraints in cervix model development.

TABLE II: Performance evaluation

Initialization | NHS | ALTS
---|---|---
Method | ACC | Recall | Precision | F1_Score | ACC | Recall | Precision | F1_Score
Random | 79.73 | 0.5306 | 0.7879 | 0.6341 | 76.77 | 0.5735 | 0.6964 | 0.6290
ImageNet | 80.41 | 0.4898 | 0.8571 | 0.6234 | 77.78 | 0.5588 | 0.7308 | 0.6333
SSL (N 8) | 83.11 | 0.6327 | 0.8158 | 0.7126 | 79.80 | 0.7206 | 0.7000 | 0.7101
SSL (N 16) | 83.11 | 0.6939 | 0.7727 | 0.7312 | 80.30 | 0.6029 | 0.7736 | 0.6777
SSL (N 32) | 83.11 | 0.6531 | 0.8000 | 0.7191 | 79.29 | 0.6176 | 0.7368 | 0.6720
CSSL (N 8) | 86.49 | 0.7143 | 0.8537 | 0.7778 | 81.82 | 0.7794 | 0.7162 | 0.7465
CSSL (N 16) | 85.81 | 0.7551 | 0.8043 | 0.7789 | 81.31 | 0.6324 | 0.7818 | 0.6992
CSSL (N 32) | 85.14 | 0.7959 | 0.7647 | 0.7800 | 81.31 | 0.6471 | 0.7719 | 0.7040
CSFSSL | 84.46 | 0.6735 | 0.8250 | 0.7416 | 79.80 | 0.6912 | 0.7121 | 0.7015
PPFSSL$\\_$NHS | 84.46 | 0.7347 | 0.7826 | 0.7579 | 80.81 | 0.6618 | 0.7500 | 0.7031
PPFSSL$\\_$ALTS | 84.46 | 0.7755 | 0.7600 | 0.7677 | 80.30 | 0.6618 | 0.7377 | 0.6977

##  V Conclusion and scope of future work

This paper discusses the challenges behind cervical image analysis due to
labeling unavailability and variability and presents a novel direction to
address them. Experimental results shows that the self-supervised learning
algorithm is proved to be an efficient and effective candidate to deal with
the label scarcity as well as the labeling variability. In addition, the
experimentations on Federated Self-Supervised Learning shed into light to deal
with the data-sharing restrictions. To the best of our knowledge, this is the
first attempt to develop a cervix model in light of a domain-specific pre-
trained model for task-specific fine-tuning.

The development of an improved cervix model by uniting larger datasets from
different sources independent of labeling criteria, unavailability of labels,
and different imaging devices is the immediate future scope of this work. The
cervical image datasets for which the labeling efforts made by researchers at
NCI are in progress will be used for this research. The engineering
implementation of federated learning to work in a real scenario is another
important future work. The presented idea can be employed by other medical
image analysis tasks for utilizing unlabeled data and providing the same
domain transfer learning.

##  VI Acknowledgement

We are very much grateful to Dr. Mark Schiffman of the Division of Cancer
Epidemiology and Genetics, National Cancer Institute, National Institutes of
Health, and his team for providing us the images and labels used in this
paper.

## References

  * [1] L. Hu, D. Bell, S. Antani, Z. Xue, K. Yu, M. P. Horning, N. Gachuhi, B. Wilson, M. S. Jaiswal, B. Befano, L. R. Long, R. Herrero, M. H. Einstein, R. D. Burk, M. Demarco, J. C. Gage, A. C. Rodriguez, N. Wentzensen, and M. Schiffman, “An Observational Study of Deep Learning and Automated Evaluation of Cervical Images for Cancer Screening,” _JNCI: Journal of the National Cancer Institute_ , vol. 111, no. 9, pp. 923–932, 01 2019.
  * [2] Y. Srinivasan, B. Nutter, S. Mitra, B. Phillips, and E. Sinzinger, “Classification of cervix lesions using filter bank-based texture mode,” in _19th IEEE Symposium on Computer-Based Medical Systems (CBMS’06)_ , 2006, pp. 832–840.
  * [3] E. Kim and X. Huang, _A Data Driven Approach to Cervigram Image Analysis and Classification_. Dordrecht: Springer Netherlands, 2013, pp. 1–13.
  * [4] K. Fernandes, J. S. Cardoso, and J. Fernandes, “Automated methods for the decision support of cervical cancer screening using digital colposcopies,” _IEEE Access_ , vol. 6, pp. 33 910–33 927, 2018.
  * [5] G. Litjens, T. Kooi, B. E. Bejnordi, A. A. A. Setio, F. Ciompi, M. Ghafoorian, J. A. van der Laak, B. van Ginneken, and C. I. Sánchez, “A survey on deep learning in medical image analysis,” _Medical Image Analysis_ , vol. 42, pp. 60 – 88, 2017. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S1361841517301135
  * [6] T. Ching, D. S. Himmelstein, B. K. Beaulieu-Jones, A. A. Kalinin, B. T. Do, G. P. Way, E. Ferrero, P.-M. Agapow, M. Zietz, M. M. Hoffman, W. Xie, G. L. Rosen, B. J. Lengerich, J. Israeli, J. Lanchantin, S. Woloszynek, A. E. Carpenter, A. Shrikumar, J. Xu, E. M. Cofer, C. A. Lavender, S. C. Turaga, A. M. Alexandari, Z. Lu, D. J. Harris, D. DeCaprio, Y. Qi, A. Kundaje, Y. Peng, L. K. Wiley, M. H. S. Segler, S. M. Boca, S. J. Swamidass, A. Huang, A. Gitter, and C. S. Greene, “Opportunities and obstacles for deep learning in biology and medicine.” _Journal of the Royal Society, Interface_ , April, 2018.
  * [7] Q. Yang, Y. Zhang, W. Dai, and S. J. Pan, _Transfer Learning_. Cambridge University Press, 2020.
  * [8] Q. Yang, Y. Liu, Y. Cheng, Y. Kang, T. Chen, and H. Yu, 2019.
  * [9] M. Yang, A. Wong, H. Zhu, H. Wang, and H. Qian, “Federated learning with class imbalance reduction,” 2020.
  * [10] Z. Zhou, V. Sodha, J. Pang, M. B. Gotway, and J. Liang, “Models genesis,” _Medical Image Analysis_ , vol. 67, p. 101840, 2021.
  * [11] T. Chen, S. Kornblith, M. Norouzi, and G. Hinton, “A simple framework for contrastive learning of visual representations,” in _Proceedings of the 37th International Conference on Machine Learning_ , ser. Proceedings of Machine Learning Research, H. D. III and A. Singh, Eds., vol. 119. Virtual: PMLR, 13–18 Jul 2020, pp. 1597–1607.
  * [12] M. C. Bratti, A. C. Rodríguez, M. Schiffman, A. Hildesheim, J. Morales, M. Alfaro, D. Guillén, M. Hutchinson, M. E. Sherman, C. Eklund, J. Schussler, J. Buckland, L. A Morera, F. Cárdenas, M. Barrantes, E. Pérez, T. J Cox, R. D Burk, and R. Herrero, “Description of a seven-year prospective study of human papillomavirus infection and cervical neoplasia among 10000 women in guanacaste, costa rica,” _Pan American journal of public health_ , vol. 2, no. 15, pp. 75–89, Feb 2004.
  * [13] M. Schiffman and M. E. Adrianza, “Ascus-lsil triage study. design, methods and characteristics of trial participants,” _Acta Cytol_ , vol. 44(5), pp. 726–742, Sept-Oct 2000.
  * [14] F. Chollet _et al._ , “Keras,” https://keras.io, 2015.


the right hand side with respect to $\sigma^{2}_{i}$ yields the marginal
likelihood under the fractional prior:

$\displaystyle
p(\tilde{{\mathbf{y}}}_{i}|\boldsymbol{\delta}_{r},{\mathbf{f}}_{r})=\frac{b^{q_{i}/2}\Gamma(c_{T})(C_{i0})^{c_{0}}}{(2\pi)^{T(1-b)/2}\Gamma(c_{0})(C_{iT}^{\boldsymbol{\delta}})^{c_{T}}}.$
(B.4)

###  B.2 Multimove sampling of a set of indicators in a column

Another important building block of MCMC inference for sparse Bayesian factor
models is sampling all indicators $\delta_{ij}$ in column ${j}$ for a set of
rows $i\in I_{j}\subseteq\\{1,\ldots,m\\}$, conditional on the factors
${\mathbf{f}}_{r}=({\mathbf{f}}^{r}_{1},\ldots,{\mathbf{f}}^{r}_{T})$, the
remaining indicator $\boldsymbol{\delta}_{r}^{u}$ and the hyperparameter
$\tau_{j}$, without conditioning on the model parameters
$\boldsymbol{\beta}_{r}$ and $\sigma^{2}_{1},\ldots,\sigma^{2}_{m}$, see Step
(D) of Algorithm 1.

According to the prior (3.2), the indicators $\delta_{ij}$ are independent
apriori conditional on the hyperparameter $\tau_{j}$, with the log prior odds
$O^{\mbox{\rm\tiny pr}}_{ij}$ of $\delta_{ij}=1$ versus $\delta_{ij}=0$ being
given by:

$\displaystyle O^{\mbox{\rm\tiny pr}}_{ij}=\log\frac{\mbox{\rm
Pr}(\delta_{ij}=1|\tau_{j})}{\mbox{\rm
Pr}(\delta_{ij}=0|\tau_{j})}=\log\frac{\tau_{j}}{1-\tau_{j}}.$ (B.5)

To sample $\delta_{ij}$ conditional $\boldsymbol{\delta}_{r}^{u}$ and
${\mathbf{f}}_{r}$, without conditioning on $\boldsymbol{\beta}_{r}$ and
$(\sigma^{2}_{1},\ldots,\sigma^{2}_{m})$, the log posterior odds
$O^{\mbox{\rm\tiny post}}_{ij}$, given by

$\displaystyle O^{\mbox{\rm\tiny post}}_{ij}$ $\displaystyle=$
$\displaystyle\log\frac{\mbox{\rm
Pr}(\delta_{ij}=1|\boldsymbol{\delta}_{r}^{u},\tau_{j},\tilde{{\mathbf{y}}}_{i},{\mathbf{f}}_{r})}{\mbox{\rm
Pr}(\delta_{ij}=0|\boldsymbol{\delta}_{r}^{u},\tau_{j},\tilde{{\mathbf{y}}}_{i},{\mathbf{f}}_{r})}=\log\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}+\log\frac{\mbox{\rm
Pr}(\delta_{ij}=1|\tau_{j})}{\mbox{\rm Pr}(\delta_{ij}=0|\tau_{j})}$ (B.6)
$\displaystyle=$ $\displaystyle O_{ij}+O^{\mbox{\rm\tiny pr}}_{ij},$

is required which combines the log prior odds $O^{\mbox{\rm\tiny pr}}_{ij}$
given in (B.5) with the log likelihood ratio $O_{ij}$, given by:

$\displaystyle
O_{ij}=\log\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}.$
(B.7)

The likelihood ratio $O_{ij}$ is easily computed from the marginal likelihoods
$p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij},\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})$
where, respectively, $\delta_{ij}=1$ and $\delta_{ij}=0$. As discussed in
Section B.1, these marginal likelihoods are available in closed form and
marginal likelihood computation can be done individually for each row $i\in
I_{j}$, separately for $\delta_{ij}=0$ and $\delta_{ij}=1$.

However, this procedure is likely to be inefficient, in particular, if the set
$I_{j}$ is large. To achieve greater efficiency, Algorithm 4 outlined below
provides a technique to compute directly the log likelihood ratio $O_{ij}$
(rather than the individual marginal likelihoods) simultaneously for all rows
$i\in I_{j}$. This allows joint sampling of all indicators $\delta_{ij}$ in
column $j$ for all rows $i\in I_{j}$.

The precise form of the log likelihood ratio $O_{ij}$ of $\delta_{ij}=1$
versus $\delta_{ij}=0$ defined in (B.7) depends on the remaining indicators
$\boldsymbol{\delta}_{i,-j}$ in row $i$. The computation of $O_{ij}$ is easily
vectorized for all rows $i\in I_{j}$ where all elements of
$\boldsymbol{\delta}_{i,-j}$ are zero. In this case, a model where observation
$y_{it}$ is dedicated to factor $j$ ($\delta_{ij}=1$) is compared to a model
where $y_{it}$ is uncorrelated with all remaining observations
($\delta_{ij}=0$). In this case, $O_{ij}$ is easily obtained from the marginal
likelihood of a dedicated model with ${j_{i}}=j$ and the “null” model. As
shown in Algorithm 4, it is possible (but less straightforward) to vectorize
the computation of the log likelihood ratio $O_{ij}$ also for the remaining
rows $i\in I_{j}$ where at least one element of $\boldsymbol{\delta}_{i,-j}$
is different from zero.

######  Algorithm 4.

Multimove sampling of indicators in a column. Sample all indicators
$\delta_{ij}$ in column $\boldsymbol{\delta}_{\cdot,j}$ jointly for all rows
$i\in I_{j}\subseteq\\{1,\ldots,m\\}$ conditional on the factors
${\mathbf{f}}_{r}=({\mathbf{f}}^{r}_{1},\ldots,{\mathbf{f}}^{r}_{T})$, the
remaining indicators $\boldsymbol{\delta}_{r}^{u}$ and the hyperparameter
$\tau_{j}$, without conditioning on the model parameters
$\boldsymbol{\beta}_{r}$ and $\sigma^{2}_{1},\ldots,\sigma^{2}_{m}$ using the
following steps:

  * (I-a)

Compute the log likelihood ratio $O_{ij}$ for all rows $i\in I_{j}$ where all
elements of $\boldsymbol{\delta}_{r}^{u}$ are zero as

$\displaystyle
O_{ij}=\log\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}=\log\frac{\Gamma(c_{T})(C_{iT}^{n})^{c_{T}^{n}}}{\Gamma(c_{T}^{n})(C_{iT})^{c_{T}}}+D_{ij}.$
(B.8)

$c_{T}^{n}$ and $C_{iT}^{n}$ are the posterior moments of the null model given
in (A.3). $c_{T}$ and $C_{iT}$ are the posterior moments of $\sigma_{i}^{2}$
for a dedicated measurement with ${j_{i}}=j$, given in (A.8) and (A.9),
respectively, for a hierarchical Gaussian and a fractional prior. For a
hierarchical Gaussian prior, $D_{ij}=0.5\log(B_{iT}/B_{i0,jj})$, where
$B_{i0,jj}$ is $j$th diagonal element of the prior variance
${\mathbf{B}}_{i0}$ and $B_{iT}$ is the posterior scale factor for a dedicated
measurement with ${j_{i}}=j$, given in (A.8). For a fractional prior,
$D_{ij}=0.5\log(b(2\pi)^{bT})$.

This step is trivial to vectorize.

  * (I-b)

For all rows $i\in\\{i_{1},\ldots,i_{n}\\}\subseteq I_{j}$ where
$\boldsymbol{\delta}_{i,-j}$ is non-zero, compute

$\displaystyle
O_{ij}=\log\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}=c_{T}\log\frac{C_{iT}^{0}}{C_{iT}^{1}}+D_{ij},$
(B.9)

where $c_{T}$ and $C_{iT}^{\delta_{ij}}$ are the posterior moments of
$\sigma^{2}_{i}|\delta_{ij},\cdot$ given in (A.4) and $C_{iT}^{0}$ refers to a
model with $\delta_{ij}=0$, while $C_{iT}^{1}$ refer to a model with
$\delta_{ij}=1$. For a hierarchical Gaussian prior,

$\displaystyle
D_{ij}=0.5\log(|{\mathbf{B}}_{iT}^{1}|/|{\mathbf{B}}_{iT}^{0}|)-0.5\log(|{\mathbf{B}}_{i0}^{1}|/|{\mathbf{B}}_{i0}^{0}|),$
(B.10)

where ${\mathbf{B}}_{i0}^{\delta_{ij}}$ and ${\mathbf{B}}_{iT}^{\delta_{ij}}$
refer to the prior and posterior moments of
$\boldsymbol{\beta}_{i\cdot}^{\boldsymbol{\delta}}|\delta_{ij},\cdot$ given in
(A.5). ${\mathbf{B}}_{i0}^{1}$ and ${\mathbf{B}}_{iT}^{1}$ refer to the prior
and posterior moments for a model where $\delta_{ij}=1$, while
${\mathbf{B}}_{i0}^{0}$ and ${\mathbf{B}}_{iT}^{0}$ refer to the prior and
posterior moments for a model where $\delta_{ij}=0$. For a fractional prior,
$D_{ij}\equiv 0.5\log b$.

Use Algorithm 5 to determine $C_{iT}^{1}$, $C_{iT}^{0}$, as well as $D_{ij}$
for a hierarchical Gaussian prior, simultaneously for all rows
$i\in\\{i_{1},\ldots,i_{n}\\}\subseteq I_{j}$.

  * (I-c)

Determine the vector of the log posterior odds $O^{\mbox{\rm\tiny
post}}_{ij}=O_{ij}+O^{\mbox{\rm\tiny pr}}_{ij}$ for all rows $i\in I_{j}$.
Joint sampling of $\delta_{ij}|\tau_{j},\cdot$ is easily vectorized:

    * (I-c1)

Propose $\delta_{ij}^{\mbox{\rm\tiny new}}=1-\delta_{ij}$ for $i\in I_{j}$.

    * (I-c2)

Draw a vector of $|I_{j}|$ random variables
$U_{i}\sim\mathcal{U}\left[0,1\right]$, indexed by $i\in I_{j}$.

    * (I-c3)

For all rows $i\in I_{j}$, where $\delta_{ij}=0$, accept the proposal
$\delta_{ij}^{\mbox{\rm\tiny new}}=1$, iff $\log U_{i}\leq O^{\mbox{\rm\tiny
post}}_{ij}$;

    * (I-c4)

For all rows $i\in I_{j}$, where $\delta_{ij}=1$, accept the proposal
$\delta_{ij}^{\mbox{\rm\tiny new}}=0$, iff $\log U_{i}\leq-O^{\mbox{\rm\tiny
post}}_{ij}$.

Using, respectively, (B.2) and (B.4), the expression for $O_{ij}$ in (B.9) is
easily derived. Since the indicators in column $j$ are independent given
$\tau_{j}$, Step (I-c) is based on $|I_{j}|$ independent Metropolis-Hastings
(MH) steps each of which proposes to update $\delta_{ij}$ by flipping the
indicator, i.e. $\delta_{ij}^{\mbox{\rm\tiny
new}}=1-\delta_{ij}$.101010Alternatively, a Gibbs step may be used, i.e. set
$\delta_{ij}^{\mbox{\rm\tiny new}}=1$, iff $\log(U_{i}/(1-U_{i}))\leq
O^{\mbox{\rm\tiny post}}_{ij}$, otherwise $\delta_{ij}^{\mbox{\rm\tiny
new}}=0$. However, simulation experiments indicate that the MH step is more
efficient. It easy to verify that the acceptance rules formulated in Step
(I-c3) and (I-c4) are equivalent to the more convential rule to accept
$\delta_{ij}^{\mbox{\rm\tiny new}}$ with probability

$\min\left\\{1,\frac{\mbox{\rm Pr}(\delta_{ij}^{\mbox{\rm\tiny
new}}|\boldsymbol{\delta}_{r}^{u},\tau_{j},\tilde{{\mathbf{y}}}_{i},{\mathbf{f}}_{r})}{\mbox{\rm
Pr}(\delta_{ij}|\boldsymbol{\delta}_{r}^{u},\tau_{j},\tilde{{\mathbf{y}}}_{i},{\mathbf{f}}_{r})}\right\\}=\min\left\\{1,\exp(O^{\mbox{\rm\tiny
post}}_{ij})\right\\}.$

######  Algorithm 5.

To compute all relevant posterior moments in (B.9) simultaneously for all rows
$\\{i_{1},\ldots,i_{n}\\}$, proceed as follows:

  * (a)

Set the indicator $\delta_{i_{l},j}=1$ in each row
$i_{l}\in\\{i_{1},\ldots,i_{n}\\}$. Reorder the columns of the factor loading
matrix in such a way, that the $j$th column appears last. This is simply done
by permuting the column of ${\mathbf{F}}$ appropriately before defining
${\mathbf{X}}_{i_{l}}^{\boldsymbol{\delta}}$.111111While the fractional prior
is not affected by this, it might be necessary to reorder the prior mean and
the prior covariance matrix for a hierarchical Gaussian prior.

  * (b)

Set up the information matrix ${\mathbf{P}}$ and the covector ${\mathbf{m}}$
of the corresponding joint posterior of all nonzero factor loadings in the
rows $i_{1},\ldots,i_{n}$ as described in Algorithm 3. Compute the Cholesky
decomposition ${\mathbf{L}}$ of ${\mathbf{P}}$ and the corresponding vector
${\mathbf{x}}$ solving ${\mathbf{L}}{\mathbf{x}}={\mathbf{m}}$.

  * (c)

Knowing ${\mathbf{L}}$ and ${\mathbf{x}}$, a vectorized computation of the log
likelihood ratio (B.8) for all rows $i_{l}\in\\{i_{1},\ldots,i_{n}\\}$ is
possible. The posterior moments $C_{i_{l},T}^{1}$ are directly available from
the appropriate sub vectors ${\mathbf{x}}_{i_{l}}$ of ${\mathbf{x}}$, defined
in (A.11). When we switch from $\delta_{i_{l},j}=1$ to a model where
$\delta_{i_{l},j}=0$, then for a hierarchical Gaussian prior,

$\displaystyle
C_{i_{l},T}^{0}=C_{i_{l},T}^{1}+\frac{1}{2}(x^{\star}_{i_{l}})^{2}.$ (B.11)

Furthermore,

$\displaystyle
0.5\log(|{\mathbf{B}}_{i_{l},T}^{1}|/|{\mathbf{B}}_{i_{l},T}^{0}|)=-\log
L^{\star}_{i_{l}},$ (B.12)

where $L^{\star}_{i_{l}}=({\mathbf{L}}_{i})_{q_{i_{l}},q_{i_{l}}}$ is the last
diagonal element of the submatrix ${\mathbf{L}}_{i_{l}}$. Therefore,

$D_{ij}=-\log L^{\star}_{i_{l}}-0.5\log B_{i0,jj}.$

For the fractional prior,

$\displaystyle
C_{i_{l},T}^{0}=C_{i_{l},T}^{1}+\frac{1-b}{2}(x^{\star}_{i_{l}})^{2},$ (B.13)

where $x^{\star}_{i_{l}}=({\mathbf{x}}_{i_{l}})_{q_{i_{l}}}$ is the last
element of ${\mathbf{x}}_{i_{l}}$.

Derivation of Step (c). When we switch from a model where all indicator
$\delta_{i_{1},j}=\ldots=\delta_{i_{n},j}=1$ are equal to one to a model where
all indicators $\delta_{i_{1},j}=\ldots=\delta_{i_{n},j}=0$ are zero, then the
information matrix ${\mathbf{P}}^{0}$ and the covector ${\mathbf{m}}^{0}$ of
the joint posterior of the remaining nonzero factor loadings is obtained from
${\mathbf{P}}$ and ${\mathbf{m}}$ simply by deleting all rows and columns
corresponding to $\delta_{i_{1},j},\ldots,\delta_{i_{n},j}$, and the Cholesky
decomposition ${\mathbf{L}}^{0}$ of ${\mathbf{P}}^{0}$ is obtained from
${\mathbf{L}}$ in the same way. Also the vector ${\mathbf{x}}^{0}$ solving
${\mathbf{L}}^{0}{\mathbf{x}}^{0}={\mathbf{m}}^{0}$ is obtained from
${\mathbf{x}}$ simply by deleting the rows corresponding to
$\delta_{i_{1},j},\ldots,\delta_{i_{n},j}$. This last result is easily seen by
considering the subsystem
${\mathbf{L}}_{i_{l}}{\mathbf{x}}_{i_{l}}={\mathbf{m}}_{i_{l},T}^{\boldsymbol{\delta}}$
corresponding to the $i_{l}$th row. Because

$\displaystyle{\mathbf{L}}_{i_{l}}=\left(\begin{array}[]{cc}{\mathbf{L}}_{i_{l}}^{0}&{\mathbf{O}}\\\
\mathbf{l}_{i_{l}}&({\mathbf{L}}_{i})_{q_{i_{l}},q_{i_{l}}}\\\
\end{array}\right)=\left(\begin{array}[]{cc}{\mathbf{L}}_{i_{l}}^{0}&{\mathbf{O}}\\\
\mathbf{l}_{i_{l}}&L^{\star}_{i_{l}}\\\ \end{array}\right),\qquad$ (B.18)

we obtain
${\mathbf{L}}_{i_{l}}^{0}{\mathbf{x}}_{i_{l}}^{0}={\mathbf{m}}^{0}_{i_{l}}$,
where ${\mathbf{x}}_{i_{l}}^{0}$ is obtained from ${\mathbf{x}}_{i_{l}}$ by
deleting the $q_{i_{l}}$th element
$x^{\star}_{i_{l}}=({\mathbf{x}}_{i_{l}})_{q_{i_{l}}}$. Hence,
${\mathbf{x}}_{i_{l}}^{0}$ defines the desired subvector of ${\mathbf{x}}^{0}$
to compute $C_{i_{l},T}^{0}$ as in (A.11). Since
$({\mathbf{x}}_{{i_{l}}}^{0})^{{}^{\prime}}{\mathbf{x}}_{i_{l}}^{0}={\mathbf{x}}_{i_{l}}^{{}^{\prime}}{\mathbf{x}}_{i_{l}}-({\mathbf{x}}_{i_{l}})_{q_{i_{l}}}^{2}$
we obtain from (A.4) that (B.13) holds. Note, however, that this simple
relationship would not hold without reordering the columns as described above.

Finally, to compute the log likelihood ratio for a hierarchical Gaussian
prior, the ratio of the determinants
$|{\mathbf{B}}^{1}_{i_{l},T}|/|{\mathbf{B}}_{i_{l},T}^{0}|$ is required. Since
the lower triangular matrices ${\mathbf{L}}_{i_{l}}$ and
${\mathbf{L}}_{i_{l}}^{0}$ are, respectively, the Cholesky decomposition of
$({\mathbf{B}}_{i_{l},T}^{1})^{-1}$ and $({\mathbf{B}}_{i_{l},T}^{0})^{-1}$,
we obtain:

$\displaystyle
1/|{\mathbf{B}}_{i_{l},T}^{1}|^{1/2}=|({\mathbf{B}}_{i_{l},T}^{1})^{-1}|^{1/2}=|{\mathbf{L}}_{i_{l}}|,$
(B.19)

where $|{\mathbf{L}}_{i_{l}}|$ is the product of the diagonal elements of
${\mathbf{L}}_{i_{l}}$. Computing $|{\mathbf{B}}_{i_{l},T}^{0}|$ in the same
way and using (B.18) proves (B.12).

##  Appendix C Updating shrinkage parameters in Step (S)

Under a hierarchical Gaussian prior on the factor loadings $\beta_{ij}$, Steps
(P), (D) and (L) are performed conditional on all hyperparameters of this
prior. All local, column specific and global scaling parameters in the prior
of the loadings $\beta_{ij}$ are updated in Step (S) of Algorithm 1. This step
relies on a representation of an F-distribution as a Gamma scale mixture of
inverse gamma distributions and avoids the GIG-distribution which was used in
Cadonna et al. (2020).121212A r.v. $X\sim\mbox{\rm F}\left(2a,2c\right)$ has a
representation as a Gamma scale mixture of inverse gamma distributions:
$\displaystyle X|b\sim\mathcal{G}^{-1}\left(c,b\right),\quad
b\sim\mathcal{G}\left(a,\frac{a}{c}\right).$ (C.1) Hence,  $\displaystyle
b|X\sim\mathcal{G}\left(a+c,\frac{a}{c}+\frac{1}{X}\right).$ (C.2) If
$Y|X\sim\mathcal{N}\left(0,dX\right)$, then  $\displaystyle
X|b,Y\sim\mathcal{G}^{-1}\left(c+\frac{1}{2},b+\frac{Y^{2}}{2d}\right).$ (C.3)

#### Step (S-a).

Under prior (3.15), sample the local shrinkage parameters $\omega_{ij}$ using
the representation of the F-distribution as a scale mixture of inverse gamma
distributions:

$\displaystyle\omega_{ij}|b^{\omega}_{ij}\sim\mathcal{G}^{-1}\left(c^{\omega},b^{\omega}_{ij}\right),\quad
b^{\omega}_{ij}|\omega_{ij}\sim\mathcal{G}\left(a^{\omega},\frac{a^{\omega}}{c^{\omega}}\right),$

This yields a two step sampler, where $b^{\omega}_{ij}|\omega_{ij}$ is imputed
from

$\displaystyle
b^{\omega}_{ij}|\omega_{ij}\sim\mathcal{G}\left(a^{\omega}+c^{\omega},\frac{a^{\omega}}{c^{\omega}}+\frac{1}{\omega_{ij}}\right),$
(C.4)

and $\omega_{ij}|b^{\omega}_{ij},{\beta}_{ij},\delta_{ij}$ given
$b^{\omega}_{ij}$ is sampled from

$\displaystyle\omega_{ij}|b^{\omega}_{ij},{\beta}_{ij},\delta_{ij}\sim\mathcal{G}^{-1}\left(c^{\omega}+\delta_{ij}\frac{1}{2},b^{\omega}_{ij}+\delta_{ij}\frac{{\beta}_{ij}^{2}}{2\kappa\theta_{j}\sigma^{2}_{i}}\right).$
(C.5)

This has to be done for the entire matrix $\boldsymbol{\omega}_{r}$, as
$\omega_{ij}$ is needed to compute the odd likelihood ratio $O_{ij}$ of
$\delta_{ij}=1$ versus $\delta_{ij}=0$ (see (B.7)) also above the current
pivot.

#### Step (S-b).

Sample the column specific shrinkage parameters $\boldsymbol{\theta}_{r}$ (if
any). For
$\theta_{j}\sim\mathcal{G}^{-1}\left(c^{\theta},b^{\theta}_{j}\right)$, this
yields for all $j=1$ to $r$:

$\displaystyle\theta_{j}|\kappa,\boldsymbol{\beta}_{r},b^{\theta}_{j}\sim\mathcal{G}^{-1}\left(c^{\theta}+\frac{d_{j}}{2},b^{\theta}_{j}+\frac{1}{2\kappa}\sum_{i:\delta_{ij}=1}^{m}\frac{{\beta}_{ij}^{2}}{\sigma^{2}_{i}\omega_{ij}}\right).$
(C.6)

For $\theta_{j}\sim\mbox{\rm F}\left(2a^{\theta},2c^{\theta}\right)$, we use
again the representation of the F-distribution as a scale mixture of inverse
gamma distributions:

$\displaystyle\theta_{j}|b^{\theta}_{j}\sim\mathcal{G}^{-1}\left(c^{\theta},b^{\theta}_{j}\right),\quad
b^{\theta}_{j}\sim\mathcal{G}\left(a^{\theta},\frac{a^{\theta}}{c^{\theta}}\right),$

and impute $b^{\theta}_{j}|\theta_{j}$ from

$\displaystyle
b^{\theta}_{j}|\theta_{j}\sim\mathcal{G}\left(a^{\theta}+c^{\theta},\frac{a^{\theta}}{c^{\theta}}+\frac{1}{\theta_{j}}\right),$
(C.7)

before we update $\theta_{j}|\kappa,\boldsymbol{\beta}_{r}$ from (C.6).

#### Step (S-c).

If $\kappa$ is random with prior
$\kappa\sim\mathcal{G}^{-1}\left(c^{\kappa},b^{\kappa}\right)$, then
$\kappa|{\boldsymbol{\beta}}_{r},{\mathbf{\Sigma}}_{r},\boldsymbol{\omega}_{r}$
is updated from

$\displaystyle\kappa|{\boldsymbol{\beta}}_{r},{\mathbf{\Sigma}}_{r},\boldsymbol{\omega}_{r}\sim\mathcal{G}^{-1}\left(c^{\kappa}+\frac{d}{2},b^{\kappa}+1/2\sum_{i=1}^{m}\frac{1}{\sigma^{2}_{i}}\sum_{j=1}^{r}\delta_{ij}\frac{{\beta}_{ij}^{2}}{\theta_{j}\omega_{ij}}\right),$
(C.8)

where $d=\sum_{i=1}^{m}\sum_{j=1}^{r}\delta_{ij}$ is the total number of non-
zero loadings in the CFA model. If $\kappa\sim\mbox{\rm
F}\left(2a^{\kappa},2c^{\kappa}\right)$, then using the same representation of
the F-distribution, we impute $b^{\kappa}|\kappa$ from

$\displaystyle
b^{\kappa}|\kappa\sim\mathcal{G}\left(a^{\kappa}+c^{\kappa},\frac{a^{\theta}}{c^{\theta}}+\frac{1}{\kappa}\right),$
(C.9)

before we update
$\kappa|{\boldsymbol{\beta}}_{r},{\mathbf{\Sigma}}_{r},\boldsymbol{\omega}_{r}$
from (C.8).

##  Appendix D Details on split and merge in Step (R)

#### Proof of (4.3).

For any zero column $j_{\footnotesize\rm sp}$ in a sparsity matrix
$\boldsymbol{\delta}_{k}$ with $r$ active and $r_{\footnotesize sp}$ spurious
columns, there are $m-(r+r_{\footnotesize sp})$ unrestricted elements that can
be subjected to variable selection and take the value one with probability
$\tau_{j_{\footnotesize\rm sp}}|k\sim\mathcal{B}\left(a_{k},b_{k}\right)$
following the prior (3.2). Hence,

$\displaystyle\mbox{\rm Pr}(d_{j_{\footnotesize\rm sp}}=0|r,r_{\footnotesize
sp})=\int(1-\tau_{j_{\footnotesize\rm sp}})^{m-(r+r_{\footnotesize
sp})}p(\tau_{j_{\footnotesize\rm sp}})d\,\tau_{j_{\footnotesize\rm sp}}=$
(D.1) $\displaystyle\displaystyle\frac{B(a_{k},b_{k}+m-(r+r_{\footnotesize
sp}))}{B(a_{k},b_{k})}$

whereas

$\displaystyle\mbox{\rm Pr}(d_{j_{\footnotesize\rm sp}}=1|r,r_{\footnotesize
sp})=\displaystyle(m-(r+r_{\footnotesize sp}))\int\tau_{j_{\footnotesize\rm
sp}}(1-\tau_{j_{\footnotesize\rm sp}})^{m-(r+r_{\footnotesize
sp})-1}p(\tau_{j_{\footnotesize\rm sp}})d\,\tau_{j_{\footnotesize\rm sp}}=$
$\displaystyle\displaystyle\frac{(m-(r+r_{\footnotesize
sp}))B(a_{k}+1,b_{k}+m-(r+r_{\footnotesize sp})-1)}{B(a_{k},b_{k})}.$ (D.2)

Therefore,

$\displaystyle\frac{\mbox{\rm Pr}(d_{j_{\footnotesize\rm
sp}}=1|r,r_{\footnotesize sp})}{\mbox{\rm Pr}(d_{j_{\footnotesize\rm
sp}}=0|r,r_{\footnotesize sp})}=\frac{a_{k}(m-r-r_{\footnotesize
sp})}{b_{k}+m-r-r_{\footnotesize sp}-1}.$

#### Proposing factors in a spurious column - Proof of (4.2).

Whenever a zero column is turned into a spurious column (with columns index
$j_{\footnotesize\rm sp}$), factors ${\mathbf{f}}_{j_{\footnotesize\rm
sp}}=(f_{j_{\footnotesize\rm sp},1},\ldots,f_{j_{\footnotesize\rm sp},T})$ are
proposed, while holding the factors ${\mathbf{f}}^{r}_{t},t=1,\ldots,T,$ in
all active columns fixed. Draws of $f_{j_{\footnotesize\rm sp},t}$ could be
proposed from the prior $f_{j_{\footnotesize\rm
sp},t}\sim\mathcal{N}\left(0,1\right)$, since column $j_{\footnotesize\rm sp}$
was a zero column before splitting. However, with $l_{\footnotesize\rm sp}$
being the pivot row, $y_{l_{\footnotesize\rm sp},t}$ is a measurement that
contains information about $f_{j_{\footnotesize\rm sp},t}$ in a spurious
column and its likelihood can be combined with the prior to define the
conditional posterior density $p(f_{j_{\footnotesize\rm
sp},t}|{\mathbf{f}}^{r}_{t},\boldsymbol{\beta}_{r},\sigma^{2}_{l_{\footnotesize\rm
sp}},y_{l_{\footnotesize\rm sp},t})$ given $y_{l_{\footnotesize\rm sp},t}$ in
addition to
${\mathbf{f}}^{r}_{t},\boldsymbol{\beta}_{r},\sigma^{2}_{l_{\footnotesize\rm
sp}},y_{l_{\footnotesize\rm sp},t}$.

It is easy to verify from the filter given in (4.1) that for a spurious column
$j_{\footnotesize\rm sp}$ with leading element $\Xi_{l_{\footnotesize\rm
sp}}$, the conditional density $p(f_{j_{\footnotesize\rm
sp},t}|{\mathbf{f}}^{r}_{t},\boldsymbol{\beta}_{r},\sigma^{2}_{l_{\footnotesize\rm
sp}},y_{l_{\footnotesize\rm sp},t})$ of $f_{j_{\footnotesize\rm sp},t}$ is
given by:

$\displaystyle f_{j_{\footnotesize\rm
sp},t}|{\mathbf{f}}^{r}_{t},\boldsymbol{\beta}_{r},\sigma^{2}_{l_{\footnotesize\rm
sp}},y_{l_{\footnotesize\rm sp},t}\sim\mathcal{N}\left(E_{j_{\footnotesize\rm
sp},t},V_{j_{\footnotesize\rm sp}}\right),$ $\displaystyle\displaystyle
V_{j_{\footnotesize\rm sp}}=\left(1+\frac{\Xi_{l_{\footnotesize\rm
sp}}^{2}}{(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny
new}}}\right)^{-1}=\frac{(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}}{(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}+\Xi_{l_{\footnotesize\rm sp}}^{2}},\quad
E_{j_{\footnotesize\rm sp},t}=\frac{V_{j_{\footnotesize\rm
sp}}\Xi_{l_{\footnotesize\rm sp}}}{(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}}\tilde{y}_{l_{\footnotesize\rm
sp},t}=\frac{\Xi_{l_{\footnotesize\rm sp}}}{(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}+\Xi_{l_{\footnotesize\rm
sp}}^{2}}\tilde{y}_{l_{\footnotesize\rm sp},t},$

where the pseudo outcome $\tilde{y}_{l_{\footnotesize\rm sp},t}$ is given by
$\tilde{y}_{l_{\footnotesize\rm sp},t}=y_{l_{\footnotesize\rm
sp},t}-\boldsymbol{\beta}_{r,l_{\footnotesize\rm
sp},\cdot}{\mathbf{f}}^{r}_{t}$. Using the definition of
$\Xi_{l_{\footnotesize\rm sp}}$ and $(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}$ in terms of $\sigma^{2}_{l_{\footnotesize\rm
sp}}$ and a uniform random variable $U_{j_{\footnotesize\rm sp}}\in[-1,1]$
given in (4.6), we obtain the simple expressions for the posterior moments
given in (4.2):

$\displaystyle V_{j_{\footnotesize\rm
sp}}=\frac{(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny
new}}}{(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny
new}}+\Xi_{l_{\footnotesize\rm sp}}^{2}}=1-U_{j_{\footnotesize\rm
sp}}^{2},\quad E_{j_{\footnotesize\rm sp},t}=\frac{\Xi_{l_{\footnotesize\rm
sp}}}{(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny
new}}+\Xi_{l_{\footnotesize\rm sp}}^{2}}\tilde{y}_{l_{\footnotesize\rm
sp},t}=\frac{U_{j_{\footnotesize\rm
sp}}}{\sqrt{(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny
new}}}}\tilde{y}_{l_{\footnotesize\rm sp},t}.$

#### Turning spurious into active columns.

In Step (R-D), variable selection is performed on the spurious columns of the
current EFA model by marginalizing over the idiosyncratic variances and the
factor loading matrix $(\boldsymbol{\beta}_{r},\,\boldsymbol{\Xi})$ as in Step
(D) of Algorithm 1. By starting with the last spurious column, we ensure that
the application of Step (D) is valid. Indeed, ordering the spurious pivots in
Step (L) by size guarantees that all rows of the current EFA model below
$l_{\footnotesize\rm sp}$ are equal to the corresponding rows of the current
CFA model in the first $r$ columns and zero in all other columns. Hence,
updating the zero column $r+1$ below $l_{\footnotesize\rm sp}$ can be regarded
as an attempt to introduce an additional active column in the current CFA
model and to increase $r$, while $r_{\footnotesize sp}$ decreases. If a column
remains spurious, then we integrate in the current EFA model immediately over
the corresponding factor, while the current number of spurious columns
$r_{\footnotesize sp}$ is preserved.

#### Step (R-Da) for hierarchical Gaussian priors with unknown shrinkage
factors.

Step (R-Da) has to be adjusted slightly for hierarchical Gaussian priors with
unknown shrinkage factors. A column specific shrinkage parameter
$\theta_{j_{\footnotesize\rm sp}}$ is needed to compute the odds of turning
zero loadings in the spurious column into non-zero ones. In addition, under
the shrinkage prior (3.15), local shrinkage parameters
$\omega_{i,j_{\footnotesize\rm sp}}$ are needed for all rows $i=1,\ldots,m$ in
the spurious column $j_{\footnotesize\rm sp}$ (note that
$\omega_{i,j_{\footnotesize\rm sp}}$ will be needed later in step (L), even if
$i<l_{\footnotesize\rm sp}$ lies above the spurious pivot element). All these
quantities could be simply simulated from the prior. However, we found
considerably gain in exploited the information contained in
$\Xi_{l_{\footnotesize\rm sp}}$ and $(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}$, when sampling $\theta_{j_{\footnotesize\rm sp}}$
and $\omega_{l_{\footnotesize\rm sp},j_{\footnotesize\rm sp}}$, which
basically results from conditioning on the r.v. $U_{j_{\footnotesize\rm sp}}$.

For the Student-$t$ slab prior (3.14), the information contained in
$\Xi_{l_{\footnotesize\rm sp}}$ and $(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}}$ is exploited by sampling
$\theta_{j_{\footnotesize\rm sp}}|\Xi_{l_{\footnotesize\rm
sp}},(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny new}}$ from (C.6):

$\displaystyle\theta_{j_{\footnotesize\rm sp}}|\kappa,U_{j_{\footnotesize\rm
sp}}\sim\mathcal{G}^{-1}\left(c^{\theta}+\frac{1}{2},b^{\theta}_{j_{\footnotesize\rm
sp}}+\frac{U_{j_{\footnotesize\rm sp}}^{2}}{2\kappa(1-U_{j_{\footnotesize\rm
sp}}^{2})}\right),$ (D.3)

where the scale parameter is obtained from combining (4.6) and (C.6):

$\displaystyle\frac{\beta_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}}^{2}}{({\mathbf{\Sigma}}_{r+1})_{l_{\footnotesize\rm
sp},l_{\footnotesize\rm sp}}}=\frac{\Xi_{l_{\footnotesize\rm
sp}}^{2}}{(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny
new}}}=\frac{U_{j_{\footnotesize\rm sp}}^{2}}{1-U_{j_{\footnotesize\rm
sp}}^{2}}.$

Under the prior $\theta_{j_{\footnotesize\rm
sp}}\sim\mathcal{G}^{-1}\left(c^{\theta},b^{\theta}\right)$, the scale
$b^{\theta}_{j_{\footnotesize\rm sp}}=b^{\theta}$, whereas under the F-prior
on $\theta_{j_{\footnotesize\rm sp}}$, $b^{\theta}_{j_{\footnotesize\rm sp}}$
is imputed from the prior. This is achieved by sampling
$\theta_{j_{\footnotesize\rm sp}}\sim\mbox{\rm
F}\left(2a^{\theta},2c^{\theta}\right)$ and $b^{\theta}_{j_{\footnotesize\rm
sp}}|\theta_{j_{\footnotesize\rm sp}}$ from (C.7):

$\displaystyle b^{\theta}_{j_{\footnotesize\rm
sp}}|\theta_{j_{\footnotesize\rm
sp}}\sim\mathcal{G}\left(a^{\theta}+c^{\theta},\frac{a^{\theta}}{c^{\theta}}+\frac{1}{\theta_{j_{\footnotesize\rm
sp}}}\right).$

For the shrinkage prior (3.15), first local shrinkage parameter
$\omega_{i,j_{\footnotesize\rm sp}}$ are sampled from the prior for all rows,
i.e. $\omega_{i,j_{\footnotesize\rm sp}}\sim\mbox{\rm
F}\left(2a^{\omega},2c^{\omega}\right)$, $i=1,\ldots,m$. The information
contained in $\Xi_{l_{\footnotesize\rm sp}}$ and
$(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny new}}$ is then
exploited by sampling $\theta_{j_{\footnotesize\rm
sp}}|\Xi_{l_{\footnotesize\rm sp}},(\sigma^{2}_{l_{\footnotesize\rm
sp}})^{\mbox{\rm\tiny new}},\omega_{l_{\footnotesize\rm
sp},j_{\footnotesize\rm sp}}$ from (C.6):

$\displaystyle\theta_{j_{\footnotesize\rm sp}}|\kappa,U_{j_{\footnotesize\rm
sp}},\omega_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}}\sim\mathcal{G}^{-1}\left(c^{\theta}+\frac{1}{2},b^{\theta}_{j_{\footnotesize\rm
sp}}+\frac{U_{j_{\footnotesize\rm sp}}^{2}}{2\kappa\omega_{l_{\footnotesize\rm
sp},j_{\footnotesize\rm sp}}(1-U_{j_{\footnotesize\rm sp}}^{2})}\right),$

where $b^{\theta}_{j_{\footnotesize\rm sp}}$ is defined as for the Student-$t$
slab prior (3.14). Finally, we resample $\omega_{l_{\footnotesize\rm
sp},j_{\footnotesize\rm sp}}$ by first imputing
$b^{\omega}_{l_{\footnotesize\rm sp},j_{\footnotesize\rm sp}}$ from the prior
by sampling $b^{\omega}_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}}|\omega_{l_{\footnotesize\rm sp},j_{\footnotesize\rm sp}}$ from (C.4):

$\displaystyle b^{\omega}_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}}|\omega_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}}\sim\mathcal{G}\left(a^{\omega}+c^{\omega},\frac{a^{\omega}}{c^{\omega}}+\frac{1}{\omega_{l_{\footnotesize\rm
sp},j_{\footnotesize\rm sp}}}\right),$

(where $\omega_{l_{\footnotesize\rm sp},j_{\footnotesize\rm sp}}$ has been
sampled from the prior in the first step) and then sampling from (C.5):

$\displaystyle\omega_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}}|b^{\omega}_{l_{\footnotesize\rm sp},j_{\footnotesize\rm
sp}},\theta_{j_{\footnotesize\rm sp}},\kappa,U_{j_{\footnotesize\rm
sp}}\sim\mathcal{G}^{-1}\left(c^{\omega}+\frac{1}{2},b^{\omega}_{l_{\footnotesize\rm
sp},j_{\footnotesize\rm sp}}+\frac{U_{j_{\footnotesize\rm
sp}}^{2}}{2\kappa\theta_{j_{\footnotesize\rm sp}}(1-U_{j_{\footnotesize\rm
sp}}^{2})}\right).$

In this way, we exploit the information contained in $\Xi_{l_{\footnotesize\rm
sp}}$ and $(\sigma^{2}_{l_{\footnotesize\rm sp}})^{\mbox{\rm\tiny new}}$ for
both shrinkage parameters.

##  Appendix E Step (L) - Updating the pivots

This subsection provides details on Step (L) of Algorithm 1, which was briefly
discussed in Section 4.2. Four local moves are applied which were illustrated
in Figure 1 in Section 4.2.

#### Shifting the pivot.

A shift move is selected with probability $p_{\mbox{\rm\footnotesize shift}}$.
Let $l_{\star}$ denote the index of the first nonzero row below $l_{j}$ (i.e.
$\delta_{l_{\star},j}=1$, $\delta_{ij}=0,l_{j}<i<l_{\star}$). If
$l_{\star}>2$, then a new pivot $l^{\mbox{\rm\tiny new}}_{j}$ is proposed by
sampling uniformly from the set ${\cal
M}(l_{\star},{\mathbf{l}}_{r,-j})=\\{1,\ldots,l_{\star}-1\\}\cap{\cal
L}({\mathbf{l}}_{r,-j})$, where ${\cal L}({\mathbf{l}}_{r,-j})$ is the set of
pivots outside of column $j$. If ${\cal M}(l_{\star},{\mathbf{l}}_{r,-j})$ is
empty, then no shift move is performed. Otherwise, two indicators in column
$j$ are changed, namely $\delta_{l_{j}^{\mbox{\rm\tiny new}},j}$ from zero to
one and $\delta^{\mbox{\rm\tiny new}}_{l_{j},j}$ from one to zero, while the
remaining elements of $\boldsymbol{\delta}_{r}$ are unchanged. The new
indicator matrix $\boldsymbol{\delta}_{r}^{\mbox{\rm\tiny new}}$ is accepted
with probability $\min(1,\alpha_{\mbox{\rm\footnotesize shift}})$, where

$\displaystyle\alpha_{\mbox{\rm\footnotesize
shift}}=\exp(O_{l_{j}^{\mbox{\rm\tiny
new}},j}-O_{l_{j},j})R_{\mbox{\rm\footnotesize shift}},$

with $O_{ij}$ being the log likelihood ratio of $\delta_{ij}=1$ versus
$\delta_{ij}=0$ defined in (B.7). A shift move does not change the number of
nonzero elements $d_{j}$ in column $j$ and the prior ratio (4.7) of this move
simplifies to:

$\displaystyle R_{\mbox{\rm\footnotesize
shift}}=\frac{B(a_{k}+d_{j}-1,b_{k}+m-l_{j}^{\mbox{\rm\tiny
new}}-d_{j}+1)}{B(a_{k}+d_{j}-1,b_{k}+m-l_{j}-d_{j}+1)}.$ (E.1)

Since $l^{\mbox{\rm\tiny new}}_{j}$ is sampled from a set ${\cal
M}(l_{\star},{\mathbf{l}}_{r,-j})$ that does not depend on $l_{j}$, the
proposal density is symmetric and the proposal ratio $q(l^{\mbox{\rm\tiny
new}}_{j}|l_{\star},{\mathbf{l}}_{r,-j})/q(l_{j}|l_{\star},{\mathbf{l}}_{r,-j})$
cancels from $\alpha_{\mbox{\rm\footnotesize shift}}$.

#### Switching pivots.

This move is selected with probability $p_{\mbox{\rm\footnotesize switch}}$
and is performed only if $r>1$. A nonzero column $\ell\neq j$ is selected
randomly and all indicators between (and including) row $l_{j}$ and $l_{\ell}$
that are different are switched between the two columns, i.e.
$\delta_{ij}^{\mbox{\rm\tiny new}}=1-\delta_{ij}$ and
$\delta_{i\ell}^{\mbox{\rm\tiny new}}=1-\delta_{i\ell}$ for all $i\in{\cal
S}_{j,\ell}=\\{i:\min(l_{\ell},l_{j})\leq
i\leq\max(l_{\ell},l_{j}),\delta_{ij}\neq\delta_{i\ell}\\}$. Evidently, this
move switches pivots between the two columns. Since the corresponding proposal
density satisfies $q(\boldsymbol{\delta}_{r}^{\mbox{\rm\tiny
new}}|\boldsymbol{\delta}_{r})=q(\boldsymbol{\delta}_{r}|\boldsymbol{\delta}_{r}^{\mbox{\rm\tiny
new}})$, $\boldsymbol{\delta}_{r}^{\mbox{\rm\tiny new}}$ is accepted with
probability $\min(1,\alpha_{\mbox{\rm\footnotesize switch}})$, where

$\displaystyle\alpha_{\mbox{\rm\footnotesize switch}}=\prod_{i\in{\cal
S}_{j,\ell}}\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}},\delta_{i\ell}^{\mbox{\rm\tiny
new}},\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij},\delta_{i\ell},\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}R_{\mbox{\rm\footnotesize
switch}},$

where the prior odds ratio of this move is derived from (4.7):

$\displaystyle R_{\mbox{\rm\footnotesize
switch}}=\prod_{{\tilde{l}}=j,l}\frac{B(a_{k}+d_{\tilde{l}}^{\mbox{\rm\tiny
new}}-1,b_{k}+m-l_{\tilde{l}}^{\mbox{\rm\tiny
new}}-d_{\tilde{l}}^{\mbox{\rm\tiny
new}}+1)}{B(a_{k}+d_{\tilde{l}}-1,b_{k}+m-l_{\tilde{l}}-d_{\tilde{l}}+1)}.$
(E.2)

If $\delta_{ij}=0$ (and consequently $\delta_{i\ell}=1$), then we obtain:

$\displaystyle\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}}=1,\delta_{i\ell}^{\mbox{\rm\tiny
new}}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\delta_{i\ell}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}=$
$\displaystyle\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}}=1,\delta_{i\ell}^{\mbox{\rm\tiny
new}}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\delta_{i\ell}^{\mbox{\rm\tiny
new}}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\delta_{i\ell}^{\mbox{\rm\tiny
new}}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=0,\delta_{i\ell}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}=\exp(O_{ij|\ell}-O_{i\ell|j}),$

where $O_{i,j_{1}|j_{2}}$ is the log likelihood ratio of $\delta_{i,j_{1}}=1$
versus $\delta_{i,j_{1}}=0$ provided that the indicator $\delta_{i,j_{2}}=0$.
It can be obtained as the likelihood ratio $O_{i,j_{1}}$ given in (B.7), with
$\delta_{i,j_{2}}=0$ for both models. On the other hand, if $\delta_{ij}=1$
(and consequently $\delta_{i\ell}=0$), then

$\displaystyle\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}}=0,\delta_{i\ell}^{\mbox{\rm\tiny
new}}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=1,\delta_{i\ell}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}=$
$\displaystyle\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}}=0,\delta_{i\ell}^{\mbox{\rm\tiny
new}}=1,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}}=0,\delta_{i\ell}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}\frac{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}^{\mbox{\rm\tiny
new}}=0,\delta_{i\ell}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}{p(\tilde{{\mathbf{y}}}_{i}|\delta_{ij}=1,\delta_{i\ell}=0,\boldsymbol{\delta}_{r}^{u},{\mathbf{f}}_{r})}=\exp(O_{i\ell|j}-O_{ij|\ell}).$

Therefore

$\displaystyle\alpha_{\mbox{\rm\footnotesize
switch}}=R_{\mbox{\rm\footnotesize switch}}\exp\left(\sum_{i\in{\cal
S}_{j,\ell}:\delta_{ij}=0}(O_{ij|\ell}-O_{i\ell|j})+\sum_{i\in{\cal
S}_{j,\ell}:\delta_{ij}=1}(O_{i\ell|j}-O_{ij|\ell})\right).$ (E.3)

This move allows changes in $d_{j}$ and $d_{\ell}$, but leaves the overall
number $d$ of nonzero elements unchanged.

#### Adding or deleting a pivot.

Finally, a reversible pair of moves is selected with probability
$1-p_{\mbox{\rm\footnotesize shift}}-p_{\mbox{\rm\footnotesize switch}}$. The
add move introduces a new pivot $l_{j}^{\mbox{\rm\tiny new}}$ in a row above
the current pivot $l_{j}$ which is not occupied by the pivots of the other
columns. Hence, $l^{\mbox{\rm\tiny new}}_{j}$ is selected randomly from the
set ${\cal A}(l_{j},{\mathbf{l}}_{r,-j})=\\{1,\ldots,l_{j}-1\\}\cap{\cal
L}({\mathbf{l}}_{r,-j})$, i.e. $\delta_{l_{j}^{\mbox{\rm\tiny
new}},j}^{\mbox{\rm\tiny new}}=1$, while the remaining elements of
$\boldsymbol{\delta}_{r}$ are unchanged (in particular
$\delta_{l_{j},j}^{\mbox{\rm\tiny new}}=\delta_{l_{j},j}=1$). An add move is
only possible, if $|{\cal A}(l_{j},{\mathbf{l}}_{r,-j})|>0$.131313The number
of rows in ${\cal A}(l_{j},{\mathbf{l}}_{r,-j})$ is equal to $l_{j}-z_{j}$,
where $z_{j}=\\#\\{l_{j^{\prime}}\in{\mathbf{l}}:l_{j^{\prime}}<l_{j}\\}$ is
the rank of $l_{j}$ among the leading indices. Hence, an add move is possible,
whenever $l_{j}>z_{j}$.

The corresponding reverse move is deterministic and deletes the current pivot
$l_{j}$, turning $l_{j}^{\mbox{\rm\tiny new}}=l_{\star}$ into the new pivot
where $l_{\star}$ is the row index of the first nonzero element in
$\boldsymbol{\delta}_{r}$ below $l_{j}$. Hence,
$\delta_{l_{j},j}^{\mbox{\rm\tiny new}}=0$, while all other elements of
$\boldsymbol{\delta}_{r}$ remain unchanged. A delete move is only performed,
if $l_{\star}$ is not a pivot in any other column of $\boldsymbol{\delta}_{r}$
(that is $l_{\star}\in{\cal L}({\mathbf{l}}_{r,-j})$).

If neither an add nor a delete move is possible, then $l_{j}$ remains
unchanged. Otherwise, either an add or a delete move is selected with
probability $p_{\mbox{\rm\footnotesize add}}(\boldsymbol{\delta}_{r})$ and
$1-p_{\mbox{\rm\footnotesize add}}(\boldsymbol{\delta}_{r})$. If both add and
delete moves are possible, then $p_{\mbox{\rm\footnotesize
add}}(\boldsymbol{\delta}_{r})=p_{a}$, with $p_{a}$ being a tuning parameter;
if only an add move is possible, then $p_{\mbox{\rm\footnotesize
add}}(\boldsymbol{\delta}_{r})=1$, whereas $p_{\mbox{\rm\footnotesize
add}}(\boldsymbol{\delta}_{r})=0$, if only a delete move is possible. Note
that whenever an add move is selected, the reverse delete move is always
possible; and similarly, the reverse add move is always possible, whenever a
delete move is selected. This move changes $d_{j}$ and increases or decreases
the overall number of nonzero elements by one.

The acceptance probability for an add move is equal to
$\min(1,\alpha_{\mbox{\rm\footnotesize add}})$, with

$\displaystyle\alpha_{\mbox{\rm\footnotesize
add}}=\exp(O_{l_{j}^{\mbox{\rm\tiny new}},j})R_{\mbox{\rm\footnotesize
add}}\frac{|{\cal A}(l_{j},{\mathbf{l}}_{r,-j})|(1-p_{\mbox{\rm\footnotesize
add}}(\boldsymbol{\delta}_{r}^{\mbox{\rm\tiny
new}}))}{p_{\mbox{\rm\footnotesize add}}(\boldsymbol{\delta}_{r})},$

where $O_{ij}$ is the log likelihood ratio of $\delta_{ij}=1$ versus
$\delta_{ij}=0$ defined in (B.7). Since $d_{j}^{\mbox{\rm\tiny new}}=d_{j}+1$,
the prior ratio (4.7) of this move simplifies to:

$\displaystyle R_{\mbox{\rm\footnotesize
add}}=\frac{B(a_{k}+d_{j},b_{k}+m-l_{j}^{\mbox{\rm\tiny
new}}-d_{j})}{B(a_{k}+d_{j}-1,b_{k}+m-l_{j}-d_{j}+1)}.$ (E.4)

The acceptance probability for a delete move is equal to
$\min(1,\alpha_{\mbox{\rm\footnotesize del}})$, with

$\displaystyle\alpha_{\mbox{\rm\footnotesize
del}}=\exp(-O_{l_{j},j})R_{\mbox{\rm\footnotesize
del}}\frac{p_{\mbox{\rm\footnotesize
add}}(\boldsymbol{\delta}_{r}^{\mbox{\rm\tiny new}})}{|{\cal
A}(l_{j}^{\mbox{\rm\tiny
new}},{\mathbf{l}}_{r,-j})|(1-p_{\mbox{\rm\footnotesize
add}}(\boldsymbol{\delta}_{r}))}.$

Since $d_{j}^{\mbox{\rm\tiny new}}=d_{j}-1$, the prior ratio (4.7) of this
move simplifies to:

$\displaystyle R_{\mbox{\rm\footnotesize
del}}=\frac{B(a_{k}+d_{j}-2,b_{k}+m-l_{j}^{\mbox{\rm\tiny
new}}-d_{j}+2)}{B(a_{k}+d_{j}-1,b_{k}+m-l_{j}-d_{j}+1)}.$ (E.5)

#### Tuning parameters.

These four moves involve three tuning probabilities, namely
$p_{\mbox{\rm\footnotesize shift}}$, $p_{\mbox{\rm\footnotesize switch}}$, and
$p_{a}$, with $1-p_{\mbox{\rm\footnotesize shift}}-p_{\mbox{\rm\footnotesize
switch}}>0$ and $0<p_{a}<1$.

##  Appendix F Details on Step (H)

In this section, we provide details on updating the hyperparameters $\alpha$
and $\gamma$ for the 2PB prior (3.4) without conditioning on the slab
probabilities $\tau_{1},\ldots,\tau_{k}$. Inference is performed conditionally
on the current columns sizes $d_{1},\ldots,d_{k}$ and the number of spurious
columns $r_{\footnotesize sp}$ in the EFA model and the current pivots in the
CFA model:

$\displaystyle p(\alpha,\gamma|\cdot)\propto
p(\alpha)p(\gamma)\prod_{j:d_{j}>1}\mbox{\rm
Pr}(\boldsymbol{\delta}^{r}_{\cdot,j}|\alpha,\gamma,l_{j})\,\prod_{j:d_{j}\leq
1}\mbox{\rm Pr}(d_{j}|\alpha,\gamma),$

where $\alpha\sim\mathcal{G}\left(a^{\alpha},b^{\alpha}\right)$ and
$\gamma\sim\mathcal{G}\left(a^{\gamma},b^{\gamma}\right)$.

The conditional posterior $p(\alpha,\gamma|\cdot)$ is based on the likelihood
(3.11) for active columns, the likelihood (D.2) for spurious columns and the
likelihood (D.1) for zero columns:141414Note that in (D.2), $r_{\footnotesize
sp}$ is the number of spurious columns excluding column $j_{\footnotesize\rm
sp}$. Hence, the likelihood is computed sequentially from the first to the
last spurious column.

$\displaystyle\displaystyle
p(\alpha,\gamma|\cdot)\propto\frac{p(\alpha)p(\gamma)}{B(\frac{\alpha\gamma}{k},\gamma)^{k}}\cdot
B(\frac{\alpha\gamma}{k},\gamma+m-r-r_{\footnotesize sp}))^{k-r-
r_{\footnotesize sp}}$
$\displaystyle\displaystyle\cdot\prod_{j:d_{j}>1}B(\frac{\alpha\gamma}{k}+d_{j}-1,\gamma+m-l_{j}-d_{j}+1)\cdot\prod_{j_{\footnotesize\rm
sp}=1}^{r_{\footnotesize sp}}B(\frac{\alpha\gamma}{k}+1,\gamma+m-r-
j_{\footnotesize\rm sp}).$

An MH-step based on a random walk proposal for, respectively, $\log\alpha$ and
$\log\gamma$ is implemented to sample $p(\alpha|\gamma,\cdot)$ and
$p(\gamma|\alpha,\cdot)$,.

##  Appendix G More about boosting in Step (A)

Algorithm 6 described a generic boosting strategy in Step (A) in Algorithm 1.

######  Algorithm 6 (Boosting MCMC).

Step (A) in Algorithm 1 is implemented in three steps:

  * (A-a)

Choose a (current) value $\boldsymbol{\Psi}$ and move from the CFA model (2.8)
to the expanded model (4.8) using the following transformation:

$\displaystyle\tilde{\mathbf{f}}^{r}_{t}=(\boldsymbol{\Psi})^{1/2}{\mathbf{f}}^{r}_{t},\quad\tilde{\boldsymbol{\beta}}_{r}=\boldsymbol{\beta}_{r}(\boldsymbol{\Psi})^{-1/2}.$
(G.1)

  * (A-b)

Sample a new value $\boldsymbol{\Psi}^{\mbox{\rm\tiny new}}$ in model (4.8)
conditional on
$\tilde{{\mathbf{f}}}_{r}=(\tilde{\mathbf{f}}^{r}_{1},\ldots,\tilde{\mathbf{f}}^{r}_{T})$
and $\tilde{\boldsymbol{\beta}}_{r}$ from the conditional posterior
$p(\boldsymbol{\Psi}|\tilde{{\mathbf{f}}}_{r},\tilde{\boldsymbol{\beta}}_{r})$
given by

$\displaystyle
p(\boldsymbol{\Psi}|\tilde{{\mathbf{f}}}_{r},\tilde{\boldsymbol{\beta}}_{r})\propto
p(\boldsymbol{\Psi})p(\tilde{\boldsymbol{\beta}}_{r}|\boldsymbol{\Psi})\prod_{j=1}^{r}\Psi_{j}^{-T/2}\exp\left\\{-\frac{1}{2\Psi_{j}}\sum_{t=1}^{T}\tilde{f}_{jt}^{2}\right\\}.$
(G.2)

  * (A-c)

Move back from model (4.8) to the CFA model (2.8) using the inverse of
transformation (G.1) for $j=1,\ldots,r$:

$\displaystyle\beta_{ij}^{\mbox{\rm\tiny
new}}=\beta_{ij}\sqrt{\Psi_{j}^{\mbox{\rm\tiny new}}/\Psi_{j}},\quad
i=1,\ldots,m,\quad f_{jt}^{\mbox{\rm\tiny
new}}=f_{jt}\sqrt{\Psi_{j}/\Psi_{j}^{\mbox{\rm\tiny new}}}\quad t=1,\ldots,T.$

Boosting affects the factor loading matrix $\boldsymbol{\beta}_{r}$ and all
factors ${\mathbf{f}}_{r}=({\mathbf{f}}^{r}_{1},\ldots,{\mathbf{f}}^{r}_{T})$,
see Step (A-c). The main difference between ASIS and MDA lies in the choice of
the current value of $\boldsymbol{\Psi}$ in Step (A-a). While $\Psi_{j}$ is
chosen in a deterministic fashion for ASIS, typically involving a specific
factor loading $\beta_{n_{j},j}$ in each column, $\Psi_{j}$ is sampled from a
working prior for MDA. This leads to different priors $p(\boldsymbol{\Psi})$
and $p(\tilde{\boldsymbol{\beta}}_{r}|\boldsymbol{\Psi})$ in (G.2).

#### ASIS for fractional priors.

For boosting based on ASIS, a nonzero factor loading $\beta_{n_{j},j}$ is
chosen in column $j$ to define the current value of $\Psi_{j}$ as
$\sqrt{\Psi_{j}}=\beta_{n_{j},j}$. This creates a factor loading matrix
$\tilde{\boldsymbol{\beta}}_{r}$ in the expanded model (4.8) where
$\tilde{\beta}_{n_{j},j}=1$ whereas
$\tilde{\beta}_{i,j}=\beta_{ij}/\beta_{n_{j},j}$ for $i\neq n_{j}$.
$\beta_{n_{j},j}$ can be chosen as the pivot element in each column, i.e.
$n_{j}=l_{j}$, or such that $|\beta_{n_{j},j}|$ is maximized for all loadings
in column $j$. Apart from this choice, ASIS requires no further tuning. The
implied prior of $\Psi_{j}=\beta_{n_{j},j}^{2}$ is given by
$p(\Psi_{j})\propto\Psi_{j}^{-1/2}$, whereas
$p(\tilde{\boldsymbol{\beta}}_{r}|\boldsymbol{\Psi})\propto\prod_{j=1}^{r}\Psi_{j}^{(d_{j}-1)/2}$,
since
$p(\tilde{\beta_{ij}}|\sigma^{2}_{i},\Psi_{j},\delta_{ij}=1)\propto\Psi_{j}^{1/2}$
for all $i\neq n_{j}$ with $\delta_{ij}=1$. Hence,
$p(\boldsymbol{\Psi}|\tilde{{\mathbf{f}}}_{r},\tilde{\boldsymbol{\beta}}_{r})$
defined in (G.2) reduces to a product of inverted Gamma distribution where for
all $j=1,\ldots,r$:

$\displaystyle\Psi_{j}^{\mbox{\rm\tiny
new}}|{\mathbf{f}}_{r},\boldsymbol{\beta}_{r}\sim\mathcal{G}^{-1}\left(\frac{T-d_{j}}{2},\frac{\beta_{n_{j},j}^{2}}{2}\sum_{t=1}^{T}f_{jt}^{2}\right),$

where $d_{j}=\sum_{i=1}^{m}\delta_{ij}$ counts the nonzero elements in column
$j$ (which is at least equal to 2 by definition of $\boldsymbol{\delta}_{r}$).

#### MDA for fractional priors.

In MDA, the current value $\boldsymbol{\Psi}$ is drawn from a working prior
$p(\boldsymbol{\Psi})$ which is independent both of $\boldsymbol{\beta}_{r}$
and ${\mathbf{\Sigma}}_{r}$. This guarantees that the prior distribution of
$\boldsymbol{\beta}_{r}$ remains unchanged, despite moving between the two
models. Below, an inverted Gamma working prior
$\Psi_{j}\sim\mathcal{G}^{-1}\left(\nu_{j},q_{j}\right)$ is applied. The prior
in the expanded model reads
$p(\tilde{\boldsymbol{\beta}}_{r}|\boldsymbol{\Psi})\propto\prod_{j=1}^{r}\Psi_{j}^{d_{j}/2}$,
since
$p(\tilde{\beta_{ij}}|\sigma^{2}_{i},\Psi_{j},\delta_{ij}=1)\propto\Psi_{j}^{1/2}$
for all $i$ with $\delta_{ij}=1$. The two last terms in (G.2) factor into a
product of independent inverted Gamma distributions which leads to an inverted
Gamma posterior for each $\Psi_{j}$:

$\displaystyle\Psi_{j}^{\mbox{\rm\tiny
new}}|\Psi_{j},{\mathbf{f}}_{r}\sim\mathcal{G}^{-1}\left(\nu_{j}-d_{j}/2+\frac{T}{2},q_{j}+\Psi_{j}/2\sum_{t=1}^{T}f_{tj}^{2}\right).$

#### Boosting for hierarchical Gaussian priors.

Column boosting is based on interweaving $\theta_{j}$ into the state equation
by choosing $\Psi_{j}=\theta_{j}$. This leads to a prior for the loadings
$\tilde{\beta}_{ij}$ which is independent of $\theta_{j}$, while $\theta_{j}$
acts as variance of the $j$th factors in the expanded model,
$\tilde{f}_{jt}|\theta_{j}\sim\mathcal{N}\left(0,\theta_{j}\right)$. Hence, in
Step (A), under the inverse gamma prior
$\theta_{j}\sim\mathcal{G}^{-1}\left(c^{\theta},b^{\theta}_{j}\right)$,
$\theta_{j}$ is updated in the expanded model,

$\displaystyle\theta_{j}^{\mbox{\rm\tiny
new}}|\theta_{j},{\mathbf{f}}_{r}\sim\mathcal{G}^{-1}\left(c^{\theta}+\frac{T}{2},B_{j}\right),\quad
B_{j}=b^{\theta}_{j}+1/2\sum_{t=1}^{T}\tilde{f}_{jt}^{2}=b^{\theta}_{j}+\theta_{j}/2\sum_{t=1}^{T}f_{tj}^{2}.$

Moving back, all factors and factor loadings in the CFA model are updated:

$\displaystyle\beta_{ij}^{\mbox{\rm\tiny
new}}=\sqrt{\theta_{j}^{\mbox{\rm\tiny
new}}}\tilde{\beta}_{ij}=\frac{\sqrt{\theta_{j}^{\mbox{\rm\tiny
new}}}}{\sqrt{\theta_{j}}}\beta_{ij},\quad f_{jt}^{\mbox{\rm\tiny
new}}=\tilde{f}_{jt}/\sqrt{\theta_{j}^{\mbox{\rm\tiny
new}}}=\frac{\sqrt{\theta_{j}}}{\sqrt{\theta_{j}^{\mbox{\rm\tiny
new}}}}f_{jt}.$

A boosting step that we found to be useful for hierarchical Gaussian priors in
addition to column boosting is interweaving the global shrinkage parameter
$\kappa$ into the prior
$\theta_{j}\sim\mathcal{G}^{-1}\left(c^{\theta},b^{\theta}_{j}\right)$ of the
shrinkage parameter $\theta_{j}$ through the transformation
$\tilde{\theta}_{j}=\kappa\theta_{j}$:

$\displaystyle\beta_{ij}|\delta_{ij}=1,\omega_{ij},\tilde{\theta}_{j},\sigma^{2}_{i}\sim\mathcal{N}\left(0,\tilde{\theta}_{j}\omega_{ij}\sigma^{2}_{i}\right),$
$\displaystyle\tilde{\theta}_{j}|\kappa\sim\mathcal{G}^{-1}\left(c^{\theta},\kappa
b^{\theta}_{j}\right),\quad\kappa\sim\mathcal{G}^{-1}\left(c^{\kappa},b^{\kappa}\right).$

The global shrinkage parameter $\kappa$ is then updated in Step (A) under this
representation, combining the likelihood

$\displaystyle
p(\tilde{\theta}_{1},\ldots,\tilde{\theta}_{r}|\kappa)\propto\kappa^{rc^{\theta}}\exp\left(-\kappa\sum_{j=1}^{r}\frac{b^{\theta}_{j}}{\tilde{\theta}_{j}}\right),$

which is the kernel of a gamma density in $\kappa$ with the inverse gamma
prior $p(\kappa)$. This yields a generalized inverse Gaussian
distribution151515The generalized inverse Gaussian distribution,
$Y\sim\mathcal{GIG}\left(p,a,b\right)$, is a three-parameter family of
probability distribution with support $y\in\mathbb{R}^{+}$. The density is
given by  $\displaystyle\displaystyle
f(y)=\frac{(a/b)^{p/2}}{2K_{p}(\sqrt{ab})}y^{p-1}e^{-(a/2)y}e^{-b/(2y)},$
where $K_{p}(z)$ is the modified Bessel function of the second kind, $a>0$,
$b>0$ and $p$ is a real parameter. to update $\kappa^{\mbox{\rm\tiny new}}$:

$\displaystyle\kappa^{\mbox{\rm\tiny
new}}|\theta_{1},\ldots,\theta_{r},\kappa\sim\mathcal{GIG}\left(rc^{\theta}-c^{\kappa},\frac{2}{\kappa}\sum_{j=1}^{r}\frac{b^{\theta}_{j}}{\theta_{j}},2b^{\kappa}\right).$

Given $\kappa^{\mbox{\rm\tiny new}}$, the column specific parameters are
updated as

$\displaystyle\theta_{j}^{\mbox{\rm\tiny
new}}=\tilde{\theta}_{j}/\kappa^{\mbox{\rm\tiny
new}}=\frac{\kappa}{\kappa^{\mbox{\rm\tiny new}}}\theta_{j},\qquad
j=1,\ldots,r.$

##  Appendix H Initialising Algorithm 1

To check the mixing of the MCMC chain, two (or more) independent runs with
different initial dimensions $r$ of the CFA model and a positive number
$r_{\footnotesize sp}$ of spurious columns are performed. On the one hand, $r$
is chosen to be smaller than the expected number of factors, on the other hand
a large value close to $k-r_{\footnotesize sp}$ is chosen.

Initial values for the leading indices ${\mathbf{l}}_{r}=(l_{1},\ldots,l_{r})$
of the sparsity matrix $\boldsymbol{\delta}_{r}$ are chosen by first sampling
$l_{1}$ from $\\{1,\ldots,u_{1}\\}$, where $u_{1}$ is a small number, e.g. 5.
Then for $j=2,\ldots,r$, we sample $l_{j}$ uniformly from the set ${\cal
L}({\mathbf{l}}_{-j})$ with ${\mathbf{l}}_{-j}=\\{l_{1},\ldots,l_{j-1}\\}$.
Factor loadings below the leading elements are allowed to be zero with
positive probability $p_{0}$, e.g. $p_{0}=0.5$. We draw at most 100 initial
values $\boldsymbol{\delta}_{r}$ (including the leading indices) in this way,
until a matrix $\boldsymbol{\delta}_{r}$ is obtained which satisfies the 3579
counting rule. If no such indicator matrix is found, then we add enough
nonzero elements in each column of $\boldsymbol{\delta}_{r}$ (e.g., by setting
$\delta_{l_{j}+1,j}=1,\ldots,\delta_{l_{j}+3,j}=1$) to ensure variance
identification for the initial value.

To obtain starting values for the factors ${\mathbf{f}}_{r}$, we perform a few
(say 100) MCMC iterations in the confirmatory factor model corresponding to
$\boldsymbol{\delta}_{r}$, which is initialized by sampling the factors
${\mathbf{f}}_{r}=({\mathbf{f}}^{r}_{1},\ldots,{\mathbf{f}}^{r}_{T})$ from the
prior: $f_{jt}\sim\mathcal{N}\left(0,1\right)$, $j=1,\ldots,r$,
$t=1,\ldots,T$. While holding $\boldsymbol{\delta}_{r}$ fixed, we iterate
between sampling the model parameters $\boldsymbol{\beta}_{r}$ and
$\sigma^{2}_{1},\ldots,\sigma^{2}_{m}$ as in Step (P) and sampling the factors
${\mathbf{f}}_{r}$ as in Step (F) of Algorithm 1. The resulting factors
${\mathbf{f}}_{r}$ serve as starting values for the full-blown MCMC scheme
described in Algorithm 1.

##  Appendix I More details for the exchange rate data

#### Data.

The 22 exchange rates analysed in Section 6.1 are listed in Table I.1.

Table I.1: Exchange rate data; currency abbreviations. | 1 | AUD | Australia dollar
---|---|---
2 | CAD | Canada dollar
3 | CHF | Switzerland franc
4 | CZK | Czech R. koruna
5 | DKK | Denmark krone
6 | GBP | UK pound
7 | HKD | Hong Kong dollar
8 | IDR | Indonesia rupiah
9 | JPY | Japan yen
10 | KRW | South Korea won
11 | MXN | Mexican Peso
|  | 12 | MYR | Malaysia ringgit
---|---|---
13 | NOK | Norway krone
14 | NZD | New Zealand dollar
15 | PHP | Philippines peso
16 | PLN | Poland zloty
17 | RON | Romania fourth leu
18 | RUB | Russian ruble
19 | SEK | Sweden krona
20 | SGD | Singapore dollar
21 | THB | Thailand baht
22 | USD | US dollar

#### Running MCMC.

All computations are based on Algorithm 1 and 2. For tuning, we choose
$p_{s}=0.5$ in Step (R) and $p_{\mbox{\rm\footnotesize
shift}}=p_{\mbox{\rm\footnotesize switch}}=1/3$ and $p_{a}=0.5$ in Step (L).
Boosting in Step (A) relies on ASIS with $\sqrt{\Psi_{j}}$ being the largest
loading (in absolute value) in each non-zero column, see Appendix G.

Table I.2: Exchange rate data; inclusion probabilities for the sparsity matrix $\boldsymbol{\delta}_{4}$ for the 2PB prior, averaged over all ordered GLT draws with pivots equal to ${\mathbf{l}}^{\star}=(1,2,5,7)$. Currency | Factor 1 | Factor 2 | Factor 3 | Factor 4
---|---|---|---|---
AUD | 1 | 0 | 0 | 0
CAD | 1 | 1 | 0 | 0
CHF | 0.01 | 0.12 | 0 | 0
CZK | 0 | 0.2 | 0 | 0
DKK | 0.01 | 1 | 1 | 0
GBP | 0.04 | 1 | 0.05 | 0
HKD | 0.01 | 1 | 0.97 | 1
IDR | 0.02 | 1 | 0.03 | 1
JPY | 0.09 | 1 | 0.01 | 0.02
KRW | 0 | 1 | 0.07 | 0.02
MXN | 0 | 0.16 | 0.01 | 0.02
MYR | 1 | 0.06 | 0.01 | 0.01
NOK | 0 | 1 | 0.01 | 0.03
NZD | 0.06 | 0.45 | 0.04 | 0.02
PHP | 0 | 1 | 0.94 | 0.04
PLN | 0 | 1 | 0.02 | 0.75
RON | 0.1 | 0.07 | 0.25 | 0.01
RUB | 0.21 | 0.1 | 0.1 | 0.16
SEK | 0 | 1 | 0.01 | 1
SGD | 0.02 | 1 | 0.03 | 0.99
THB | 0 | 1 | 0.01 | 0.02
USD | 0.01 | 1 | 1 | 0.01

#### Further results.

Table I.2 reports $\mbox{\rm
Pr}(\delta_{ij}=1|{\mathbf{y}},{\mathbf{l}}^{\star})$ for the 2PB prior,
averaged over all ordered GLT draws with pivots equal to
${\mathbf{l}}^{\star}=(1,2,5,7)$. This table is used to derive the MPM. Table
I.3 shows the posterior mean of the factor loading matrix, the idiosyncratic
variances and the communalities, obtained by averaging over all draws where
the pivots of $\boldsymbol{\delta}_{4}$ coincide with ${\mathbf{l}}^{\star}$.
As expected, nonzero factor loadings have relatively high communalities for
the different currencies, whereas for zero rows the communalities are
practically equal to zero percent.

Table I.3: Exchange rate data; posterior mean of the factor loadings
$\Lambda_{ij}$, the communalities $R^{2}_{ij}$ (in percent) and the
idiosyncratic variances $\sigma_{i}^{2}$ (2PB prior) for a 4-factor model with
the GLT constraint ${\mathbf{l}}^{\star}=(1,3,5,7)$.

| Factor loadings | Communalities |
---|---|---|---
Currency | $\Lambda_{i1}$ | $\Lambda_{i2}$ | $\Lambda_{i3}$ | $\Lambda_{i4}$ | $R^{2}_{i1}$ | $R^{2}_{i2}$ | $R^{2}_{i3}$ | $R^{2}_{i4}$ | $\sigma_{i}^{2}$
AUD | 0.962 | 0 | 0 | 0 | 88.4 | 0 | 0 | 0 | 0.12
CAD | 0.391 | 0.601 | 0 | 0 | 16.9 | 38.6 | 0 | 0 | 0.42
CHF | $\approx 0$ | -0.014 | 0 | 0 | 0.01 | 0.3 | 0 | 0 | 0.98
CZK | $\approx 0$ | 0.034 | 0 | 0 | 0.01 | 0.8 | 0 | 0 | 0.98
DKK | $\approx 0$ | 1.06 | 0.22 | 0 | $\approx 0$ | 95.3 | 4.1 | 0 | 0.01
GBP | 0.01 | 0.569 | -0.01 | 0 | 0.26 | 31.4 | 0.27 | 0 | 0.71
HKD | $\approx 0$ | 0.502 | 0.398 | 0.762 | $\approx 0$ | 21.8 | 14.2 | 49.4 | 0.17
IDR | $\approx 0$ | 0.8 | $\approx 0$ | 0.427 | 0.05 | 57.2 | 0.06 | 16.9 | 0.28
JPY | 0.013 | 0.924 | $\approx 0$ | $\approx 0$ | 0.23 | 75.6 | $\approx 0$ | 0.01 | 0.27
KRW | $\approx 0$ | 1.06 | $\approx 0$ | $\approx 0$ | $\approx 0$ | 95.9 | 0.04 | $\approx 0$ | 0.05
MXN | $\approx 0$ | 0.025 | $\approx 0$ | $\approx 0$ | $\approx 0$ | 0.58 | 0.04 | 0.04 | 0.98
MYR | 0.793 | $\approx 0$ | $\approx 0$ | $\approx 0$ | 61.2 | 0.04 | 0.01 | 0.01 | 0.40
NOK | $\approx 0$ | 0.88 | $\approx 0$ | $\approx 0$ | $\approx 0$ | 69.7 | 0.01 | 0.05 | 0.33
NZD | 0.016 | 0.108 | -0.01 | $\approx 0$ | 0.50 | 3.17 | 0.29 | 0.04 | 0.95
PHP | $\approx 0$ | 0.549 | -0.415 | $\approx 0$ | 0.01 | 28.7 | 17.9 | 0.17 | 0.56
PLN | $\approx 0$ | 0.995 | $\approx 0$ | 0.128 | $\approx 0$ | 85.6 | 0.01 | 2.02 | 0.14
RON | 0.027 | $\approx 0$ | -0.08 | $\approx 0$ | 0.87 | 0.09 | 3.08 | 0.02 | 0.95
RUB | -0.067 | 0.012 | 0.029 | 0.043 | 2.37 | 0.28 | 0.95 | 1.37 | 0.94
SEK | $\approx 0$ | 0.974 | $\approx 0$ | 0.31 | $\approx 0$ | 81.8 | 0.01 | 8.7 | 0.11
SGD | $\approx 0$ | 0.748 | $\approx 0$ | 0.392 | 0.05 | 50.9 | 0.05 | 14.6 | 0.38
THB | $\approx 0$ | 0.583 | $\approx 0$ | $\approx 0$ | 0.01 | 32.7 | 0.01 | 0.05 | 0.70
USD | $\approx 0$ | 1.06 | 0.215 | $\approx 0$ | $\approx 0$ | 95.3 | 4.1 | $\approx 0$ | 0.01

Note: entries with $|\Lambda_{ij}|<0.01$ and entries with $R^{2}_{ij}<0.1$ are
indicated by $\approx 0$.