O90 O1em m #3

# Creating area level indices of behaviours impacting cancer in Australia with
a Bayesian generalised shared component model

James Hogg1 Susanna Cramb3 Jessica Cameron1,2 Peter Baade1,2 Kerrie Mengersen1
Corresponding author<EMAIL_ADDRESS>
(1Centre for Data Science, Queensland University of Technology (QUT), 2 George
St, Brisbane City, 4000, Queensland, Australia
2Viertel Cancer Research Centre, Cancer Council Queensland (CCQ), 553 Gregory
Terrace, Fortitude Valley, 4006, Queensland, Australia
3Australian Centre for Health Services Innovation, School of Public Health and
Social Work, Queensland University of Technology (QUT), 2 George St, Brisbane
City, 4000, Queensland, Australia
)

###### Abstract

This study develops a model-based index creation approach called the
Generalised Shared Component Model (GSCM) by drawing on the large field of
factor models. The proposed fully Bayesian approach accommodates
heteroscedastic model error, multiple shared factors and flexible spatial
priors. Moreover, unlike previous index approaches, our model provides indices
with uncertainty. Focusing on unhealthy behaviors that lead to cancer, the
proposed GSCM is used to develop the Area Indices of Behaviors Impacting
Cancer product — representing the first area level cancer risk factor index in
Australia. This advancement aids in identifying communities with elevated
cancer risk, facilitating targeted health interventions.

_Keywords_ Bayesian spatial statistics, factor models, unhealthy behaviors,
indices, shared component models

##  1 Introduction

In the past two decades, the resolution and accessibility of health data for
small geographic areas have improved considerably. While these improvements
allow for more granular insights into spatial disparities in specific diseases
or health characteristics, the increasing quantity and complexity of such data
can make it more challenging to collate, draw meaningful conclusions and
formulate effective policy decisions [1]. Recently, interactive and publicly
available health atlases (e.g., Australian Cancer Atlas [2] or PLACES platform
[3]) have gained popularity as valuable tools for disseminating health data
for small areas [2, 4, 3]. These atlases often provide estimates and measures
of uncertainty (e.g., measurement error) for multiple health measures, opening
up avenues for further spatial analysis [5, 6].

Multidimensional data often contain redundancy; that is, some dimensions may
not provide substantive additional information. In this case, the data could
be reasonably summarised by a lower-dimensional latent space, called an index
[1]. An index should measure a multifaceted concept that cannot be captured by
a single feature. Indices often offer a simple interpretation of a complex
phenomenon, can help determine policy priorities, and can be used to benchmark
or monitor the performance of current policies [7].

There are a plethora of methods routinely used to derive indices [8, 7]
including rank sum, standardisation [7] and principal component analysis (PCA)
[9]. However, these methods do not naturally accommodate spatial
autocorrelation and area-specific uncertainty, which are common in atlas data.
Since atlas data are primarily derived from models, their uncertainty is
largely model-based. Like other work [10], we use the generic term measurement
error to denote this model-based uncertainty.

While spatial PCA can accommodate spatial autocorrelation [11, 12], model-
based approaches are gaining popularity in the creation of indices [1, 13, 14,
15] as these approaches can also accommodate measurement error and provide
uncertainty measures. Furthermore, previous work highlights that by fitting
Bayesian models [13, 14, 16], researchers can provide probabilistic statements
about index values. This facilitates a sophisticated decision-making process
that uses the distributions of index values rather than point estimates. With
access to distributions, researchers can more effectively determine whether
variations in index values are truly meaningful.

Shared component models (SCMs) [17] have been used to create indices from
spatial data [14]. These models assume that the input data originate from a
linear combination of a _single_ shared factor and residual errors [18, 19].
However, conventional SCMs do not naturally accommodate measurement error [20,
21, 22, 23, 24] and can be limiting when the input data lack strong
correlations [20].

This study introduces a Bayesian generalised SCM (GSCM) capable of
accommodating heteroscedastic measurement error, multiple shared factors and
flexible spatial priors [25]. Innovations are introduced with respect to the
modelling, computation, and application. In terms of the model, while Bayesian
spatial factor models have been employed previously to generate health indices
with uncertainty measures [13, 15], the GSCM bridges the mathematical and
intuitive gaps between SCMs and Bayesian spatial factor models (BSFM); a
connection that, to our knowledge, has yet to be acknowledged. To help
distinguish between SCMs and BSFMs, we define a SCM as a model with a single
shared factor and spatial dependencies in both the shared and residual terms.
A BSFM, on the other hand, has at least one shared factor and spatial
dependency in at least one component of the model.

On the computational front, the model was fitted using Stan [26], a flexible
tool that facilitates efficient analysis of diverse Bayesian models. As a side
effect of this work, we also present an efficient implementation of the Leroux
prior [27] in Stan, which, to our knowledge, represents the first feasible
implementation of this increasingly common spatial prior using Hamiltonian
Monte Carlo.

The motivation for this work is a significant public health challenge; cancer.
In light of the substantial global health burden [28], cancer prevention
strategies focus on reducing the prevalence of unhealthy behaviors (e.g., risk
factors) that can lead to cancer such as drinking, smoking, obesity, poor
diet, and inadequate activity [29]. Previous research has found spatial
disparities in unhealthy behaviors that can lead to cancer [4, 3, 30, 31],
such as smoking [32]. Yet, while area-level health indices are commonplace
[33], cancer-specific indices have primarily been developed at the individual
level [34, 35]. Area-level cancer risk indices have remained unexplored.

To address this gap and showcase the potential of the GSCM approach, we
conduct a comprehensive case study using publicly available modelled data on
unhealthy behaviors in Australia. The aim of this study is threefold:

  * (1)

generate a summary measure of unhealthy behaviors in Australia to help
understand any underlying spatial patterns;

  * (2)

identify areas with elevated prevalence of unhealthy behaviors that may cause
cancer;

  * (3)

provide decision support to help policymakers focus their cancer prevention
strategies.

To do this, we developed the first area-level cancer risk index in Australia
called the Area Indices of Behaviors Impacting Cancer (AIBIC), which provides
a relative summary measure to compare the prevalence of unhealthy behaviors
between areas. Unlike other applications of model-based index creation, which
provide a single index only [13, 15, 14], we provide four indices in the AIBIC
and explore their distinct properties. The first two indices stem directly
from the fitted GSCM, while the latter two are unique transformations of the
former pair.

##  2 Methods

###  2.1 Data

The presented analyses utilize data from our previous Bayesian spatial small
area estimation study [31] which used data from the 2017-2018 National Health
Survey (NHS). The NHS, conducted by the Australian Bureau of Statistics (ABS),
is an Australia-wide population-level health survey conducted every 3-4 years
[36, 37]. In our previous study [31] we generated modelled estimates and
measures of uncertainty for the proportion of persons aged 15 years and older
in 2017-18 engaging in behaviors that increase the risk of being diagnosed
with cancer [29]. Modelled results were provided for 2221 small areas, known
as Statistical Area Level 2 or “SA2s”, as defined by the 2016 Australian
Statistical Geography Standard (ASGS) [38].

For the AIBIC product, we used five features from Hogg _et al._ [31], namely,
the proportion of the following unhealthy behaviors: current smokers, risky
alcohol consumption, overweight/obese, inadequate diet and inadequate physical
activity based on leisure only. The selection of these unhealthy behaviors
followed the same rationale as Hogg _et al._ [31], involving consultations
with experts, a review of relevant literature [39, 40, 41, 29] and the
evaluation of data availability in the 2017-2018 NHS.

These data on unhealthy behaviors were available as point estimates and
standard deviations of proportions. For all the unhealthy behaviors, each was
defined so that a higher proportion indicated a higher prevalence of poorer
health behavior [31]. To enable the use of a Gaussian likelihood in our model
we first applied the inverse logit transformation, resulting in
$\mathbf{Y}_{k}$ and $\mathbf{S}_{k}$ being $N$-dimensional vectors of
estimates and corresponding standard deviations for risk factor $k=1,\dots,5$,
where $N=2221$ is the total number of small areas (i.e., SA2s). Finally, these
logit-transformed point estimates were centered and scaled, with the
corresponding standard deviations appropriately adjusted (see Section A of the
Additional File). The final set of input features is available on Github [42].

Alongside the data on unhealthy behaviors, in this study, we utilised
population data for persons 15 years and older at the SA2 level. The vector of
population counts, $\mathbf{P}=\left(P_{1},\dots,P_{N}\right)$, was sourced
from the ABS’s Estimated Resident Population product and calculated as the
average annual count for 2017 and 2018 [43].

####  2.1.1 Risk factors

Current smoking was defined as those who reported to be current smokers
(including daily, weekly or less than weekly), and had smoked at least 100
cigarettes in their lifetime.

Risky alcohol consumption was defined as individuals who exceeded the revised
2020 National Health and Medical Research Council (NHMRC) guidelines [44] of
up to 10 standard drinks/week and no more than 4 standard drinks in any day.
Compliance with the guidelines was assessed using self-reported alcohol
consumption during the last three drinking days from the preceding seven days.

Inadequate diet, which was based on self-reported diet, was defined as those
who did not meet both the fruit (2 serves/day) and vegetable (5 serves/day)
2013 NHMRC Australian Dietary guidelines [45]. Participants of the 2017-18 NHS
were asked to report the number of servings of fruit and serves of vegetables
they usually ate each day.

Overweight/obese was defined as persons with a Body Mass Index greater or
equal to 25. The NHS data collects both measured and self-reported height and
weight. Measurements were voluntary, so 40% of the measured values in the NHS
were imputed by the ABS using the hot-decking method [46].

