# Heterogeneous Effects in the Built Environment

Adam T. Peterson
Department of Biostatistics
University of Michigan
Ann Arbor, MI 48104
<EMAIL_ADDRESS>
&Brisa N. Sánchez
Department of Biostatistics and Epidemiology
Drexel University
Philadelphia, PA
<EMAIL_ADDRESS>
Emma Sanchez-Vaznaugh
Department of Public Health
University of California San-Francisco
San Francisco, CA

###### Abstract

We present an approach to estimate distance-dependent heterogeneous
associations between point-referenced exposures to built environment
characteristics and health outcomes. By estimating associations that depend
non-linearly on distance between subjects and point-referenced exposures, this
method addresses the modifiable area-unit problem that is pervasive in the
built environment literature. Additionally, by estimating heterogeneous
effects, the method also addresses the uncertain geographic context problem.
The key innovation of our method is to combine ideas from the non-parametric
function estimation literature and the Bayesian Dirichlet process literature.
The former is used to estimate nonlinear associations between subject’s
outcomes and proximate built environment features, and the latter identifies
clusters within the population that have different effects. We study this
method in simulations and apply our model to study heterogeneity in the
association between fast food restaurant availability and weight status of
children attending schools in Los Angeles, California.

_Keywords_ Bayesian Non-Parametric $\cdot$ Built Environment $\cdot$
Heterogeneous Effects.

##  1 Introduction

The relationship between amenities in or near residential, work or
school–neighborhood environments and health is receiving increasing attention,
given that these environments can influence health-related behaviors and
subsequent outcomes. Where spatial proximity to supermarkets is associated
with diet, so too are recreational facilities associated with physical
activity and fast food restaurants near schools associated with child obesity
(Baek et al., 2016a, 2017; Kaufman et al., 2019; Kern et al., 2017). Work in
this area has been limited by the lack of knowledge of what geographic units
are most relevant for exposure assessment, i.e. the well known modifiable unit
areal problem (MAUP) (Fotheringham and Wong, 1991; Spielman and Yoo, 2009;
Wong, 2009; Guo and Bhat, 2004; Ji et al., 2009; James et al., 2014).
Additionally, there may also be measured or unmeasured person-level behaviors
or characteristics that give rise to the “uncertain geographic context
problem” (UGCP) (Macintyre et al., 2002; Kwan, 2013, 2018). Whereas the former
establishes that using different spatial units or spatial scales to define
exposure measures will yield different estimates of association, the latter
acknowledges that the most relevant spatial unit may differ from place to
place or subject to subject due to place or person characteristics such as
predominant transport modes in a given area or vehicle ownership, among
others.

Recent work has begun addressing these issues by foregoing the pre-
specification of the spatial unit used to construct exposure metrics (Baek et
al., 2016a; Peterson and Sanchez, 2018). Instead, the association between
proximity to amenities of interest, broadly referred to as built environment
features (BEFs), and subjects’ outcomes is estimated as a continuous function
of distance between subjects and amenities. Whereas typical models regress the
outcome on a BEF metric that depends on a pre-defined scale, these new methods
use all the pair-wise distances between subjects and BEFs as inputs to the
model. Specifically, in order to address the MAUP, an idealized smooth
function $f(d)$ is used to represent the association between the health
outcome of interest and a single BEF placed at distance $d>0$ from the
subject. Having $f(d)$ as the objective of inference enables the visualization
of whether and how the association between availability of amenities and
outcomes dissipates with distance, as well as estimation of the spatial scale,
defined as the distance at which the association is negligible, i.e.
$d:f(d)=0$.

The function $f(d)$ has been modeled in different ways: Baek et al. (2016a)
estimated $f(d)$ non-parametrically by first discretizing the distances into a
grid, and using the count of distances within bins defined by the grid as
predictors in a Distributed Lag Model (DLM), i.e., the count of distances
within each bin are conceptualized as distributed lag predictors, indexed by
the corresponding value of the grid. The coefficients corresponding to each
distributed lag predictor are smoothed using splines, yielding estimates of
$f(d)$ at the values of $d$ used to construct the grid. Alternatively,
Peterson et al. (2021) modeled $f(d)$ parametrically, typically using
exponential functions to enforce the substantive belief that the association
between health outcomes and spatial availability of amenities monotonically
decays across distance, e.g. $f(d)\propto\exp(-\frac{d}{\theta})$. However,
the estimation of $f(d)$ at the population level, as the previous methods
propose, fails to account for the concerns the UGCP raises regarding
unmeasured person-level behaviors or place-level factors that may determine
subject – or location – specific spatial association.

Building upon their work in DLMs, Baek et al. (2016b, 2017) constructed a
hierarchical DLM (HDLM) allowing for the estimated $f(d)$ to vary between
subjects and or locations, according to pre-specified groups (e.g., different
$f(d)$ by sex), as well as unexplained variation in the association (i.e.,
using the idea of random coefficients to estimate $f(d)$ for individual
subjects). However, the HDLM approach, has some disadvantages: (1) it uses
discretized distances to estimate association across space, unnecessarily
coarsening the exposure information; (2) it requires pre-specifying the groups
where heterogeneity in the association may occur (covariates and or subjects);
and (3) by enforcing that heterogeneity in the association estimates to occur
at the subject-level through random effects, it loses possible gains in
precision that could result from pooling subjects with similar levels of
association.

Motivated by the desire to identify schools where pupils may be at greater
risk of obesity related to the proximity of fast food restaurants (FFRs), we
propose a model that clusters schools-specific association curves, $f(d)$
according to the strength of association between the spatial proximity of
nearby FFRs and child obesity. Clustering provides investigators and
policymakers with a greater understanding of the kinds of relationships that
exist between students and their environment as well as identifies schools
where students may be at greater risk, as identifying risk groups may help
prioritize population level interventions. The data for this motivating study
consists of body weight status of children nested within schools across Los
Angeles County during academic years 2001-2008. Distances between schools and
FFRs are calculated from geocoded school addresses, supplied by the California
Department of Education, and geocoded FFR business addresses from the National
Establishment Time Series Database(Walls, 2013)

Our method uses the Dirichlet Process Mixture (DPM) prior and a spline basis
function expansion to non-parametrically estimate both the number of cluster-
BEF effects, and the nonlinear association functions across space,
respectively. We name our method the Spatial Aggregated Predictor - Dirichlet
Process, to reflect this dual non parametric estimation, but refer to it more
generally as STAP-DP, given its potential for also modeling temporal
exposures. Our approach is inspired by the work of Rodriguez et al. (2014) and
Ray and Mallick (2006) on clustering functions using DPM family priors. We use
the penalized spline approach developed by O’Sullivan (1986) and further
popularized by (Wahba, 1990; Wood, 2017) to construct the estimates of the
association functions, and use the DPM to cluster the spline coefficients.

Section 2 describes the model that estimates homo- and heterogeneous BEF
effects. Section 3 studies the performance of the STAP-DP model in a variety
of simulated data settings and discusses how the results may inform normative
practice. Section 4 describes the application of the STAP-DP model to the
motivating study on child obesityin Los Angeles. We conclude our work with a
discussion of the model and future directions to explore.

##  2 Model

We now introduce the STAP-DP framework, describing how we incorporate the
estimation of heterogeneous BEF effects into a regression framework. We limit
our discussion to the estimation of only one BEF’s effects in space, FFRs for
example, as the extension to multiple BEFs is straightforward. We organize our
discussion into four parts. First, we build intuition for our approach by
defining the STAP estimated via spline basis functions at the population
level, i.e., homogeneous effect. Then, we define how to extend the STAP model
to estimate heterogeneous effects – at the latent cluster level – for a
univariate outcome. In the final two sections we generalize the clustering
framework for repeated outcome measures and discuss estimation.

###  2.1 The STAP Model

Suppose a continuous outcome $Y_{i}$ ($i=1,...N$) and corresponding covariates
$\bm{X}_{i}\in\mathbb{R}^{n\times p}$ are observed for a sample of $N$
subjects. Additionally, spatial data, $\mathcal{D}_{i}$ which contains
distances, $d$, between subject $i$ and all FFRs within some substantively
determined radius $R$, are also measured. The inferential objective is to
estimate function $f(d)$, which represents the expected difference in the
outcome associated with placing a single FFR at distance $d$ after adjusting
for covariates $\bm{X}_{i}$. Defining
$F(\mathcal{D}_{i}):=\sum_{d\in\mathcal{D}_{i}}f(d)$, as the aggregated FFR
effect under the assumption of additivity, we complete the initial STAP model
formulation:

$\displaystyle Y_{i}$
$\displaystyle=\bm{X}_{i}^{T}\bm{\gamma}+F(\mathcal{D}_{i})+\epsilon_{i},$ (1)
$\displaystyle\epsilon_{i}$ $\displaystyle\stackrel{{\scriptstyle
iid}}{{\sim}}N(0,\sigma^{2}),$

where $\epsilon_{i}$ is the residual error, with variance $\sigma^{2}$.

As mentioned in Section 1, there are a number of approaches to model $f(d)$.
In this work, we propose to model $f(d)$ as a linear combination of basis
functions, $\\{\phi\\}_{l=1}^{L}$, which allows us to rewrite
$F(\mathcal{D}_{i})$ as follows:

$\displaystyle F(\mathcal{D}_{i})$
$\displaystyle=\sum_{d\in\mathcal{D}_{i}}f(d)=\sum_{d\in\mathcal{D}_{i}}\sum_{l=1}^{L}\beta_{l}\phi_{l}(d),$
(2)

where $\phi_{l}(d)$ is the evaluation of the distance through the $l$th basis
function and $\beta_{l}$ is the corresponding regression coefficient. In this
work we use $L$ spline basis functions defined across a set of equally spaced
knots, though other knot placements or basis functions could be used. In order
to avoid over fitting when $L$ is large, the regression coefficients are
regularized through the use of a quadratic penalty on $\bm{\beta}$ implemented
through a smoothing matrix $S$ and tuned by penalty parameter $\tau$. We use
the difference penalty matrices of Eilers and Marx (1996), a widely used
spline penalty formulation.

Within a Bayesian paradigm, this penalty is equivalent to specifying a
multivariate normal prior with improper precision matrix $\tau S$. We adopt a
variant of this Bayesian approach and, to improve computational efficiency in
our more complex model formulations discussed in the next subsection, we first
transform the spline basis function expansion matrix, $\Phi(d)$, such that the
transformed coefficients can have independent normal priors (Wood, 2004,
2016). While centering constraints are often imposed on $\Phi(d)$ to avoid
collinearity with the intercept in $\bm{X}$, this constraint is not needed in
our model (see supplementary material). Given that $r_{S}=$rank($S)<L$, two
precision parameters for the priors are used, one for the first $r_{S}$
coefficients and a second for the last $L-r_{S}$ coefficients:

