# Flexible hierarchical risk modeling for large insurance data via NumPyro

Christopher Krapu and Mark Borsuk
Department of Civil & Environmental Engineering
Duke University

(November 2023)

## Abstract

Data analysis and individual policy-level modeling for insurance involves
handling large data sets with strong spatiotemporal correlations, non-Gaussian
distributions, and complex hierarchical structures. In this research, we
demonstrate that by utilizing gradient-based Markov chain Monte Carlo (MCMC)
techniques accelerated by graphics processing units, the trade-off between
complex model structures and scalability for inference is overcome at the
million-record size. By implementing our model in NumPyro, we leverage its
built-in MCMC capabilities to fit a model with multiple sophisticated
components such as latent conditional autoregression and spline-based exposure
adjustment, achieving an 8.8x speedup compared to CPU-based implementations.
We apply this model to a case study of 2.6 million individual policy-level
claim count records for automobile insurance from Brazil in 2011. We
illustrate how this modeling approach significantly advances current risk
assessment processes for numerous, closely related outcomes. The code and data
are available at https://github.com/ckrapu/bayes-at-scale.

##  1 Introduction

Modeling for risk assessment and prediction in insurance is challenging due to
the tension between (1) the desire for accurate forecasts of expected claim
counts and/or severity, and (2) the requirement for reliable uncertainty
quantification. As total policy value typically exceeds available liquid
assets for most types of modern insurance, a key objective in this setting is
to forecast the probability of ruinous losses under a chosen ratemaking and
loss adjustment scheme. Unfortunately, this probability does not allow a
straightforward factorization into independent distribution functions because
of strong correlations across policies stemming from similarities in
geography, demographics, property type, and other factors. Consequently, an
effective modeling workflow must account for as many of these correlation
sources as possible. A vitally important output of such statistical modeling
is the predictive covariance matrix which, for new data points, quantifies
their cross-correlation structure. Although a machine learning-based approach
can utilize large data sets for accurate record-level predictions, which may
be highly useful operationally [Ding et al., 2020], it is uncertain whether
any current research for calibrated uncertainty quantification from machine
learning [Kuleshov et al., 2018, Angelopoulos and Bates, 2022] applies to the
multivariate context induced by multiple policies or outcomes from the same
policy. Moreover, it is often desirable for models to be useful in both small-
and large-data situations such as in [Zhang, 2017], and to be structured in a
way that allows for some interpretability. These requirements are not yet met
by the machine learning models frequently used in insurance modeling [Lupton,
2022].

### Parameter estimation and uncertainty quantification

With regard to the statistical approach, in scenarios where the data
generating process, synonymous with _forward model_ , can be written down
easily, probabilistic inference might not be straightforward or
computationally feasible for a maximum likelihood or restricted maximum
likelihood approach (see [Wahl et al., 2022] for a recent example). Generally,
these methods will require a bespoke algorithm for model fitting and producing
parameter estimates, restricting their usability to experts in statistical
subfields. Packages designed for a broader range of probabilistic models [Lunn
et al., 2009, Carpenter et al., 2017, Salvatier et al., 2016] using sampling-
based Monte Carlo inference methods like Markov chain Monte Carlo (MCMC)
[Hastings, 1970, Brooks, 1998, Geyer, 1992] have found widespread acceptance.
A major challenge of MCMC methods, including Gibbs sampling and the
Metropolis-Hastings algorithm, is the substantial time complexity of
$\mathcal{O}(d^{2})$ where $d$ denotes the dimension of the target density,
generally equivalent to the combined number of latent variables and/or
parameters in a model. This quadratic scaling arises from a general
$\mathcal{O}(d)$ time to evaluate the model log posterior, scaling with the
number of parameters, and a theoretical $\mathcal{O}(d)$ mixing time [Roberts
and Rosenthal, 1998, Beskos and Stuart, 2009] for a broad class of
distributions. Taken together, these lead to the $\mathcal{O}(d^{2})$ scaling.
Alternatives to MCMC for Bayesian model fitting include stochastic variational
inference [Hoffman et al., 2013] and the integrated nested Laplace
approximation reviewed [Rue et al., 2016] and benchmarked against MALA in
[Taylor and Diggle, 2014]. A more comprehensive review of approximate Bayes
methods is given in [Martin et al., 2023].

However, we now have MCMC methods that incorporate proposal generation using
gradients of the log-posterior with regard to latent variables and parameters
[Duane et al., 1987, Roberts and Tweedie, 1996, Neal, 2011]. For certain
distributions, these methods can have $\mathcal{O}(d^{4/3})$ or even
$\mathcal{O}(d^{5/4})$ time complexity for the Metropolis-adjusted Langevin
algorithm (MALA) and Hamiltonian Monte Carlo, respectively. These algorithms
are a substantial field of research and are explained thoroughly in
[Betancourt, 2017]. An obstacle to the general use of these MCMC algorithms is
the need for gradients of the log-posterior density. However, automatic
differentiation software like Torch [Paszke et al., 2019], TensorFlow [Abadi
et al., 2015], Jax [Frostig et al., 2018], and Theano [Bergstra et al., 2010]
enable this and integrate with probabilistic programming frameworks like Stan
[Carpenter et al., 2017], PyMC [Salvatier et al., 2016, Abril-Pla et al.,
2023], Pyro [Bingham et al., 2019], and NumPyro [Phan et al., 2019] for
statistical modeling. Techniques for scaling differentiation to large models
in machine learning also apply to MCMC, benefiting from GPU-based parallel
computing for rapid model fitting. The No-U-Turn sampler [Hoffman and Gelman,
2014] has become a general-purpose method due to few tuning parameters and
diverse applicability. Initial efforts [Chin and McNulty, 2023] have exposed
the insurance modeling community to a probabilistic programming approach, yet
there is an unmet need in the literature to communicate the workflow’s pros
and cons against the unique requirements of risk assessment and insurance
modeling. In this work, we use NumPyro, recognizing several promising features
that could benefit the modeling community.

To summarize, by writing models in probabilistic programming frameworks, we
can achieve both flexibility and scalability. A detailed discussion of
probabilistic programming, particularly in environmental statistics, is
provided in [Krapu and Borsuk, 2019]. By employing general, model-agnostic
inference methods, we can incorporate various model components; if we code the
forward model in such a framework, we then automatically leverage GPU with
minimal effort. We illustrate this by developing a log-additive spatial count
regression model with complex submodel components. In the Methods section, we
describe our statistical model and its intended data. The Results section
presents summaries and analyses of the model’s inferential outputs, while the
Discussion covers limitations and further improvement opportunities.

##  2 Methods

To provide context to the material from the previous section, we designed a
case study to understand the capabilities and constraints of the proposed
approaches. This study involves an analysis of 2.6 million auto insurance
policies from _SUSEP_ , the Brazilian federal insurance supervisory agency.
The records, originally compiled for the supplemental material of the book
_Computational Actuarial Science with R_ [Charpentier, 2015, Dutang and
Charpentier, 2022], detail individual auto policies from 2011, which may
include claims. Our focus, however, is solely on the collision data. Each
record is linked to a municipality name and contains the automobile brand
(e.g., Volkswagen, Kia, Iveco) and the model (e.g., Taurus, Camry, F-150). We
preserved the original brand categories but simplified the model
classification into 15 vehicle types depicted in Figure 7. Each record also
specifies the vehicle year and an exposure value.

### Model

The general form of the model considered in this work is a log-additive
regression for count data, which records the number of collision claims per
policy. We use a Poisson likelihood, expressing the $i$-th policy’s Poisson
rate parameter $\lambda_{i}$ as

$\displaystyle y_{i}$ $\displaystyle\sim\mathrm{Poisson}(\lambda_{i})$ (1)
$\displaystyle\log\lambda_{i}$
$\displaystyle=\underbrace{\log{\alpha_{i}}}_{\textbf{Base}\atop\textbf{exposure}}+\underbrace{g(\alpha_{i})}_{\textbf{Exposure}\atop\textbf{adjustment}}+\underbrace{\sum_{k=1}^{P}v_{x_{ik}}}_{\textbf{Categorical}\atop\textbf{predictors}}+\underbrace{u_{j[i]}}_{\textbf{City}\atop\textbf{effect}}+\underbrace{S_{t_{i}}}_{\textbf{Time}\atop\textbf{effect}}$
(2)

where $\alpha_{i}$ denotes the base exposure, $g(\alpha_{i})$ denotes an
exposure adjustment using splines, $\sum_{k=1}^{P}v_{x_{ik}}$ captures
contributions from policy-level categorical variables related to the vehicle’s
brand and category, $u_{j[i]}$ is a city-level random effect associated with
the city $j$ of policy $i$, and $S(t_{i})$ represents a Gaussian random walk
over time with $t_{i}$ indicating the vehicle’s manufacturing year for policy
$i$. We will provide further details on these terms below. There are 7,756
free parameters and latent variables in this model, with most arising from the
latent city-level spatial effect. Figure 1 presents a visual representation of
this model in plate notation.