Inadequate physical activity based on leisure only was defined as those who
did not meet the 2014 Department of Health Physical Activity guidelines [47].
This guideline stipulated that each week, adults aged 18-64 should either do 2
1/2 to 5 hours of moderate-intensity physical activity or 1 1/4 to 2 1/2 hours
of vigorous-intensity physical activity or an equivalent combination of both,
plus muscle-strengthening activities at least 2 days each week.

###  2.2 Statistical Model

Shared component models (SCM) are a multivariate statistical tool developed in
the context of disease mapping. The traditional assumption underpinning the
SCM is that the input features share a single common factor [20, 21] or latent
field that “manifests” in the separate features. To accommodate residual
heterogeneity not captured by the shared factor, SCMs also include feature-
specific residual errors [17]. In disease mapping, these models assume spatial
dependence in both the shared factor and the feature-specific residual errors.
Borrowing from the massive literature on factor models, incorporating multiple
factors into the SCM enhances its ability to explain variation in the
multivariate data, improves the approximation of the covariance structure, and
facilitates more informative inferences.

Given the expressed interest in accommodating the uncertainty of the input
features, there is a need to include area-specific measurement error.
Following the classical measurement error framework [48], which is abundant in
spatial models [49, 50], we adopt a similar method to that of Jahan _et al._
[5]. With the constraint that the number of shared factors ($L$) was much
fewer than the number of features (i.e., $L<<K$), the proposed GSCM is written
as a $N$-dimensional multivariate normal distribution for the $k$th feature,

$\displaystyle\mathbf{Y}_{k}$ $\displaystyle\sim$
$\displaystyle\text{MVN}_{N}\left(\boldsymbol{\mu}_{k},\text{diag}\left(\boldsymbol{\sigma}_{k}\right)\right)$
(1) $\displaystyle\boldsymbol{\mu}_{k}$ $\displaystyle=$
$\displaystyle\mathbf{z}\left(\boldsymbol{\Lambda}_{k}\right)^{T}+\boldsymbol{\epsilon}_{k},$

where $\mathbf{Y}_{k}$, $\boldsymbol{\sigma}_{k}$ and $\boldsymbol{\mu}_{k}$
are $N$-dimensional vectors of estimates, corresponding standard deviations
and true values for the $k$th risk factor. Unlike Jahan _et al._ [5], who
placed a prior on their standard deviation estimates, in this work
$\boldsymbol{\sigma}_{k}$ is taken to be equal to $\mathbf{S}_{k}$ (i.e., the
standard deviations from Hogg _et al._ [31]).

In line with the classical measurement error framework [48], while the
estimates are assumed independent, the true values are assumed to exhibit a
complex correlation structure, captured in the $\boldsymbol{\mu}_{k}$’s.
$\mathbf{z}$ is a $N\times L$ matrix of shared factor scores, with independent
columns and $\boldsymbol{\epsilon}_{k}$ a $N$-dimensional vector of feature-
specific residual errors for the $k$th feature with associated variance
$\tau_{k}^{2}$. $\boldsymbol{\Lambda}_{k}$ is the $k$th row of the $K$ by $L$
matrix of factor loadings. We assume that the $\boldsymbol{\epsilon}_{k}$’s
are independent of each other and also from $\mathbf{z}$. The shared factors,
$\mathbf{z}$, are assumed homoscedastic, where the columns have fixed
variances equal to 1 [1, 51]. In this work, these shared factors are used to
simplify our understanding of the spatial patterns of unhealthy behaviors and
represent the building blocks for the AIBIC product described later. The
covariance structure of the GSCM is described in Section C of the Additional
File.

Since factor models are notorious for being unidentifiable, there is an
abundance of research in this area [52, 53, 54]. In this work, we use the
hierarchical structural constraint [55] where the factor loading matrix,
$\boldsymbol{\Lambda}$, is constrained to be lower triangular with strictly
positive diagonal entries $\lambda_{ii}>0$ [1]. This constraint enforces
bounds on the number of free parameters in $\boldsymbol{\Lambda}$ and the
maximum number of shared factors. In this work, where $K=5$, $L$ was set to
$2$, requiring the estimation of nine parameters in $\boldsymbol{\Lambda}$.

####  2.2.1 Spatial dependence

Although the $N$-dimensional vectors $\boldsymbol{\epsilon}_{k}$ and
$\mathbf{z}$ are independent of each other, we assume spatial dependence
within each vector. Whilst previous research has used spatial priors on the
loading matrix [56, 57], shared factors [1, 13], or feature-specific residual
errors only [58], traditional SCMs assume spatial dependence in the shared
factors and feature-specific residual errors [14, 17]. Although marginal
spatial priors have been explored, particularly in geostatistical applications
[13, 16, 51], conditional autoregressive priors [59] are more common in SCMs
[20, 21, 15] and disease mapping applications in general [25].

In this work, we fit independent $N$-dimensional conditional spatial priors on
the columns of $\mathbf{z}$ and the $\boldsymbol{\epsilon}_{k}$’s. To allow
for a blend of spatially structured and unstructured variation, we used the
Leroux prior [27], hereafter called the LCAR (see Section B of the Additional
File).

The LCAR represents a generalisation of two priors due to the inclusion of a
spatial autocorrelation parameter ranging from 0 to 1. When the spatial
autocorrelation parameter is equal to 1, the LCAR becomes the common intrinsic
conditional autoregressive (ICAR) prior [22], which has been used in SCMs and
BSFMs before [14, 56]. When the spatial autocorrelation parameter is equal to
0, the LCAR collapses to an unstructured independent and identically
distributed (IID) normal prior. The spatial autocorrelation parameters from
the shared factors and feature-specific residual errors are denoted as
$\rho_{l}$ and $\kappa_{k}$, respectively.

The GSCM represents several generic models as special cases. The conventional
(non-measurement error) SCM is recovered when uncertainty of the input data is
ignored (e.g., setting $\sigma_{nk}=0,\forall n,k$) and $L=1$. Continuing to
ignore measurement error, a BSFM is recovered similar to that of Nethery _et
al._ [1] when all $\kappa_{k}=0$ or that of Mezzetti _et al._ [58] when all
$\rho_{k}=0$. Furthermore, when all $\rho_{k}=0,\kappa_{k}=0$ the GSCM becomes
a generic Bayesian factor model.

####  2.2.2 Priors

To complete the Bayesian GSCM, prior distributions are chosen for the
parameters
$\boldsymbol{\Lambda},\boldsymbol{\tau},\boldsymbol{\rho},\boldsymbol{\kappa}$.
We use weakly informative $\text{N}\left(0,1\right)$ priors on the free
parameters of $\boldsymbol{\Lambda}$. Mildly informative
$\text{Gamma}\left(2,3\right)$ priors, which have a density of 0.04 below 0.1,
were used for $\boldsymbol{\tau}=\left(\tau_{1},\dots,\tau_{5}\right)$. Mildly
informative $\text{Beta}\left(6,2\right)$ priors, which have a density of 0.94
above 0.5, were used for
$\boldsymbol{\rho}=\left(\rho_{1},\rho_{2}\right),\boldsymbol{\kappa}=\left(\kappa_{1},\dots,\kappa_{5}\right)$.
These mildly informative priors were used to improve the convergence and
identifiability of the model and avoid degenerate feature-specific residual
errors. A sensitivity analysis found that these priors had minimal effect on
the inference of $\mathbf{z}$.

####  2.2.3 Computation

We used fully Bayesian inference with MCMC via the R package rstan Version
2.26.11 [26]. Where possible we used the non-mean centered parameterisation
for the shared factors and feature-specific residual errors. Although
efficient implementations of the proper CAR [60] and ICAR [61] are available,
to the best of our knowledge this is the first implementation of the LCAR in
Stan. This computational novelty is widely applicable to any spatial model
fitted using Stan, not just the model proposed in the current paper. Code and
details can be found on Github [42] or in Section B.1 of the Additional File.

We used 4000 warmup and 6000 post-warmup draws for each of the four chains.
For storage reasons we thinned the final posterior draws by three, resulting
in 8000 draws for inference. Convergence of the models was assessed using
trace and autocorrelation plots, effective sample size, and $\hat{R}$ [62].

####  2.2.4 Model fit metrics

To formally compare models we used the Deviance Information Criterion (DIC),
Widely Applicable Information Criterion (WAIC) [63], and mean absolute error
(MAB) between the estimates, $\mathbf{Y}_{k}$, and the posterior median of the
true values, $\boldsymbol{\mu}_{k}$. For these metrics, smaller values were
preferred.

####  2.2.5 Summaries

To address the first aim of this study, which was to generate a summary
measure of unhealthy behaviors in Australia, the results from the GSCM are
used to acquire draws for four $N$-dimensional vectors. Each vector represents
the underlying scores for the four indices in the AIBIC product described
later.

The first two vectors are the columns of the shared factor scores,
$\mathbf{z}$. The third is a weighted sum of these [7], where the weights,
$w_{1},w_{2}$, are based on the squared factor loadings (see the Additional
File for details), via

$w_{l}=\left(\sum_{k=1}^{5}\lambda_{kl}^{2}\right)/\left(\sum_{k=1}^{5}\sum_{l=1}^{2}\lambda_{kl}^{2}\right),l=1,2.$

The fourth version is obtained by weighting the third version by the SA2
population.

After acquiring these four vectors, the posterior draws were converted into
posterior ranks and posterior percentiles [15, 9, 7]. Posterior ranks were
derived by ranking each draw, while posterior percentiles were derived by
assigning each draw a percentile value within the range of 1 to 100. Although
the posterior percentiles and ranks provide similar inference, the ranks
provide slightly more granular information, whilst percentiles are slightly
easier to interpret.