$\displaystyle\bm{\beta}_{1}$ $\displaystyle\sim
MVN_{r_{S}}\left(\bm{0},\sigma^{2}\tau_{1}^{-1}\bm{I}_{r_{S}}\right)\quad\bm{\beta}_{2}\sim
MVN_{L-r_{S}}(\bm{0},\sigma^{2}\tau_{2}^{-1}\bm{I}_{L-r_{S}})$ (3)
$\displaystyle\tau_{z}$ $\displaystyle\stackrel{{\scriptstyle
iid}}{{\sim}}\text{Gamma}(a_{\tau},b_{\tau})\quad z=1,2.$

In $(\ref{eqn:prior})$ we denote $\bm{\beta}_{z}$, $z=1,2$, as the regression
coefficients in the penalty range and null space, respectively.
Correspondingly, $\tau_{1}$ and $\tau_{2}$ are the respective precisions for
these separate subsets of $\bm{\beta}$. For ease of further exposition we
define $\bm{\Lambda}$ as the diagonal covariance matrix which has
$\tau_{1}^{-1}$ as the first $r_{S}$ diagonal elements and $\tau_{2}^{-1}$ as
the last $L-r_{S}$ diagonal elements, so that the prior in (3) can be written
simply as $\bm{\beta}\sim MVN_{L}(\bm{0},\sigma^{2}\bm{\Lambda})$. We place
independent conjugate Gamma priors on $\bm{\tau}=(\tau_{1},\tau_{2})$ so that
both $\bm{\beta}$’s and $\bm{\tau}$’s conditional posterior distributions are
available in closed form.

###  2.2 STAP-DP with Univariate Outcomes

In alignment with this work’s goal to estimate heterogeneous effects, we
replace $F_{i}(\mathcal{D}_{i})$ with $F_{i}(\mathcal{D}_{i})$ in (1) while
allowing for clustering in the $f_{i}(d)$. Given that $f_{i}(d)$ is
represented by the fixed spline functions and random coefficients
$\bm{\beta}$, we implement this clustering goal by placing a DP prior on the
vector of regression coefficients, $\bm{\beta}$, and associated penalty
parameter, $\bm{\tau}$:

$\displaystyle(\bm{\beta},\bm{\tau})$ $\displaystyle\sim P$ (4) $\displaystyle
P$ $\displaystyle\sim DP(\alpha,P_{0})$ $\displaystyle P_{0}$
$\displaystyle\equiv
MVN_{L}\left(\bm{0},\sigma^{2}\bm{\Lambda}\right)\times\prod_{z=1}^{2}\text{Gamma}(a_{\tau},b_{\tau}).$

In (4), $P$ is a random measure drawn from Dirichlet Process
$DP(\alpha,P_{0})$, where $\alpha>0$ is a concentration parameter reflecting
the variability of distribution $P$ around base measure $P_{0}$ (Ferguson,
1973; Gelman et al., 2013). $P_{0}$ is chosen to retain the prior previously
discussed in (3).

By placing the DP prior on $(\bm{\beta},\bm{\tau})$, clustering is induced on
the $f_{i}(d)$ as can be seen from the stick breaking construction of the DP:
$P=\sum_{k=1}^{\infty}\pi_{k}\delta_{(\bm{\beta}^{\star},\bm{\tau}^{\star})}(\cdot)$.
In this representation $\pi_{k}$ represents the probability the $i$th
observation is assigned to the $k$th exposure function and $\delta(\cdot)$ is
the dirac-delta function. Each $\pi_{k}$, itself is composed of the “broken
sticks” created from variables drawn from a Beta distribution:
$\pi_{k}=v_{k}\prod_{u<k}(1-v_{u});v_{k}\sim\text{Beta}(1,\alpha)$.

Combining all these pieces together, our proposed STAP-DP model for univariate
outcome $Y_{i}$ takes the following form:

$\displaystyle Y_{i}$
$\displaystyle=\bm{X}_{i}\bm{\gamma}+\sum_{d\in\mathcal{D}_{i}}\sum_{l=1}^{L}\beta_{il}\phi_{l}(d)+\epsilon_{i}$
(5) $\displaystyle\epsilon_{i}$ $\displaystyle\stackrel{{\scriptstyle
iid}}{{\sim}}N(0,\sigma^{2})$ $\displaystyle(\bm{\beta},\bm{\tau})$
$\displaystyle\sim P$ $\displaystyle P$ $\displaystyle\sim DP(\alpha,P_{0})$
$\displaystyle P_{0}$ $\displaystyle\equiv
MVN_{L}\left(\bm{0},\sigma^{2}\bm{\Lambda}\right)\times\prod_{z=1}^{2}\text{Gamma}(a_{\tau},b_{\tau}).$

A final comment is warranted regarding the choice of the number and placement
of the $L$ knots in constructing the splines. While our approach follows
previous work in placing a sufficient number of knots equally across the
domain of observed distances, deciding what number of knots is “sufficient”
requires greater statistical judgement than in standard applications. Clusters
may be more difficult to detect when the dimension on which clusters are
formed (i.e., number of coefficients) is large and the between-cluster
differences are small (low signal effects). Conversely, more clusters may be
identified in a setting with a stronger signal and greater number of knots.
Thus, $L$ must be chosen to balance accuracy in both function estimation and
cluster discrimination.

###  2.3 STAP-DP with Repeated Measurements

Extending (5) to correlated outcomes, we consider the setting in which
subjects are measured repeatedly over time, for $j=1,...,n_{i}$ occasions.
This results in outcome $Y_{ij}$ ($i=1,...,N,j=1,...,n_{i})$ modeled as a
function of covariates $\bm{X}_{ij}$, and their corresponding coefficients
$\bm{\gamma}$. The distance set adopts the new visit-specific index as well,
i.e., $\mathcal{D}_{ij}$, indicating it may vary over time; for instance FFRs
may open and close between measurement occasions. Finally, a subset of
$\bm{X}_{ij}$, $\bm{Z}_{ij}$, is included in the model, along with subject-
specific coefficients $\bm{b}_{i}\sim MVN(0,\Sigma)$ to account for within
subject variability in standard fashion (Fitzmaurice et al., 2008). Augmenting
(4) accordingly, we arrive at our final model:

$\displaystyle Y_{ij}$
$\displaystyle=\bm{X}_{ij}^{T}\bm{\gamma}+\sum_{d\in\mathcal{D}_{ij}}\sum_{l=1}^{L}\beta_{il}\phi_{l}(d)+\bm{Z}_{ij}^{T}\bm{b}_{i}+\epsilon_{i}$
(6) $\displaystyle\bm{b}_{i}$ $\displaystyle\stackrel{{\scriptstyle
iid}}{{\sim}}N(\bm{0},\Sigma)$ $\displaystyle\epsilon_{i}$
$\displaystyle\stackrel{{\scriptstyle iid}}{{\sim}}N(0,\sigma^{2})$
$\displaystyle(\bm{\beta},\bm{\tau})$ $\displaystyle\sim P$ $\displaystyle P$
$\displaystyle\sim DP(\alpha,P_{0})$ $\displaystyle P_{0}$
$\displaystyle\equiv
MVN_{L}\left(\bm{0},\sigma^{2}\bm{\Lambda}\right)\times\prod_{z=1}^{2}\text{Gamma}(a,b).$

###  2.4 Estimation

In order to fit models of the form described in (4) and (6), we truncate the
DP so that a blocked Gibbs sampler can be used to draw samples from the
posterior (Gelman et al., 2013). While this sampler is fairly straightforward,
it bears mentioning that $\Phi(d)$ has to be adjusted at each iteration of
sampling so that any DP components associated with 0 or some small number of
observations are not included in the usual matrix inversion used to estimate
the mean of the conditional posterior distribution for the regression
coefficients, $\bm{\beta}^{*}=[\bm{\gamma},\bm{\beta}]^{T}$. Instead,
coefficients for those low-member cluster components are sampled with draws
from the prior. For example, if on the $m$th iteration, none of the $N$
observations are assigned to the $k$th DP component, then the samples of the
spline regression coefficients for that iteration, $\bm{\beta}_{k}^{(m)}$, are
drawn from a $MVN_{L}(\bm{0},\sigma^{2}\bm{\Lambda}_{k})$ prior, where
$\bm{\Lambda}_{k}$ is the cluster specific covariance matrix, and the columns
of zeros that would otherwise be included in $\Phi(d)$ are omitted.

We present the closed form conditional posteriors and associated algorithm in
the Supplementary Material. Our algorithm is implemented in C++ which can be
called from our R package rstapDP (Peterson, 2020a; Peterson and Sánchez, ).

For both our simulations and California data we use rstapDP to fit the STAP-DP
in R (v.4.0.2) R Core Team (2019) on a MacOS Catalina operating system with a
2.8 GhZ Quad-Core Intel Core i7 processor.

##  3 Simulations

###  3.1 Simulation Design

For a given sample size, the ability of the STAP-DP model to correctly
classify subjects depends on (a) the proportion of subjects belonging to that
cluster, (b) the difference in the $f_{i}(d)$ functional forms, and (c) the
distribution of distances (i.e., exposure information) present within each
cluster. As the first of these three principles follows straightforward
sample-size intuitions, in this section we study the STAP-DP’s ability to
correctly recover cluster specific functions,$f_{i}(d)$, and cluster
partitions in the latter two settings. Using simulated data we vary: (i)
cluster effect size and (ii) distance distributions in order to see how these
may impact correct cluster classification. We focus on evaluating cluster
classification accuracy as it is the upstream predictor of all remaining model
components, like the estimation of the $f_{i}(d)$, which are all standard
Bayes estimators conditional on the correct cluster classification.

We evaluate our method’s ability to correctly classify subjects using a
partition loss function developed by (Binder, 1978) and used regularly in DP
and other mixture model applications where label-switching may be of concern
(Lau and Green, 2007; Wade et al., 2018; Rodriguez et al., 2008). Our
employment of the loss function equally weights correct and incorrect
classification, using the subjects’ true and estimated class indicators,
$\zeta_{i},\hat{\zeta}_{i}$, respectively:

$\displaystyle\mathcal{\psi}(\bm{\zeta},\hat{\bm{\zeta}})=\sum_{(i,i^{\prime});i<i^{\prime}<N}I(\zeta_{i}=\zeta_{i^{\prime}},\hat{\zeta}_{i}\neq\hat{\zeta}_{i^{\prime}})+I(\zeta_{i}\neq\zeta_{i^{\prime}},\hat{\zeta}_{i}=\hat{\zeta}_{i^{\prime}}).$
(7)