Figure 1: Graphical representation of the statistical model. Shaded circles
denote observed variables, while unshaded circles denote latent variables and
parameters. Rectangular plates indicate replication of the variables within
the plate.

#### Exposure adjustment

While the exposure score is highly correlated with observed claim counts and
provides substantial predictive power, early analysis indicated a monotonic
relationship between log exposure and log claim count that strongly deviated
from linearity, particularly at the higher percentiles of the exposure
distribution. Considering the data includes various claim types beyond
collision, the exposure score likely represents a broader risk measure, rather
than being finely tuned for collision claims alone. To adjust this exposure
score, we modeled an additive effect on log exposure $\log\alpha_{i}$ using
basis splines to construct a nonlinear function of $\log\alpha_{i}$. The
nonlinear function $g(\alpha_{i})$ is defined as:

$\displaystyle\sigma_{g}$ $\displaystyle\sim\mathrm{HalfNormal}(1)$ (3)
$\displaystyle c_{l}$
$\displaystyle\overset{\text{iid}}{\sim}\mathrm{N}(0,\sigma_{g})$ (4)
$\displaystyle g(\alpha_{i})$
$\displaystyle=\sum_{l=1}^{L}b_{l}(\log\alpha_{i})\cdot c_{l}$ (5)

where $b_{l}(\cdot)$ denotes the basis spline functions, $L$ is the number of
spline coefficients, and $c_{l}$ are the spline coefficients estimated by the
model. We positioned the spline knots at evenly spaced percentiles of the log-
transformed exposure values to ensure a uniform distribution throughout the
range of the data. This technique captures subtle variations in the exposure-
claim count relationship, especially in the distribution’s tails. We employed
the interpolate module from the software package SciPy [Virtanen et al., 2020]
to determine $b_{l}$. We acknowledge that [Basile et al., 2014] combines
nonlinear modeling with spatially-correlated latent variables, similarly to
our approach.

#### Categorical predictors

Each policy record is associated with a vehicle brand and category,
represented as discrete, unordered variables. Consequently, in equation (2),
$P=2$, but we could include additional covariates that are continuous or
discrete. Letting $n_{k}$ denote the number of classes in the $k$-th
categorical variable, we have

$\displaystyle\sigma^{(k)}_{v}$ $\displaystyle\sim\mathrm{HalfNormal}(1)$ (6)
$\displaystyle v^{(a)}_{1},...,v^{(k)}_{n_{k}}$
$\displaystyle\overset{\text{iid}}{\sim}\mathrm{N}(0,\sigma^{(k)}_{v})$ (7)

The interpretation of the $v^{(k)}_{1},...,v^{(k)}_{n_{k}}$ variables is that
they adjust the risk according to the vehicle’s brand and category, such as a
motorcycle or pickup. The priors used here could be used with policy-level
continuous covariates, although none are present in the current data set. For
an example of a model that incorporates per-category effects for auto
insurance modeling, see [Frees and Valdez, 2008]

#### City effect

Most insured assets are tied to specific locations or areas on Earth, so it is
common to address unmodeled variation using spatial statistical techniques.
Examples include kriging to account for spatial correlation in drought
severity [Paulson et al., 2010], employing spatially autocorrelated indicator
variables for feature selection [Paulson et al., 2010], and using a Gaussian
line process to represent hail damage Miralles et al. [2023]. For our study,
policy records are categorized by cities; therefore, we treat these as _areal_
data suitable for modeling with conditional autoregressions. Modeling spatial
correlation at the level of individual policies is impractical because we do
not have spatial coordinates more precise than the policyholder’s city. As a
result, we opt to model at the city level, utilizing spatial priors to account
for correlations across the $J$ components of the city effect vector $\bm{u}$,
which we model as a latent conditional autoregression (CAR) [Cressie and Chan,
1989]. There are $J=3785$ cities and municipalities, and we organize the
autoregressive prior for $u_{j}$, the effect from city $j$, around a modified
Besag-York-Mollie (BYM) prior [Besag et al., 1991], reparameterized to enhance
identifiability [Simpson et al., 2015], with city-level covariates

$\displaystyle\sigma_{u},\sigma_{\epsilon}$
$\displaystyle\sim\mathrm{HalfNormal}(1)$ $\displaystyle\delta_{j}$
$\displaystyle\overset{\text{iid}}{\sim}\mathrm{N}(0,1)$ $\displaystyle u_{j}$
$\displaystyle=\sum_{m=1}^{M}\gamma_{m}z_{mj}+\epsilon_{j}$
$\displaystyle\phi$ $\displaystyle\sim\mathrm{Beta}(1,1)$ $\displaystyle\rho$
$\displaystyle\sim\mathrm{Beta}(2,2)$ $\displaystyle\bm{\epsilon}$
$\displaystyle=\sigma_{\epsilon}\left(\sqrt{1-\phi}\bm{\delta}+\sqrt{\phi}\bm{\eta}\right)$
$\displaystyle\bm{\eta}$
$\displaystyle\sim\mathrm{N}\left(\bm{0},\bm{D}(\bm{I}_{J}-\rho\bm{D^{-1}W})^{-1}\right)$

Here, the city effect is represented by city-level covariates $z_{mj}$,
independent Gaussian noise $\delta_{j}$, and spatially correlated terms
$\eta_{j}$. $W$ is the spatial adjacency matrix where $W_{ab}=1$ if $a$ is a
neighbor of $b$. This matrix is created by finding the five nearest
neighboring cities for each city and setting the corresponding entries in $D$
to 1. We then make $W$ symmetric by ensuring that $W_{ab}=1$ if $W_{ba}=1$ for
all distinct $a$ and $b$. Non-neighboring entries in $D$ remain zero.
Therefore, $D$ is a $J$-dimensional diagonal matrix where $D_{jj}$ is the
count of neighbors for city $j$, which must be at least five. For
illustrations of latent Gaussian Markov random field applications, see [Jin et
al., 2005, Gelfand and Vounatsou, 2003]. The parameter $\rho$ represents the
spatial autocorrelation level; the $\mathrm{Beta}(2,2)$ prior promotes edge-
avoiding behavior away from 0 or 1 fosters more effective MCMC convergence.
The $\phi$ parameter dictates the reliance of the city effect on uncorrelated
($\delta_{j}$) versus correlated terms ($\eta_{j}$) when assessing the city’s
influence on policyholder risk. We apply the sparse CAR likelihood in NumPyro,
leading to $\mathcal{O}(J)$ computation of $\log p(\bm{\eta})$, as detailed in
[Jin et al., 2005].

Figure 2: Adjacency network for conditional autoregressive city-level effect.
Blue circles denote cities as identified by the Google Maps Geocoding API.
Gray lines depict the city adjacency graph $W$ edges established for our
model’s conditional autoregression. Edges are symmetric, and the number of
neighboring cities can vary.

#### Temporal effect

To include vehicle age in the model structure in a way that is flexible and
possibly nonlinear, we created a latent time series. This series is a Gaussian
random walk spanning from years 1971 to 2011, with yearly increments. The
purpose of this approach is to construct a time-varying random effect [Boucher
et al., 2008] that ensures modest year-to-year variations in the effect but
allows for significant overall changes over the 40 years without prescribing a
specific parametric form.

$\displaystyle\sigma_{\xi}$ $\displaystyle\sim\mathrm{HalfNormal}(1)$
$\displaystyle\xi_{t}$
$\displaystyle\overset{\text{iid}}{\sim}\mathrm{N}(0,\sigma_{\xi})$
$\displaystyle S_{t^{\prime}}$ $\displaystyle=\sum_{t<t^{\prime}}\xi_{t}$

The variables $\xi_{t}$ represent increments in the latent time series, and
$S_{t^{\prime}}$ is the corresponding cumulative sum. This model is not suited
for continuous temporal coordinates; for such cases, a latent Gaussian process
[Jia et al., 2023] or changepoint model [Chib, 1998] would be preferable.

#### Scale priors

To fully specify a Bayesian model, we must assign prior distributions to scale
or variance parameters such as $\sigma_{x}i$ and $\sigma_{g}$, ensuring these
priors are defined only for nonnegative values. We prefer the half-normal or
folded-normal prior for these parameters as they impart minimal information on
the log scale; the $\mathrm{HalfNormal}(1)$ prior places considerable
probability mass on values ranging from 0 to 3. In terms of categorical
effects, a log scale value of 3 indicates that expected claim counts are
$\exp(3)\approx 20$ times greater on the original scale of the data. This
concept is elaborated in [Gelman and Carpenter, 2020]. Other priors such as
the Half-Cauchy [Gelman, 2006] or horseshoe [Carvalho et al., 2009] are
alternatives for stronger regularization.