The final step to produce the AIBIC product was to appropriately summarise the
draws of these vectors. We use the posterior median and 95% highest posterior
density interval (HPDI). To address the second and third aims of this study,
which were to identify areas with elevated cancer risk and provide decision
support to easily select appropriate areas, we calculated the posterior
probabilities of the vectors exceeding the 80th, 95th, and 99th percentiles
[13, 14]. Additionally, we derived the posterior probability of the vectors
being ranked within the top 10, 20, and 100 areas.

The point estimates, HDPIs, and posterior probabilities for the vectors,
percentiles, and ranks are available on Github [42]. These summaries are
provided for each of the four indices, with percentiles and ranks computed
both nationally and for each of the eight states and territories of Australia.

##  3 Results

###  3.1 Exploratory Analyses

This section explores the unhealthy behaviors in terms of spatial dependency
and correlation structure (Figures 1 and 2). This exploratory analysis is used
to inform the fitting of the GSCM in the following sections. Section A in the
Additional File provides descriptive statistics for the input features.

The correlation coefficients of each pair of features are illustrated in
Figure 1. It indicates a strong positive relationship between current smoking,
inadequate physical activity and overweight/obese, and a negative correlation
between risky alcohol consumption and inadequate physical activity. Exploring
the standard errors of the input data reveals that risky alcohol consumption
has the lowest error (median of 0.32), while inadequate diet (median of 0.56)
and overweight (median of 0.51) have the highest. Univariate Moran’s $I$ tests
[64] indicate that all the risk factors exhibit a high level of spatial
autocorrelation.

Figure 1: Visualisation of the correlation matrix for the point estimates of
the five unhealthy behaviors. Figure 2: Visualisation of the loadings from a
Principal Components Analysis (PCA) applied to the point estimates for the
five unhealthy behaviors. The horizontal bars are colored according to the
direction of the loading with red indicating a negative loading. Above the
plots, the percentage of variation explained by Principal Component 1 (PC1)
and PC2 are given.

Initial exploratory analysis also involved a Principal Components Analysis
(PCA). The first two components accounted for 72% of the total variation
(Figure 2). Notably, the first principal component captured less than half of
the total variation (42%) suggesting a single shared factor would be
insufficient. Furthermore, risky alcohol consumption exhibited negative or
weak correlations with the majority of the other unhealthy behaviors and had a
very low loading with the first principal component. Thus, using a second
factor allowed us to effectively capture the distinct patterns associated with
risky alcohol consumption, ensuring that these patterns were not diluted
within a single shared factor or lost within the residuals.

The results from the PCA were used to inform the ordering of the features in
$\boldsymbol{\Lambda}$ in the GSCM [56]. This ordering does not impact the
theoretical model or its predictions [65, 55], but can influence the shared
factors and consequently, the interpretation of the indices. Given that a
maximum of two shared factors is possible, only the ordering for the first two
features matter in this case. We used risky alcohol consumption followed by
current smoking; a choice driven by the negligible loading on risky alcohol
consumption with principal component 1 (PC1) and strong positive loading with
PC2. Moreover, current smoking was chosen as the second feature due to its
strong positive loading with both PC1 and PC2. This ordering aimed to
decorrelate the two shared factors.

###  3.2 GSCM model results

To identify the best specification, we initially fit a range of GSCMs with
increasing complexity, as described in the following paragraphs. This process
allowed us to assess the impact of measurement error, the number of shared
factors and spatial structures, with a small section dedicated to each.
Following this assessment, we select a specific GSCM specification, which
serves as the basis for developing the AIBIC product.

Note that to improve identifiability and convergence, we used a non-spatial
(e.g., IID) prior for the first feature-specific residual error in all models.

####  3.2.1 Impact of measurement error

To assess the impact of ignoring measurement error, we compare the results
from two-factor GSCMs where one model accommodates the measurement error and
the other does not. Both models use LCARs for the shared factors and IID
priors for the feature-specific residual errors. Figure 3, which compares the
posterior median and 95% HDPIs, confirms the expected result; the GSCM with
measurement error provides shared factor 1 scores with distinct patterns (plot
(a)) and considerably increased uncertainty (plot (b)). The outliers on the
right of plot (a) are mostly lower socioeconomic areas within major cities.
Figure S1 in Section E of the Additional File provides a similar plot for the
shared factor 2 scores, illustrating slightly more agreement between the two
models, but a similar increase in uncertainty.

The factor loadings are considerably impacted when the measurement error is
accommodated (Table 1). The large relative changes, particularly for
overweight/obese and inadequate diet, are an artefact of the higher standard
errors of these input features. The largest change is for overweight/obese
which loads with 0.40 for factor 1 when measurement error is ignored, but 0.07
when measurement error is accommodated. This represents both an over 5-fold
decrease in central tendency and an over 4-fold increase in relative standard
error.

Note the 46% decrease in the loading of risky alcohol consumption with factor
1 (Table 1). This is noteworthy, particularly since this feature loads the
strongest with this factor and has the smallest standard errors among the
features, potentially explaining the marked differences in factor loadings
when the measurement error is accommodated.

Considering the significant impact of measurement error, all subsequent models
will accommodate the known measurement error.

Figure 3: Comparison of the first shared factor when estimated using a two-factor GSCM with and without accommodating the measurement error (ME). Scatter plot (a) gives the posterior median and 95% HPDIs from the GSCM with ($y$-axis) and without ($x$-axis) ME. The red diagonal line represents equivalence of the $x$ and $y$ axes. Caterpillar plot (b) compares the posterior standard deviations ranked by the “no ME” model. Table 1: Posterior medians (95% highest posterior density confidence intervals (HPDI)) of the factor loadings from two-factor GSCMs where one accommodates the measurement error (ME), while the other does not. Both models use LCARs for the shared factors and non-spatial priors for the feature-specific residual errors.  | No ME | With ME
---|---|---
| Factor 1 | Factor 2 | Factor 1 | Factor 2
Risky alcohol consumption | 1.17 (1.07, 1.27) | 0 | 0.76 (0.73, 0.80) | 0
Current smoking | 0.04 (-0.04, 0.14) | 0.62 (0.57, 0.66) | -0.12 (-0.16, -0.08) | 0.54 (0.51, 0.57)
Inadequate physical activity | -0.77 (-0.84, -0.69) | 0.71 (0.67, 0.75) | -0.48 (-0.52, -0.43) | 0.54 (0.51, 0.58)
Inadequate diet | -0.30 (-0.37, -0.23) | 0.28 (0.24, 0.32) | -0.09 (-0.13, -0.06) | 0.30 (0.27, 0.32)
Overweight /obese | 0.40 (0.34, 0.46) | 0.68 (0.64, 0.72) | 0.07 (0.02, 0.12) | 0.66 (0.62, 0.69)

####  3.2.2 Impact of multiple shared factors

To assess the impact of using $L>1$ in the GSCM, we compared the first shared
factor when fitting a one-factor and two-factor GSCM with LCARs for both the
shared factors and feature-specific residual errors. Figure 4, which compares
the posterior median and 95% HDPIs for the scores, confirms our expectation:
the inclusion of the second factor in the GSCM has a negligible effect on both
the central tendency and uncertainty of the factor 1 scores.

More broadly though, the two-factor model offers a more accurate
representation of the data’s covariance structure, leading to the WAIC
favoring it as the superior model. The WAIC for the one-factor and two-factor
models are 16744 and 15591, respectively.

Furthermore, Table 2 illustrates the large loadings of all the features with
the second factor, supporting its utility in this work. The factor 1 loadings
are only minimally affected when factor 2 is included, with the greatest
change occurring for current smoking and inadequate diet. The difficulty in
interpreting factor 1, given the negative loadings, is discussed in Sections
3.3 and 4.

Table 2: Posterior medians (95% highest posterior density confidence intervals (HPDI)) of the factor loadings from a one-factor and two-factor GSCM with LCARs for both the shared factors and feature-specific residual errors. | One-factor GSCM | Two-factor GSCM
---|---|---
| Factor 1 | Factor 1 | Factor 2
Risky alcohol consumption | 0.78 (0.74, 0.81) | 0.77 (0.73, 0.81) | 0
Current smoking | -0.20 (-0.25, -0.16) | -0.14 (-0.19, -0.10) | 0.55 (0.51, 0.58)
Inadequate physical activity | -0.46 (-0.51, -0.41) | -0.47 (-0.52, -0.42) | 0.51 (0.47, 0.56)
Inadequate diet | -0.29 (-0.34, -0.24) | -0.23 (-0.28, -0.18) | 0.40 (0.36, 0.45)
Overweight /obese | 0.04 (-0.01, 0.10) | 0.08 (0.03, 0.13) | 0.63 (0.58, 0.67)
Figure 4: Comparison of the first shared factor when estimated using the
1-factor GSCM and the 2-factor GSCM where both have LCARs for the shared
factors and feature-specific residual errors. Scatter plot (a) gives the
posterior median and 95% HPDIs from the 1-factor ($y$-axis) and the 2-factor
($x$-axis) models. The red diagonal line represents the equivalence of the $x$
and $y$ axes. Caterpillar plot (b) compares the posterior standard deviations
ranked by the 2-factor model.

####  3.2.3 Impact of spatial priors

We explored which combination of spatial priors provided the best fitting two-
factor GSCM to these data. For the shared factors and feature-specific
residual errors, three priors were possible (IID, LCAR, and ICAR) resulting in
nine possible models. The results (Table 3) show that the models with spatial
priors for both the shared factors _and_ feature-specific residual errors
provided a superior fit to the multivariate data.