Conceptually, (7) tallies the number of times that observations $i$ and
$i^{\prime}$ are incorrectly assigned to different clusters, when they in fact
belong in the same cluster, as well as tallies when they are incorrectly
assigned to the same cluster.

In each simulation setting discussed below, we generate 25 datasets and then
fit the STAP-DP model shown in (4), truncating the DP at K= 50 and using
weakly informative Gamma(1,1) priors on $\sigma^{-2},\alpha,\tau_{1}$ and
$\tau_{2}$, respectively. We draw 2000 samples from the posterior distribution
for inference via Gibbs Sampling using rstapDP after discarding 2000 initial
samples for burn-in. Across all 25 simulations we evaluate the loss (7) across
all $M=2000$ iterations of the posterior samples drawn via Gibbs sampling.
Given that the loss function does not have a standard range, we normalize the
loss results by the maximum loss across all simulation settings, so as to make
the results more interpretable relative to one another. We have organized the
files used to run the simulations in the STAPDPSimulations R package available
via Github.

###  3.2 Cluster Effect Size

Our first simulation study focuses on model performance as a function of the
difference between two clusters’ $f(d)$, defined in (8) below, with each
observation having a 50% probability of being assigned to either of these
clusters. The cluster function set-up is intended to mimic a hypothetical high
and low risk population scenario, in which subjects with equivalent exposure
to the same BEFs experience different effects according to which risk
population they belong. For each subject we generate a random number of
distances uniformly so that the average number of BEFs is 15. Conditional on
the number of distances, the distances themselves are then generated according
to the “Skew” distribution shown in Figure S1 in the Supplementary
Information. This distribution was selected in order to test our model’s
performance under a “worst case scenario”, given that with this distribution
there is relatively less exposure information at the distances where the BEF
effects are non-zero. Specifically, the generative model takes the following
form:

$\displaystyle Y_{i}$
$\displaystyle=26+.5Z_{i}+\sum_{d\in\mathcal{D}_{i}}f_{\zeta_{i}}(d)+\epsilon_{i}$
(8) $\displaystyle\epsilon_{i}$ $\displaystyle\sim N(0,\sigma^{2}=1)\qquad
i=1,...,200$ $\displaystyle f_{1}(d)$
$\displaystyle=\exp\left\\{\left(\frac{-d}{.5}\right)^{5}\right\\}\qquad$
$\displaystyle f_{2}(d)$
$\displaystyle=\nu\exp\left\\{\left(\frac{-d}{.5}\right)^{5}\right\\}\qquad\nu=(0,.25,.5,.75)$
$\displaystyle P(\zeta_{i}=1)$ $\displaystyle=P(\zeta_{i}=2)=.5;$

where $Z_{i}$ is a covariate generated as a fair Bernoulli random variable,
$\zeta_{i}$ is the subject specific cluster label indicating the true BEF
effect, $f(d)$, for the $i$th observation, and $\nu$ represents the varying
effect size at $d=0$.

The relative loss as a function of the effect size $\nu$ is shown in Figure 1.
As expected, one can see a decrease in relative loss and consequently improved
classification as the effect size increases.

Figure 1: Relative loss as a function of the difference in effect size:
$(1-\nu)$; see (8) for more details. Point estimates and error lines represent
median, 2.5 and 97.5 quantiles of loss across simulations, respectively.

###  3.3 Distance Distributions

As our method non-parametrically estimates cluster functions
$f_{i}(\mathcal{D}_{i})$ across continuously measured space using a basis
function expansion, correct estimation of the function requires there to be
BEFs observed at the relevant distances, $d:f(d)\neq 0$, within the study area
of interest. Of course these “relevant” distances are not known a priori and
so it is to the benefit of the investigator to err on the side of caution in
specifying a larger study area if possible. However, despite any preparatory
work that may be done to ensure an adequate area is included at the level of
the sample study, it is not clear how differing distributions of distances at
the latent cluster level may impact inference. For example, will suburbanites’
lower exposure to proximate FFRs impact the ability of the stapDP model to
discern the impact of FFRs on their health relative to their more exposure
rich urban counterparts? For this reason our second simulation study examines
how exposure to different distance distributions may impact classification.

We study this problem by considering three different generative distance
distributions which we label “Uniform”,“CA” and “Skew”. The first,
straightforwardly, refers to the idealistic - but unrealistic - scenario in
which there is equivalent information available at all distances within the
study area. The second two cases refer to more realistic situations in which
there are more likely to be a higher number of BEFs found further away from
the subject than close by – a consequence of area’s quadratic growth as a
function of distance. We create the first of these skewed distributions, “CA”,
by using maximum likelihood to fit a beta distribution to the distribution of
distances in our motivating California data distance distribution and the
second by altering a beta distribution to be a more extreme version of the
first. We generate distances under each distribution for each cluster in order
to examine how differing exposure patterns between clusters impact cluster
classification. The densities of each of these distributions can be found in
Figure S1 in the supplementary material.

Since the exposure information depends both on the distribution of distances
and the number of BEFs, we generate scenarios where the amount of information
increases as a function of the number of BEFs within the same distribution. We
simulate data under the same model as proposed in (8), with $\nu=0.25$, to
illustrate how a substantial, but not obvious, difference in cluster functions
manifest across the varying distance distribution settings. Fitting our STAP-
DP model under the priors and sampler settings previously described, we plot
the results below in Figure 2.

Figure 2: Relative loss as a function of different distance distributions.
Points and lines represent median, 2.5 and 97.5 quantiles of loss,
respectively. Row labels represent the distance distribution of the lower
effect size cluster and columns that of the higher effect size cluster.

Figure 2 shows a number of patterns worth highlighting. First, across all
distance distributions we observe a decrease in loss as the number of built
environment features increases. This is as expected – more information or
exposure results in a more easily detectable signal. Further, distance
distribution combinations that include more information result in lower levels
of relative loss compared to more skewed distance distributions. There are
number of cases where this can be seen in Figure 2, the most obvious being the
top-left diagonal panel where both clusters have Uniform distance
distributions; this has the lowest loss values across all panels due to the
relative abundance in exposure information. This pattern holds when comparing
to the more skewed distributions: The uniform-CA combination has higher error
than the uniform-uniform combination when there are a relatively small number
of built environment features present.

##  4 Fast food restaurants near schools and child obesity among public school
students in Los Angeles

There is a pressing need to understand contextual determinants of child
obesity, in order to implement population level strategies to reduce and
prevent it(McGuire, 2012). The food environment near schools has been proposed
and studied as a contextual factor that influences children’s diet, and thus
obesity (Currie et al., 2010; Davis and Carpenter, 2009; Sánchez et al., 2012;
Baek et al., 2016a). We use data on body weight of children attending public
schools in Los Angeles, CA, along with data on the locations of FFRs as a
marker of the food environment near schools, and apply our proposed method to
identify schools where children may be at higher risk of obesity, related to
food environment exposures. Identifying these schools may help prioritize or
tailor population-level interventions to address child obesity.

###  4.1 Data Description

Every year public schools in the State of California collect data on the
fitness status of pupils in 5th, 7th and 9th grade, as part of a state
mandate, using the Fitnessgram battery of tests (Institute, 2001). Child-level
data available for our analysis were collected during academic years 2001-2008
on 5th and 7th graders, and consists of children’s weight (Kg) and height (m)
and, and the following categorical covariates: sex, race-ethnicity, fitness
status (unfit, fit, fit above standard), and grade level. Weight and height
are transformed to body mass index (Kg/m2), and standardized to BMIz scores
according to age- and sex-specific growth curves published by the United
States’ Center for Disease Control. In contrast to adults, standardization is
needed when analyzing data from children of different ages, given children’s
rapid growth. To aid in managing the large database, and given that all child-
level covariates are categorical, the dataset is “collapsed” so that each row
represents a group of children within each school defined by the cross-
classification of categorical child-level characteristics described above. The
average BMIz of children in the group is the outcome of interest. Given the
categorical nature of the covariates and our use of weighting by the size of
the group represented by each row (below), this approach yields exactly the
same results as would be obtained if data had not been grouped, thus avoiding
biases in an ecological analysis (Schoenborn, 2002).

Data on school-level characteristics are also available from the CDE website
(see Table 1), and, importantly, so is the geocode of the school. School
geocodes were used for two purposes. First, the geocodes were used to link
schools to census tract level covariates. Second, the school geocodes were
used to calculate the distances between the school and the geocoded location
of each FFR in the LA area. FFRs were identified from the National
Establishment Time Series (NETS) database (Walls, 2013), using a published
algorithm that classifies specific food establishments as FFRs (Auchincloss et
al., 2012). Only FFRs within five miles of schools were kept for this
analysis. This distance was chosen to be a conservative as previous work
estimated that the distance at which FFRs cease to have an effect on childhood
obesity is approximately one mile (Baek et al., 2016a).

###  4.2 Los Angeles STAP-DP Model

We fit models estimating both the population-level and latent cluster-level
effects – STAP and STAP-DP models, respectively. Given the available data
consisting of subgroups of children defined by the cross-classification of
categorical covariates, we use the standardized average BMI within the
subgroup as the outcome. The models adjust for the student group and school-
level covariates listed in Table 1. Denoting these covariates as
$\bm{X}_{ijq}$ for student group $q=1,...,n_{ij}$, measured at year
$j=2001,...,2008$, attending school $i=1,...,N$, and using notation as
described in (6) our model for analyzing the Los Angeles data is:

$\displaystyle\text{BMIz}_{ijq}$
$\displaystyle=\bm{X}_{ijq}^{T}\bm{\gamma}+F_{i}(\mathcal{D}_{ij})+b_{i1}+b_{i2}\frac{\text{year}_{ij}}{10}+\epsilon_{ijq},$
(9) $\displaystyle\epsilon_{iqj}$ $\displaystyle\sim
N\left(0,\frac{\sigma^{2}}{n_{iqj}}\right),$ $\displaystyle\bm{b}_{i}$
$\displaystyle\sim MVN_{2}(\bm{0},\Sigma),$