All modeling components were executed in NumPyro and Jax, utilizing Jax’s
alternative module for NumPy, jax.numpy. The following section delves into
data preparation and summary statistics.

##  3 Data

To fit the model from the previous section, we used 2,658,372 records of auto
insurance policies from 2011 in Brazil, as collected by SUSEP, the Brazilian
insurance oversight agency, and collated as the brvehins2 data set maintained
by [Dutang and Charpentier, 2022] as supporting material for the textbook
_Computational Actuarial Science in R_ [Charpentier, 2015].

We removed 944 rows with vehicle manufacture dates prior to 1971 to avoid
introducing extra parameters due to the time-varying effect corresponding to a
small sample size. We also discarded an additional 8,346 rows that had
exposure values set exactly to zero. 141,004 records lacked information on
vehicle type; these were also removed from the data set. From the vehicle
model information, we extracted the vehicle brand ($n_{brand}=104$) and the
vehicle category ($n_{category}=15$); both are listed in Figure 7. This
classification data from SUSEP has proven to be a strong predictor of the risk
of theft in previous studies [de Azevêdo et al., 2023]. [Peres et al., 2019]
provides a comprehensive analysis of the Brazilian auto insurance market.

Nonzero claim counts are relatively rare in this data set, with 92% of
policies reporting no collision claims during 2011. The marginal distributions
of claim counts and exposure values are shown in Figure 3. We observe a
distinct peak in the histogram of exposure values at the median value, which
may indicate flawed reporting, the application of a standard exposure scoring
formula to a large number of similar inputs, or another cause. We chose not to
intervene or further investigate this anomaly.

Figure 3: Histograms for the observed claim counts and exposure values for all
policies considered.

Originally, this data set did not include coordinates linking each city or
municipality with its geographic location. To address this, we employed the
Google Maps Geocoding API to map each of 3,785 cities to a city/state/country
string. For instance, ’Imbuia, Santa Catarina, Brazil’ refers to the city of
Imbuia in the state of Santa Catarina. This method successfully linked all
unique cities to latitude and longitude. 1,404 rows listed the state as
‘’Santa Catarina‘’, but had no city or city code and were listed as Area value
of ’Oeste’. These were renamed to ’Oeste Catarinense’. The states with the
most policies were Sao Paulo (860,584), Minas Gerais (346,948), and Rio Grande
do Sul (339,644). The states with the fewest policies were Piaui (5,112), Acre
(2,632), and Rondonia (1,928).

Having identified coordinates for each city, we performed a spatial join with
a data layer for the municipal boundaries. Although the linked municipalities
nearly cover the more populated areas of Brazil, many are geographically
isolated with no adjacent municipalities, preventing us from creating an
adjacency matrix based on the municipal boundaries.

To understand the spatial structure of collision risk, we included city-level
covariates. For this, we calculated zonal averages of several globally-
available raster data layers using the municipal boundaries for each city
through Google Earth Engine [Gorelick et al., 2017]. The layers used for this
covariate extraction are listed below. We calculated the mean for each pixel
touching a municipal polygon in each case.

Water occurrence: For each pixel, the Global Surface Water project [Pekel et
al., 2016] monitors the proportion of time water was observed from 1984 - 2022
using Landsat data at a 30 m spatial resolution.

Precipitation: The WorldClim project [Fick and Hijmans, 2017] provides long-
term monthly climate data interpolated from weather stations from 1970-2000 at
a spatial resolution of 1 km. We selected average rainfall for December,
typically the wettest month.

Population density: To reflect the increased need for collision insurance with
higher population density [Sherden, 1984], we utilized the WorldPop project’s
[Tatem, 2017] estimate of population density for 2011 at a 100 m resolution.

Enhanced vegetation index: The MODIS satellite data set provides an average
EVI for January 2011, used as an indicator of urban development and
agricultural land use [Huete et al., 2002].

Elevation: The elevation data from the Shuttle Radar Topography Mission, with
a 90 m resolution, serve as a simple proxy for distance from the coast and are
correlated with development patterns.

Topographic diversity: The topographic diversity index from [Theobald et al.,
2015], measured at a 270 m resolution, helps identify areas with considerable
variations in elevation or landforms.

Forest cover: Including the global forest cover data at a 270 m resolution
from [Shimada et al., 2014] provides additional land use information.

Travel friction index: The Global Friction Surface data sets provide indices
for the difficulty of motorized and nonmotorized travel at a 1 km resolution,
relevant for average trip velocity and collision risk [Weiss et al., 2018].

The specific image IDs and bands applied are detailed in Table 1. All
covariates were standardized to have a mean of zero and unit variance before
modeling. Scatter plots illustrating cross-variable correlations and marginal
distributions are presented in Figure 4.

Figure 4: City-level covariates obtained through zonal averaging. Correlations for each pair of variables are displayed in the upper right of the subplot. Covariate | Image ID | Band name | Spatial resolution
---|---|---|---
Water | JRC/GSW1_4/GlobalSurfaceWater | occurrence | 30 meters
Population | WorldPop/GP/100m/pop/BRA_2011 | population | 100 meters
Vegetation | MODIS/061/MOD13A3/2011_01_01 | EVI | 1000 meters
Elevation | CGIAR/SRTM90_V4 | elevation | 90 meters
Topography | CSP/ERGo/1_0/Global/ALOS_topoDiversity | constant | 270 meters
Forest | JAXA/ALOS/PALSAR/YEARLY/FNF4/2017 | fnf | 270 meters
Motorized travel | Oxford/MAP/friction_surface_2019 | friction | 1000 meters
Nonmotorized travel | Oxford/MAP/friction_surface_2019 | friction_walking_only | 1000 meters
Rainfall | WORLDCLIM/V1/MONTHLY/12 | prec | 1000 meters
Table 1: Earth Engine images used for city-level covariates. All images were
processed using the zonal mean to produce a city-level covariate matrix with
dimensions of 3785 by 9

### Model Fitting and Benchmarking

#### Benchmarking GPU-based acceleration

We compared the time required to evaluate the model’s log-posterior density
using the full data set mentioned previously. This comparison is more suitable
than measures like minimal effective sample size per second because it avoids
bias from inefficiencies in specific implementations of Hamiltonian Monte
Carlo or the No-U-Turn sampler, and directly targets the fundamental
computational challenges in performing gradient-based MCMC for our purpose.
Calculating the gradient 100 times on both GPU and CPU, we recorded a median
evaluation time of 5 milliseconds on the GPU and 44 milliseconds on the CPU.
This indicates approximately an 8.8 times speed increase when moving from CPU
to GPU. All experiments were conducted on a Lambda Labs virtual machine
outfitted with an Nvidia A10 GPU with Nvidia driver version 525.85.12 and CUDA
version 12.0. The Jax and NumPyro versions installed were 0.4.13 and 0.12.1,
respectively. The outcomes from this benchmarking are consistent with the
runtime for the entire model fitting process detailed in the next paragraph.

#### Model fitting

To obtain parameter estimates under the specified model, we ran the No-U-Turn
sampler in NumPyro for 2000 warmup iterations and 2000 sampling iterations.
The warmup iterations usually take place before the Markov chain has converged
to the stationary distribution and are discarded. We retained every 20th
sample due to GPU memory constraints. For finding a starting point for
sampling, we used the built-in stochastic variational inference (SVI) [Hoffman
et al., 2013] optimization routine for 200,000 iterations and the Adam
optimizer with a learning rate of $10^{-3}$. The chains were run sequentially;
the total runtime of the entire model fitting procedure was 6.7 hours, with
chain 1 requiring 194 minutes, chain 2 requiring 189 minutes, and the SVI
initialization requiring 19 minutes. To assess MCMC convergence, we calculated
$\hat{R}$ [Gelman and Rubin, 1992] and effective sample size, noting all
$\hat{R}$ values were under 1.10 and all effective sample sizes were over 35.
The variables with the largest $\hat{R}$ values and smallest effective sample
sizes were elements of $\bm{\eta}$, the spatially-correlated contributions to
the city-level random effect. Based on these statistics, we found no evidence
of non-convergence of the Markov chains.

##  4 Results

### Model fit