Table 3: Comparison of the two-factor GSCM when altering the prior structure for the shared factors and feature-specific residual errors. Model comparison is carried out using the Deviance Information Criterion (DIC), Widely Applicable Information Criterion (WAIC), and mean absolute error (MAB). The bolded quantities represent the lowest value in each column where appropriate. Shared | Specific | DIC | WAIC | MAB
---|---|---|---|---
LCAR | ICAR | 16457 | 15583 | 0.329
ICAR | ICAR | 16443 | 15589 | 0.330
LCAR | LCAR | 16500 | 15591 | 0.327
ICAR | LCAR | 16489 | 15601 | 0.328
IID | ICAR | 17740 | 16427 | 0.295
IID | LCAR | 17827 | 16476 | 0.293
LCAR | IID | 18209 | 17987 | 0.388
ICAR | IID | 18201 | 17990 | 0.389
IID | IID | 19634 | 19299 | 0.363

###  3.3 Area Indices of Behaviors Impacting Cancer (AIBIC)

By considering the model selection criteria together (Table 3), the two-factor
GSCM with LCAR priors for the shared factors and ICAR priors for the feature-
specific residual errors was considered the optimal model. The factor loadings
are given in Table 4, and summaries of the other model parameters are provided
in Section D of the Additional File.

Table 4: Posterior medians (and 95% highest posterior density confidence intervals (HPDI)) of the factor loadings from the selected GSCM. | Factor 1 | Factor 2
---|---|---
Risky alcohol consumption | 0.77 (0.73, 0.81) | 0
Current smoking | -0.15 (-0.19, -0.10) | 0.54 (0.50, 0.58)
Inadequate physical activity | -0.47 (-0.52, -0.42) | 0.51 (0.47, 0.55)
Inadequate diet | -0.24 (-0.29, -0.19) | 0.40 (0.35, 0.45)
Overweight /obese | 0.09 (0.03, 0.14) | 0.64 (0.60, 0.69)

The shared factors from this model are used to construct the AIBIC product,
which provides a relative summary measure to compare the prevalence of
unhealthy behaviors between SA2s. A low index value for an area is likely to
be an area where a higher proportion of the residents engage in unhealthy
behaviors that can lead to cancer.

The AIBIC product helps to address the three aims of this study: generating a
summary measure of unhealthy behaviors, identifying areas with a higher
prevalence of unhealthy behaviors, and providing decision support for cancer
prevention in Australia. The product is tailored for broad cancer prevention
strategies that aim to simultaneously address multiple unhealthy behaviors
that can lead to cancer.

In this section, we propose and compare the four AIBIC indices, each with a
unique purpose and interpretation (Table 5). While the first two are derived
directly from the selected GSCM, the second two are unique summations of the
first two. The AIBIC product can be found on Github [42].

Table 5: Terminology for the four indices of the AIBIC product. In the table, $z_{n1}$ is the first shared factor raw score for the $n$th area and $P_{n}$ is the corresponding population count. | Name | Source
---|---|---
1 | Factor 1 | $z_{n1}$
2 | Factor 2 | $z_{n2}$
3 |  Health Behaviors Index (HBI) | $w_{1}z_{n1}+w_{2}z_{n2}$
4 |  Population Adjusted Health Behaviors Index (PAHBI) | $P_{n}\left(w_{1}z_{n1}+w_{2}z_{n2}\right)$

####  3.3.1 Factor 1

This first shared factor ($\mathbf{z}_{1}$), mapped in Figure S4 of the
Additional File, could be considered an Alcohol, Activity Index as the
loadings with the highest absolute magnitude appear on risky alcohol
consumption and inadequate physical activity (Table 4). The direction of the
loadings suggests that areas with a high factor 1 score are areas with a high
prevalence of risky alcohol consumption, but a low prevalence of inadequate
physical activity. We discuss the implications of the negative loadings later
(Section 4).

####  3.3.2 Factor 2

This second shared factor ($\mathbf{z}_{2}$), mapped in Figure S5 in the
Additional File, could be considered a Smoking, Activity and Weight Index,
given its strong loadings with these unhealthy behaviors. Thus, factor 2
represents an approximate average of the unhealthy behaviors, excluding risky
alcohol consumption.

####  3.3.3 Health Behaviors Index (HBI)

Although factor 1 or 2 could feasibly be used to address the first and second
aims of this study, justifying the sole use of either is challenging. While
factor 1 can be interpreted as risky alcohol consumption minus the average of
the other features, factor 2 averages all the features, except risky alcohol
consumption.

Rather than select one, we followed previous advice in the Handbook of
Constructing Composite Indicators [7], and derived a weighted average of
factors 1 and 2 instead. This weighted average is the Health Behaviors Index
(HBI) (Table 5). The Pearson correlations are 0.66 between factor 1 and the
HBI, and 0.90 between factor 2 and the HBI (Figure 5); reflecting the uneven
weights of 0.45 and 0.55 (Section 2.2.5).

The HBI provides a single value for each area that is more representative of
the five unhealthy behaviors, both in terms of point estimates and
uncertainty. Figure 6 illustrates the spatial disparities in the HBI with
lower values occurring in the major cities, indicating a generally lower
prevalence of unhealthy behaviors in these areas. Moreover, within each
remoteness category, HBI values are generally lower for less socioeconomically
disadvantaged areas (Figure S6 in the Additional File). These results
generally agree with the spatial patterns of the underlying unhealthy
behaviors [31].

Figure 5: Pairs plot of the four indices in the AIBIC product. Posterior
medians of the percentiles are used. The plots in the lower triangular squares
are scatter plots, the diagonal elements are density representations and the
empirical correlations between the four vectors are given in the upper
triangular. Figure 6: Choropleth maps displaying the HBI, which is based on
the summation of the two shared factor scores from the two-factor GSCM. The
results are the posterior median of the posterior percentiles for 2221 SA2s
across Australia. Higher values indicate areas with a higher prevalence of
unhealthy behaviors. The map includes insets for the eight capital cities for
each state and territory, with white boxes on the main map indicating the
location of the inset. White lines represent the boundaries of the eight
states and territories of Australia. Figure 7: Choropleth maps displaying the
PAHBI, which is based on the summation of the two shared factor scores from
the two-factor GSCM. The results are the posterior median of the posterior
percentiles for 2221 SA2s across Australia. Higher percentiles indicate areas
with a higher prevalence of unhealthy behaviors. The map includes insets for
the eight capital cities for each state and territory, with white boxes on the
main map indicating the location of the inset. White lines represent the
boundaries of the eight states and territories of Australia.

####  3.3.4 Population Adjusted Health Behaviors Index (PAHBI)

Considering the third aim of this study (e.g., decision support), areas with
high HBI values do not necessarily correspond with areas with large
populations, where the impact may be more pronounced. To address this aim, the
Population Adjusted HBI (PAHBI) is a population-weighted version of the HBI
(Table 5).

In the PAHBI, areas with high (or low) HBI values but small populations are
shrunk. These areas can be seen in the PAHBI versus HBI scatter plot of Figure
5 as the line of observations near the 50th percentile. On the other hand,
areas with large populations have their index value made more extreme; ranking
the area higher (or lower).

Similar to the HBI scores, the PAHBI scores (Figure 7) illustrate strong
spatial disparities of unhealthy behaviors. Although remote areas (central
Australia) do have lower PAHBI values than HBI values, even after accounting
for population these regions continue to represent areas with the greatest
potential for impact.

###  3.4 Comparison of indices

The motivation for the creation of an index is to enable policymakers to
prioritize health interventions in a geographically targeted manner. Figure S3
in the Additional File shows the correlations between the AIBIC indices and
the individual unhealthy behaviors.

In this section, we compare how the different indices result in unique
decisions when a policymaker can only, due to feasibility or monetary
constraints, select four areas to apply a new health campaign aimed at
reducing unhealthy behaviors. LABEL:tab:top4 provides these decisions based on
the AIBIC derived from the GSCM and those from the non-model-based methods of
rank sum and PCA applied to the input estimates.

For the non-model-based methods, areas were selected according to the point
estimates only. For the four indices from the AIBIC, the selected areas are
those with the highest posterior probability of having percentiles above the
95th. Only areas with populations greater than 99 were included. Table S5 in
the Additional File provides point estimates and standard errors for the
selected areas in LABEL:tab:top4.

While the non-model-based methods can be used to efficiently select the areas
with the highest prevalence of unhealthy behaviors, these methods do not
accommodate the area-specific uncertainty or population size (LABEL:tab:top4).
Considering PC1 it selects four areas with very high current smoking and
inadequate activity rates, with the standard errors of the estimates for these
unhealthy behaviors exceeding the 98th percentile of the standard errors.
Similarly, PC2 selects four areas with very high current smoking and risky
alcohol consumption rates, coupled with small populations (some under 200).
The standard errors for these areas are all above the 97th percentile of the
standard errors.

Considering the AIBIC, specifically the PAHBI, the standard errors for the
unhealthy behaviors of the areas selected are all below the 92nd percentile
(median 70th) of the standard errors. The PAHBI helps to target over 71100
people, compared to around 11300 and 5770 with PC1 and PC2, respectively. Note
that the selection of which index to use is contingent on the specific
requirements and context of the use case.