where $n_{ijq}$ represents the number of students in student group $q$ during
year $j$ at school $i$. Given that FFRs may open or close during the study
period, $\mathcal{D}_{ij}$ represents the distances between school $i$ and
FFRs available within 5 miles during year $j$. Similar to our simulations, we
place a weakly informative Gamma(1,1) prior on each of the penalty parameters
in $\bm{\tau}$, associated with each cluster regression coefficients, the
residual precision $\sigma^{-2}$, and the concentration parameter $\alpha$.
The Gamma(1,1) prior on the concentration parameter is a common prior setting
in the DP literature, reflecting the a priori expectation that the
concentration parameter is $1$, so that fewer clusters are favored (Rodriguez
et al., 2008; Gelman et al., 2013). Additionally, we place a non-informative
Jeffrey’s prior on the covariance matrix for the school specific $\bm{b}_{i}$
vector: $p(\Sigma^{-1})\propto\mid\Sigma\mid^{\frac{3}{2}}$. Estimation is
conducted through rstapDP, drawing 2000 samples from each of 2 independent
MCMC chains after 8000 samples have been iterated as “burn-in” on each chain.
We check convergence via $\hat{R}$ diagnostic Vehtari et al. (2020) and
visually inspecting traceplots. We use $L=7$ coefficients in our spline basis
function expansion, and similarly use this basis to estimate the $f_{i}(d)$ on
a grid of values, calculating the 95% point-wise credible interval at each
distance grid point. We also calculate the posterior probability of co-
clustering which can be arranged in a matrix $\bm{P}\in\mathbb{R}^{N\times N}$
so that $\bm{P}_{i,i^{\prime}}=P($ school $i$ is co-clustered with school
$i^{\prime}$ across post burn-in iterations$)$. School cluster characteristics
are tabulated using the cluster mode school assignment calculated using (7) as
implemented in the rstapDP package via the assign_mode function.

For comparative purposes, we fit a model similar to (9) in all ways save for
restricting the $f_{i}(d)$ to be estimated at the population level - $f(d)$.
We fit this model using Hamiltonian Monte Carlo via the rsstap R package
(Peterson, 2020b), drawing 1000 samples after 1000 warm-up across 4
independent MCMC chains. Convergence is assessed via $\hat{R}$ diagnostic and
we calculate the analogous posterior estimate for $f(d)$ across the same grid
of distance values.

###  4.3 Los Angeles Results

Figure 3 shows four functions corresponding to the 3 estimated cluster
functions from the STAP-DP model, as well as the 1 homogeneous effect
estimate. We name the three cluster effects “Majority”, “High Risk” and “Low
Risk”. These names derive from the proportion of schools assigned to the
cluster as well as the relative effect size associated with the function at
and around distance 0 mi from the school: In the cluster labeled “High Risk”,
one additional FFR placed at distance 0 from a school is associated with an
expected 0.46 higher BMIz among students attending those schools (95% CI:
0.36, 0.58), all else equal. In contrast, placing one FFR at distance 0 from
the schools assigned to the “Low Risk” cluster is associated with lower BMIz
score, by -0.15 (95% CI: -0.17,-0.12), all else equal. The analogous values
for the “Majority” and homogeneous function estimates, are 0.01 (95% CI:
0.00,0.02) and 0.01 (95% CI: 0.01,0.0132), respectively. In all clusters, the
estimated associations rapidly decay with increasing distance, with all
association estimates effectively zero by 1 mile.

Figure 3: Changes in student BMI associated with FFR exposure across 5 mi.
Line and band represents median and 95% posterior credible interval. The
number following Each cluster label represents the median proportion of
schools that are assigned to that cluster. Dotted line (colored red online)
represents 0 reference line. Please note that the y-axis is different between
the two rows.

We now turn our attention to the matrix of co-clustering probabilities
$\bm{P}$ which we visualize using a heat-map in Figure S2 in the Supplementary
material, after applying Rodriguez et al. (2008)’s hierarchical sorting
algorithm to group schools with similar co-clustering probabilities together.
The estimate of co-clustering probabilities shows approximately 550 schools
(90%) are consistently co-clustered within one of the three clusters,
reflecting a high degree of model certainty in cluster configurations for
these schools. The remaining $\approx$ 50 schools show a greater uncertainty
between being classified in the “Majority” or “Low Risk” cluster. This
uncertainty likely stems from an insufficient number of FFRs present within
the relevant $\approx$ 1 mile distance from schools, where the cluster effects
are most discernibly different.

Further examination of the school characteristics associated with each cluster
details several suggestive, though not conclusive between-cluster differences
(Table 1). Differences between the three clusters of schools are fairly muted,
with summary statistics across student- and school-level measures describing
similar student populations and levels of median household income and
education amongst the neighborhoods of schools in each cluster. The most
noteworthy differences amongst the clusters are in the number of FFRs within
$\frac{1}{2}$ mile of the school – lower in the high risk group as compared to
the other two clusters – and the total enrollment – higher in the low and high
risk groups as compared to the majority cluster.

| Overall | Majority | Low Risk | High Risk
---|---|---|---|---
# Students | 752,529 | 655,017 | 62,573 | 34,939
% Obese | 52 | 52 | 51 | 46
Average BMIz | 0.84 | 0.85 | 0.78 | 0.71
% Female | 49 | 49 | 49 | 49
Race/Ethnicity |  |  |  |
% Asian | 3 | 2 | 2 | 7
% Black | 10 | 10 | 14 | 13
% Hispanic | 79 | 80 | 73 | 62
% White | 8 | 8 | 10 | 18
School Characteristic1 | N = 593 | N = 535 | N = 36 | N = 22
Total Enrollment (100’s of students) | 7.7 (5.0, 13.0) | 7.5 (4.9, 12.8) | 10.7 (6.3, 15.4) | 10.7 (6.2, 12.7)
# FFRs within 1/2 mile | 23 (8, 45) | 23 (8, 46) | 24 (8, 48) | 14 (5, 24)
# FFRs within 1 mile | 101 (69, 143) | 103 (68, 144) | 92 (70, 126) | 89 (74, 115)
% Free or Reduced Price Meals | 0.86 (0.69, 0.94) | 0.86 (0.69, 0.94) | 0.85 (0.75, 0.93) | 0.80 (0.58, 0.95)
Education2 | 14 (5, 27) | 14 (5, 28) | 13 (7, 25) | 13 (5, 23)
Income3 (1000 USD) | 34 (25, 49) | 34 (25, 49) | 33 (28, 52) | 32 (25, 55)
School Type |  |  |  |
Elementary | 468 (79%) | 427 (80%) | 26 (72%) | 15 (68%)
K-12 | 7 (1.2%) | 7 (1.3%) | 0 (0%) | 0 (0%)
Middle | 90 (15%) | 75 (14%) | 9 (25%) | 6 (27%)
High School | 12 (2.0%) | 11 (2.1%) | 1 (2.8%) | 0 (0%)
Other | 16 (2.7%) | 15 (2.8%) | 0 (0%) | 1 (4.5%)
Urbanicity |  |  |  |
Suburban | 75 (13%) | 66 (12%) | 6 (17%) | 3 (14%)
Urban | 518 (87%) | 469 (88%) | 30 (83%) | 19 (86%)

Table 1: Characteristics of children and schools in each cluster, assigned
using the mode cluster
1Statistics presented: Median (IQR); n (%)
2 Percent of individuals 25 years or older within the school’s census tract
with at least a bachelors degree.
3 Median Household Income among residents in the census tract where the school
is located.

The lack of substantial differences in measured characteristics between these
two clusters is noteworthy, suggesting that none of the observed
characteristics appear to modify the obesity risk associated with FFR
exposure. Although prior research shows that area-level education and income
modify the effects of child obesity interventions, for instance, the
protective effects of socioeconomic factors do not appear to extend to
children’s risk of obesity as due to proximate FFR exposure – at least within
this population. This lack of difference in socioeconomic characteristics
suggests that there are unmeasured variables that account for the
heterogeneous FFR effects. One potential unmeasured factor could be the type
of FFRs proximal to the schools, e.g. chains vs non-chain FFR’s. Another
possibility could be unmeasured student-level measures of wealth – which may
modify obesity risk – but we are not able to account for in our analysis.

##  5 Discussion

This work proposed a modeling approach to identify heterogeneity in distance-
dependent BEF effects. By allowing flexibility both across space and
identifying subgroups of subjects with different effects, this modeling
framework addresses two problems raised in the built environment literature,
namely the MAUP and the UGCP, respectively. The modeling approach was shown to
work well in both simulated data, as well as the data that motivated this
work, concerning children’s BMI and proximity to FFRs near their schools.
While spatial point pattern built environment data are the primary motivation
for this methodology, it could be also be applied to temporal or spatio-
temporal data, the latter which we discuss in greater detail below.

Similar to the HDLM proposed by Baek et al. (2016b), we seek to allow for
differences across subjects, or other substantively defined groups like
schools, in the BEF associations across space. In contrast to that work, we
pool subjects with similar association effects through the DPM, allowing us to
identify latent risk subject groups.

In simulations our model demonstrated classification robustness to differing
distributions of distances and expected improvement in classification due to
increased information through BEF exposure or effect size. Our analysis of
Fitnessgram data illustrated how one can analyze these data in terms of the
spatial effects estimated as well as the characteristics associated with each
latent cluster. The software to fit this model and perform the necessary
auxiliary functions is freely available through our R package rstapDP
(Peterson, 2020a).

There are a number of future directions with which to take this work. One
obvious direction would be to extend the modeling framework for more general
exponential family error distributions, though this makes estimation more
difficult, as the posterior distribution of $\bm{\beta}$ is no longer
available in closed form. Work by Ferrari (2020) has used a Riemann
Hamiltonian Monte Carlo sampler in this context for models similar to ours,
without smooth functional terms. This could provide one avenue to pursue.
Another direction to explore would be to incorporate temporally indexed BEF
data to enable spatio-temporal function estimation via tensor product of the
spline basis function expansion used here. This approach would allow for
cluster estimates across space and time, increasing the dimensionality and
consequently, relevancy, of this work to more precisely target and understand
how environments shape health and health behaviors across both time and space.

Finally, while there have been numerous methods to identify associations
between subjects and BEF exposure we believe this to be the first to utilize
techniques in both the Bayesian and functional non-parametric literature to
identify heterogeneous BEF effects across a population.

## Acknowledgements

This research was partially supported by NIH grants R01-HL131610 (PI: Sánchez)
and R01-HL136718 (MPIs: Sanchez-Vaznaugh and Sánchez).