Figure 5 displays several visual indicators of in-sample model fit. Since a
Poisson likelihood was adopted, the relative frequency of predicted claim
counts falling outside the predictive $2\sigma$ credible intervals implied by
the Poisson model ($\lambda_{i}\pm 2\sqrt{\lambda_{i}})$ is a critical measure
of overdispersion. Although Figure 5A illustrates this graphically, an
overdispersion test was also conducted by calculating the ratio of the sum of
squared Pearson residuals to the degrees of freedom, yielding a test statistic
of 1.107 which indicates minimal evidence of overdispersion. Regarding
predictive means, the observed and predicted quantiles match well, as shown in
Figure 5D. Figure 5B illustrates a strong correlation between predicted and
true claim counts which increases with larger $\lambda_{i}$, consistent with
the $\sqrt{\lambda_{i}}$ scaling of expected deviation for a Poisson
distribution. The posterior mean values of $\lambda_{i}$ were used to
determine the per-record log-likelihood and absolute error. On average, a mean
log-likelihood of $-0.27$ and mean absolute error of $0.16$ were produced. The
coefficient of determination specific to Poisson-distributed data
$R^{2}_{P,P}$ was computed using Pearson residuals [Cameron and Windmeijer,
1996], resulting in $R^{2}_{P,P}=0.94$. It should be noted that this model
incorporates the exposure score as a predictive element, which could be the
result of a prior modeling effort.

Figure 5: Relationship between true and in-sample posterior predictive means.

### Nonlinear exposure adjustment

The posterior estimate of $g(\alpha_{i})$ displayed in Figure 6 indicates a
relatively stable relation between exposure and $g$ for exposures in the range
$10^{-2}$ to $1.0$. From $1.0$ upwards, there is an increasing linear trend.
Since this model lacks a global intercept and regularizes the parameters
towards zero through scale priors with modes at zero, the posterior mean of
$g(\alpha_{i})$ is not centered at zero. In summary, these estimates indicate
that the optimal adjustment to the exposure score remains fairly constant for
exposure values near the modeled range from $10^{-1}$ to $1$, with a more
pronounced adjustment for policies with high exposure values.

Figure 6: Posterior estimates of nonlinear exposure adjustment

### Categorical predictors

As demonstrated in Figure 7, the posterior estimates for the categorical
predictors vary significantly, with some categories showing notable increases
in expected claim counts and others decreases. Brands with large impacts
include _Engresa_ , _Subaru_ , and _Ducati_. Approximately one quarter of the
brands showed posterior coefficient distributions markedly different from
zero, as determined by 95% credible intervals. In terms of vehicle type,
sports cars tend to be associated with higher claim counts, whereas trucks,
buses, compact cars, and motorcycles typically have lower claim counts.
Nonetheless, only compact cars and motorcycles display significantly negative
coefficients, according to 95% credible intervals.

Figure 7: Posterior summaries for parameters linked to vehicle brand and
category. Due to space constraints, we show only every fourth brand. All
categories are displayed. Lines indicate 2$\sigma$ credible interval bounds
and dots represent posterior mean estimates of model parameters.

### Spatial effect

Despite the inclusion of several city-level covariates, we find that none of
these covariates’ coefficients have values significantly different from zero,
assessed via the 95% posterior credible intervals (Figure 8). The posterior
mean values for these coefficients suggest that population density and
vegetation may be positively associated with claim counts, while rainfall and
forest cover appear to be negatively associated with claim counts. There is
also evidence of spatial autocorrelation in $\bm{\eta}$ values, demonstrating
the utility of this model component in identifying patterns of risk which are
not fully explained by covariates, as depicted in Figure 9.

Figure 8: Posterior summaries for parameters relating city-level covariates to
predicted counts. Lines indicate 2$\sigma$ credible interval bounds and dots
indicate posterior means. Figure 9: Posterior mean estimates of spatially
correlated values of $\eta$ in the state of Rio Grande do Sul.

### Temporal effect

The impact of car manufacture date on expected collision claim counts is shown
in Figure 10. The relationship between car age and claim count is relatively
level for vehicles aged 20 to 30 years. We also noted a nearly linear increase
in this effect for cars younger than 20 years, with the trend intensifying as
vehicle age decreases. In summary, the risk assessed is stable for vehicles
ranging from 20 to 30 years old and decreases with the age for cars under 20
years old.

Figure 10: Posterior summary for time-varying effect.

##  5 Discussion

The case study presented in this work has integrated several modeling
subcomponents to address different aspects of risk modeling for insurance. We
note while elements such as the spline-based exposure adjustment, vehicle
category, and time-varying effect demonstrated strong signals manifesting as
coefficients or parameters with posterior credible intervals excluding zero,
the city-level covariates did not. We also observed a relatively large number
of model parameters and coefficients with wide posterior credible intervals
crossing zero, suggesting that the data is not strong enough to constrain
these values to useful ranges. A major limitation of this work is the absence
of shrinkage priors, which are designed specifically for sparsification
[Carvalho et al., 2009], preventing us from applying substantial
regularization to push parameter estimates towards zero when needed. Employing
such priors would also allow for more model components, like feature-brand
interactions or varying effects for the categorical predictors. A more
comprehensive analysis could include policy-level covariates like driver age,
driving record, and other risk or exposure indicators. This was intended as a
demonstration for insurance-relevant models using large ($>1000$) parameter
and data sets in a general-purpose statistical framework like NumPyro.
Nonetheless, there is no barrier to creating models with significantly more
parameters and data, given enough computational resources. We anticipate
future studies using multi-GPU setups, enabling models with $10^{5}$ or more
parameters and $>10^{8}$ observations.

### Future work

We anticipate further interest in statistical modeling for insurance with
large datasets and substantial models, and we view several avenues as
particularly relevant directions for future work. There exists considerable
potential for extending the model to include multiple types of counts with
interesting correlation structures [Bermúdez and Karlis, 2011] to help share
information from different observation types [Anastasiadis and Chukova, 2012].
For example, we can imagine a scenario where the risk due to robbery and the
risk due to fire are highly correlated because of underlying factors related
to emergency services and population density. Individual models of each type
may be fit on data without a strong enough signal to constrain parameter
estimates; together, however, a multivariate model may share information and
consequently enable more insights into joint risk factors. To accommodate data
with varying degrees of dispersion reflecting different amounts of underlying
heterogeneity, using likelihoods suitable for multiple dispersion types
Guikema and Coffelt [2008] may be productive. For a more practical and useful
modeling exercise of direct relevance to insurance operations, we may want to
model claim frequency and severity together [Frees and Valdez, 2008, Chin and
McNulty, 2023] since both types of outcomes are essential for pricing and
reserving.

In early modeling iterations, we identified prominent disparities in risk
between municipalities close to major urban centers and those farther away.
While the model version ultimately presented explains a significant amount of
this variation using covariates and the spatial random effect, we are
interested in using a latent mixture model to cluster together cities with
similar excess risk patterns. Generally, this type of model is hard to
accommodate in a gradient-based MCMC framework as the latent class label is
usually discrete and thus does not admit a log-posterior which is
differentiable in all latent variables. However, both Stan and NumPyro have
the ability to either enumerate or automatically marginalize out latent
variables of modest cardinality, thus avoiding this issue and allowing for the
use of gradient-based MCMC. We would also be interested in using
computationally less expensive approximations to the Poisson likelihood
[Cameron and Johansson, 1997] to speed up model fitting.

##  6 Conclusion

This study demonstrates the value of state-of-the-art statistical modeling
frameworks like NumPyro that work seamlessly on GPUs to speed up model fitting
without sacrificing any flexibility in the selection of model components. We
developed a log-additive model for collision claim frequencies of auto
insurance policies in Brazil. This model includes a nonlinear exposure
adjustment, categorical predictors, a city-level random effect with spatial
correlation, and a dynamic temporal effect. The model succeeded in explaining
much of the data variance, though the city-level covariates showed weak
effects. Moreover, we discovered that NumPyro’s GPU-based implementation made
the model fitting process nine times faster than using a CPU, as shown in the
log posterior gradient calculations, enabling completion in a practical
timeframe. This study should interest both practitioners and researchers in
the insurance sector, as well as those who employ Bayesian methods for
analyzing extensive spatial datasets.

##  7 Data Availability

All code and data are available for replication at
https://github.com/ckrapu/bayes-at-scale.