Table 6: Comparison of the selected top four areas when using a variety of indices. Only areas with populations greater than 99 were considered. For rank sum, PC1 and PC2, the four areas with the highest value are chosen. For the AIBIC indices, the selected areas have the highest posterior probability of having percentiles above the 95th. Areas that are _italicised_ are also in the top four areas according to their posterior probability of also being ranked in the top 10 areas across Australia for the same index. For each selected area and unhealthy behavior, we provide the actual proportion estimates (and standard errors) as percentages on the left and the nationwide percentiles of the same proportion estimates (and standard errors) on the right. A percentile of 100 should be interpreted as a very high prevalence estimate (or high standard error) of that unhealthy behavior across Australia. |  |  |  |  Inadequate physical activity |  Risky alcohol consumption |  Inadequate diet |  Overweight /obese |  Current smoking |  Inadequate physical activity |  Risky alcohol consumption |  Inadequate diet |  Overweight /obese |  Current smoking
---|---|---|---|---|---|---|---|---|---|---|---|---|---
Index |  Area Name | State | Population | Percentages (%) | Percentiles
Rank Sum | Tuncurry | NSW | 5,684 | 90 (3.0) | 42 (5.7) | 49 (5.5) | 77 (4.4) | 28 (6.5) | 92 (82) | 97 (95) | 78 (94) | 99 (89) | 97 (96)
|  Sussex Inlet - Berrara | NSW | 3,971 | 89 (2.7) | 41 (5.0) | 50 (4.8) | 78 (3.6) | 25 (5.2) | 89 (71) | 96 (93) | 80 (90) | 100 (79) | 95 (94)
| Mannum | SA | 5,548 | 90 (2.5) | 38 (4.6) | 51 (4.4) | 76 (3.5) | 22 (4.4) | 93 (66) | 91 (91) | 86 (88) | 97 (75) | 88 (92)
|  Far South West | QLD | 2,244 | 88 (2.4) | 36 (3.9) | 54 (3.8) | 76 (3.0) | 25 (4.0) | 85 (62) | 83 (80) | 96 (78) | 96 (51) | 95 (89)
PC1 (42% of variation) |  Palm Island | QLD | 1,618 | 94 (16.2) | 33 (23.4) | 61 (22.8) | 64 (24.3) | 76 (24.0) | 100 (100) | 71 (100) | 100 (100) | 34 (100) | 100 (100)
| Yarrabah | QLD | 1,776 | 95 (18.7) | 28 (25.6) | 49 (26.0) | 60 (27.6) | 80 (25.7) | 100 (100) | 41 (100) | 77 (100) | 15 (100) | 100 (100)
|  East Arnhem | NT | 6,085 | 96 (17.7) | 26 (25.3) | 33 (24.7) | 68 (26.6) | 80 (26.1) | 100 (100) | 29 (100) | 2 (100) | 57 (100) | 100 (100)
|  Northern Peninsula | QLD | 1,864 | 93 (9.3) | 30 (15.8) | 62 (16.0) | 67 (16.7) | 58 (19.7) | 100 (98) | 51 (99) | 100 (99) | 48 (99) | 100 (99)
PC2 (30% of variation) | Hume | ACT | 423 | 62 (26.9) | 86 (21.0) | 14 (16.0) | 69 (22.9) | 78 (28.0) | 1 (100) | 100 (100) | 1 (99) | 62 (100) | 100 (100)
|  Melbourne Airport | VIC | 168 | 68 (16.5) | 72 (15.3) | 24 (11.0) | 63 (15.3) | 42 (20.7) | 1 (100) | 100 (99) | 1 (98) | 29 (99) | 99 (99)
|  Brisbane Airport | QLD | 209 | 68 (15.1) | 65 (15.2) | 25 (10.8) | 60 (14.2) | 26 (16.9) | 1 (99) | 100 (99) | 1 (98) | 16 (99) | 96 (99)
|  Rockbank - Mount Cottrell | VIC | 4,970 | 82 (9.3) | 46 (13.3) | 33 (9.7) | 70 (10.7) | 44 (16.3) | 28 (98) | 99 (98) | 2 (97) | 68 (98) | 99 (99)
Factor 1 |  _Main_ _Beach_ | _QLD_ | _3,751_ | 76 (4.1) | 50 (4.9) | 36 (3.8) | 60 (4.5) | 14 (3.7) | 5 (93) | 99 (92) | 3 (78) | 14 (90) | 49 (83)
| Barton | ACT | 1,506 | 67 (4.0) | 43 (4.0) | 49 (3.5) | 59 (3.7) | 11 (2.7) | 1 (93) | 98 (83) | 72 (71) | 13 (80) | 23 (56)
|  _Stockton -_ _Fullerton_ _Cove_ | _NSW_ | _6,764_ | 78 (3.4) | 45 (4.2) | 39 (3.4) | 68 (3.7) | 16 (3.3) | 9 (88) | 99 (88) | 5 (67) | 56 (80) | 59 (75)
| _Coolangatta_ | _QLD_ | _5,738_ | 76 (3.5) | 43 (4.1) | 39 (3.4) | 66 (3.6) | 15 (3.3) | 6 (89) | 97 (84) | 6 (69) | 43 (78) | 55 (75)
Factor 2 | Wingham | NSW | 4,474 | 90 (2.1) | 30 (3.7) | 53 (3.8) | 76 (3.2) | 26 (4.5) | 94 (42) | 53 (76) | 93 (79) | 96 (64) | 97 (92)
| _Tuncurry_ | _NSW_ | _5,684_ | 90 (3.0) | 42 (5.7) | 49 (5.5) | 77 (4.4) | 28 (6.5) | 92 (82) | 97 (95) | 78 (94) | 99 (89) | 97 (96)
| Casino | NSW | 9,958 | 90 (2.0) | 29 (3.3) | 54 (3.5) | 75 (3.0) | 27 (3.8) | 95 (32) | 44 (63) | 96 (70) | 95 (54) | 97 (86)
|  _Torres_ _Strait_ _Islands_ | _QLD_ | _3,045_ | 94 (10.6) | 29 (17.7) | 53 (18.4) | 67 (18.7) | 62 (21.3) | 100 (98) | 48 (100) | 93 (100) | 52 (100) | 100 (99)
HBI | _Tuncurry_ | _NSW_ | _5,684_ | 90 (3.0) | 42 (5.7) | 49 (5.5) | 77 (4.4) | 28 (6.5) | 92 (82) | 97 (95) | 78 (94) | 99 (89) | 97 (96)
| _Forster_ | _NSW_ | _12,592_ | 88 (2.5) | 38 (4.2) | 49 (4.0) | 76 (3.4) | 24 (4.4) | 81 (68) | 90 (87) | 75 (82) | 97 (71) | 93 (92)
|  Forster- Tuncurry Region | NSW | 5,319 | 85 (2.1) | 35 (3.3) | 49 (3.1) | 76 (2.6) | 24 (3.4) | 61 (43) | 82 (63) | 73 (55) | 97 (21) | 93 (78)
|  Laurieton - Bonny Hills | NSW | 14,976 | 84 (2.8) | 39 (4.1) | 50 (3.9) | 78 (2.9) | 22 (3.7) | 51 (77) | 93 (85) | 79 (81) | 100 (49) | 89 (85)
PAHBI |  _Laurieton -_ _Bonny Hills_ | _NSW_ | _14,976_ | 84 (2.8) | 39 (4.1) | 50 (3.9) | 78 (2.9) | 22 (3.7) | 51 (77) | 93 (85) | 79 (81) | 100 (49) | 89 (85)
|  _Wonthaggi -_ _Inverloch_ | _VIC_ | _19,804_ | 83 (2.5) | 37 (3.4) | 49 (3.2) | 77 (2.6) | 22 (3.1) | 38 (65) | 86 (66) | 72 (58) | 98 (19) | 88 (71)
|  _Port_ _Macquarie -_ _East_ | _NSW_ | _23,748_ | 81 (2.7) | 39 (3.5) | 46 (3.2) | 73 (2.9) | 17 (2.9) | 24 (73) | 92 (69) | 51 (58) | 84 (46) | 65 (63)
| Forster | NSW | 12,592 | 88 (2.5) | 38 (4.2) | 49 (4.0) | 76 (3.4) | 24 (4.4) | 81 (68) | 90 (87) | 75 (82) | 97 (71) | 93 (92)

##  4 Discussion

Our results indicate that for the data utilised, ignoring the heteroscedastic
measurement error, relying only on a single shared factor, or accommodating
spatial autocorrelation in only the shared factors or feature-specific
residual errors could lead to incomplete or inefficient policy decisions. Our
solution, the GSCM, is a flexible model-based approach to index creation.
Bolstering its applicability, the GSCM encompasses several previous models [1,
58] used for index creation as special cases. With the computational benefits
of Stan [26] and a fast implementation of the LCAR [27], the GSCM can feasibly
be used to model a wide range of input data with varying degrees of spatial
autocorrelation and measurement error. Although we illustrated our approach
using atlas data from one data source, the GSCM can feasibly be applied to any
geographical health dataset or datasets, with or without uncertainty. It would
be feasible to fit a GSCM where some features have heteroscedastic measurement
error whilst others remain free of such error. Moreover, our model could be
easily incorporated into larger, more complex models such as spatial
structural equation models [16, 66].

As anticipated, the incorporation of the heteroscedastic measurement error
considerably altered the shared factor scores (Figure 3), but the inclusion of
a second factor had minimal influence on the scores of the first factor
(Figure 4). However, the latter comparison was made using a GSCM with spatial
shared factors and spatial feature-specific residual errors. If instead, IID
priors were used for the feature-specific residual errors, the shared factor 1
scores exhibited very different patterns, with the one-factor model producing
scores with higher uncertainty (see Figure S2 in the Additional File). This
result supports the use of traditional SCMs (which assume spatial dependence
in both the shared and feature-specific parts) over BSFMs (which assume
spatial dependence in only one of these). By allowing for spatial dependence
to be captured by a greater number of model terms, our results became more
robust. Further, GSCMs with IID feature-specific residual errors, akin to
Nethery _et al._ [1] BSFM, underperformed irrespective of whether the shared
factors were spatial (Table 3). Future research could investigate whether the
BSFM could be improved by using alternative spatial priors like the BYM [59]
and BYM2 [67]. Consistent with previous research [1, 15], the worst-performing
model was the standard Bayesian factor model which excluded spatial dependence
entirely (Table 3).