## References

  * Baek et al. [2016a] Jonggyu Baek, Brisa N Sánchez, Veronica J Berrocal, and Emma V Sanchez-Vaznaugh.  Distributed lag models: examining associations between the built environment and health.  _Epidemiology (Cambridge, Mass.)_ , 27(1):116, 2016a.
  * Baek et al. [2017] Jonggyu Baek, Jana A Hirsch, Kari Moore, Loni Philip Tabb, Tonatiuh Barrientos-Gutierrez, Lynda D Lisabeth, Ana V Diez-Roux, and Brisa N Sánchez.  Methods to study variation in the associations between food store availability and body mass in the multi-ethnic study of atherosclerosis.  _Epidemiology (Cambridge, Mass.)_ , 28(3):403, 2017.
  * Kaufman et al. [2019] Tanya K Kaufman, Andrew Rundle, Kathryn M Neckerman, Daniel M Sheehan, Gina S Lovasi, and Jana A Hirsch.  Neighborhood recreation facilities and facility membership are jointly associated with objectively measured physical activity.  _Journal of urban health_ , 96(4):570–582, 2019\.
  * Kern et al. [2017] David M Kern, Amy H Auchincloss, Mark F Stehr, Ana V Diez Roux, Latetia V Moore, Genevieve P Kanter, and Lucy F Robinson.  Neighborhood prices of healthier and unhealthier foods and associations with diet quality: Evidence from the multi-ethnic study of atherosclerosis.  _International journal of environmental research and public health_ , 14(11):1394, 2017.
  * Fotheringham and Wong [1991] A Stewart Fotheringham and David WS Wong.  The modifiable areal unit problem in multivariate statistical analysis.  _Environment and planning A_ , 23(7):1025–1044, 1991.
  * Spielman and Yoo [2009] Seth E Spielman and Eun-hye Yoo.  The spatial dimensions of neighborhood effects.  _Social science & medicine_, 68(6):1098–1105, 2009.
  * Wong [2009] David Wong.  The modifiable areal unit problem (MAUP).  _The SAGE handbook of spatial analysis_ , 105(23):2, 2009.
  * Guo and Bhat [2004] Jessica Y Guo and Chandra R Bhat.  Modifiable areal units: problem or perception in modeling of residential location choice?  _Transportation Research Record_ , 1898(1):138–147, 2004.
  * Ji et al. [2009] Chunlin Ji, Daniel Merl, Thomas B Kepler, and Mike West.  Spatial mixture modelling for unobserved point processes: Examples in immunofluorescence histology.  _Bayesian analysis (Online)_ , 4(2):297, 2009\.
  * James et al. [2014] Peter James, David Berrigan, Jaime E Hart, J Aaron Hipp, Christine M Hoehner, Jacqueline Kerr, Jacqueline M Major, Masayoshi Oka, and Francine Laden.  Effects of buffer size and shape on associations between the built environment and energy balance.  _Health & place_, 27:162–170, 2014.
  * Macintyre et al. [2002] Sally Macintyre, Anne Ellaway, and Steven Cummins.  Place effects on health: how can we conceptualise, operationalise and measure them?  _Social Science & Medicine_, 55(1):125–139, 2002.  ISSN 0277-9536.  doi:https://doi.org/10.1016/S0277-9536(01)00214-3.  URL https://www.sciencedirect.com/science/article/pii/S0277953601002143.  Selected papers from the 9th International Symposium on Medical G eography.
  * Kwan [2013] Mei-Po Kwan.  Beyond space (as we knew it): Toward temporally integrated geographies of segregation, health, and accessibility: Space–time integration in geography and giscience.  _Annals of the Association of American Geographers_ , 103(5):1078–1086, 2013.
  * Kwan [2018] Mei-Po Kwan.  The limits of the neighborhood effect: Contextual uncertainties in geographic, environmental health, and social science research.  _Annals of the American Association of Geographers_ , 108(6):1482–1490, 2018.
  * Peterson and Sanchez [2018] Adam Peterson and Brisa Sanchez.  rstap: An r package for spatial temporal aggregated predictor models.  _arXiv preprint arXiv:1812.10208_ , 2018.
  * Peterson et al. [2021] Adam Peterson, Jana Hirsch, and Brisa Sanchez.  Spatial temporal aggregated predictors to examine built environment health effects.  2021\.
  * Baek et al. [2016b] Jonggyu Baek, Emma V Sanchez-Vaznaugh, and Brisa N Sánchez.  Hierarchical distributed-lag models: exploring varying geographic scale and magnitude in associations between the built environment and health.  _American journal of epidemiology_ , 183(6):583–592, 2016b.
  * Walls [2013] D Walls.  National establishment time-series (NETS) database: 2012 database description.  _Oakland: Walls & Associates_, 2013.
  * Rodriguez et al. [2014] Abel Rodriguez, David B Dunson, et al.  Functional clustering in nested designs: modeling variability in reproductive epidemiology studies.  _Annals of Applied Statistics_ , 8(3):1416–1442, 2014.
  * Ray and Mallick [2006] Shubhankar Ray and Bani Mallick.  Functional clustering by bayesian wavelet methods.  _Journal of the Royal Statistical Society: Series B (Statistical Methodology)_ , 68(2):305–332, 2006.
  * O’Sullivan [1986] Finbarr O’Sullivan.  A statistical perspective on ill-posed inverse problems.  _Statistical science_ , pages 502–518, 1986.
  * Wahba [1990] Grace Wahba.  _Spline models for observational data_.  SIAM, 1990.
  * Wood [2017] S.N Wood.  _Generalized Additive Models: An Introduction with R_.  Chapman and Hall/CRC, 2 edition, 2017.
  * Eilers and Marx [1996] Paul HC Eilers and Brian D Marx.  Flexible smoothing with b-splines and penalties.  _Statistical science_ , pages 89–102, 1996.
  * Wood [2004] Simon N Wood.  Stable and efficient multiple smoothing parameter estimation for generalized additive models.  _Journal of the American Statistical Association_ , 99(467):673–686, 2004.
  * Wood [2016] Simon N Wood.  Just another gibbs additive modeller: interfacing jags and mgcv.  _arXiv preprint arXiv:1602.02539_ , 2016.
  * Ferguson [1973] Thomas S Ferguson.  A Bayesian analysis of some nonparametric problems.  _The annals of statistics_ , pages 209–230, 1973.
  * Gelman et al. [2013] Andrew Gelman, John B Carlin, Hal S Stern, David B Dunson, Aki Vehtari, and Donald B Rubin.  _Bayesian data analysis_.  Chapman and Hall/CRC, 2013.
  * Fitzmaurice et al. [2008] Garrett Fitzmaurice, Marie Davidian, Geert Verbeke, and Geert Molenberghs.  _Longitudinal data analysis_.  CRC press, 2008.
  * Peterson [2020a] Adam Peterson.  _rstapDP: Functional Dirichlet Process Spatial Temporal Aggregated Predictors in R_ , 2020a.  version 0.1.1.
  * [30] Adam Peterson and Brisa Sánchez.  Supplementary material for “heterogeneous effects in the built environment”.
  * R Core Team [2019] R Core Team.  _R: A Language and Environment for Statistical Computing_.  R Foundation for Statistical Computing, Vienna, Austria, 2019.  URL https://www.R-project.org/.
  * Binder [1978] David A Binder.  Bayesian cluster analysis.  _Biometrika_ , 65(1):31–38, 1978.
  * Lau and Green [2007] John W Lau and Peter J Green.  Bayesian model-based clustering procedures.  _Journal of Computational and Graphical Statistics_ , 16(3):526–558, 2007.
  * Wade et al. [2018] Sara Wade, Zoubin Ghahramani, et al.  Bayesian cluster analysis: Point estimation and credible balls (with discussion).  _Bayesian Analysis_ , 13(2):559–626, 2018.
  * Rodriguez et al. [2008] Abel Rodriguez, David B Dunson, and Alan E Gelfand.  The nested Dirichlet process.  _Journal of the American Statistical Association_ , 103(483):1131–1154, 2008.
  * McGuire [2012] Shelley McGuire.  Institute of Medicine (IOM) Early Childhood Obesity Prevention Policies. Washington, DC: The National Academies Press; 2011.  _Advances in Nutrition_ , 3(1):56–57, 01 2012\.  ISSN 2161-8313.  doi:10.3945/an.111.001347.  URL https://doi.org/10.3945/an.111.001347.
  * Currie et al. [2010] Janet Currie, Stefano DellaVigna, Enrico Moretti, and Vikram Pathania.  The effect of fast food restaurants on obesity and weight gain.  _American Economic Journal: Economic Policy_ , 2(3):32–63, 2010.
  * Davis and Carpenter [2009] Brennan Davis and Christopher Carpenter.  Proximity of fast-food restaurants to schools and adolescent obesity.  _American Journal of Public Health_ , 99(3):505–510, 2009.
  * Sánchez et al. [2012] Brisa N Sánchez, Emma V Sanchez-Vaznaugh, Ali Uscilka, Jonggyu Baek, and Lindy Zhang.  Differential associations between the food environment near schools and childhood overweight across race/ethnicity, gender, and grade.  _American journal of epidemiology_ , 175(12):1284–1293, 2012.
  * Institute [2001] The Cooper Institute.  The fitnessgram assessment, 2001.  URL https://fitnessgram.net/assessment/.  accessed 2020.
  * Schoenborn [2002] Charlotte A Schoenborn.  _Body weight status of adults: United States, 1997-98_.  Number 330. US Department of Health and Human Services, Centers for Disease Control and …, 2002.
  * Auchincloss et al. [2012] Amy H Auchincloss, Kari AB Moore, Latetia V Moore, and Ana V Diez Roux.  Improving retrospective characterization of the food environment for a large region in the United States during a historic time period.  _Health & place_, 18(6):1341–1347, 2012.
  * Vehtari et al. [2020] Aki Vehtari, Andrew Gelman, Daniel Simpson, Bob Carpenter, Paul-Christian Bürkner, et al.  Rank-normalization, folding, and localization: An improved $\hat{R}$ for assessing convergence of mcmc.  _Bayesian Analysis_ , 2020.
  * Peterson [2020b] Adam Peterson.  _rsstap: Spline Spatial Temporal Aggregated Predictors in R_ , 2020b.  R package version 0.1.0.
  * Ferrari [2020] Diogo Ferrari.  Modeling context-dependent latent effect heterogeneity.  _Political Analysis_ , 28(1):20–46, 2020.


# Causal Video Summarizer for Video Exploration

###### Abstract

Recently, video summarization has been proposed as a method to help video
exploration. However, traditional video summarization models only generate a
fixed video summary which is usually independent of user-specific needs and
hence limits the effectiveness of video exploration. Multi-modal video
summarization is one of the approaches utilized to address this issue. Multi-
modal video summarization has a video input and a text-based query input.
Hence, effective modeling of the interaction between a video input and text-
based query is essential to multi-modal video summarization. In this work, a
new causality-based method named Causal Video Summarizer (CVS) is proposed to
effectively capture the interactive information between the video and query to
tackle the task of multi-modal video summarization. The proposed method
consists of a probabilistic encoder and a probabilistic decoder. Based on the
evaluation of the existing multi-modal video summarization dataset,
experimental results show that the proposed approach is effective with the
increase of $+5.4$% in accuracy and $+4.92$% increase of $F1$-score, compared
with the state-of-the-art method.