## References

  * Abadi et al. [2015] Martın Abadi, Ashish Agarwal, Paul Barham, Eugene Brevdo, Zhifeng Chen, Craig Citro, Greg S Corrado, Andy Davis, Jeffrey Dean, Matthieu Devin, Sanjay Ghemawat, Ian Goodfellow, Andrew Harp, Geoffrey Irving, Michael Isard, Yangqing Jia, Rafal Jozefowicz, Lukasz Kaiser, Manjunath Kudlur, Josh Levenberg, Dan Mane, Rajat Monga, Sherry Moore, Derek Murray, Chris Olah, Mike Schuster, Jonathon Shlens, Benoit Steiner, Ilya Sutskever, Kunal Talwar, Paul Tucker, Vincent Vanhoucke, Vijay Vasudevan, Fernanda Viegas, Oriol Vinyals, Pete Warden, Martin Wattenberg, Martin Wicke, Yuan Yu, and Xiaoqiang Zheng.  TensorFlow: Large-Scale Machine Learning on Heterogeneous Distributed Systems.  page 19, 2015.
  * Abril-Pla et al. [2023] Oriol Abril-Pla, Virgile Andreani, Colin Carroll, Larry Dong, Christopher J. Fonnesbeck, Maxim Kochurov, Ravin Kumar, Junpeng Lao, Christian C. Luhmann, Osvaldo A. Martin, Michael Osthege, Ricardo Vieira, Thomas Wiecki, and Robert Zinkov.  PyMC: A modern, and comprehensive probabilistic programming framework in Python.  _PeerJ Computer Science_ , 9:e1516, September 2023.  ISSN 2376-5992.  doi: 10.7717/peerj-cs.1516.
  * Anastasiadis and Chukova [2012] Simon Anastasiadis and Stefanka Chukova.  Multivariate insurance models: An overview.  _Insurance: Mathematics and Economics_ , 51(1):222–227, July 2012.  ISSN 0167-6687.  doi: 10.1016/j.insmatheco.2011.01.013.
  * Angelopoulos and Bates [2022] Anastasios N. Angelopoulos and Stephen Bates.  A Gentle Introduction to Conformal Prediction and Distribution-Free Uncertainty Quantification, December 2022.
  * Basile et al. [2014] Roberto Basile, María Durbán, Román Mínguez, Jose María Montero, and Jesús Mur.  Modeling regional economic dynamics: Spatial dependence, spatial heterogeneity and nonlinearities.  _Journal of Economic Dynamics and Control_ , 48:229–245, November 2014.  ISSN 0165-1889.  doi: 10.1016/j.jedc.2014.06.011.
  * Bergstra et al. [2010] James Bergstra, Olivier Breuleux, Frédéric Bastien, Pascal Lamblin, Razvan Pascanu, Guillaume Desjardins, Joseph Turian, David Warde-Farley, and Yoshua Bengio.  Theano: A CPU and GPU Math Compiler in Python.  page 7, 2010.
  * Bermúdez and Karlis [2011] Lluís Bermúdez and Dimitris Karlis.  Bayesian multivariate Poisson models for insurance ratemaking.  _Insurance: Mathematics and Economics_ , 48(2):226–236, March 2011.  ISSN 0167-6687.  doi: 10.1016/j.insmatheco.2010.11.001.
  * 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, March 1991.  ISSN 1572-9052.  doi: 10.1007/BF00116466.
  * Beskos and Stuart [2009] Alexandros Beskos and Andrew Stuart.  Computational Complexity of Metropolis-Hastings Methods in High Dimensions.  In Pierre L’ Ecuyer and Art B. Owen, editors, _Monte Carlo and Quasi-Monte Carlo Methods 2008_ , pages 61–71, Berlin, Heidelberg, 2009\. Springer.  ISBN 978-3-642-04107-5.  doi: 10.1007/978-3-642-04107-5˙4.
  * Betancourt [2017] Michael Betancourt.  A Conceptual Introduction to Hamiltonian Monte Carlo.  _arXiv preprint:1701.02434_ , page 60, 2017.
  * Bingham et al. [2019] Eli Bingham, Jonathan P. Chen, Martin Jankowiak, Fritz Obermeyer, Neeraj Pradhan, Theofanis Karaletsos, Rohit Singh, Paul Szerlip, Paul Horsfall, and Noah D. Goodman.  Pyro: Deep Universal Probabilistic Programming.  _Journal of Machine Learning Research_ , 20(28):1–6, 2019.  ISSN 1533-7928.
  * Boucher et al. [2008] Jean-Philippe Boucher, Michel Denuit, and Montserrat Guillén.  Models of Insurance Claim Counts with Time Dependence Based on Generalization of Poisson and Negative Binomial Distributions.  _Variance_ , 2(1), 2008.
  * Brooks [1998] Stephen P. Brooks.  Markov Chain Monte Carlo Method and Its Application.  _Journal of the Royal Statistical Society. Series D (The Statistician)_ , 47(1):69–100, 1998.  ISSN 0039-0526.
  * Cameron and Johansson [1997] A. Colin Cameron and Per Johansson.  Count Data Regression Using Series Expansions: With Applications.  _Journal of Applied Econometrics_ , 12(3):203–223, 1997.  ISSN 0883-7252.
  * Cameron and Windmeijer [1996] A. Colin Cameron and Frank A. G. Windmeijer.  R-Squared Measures for Count Data Regression Models with Applications to Health-Care Utilization.  _Journal of Business & Economic Statistics_, 14(2):209–220, 1996.  ISSN 0735-0015.  doi: 10.2307/1392433.
  * Carpenter et al. [2017] Bob Carpenter, Andrew Gelman, Matthew D. Hoffman, Daniel Lee, Ben Goodrich, Michael Betancourt, Marcus Brubaker, Jiqiang Guo, Peter Li, and Allen Riddell.  Stan: A Probabilistic Programming Language.  _Journal of Statistical Software_ , 76(1), 2017.  ISSN 1548-7660.  doi: 10.18637/jss.v076.i01.
  * Carvalho et al. [2009] Carlos M. Carvalho, Nicholas G. Polson, and James G. Scott.  Handling Sparsity via the Horseshoe.  In _Proceedings of the Twelth International Conference on Artificial Intelligence and Statistics_ , pages 73–80. PMLR, April 2009\.
  * Charpentier [2015] Arthur Charpentier.  _Computational Actuarial Science with R_.  1st edition, 2015.  ISBN 978-1-138-03378-8.
  * Chib [1998] Siddhartha Chib.  Estimation and comparison of multiple change-point models.  _Journal of Econometrics_ , 86(2):221–241, October 1998.  ISSN 0304-4076.  doi: 10.1016/S0304-4076(97)00115-2.
  * Chin and McNulty [2023] Stephanie Chin and Greg McNulty.  A Bayesian Model for Estimating Long Tailed Excess of Loss Reinsurance Loss Costs.  _CAS E-Forum_ , Summer, August 2023.
  * Cressie and Chan [1989] Noel Cressie and Ngai H. Chan.  Spatial Modeling of Regional Variables.  _Journal of the American Statistical Association_ , 84(406):393–401, 1989.  ISSN 0162-1459.  doi: 10.2307/2289922.
  * de Azevêdo et al. [2023] Beatriz Margarida Zanotto de Azevêdo, Luiz Augusto Finger França Maluf, Joelson Oliveira Sampaio, and Vinicius Augusto Brunassi Silva.  Estimating vehicle types effects on auto insurance premiums in São Paulo City: A GAMLSS Approach.  _CONTRIBUCIONES A LAS CIENCIAS SOCIALES_ , 16(10):21008–21044, October 2023.  ISSN 1988-7833.  doi: 10.55905/revconv.16n.10-140.
  * Ding et al. [2020] Kexing Ding, Baruch Lev, Xuan Peng, Ting Sun, and Miklos A. Vasarhelyi.  Machine learning improves accounting estimates: Evidence from insurance payments.  _Review of Accounting Studies_ , 25(3):1098–1134, September 2020.  ISSN 1573-7136.  doi: 10.1007/s11142-020-09546-9.
  * Duane et al. [1987] Simon Duane, A.D. Kennedy, Brian J. Pendleton, and Duncan Roweth.  Hybrid Monte Carlo.  _Physics Letters B_ , 195(2), September 1987.
  * Dutang and Charpentier [2022] Christophe Dutang and Arthur Charpentier.  CASdatasets, 2022.
  * Fick and Hijmans [2017] Stephen E. Fick and Robert J. Hijmans.  WorldClim 2: New 1-km spatial resolution climate surfaces for global land areas.  _International Journal of Climatology_ , 37(12):4302–4315, October 2017.  ISSN 08998418.  doi: 10.1002/joc.5086.
  * Frees and Valdez [2008] Edward W. Frees and Emiliano A. Valdez.  Hierarchical Insurance Claims Modeling.  _Journal of the American Statistical Association_ , 103(484):1457–1469, 2008.  ISSN 0162-1459.
  * Frostig et al. [2018] Roy Frostig, Matthew James Johnson, and Chris Leary.  Compiling machine learning programs via high-level tracing.  In _SysML_ , Stanford, CA, February 2018.
  * Gelfand and Vounatsou [2003] Alan E. Gelfand and Penelope Vounatsou.  Proper multivariate conditional autoregressive models for spatial data analysis.  _Biostatistics (Oxford, England)_ , 4(1):11–25, January 2003.  ISSN 1465-4644.  doi: 10.1093/biostatistics/4.1.11.
  * Gelman [2006] Andrew Gelman.  Prior distributions for variance parameters in hierarchical models (comment on article by Browne and Draper).  _Bayesian Analysis_ , 1(3):515–534, September 2006.  ISSN 1936-0975, 1931-6690.  doi: 10.1214/06-BA117A.
  * Gelman and Carpenter [2020] Andrew Gelman and Bob Carpenter.  Bayesian Analysis of Tests with Unknown Specificity and Sensitivity.  _Journal of the Royal Statistical Society. Series C, Applied Statistics_ , 69(5):1269–1283, August 2020.  ISSN 0035-9254.  doi: 10.1111/rssc.12435.
  * Gelman and Rubin [1992] Andrew Gelman and Donald B. Rubin.  Inference from Iterative Simulation Using Multiple Sequences.  _Statistical Science_ , 7(4), 1992.
  * Geyer [1992] Charles J. Geyer.  Practical Markov Chain Monte Carlo.  _Statistical Science_ , 7(4):473–483, 1992.  ISSN 0883-4237.
  * Gorelick et al. [2017] Noel Gorelick, Matt Hancher, Mike Dixon, Simon Ilyushchenko, David Thau, and Rebecca Moore.  Google Earth Engine: Planetary-scale geospatial analysis for everyone.  _Remote Sensing of Environment_ , 202:18–27, December 2017\.  ISSN 0034-4257.  doi: 10.1016/j.rse.2017.06.031.
  * Guikema and Coffelt [2008] Seth D. Guikema and Jeremy P. Coffelt.  A flexible count data regression model for risk analysis.  _Risk Analysis: An Official Publication of the Society for Risk Analysis_ , 28(1):213–223, February 2008.  ISSN 1539-6924.  doi: 10.1111/j.1539-6924.2008.01014.x.
  * Hastings [1970] W. K. Hastings.  Monte Carlo Sampling Methods Using Markov Chains and Their Applications.  _Biometrika_ , 57(1):97–109, 1970.  ISSN 0006-3444.  doi: 10.2307/2334940.
  * Hoffman and Gelman [2014] Matthew D Hoffman and Andrew Gelman.  The No-U-Turn Sampler: Adaptively Setting Path Lengths in Hamiltonian Monte Carlo.  _Journal of Machine Learning Research_ , page 31, 2014.
  * Hoffman et al. [2013] Matthew D Hoffman, David Blei, Chong Wang, and John Paisley.  Stochastic Variational Inference.  _Journal of Machine Learning Research_ , 2013.
  * Huete et al. [2002] A Huete, K Didan, T Miura, E. P Rodriguez, X Gao, and L. G Ferreira.  Overview of the radiometric and biophysical performance of the MODIS vegetation indices.  _Remote Sensing of Environment_ , 83(1):195–213, November 2002.  ISSN 0034-4257.  doi: 10.1016/S0034-4257(02)00096-2.
  * Jia et al. [2023] Yisu Jia, Stefanos Kechagias, James Livsey, Robert Lund, and Vladas Pipiras.  Latent Gaussian Count Time Series.  _Journal of the American Statistical Association_ , 118(541):596–606, January 2023.  ISSN 0162-1459.  doi: 10.1080/01621459.2021.1944874.
  * Jin et al. [2005] Xiaoping Jin, Bradley P. Carlin, and Sudipto Banerjee.  Generalized Hierarchical Multivariate CAR Models for Areal Data.  _Biometrics_ , 61(4):950–961, December 2005.  ISSN 0006341X.  doi: 10.1111/j.1541-0420.2005.00359.x.
  * Krapu and Borsuk [2019] Christopher Krapu and Mark Borsuk.  Probabilistic programming: A review for environmental modellers.  _Environmental Modelling & Software_, 114:40–48, April 2019.  ISSN 1364-8152.  doi: 10.1016/j.envsoft.2019.01.014.
  * Kuleshov et al. [2018] Volodymyr Kuleshov, Nathan Fenner, and Stefano Ermon.  Accurate Uncertainties for Deep Learning Using Calibrated Regression.  In _Proceedings of the 35th International Conference on Machine Learning_ , pages 2796–2804. PMLR, July 2018.
  * Lunn et al. [2009] David Lunn, David Spiegelhalter, Andrew Thomas, and Nicky Best.  The BUGS project: Evolution, critique and future directions.  _Statistics in Medicine_ , 28(25):3049–3067, 2009\.  ISSN 1097-0258.  doi: 10.1002/sim.3680.
  * Lupton [2022] Daniel Lupton.  Machine Learning in Insurance.  2022\.
  * Martin et al. [2023] Gael M. Martin, David T. Frazier, and Christian P. Robert.  Approximating Bayes in the 21st Century.  _Statistical Science_ , -1(-1), January 2023.  ISSN 0883-4237.  doi: 10.1214/22-STS875.
  * Miralles et al. [2023] Ophélia Miralles, Anthony C. Davison, and Timo Schmid.  Bayesian modeling of insurance claims for hail damage, August 2023.
  * Neal [2011] Radford Neal.  MCMC Using Hamiltonian Dynamics.  In Steve Brooks, Andrew Gelman, Galin Jones, and Xiao-Li Meng, editors, _Handbook of Markov Chain Monte Carlo_ , volume 20116022. Chapman and Hall/CRC, May 2011.  ISBN 978-1-4200-7941-8 978-1-4200-7942-5.  doi: 10.1201/b10905-6.
  * Paszke et al. [2019] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zach DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.  PyTorch: An Imperative Style, High-Performance Deep Learning Library.  https://arxiv.org/abs/1912.01703v1, December 2019.
  * Paulson et al. [2010] Nicholas D. Paulson, Chad E. Hart, and Dermot J. Hayes.  A spatial Bayesian approach to weather derivatives.  _Agricultural Finance Review_ , 70(1):79–96, May 2010.  ISSN 0002-1466.  doi: 10.1108/00021461011042657.
  * Pekel et al. [2016] Jean-François Pekel, Andrew Cottam, Noel Gorelick, and Alan S. Belward.  High-resolution mapping of global surface water and its long-term changes.  _Nature_ , 540(7633):418–422, December 2016.  ISSN 0028-0836, 1476-4687.  doi: 10.1038/nature20584.
  * Peres et al. [2019] Vivileine Maria Peres, Wilfredo L. Maldonado, and Osvaldo Candido.  Automobile insurance in Brazil: Market concentration and demand.  _Revista Contabilidade & Finanças_, 30:396–408, May 2019.  ISSN 1519-7077, 1808-057X.  doi: 10.1590/1808-057x201808300.
  * Phan et al. [2019] Du Phan, Neeraj Pradhan, and Martin Jankowiak.  Composable Effects for Flexible and Accelerated Probabilistic Programming in NumPyro.  _arXiv:1912.11554 [cs, stat]_ , December 2019.
  * Roberts and Rosenthal [1998] Gareth O. Roberts and Jeffrey S. Rosenthal.  Optimal scaling of discrete approximations to Langevin diffusions.  _Journal of the Royal Statistical Society: Series B (Statistical Methodology)_ , 60(1):255–268, February 1998.  ISSN 1369-7412, 1467-9868.  doi: 10.1111/1467-9868.00123.
  * Roberts and Tweedie [1996] Gareth O. Roberts and Richard L. Tweedie.  Exponential convergence of Langevin distributions and their discrete approximations.  _Bernoulli_ , 2(4):341–363, December 1996.  ISSN 1350-7265.
  * Rue et al. [2016] Håvard Rue, Andrea Riebler, Sigrunn H. Sørbye, Janine B. Illian, Daniel P. Simpson, and Finn K. Lindgren.  Bayesian Computing with INLA: A Review, September 2016.
  * Salvatier et al. [2016] John Salvatier, Thomas V. Wiecki, and Christopher Fonnesbeck.  Probabilistic programming in Python using PyMC3.  _PeerJ Computer Science_ , 2:e55, April 2016.  ISSN 2376-5992.  doi: 10.7717/peerj-cs.55.
  * Sherden [1984] William A. Sherden.  An Analysis of the Determinants of the Demand for Automobile Insurance.  _The Journal of Risk and Insurance_ , 51(1):49–62, 1984.  ISSN 0022-4367.  doi: 10.2307/252800.
  * Shimada et al. [2014] Masanobu Shimada, Takuya Itoh, Takeshi Motooka, Manabu Watanabe, Tomohiro Shiraishi, Rajesh Thapa, and Richard Lucas.  New global forest/non-forest maps from ALOS PALSAR data (2007–2010).  _Remote Sensing of Environment_ , 155:13–31, December 2014\.  ISSN 0034-4257.  doi: 10.1016/j.rse.2014.04.014.
  * Simpson et al. [2015] Daniel P. Simpson, Håvard Rue, Thiago G. Martins, Andrea Riebler, and Sigrunn H. Sørbye.  Penalising model component complexity: A principled, practical approach to constructing priors, August 2015.
  * Tatem [2017] Andrew J. Tatem.  WorldPop, open data for spatial demography.  _Scientific Data_ , 4(1):170004, January 2017\.  ISSN 2052-4463.  doi: 10.1038/sdata.2017.4.
  * Taylor and Diggle [2014] Benjamin M. Taylor and Peter J. Diggle.  INLA or MCMC? A tutorial and comparative evaluation for spatial prediction in log-Gaussian Cox processes.  _Journal of Statistical Computation and Simulation_ , 84(10):2266–2284, October 2014.  ISSN 0094-9655.  doi: 10.1080/00949655.2013.788653.
  * Theobald et al. [2015] David M. Theobald, Dylan Harrison-Atlas, William B. Monahan, and Christine M. Albano.  Ecologically-Relevant Maps of Landforms and Physiographic Diversity for Climate Adaptation Planning.  _PLOS ONE_ , 10(12):e0143619, December 2015.  ISSN 1932-6203.  doi: 10.1371/journal.pone.0143619.
  * Virtanen et al. [2020] Pauli Virtanen, Ralf Gommers, Travis E. Oliphant, Matt Haberland, Tyler Reddy, David Cournapeau, Evgeni Burovski, Pearu Peterson, Warren Weckesser, Jonathan Bright, Stéfan J. van der Walt, Matthew Brett, Joshua Wilson, K. Jarrod Millman, Nikolay Mayorov, Andrew R. J. Nelson, Eric Jones, Robert Kern, Eric Larson, C. J. Carey, İlhan Polat, Yu Feng, Eric W. Moore, Jake VanderPlas, Denis Laxalde, Josef Perktold, Robert Cimrman, Ian Henriksen, E. A. Quintero, Charles R. Harris, Anne M. Archibald, Antônio H. Ribeiro, Fabian Pedregosa, and Paul van Mulbregt.  SciPy 1.0: Fundamental algorithms for scientific computing in Python.  _Nature Methods_ , 17(3):261–272, March 2020\.  ISSN 1548-7105.  doi: 10.1038/s41592-019-0686-2.
  * Wahl et al. [2022] Jens Christian Wahl, Fredrik Lohne Aanes, Kjersti Aas, Sindre Froyn, and Daniel Piacek.  Spatial modelling of risk premiums for water damage insurance.  _Scandinavian Actuarial Journal_ , 2022(3):216–233, March 2022.  ISSN 0346-1238.  doi: 10.1080/03461238.2021.1951346.
  * Weiss et al. [2018] D. J. Weiss, A. Nelson, H. S. Gibson, W. Temperley, S. Peedell, A. Lieber, M. Hancher, E. Poyart, S. Belchior, N. Fullman, B. Mappin, U. Dalrymple, J. Rozier, T. C. D. Lucas, R. E. Howes, L. S. Tusting, S. Y. Kang, E. Cameron, D. Bisanzio, K. E. Battle, S. Bhatt, and P. W. Gething.  A global map of travel time to cities to assess inequalities in accessibility in 2015.  _Nature_ , 553(7688):333–336, January 2018.  ISSN 1476-4687.  doi: 10.1038/nature25181.
  * Zhang [2017] Yanwei Zhang.  Bayesian Analysis of Big Data in Insurance Predictive Modeling Using Distributed Computing, July 2017.