Given the negative loadings in factor 1 (Table 4), interpreting the overall
cancer risk profile based on factor 1 or factor 2 alone posed a challenge [9].
This partly motivated the combined indices. Fortunately, the negative loadings
were anticipated based on previous work and our exploratory analysis of the
features (Figure 1). Previous work indicates an inverse relationship between
alcohol consumption and area level socioeconomic disadvantage [68, 69, 4].

The three aims of this work were met by creating the Area Indices of Behaviors
Impacting Cancer (AIBIC). The AIBIC, developed from the output of the GSCM,
provides four indices designed to help identify areas with a high (or low)
prevalence of persons who engage in behaviors that can lead to cancer,
including current smoking, risky alcohol consumption, overweight/obese,
inadequate diet and inadequate physical activity. While factor 1 and 2 are
specific to certain unhealthy behaviors, the HBI is more general and thus more
suitable for future spatial regression models or analyses. The PAHBI could be
used for direct policy decisions, particularly in contexts where population
density is relevant. Our analyses shed light on the different policy decisions
that could be made by examining indices crafted with different methods.
Indeed, indices inherently involve a loss of information, rendering them less
suitable for decisions relating to a single feature (e.g., current smoking
only). In such cases, a direct examination of the specific feature is always
preferable.

The AIBIC (Figures 6 and 7) reveals a pronounced spatial autocorrelation in
unhealthy behaviors and appears to reflect the remoteness [70] and
socioeconomic [9] patterns of Australia (see Figures S6 and S7 in the
Additional File). This pattern could be an artifact of the modelling [31] that
produced the input data for this work. Although this initially raises concerns
regarding the novelty of the indices, subsequent analyses demonstrated that
the AIBIC does not accurately predict the Index of Relative Socio-Economic
Disadvantage from SEIFA [9], with the highest accuracy being a mere 16%. On
the other hand, remoteness [70] proves to be an adequate predictor of AIBIC
scores, evidenced by a maximum adjusted $R^{2}$ of 45%. Thus, future research
could explore the benefits or limitations of using the AIBIC, SEIFA, and
remoteness variables in the same ecological model.

Several limitations were encountered in this study. Firstly, the absence of an
established area-level cancer risk index posed challenges for validation [71].
We sought to address this by directly comparing policy decisions derived from
our indices to the underlying features. Secondly, we faced challenges with
convergence and identifiability in our case study, primarily due to the strong
association of risky alcohol consumption with shared factor 1, which led to
the feature-specific factor for this feature being weakly identified. To
mitigate this, we used more informative hyperpriors and replaced the spatial
prior with an IID prior.

The third limitation is the use of only five features and two shared factors.
This may limit the transferability of the GSCM to more extensive dimension
reduction tasks. Our implementation of the LCAR, although efficient, still
encounters computational challenges when fitting many spatial priors.
Subsequent development of these models could explore the use of direct
sampling [72] or approximation methods such as variational inference [73].

The usefulness of an index is constrained by the quality and breadth of the
input data. Unfortunately, we lacked small area level data on factors like sun
exposure, a significant determinant of melanoma, and dietary factors,
including red meat and fibre consumption, which are closely linked to
colorectal cancer [29, 39]. As these data become available, future work could
update the AIBIC.

Future methodological work could continue to tailor or generalize the GSCM to
include dynamic factor loadings. We experimented with using a unique factor
loading matrix for each of the five remoteness categories across Australia.
Although this significantly improved model fit, it made interpretation
difficult and so wasn’t pursued. Another avenue of future work could be the
use of non-Gaussian distributions for the likelihood, such as the Beta
distribution.

Unlike Nethery _et al._ [1], we did not include population as a feature in our
model. Our method of incorporating population into the PAHBI represents just
one approach, with alternative methods warranting future exploration.
Moreover, there is a lack of development in alternative methods of combining
numerous factors into a single index [7]. Finally, our methodology, like other
factor models, is impacted by the ordering of the features, the chosen
approach for ensuring identifiability and the aggregation level of the spatial
data [74]. Different choices can result in marked variability in inferences,
highlighting possible avenues to improve the robustness of model-based methods
of index creation.

In conclusion, in the current paper, we have proposed a generalised shared
component model that can effectively model complex correlations in
multidimensional health data, which is often derived from atlas platforms such
as the Australian Cancer Atlas [2] and Social Health Atlases of Australia [4].
The GSCM allows for multiple shared factors and heteroscedastic measurement
error and provides uncertainty measures for improved decision-making.

We applied this model to small area-level data on unhealthy behaviors that can
lead to cancer, to develop the first area-level cancer risk index in
Australia; the Area Indices of Behaviors Impacting Cancer (AIBIC).
Successfully addressing the three aims of this study, the AIBIC product
provides a relative summary measure to compare the prevalence of unhealthy
behaviors between areas and is predominantly relevant for analyses and
decision-making for cancer prevention. But, the AIBIC is also applicable to
other disease modelling tasks and depending on the research question, could be
used in conjunction with other indices such as SEIFA. Thus, not only is the
proposed model a useful generalisation of a variety of commonly used model-
based methods for creating indices, but the set of indices created in this
work have wide applicability in public health.

## Backmatter

### Abbreviations

ABS: Australian Bureau of Statistics; AIBIC: Area Indices of Behaviors
Impacting Cancer; BSFM: Bayesian spatial factor model; GSCM: Generalised
shared component model; HPDI: Highest posterior density interval; ICAR:
Intrinsic conditional autoregressive; IID: Independent and identically
distributed; LCAR: Leroux conditional autoregressive; NHMRC: National Health
and Medical Research Council; PC1: Principal Component 1; PC2: Principal
Component 2; PCA: Principal components analysis; SCM: Shared component model;
SEIFA: Socio-Economic Indices for Areas

### Acknowledgments

We thank the Australian Bureau of Statistics (ABS) for designing and
collecting the National Health Survey data. The views expressed in this paper
are those of the authors and do not necessarily reflect the policy of QUT or
CCQ.

### Competing interests

The authors declare that they have no competing interests.

### Declaration of Generative AI and AI-assisted technologies in the writing
process

During the preparation of this work, the authors used ChatGPT4.0 to improve
the clarity and flow of the text. After using this tool, the authors
extensively reviewed and edited the content. The authors take full
responsibility for the content of the publication.

### Author’s contributions

JH: Conceptualisation, Methodology, Software, Writing-Original Draft, Writing-
Review and Editing, Visualisation, Project Administration, Resources, Formal
Analysis, Validation, Investigation, Data Curation. SC: Supervision, Writing-
Review and Editing, Investigation. JC: Supervision, Writing-Review and
Editing, Investigation. PB: Supervision, Writing-Review and Editing,
Investigation. KM: Methodology, Supervision, Writing-Review and Editing,
Investigation.

### Funding