##  1 Introduction

Video contents are growing at an ever-increasing speed and beyond the capacity
of an individual for full comprehension. According to [1], more than 18,000
hours of videos are uploaded to YouTube per minute. Exploring this quantity of
video data is a daunting task, and video summarization is gaining traction as
the ideal solution [2, 3]. The main idea of video summarization is to
automatically generate a short video clip that summarizes the content of an
original, longer video by capturing its important parts [4, 5, 6]. However,
traditional video summarization approaches, e.g., [7, 8], only create a fixed
video summary for a given video. Given that viewers may have different
preferences and videos have different ideal summarization, this issue reduces
the effectiveness of video exploration.

Fig. 1: Multi-modal video summarization. The input video is summarised taking
into account text-based queries. “Input query-1: Sport of running” and “Input
query-2: Sport of cycling’ independently drive the model to generate video
summaries that contain running-related content and cycling-related content,
respectively.

Multi-modal video summarization has been proposed as an approach to improve
the effectiveness of video exploration [9, 10]. It generates video summaries
for a given video based on the text-based query provided by the user,
illustrated in Figure 1. Conventional video summarization only has video input
modality, while an efficient choice for multi-modal video summarization is a
text-based query, in addition to video [9, 10]. When multi-modal video
summarization is used to help video exploration, effectively modeling the
implicit relation/interaction between the text-based query and the video is
important [9]. In [10], the proposed model exploits a joint representation of
vision and language to perform multi-modal video summarization. However, the
implicit interaction between the query and the video is not properly modeled
because the simple average-based method used in [10] probably is not effective
enough.

As stated in [11, 12, 13], causal effect modeling, visualized in Figure 2, is
helpful for a machine learning task, e.g., image/frame classification, and
affects model performance in a positive way. In this work, a new causal video
summarizer (CVS) is proposed that tackles the aforementioned issue to improve
the performance of a multi-modal video summarization model. The proposed
causality-based model consists of a multi-modal feature processing module
(MFPM), probabilistic encoding module (PEM), and a probabilistic decoding
module (PDM), referring to Figure 3. Studying multi-modal video summarization
from the causality perspective [11] eliminates the need for a priori
definition of objectives [9] based on high-level concepts. For the implicit
interaction between the video and query, an attention mechanism is applied to
better capture the interactive information.

According to [14], typically the generated video summary by a video
summarization algorithm is composed of a set of representative video frames or
video fragments. Frame-based video summaries are not restricted by timing or
synchronization issues and, therefore, they provide more flexibility in terms
of data organization for video exploration purpose [15, 14]. In this work, the
proposed CVS is validated on the frame-based multi-modal video summarization
dataset [10]. Experimental results show that the proposed method is effective
and significantly increases both the accuracy and $F1$-score, compared with
the state-of-the-art method.

##  2 Related Work

2.1 Video Summarization with A Single Modality

Recently, a number of methods with different techniques have been proposed for
video summarization, such as [4, 16, 17, 1, 3]. In [4, 16], human-explainable
concepts, e.g., interestingness and representativeness, are used to
characterize a good video summary. Then, objectives are built based on those
concepts. Visual and textual information related to a video can then be
captured by the targeted objectives. However, a notable limitation is that it
is hard to define/measure all possible confounders [11] a priori, such as
compactness [17], uniformity, and representativeness [4]. In [1, 17], an
attention mechanism, a generative adversarial network, and reinforcement
learning are used to perform video summarization. Since video summarization
with a single modality cannot effectively help video exploration [9, 10], in
this work a new multi-modal video summarization model is proposed to address
this problem.

Fig. 2: Explanation of the causal effect modeling concept, i.e., causal
graphical model [11], in multi-modal video summarization. t is a treatment,
e.g., visual or textual perturbation [12]. y is an outcome, e.g., an
importance score of a video frame or a relevance score between the input text-
based query and video. Z is an unobserved confounder, e.g., representativeness
or interestingness [4, 16, 9]. X is a noisy view [11] on the hidden confounder
Z, say the input text query and video.

(a) Acting as an encoder.

(b) Acting as a decoder.

Fig. 3: This figure shows the proposed end-to-end causal attention model for
multi-modal video summarization, the practical implementation of the causal
graphical model in Figure 2. In (a), the input video and input text query are
processed by MFPM for generating the feature map x which is the input of PEM.
Then, PEM generates the corresponding probabilistic encoding
$q(\textbf{z}|\textbf{x},y,t)$ of x with the operations of the introduced
causal attention. In (b), PDM takes $p(\textbf{z})$, estimated by
$q(\textbf{z}|\textbf{x},y,t)$, to approximately reconstruct x based on
$p(\textbf{x}|\textbf{z})$ and generate the prediction score labels $y$ based
on $p(y|\textbf{z},t)$ at the same time. Please refer to Methodology section
for details.

2.2 Multi-modal Video Summarization

Instead of only considering the visual input, several works have investigated
the potential of using additional modalities, e.g., viewers’ comments, video
captions, or any other available contextual data, to help models’ performance.
[18, 19, 10, 20, 21, 9, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34,
35, 36]. In [25], a multi-modal video summarization method is introduced for
key-frame extraction from first-person videos. In [18], a multi-modal deep-
learning-based model is proposed to summarize videos of soccer games. In [19],
a semantic-based video fragment selection and a visual-to-text mapping are
applied based on the relevance between the original and the automatically-
generated video descriptions, with the help of semantic attended networks.
Existing multi-modal video summarization models do not focus on causal effect
modeling [11]. Hence, the interactive information between the video and query
is difficult to capture effectively. In this work, a new causality-based
multi-modal video summarization method is introduced, based on a probabilistic
encoder-decoder framework and attention mechanism, to better capture the
interaction between the video and query and make the video exploration more
effective.

2.3 Causality and Variational Autoencoders

Proxy variables, e.g., input data X in Figure 2, and the challenges about how
to use them correctly have been considered in the literature of causal
inference [37]. In many observational studies [4, 16, 38], understanding the
best way to derive and measure possible proxy variables is crucial. Building
on the previous work [39, 40], the authors of [11] and recent works [41, 38]
for vision applications have tried to exploit proxy variables to study
conditions for causal identifiability. In many cases, the general idea is that
one should first attempt to infer the joint distribution
$p(\textbf{X},\textbf{Z})$ between the hidden confounders and the proxy, and
then use that knowledge to adjust for the hidden confounders [42, 43, 44].
Take the causal graphical model in Figure 2 as an example. The authors of [43]
have shown that if X and Z are categorical, with X having at least as many
categories as Z, and with the matrix $p(\textbf{X},\textbf{Z})$ being full-
rank, one could use a matrix inversion formula, an approach called “effect
restoration” [44], to identify the causal effect of t on y, referring to
Methodology section for details. Recently, the authors of [42] have given the
conditions under which one could identify more complicated and general proxy
models. The proposed causality-based method for multi-modal video
summarization is mainly inspired by [11, 45, 43, 37].

##  3 Methodology

3.1 Causal Video Summarizer (CVS)

In this section, details of the proposed CVS are described. Note that causal
effect modeling for real-world multi-modal video summarization is very
complicated [11]. Hence, in this work, the assumptions mentioned in [11] are
followed to model the problem of multi-modal video summarization based on the
concept of causal effect inference, i.e., causal graphical model. The proposed
CVS is mainly composed of a multi-modal feature processing module (MFPM), a
probabilistic encoding module (PEM), and a probabilistic decoding module
(PDM), referring to Figure 3. The generation of a good video summary is
affected by latent factors. In this work, the latent factors are considered as
the causal effect. Specifically, the concept of causal effect inference [11],
i.e., a causal graphical model, is used to model the multi-modal video
summarization problem, referring to Figure 2. Figure 2 contains the four key
components of the causal graphical model: “input data (X)”, “latent confounder
(Z)”, “treatment (t)”, and “output result (y)”. From the modeling perspective
of multi-modal video summarization, X is “an input video with a text query”. t
is “a visual or textual treatment”. Note that treatment in causal effect
modeling is a way to make an input characteristic more salient and help a
model learn better [11]. y is “a relevance score between the input text-based
query and video frame or an importance score of a video frame”. Z is “a
variational latent representation” which the proposed model aims to learn from
X for reconstruction [45, 11], referring to Figure 3 for the practical
implementation. The proposed causal model simultaneously generates outcome
score labels y when X is reconstructed, illustrated in Figure 3. Thereafter,
video summaries can be created based on the generated outcome score labels.

Objective function for training. According to [45, 11, 46], variational
autoencoder (VAE) aims to learn a variational latent representation Z from
data X for reconstruction; it is capable of learning the latent variables when
used in a causal graphical model. Hence, the proposed CVS is built on top of
VAE. Now, it is time to form a single objective for the encoder, in Figure
3-(a), and the decoder, in Figure 3-(b), to learn meaningful causal
representations in the latent space and generate video summaries. Based on
Figure 2, we know that the true posterior over Z depends on not just X, but
also on t and y. We have to know the treatment assignment t along with its
outcome y before inferring the distribution over Z. Hence, the following two
auxiliary distributions are introduced for the treatment assignment t and the
outcome y, referring to Equations (1), and (2).

$q(t_{i}|\textbf{x}_{i})=Bern(\sigma(g_{\phi_{5}}(\textbf{x}_{i})))$ (1)
$q(y_{i}|\textbf{x}_{i},t_{i})=\sigma(t_{i}g_{\phi_{6}}(\textbf{x}_{i})+(1-t_{i})g_{\phi_{7}}(\textbf{x}_{i})),$
(2)

where $g_{\phi_{k}}(\cdot)$ is a neural network with variational parameters
$\phi_{k}$ for $k=5,6,7$. Note that, unlike traditional VAEs simply passing
the feature map directly to the latent space, i.e., the upper path of PEM in
Figure 3-(a)), the feature map is also sent to the other two paths, i.e., the
middle and the lower paths in PEM in Figure 3-(a)), for the posterior
estimations of the treatment $t_{i}$ and the outcome $y_{i}$.

The introduced auxiliary distributions help us predict $t_{i}$ and $y_{i}$ for
new samples. To estimate the parameters of these two distributions,
$q(t_{i}|\textbf{x}_{i})$ and $q(y_{i}|\textbf{x}_{i},t_{i})$, we add an
auxiliary objective, referring to Equation (3), to the model training
objective over $N$ data samples.