#  LaTeX Author Guidelines for CVPR Proceedings

First Author
Institution1
Institution1 address
<EMAIL_ADDRESS>Second Author
Institution2
First line of institution2 address
<EMAIL_ADDRESS>
###### Abstract

The ABSTRACT is to be in fully justified italicized text, at the top of the
left-hand column, below the author and affiliation information. Use the word
“Abstract” as the title, in 12-point Times, boldface type, centered relative
to the column, initially capitalized. The abstract is to be in 10-point,
single-spaced type. Leave two blank lines after the Abstract, then begin the
main text. Look at previous CVPR abstracts to get a feel for style and length.

##  1 Introduction

Please follow the steps outlined below when submitting your manuscript to the
IEEE Computer Society Press. This style guide now has several important
modifications (for example, you are no longer warned against the use of sticky
tape to attach your artwork to the paper), so all authors should read this new
version.

###  1.1 Language

All manuscripts must be in English.

###  1.2 Dual submission

Please refer to the author guidelines on the CVPR 2022 web page for a
discussion of the policy on dual submissions.

###  1.3 Paper length

Papers, excluding the references section, must be no longer than eight pages
in length. The references section will not be included in the page count, and
there is no limit on the length of the references section. For example, a
paper of eight pages with two pages of references would have a total length of
10 pages. There will be no extra page charges for CVPR 2022.