JH was supported by the Queensland University of Technology (QUT) Centre for
Data Science and Cancer Council QLD (CCQ) Scholarship. SC receives salary and
research support from a National Health and Medical Research Council
Investigator Grant (#2008313).

### Additional File

The Additional File provides comprehensive details on the input data, the
proposed GSCM, and the efficient stan implementation of the LCAR prior. It
also includes details on the model results and further plots visualising the
AIBIC product.

## References

  * Nethery et al. [2015] R. C. Nethery, J. L. Warren, A. H. Herring, K. A. Moore, K. R. Evenson, and A. V. Diez-Roux.  A common spatial factor analysis model for measured neighborhood-level characteristics: The Multi-Ethnic Study of Atherosclerosis.  _Health and Place_ , 36:35–46, 2015.  doi:10.1016/j.healthplace.2015.08.009.
  * Duncan et al. [2019] E. W. Duncan, S. M. Cramb, J. F. Aitken, K. L. Mengersen, and P. D. Baade.  Development of the Australian Cancer Atlas: Spatial modelling, visualisation, and reporting of estimates.  _International Journal of Health Geographics_ , 18(1), 2019\.  doi:10.1186/s12942-019-0185-9.
  * Centre for Disease Control [2021] Centre for Disease Control.  PLACES: Local data for better health, 2021.  URL https://www.cdc.gov/places.
  * Public Health Information Development Unit [2018] Public Health Information Development Unit.  Social Health Atlases of Australia, 2018.  URL https://phidu.torrens.edu.au/social-health-atlases.
  * Jahan et al. [2020a] Farzana Jahan, Earl W Duncan, Susanna M Cramb, Peter D Baade, and Kerrie L Mengersen.  Augmenting disease maps: A Bayesian meta-analysis approach.  _Royal Society Open Science_ , 7(8):192151, 2020a.  doi:10.1098/rsos.192151.
  * Jahan et al. [2020b] F. Jahan, E. W. Duncan, S. M. Cramb, P. D. Baade, and K. L. Mengersen.  Multivariate Bayesian meta-analysis: Joint modelling of multiple cancer types using summary statistics.  _International Journal of Health Geographics_ , 19(1), 2020b.  doi:10.1186/s12942-020-00234-0.
  * Nardo et al. [2008] Michela Nardo, Michaela Saisana, Andrea Saltelli, Stefano Tarantola, Anders Hoffmann, Enrico Giovannini, Candice Stevens, Gunseli Baygan, Karsten Olsen, and Sarah Moore.  Handbook on constructing composite indicators.  Report, European Commission, 2008.
  * Ceely [2020] Greg Ceely.  Methods used to develop the Health Index for England: 2015 to 2018.  Report, Office for National Statistics, United Kingdom, 2020.
  * Australian Bureau of Statistics [2016a] Australian Bureau of Statistics.  Technical paper: Socio-Economic Indexes for Areas (SEIFA).  Report, Australian Bureau of Statistics, 2016a.
  * Honaker and Plutzer [2011] James Honaker and Eric Plutzer.  Small area estimation with multiple overimputation.  _Midwest Political Science Association, Chicago_ , 2011.
  * Demsar et al. [2013] Urska Demsar, Paul Harris, Chris Brunsdon, A. Stewart Fotheringham, and Sean McLoone.  Principal component analysis on spatial data: An overview.  _Annals of the Association of American Geographers_ , 103(1):106–128, 2013.  doi:10.1080/00045608.2012.689236.
  * de Asís López et al. [2024] Francisco de Asís López, Celestino Ordóñez, and Javier Roca-Pardiñas.  A generalized additive model (GAM) approach to principal component analysis of geographic data.  _Spatial Statistics_ , 59:100806, 2024.  doi:https://doi.org/10.1016/j.spasta.2023.100806.
  * Hogan and Tchernis [2004] Joseph W. Hogan and Rusty Tchernis.  Bayesian factor analysis for spatially correlated data, with application to summarizing area-level material deprivation from census data.  _Journal of the American Statistical Association_ , 99(466):314–324, 2004.  doi:10.1198/016214504000000296.
  * Norwood et al. [2020] T. A. Norwood, C. Encisa, X. Wang, L. Seliske, J. Cunningham, and P. De.  A Bayesian shared components modeling approach to develop small area indicators of social determinants of health with measures of uncertainty.  _Canadian Journal of Public Health_ , 111(3):342–357, 2020.  doi:10.17269/s41997-020-00321-6.
  * Hedibert et al. [2012] F. Lopes Hedibert, M. Schmidt Alexandra, Salazar Esther, Gómez Mariana, and Achkar Marcel.  Measuring the vulnerability of the Uruguayan population to vector-borne diseases via spatially hierarchical factor models.  _The Annals of Applied Statistics_ , 6(1):284–303, 2012.  doi:10.1214/11-AOAS497.
  * Liu et al. [2005] Xuan Liu, Melanie M. Wall, and James S. Hodges.  Generalized spatial structural equation models.  _Biostatistics_ , 6(4):539–557, 2005.  doi:10.1093/biostatistics/kxi026.
  * Knorr-Held and Best [2001] Leonhard Knorr-Held and Nicola G. Best.  A shared component model for detecting joint and selective clustering of two diseases.  _Journal of the Royal Statistical Society Series A: Statistics in Society_ , 164(1):73–85, 2001.  doi:10.1111/1467-985X.00187.
  * Chidumwa et al. [2021] Glory Chidumwa, Innocent Maposa, Paul Kowal, Lisa K Micklesfield, and Lisa J Ware.  Bivariate joint spatial modeling to identify shared risk patterns of hypertension and diabetes in South Africa: Evidence from WHO SAGE South Africa wave 2.  _International Journal of Environmental Research and Public Health_ , 18(1):359, 2021.  doi:10.3390/ijerph18010359.
  * Law et al. [2020] Jane Law, Matthew Quick, and Afraaz Jadavji.  A Bayesian spatial shared component model for identifying crime-general and crime-specific hotspots.  _Annals of GIS_ , 26(1):65–79, 2020.  doi:10.1080/19475683.2020.1720290.
  * MacNab [2010] Ying C MacNab.  On Bayesian shared component disease mapping and ecological regression with errors in covariates.  _Statistics in Medicine_ , 29(11):1239–1249, 2010\.  doi:10.1002/sim.3875.
  * Cramb et al. [2015] Susanna M. Cramb, Peter D. Baade, Nicole M. White, Louise M. Ryan, and Kerrie L. Mengersen.  Inferring lung cancer risk factor patterns through joint Bayesian spatio-temporal analysis.  _Cancer Epidemiology_ , 39(3):430–439, 2015.  doi:https://doi.org/10.1016/j.canep.2015.03.001.
  * Best et al. [2005] Nicky Best, Sylvia Richardson, and Andrew Thomson.  A comparison of Bayesian spatial models for disease mapping.  _Statistical Methods in Medical Research_ , 14(1):35–59, 2005.  doi:https://doi.org/10.1191/0962280205sm388oa.
  * Stoppa et al. [2022] Giorgia Stoppa, Carolina Mensi, Lucia Fazzo, Giada Minelli, Valerio Manno, Dario Consonni, Annibale Biggeri, and Dolores Catelan.  Spatial analysis of shared risk factors between pleural and ovarian cancer mortality in Lombardy (Italy).  _International Journal of Environmental Research and Public Health_ , 19(6), 2022.  doi:10.3390/ijerph19063467.
  * Mahaki et al. [2011] Behzad Mahaki, Yadollah Mehrabi, Amir Kavousi, Mohammad E Akbari, Thomas Waldhoer, Volker J Schmid, and Mehdi Yaseri.  Multivariate disease mapping of seven prevalent cancers in Iran using a shared component model.  _Asian Pacific Journal of Cancer Prevention_ , 12(9):2353–2358, 2011.
  * Cramb et al. [2020] Susanna Cramb, Earl Duncan, Peter Baade, and Kerrie L. Mengersen.  _A Comparison of Bayesian Spatial Models for Cancer Incidence at a Small Area Level: Theory and Performance_ , pages 245–274.  Springer International Publishing, 2020.  ISBN 978-3-030-42553-1.  doi:10.1007/978-3-030-42553-1_10.
  * Stan Development Team [2023] Stan Development Team.  Stan, 2023.  URL https://mc-stan.org.
  * Leroux et al. [2000] Brian Leroux, Xingye Lei, and Norman Breslow.  _Estimation of Disease Rates in Small Areas: A new Mixed Model for Spatial Dependence_ , pages 179–191.  Springer, 2000.  ISBN 978-1-4612-1284-3.
  * Sung et al. [2021] Hyuna Sung, Jacques Ferlay, Rebecca L Siegel, Mathieu Laversanne, Isabelle Soerjomataram, Ahmedin Jemal, and Freddie Bray.  Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries.  _CA: A Cancer Journal for Clinicians_ , 71(3):209–249, 2021.  doi:10.3322/caac.21660.
  * Whiteman et al. [2015] D. C. Whiteman, P. M. Webb, A. C. Green, R. E. Neale, L. Fritschi, C. J. Bain, D. M. Parkin, L. F. Wilson, C. M. Olsen, C. M. Nagle, N. Pandeya, S. J. Jordan, A. Antonsson, B. J. Kendall, M. C. B. Hughes, T. I. Ibiebele, K. Miura, S. Peters, and R. N. Carey.  Cancers in Australia in 2010 attributable to modifiable factors: Introduction and overview.  _Australian and New Zealand Journal of Public Health_ , 39(5):403–407, 2015.  doi:10.1111/1753-6405.12471.
  * Hermes and Poulsen [2012] K. Hermes and M. Poulsen.  Small area estimates of smoking prevalence in London. Testing the effect of input data.  _Health and Place_ , 18(3):630–638, 2012.  doi:10.1016/j.healthplace.2011.12.010.
  * Hogg et al. [2023] James Hogg, Jessica Cameron, Susanna Cramb, Peter Baade, and Kerrie Mengersen.  Mapping the prevalence of cancer risk factors at the small area level in Australia.  _International Journal of Health Geographics_ , 22(1):37, 2023.  doi:10.1186/s12942-023-00352-5.
  * Das et al. [2023] Sumonkanti Das, Bernard Baffour, Alice Richardson, Susanna Cramb, and Stephen Haslett.  Daily smoking prevalence for small domains in Australia.  _Research Square preprint_ , 2023.  doi:10.21203/rs.3.rs-2829471/v1.
  * Ashraf et al. [2019] Khalid Ashraf, Chirk Jenn Ng, Chin Hai Teo, and Kim Leng Goh.  Population indices measuring health outcomes: A scoping review.  _Journal of Global Health_ , 9(1):010405–010405, 2019.  doi:10.7189/jogh.09.010405.
  * Colditz et al. [2000] G. A. Colditz, K. A. Atwood, K. Emmons, R. R. Monson, W. C. Willett, D. Trichopoulos, and D. J. Hunter.  Harvard report on cancer prevention volume 4: Harvard Cancer Risk Index.  _Cancer Causes and Control_ , 11(6):477–488, 2000\.  doi:10.1023/A:1008984432272.
  * Krishna Rao et al. [2016] S. Krishna Rao, G. C. Mejia, R. M. Logan, M. Kulkarni, V. Kamath, D. J. Fernandes, S. Ray, and K. Roberts-Thomson.  A screening model for oral cancer using risk scores: Development and validation.  _Community Dentistry and Oral Epidemiology_ , 44(1):76–84, 2016.  doi:10.1111/cdoe.12192.
  * Australian Bureau of Statistics [2018] Australian Bureau of Statistics.  National Health Survey: First results methodology, 2018.  URL https://www.abs.gov.au/methodologies/national-health-survey-first-results-methodology/2017-18.
  * Australian Bureau of Statistics [2017a] Australian Bureau of Statistics.  Microdata: National Health Survey 2017-18 [DataLab], 2017a.
  * Australian Bureau of Statistics [2011] Australian Bureau of Statistics.  Australian Statistical Geography Standard (ASGS), 2011.  URL https://www.abs.gov.au/websitedbs/d3310114.nsf/home/australian+statistical+geography+standard+(asgs).
  * World Cancer Research Fund and American Institute for Cancer Research [2018] World Cancer Research Fund and American Institute for Cancer Research.  Exposures, risk factors and cancer, 2018.  URL https://www.wcrf.org/diet-and-cancer/exposures/.
  * Australian Institute of Health Welfare [2019] Australian Institute of Health Welfare.  Australian burden of disease study: impact and causes of illness and death in australia 2015.  Report, Australian Institute of Health Welfare, 2019.  URL https://www.aihw.gov.au/reports/burden-of-disease/burden-disease-study-illness-death-2015.
  * Rezende et al. [2020] Leandro F. M. Rezende, Eliana Murata, Beatriz Giannichi, Luciana Yuki Tomita, Gabriela Arantes Wagner, Zila M. Sanchez, Carlos Celis-Morales, and Gerson Ferrari.  Cancer cases and deaths attributable to lifestyle risk factors in chile.  _BMC Cancer_ , 20(1):693, 2020.  doi:10.1186/s12885-020-07187-4.
  * Hogg [2024] James Hogg.  gscm, 2024.  URL https://github.com/JamieHogg-depo/gscm.
  * Australian Bureau of Statistics [2023] Australian Bureau of Statistics.  ERP by SA2 (ASGS 2016), Age and Sex, 2001 Onwards, 2023.  URL https://www.abs.gov.au/about/data-services/data-explorer.
  * National Health and Medical Research Council [2020] National Health and Medical Research Council.  Australian guidelines to reduce health risks from drinking alcohol.  Report, National Health and Medical Research Council, 2020.
  * Health and Council [2013] National Health and Medical Research Council.  Australian dietary guidelines.  Report, National Health and Medical Research Council, 2013.
  * Australian Bureau of Statistics [2017b] Australian Bureau of Statistics.  4363.0 - National Health Survey: Users’ Guide, 2017-18, 2017b.
  * Department of Health [2014] Department of Health.  Physical activity and exercise guidelines for all Australians, 2014\.  URL https://www.health.gov.au/resources/collections/collection-of-physical-activity-and-sedentary-behaviour-guidelines-for-all-ages.
  * Carroll et al. [2006] Raymond J. Carroll, David Ruppert, Leonard A. Stefanski, and Ciprian M. Crainiceanu.  _Measurement Error in Nonlinear Models : A Modern Perspective, Second Edition_.  CRC Press LLC, London, United Kingdon, 2006.  ISBN 9781420010138.
  * Li et al. [2009] Yi Li, Haicheng Tang, and Xihong Lin.  Spatial linear mixed models with covariate measurement errors.  _Statistica Sinica_ , 19(3):1077–1093, 2009.  doi:https://www.jstor.org/stable/24308946.
  * Huque et al. [2014] Md Hamidul Huque, Howard Bondell, and Louise Ryan.  On the impact of covariate measurement error on spatial regression modelling.  _Environmetrics_ , 25(8):560–570, 2014.  doi:10.1002/env.2305.
  * Wang and Wall [2003] Fujun Wang and Melanie M. Wall.  Generalized common spatial factor model.  _Biostatistics_ , 4(4):569–582, 2003.  doi:10.1093/biostatistics/4.4.569.
  * Papastamoulis and Ntzoufras [2022] Panagiotis Papastamoulis and Ioannis Ntzoufras.  On the identifiability of Bayesian factor analytic models.  _Statistics and Computing_ , 32(2):23, 2022.  doi:10.1007/s11222-022-10084-4.
  * Ma and Liu [2021] Yucong Ma and Jun S. Liu.  On posterior consistency of Bayesian factor models in high dimensions.  _arXiv preprint arXiv:2006.01055v3_ , 2021.
  * Ghosh and Dunson [2009] J. Ghosh and D. B. Dunson.  Default prior distributions and efficient posterior computation in Bayesian factor analysis.  _The Journal of Computational and Graphical Statistics_ , 18(2):306–320, 2009.  doi:10.1198/jcgs.2009.07145.
  * Aguilar and West [2000] Omar Aguilar and Mike West.  Bayesian dynamic factor models and portfolio allocation.  _Journal of Business and Economic Statistics_ , 18(3):338–357, 2000.
  * Shin and Ferreira [2023] Hwasoo Shin and Marco A. R. Ferreira.  Dynamic ICAR spatiotemporal factor models.  _Spatial Statistics_ , 56:100763, 2023.  doi:https://doi.org/10.1016/j.spasta.2023.100763.
  * Lopes et al. [2011] Hedibert Freitas Lopes, Dani Gamerman, and Esther Salazar.  Generalized spatial dynamic factor models.  _Computational Statistics and Data Analysis_ , 55(3):1319–1330, 2011.  doi:https://doi.org/10.1016/j.csda.2010.09.020.
  * Mezzetti [2012] Maura Mezzetti.  Bayesian factor analysis for spatially correlated data: Application to cancer incidence data in Scotland.  _Statistical Methods and Applications_ , 21:49–74, 2012\.  doi:https://doi.org/10.1007/s10260-011-0177-9.
  * Besag et al. [1991] Julian Besag, Jeremy York, and Annie Mollié.  Bayesian image restoration, with two applications in spatial statistics.  _Annals of the Institute of Statistical Mathematics_ , 43(1):1–20, 1991.  doi:10.1007/BF00116466.
  * Donegan [2021] C. Donegan.  Building spatial Conditional Autoregressive (CAR) models in the Stan programming language.  _OSF Preprints_ , 2021.  doi:https://doi.org/10.31219/osf.io/3ey65.
  * Morris et al. [2019] M. Morris, K. Wheeler-Martin, D. Simpson, S. J. Mooney, A. Gelman, and C. DiMaggio.  Bayesian hierarchical spatial models: Implementing the Besag York Mollié model in Stan.  _Spatial and Spatio-temporal Epidemiology_ , 31(1), 2019\.  doi:10.1016/j.sste.2019.100301.
  * Aki et al. [2021] Vehtari Aki, Gelman Andrew, Simpson Daniel, Carpenter Bob, and Bürkner Paul-Christian.  Rank-normalization, folding, and localization: An improved $\widehat{R}$ for assessing convergence of MCMC (with discussion).  _Bayesian Analysis_ , 16(2):667–718, 2021.  doi:10.1214/20-BA1221.
  * Vehtari et al. [2017] A. Vehtari, A. Gelman, and J. Gabry.  Practical Bayesian model evaluation using leave-one-out cross-validation and WAIC.  _Statistics and Computing_ , 27(5):1413–1432, 2017.  doi:10.1007/s11222-016-9696-4.
  * Moran [1950] Patrick Alfred Pierce Moran.  Notes on continuous stochastic phenomena.  _Biometrika_ , 37(1):17–23, 1950.  doi:10.2307/2332142.
  * Lopes and West [2004] Hedibert Freitas Lopes and Mike West.  Bayesian model assessment in factor analysis.  _Statistica Sinica_ , pages 41–67, 2004.  doi:https://www.jstor.org/stable/24307179.
  * Roman and Brandt [2021] Zachary J. Roman and Holger Brandt.  A latent auto-regressive approach for Bayesian structural equation modeling of spatially or socially dependent data.  _Multivariate Behavioral Research_ , pages 1–25, 2021.  doi:10.1080/00273171.2021.1957663.
  * Riebler et al. [2016] Andrea Riebler, Sigrunn H Sørbye, Daniel Simpson, and Håvard Rue.  An intuitive Bayesian spatial model for disease mapping that accounts for scaling.  _Statistical Methods in Medical Research_ , 25(4):1145–1165, 2016.  doi:10.1177/0962280216660421.
  * Australian Institute of Health and Welfare [2020] Australian Institute of Health and Welfare.  National drug strategy household survey 2019.  Report, AIHW, 2020.
  * Collins [2016] S. E. Collins.  Associations between socioeconomic factors and alcohol outcomes.  _Alcohol Research_ , 38(1):83–94, 2016.  doi:https://doi.org/10.1007/s10462-011-9236-8.
  * Australian Bureau of Statistics [2016b] Australian Bureau of Statistics.  1270.0.55.005 - Australian Statistical Geography Standard (ASGS): Volume 5 - remoteness structure, july 2016, 2016b.  URL https://www.abs.gov.au/ausstats/abs@.nsf/Latestproducts/1270.0.55.005Main%20Features15July%202016?opendocument&tabname=Summary&prodno=1270.0.55.005&issue=July%202016&num=&view=.
  * Kaltenthaler et al. [2004] E. Kaltenthaler, R. Maheswaran, and C. Beverley.  Population-based health indexes: A systematic review.  _Health Policy_ , 68(2):245–55, 2004.  doi:10.1016/j.healthpol.2003.10.005.
  * White et al. [2019] Gentry White, Dongchu Sun, and Paul Speckman.  Direct sampling of bayesian thin-plate splines for spatial smoothing.  _arXiv preprint arXiv:1906.05575_ , 2019.  doi:arXiv:1906.05575.
  * Fox and Roberts [2012] Charles W. Fox and Stephen J. Roberts.  A tutorial on variational Bayesian inference.  _Artificial Intelligence Review_ , 38:85–95, 2012.
  * Buckingham et al. [2024] William R. Buckingham, W. Ryan Powell, Sarah Anne Keller, Kelia J. Hansmann, and Amy J. H. Kind.  Bigger isn’t better: Why small area geographies are best for actionable index development.  _Papers in Applied Geography_ , pages 1–7, 2024.  doi:10.1080/23754931.2024.2312192.