$\mathcal{L}_{auxiliary}=$ $\sum^{N}_{i=1}(\log
q(t_{i}=t_{i}^{*}|\textbf{x}_{i}^{*})+\log
q(y_{i}=y_{i}^{*}|\textbf{x}_{i}^{*},t_{i}^{*})),$ (3)

where $\textbf{x}_{i}^{*}$, $t_{i}^{*}$ and $y_{i}^{*}$ are observed values in
the training set.

Finally, we have the following overall training objective for the encoder and
decoder networks. See Equation (4).

$\mathcal{L}_{causal}=\mathcal{L}_{auxiliary}~{}+$
$\sum^{N}_{i=1}\mathbb{E}_{q(\textbf{z}_{i}|\textbf{x}_{i},t_{i},y_{i})}[\log
p(\textbf{x}_{i},t_{i}|\textbf{z}_{i})+\log p(y_{i}|t_{i},\textbf{z}_{i})~{}+$
$\log p(\textbf{z}_{i})-\log q(\textbf{z}_{i}|\textbf{x}_{i},t_{i},y_{i})].$
(4)

See [11] for the detailed derivation and explanation of Equations (1) to (4),
$\log p(\textbf{x}_{i},t_{i}|\textbf{z}_{i})$, $\log
p(y_{i}|t_{i},\textbf{z}_{i})$, and $\log p(\textbf{z}_{i})$.

3.2 Spatial and Channel-wise Attentions

A dual attention network, based on self-attention, has been proposed to
adaptively integrate local features with their global dependencies [47, 48].
The dual attention network consists of spatial and channel-wise attention
modules. These two modules are capable of modeling the interdependencies in
spatial and channel dimensions. To better capture the implicit interactive
information, a causal attention mechanism, based on the dual attention
network, is introduced to reinforce the probabilistic encoder in Figure 3-(a).
In the PEM, the channel-wise attention selectively emphasizes interdependent
channel maps by integrating associated features among all channel maps. The
spatial attention selectively aggregates the feature at each position by a
weighted sum of the features at all positions.

3.3 Video Summary Generation

In Figure 3-(a), during the inference phase, the input video and input text-
based query go through the MFPM for generating feature maps x, and the PEM for
generating the probabilistic encoding of the feature maps
$q(\textbf{z}|\textbf{x},y,t)$. Then, the generated probabilistic encoding of
the feature maps is sent to the probabilistic decoder, referring to Figure
3-(b), to reconstruct x, based on $p(\textbf{x}|\textbf{z})$, and generate
prediction score labels $y$, based on $p(y|\textbf{z},t)$, for the input video
and query simultaneously. Finally, based on these generated score labels, a
set of video frames is selected from the original input video to create the
final video summary. Note that the video summary budget is considered as a
user-defined hyperparameter in video exploration [10].

##  4 Experiments

Fig. 4: This figure shows the randomly selected qualitative results of the
proposed causal model and the baseline model. According to this figure, we see
that the proposed model (in green) prediction result frame pattern is more
similar to the ground truth (in red) pattern than the baseline (in blue). Note
that in each frame pattern, orange denotes “not selected frames” and gray
denotes “selected frames”. Each video here has 199 frames and we also show the
corresponding indices of the selected frames in this figure.  Table 1:
Ablation study of the introduced attention mechanism and comparison with the
state-of-the-art multi-modal video summarization model “QueryVS” [10]. The
results are based on the introduced CVSD and the metric of accuracy from [10].
The attention mechanism is effective and the proposed CVS outperforms the
“QueryVS” model in the metric of accuracy under six different situations of
inputs. Please refer to the Effectiveness Analysis section for the definitions
of notations in the table.

Model | $V_{s\&p}$ | $V_{blur}$ | $T+V_{s\&p}$ | $T+V_{blur}$ | $T_{k}+V_{s\&p}$ | $T_{k}+V_{blur}$
---|---|---|---|---|---|---
QueryVS [10] | 0.629 | 0.507 | 0.547 | 0.533 | 0.501 | 0.488
mygray CVS with attention | 0.687 | 0.596 | 0.554 | 0.535 | 0.532 | 0.524
CVS without attention | 0.676 | 0.423 | 0.442 | 0.530 | 0.514 | 0.500

Table 2: Comparison with the state-of-the-art model “QueryVS” [10] validated on [10]’s dataset. The results are based on the metric of accuracy [10] and $F_{1}$-score [49]. Proposed CVS outperforms the model in [10] by $+5.4$% in accuracy and $+4.92$% in $F_{1}$-score. Evaluation Metric | Accuracy |  $F_{1}$-score
---|---|---
QueryVS [10] | 0.704 | 0.508
mygray CVS | 0.742 | 0.533

4.1 Dataset and Evaluation Metrics

In the experiments, [10]’s multi-modal video summarization dataset and the
introduced causal video summarization dataset (CVSD) with treatment labels are
used to verify the proposed method. The introduced CVSD is constructed based
on [10]’s dataset. [10]’s dataset is composed of $190$ videos with a duration
of two to three minutes for each video. Each video in their dataset is
retrieved based on a given text-based query. The entire dataset is divided
into splits of 60%/20%/20% for training/validation/testing, respectively.
Annotations from human experts are necessary for the automatic evaluation of
multi-modal video summarization. Hence, the authors of [10] sample all of the
$190$ videos at one frame per second (fps), and Amazon Mechanical Turk (AMT)
is used to annotate every frame with its relevance score with respect to the
given text-based query. A single ground truth relevance score is created for
each query-video pair by merging the corresponding relevance annotations from
AMT workers. In [10]’s dataset, the maximum number of words of a text query is
$8$ and the maximum number of frames of video is $199$.

For the evaluation metric, the authors of [10] propose to use accuracy as the
evaluation metric. That is, a predicted relevance score is considered correct
when the predicted score is the same as the majority of human experts’
provided score. In addition, motivated by [49, 3], in this work $F_{1}$-score
is also used to quantify the performance of the proposed model by measuring
the agreement between the gold standard and predicted scores provided by the
human experts.

4.2 Experimental Setup

In this work, [10]’s dataset and the introduced CVSD have the same
preprocessing as follows. Since the video lengths in the CVSD and [10]’s
dataset are varied, the number of frames for each video is different based on
fps $=1$. The maximum number of frames of a video is $199$ in the used
dataset. Hence, the frame-repeating video preprocessing technique [10] is
adopted to make all the videos have the same length of $199$. ResNet [50] pre-
trained on ImageNet [51] is used to extract $199$ frame-based features for
each video, and the feature used is located in the visual layer one layer
below the classification layer. The input frame size of the ResNet is $224$ by
$224$ with red, green, and blue channels. Each image channel is normalized by
mean $=(0.4280,0.4106,0.3589)$ and standard deviation
$=(0.2737,0.2631,0.2601)$. PyTorch is used for the implementation and to train
models with $50$ epochs, $0.01$ learning rate, and Adam optimizer [52]. For
the hyperparameters of the Adam, $\beta_{1}=0.9$ and $\beta_{2}=0.999$ are the
coefficients used for computing moving averages of gradient and its square.
The term added to the denominator to improve numerical stability is
$\epsilon=10^{-8}$.

4.3 Effectiveness Analysis

Effectiveness of Causal Attention Mechanism. In the experiments, the CVSD is
used to verify the effectiveness of the introduced attention mechanism under
the following six different situations of inputs: “visual input with salt and
pepper treatment ($V_{s\&p}$)”, “visual input with blurring treatment
($V_{blur}$)”, “clean text-based query input and visual input with salt and
pepper treatment ($T+V_{s\&p}$)”, “clean text-based query input and visual
input with blurring treatment ($T+V_{blur}$)”, “disturbed text-based query
input and visual input with salt and pepper treatment ($T_{k}+V_{s\&p}$)”, and
“disturbed text-based query input and visual input with blurring treatment
($T_{k}+V_{blur}$)”. Note that $T_{k}$ denotes randomly removing $k$ words,
e.g., $k=2$, from a text-based query input. According to Table 4, the results
show that the introduced attention mechanism results in performance
improvement in the causal model. The main reason is that the introduced
attention mechanism helps to better capture the implicit interactive
information between the visual and textual inputs.

Effectiveness of Causal Effect Modeling. The main difference between the
proposed causal model and the existing state-of-the-art multi-modal video
summarization method, e.g., [10], is the causal effect modeling. According to
[11], proper causal effect modeling improves the performance of a machine
learning model. In this work, we claim that if the causal effect modeling is
effective, the proposed model’s performance can be improved. Based on Table 4
and Table 2, the results show that the proposed causal model for multi-modal
video summarization beats the state-of-the-art [10], with the increase of
$+5.4$% in accuracy and $+4.92$% increase of $F1$-score. The main reason is
the proposed model is a causality-based model which properly considers the
causal effect modeling. Hence, the claim is well proved by the experimental
results. Qualitative results are presented in Figure 4.

##  5 Conclusion

In this work, a new perspective is introduced to build an end-to-end deep
causal model for multi-modal video summarization. The proposed Causal Video
Summarizer is based on a probabilistic encoder-decoder architecture.
Experimental results show that the proposed causal model is effective and
achieves state-of-the-art performance, $+5.4$% in accuracy and $+4.92$%
increase of $F1$-score.

##  6 Acknowledgments

This project has received funding from the European Union’s Horizon 2020
research and innovation programme under the Marie Skłodowska-Curie grant
agreement No 765140.