Overlength papers will simply not be reviewed. This includes papers where the
margins and formatting are deemed to have been significantly altered from
those laid down by this style guide. Note that this LaTeX guide already sets
figure captions and references in a smaller font. The reason such papers will
not be reviewed is that there is no provision for supervised revisions of
manuscripts. The reviewing process cannot determine the suitability of the
paper for presentation in eight pages if it is reviewed in eleven.

###  1.4 The ruler

The LaTeX style defines a printed ruler which should be present in the version
submitted for review. The ruler is provided in order that reviewers may
comment on particular lines in the paper without circumlocution. If you are
preparing a document using a non-LaTeX document preparation system, please
arrange for an equivalent ruler to appear on the final output pages. The
presence or absence of the ruler should not change the appearance of any other
content on the page. The camera-ready copy should not contain a ruler. (LaTeX
users may use options of cvpr.sty to switch between different versions.)

Reviewers: note that the ruler measurements do not align well with lines in
the paper — this turns out to be very difficult to do well when the paper
contains many figures and equations, and, when done, looks ugly. Just use
fractional references (_e.g_., this line is $087.5$), although in most cases
one would expect that the approximate location will be adequate.

###  1.5 Paper ID

Make sure that the Paper ID from the submission system is visible in the
version submitted for review (replacing the “*****” you see in this document).
If you are using the LaTeX template, make sure to update paper ID in the
appropriate place in the tex file.

###  1.6 Mathematics

Please number all of your sections and displayed equations as in these
examples:

$E=m\cdot c^{2}$ (1)

and

$v=a\cdot t.$ (2)

It is important for readers to be able to refer to any particular equation.
Just because you did not refer to it in the text does not mean some future
reader might not need to refer to it. It is cumbersome to have to use
circumlocutions like “the equation second from the top of page 3 column 1”.
(Note that the ruler will not be present in the final copy, so is not an
alternative to equation numbers). All authors will benefit from reading
Mermin’s description of how to write mathematics:
http://www.pamitc.org/documents/mermin.pdf.

###  1.7 Blind review

Many authors misunderstand the concept of anonymizing for blind review. Blind
review does not mean that one must remove citations to one’s own work—in fact
it is often impossible to review a paper unless the previous citations are
known and available.

Blind review means that you do not use the words “my” or “our” when citing
previous work. That is all. (But see below for tech reports.)

Saying “this builds on the work of Lucy Smith [1]” does not say that you are
Lucy Smith; it says that you are building on her work. If you are Smith and
Jones, do not say “as we show in [7]”, say “as Smith and Jones show in [7]”
and at the end of the paper, include reference 7 as you would any other cited
work.

An example of a bad paper just asking to be rejected:

> An analysis of the frobnicatable foo filter.
>
> In this paper we present a performance analysis of our previous paper [1],
> and show it to be inferior to all previously known methods. Why the previous
> paper was accepted without this analysis is beyond me.
>
> [1] Removed for blind review

An example of an acceptable paper:

> An analysis of the frobnicatable foo filter.
>
> In this paper we present a performance analysis of the paper of Smith _et
> al_. [1], and show it to be inferior to all previously known methods. Why
> the previous paper was accepted without this analysis is beyond me.
>
> [1] Smith, L and Jones, C. “The frobnicatable foo filter, a fundamental
> contribution to human knowledge”. Nature 381(12), 1-213.

If you are making a submission to another conference at the same time, which
covers similar or overlapping material, you may need to refer to that
submission in order to explain the differences, just as you would if you had
previously published related work. In such cases, include the anonymized
parallel submission [Authors14] as supplemental material and cite it as

> [1] Authors. “The frobnicatable foo filter”, F&G 2014 Submission ID 324,
> Supplied as supplemental material fg324.pdf.

Finally, you may feel you need to tell the reader that more details can be
found elsewhere, and refer them to a technical report. For conference
submissions, the paper must stand on its own, and not require the reviewer to
go to a tech report for further details. Thus, you may say in the body of the
paper “further details may be found in [Authors14b]”. Then submit the tech
report as supplemental material. Again, you may not assume the reviewers will
read this material.

Sometimes your paper is about a problem which you tested using a tool that is
widely known to be restricted to a single institution. For example, let’s say
it’s 1969, you have solved a key problem on the Apollo lander, and you believe
that the CVPR70 audience would like to hear about your solution. The work is a
development of your celebrated 1968 paper entitled “Zero-g frobnication: How
being the only people in the world with access to the Apollo lander source
code makes us a wow at parties”, by Zeus _et al_.