## References

  * [1] Zhong Ji, Kailin Xiong, Yanwei Pang, and Xuelong Li,  “Video summarization with attention-based encoder-decoder networks,”  IEEE Transactions on Circuits and Systems for Video Technology, 2019\.
  * [2] Danila Potapov, Matthijs Douze, Zaid Harchaoui, and Cordelia Schmid,  “Category-specific video summarization,”  in ECCV, 2014, pp. 540–555.
  * [3] Yale Song, Jordi Vallmitjana, Amanda Stent, and Alejandro Jaimes,  “Tvsum: Summarizing web videos using titles,”  in Conference on Computer Vision and Pattern Recognition, 2015, pp. 5179–5187.
  * [4] Michael Gygli, Helmut Grabner, and Luc Van Gool,  “Video summarization by learning submodular mixtures of objectives,”  in Conference on Computer Vision and Pattern Recognition, 2015, pp. 3090–3098.
  * [5] Ke Zhang, Wei-Lun Chao, Fei Sha, and Kristen Grauman,  “Summary transfer: Exemplar-based subset selection for video summarization,”  in Proceedings of the IEEE conference on computer vision and pattern recognition, 2016, pp. 1059–1067.
  * [6] Yujia Zhang, Michael Kampffmeyer, Xiaoguang Zhao, and Min Tan,  “Dtr-gan: Dilated temporal relational adversarial network for video summarization,”  in Proceedings of the ACM Turing Celebration Conference-China, 2019, pp. 1–6.
  * [7] Sandra Eliza Fontes De Avila, Ana Paula Brandão Lopes, Antonio da Luz Jr, and Arnaldo de Albuquerque Araújo,  “Vsumm: A mechanism designed to produce static video summaries and a novel evaluation method,”  Pattern Recognition Letters, vol. 32, no. 1, pp. 56–68, 2011.
  * [8] Aidean Sharghi, Ali Borji, Chengtao Li, Tianbao Yang, and Boqing Gong,  “Improving sequential determinantal point processes for supervised video summarization,”  in European Conference on Computer Vision, 2018, pp. 517–533.
  * [9] Arun Balajee Vasudevan, Michael Gygli, Anna Volokitin, and Luc Van Gool,  “Query-adaptive video summarization via quality-aware relevance estimation,”  in Proceedings ACM MM, 2017, pp. 582–590.
  * [10] Jia-Hong Huang and Marcel Worring,  “Query-controllable video summarization,”  in ICMR, 2020, pp. 242–250.
  * [11] Christos Louizos, Uri Shalit, Joris M Mooij, David Sontag, Richard Zemel, and Max Welling,  “Causal effect inference with deep latent-variable models,”  in NIPS, 2017, pp. 6446–6456.
  * [12] Chao-Han Huck Yang, I Hung, Te Danny, Yi Ouyang, and Pin-Yu Chen,  “Training a resilient q-network against observational interference,”  AAAI, 2022.
  * [13] Chao-Han Huck Yang, I Hung, Te Danny, Yi-Chieh Liu, and Pin-Yu Chen,  “Treatment learning transformer for noisy image classification,”  arXiv preprint arXiv:2203.15529, 2022.
  * [14] Evlampios Apostolidis, Eleni Adamantidou, Alexandros I Metsai, Vasileios Mezaris, and Ioannis Patras,  “Video summarization using deep neural networks: A survey,”  arXiv preprint arXiv:2101.06072, 2021.
  * [15] Janko Calic, David P Gibson, and Neill W Campbell,  “Efficient layout of comic-like video summaries,”  IEEE Transactions on Circuits and Systems for Video Technology, vol. 17, no. 7, pp. 931–936, 2007.
  * [16] Bryan A Plummer, Matthew Brown, and Svetlana Lazebnik,  “Enhancing video summarization via vision-language embedding,”  in Conference on Computer Vision and Pattern Recognition, 2017, pp. 5781–5789.
  * [17] Yujia Zhang, Michael Kampffmeyer, Xiaodan Liang, Dingwen Zhang, Min Tan, and Eric P Xing,  “Dilated temporal relational adversarial network for generic video summarization,”  Multimedia Tools and Applications, vol. 78, no. 24, pp. 35237–35261, 2019.
  * [18] Melissa Sanabria, Frédéric Precioso, and Thomas Menguy,  “A deep architecture for multimodal summarization of soccer games,”  in International Workshop on Multimedia Content Analysis in Sports, 2019, pp. 16–24.
  * [19] Huawei Wei, Bingbing Ni, Yichao Yan, Huanyu Yu, Xiaokang Yang, and Chen Yao,  “Video summarization via semantic attended networks,”  in AAAI, 2018, vol. 32.
  * [20] Jia-Hong Huang, Luka Murn, Marta Mrak, and Marcel Worring,  “Gpt2mvs: Generative pre-trained transformer-2 for multi-modal video summarization,”  in ICMR, 2021, pp. 580–589.
  * [21] Jia-Hong Huang, Modar Alfadly, Bernard Ghanem, and Marcel Worring,  “Assessing the robustness of visual question answering,”  arXiv:1912.01452, 2019.
  * [22] Jia-Hong Huang, Cuong Duc Dao, Modar Alfadly, and Bernard Ghanem,  “A novel framework for robustness analysis of visual qa models,”  in AAAI, 2019, vol. 33, pp. 8449–8456.
  * [23] Jia-Hong Huang, Modar Alfadly, and Bernard Ghanem,  “Vqabq: Visual question answering by basic questions,”  VQA Challenge Workshop, CVPR, 2017.
  * [24] Tao Hu, Pascal Mettes, Jia-Hong Huang, and Cees GM Snoek,  “Silco: Show a few images, localize the common object,”  in ICCV, 2019, pp. 5067–5076.
  * [25] Yujie Li, Atsunori Kanemura, Hideki Asoh, Taiki Miyanishi, and Motoaki Kawanabe,  “Extracting key frames from first-person videos in the common space of multiple sensors,”  in ICIP. IEEE, 2017, pp. 3993–3997.
  * [26] Jia-Hong Huang, Modar Alfadly, and Bernard Ghanem,  “Robustness analysis of visual qa models by basic questions,”  VQA Challenge and Visual Dialog Workshop, CVPR, 2018.
  * [27] Jia-Hong Huang, Ting-Wei Wu, and Marcel Worring,  “Contextualized keyword representations for multi-modal retinal image captioning,”  in ICMR, 2021, pp. 645–652.
  * [28] Jia-Hong Huang, Ting-Wei Wu, Chao-Han Huck Yang, and Marcel Worring,  “Deep context-encoding network for retinal image captioning,”  in ICIP. IEEE, 2021, pp. 3762–3766.
  * [29] Jia-Hong Huang, Ting-Wei Wu, Chao-Han Huck Yang, and Marcel Worring,  “Longer version for” deep context-encoding network for retinal image captioning”,”  arXiv preprint arXiv:2105.14538, 2021.
  * [30] C-H Huck Yang, Fangyu Liu, Jia-Hong Huang, Meng Tian, I-Hung Lin, Yi Chieh Liu, Hiromasa Morikawa, Hao-Hsiang Yang, and Jesper Tegner,  “Auto-classification of retinal diseases in the limit of sparse data using a two-streams machine learning model,”  in ACCV. Springer, 2018, pp. 323–338.
  * [31] Yi-Chieh Liu, Hao-Hsiang Yang, C-H Huck Yang, Jia-Hong Huang, Meng Tian, Hiromasa Morikawa, Yi-Chang James Tsai, and Jesper Tegner,  “Synthesizing new retinal symptom images by multiple generative models,”  in ACCV. Springer, 2018, pp. 235–250.
  * [32] C-H Huck Yang, Jia-Hong Huang, Fangyu Liu, Fang-Yi Chiu, Mengya Gao, Weifeng Lyu, Jesper Tegner, et al.,  “A novel hybrid machine learning model for auto-classification of retinal diseases,”  Workshop on Computational Biology, ICML, 2018.
  * [33] Riccardo Di Sipio, Jia-Hong Huang, Samuel Yen-Chi Chen, Stefano Mangini, and Marcel Worring,  “The dawn of quantum natural language processing,”  ICASSP, 2022.
  * [34] Jia-Hong Huang,  “Robustness analysis of visual question answering models by basic questions,”  King Abdullah University of Science and Technology, Master Thesis, 2017.
  * [35] Jia-Hong Huang, C-H Huck Yang, Fangyu Liu, Meng Tian, Yi-Chieh Liu, Ting-Wei Wu, I Lin, Kang Wang, Hiromasa Morikawa, Hernghua Chang, et al.,  “Deepopht: medical report generation for retinal images via deep models and visual explanation,”  in WACV, 2021, pp. 2442–2452.
  * [36] Jia-Hong Huang, Ting-Wei Wu, C-H Huck Yang, Zenglin Shi, I Lin, Jesper Tegner, Marcel Worring, et al.,  “Non-local attention improves description generation for retinal images,”  in WACV, 2022, pp. 1606–1615.
  * [37] Peter A Frost,  “Proxy variables and specification bias,”  The review of economics and Statistics, p. 323, 1979.
  * [38] Dong Zhang, Hanwang Zhang, Jinhui Tang, Xiansheng Hua, and Qianru Sun,  “Causal intervention for weakly-supervised semantic segmentation,”  NIPS, 2020.
  * [39] Sander Greenland and DAVID G KLEINBAUM,  “Correcting for misclassification in two-way tables and matched-pair studies,”  International Journal of Epidemiology, vol. 12, no. 1, pp. 93–97, 1983.
  * [40] Jan Selén,  “Adjusting for errors in classification and measurement in the analysis of partly and purely categorical data,”  Journal of the American Statistical Association, vol. 81, no. 393, pp. 75–81, 1986.
  * [41] Ehsan Abbasnejad, Damien Teney, Amin Parvaneh, Javen Shi, and Anton van den Hengel,  “Counterfactual vision and language learning,”  in CVPR, 2020, pp. 10044–10054.
  * [42] Wang Miao, Zhi Geng, and Eric J Tchetgen Tchetgen,  “Identifying causal effects with proxy variables of an unmeasured confounder,”  Biometrika, vol. 105, no. 4, pp. 987–993, 2018.
  * [43] Judea Pearl,  “On measurement bias in causal inference,”  arXiv preprint arXiv:1203.3504, 2012.
  * [44] Manabu Kuroki and Judea Pearl,  “Measurement bias and effect restoration in causal inference,”  Biometrika, vol. 101, no. 2, pp. 423–437, 2014.
  * [45] Diederik P Kingma and Max Welling,  “Auto-encoding variational bayes,”  arXiv:1312.6114, 2013.
  * [46] Uri Shalit, Fredrik D Johansson, and David Sontag,  “Estimating individual treatment effect: generalization bounds and algorithms,”  in ICML, 2017, p. 3076.
  * [47] Jun Fu, Jing Liu, Haijie Tian, Yong Li, Yongjun Bao, Zhiwei Fang, and Hanqing Lu,  “Dual attention network for scene segmentation,”  in CVPR, 2019, p. 3146.
  * [48] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Lukasz Kaiser, and Illia Polosukhin,  “Attention is all you need,”  arXiv preprint arXiv:1706.03762, 2017.
  * [49] George Hripcsak and Adam S Rothschild,  “Agreement, the f-measure, and reliability in information retrieval,”  Journal of the American medical informatics association, vol. 12, no. 3, pp. 296–298, 2005.
  * [50] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun,  “Deep residual learning for image recognition,”  in CVPR, 2016, pp. 770–778.
  * [51] Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei,  “Imagenet: A large-scale hierarchical image database,”  in CVPR. IEEE, 2009, pp. 248–255.
  * [52] Diederik P Kingma and Jimmy Ba,  “Adam: A method for stochastic optimization,”  arXiv:1412.6980, 2014.