You can handle this paper like any other. Do not write “We show how to improve
our previous work [Anonymous, 1968]. This time we tested the algorithm on a
lunar lander [name of lander removed for blind review]”. That would be silly,
and would immediately identify the authors. Instead write the following:

> We describe a system for zero-g frobnication. This system is new because it
> handles the following cases: A, B. Previous systems [Zeus et al. 1968] did
> not handle case B properly. Ours handles it by including a foo term in the
> bar integral.
>
> …
>
> The proposed system was integrated with the Apollo lunar lander, and went
> all the way to the moon, don’t you know. It displayed the following
> behaviours, which show how well we solved cases A and B: …

As you can see, the above text follows standard scientific convention, reads
better than the first version, and does not explicitly name you as the
authors. A reviewer might think it likely that the new paper was written by
Zeus _et al_., but cannot make any decision based on that guess. He or she
would have to be sure that no other authors could have been contracted to
solve problem B.

FAQ
Q: Are acknowledgements OK?
A: No. Leave them for the final copy.
Q: How do I cite my results reported in open challenges? A: To conform with
the double-blind review policy, you can report results of other challenge
participants together with your results in your paper. For your results,
however, you should not identify yourself and should not mention your
participation in the challenge. Instead present your results referring to the
method proposed in your paper and draw conclusions based on the experimental
comparison to other results.

Figure 1: Example of caption. It is set in Roman so that mathematics (always
set in Roman: $B\sin A=A\sin B$) may be included without an ugly clash.

###  1.8 Miscellaneous

Compare the following:
`$conf_a$` $conf_{a}$ `$\mathit{conf}_a$` $\mathit{conf}_{a}$
See The TeXbook, p165.

The space after _e.g_., meaning “for example”, should not be a sentence-ending
space. So _e.g_. is correct, e.g. is not. The provided `\eg` macro takes care
of this.

When citing a multi-author paper, you may save space by using “et alia”,
shortened to “ _et al_.” (not “et. al.” as “et” is a complete word). If you
use the `\etal` macro provided, then you need not worry about double periods
when used at the end of a sentence as in Alpher _et al_. However, use it only
when there are three or more authors. Thus, the following is correct:
“Frobnication has been trendy lately. It was introduced by Alpher [Alpher02],
and subsequently developed by Alpher and Fotheringham-Smythe [Alpher03], and
Alpher _et al_. [Alpher04].”

This is incorrect: “… subsequently developed by Alpher _et al_. [Alpher03] …”
because reference [Alpher03] has just two authors.

(a) An example of a subfigure.

(b) Another example of a subfigure.

Figure 2: Example of a short caption, which should be centered.

##  2 Formatting your paper

All text must be in a two-column format. The total allowable size of the text
area is $6\frac{7}{8}$ inches (17.46 cm) wide by $8\frac{7}{8}$ inches (22.54
cm) high. Columns are to be $3\frac{1}{4}$ inches (8.25 cm) wide, with a
$\frac{5}{16}$ inch (0.8 cm) space between them. The main title (on the first
page) should begin 1 inch (2.54 cm) from the top edge of the page. The second
and following pages should begin 1 inch (2.54 cm) from the top edge. On all
pages, the bottom margin should be $1\frac{1}{8}$ inches (2.86 cm) from the
bottom edge of the page for $8.5\times 11$-inch paper; for A4 paper,
approximately $1\frac{5}{8}$ inches (4.13 cm) from the bottom edge of the
page.

###  2.1 Margins and page numbering

All printed material, including text, illustrations, and charts, must be kept
within a print area $6\frac{7}{8}$ inches (17.46 cm) wide by $8\frac{7}{8}$
inches (22.54 cm) high. Page numbers should be in the footer, centered and
$\frac{3}{4}$ inches from the bottom of the page. The review version should
have page numbers, yet the final version submitted as camera ready should not
show any page numbers. The LaTeX template takes care of this when used
properly.

###  2.2 Type style and fonts

Wherever Times is specified, Times Roman may also be used. If neither is
available on your word processor, please use the font closest in appearance to
Times to which you have access.

MAIN TITLE. Center the title $1\frac{3}{8}$ inches (3.49 cm) from the top edge
of the first page. The title should be in Times 14-point, boldface type.
Capitalize the first letter of nouns, pronouns, verbs, adjectives, and
adverbs; do not capitalize articles, coordinate conjunctions, or prepositions
(unless the title begins with such a word). Leave two blank lines after the
title.

AUTHOR NAME(s) and AFFILIATION(s) are to be centered beneath the title and
printed in Times 12-point, non-boldface type. This information is to be
followed by two blank lines.

The ABSTRACT and MAIN TEXT are to be in a two-column format.

MAIN TEXT. Type main text in 10-point Times, single-spaced. Do NOT use double-
spacing. All paragraphs should be indented 1 pica (approx. $\frac{1}{6}$ inch
or 0.422 cm). Make sure your text is fully justified—that is, flush left and
flush right. Please do not place any additional blank lines between
paragraphs.

Figure and table captions should be 9-point Roman type as in Figs. 1 and 2.
Short captions should be centred.

Callouts should be 9-point Helvetica, non-boldface type. Initially capitalize
only the first word of section titles and first-, second-, and third-order
headings.

FIRST-ORDER HEADINGS. (For example, 1\. Introduction) should be Times 12-point
boldface, initially capitalized, flush left, with one blank line before, and
one blank line after.

SECOND-ORDER HEADINGS. (For example, 1.1. Database elements) should be Times
11-point boldface, initially capitalized, flush left, with one blank line
before, and one after. If you require a third-order heading (we discourage
it), use 10-point Times, boldface, initially capitalized, flush left, preceded
by one blank line, followed by a period and your text on the same line.

###  2.3 Footnotes

Please use footnotes111This is what a footnote looks like. It often distracts
the reader from the main flow of the argument. sparingly. Indeed, try to avoid
footnotes altogether and include necessary peripheral observations in the text
(within parentheses, if you prefer, as in this sentence). If you wish to use a
footnote, place it at the bottom of the column on the page on which it is
referenced. Use Times 8-point type, single-spaced.

###  2.4 Cross-references

For the benefit of author(s) and readers, please use the

      \cref{...}

command for cross-referencing to figures, tables, equations, or sections. This
will automatically insert the appropriate label alongside the cross-reference
as in this example:

> To see how our method outperforms previous work, please see Fig. 1 and Tab.
> 1. It is also possible to refer to multiple targets as once, _e.g_. to Figs.
> 1 and 2(a). You may also return to Sec. 2 or look at Eq. 2.

If you do not wish to abbreviate the label, for example at the beginning of
the sentence, you can use the

      \Cref{...}

command. Here is an example:

> Figure 1 is also quite important.

###  2.5 References

List and number all bibliographical references in 9-point Times, single-
spaced, at the end of your paper. When referenced in the text, enclose the
citation number in square brackets, for example [Authors14]. Where
appropriate, include page numbers and the name(s) of editors of referenced
books. When you cite multiple papers at once, please make sure that you cite
them in numerical order like this [Alpher02, Alpher03, Alpher05, Authors14b,
Authors14]. If you use the template as advised, this will be taken care of
automatically.

Method | Frobnability
---|---
Theirs | Frumpy
Yours | Frobbly
Ours | Makes one’s heart Frob
Table 1: Results. Ours is better.

###  2.6 Illustrations, graphs, and photographs

All graphics should be centered. In LaTeX, avoid using the center environment
for this purpose, as this adds potentially unwanted whitespace. Instead use

      \centering

at the beginning of your figure. Please ensure that any point you wish to make
is resolvable in a printed copy of the paper. Resize fonts in figures to match
the font in the body text, and choose line widths that render effectively in
print. Readers (and reviewers), even of an electronic copy, may choose to
print your paper in order to read it. You cannot insist that they do
otherwise, and therefore must not assume that they can zoom in to see tiny
details on a graphic.

When placing figures in LaTeX, it’s almost always best to use
`\includegraphics`, and to specify the figure width as a multiple of the line
width as in the example below

       \usepackage{graphicx} ...
       \includegraphics[width=0.8\linewidth]
                       {myfile.pdf}

###  2.7 Color

Please refer to the author guidelines on the CVPR 2022 web page for a
discussion of the use of color in your document.

If you use color in your plots, please keep in mind that a significant subset
of reviewers and readers may have a color vision deficiency; red-green
blindness is the most frequent kind. Hence avoid relying only on color as the
discriminative feature in plots (such as red _vs_. green lines), but add a
second discriminative feature to ease disambiguation.

##  3 Final copy

You must include your signed IEEE copyright release form when you submit your
finished paper. We MUST have this form before your paper can be published in
the proceedings.

Please direct any questions to the production editor in charge of these
proceedings at the IEEE Computer Society Press:
https://www.computer.org/about/contact.