Addr[arabic] ANote[fnsymbol]

# Estimating the timing of stillbirths in countries worldwide using a Bayesian
hierarchical penalized splines regression model

Michael Y.C. Chong  Monica Alexander

###### ?abstractname?

Reducing the global burden of stillbirths is important to improving child and
maternal health. Of interest is understanding patterns in the timing of
stillbirths — that is, whether they occur in the intra- or antepartum period —
because stillbirths that occur intrapartum are largely preventable. However,
data availability on the timing of stillbirths is highly variable across the
world, with low- and middle-income countries generally having few reliable
observations. In this paper we develop a Bayesian penalized splines regression
framework to estimate the proportion of stillbirths that are intrapartum for
all countries worldwide. The model accounts for known relationships with
neonatal mortality, pools information across geographic regions, incorporates
different errors based on data attributes, and allows for data-driven temporal
trends. A weighting procedure is proposed to account for unrepresentative
subnational data. Results suggest that the intrapartum proportion is generally
decreasing over time, but progress is slower in some regions, particularly
Sub-Saharan Africa.

<EMAIL_ADDRESS>The authors would like to acknowledge Lucia Hug,
Danzhen You, David Sharrow, Lucy Smith, Hannah Blencowe, Leontine Alkema,
Zhengfan Wang, and Jon Wakefield for their invaluble input and feedback.

##  1 Introduction

Stillbirths represent a large share of the global burden on child and maternal
health. Globally, in 2019 there were an estimated 2 million stillbirths
(babies born with no sign of life at 28 weeks of pregnancy or later) (Hug et
al. 2021). This represents a rate of around 14 stillbirths per 1000 total
births, which is of comparable magnitude to the global neonatal mortality rate
of around 17 deaths per 1000 live births (Hug et al. 2019). As such, there
have been increased calls to monitor stillbirth outcomes in addition to infant
and maternal deaths to fully understand the extent of the risks faced by
different populations during pregnancy and childbirth (Lawn et al. 2016). The
Every Newborn Action Plan (ENAP), which was endorsed by 194 WHO Member States,
calls for each country to achieve a rate of 12 stillbirths or fewer per 1000
total births by 2030 (World Health Organization 2014). While progress has been
made in reducing stillbirths, substantial inequalities across countries and
regions still persist (Hug et al. 2020).

An important part of working towards the goal of ending preventable
stillbirths is having reliable estimates of when stillbirths occur, that is,
whether a stillbirth occurs before or after the onset of labor (ante- versus
intrapartum). Stillbirths that occur intrapartum are likely to be preventable
in most situations given adequate access to healthcare (De Bernis et al.
2016), and thus are likely to be responsive to increased resources and health
policy changes. The goal of this project is thus to estimate the proportion of
stillbirths that are intrapartum for all countries and regions worldwide, over
the period 2000-2021.

Challenges in estimation exist, particularly in low- and middle-income
countries, due to data availability issues. Even when data do exist, the
quality of reporting and type of population captured varies substantially by
data source type and diagnosis method. An intrapartum stillbirth is defined as
a fetal death occurring after the onset of labor and prior to delivery. The
gold-standard classification of whether a fetus is alive after the onset of
labor is the presence of a fetal heartbeat. However, in many contexts, the
presence of a fetal heartbeat may not be documented and so diagnosis of the
presence of life occurs through other methods postpartum; for example, fetuses
who died antepartum can have skin changes consistent with maceration, tissue
injury, meconium staining, and edema (Da Silva et al. 2016). However, these
types of diagnosis methods are much more prone to misclassification compared
to heartbeat diagnoses. In addition, different countries use different
definitions of what constitutes a stillbirth, based on gestational age or
birthweight thresholds (or a mixture of both). For international
comparability, the WHO recommends using the cut-off of 1000g or more at birth
(if available), or after 28 completed weeks of gestation (Hug et al. 2020).
But different countries report stillbirths at different definitions, and
definitions may vary subnationally by jurisdiction; for example in the United
States, there are eight different definitions by combinations of gestational
age and weight (Da Silva et al. 2016). In particular, many countries report
stillbirths at 22 weeks gestation. An additional issue stems from differences
in data collection systems. While most high-income countries collect
information on stillbirths through civil registration and vital statistics
(CRVS) systems, many other countries rely on other collection systems that are
less reliable and cover only a small portion of the population, and some
countries have no information available. As such, a method to estimate the
proportion of stillbirths that are intrapartum needs to account for a wide
range of data reporting, coverage, and classification issues. However,
previous regional estimates of intrapartum stillbirths rely on median
proportion-based approaches, making no data adjustments and not accounting for
uncertainty in data, estimates, or projections (Hug et al. 2020).

To overcome these issues, we formulate a Bayesian hierarchical penalized
splines model that takes into consideration different types of data on
stillbirth timing from a wide range of reporting and diagnosis systems. The
model captures the underlying relationship between changes in the proportion
of stillbirths that are intrapartum and the neonatal mortality rate over time,
which allows for estimates and projections of the intrapartum proportion to be
made even in contexts where available data are limited. Additionally, the
model allows for data-driven trends through the use of a penalized splines
regression framework. The model is also informed by high-quality auxiliary
data to adjust for different gestational age definitions, and allows for
varying amounts of uncertainty around data from different sources. We also
propose a post-estimation weighting scheme to account for varying levels of
coverage in the observed data.

The remainder of the paper is structured as follows. Firstly, we describe data
sources and availability, before introducing the modeling framework and
presenting some results and validation. We conclude with a discussion of
possible extensions to the model.

##  2 Data

Data on the number of stillbirths by timing (intrapartum versus antepartum)
can come from a number of different sources. Data from civil registration and
vital statistics (CRVS) systems generally has national coverage and are
assumed to be relatively high quality compared to other sources. Data on
stillbirth timing also come from health management information systems (HMIS)
(such as DHIS2), which are commonly used in low- and middle-income countries
to record information from a number of health facilities (Moxon et al. 2015).
Subnational data are sourced from either health facility data or population-
based studies, most notably through the Global Network Study, which collects
data on stillbirths and neonatal deaths in multiple communities across several
different countries (Frøen et al. 2016). Data were extracted and collated
through a number of channels, including web-based searches of national
statistics’ offices, UNICEF country consultations, and literature searches.

Overall, a least one observation of the number of stillbirths that occurred
intrapartum and antepartum was available for 92 countries across the period
2000–2020. Table 1 shows the breakdown of data availability by Sustainable
Development Goal (SDG) region, suggesting that the most observations are
available in the high-income country group.

Table 1: Data availability by SDG region SDG region | Observations | Countries | Country-years
---|---|---|---
Central and Southern Asia | 163 | 7 | 65
Eastern and South-Eastern Asia | 57 | 8 | 53
Latin America and the Caribbean | 272 | 13 | 171
North America, Europe, Australia and New Zealand | 460 | 30 | 368
Northern Africa and Western Asia | 43 | 8 | 42
Oceania (exc. Australia and New Zealand) | 1 | 1 | 1
Sub-Saharan Africa | 280 | 25 | 158

Looking at the proportion of observations by data collection system (Table 2)
shows marked differences across the SDG regions. In particular, while 98% of
data in North America, Europe, Australia and New Zealand come from CRVS
systems, less than 10% of the data from Central and Southern Asia do, and
almost no data in Sub-Saharan Africa are from CRVS. The majority of data in
Central and Southern Asia are from subnational population-based studies.
Northern Africa and Western Asia have notable shares from both CRVS and HMIS
systems, while data from Sub-Saharan Africa are most from subnational
population-based studies or HMIS.

Table 2: Data source by SDG region SDG region | CRVS | Health facility | Subnat pop-based | HMIS
---|---|---|---|---
Central and Southern Asia | 0.067 | 0.110 | 0.810 | 0.012
Eastern and South-Eastern Asia | 0.684 | 0.281 | 0.035 | 0.000
Latin America and the Caribbean | 0.794 | 0.044 | 0.162 | 0.000
North America, Europe, Australia and New Zealand | 0.980 | 0.015 | 0.004 | 0.000
Northern Africa and Western Asia | 0.465 | 0.047 | 0.093 | 0.395
Oceania (exc. Australia and New Zealand) | 0.000 | 1.000 | 0.000 | 0.000
Sub-Saharan Africa | 0.007 | 0.139 | 0.375 | 0.479

Table 3 shows the proportion of observations by SDG based on gestational age
or birthweight definitions. The early stillbirth definition refers to either
stillbirth at 22 weeks or later, and/or a birthweight of 500g or above. The
late stillbirth definition refers to either stillbirth at 28 weeks or later,
and/or a birthweight of 1000g or above. For the purposes of our estimation
efforts, we are interested in estimating stillbirth timing for late
stillbirths, and so definitional differences need to be adjusted for in the
modeling framework. As shown in Table 3, observations based on the early
definition make up at least 17% of all observations across all regions. In
addition, observations where it is unclear what definition to use, or where
all stillbirths are included, make up a non-negligible proportion in Latin
America and the Caribbean, and Sub-Saharan Africa.

Table 3: Stillbirth definition by SDG region SDG region | early | late | not defined or all
---|---|---|---
Central and Southern Asia | 0.429 | 0.534 | 0.037
Eastern and South-Eastern Asia | 0.228 | 0.772 | 0.000
Latin America and the Caribbean | 0.176 | 0.570 | 0.254
North America, Europe, Australia and New Zealand | 0.465 | 0.526 | 0.009
Northern Africa and Western Asia | 0.628 | 0.326 | 0.047
Oceania (exc. Australia and New Zealand) | 1.000 | 0.000 | 0.000
Sub-Saharan Africa | 0.254 | 0.621 | 0.125

###  2.1 Additional data sources

In addition to data on the number of intrapartum and antepartum stillbirths,
we use several other data sources and estimates to help inform the model and
calculate the eventual proportion of stillbirths that are intrapartum (IPSB)
at various regional and global aggregations.

Firstly, we use estimates of the neonatal mortality rate (NMR), which is the
number of deaths in the neonatal period (first 28 days of life) per 1000 live
births, for every country and year in the estimation period of interest. NMR
estimates are produced by the UN Interagency Group on Mortality Estimation (UN
IGME) as part of annual SDG reporting efforts. Details of the estimates and
estimation can be found in UN IGME (2021). In particular, we use the full
posterior samples of NMR for each country-year in order to reflect and
propagate the uncertainty in the estimation process.

We also draw upon estimates of the total stillbirth rate (SBR), which is
defined as the number of babies born with no sign of life at 28 weeks or more
of gestation, per 1,000 total births. These estimates are used to quantify the
coverage of data sources within a country, and to weight country-level
estimates of IPSB to get regional estimates of IPSB. Similarly to the NMR, we
use UN IGME estimates (Hug et al. 2020) and use the full posterior samples to
better propagate the uncertainty in the SBR estimation process.

As shown in Table 3, data on stillbirths by timing is available at varying
gestational age definitions. In order to inform gestational age adjustments in
the estimation process, we use data from Euro-Peristat (Gissler et al. 2010).
The data provided contain intrapartum and antepartum stillbirth counts in 2015
for 17 high-income European countries corresponding to both types of
gestational age definitions. For confidentiality reasons, country names in
these data are not given. Additionally, to help inform the gestational age
adjustment for lower income countries, we draw on data from the Global Health
Network (Frøen et al. 2016), which gives information on stillbirths by timing
for eight low- and middle-income countries111Argentina, Democratic Republic of
Congo, Zambia, Guatemala, Bangladesh, India, Pakistan, and Kenya. at both
early and late gestational age definitions.

##  3 Model framework

###  3.1 Overview

The varying type, quality, and coverage of data on stillbirths by timing that
are available across different countries suggests the need for an estimation
process that accounts and adjusts for various data characteristics. Firstly,
there is a large degree of missingness in many countries across the time
period of interest (2000-2021), with many countries only having one
observation over the whole period, and more than half of the countries having
no observations at all. As such, a hierarchical model, which allows for
information exchange across countries within regions may be appropriate for
this context. The lack of available data also suggests that one or more
covariates may be useful to obtain reasonable trends over time. Secondly, the
data that do exist come from a wide range of data collection systems, which
are likely to have different levels of coverage of the overall population of
interest, and also measurement error in the classification of intrapartum
versus antepartum timing. This suggests a need for allowing for different
magnitudes of error around observations in a data (or measurement error)
model, and also potentially allow for observations to be ‘down-weighted’ if
they come from a data system that only covers a small portion of total
stillbirths. Finally, data are reported using different definitions (early or
late stillbirths). As we expect early stillbirths to have a higher proportion
of fetal deaths occurring in the intrapartum period (Getahun et al. 2007), the
estimation process should account for these definitional issues.

To address these issues, we propose a Bayesian hierarchical penalized splines
regression model to estimate the IPSB levels in each country and trends over
the period 2000-2021. In this section we describe the model for place-specific
IPSB, the weighting process to obtain country-level estimates, and how
regional estimates are derived.

###  3.2 Model for place-level IPSB

We consider data on intrapartum and antepartum stillbirths as available for a
specific ‘place.’ In many cases, a ‘place’ refers to data for a whole country
(e.g. data from a nationally-representative CRVS system), but in other
contexts, a ‘place’ may refer to a subnational region (such as a state or
province) or a single health facility. In our approach, we explicitly model
the data at the ‘place’ level, and then reweight the estimates in a second
step to obtain country-level estimates.

For data points $i=1,\dots,N$, let $y_{i}$ and $z_{i}$ denote the the number
of observed intrapartum and antepartum stillbirths respectively. Then

$y_{i}|\phi_{i}\sim\text{Binomial}(y_{i}+z_{i},\phi_{i}),$ (1)

where $y_{i}+z_{i}$ represents the total number of classified
stillbirths222Note that for some observations, the number of stillbirths of
unknown timing was also reported. For the purposes of estimation, unknown
timing stillbirths were excluded. and the parameter $\phi_{i}$ represents the
proportion of intrapartum stillbirths. The proportion $\phi_{i}$ is modeled

$\displaystyle\text{logit}(\phi_{i})=\mu_{i}+\varepsilon_{i}$ (2)
$\displaystyle\varepsilon_{i}\sim\text{Normal}(0,\sigma_{\varepsilon,s[i]}^{2})$
(3)

where $\mu_{i}$ describes the “true” inverse-logit transformed proportion for
the population, and $\varepsilon_{i}$ is an error term. The variance
$\sigma_{\varepsilon,s[i]}^{2}$ of the error term depends on the type of data
system of the observation which are expected to vary in reporting quality.
Observations from CRVS systems are expected to be the most reliable, and are
assigned $\sigma_{\varepsilon,s[i]}=0$. Observations from other types of
systems, such as those from subnational health facilities, DHIS/HMIS, or
population studies are expected to be less reliable and for these methods the
variance terms are estimated from data with half-Normal priors:

$\sigma_{\varepsilon,s}\;\begin{cases}=0&\text{if }{s=}\text{ CRVS}\\\
\sim\text{Normal}^{+}(0,1^{2})&\text{if }s=\text{health facility, DHIS/HMIS,
population study}\end{cases}$ (4)

The parameter $\mu_{i}$ represents the mean for the place-level observation
$i$, and is estimated as

$\mu_{i}=\beta_{0}+\beta_{r[i]}+\beta_{c[i]}+\beta_{p[i]}+\beta_{\text{NMR}}\log\text{NMR}_{c[i],t[i]}+\eta_{p[i],t[i]}+\gamma_{g[i],m[i]},$
(5)

where $\beta_{0}$ is a global intercept, $\beta_{r[i]}$ is a region effect for
region $r$, grouped according to the UN Sustainable Development Goals (SDG)
regional groupings, and $\beta_{c[i]}$ is a country effect for country $c$.
These are modeled as

$\displaystyle\beta_{0}\sim\text{Normal}(0,1^{2})$ (6)
$\displaystyle\beta_{r}\sim\text{Normal}(0,\sigma_{\beta_{r}}^{2})$ (7)
$\displaystyle\beta_{c}\sim\text{Normal}(0,\sigma_{\beta_{c}}^{2})$ (8)

where $\sigma_{\beta_{r}}$ and $\sigma_{\beta_{c}}$ are variance terms to be
estimated.

The parameters $\beta_{p[i]}$ represent place-level effects which are specific
to a population or subnational region $p$ as specified by descriptions in the
data. Sub-populations may vary in size depending on data sources within a
country. For instance, one sub-population may be only the births at a specific
health facility, while another may be all births at government facilities. We
assume

$\displaystyle\beta_{p}\sim\text{Normal}(0,\sigma_{\beta_{p}}^{2}),$ (9)

where $\sigma_{\beta_{p}}$ is a variance term to be estimated. For countries
where only one sub-population appears in the data, we set $\beta_{p}=0$ for
identifiability with the country effect $\beta_{c}$.

We use estimates of the (log) NMR as a covariate in the model. This is
motivated by the strong empirical relationship observed between the logit IPSB
and log NMR (as shown in Figure 1), and also by the fact that we expect the
proportion of stillbirths that are intrapartum to be positively correlated
with NMR, as both outcomes are likely to decline in response to improved
medical conditions during birth (Joyce et al. 2004). The term
$\beta_{\text{NMR}}$ represents the slope of the effect and
$\text{NMR}_{c[i],t[i]}$ refers to the UN IGME point estimate of the national
NMR of country $c$ and time $t$ of the observation.

Figure 1: Observed proportion of stillbirths that are intrapartum (IPSB)
(logit scale) versus log of the neonatal mortality rate. Colors represent
different SDG regions.

To allow for data-driven trends in addition to trends based on changes in
country-level NMR, we include a place-time specific component $\eta_{p,t}$,
which is modelled using a first-order penalized splines set up. For each place
$p$, $\eta_{p,t}$ is defined

$\eta_{p,t}=\sum_{h=1}^{H}k_{h}(t)\alpha_{h,p}$ (10)

where $k_{h}(t)$ denotes the $h^{\text{th}}$ spline function evaluated at time
$t$ and $\alpha_{h,p}$ denotes its coefficient to be estimated. We use cubic
B-splines with knots placed at integer year values. First-order differences in
the coefficients, denoted $\Delta_{h,p}$, are penalized to ensure a level of
smoothness in the resulting fit:

$\displaystyle\Delta_{h,p}=\alpha_{h,p}-\alpha_{h-1,p}$ (11)
$\displaystyle\Delta\sim\text{Normal}(0,\sigma_{\Delta}^{2})$ (12)

where $\sigma_{\Delta}$ is an estimated variance term, and the coefficients
are constrained to sum to zero

$\sum_{h=1}^{H}\alpha_{h,p}=0.$ (13)

Finally, $\gamma_{g,w}$ is an adjustment for the gestational age definition
$g$ used in the observation ($g=$ early or $g=$ late). Separate adjustment
factors are estimated for two income groups: $m=$ high income countries (HIC)
and $m=$ low-middle income countries (LMIC). Since the final desired estimate
of the intrapartum stillbirth proportion is with respect to a 28-week
gestational age definition of stillbirth, we set $\gamma_{g=\text{late},m}=0$
for observations which use a 28-week definition, or comparable “late”
definition of stillbirth, and estimate adjustment factors for
$\gamma_{g=\text{early},m}$ to accommodate observations that use a 22-week
definition or comparable “early” definition.

The adjustment factors are additionally informed by supplementary datasets for
which timing-specific stillbirth counts are available under different
gestational age definitions. In particular, the high income group is informed
by anonymized country-level data from Euro-Peristat. Country names in these
data are not available and therefore cannot be joined with the NMR covariate
to inform the region or country-level estimates directly. For the low-income
group, we incorporate data from the Global Network Maternal Newborn Health
Registry, for which timing-specific stillbirth counts according to 22-week and
28-week definitions are typically available. While these latter observations
are not anonymized, counts that correspond to the 28 week definition are
already used to inform the country intrapartum stillbirth proportions. We
therefore model these observations to inform only $\gamma$ as follows. Let
$\dot{y}_{c,g}$ and $\dot{z}_{c,g}$ denote respectively intrapartum and
antepartum stillbirth counts for some country $c$ using gestational age
definition $g$. The counts are then modeled

$\displaystyle\dot{y}_{c,g}|\rho_{c,g}\sim\text{Binomial}(\dot{y}_{c,g}+\dot{z}_{c,g},\,\rho_{c,g})$
(14) $\displaystyle\text{logit}\rho_{c,g}=\nu_{c,g}+\gamma_{g,m[c]}$ (15)
$\displaystyle\nu_{c,g}\sim\text{Normal}(0,10^{2})$ (16)

where $\nu_{c,g}$ is given a vague prior and represents the mean under the
late gestational age definition. The difference between the proportions in the
early and late definitions is therefore captured by the adjustment factor
$\gamma_{\text{early},m[c]}$, where $m[c]$ denotes the income group of country
$c$.

For high-income countries $m=$ HIC, we assign
$\gamma_{g=\text{early},m=\text{HIC}}$ the prior

$\gamma_{g=\text{early},m=\text{HIC}}\sim\text{Normal}(0,1^{2}).$ (17)

We center the prior adjustment factor for the low-middle income group on the
estimate for the high income group,

$\gamma_{g=\text{early},m=\text{LMIC}}\sim\text{Normal}(\gamma_{g=\text{early},m=\text{HIC}},\sigma_{\gamma}^{2})$
(18)

All standard deviation terms to be estimated, $\sigma_{\beta_{r}}$,
$\sigma_{\beta_{c}}$, $\sigma_{\beta_{p}}$, $\sigma_{\Delta}$, and
$\sigma_{\gamma}$, are assigned a half-Normal(0, 1) prior.

###  3.3 Weighting adjustment for country-level IPSB

####  3.3.1 Overview

Since observations may only pertain to a specific context or geography (e.g.
certain health facilities or subnational region), we do not assume that
patterns in the observations generalize to parts of the full national
population. We therefore apply a post-estimation weighting step to construct a
national estimate as a weighted average of place-level estimates.

Ideally, data pertaining to some place $p$ would cover the entire country, and
capture timing information about all known stillbirths. In this case, the
country-level estimate of the IPSB, ${\hat{\phi}}_{c,t}$, would just equal the
place-level estimate ${\hat{\phi}}_{p,t}$. In the case where the the entire
population of stillbirths is instead reported across multiple places, we would
like to know the true proportion $w_{p}$ of the country’s stillbirths
belonging to each place to use as weights. In such a scenario a reasonable
estimate of the national rate $\phi_{c,t}$ would simply be those weights
applied to the place-level rates,

$\hat{\phi}_{c,t}=\sum_{p:c(p)=c}w_{p}\hat{\phi}_{p,t}.$ (19)

In practice however, the observed sub-populations are not in general
exhaustive (i.e. $\sum w_{p}<1$), meaning that there are stillbirths that are
not captured in any of the data sources for the country. To account for this
remainder, we add an “unobserved” component so that the final estimate is

$\hat{\phi}_{c,t}=\sum_{p:c[p]=c}w_{p}\hat{\phi}^{\text{obs}}_{p,t}+\left(1-\sum_{p:c[p]=c}w_{p}\right)\hat{\phi}_{c,t}^{\text{unobs}}.$
(20)

Here, $\hat{\phi}^{\text{obs}}_{p,t}$ is the estimate specific to place $p$,
that is, using the estimate of its intercept and time trend. To construct
$\hat{\phi}_{c,t}^{\text{unobs}}$ we assume a generic place in country $c$
with unknown intercept and time trend. This assumes that the unobserved
stillbirths are centered at the prediction based on its SDG region, NMR level,
and estimated country intercept, with additional uncertainty based on the
estimated between-sub-population variation and temporal variation. Details on
these components are given in the following sections.

####  3.3.2 Construction of weights

In practice, the true weights $w_{p}$ are unknown. We instead construct an
estimate $\hat{w}_{p}$ as the ratio of the number of observed classified
stillbirths in place $p$ to the number of total stillbirths expected
nationally, based on UN IGME estimates of overall stillbirths. Under this
setup, the estimates for countries at the extremes of data quality are similar
to those which would arise from a more conventional hierarchical model where
data are modeled at the country level. For instance, if a country reports a
single data source with full coverage (e.g. high quality CRVS), then there is
only one place $p$ (which accounts for all stillbirths in the country), so
$w_{p}=1$, and the national estimate $\hat{\phi}_{c,t}$ is simply the estimate
$\hat{\phi}_{p,t}$ informed by that data source. On the other hand, if a
country has no data, then the entire estimate consists of
$\hat{\phi}_{c,t}^{\text{unobs}}$, and is informed entirely by the region and
NMR.

To construct the weights for the observed place components, first let
$s_{i}=y_{i}+z_{i}$ denote the sum of observed stillbirths classified as
intrapartum or antepartum in observation $n$. Let
$\tilde{S}_{i}=\tilde{S}_{c[i],t[i]}$ denote the estimate of total stillbirths
from the UN IGME total stillbirth rate model in the country $c$ and year $t$
corresponding to the observation. To reflect uncertainty in the number of
stillbirths, we directly use posterior samples of $\tilde{S}_{i}$ when
computing our own posterior samples.

For a place $p$, we construct its weight $\hat{w}_{p}$ as ratio of the sum of
observed classified stillbirths in that sub-population,
$\sum_{i:p[i]=p}s_{i}$, to the sum of estimated stillbirths nationally in the
country-years of those observations $\sum_{i:p[i]=p}\tilde{S}_{i}$:

$\hat{w}_{p}=\frac{\sum_{i:p[i]=p}s_{i}}{\sum_{i:p[i]=p}\tilde{S}_{i}}.$ (21)

Terms in the denominator are scaled proportionally for partial year
observations. For example, if an data point refers only to one third of the
year 2016, then we only add $S_{c[i],t[i]=2016}/3$ to the denominator for this
observation.

For a country $c$, we sum the weights of its observed sub-populations to give
the weight given to the observed population
$\hat{w}_{c}=\sum_{p:c[p]=c}\hat{w}_{p}$. If this quantity exceeds 1, we
downscale the weights to add to 1. The weight assigned to the unobserved
portion is the remainder $1-\hat{w}_{c}$.

The final estimate of the intrapartum proportion $\hat{\phi}_{c,t}$ for a
country is a weighting of observed sub-population components
$\hat{\phi}_{p,t}^{\text{obs}}$ and an unobserved component
$\hat{\phi}_{c,t}^{\text{unobs}}$:

$\hat{\phi}_{c,t}=\left(\sum_{p:c[p]=c}\hat{w}_{p}\hat{\phi}_{p,t}^{\text{obs}}\right)+(1-\hat{w}_{c})\hat{\phi}_{c,t}^{\text{unobs}}.$
(22)

####  3.3.3 Estimate for observed component

The sub-population weights $w_{p}$ are applied to sub-population estimates,
which are calculated from the estimated parameters,

$\displaystyle\hat{\mu}_{p,t}^{\text{obs}}=\hat{\beta}_{0}+\hat{\beta}_{r[c[p]])}+\hat{\beta}_{c[p]}+\hat{\beta}_{p}+\hat{\beta}_{\text{NMR}}\log\tilde{\text{NMR}}_{c[p],t}+\hat{\eta}_{p,t}$
(23)
$\displaystyle\hat{\phi}_{p,t}^{\text{obs}}=\text{logit}^{-1}(\hat{\mu}_{p,t}^{\text{obs}}).$
(24)

We directly use posterior draws of the country-year neonatal mortality
estimates from UN IGME, denoted $\tilde{\text{NMR}_{c,t}}$, to incorporate
appropriate uncertainty about the covariate.

####  3.3.4 Estimate for unobserved component

The “unobserved” component for a country is centered at the estimate given its
region and country intercepts and NMR level:

$\displaystyle\hat{\mu}_{c,t}^{\text{unobs}}=\hat{\beta}_{0}+\hat{\beta}_{r[c]}+\hat{\beta}_{c}+\tilde{\beta}_{p_{c}}+\hat{\beta}_{\text{NMR}}\log\tilde{\text{NMR}}_{c,t}+\tilde{\eta}_{c,t}$
(25)
$\displaystyle\hat{\phi}_{c,t}^{\text{unobs}}=\text{logit}^{-1}(\hat{\mu}_{c,t}^{\text{unobs}})$
(26)

where $\hat{\beta}_{0}$, $\hat{\beta}_{r[c]}$, and $\hat{\beta}_{\text{NMR}}$
are the usual parameters estimated from data, whereas $\tilde{\beta}_{p_{c}}$
and $\tilde{\eta}_{c,t}$ are new realizations of the sub-population effect and
time trend to reflect appropriate uncertainty about the unobserved population.

We assume that the unobserved component consists of a single sub-population.
For the place effect $\beta_{p}$, we use a new realization denoted
$\tilde{\beta}_{p_{c}}$,

$\tilde{\beta}_{c}\overset{RNG}{\sim}\text{Normal}(0,\hat{\sigma}_{\beta_{p}}^{2}).$
(27)

The new realization of the time trend, denoted $\tilde{\eta}_{c,t}$, is
generated according to the distribution

$\displaystyle\tilde{\eta}_{p,t}=\sum_{h=1}^{H}k_{h}(t)\tilde{\alpha}_{h,c}$
(28)
$\displaystyle\tilde{\Delta}_{h,c}=\tilde{\alpha}_{h,c}-\tilde{\alpha}_{h-1,c}$
(29)
$\displaystyle\tilde{\Delta}_{c}\overset{RNG}{\sim}\text{Normal}(0,\hat{\sigma}_{\Delta}^{2}).$
(30)

###  3.4 Obtaining regional-level estimates

The intrapartum stillbirth proportion for a region is obtained by multiplying
the country-level intrapartum stillbirth proportions by total stillbirth
counts, then aggregating at the region level and recalculating the proportion.
In particular, if $\hat{\phi}_{c,t}$ denotes the estimate for country $c$ at
time $t$, and $\tilde{S}_{c,t}$ denotes the corresponding total stillbirth
count, then $\hat{\phi}_{c,t}\tilde{S}_{c,t}$ is our estimate of the number of
intrapartum stillbirths. We obtain the intrapartum stillbirth proportion for a
region $r$ by calculating

$\hat{\phi}_{r,t}=\frac{\sum_{c:r[c]=r}\hat{\phi}_{c,t}\tilde{S}_{c,t}}{\sum_{c:r[c]=r}\tilde{S}_{c,t}},$
(31)

where the numerator represents the sum of intrapartum stillbirths in countries
in region $r$, and the denominator represents all stillbirths in countries in
region $r$. Here we again directly use posterior samples of the total
stillbirth counts $\tilde{S}_{c,t}$ when computing our own posterior samples.

###  3.5 Computation

We obtain samples of posterior distributions of parameters using Hamiltonian
Monte Carlo (HMC), a type of Markov Chain Monte Carlo sampling, implemented in
Stan using the cmdstanr interface (Gabry & Češnovar 2022). For each of 4 HMC
chains, we perform 1000 warmup iterations before taking the next 1000 post-
warmup iterations as our posterior samples. We monitor for convergence using
standard diagnostic measures. R-hat ($\hat{R}$) values for all parameters are
under 1.02 (Stan Development Team 2019), and visual inspection of traceplots
do not indicate sampling pathologies.

##  4 Results

In this section we illustrate regional estimates of the proportion of
stillbirths that are intrapartum (IPSB) over the period of interest.
Additionally, we demonstrate the impact of several key components of the
modeling process on country-level estimates of IPSB.

Figure 2 shows regional estimates of IPSB, with the shaded areas representing
90% credible intervals. There is substantial regional variation in the level
and trends of IPSB. Australia and New Zealand, and Europe, have low and
declining IPSB of around 10%, while regions such as Northern Africa and Sub-
Saharan Africa have IPSBs of 40-50% and show little evidence of decline at the
regional level over the time period. Eastern Asia shows evidence of the most
rapid decline since 2000, falling from around 40% to 12% over the 20-year
period. The largest uncertainty in estimates is in the Oceania region,
reflecting the lack of available data in this region.

Figure 2: Regional estimates of the proportion of stillbirths that are
intrapartum.

Figure 2 shows the impact of the inclusion of the gestational age adjustment
in the model. Country ‘A’ shown only has data available on the timing of
stillbirths based on the early (22 weeks gestation) period. Estimates without
the gestational age adjustment (shown in blue on the right hand side of the
graph) tend to follow the data closely. With the addition of the adjustment in
the model, the final estimates (shown in red on the left hand side) are
shifted downwards, to reflect the fact that stillbirths following the late (28
weeks) definition are likely to have a lower IPSB. Note that the the posterior
median estimate for the gestational age adjustment $\gamma_{g=\text{early},m}$
in high-income countries is 0.31 (90% CI: [0.29, 0.33]).

Figure 3: Effect of gestational age adjustment on final estimates for country
‘A’.

In Figure 4, we illustrate the effect on estimates of allowing non-sampling
error to vary by data source type for two countries with different types of
data available. In the figure, estimates shown in blue in the right-hand
column are based on a model with no non-sampling error estimated, while
estimates shown in red in the left-hand side are the final estimates from the
full model. For Country ‘B,’ which has data from a national high-quality vital
registration system, there is no effect of non-sampling error on the
estimates, as this type of error is assumed to be zero for such data sources.
In contrast, Country ‘C’ has data from several subnational sources, both from
HMIS and population-based studies. As these subnational sources have
additional error that is estimated to be non-zero, the resulting estimates
over time are smoother, and in particular are less influenced by the
particularly high observations around 2009.

Figure 4: Effect of allowing for non-sampling error to vary by source type on
final estimates for country ‘B’ (which has vital registration data) and
country ‘C’ (which has a mix of subnational data available). Points that
report on the same place or sub-population are connected by lines.

Finally we illustrate the impact of the post-estimation weighting scheme on
country-level estimates for three countries with different data availability
situations (Figure 5). In the “No weighting” setup on the right, we use a more
conventional hierarchical setup in which the $\beta_{p}$ terms are still
estimated when there are multiple data sources, but are not included when
presenting the national estimates. Instead of place-level time trends, a
single national time trend is shared across sub-populations.

Country ‘D’ has data that come from a CRVS system and capture all stillbirths
in the country. As such, the estimated weights for that data series are
approximately 1 and the post-estimation weighting scheme has no effect on the
final estimates. In contrast Countries ‘E’ and ‘F’ have one or more
subnational data sources, which capture well below the estimated total number
of stillbirths. The final estimates that include the re-weighting step thus
indicate a higher level of uncertainty around the IPSB levels.

Figure 5: Effect of weighting scheme on final estimates, for countries with a
range of different data availability situations.

###  4.1 Validation

We evaluated the model using two out-of-sample prediction exercises. In the
first exercise we fit the model using observations before 2017 as the training
set (which represent 78% of all observations), leaving out observations from
2017 and later to use as a test set. In the second exercise we perform a
10-fold cross-validation. Results from the validation exercises are given
below. We calculate mean absolute error as the average of the absolute
differences between posterior predicted medians and observed IPSB proportions
in held-out data points, and prediction interval coverage as the proportion of
observed IPSB proportions which fall within their respective posterior
prediction interval.

Table 4 shows validation metrics from the first prediction exercise, and Table
5 shows validation metrics from 10-fold cross validation. Mean absolute error
in both validation exercises is under 5%, and overall prediction interval
coverage in both cases is close to, but slightly lower than nominal, which
suggests a slight underestimation of the variability in the data.

There are some apparent regional differences in the performance of the model.
In the Oceania region, there is only one data point in 2018 which also
represents the only available observation in that country, resulting in high
error in prediction in this region. In the Northern Africa and Western Asia
region, prediction coverage is low, possibly due to relatively high between-
and within-country variability. Particularly poor performance here in the
first validation exercise is partially due to a series of data points for a
country in this region which deviate from the NMR trend beyond what is
expected from the estimated spline component. Finally, high error in the point
estimates for Central and Southern Asia region can be partially explained by a
relatively high proportion of low-quality subnational observations which are
difficult to predict. The reasonably high prediction interval coverage,
however, suggests that variability and uncertainty are well calibrated for
observations in this region.

Table 4:  Model evaluation metrics using 2000-2016 data as a training set and data from 2017 onward as a test set. Region | Mean absolute error |  95% prediction interval coverage
---|---|---
Global | 0.044 |  0.917
Central and Southern Asia | 0.095 |  0.850
Eastern and South-Eastern Asia | 0.014 |  1.000
Latin America and the Caribbean | 0.028 |  0.960
North America, Europe, Australia and New Zealand | 0.028 |  0.922
Northern Africa and Western Asia | 0.049 |  0.615
Oceania (exc. Australia and New Zealand) | 0.271 |  1.000
Sub-Saharan Africa | 0.052 |  0.979
Table 5:  Model evaluation metrics from 10-fold cross validation. Region | Mean absolute error |  95% prediction interval coverage
---|---|---
Global | 0.041 |  0.925
Central and Southern Asia | 0.102 |  0.943
Eastern and South-Eastern Asia | 0.039 |  0.980
Latin America and the Caribbean | 0.025 |  0.928
North America, Europe, Australia and New Zealand | 0.021 |  0.927
Northern Africa and Western Asia | 0.029 |  0.829
Oceania (exc. Australia and New Zealand) | 0.271 |  1.000
Sub-Saharan Africa | 0.063 |  0.916

##  5 Discussion

In this paper we proposed a Bayesian hierarchical penalized splines regression
model to estimate the proportion of stillbirths that are intrapartum in all
countries worldwide. The model sought to address a multitude of data
availability, quality, and coverage issues, which makes obtaining reliable
estimates of this indicator particularly challenging in low- and middle-income
countries. The proposed model includes a measurement error model to take into
account different data collection systems; a gestational age definition
adjustment, which is informed by high-quality data across a range of different
countries; and hierarchical pooling of information from observations within
the same country and region of the world. Trends are informed by both trends
in the neonatal mortality rate, and observed data through a penalized splines
set up, which allows for temporal patterns to be modeled in data-sparse
situations. Validation exercises suggest that the model is reasonably well
calibrated.

Results at the regional level suggest substantial variation in the levels and
trends in the proportion of stillbirths that are intrapartum globally. While
proportions are low and stable in regions such as Australia/New Zealand and
Europe (at around 10%), results suggest that more than half of all stillbirths
occur intrapartum in Southern Asia and Sub-Saharan Africa. These high levels
coupled with slow declines over time suggest there is still much progress to
be made to end preventable stillbirths worldwide.

Our approach treats each series of observations as information about the
timing of stillbirths from a particular ‘place’ or sub-population. The
coverage of a ‘place’ can range from a whole country, right down to a single
health facility in a subnational area. The subsequent country-level estimate
of IPSB is then a weighted average of the place-level effects, plus unobserved
effects, where the weights are calculated based on the relative share of
overall stillbirths observed. In comparison to assuming the every data series
is an observation of the IPSB for a whole country, this approach generally
leads to larger estimated levels of uncertainty. As the weights are estimated
based on the number of observed stillbirths alone, we are allowing larger
samples to influence final estimates more heavily than smaller samples. This
approach was partly employed to downweigh the influence of non-representative
smaller samples, as it is unlikely that ‘place’-level samples are nationally
representative, given that they may only cover government-owned health
facilities, or health facilities in a particular city or rural area, for
example. In future work, we plan to investigate the construction of these
weights in more depth. For instance, in some cases we may have detailed enough
information about the source of a subpopulation data series to infer the types
of women giving birth that were captured in the data, particularly if the
geographic location of the observations is known and can be matched with
census or survey data. In these cases, a post-stratification based approach
may be possible, to more adequately control for the non-representativeness of
different subpopulation samples.

Currently, estimates of stillbirth by timing and the total stillbirth rate are
performed separately, with the total stillbirth rate being estimated using a
Bayesian sparse regression model with temporal smoothing (Wang et al. 2022).
As detailed above, we combine our estimates of the intrapartum proportion with
posterior samples of total stillbirths to obtain regional estimates and
associated uncertainty. Future work could investigate the possibility of
estimating both quantities in the same modeling framework. As there is
generally more data available for total stillbirth rates, combining the
estimation of two approaches may help to inform estimates by timing, and would
also help to better incorporate uncertainty from different sources (Schumacher
et al. 2022).

Stillbirths have traditionally been overlooked by demographers and public
health professionals, with more emphasis being placed on improving infant and
maternal survival outcomes. But as substantial health inequalities persist
across regions of the world, there have been calls to monitor stillbirths,
working towards the goal of ending all preventable stillbirths. From a
population dynamics perspective, there has been increased work in
understanding fetal outcomes in conjunction with maternal and infant outcomes,
in order to allow for a more complete assessment of mortality conditions and
improvement (Hathi 2022). The estimation work presented here aimed to improve
estimates of stillbirths by timing, but by doing so highlighted the vast array
of data availability and quality issues in this area. Future resource efforts
should not only focus on improving access to essential healthcare during
pregnancy and childbirth, but also in helping national stakeholders to improve
data systems to collect reliable information in order to help those most in
need.

## References

reDa Silva FT, Gonik B, McMillan M, Keech C, Dellicour S, Bhange S, Tila M,
Harper DM, Woods C, Kawai AT, others (2016) Stillbirth: Case definition and
guidelines for data collection, analysis, and presentation of maternal
immunization safety data. _Vaccine_ 34:6057

preDe Bernis L, Kinney MV, Stones W, Hoope-Bender P ten, Vivio D, Leisher SH,
Bhutta ZA, Gülmezoglu M, Mathai M, Belizán JM, others (2016) Stillbirths:
Ending preventable deaths by 2030. _The lancet_ 387:703–716

preFrøen JF, Friberg IK, Lawn JE, Bhutta ZA, Pattinson RC, Allanson ER,
Flenady V, McClure EM, Franco L, Goldenberg RL, others (2016) Stillbirths:
Progress and unfinished business. _The Lancet_ 387:574–586

preGabry J, Češnovar R (2022) Cmdstanr: R interface to ’CmdStan’.

preGetahun D, Ananth CV, Kinzler WL (2007) Risk factors for antepartum and
intrapartum stillbirth: A population-based study. _American journal of
obstetrics and gynecology_ 196:499–507

preGissler M, Mohangoo AD, Blondel B, Chalmers J, Macfarlane A, Gaizauskiene
A, Gatt M, Lack N, Sakkeus L, Zeitlin J (2010) Perinatal health monitoring in
europe: Results from the EURO-PERISTAT project. _Informatics for Health and
Social Care_ 35:64–79

preHathi P (2022) Population science implications of the inclusion of
stillbirths in demographic estimates of child mortality.

preHug L, Alexander M, You D, Alkema L, Child UIG for (2019) National,
regional, and global levels and trends in neonatal mortality between 1990 and
2017, with scenario-based projections to 2030: A systematic analysis. _The
Lancet Global Health_ 7:e710–e720

preHug L, Mishra A, Lee S, You D, Moran A, Strong KL, Cao B (2020) A neglected
tragedy the global burden of stillbirths: Report of the UN inter-agency group
for child mortality estimation, 2020.

preHug L, You D, Blencowe H, Mishra A, Wang Z, Fix MJ, Wakefield J, Moran AC,
Gaigbe-Togbe V, Suzuki E, others (2021) Global, regional, and national
estimates and trends in stillbirths from 2000 to 2019: A systematic
assessment. _The Lancet_ 398:772–785

preJoyce R, Webb R, Peacock J (2004) Associations between perinatal
interventions and hospital stillbirth rates and neonatal mortality. _Archives
of Disease in Childhood-Fetal and Neonatal Edition_ 89:F51–F56

preLawn JE, Blencowe H, Waiswa P, Amouzou A, Mathers C, Hogan D, Flenady V,
Frøen JF, Qureshi ZU, Calderwood C, others (2016) Stillbirths: Rates, risk
factors, and acceleration towards 2030. _The Lancet_ 387:587–603

preMoxon SG, Ruysen H, Kerber KJ, Amouzou A, Fournier S, Grove J, Moran AC,
Vaz LM, Blencowe H, Conroy N, others (2015) Count every newborn; a measurement
improvement roadmap for coverage data. _BMC pregnancy and childbirth_ 15:1–23

preSchumacher AE, McCormick TH, Wakefield J, Chu Y, Perin J, Villavicencio F,
Simon N, Liu L (2022) A flexible bayesian framework to estimate age-and cause-
specific child mortality over time from sample registration data. _The Annals
of Applied Statistics_ 16:124–143

preStan Development Team (2019) Stan modeling language users guide and
reference manual version 2.25.

preUN IGME (2021) Levels and trends in child mortality levels and trends in
child mortality.

preWang Z, Fix MJ, Hug L, Mishra A, You D, Blencowe H, Wakefield J, Alkema L
(2022) Estimating the stillbirth rate for 195 countries using a bayesian
sparse regression model with temporal smoothing. _The Annals of Applied
Statistics_ 16:2101–2121

preWorld Health Organization (2014) Every newborn: An action plan to end
preventable deaths.

p


# AUTHOR GUIDELINES FOR IWAENC 2022 PROCEEDINGS MANUSCRIPTS

# STREAMING NOISE CONTEXT AWARE ENHANCEMENT FOR AUTOMATIC SPEECH RECOGNITION
IN MULTI-TALKER ENVIRONMENTS

###### Abstract

One of the most challenging scenarios for smart speakers is multi-talker, when
target speech from the desired speaker is mixed with interfering speech from
one or more speakers. A smart assistant needs to determine which voice to
recognize and which to ignore and it needs to do so in a streaming, low-
latency manner. This work presents two multi-microphone speech enhancement
algorithms targeted at this scenario. Targeting on-device use-cases, we assume
that the algorithm has access to the signal before the hotword, which is
referred to as the noise context. First is the Context Aware Beamformer which
uses the noise context and detected hotword to determine how to target the
desired speaker. The second is an adaptive noise cancellation algorithm called
Speech Cleaner which trains a filter using the noise context. It is
demonstrated that the two algorithms are complementary in the signal-to-noise
ratio conditions under which they work well. We also propose an algorithm to
select which one to use based on estimated SNR. When using 3 microphone
channels, the final system achieves a relative word error rate reduction of
55% at -12 dB, and 43% at 12 dB.

Index Terms—  multi-talker, array, beamformer, noise cancellation

†† We thank Alex Gruenstein and Turaj Shabestary for several useful
discussions and feedback and Jon Moeller and Brant Jameson for their help with
data collection.

##  1 Introduction

Automatic speech recognition (ASR) has shown robustness in the presence of
noise largely due to the adoption of neural network based acoustic models [1,
2, 3, 4],large scale training [5, 6, 7], and improved data augmentation
strategies [8, 9, 10]. Multiple speaker scenarios, however, still pose a
challenge [11, 12]. This is particularly true for smart speakers where it is
desired to respond to one of the speech sources and not the others. In
addition, smart speakers require streaming, low-latency solutions so the
desired speaker needs to be quickly determined and isolated from the other
sources.

Solutions aimed at separating multiple speakers have been proposed using
either a single [13, 14, 15] or multiple microphones [16, 17] are often
designed for separating the multiple voices rather than identifying one
target. Filtering based on a targeted speaker has shown success but that
targeted speaker needs to have registered with the device. [18]

A method which has shown considerable improvement isolates the desired source
using a time-frequency mask and then uses the statistics of that masked source
to steer a beamformer [19, 20]. While effective in non-speech noise, it
encounters challenges deciding between one or more voices. The entire
utterance is usually used to steer the beamformer so the streaming, low-
latency requirements of the smart speaker are not met. When operating under
such constraints, these techniques have had difficulties [21, 22].

Two techniques are studied to address speech-based noise while meeting
streaming, low-latency criteria. The device is assumed to trigger on a
wakeword or hotword phrase such as “OK Google” or “Alexa” that precedes the
user request or query. The desired speaker is, by definition, the one speaking
the hotword. These techniques function on device, which enables them to make
use of the signal that occurs directly before the hotword, referred to as the
noise context.

The first is the Context Aware Beamformer (CAB) which uses the noise context
and the detected hotword to determine the steering vector. The second is an
adaptive noise cancellation algorithm called Speech Cleaner. It has resulted
in significant improvement in hotword recognition in noisy environments[23,
24]. Here it is applied to the query as well. Both will be compared to a
beamformer steered with an ideal mask.

The rest of the paper is organized as follows. An overview of the system is
given in Section 2. In Section 3, we describe context aware beamforming.
Section 4 will describe the adaptive noise cancellation algorithm called the
Speech Cleaner. Section 5 will describe the data used and the experimental
setup. Evaluations will be provided in Section 6.

##  2 Signal Model

Each received utterance consists of a hotword followed by the query as well as
some audio before the hotword. The received short-time Fourier transform
(STFT)-processed $M$ microphone signal is written as

$\mathbf{Y}(k,n)\triangleq[Y_{0}(k,n),Y_{1}(k,n),...Y_{(M-1)}(k,n)]^{T},$ (1)

where $Y_{i}(k,n)$ is the signal from microphone $i$ at frame $n$ for
frequency subband $k$. The received signal can be written as the summation

$\mathbf{Y}(k,n)=\mathbf{X}(k,n)+\mathbf{V}(k,n),$ (2)

where $\mathbf{X}(k,n)$ is a vector of the M reverberant desired speech
signals and $\mathbf{V}(k,n)$ corresponds to the noise signals which may also
include interfering speech. Because all enhancement is conducted independently
for each frequency subband, the frequency index $n$ shall subsequently be
omitted.

##  3 Context Aware Beamforming

A beamformer is a spatial filter applied to the $M$-microphone input signal to
get an estimate of the desired signal

${\hat{X}}(n)=\mathbf{W}^{H}\mathbf{{Y}}(n),$ (3)

where $\mathbf{W}$ is a vector of the $M$ coefficients of the beamformer and
$\hat{X}(n)$ is a single channel estimate of the desired signal.

Perhaps the most challenging aspect of beamforming, especially in reverberant
and multi-talker environments, is initializing the steering vector. using the
principal eigenvector of the spatial correlation matrix of desired signal,
which will collect the most energy of the desired signal, has been shown to be
successful in reverberation[25, 26]. In mask-based beamforming[20, 16], neural
networks are used to estimating this correlation matrix. Multi-talker
conditions, where there is one desired speaker and one or more interfering
speakers, present challenges to ascertain which speaker it should be
targeting. Doing so in a streaming, low-latency manner increases those
challenges.

The Context Aware Beamformer utilizes the hotword and the noise context to
target the desired speaker. The spatial correlation matrix of the desired
speaker $\mathbf{X}(n)$ is defined as:

$\mathbf{\Phi_{XX}}=E_{n}[\mathbf{X}(n)\mathbf{X}(n)^{H}].$ (4)

This can be estimated by taking the difference between the correlation matrix
calculated when there is both noise and signal and the correlation matrix when
there is only noise:

$\hat{\mathbf{\Phi}}_{\mathbf{XX}}\approx\mathbf{\Phi}_{\mathbf{YY}}-\mathbf{\Phi}_{\mathbf{VV}},$
(5)

wherein it is assumed that the desired signal and the noise are uncorrelated
[27].

In order to make use of this approximation, time periods that are only noise
and that are signal plus noise must be identified. To that end, this work
makes uses of the hotword. When the device identifies the hotword, it is, by
definition, the desired speaker saying the hotword. Because noise will likely
also be present, the period of time that is detected as containing the hotword
can be used to estimate $\mathbf{\Phi}_{\mathbf{YY}}$.

Targeting on-device application, we assume that the model has access to the
noise context, the period of time before the hotword. This segment of time can
be assumed to contain only noise and that noise will have similar statistics
to the noise in the hotword. This can be used to estimate
$\mathbf{\Phi_{\mathbf{VV}}}$.

With these two estimates, Equation (5) can be applied to obtain an estimate of
the spatial correlation of the desired signal. The principal eigenvector of
this matrix is then used as the initial beamformer coefficients:

$\mathbf{\hat{W}}=Pr\\{\hat{\mathbf{\Phi}}_{\mathbf{XX}}\\},$ (6)

where, $Pr\\{\centerdot\\}$ signifies finding the principal eigenvector of the
matrix. The beamformer coefficients can then be adapted using least mean
squares (LMS) in accordance with the minimum variance distortionless response
(MVDR) criterion [28].

##  4 Speech Cleaner

Speech Cleaner (SC), an adaptive noise cancellation algorithm, is summarized
in this section[23, 24]. In previous works, it has been shown to result a
significant increased in robustness to noise in hotword detection. Here, it is
to applied it directly to the target query in addition to the hotword.

A finite impulse response (FIR) filter with a tapped delay line of length $L$
is applied to the signals from all microphones except for one. The summed
output of these is subtracted from the received signal at the first
microphone:

$Z_{m}(n)=Y_{0}(n)-\sum_{m=1}^{M-1}{\mathbf{U}_{m}}^{H}\tilde{\mathbf{Y}}_{m}(k,n),$
(7)

where,

$\tilde{\mathbf{Y}}_{m}(n)=[Y_{m}(n),Y_{m}(n-1),\cdots Y_{m}(n-(L-1))]^{T},$
(8)

is a vector of time delayed STFT-processed input for microphone $m$ and,

$\mathbf{U}_{m}(k)=[U_{m}(k,0),U_{m}(k,1),\cdots U_{m}(k,L-1)]^{T}$ (9)

is a vector of the filter coefficients to be applied to microphone input $m$.
The filter coefficients are specified as those that minimize the power of the
enhanced output

$\hat{\mathbf{U}}_{m}=\arg\underset{\mathbf{U}_{m}}{\min}E_{n}{[|Z_{m}(n)|^{2}]}.$
(10)

The minimization is done through adaptation during the noise context, when it
is assumed there is no desired speech present. In this way, the filter becomes
trained to cancel out the noise sources. When the hotword is detected, the
filter coefficients are frozen. The last coefficients before the hotword
detection are then applied to the query:

$\hat{{X}}(n)=Y_{0}(n)-\sum_{m=1}^{M-1}{\mathbf{\hat{U}}_{m}}^{H}\tilde{\mathbf{Y}}_{m}(k,n-l)$
(11)

to obtain a single channel estimate of the desired signal. In this way, the
filter works to cancel the noise, as it was trained, and not the desired
signal, which was not present during training. Unlike the beamformer, Speech
Cleaner has no constraint as to the amount of speech distortion introduced.

##  5 Experimental Setup

This study uses data re-rerecorded in a living-room lab with no sound
treatment. Desired speech and noise were recorded separately with a three
microphone array with two microphones on the top spaced 7.1 cm apart with a
third was on the front. The queries were played r from 7 different speaker
positions at a height of 1.5m and at 4m from the array. From each location 100
different queries prefaced with the hotword “OK Google” or “Hey Google” were
played at a volume of 40 dB over ambient room noise. Two types of noise,
speech-based from a movie and non-speech pink noise, were separately played
through a loudspeaker from the same 7 locations.

Speech and noise were mixed at several different SNR values ranging from -12
to 24 dB as well as a set with no added noise. In all cases, the recordings
contain whatever ambient noise was in the room during recording. They were
mixed such that the speech and noise did not originate from the same speaker
position and with 8s of noise recorded before the hotword commenced. The
multi-channel input was processed by the specified enhancement algorithm to
produce one channel of output which was then passed onto an ASR model.

A state-of-the-art end-to-end recurrent neural transducer model is used for
evaluations[29]. The model was trained using anonymized, hand-transcribed
English utterances from domains like VoiceSearch, Farfield, Telephony and
YouTube. Data augmentation was applied via a room simulator[8] to model
reverberation times from 0 ms and 900 ms and SNRs from 0 dB to 30 dB.
Consequently, this model is relatively robust to moderate noise conditions.

##  6 Results

Table 1: Two Channel Results SNR | Baseline | CAB | SC | SNR | Oracle
---|---|---|---|---|---
(dB) |  |  |  | Sel. | Bmfm
-12 | 100.2 | 95.4 | 64.5 | 68.0 | 95.7
-6 | 83.0 | 56.6 | 43.6 | 48.8 | 66.0
0 | 48.7 | 34.4 | 24.7 | 28.3 | 29.7
6 | 20.5 | 17.9 | 12.0 | 12.4 | 11.7
12 | 9.2 | 6.8 | 7.6 | 5.6 | 5.3
$Clean$ | 3.0 | 1.9 | 3.0 | 1.9 | 1.8
Table 2: Three Channel Results SNR | Baseline | CAB | SC | SNR | Oracle
---|---|---|---|---|---
(dB) |  |  |  | Sel. | Bmfm
-12 | 100.2 | 93.0 | 42.0 | 45.4 | 92.6
-6 | 83.0 | 66.2 | 28.9 | 33.0 | 60.4
0 | 48.7 | 31.7 | 16.6 | 20.1 | 25.3
6 | 20.5 | 11.6 | 10.9 | 11.1 | 10.1
12 | 9.2 | 4.6 | 8.0 | 5.2 | 4.3
$Clean$ | 3.0 | 1.4 | 3.9 | 1.4 | 1.6

###  6.1 Context Aware Beamformer

CAB was tested with multi-talker noise datasets using two and three microphone
configurations. No adaptation was done on the beamformer coefficients after
the initialization with principal eigenvector. The word error rate (WER)
results are presented for the two channel case in the third column of Table 1
and the three channel case in Table 2. An SNR value of $Clean$ indicates the
case where there was no added multi-talker noise. Also added for comparison is
the baseline case where the ASR model was applied directly with no
enhancement. In addition, there is an Oracle case where the spatial
correlation matrix of the desired signal is calculated directly from the
isolated desired signal. This is the ideal case for mask-based beamforming.

It can be seen that at the lowest SNR values the beamformer provided little
benefit. In these cases, where the noise power is much larger than the signal,
it can be difficult to get an accurate approximation for the desired signal
correlation from 5. In contrast, at 6 dB and above the relative WER
improvements approach 50% where it approaches the performance of the Oracle.
Going from two to three microphones provides a notable improvement.

###  6.2 Speech Cleaner

The fourth column of Tables 1 and 2 show the performance of Speech Cleaner
with the multi-talker data sets using two and three microphones, respectively.
The filter length was $L=3$ with adaption using recursive least squares (RLS)
[24] during the 8s noise context. The coefficient values at the detected
hotword start time are then used to filter the hotword and query.

The Speech Cleaner gives the most improvement at the lower SNR values where it
significantly outperforms the Oracle beamformer. At -12dB, there is a 40%
relative WER improvement with two microphones, improving to 60% with a third
mic. However, at the higher SNR values, improvement is lower and pales in
comparison to the beamformer. In the no added noise case, adding an additional
third microphone actually causes a performance degradation. Speech Cleaner,
unlike CAB, does not contain any constraints with regard to signal distortion.
In lower noise, the noise cancellation benefits are outweighed by the signal
distortion.

###  6.3 SNR Based Selection

The third and fourth columns in Tables 1 and 2 show the CAB and Speech Cleaner
to be complementary in the SNR ranges at which they best perform. Therefore,
it is desirable to select the appropriate algorithm in the corresponding
conditions.

A simple method to estimate the SNR of the utterance is used. That estimate is
used to select between the two enhancement techniques. The SNR estimate is
obtained by separately measuring the power in the hotword segment and the
power in the noise context. Because the hotword contains both signal and
noise, its power is an estimate of the power of the signal plus that of the
noise. The context power is an estimate of the noise power. Consequently,

$\widehat{SNR}=\frac{\sigma_{hotword}^{2}}{\sigma_{context}^{2}}-1=\frac{\hat{\sigma}_{X}^{2}+\hat{\sigma}_{V}^{2}}{\hat{\sigma}_{V}^{2}}-1.$
(12)

If $\widehat{SNR}$ is less than a threshold $\gamma$ then the Speech Cleaner
will be used otherwise CAB will be applied.

The fifth column of Tables 1 and 2 shows results for the two microphone array
with $\gamma$ at 6 dB. With the SNR-based selector, performance approaches
that of the Speech Cleaner for SNR values of 6 dB and below. Above 6 dB,
performance approaches that of CAB. This yields better than a 30% and 43%WER
improvement across the considered SNR range for two and three channels,
respectively.

###  6.4 Non-Speech Noise

To consider performance in non-speech noise, Table 3 presents WER for pink
noise with three channels. The SNR-based selection tracks the performance of
the Speech Cleaner past the 6dB threshold leading to some undesired results.
This suggests that either the the thresholds could benefit from more tuning,
or alternative methods of selection between the techniques should be
considered. We leave that to future work.

Table 3: Three Channel Pink Noise Results SNR | Baseline | CAB | SC | SNR
---|---|---|---|---
(dB) |  |  |  | Sel.
-12 | 100.2 | 100.0 | 42.2 | 42.2
-6 | 98.9 | 89.3 | 25.5 | 25.5
0 | 80.1 | 53.0 | 18.0 | 18.0
6 | 42.2 | 18.0 | 12.3 | 12.3
12 | 12.4 | 5.3 | 8.8 | 9.1
18 | 3.3 | 2.5 | 8.6 | 5.6
24 | 1.7 | 1.5 | 8.3 | 1.6

###  6.5 Impact of Context Length

The results thus far used all 8s of noise context. In Table 4, WERs are shown
across the SNR range for shorter noise context lengths. It can be seen that,
particularly at the higher SNRs where CAB is more effective, reducing the
noise context to 1s still captures the bulk of the improvements. Benefits are
still seen if only 0.25s of context is available.

Table 4: Impact of Noise Context Length on CAB Length | -12dB | 6dB | 0dB | 6dB | 12dB | $\infty$
---|---|---|---|---|---|---
Baseline | 101.6 | 81.6 | 46.2 | 19.3 | 7.7 | 1.9
All(8s) | 95.7 | 71.5 | 33.9 | 12.7 | 5.6 | 1.9
3s | 95.7 | 71.8 | 34.3 | 13.7 | 5.9 | 1.9
1s | 97.3 | 74.5 | 36.6 | 14.0 | 6.0 | 1.9
0.5s | 98.0 | 75.3 | 37.6 | 14.0 | 6.0 | 1.9
0.25s | 98.6 | 76.5 | 39.1 | 14.9 | 6.2 | 2.1

###  6.6 Impact of Desired Speaker in Context

It has been assumed that the desired speaker is not present during the noise
context. While that is valid often, it will not always be the case. Here, the
algorithms operate with the desired speaker talking throughout the context
period. To simulate this, part of the query is appended before the hotword. In
Table 5 the impact on the performance of CAB and Speech Cleaner for these
cases with the two channel array is shown.

These results demonstrate the importance of ensuring that this assumption
holds or allowing for the disabling of these algorithms when it does not.

Table 5: Impact of Desired Speaker in Context SNR | Baseline | CAB | Speech |
---|---|---|---|---
(dB) |  |  | Cleaner |
-12 | 101.6 | 97.9 | 77.8 |
-6 | 81.6 | 82.0 | 69.9 |
0 | 46.2 | 53.3 | 66.1 |
6 | 19.3 | 28.8 | 67.8 |
12 | 7.7 | 15.5 | 75.2 |
$Clean$ | 1.9 | 4.6 | 91.0 |

##  7 Conclusion and Future Work

Two techniques were presented to help address the multi-talker scenario on
smart speakers. The Context Aware Beamformer used the noise context and
detected hotword to determine how to target the desired speaker. CAB was shown
to enable WER reduction across the SNR range, but it was most effective at
higher SNRs. We also present an adaptive noise cancellation algorithm called
Speech Cleaner trained using the noise context. This algorithm was shown to be
very effective at SNRs below 6 dB where it outperformed an oraacle mask-based
beamformer, but could degrade performance at higher SNRs. It was demonstrated
that the two algorithms are complementary in the conditions under which they
work well. To that end, an algorithm to select which algorithm to use based on
estimated SNR estimation was demonstrated. With the SNR-based selection,
relative error rates reductions of more than 30% were demonstrated for a two
microphone array, improving to more than 43% for three microphones.

Future work will encompass more advanced techniques to select between the two
algorithms and the original raw audio. Identifying conditions when the target
speaker also appears in the noise context is another interesting direction for
the future.

## References

  * [1] R. Prabhavalkar, K. Rao, T. N. Sainath, B. Li, L. Johnson, and N. Jaitly,  “A Comparison of Sequence-to-Sequence Models for Speech Recognition,”  in Proc. Interspeech, 2017.
  * [2] E. Battenberg, J. Chen, R. Child, A. Coates, Y. G. Y. Li, H. Liu, S. Satheesh, A. Sriram, and Z. Zhu,  “Exploring Neural Transducers for End-to-end Speech Recognition,”  in Proc. ASRU, 2017.
  * [3] T. Hori, S. Watanabe, Y. Zhang, and W. Chan,  “Advances in Joint CTC-Attention Based End-to-End Speech Recognition with a Deep CNN Encoder and RNN-LM,”  in Proc. Interspeech, 2017.
  * [4] B. Li, A. Gulati, J. Yu, T. N. Sainath, C.-. Chiu, A. Narayanan, S.-Y. Chang, et al.,  “A better and faster end-to-end model for streaming asr,”  in Proc. ICASSP, 2021.
  * [5] S. Mirsamadi and J. HL Hansen,  “On multi-domain training and adaptation of end-to-end rnn acoustic models for distant speech recognition,”  in Proc. Interspeech, 2017.
  * [6] D. Hakkani-Tür, G. Tür, A. Celikyilmaz, Y-N. Chen, J. Gao, L. Deng, and Y.-Y. Wang,  “Multi-domain joint semantic frame parsing using bi-directional rnn-lstm.,”  in Proc. Interspeech, 2016.
  * [7] A. Narayanan, A. Misra, K. C. Sim, G. Pundak, A. Tripathi, M. Elfeky, P. Haghani, T. Strohman, and M. Bacchiani,  “Toward Domain-Invariant Speech Recognition via Large Scale Training,”  in Proc. of SLT, 2018.
  * [8] C. Kim, A. Misra, K. Chin, T. Hughes, A. Narayanan, T. Sainath, and M. Bacchiani,  “Generation of large-scale simulated utterances in virtual rooms to train deep-neural networks for far-field speech recognition in Google Home,”  in Proc. Interspeech, 2017.
  * [9] D. S. Park, W. Chan, Y. Zhang, C.-C. Chiu, B. Zoph, E. D. Cubuk, and Q. V. Le,  “Specaugment: A simple data augmentation method for automatic speech recognition,”  Proc. Interspeech, 2019.
  * [10] I. Medennikov, Y. Y. Khokhlov, A. Romanenko, D. Popov, N. A. Tomashenko, I. Sorokin, and A. Zatvornitskiy,  “An investigation of mixup training strategies for acoustic models in ASR,”  in Proc. Interspeech, 2018.
  * [11] J. Barker, R. Marxer, E. Vincent, and S. Watanabe,  “The third CHiME speech separation and recognition challenge: Analysis and outcomes,”  Computer Speech & Language, vol. 46, pp. 605–626, 2017.
  * [12] J. Barker, S. Watanabe, E. Vincent, and J. Trmal,  “The fifth’chime’speech separation and recognition challenge: Dataset, task and baselines,”  in Proc. Interspeech, 2018.
  * [13] J. R. Hershey, Z. Chen, J. Le Roux, and S. Watanabe,  “Deep clustering: Discriminative embeddings for segmentation and separation,”  in Proc. ICASSP. IEEE, 2016, pp. 31–35.
  * [14] Y. Luo and N. Mesgarani,  “Conv-tasnet: Surpassing ideal time–frequency magnitude masking for speech separation,”  IEEE/ACM transactions on audio, speech, and language processing, vol. 27, no. 8, pp. 1256–1266, 2019.
  * [15] Y. Luo, Z. Chen, and T. Yoshioka,  “Dual-path rnn: efficient long sequence modeling for time-domain single-channel speech separation,”  in Proc. ICASSP. IEEE, 2020, pp. 46–50.
  * [16] X. Chang, W. Zhang, Y. Qian, J. Le Roux, and S. Watanabe,  “Mimo-speech: End-to-end multi-channel multi-speaker speech recognition,”  in Proc. ASRU. IEEE, 2019, pp. 237–244.
  * [17] X. Chang, W. Zhang, Y. Qian, J. Le Roux, and S. Watanabe,  “End-to-end multi-speaker speech recognition with transformer,”  in Proc. ICASSP. IEEE, 2020, pp. 6134–6138.
  * [18] Q. Wang, H. Muckenhirn, K. Wilson, P. Sridhar, Z. Wu, J. Hershey, R. A. Saurous, R. J .Weiss, Y. Jia, and I. L. Moreno,  “Voicefilter: Targeted voice separation by speaker-conditioned spectrogram masking,”  arXiv preprint arXiv:1810.04826, 2018.
  * [19] T. Higuchi, N. Ito, T. Yoshioka, and T. Nakatani,  “Robust mvdr beamforming using time-frequency masks for online/offline asr in noise,”  in Proc. ICASSP. IEEE, 2016, pp. 5210–5214.
  * [20] J. Heymann, L. Drude, and R. Haeb-Umbach,  “Neural network based spectral mask estimation for acoustic beamforming,”  in Proc. ICASSP. IEEE, 2016, pp. 196–200.
  * [21] J. Heymann, M. Bacchiani, and T. N. Sainath,  “Performance of mask based statistical beamforming in a smart home scenario,”  in Proc. ICASSP. IEEE, 2018, pp. 6722–6726.
  * [22] B. Li, T.N. Sainath, R. J. Weiss, K. W. Wilson, and M. Bacchiani,  “Neural network adaptive beamforming for robust multichannel speech recognition,”  2016\.
  * [23] Y. Huang, T. Z. Shabestary, and A. Gruenstein,  “Hotword cleaner: dual-microphone adaptive noise cancellation with deferred filter coefficients for robust keyword spotting,”  in Proc. ICASSP. IEEE, 2019, pp. 6346–6350.
  * [24] Y. Huang, T. Z. Shabestary, A. Gruenstein, and L. Wan,  “Multi-Microphone Adaptive Noise Cancellation for Robust Hotword Detection,”  in Proc. Interspeech 2019, 2019, pp. 1233–1237.
  * [25] M. Souden, J. Benesty, and S. Affes,  “On optimal frequency-domain multichannel linear filtering for noise reduction,”  IEEE Transactions on audio, speech, and language processing, vol. 18, no. 2, pp. 260–276, 2009.
  * [26] H.L. Van Trees,  Optimum Array Processing: Part IV of Detection, Estimation, and Modulation Theory,  Detection, Estimation, and Modulation Theory. Wiley, 2004.
  * [27] H. Taherian, Z-Q Wang, and D. Wang,  “Deep learning based multi-channel speaker recognition in noisy and reverberant environments,”  in Interspeech, 2019.
  * [28] O. L. Frost,  “An algorithm for linearly constrained adaptive array processing,”  Proceedings of the IEEE, vol. 60, no. 8, pp. 926–935, 1972.
  * [29] T. N. Sainath, Y. He, B. Li, A. Narayanan, R. Pang, A. Bruguier, S. Chang, W. Li, R. Alvarez, Z. Chen, et al.,  “A streaming on-device end-to-end model surpassing server-side conventional model quality and latency,”  in Proc. ICASSP. IEEE, 2020, pp. 6059–6063.


Proton reconstruction with the Precision Proton Spectrometer (PPS) in Run 2
and the PPS at HL-LHC

Fabrizio Ferro1

on behalf of CMS Collaboration

1 INFN, Sezione di Genova (Italy)

*<EMAIL_ADDRESS>
Proceedings for the XXVIII International Workshop
on Deep-Inelastic Scattering and Related Subjects,
Stony Brook University, New York, USA, 12-16 April 2021
10.21468/SciPostPhysProc.?

## Abstract

The Precision Proton Spectrometer (PPS) started operating in 2016 and has
collected more than 110 fb-1 of data over the course of the LHC Run 2, now
fully available for physics analysis. This contribution covers the key
features of the PPS alignment and optics calibration, which have been
developed from scratch. The reconstructed proton distributions, the
performance of the PPS simulation and finally the validation of the full
reconstruction chain with physics data (dilepton events) are also illustrated.

###### Contents

  1. 1 Introduction
  2. 2 Track reconstruction
  3. 3 Alignment
  4. 4 Optics calibration
  5. 5 Proton reconstruction
  6. 6 Future of PPS
  7. 7 Conclusion

##  1 Introduction

The Precision Proton Spectrometer (PPS) [1] detector system has been installed
and integrated into the CMS [2] experiment during the LHC Run 2. It is a joint
project of the CMS and TOTEM [3] collaborations with the capability to measure
protons scattered at very small angles, and to operate at high instantaneous
luminosity. The scattered protons remain inside the beam pipe, displaced from
the central beam orbit, and can be measured by detectors placed inside movable
beam pipe insertions, called Roman Pots (RP), which approach the beam down to
a few mm. The PPS detectors have collected data corresponding to an integrated
luminosity larger than 110 fb-1 during the LHC Run 2, between 2016 and 2018.

The physics motivation behind the PPS detector is the study of central
exclusive production (CEP), i.e. the process $pp\rightarrow p^{(*)}+X+p^{(*)}$
mediated by color-singlet exchanges (photons, Pomerons, Z bosons), with at
least one of the outgoing protons detected in PPS.

Figure 1 shows the layout of the RP system installed at around 200–220 m from
the CMS interaction point (IP5), along the beam line in LHC sector 56. A
symmetric set of detectors is installed in LHC sector 45 at the opposite side
of CMS. The RPs approach the beam vertically from the top and bottom, as well
as horizontally. During standard machine operation, scattered protons undergo
a large displacement in the horizontal direction and a small vertical
displacement at the RP positions. The horizontal RPs are hence used. The
vertical RPs are used in special configurations of the machine and in low
luminosity proton-proton fills for the calibration and alignment of the
detectors.

Figure 1: Experimental setup of PPS in CMS.

Each detector arm consists of two RP stations instrumented with silicon
tracking detectors that measure the transverse displacement of protons with
respect to the beam, and one RP station with timing detectors to measure their
time-of-flight. Silicon strip sensors with a reduced insensitive region on the
edge facing the beam were initially used (10 planes per pot), but since they
cannot sustain a large radiation dose and cannot identify multiple tracks in
the same event, they have been gradually replaced by new 3D silicon pixel
sensors (6 planes per pot).

The proton reconstruction relies on the alignment of the detector planes with
respect to the LHC beam and on the knowledge of the transport matrices
parametrizing the LHC magnet lattice. The latter, referred to as the beam
optics, relate the transverse position and direction of a proton track along
the beam line to the transverse position of the proton at the nominal
interaction point ($x^{*},y^{*}$), the proton scattering angles in the
horizontal and vertical directions ($\theta_{x}^{*},\theta_{y}^{*}$) and its
momentum p, or rather its fractional momentum loss $\xi=(p_{nom}-p)/p_{nom}$,
where $p_{nom}$ is the nominal beam momentum.

##  2 Track reconstruction

The Data Acquisition System of the detectors provides as input of the
reconstruction chain the coordinates of the pixels that have collected a
charge above a certain threshold, set by a calibration procedure, and a
measurement of the collected charge in terms of ADC counts.

The set of pixel fired in the event undergoes a clusterization procedure that
collects in sub-sets, called clusters, the adjacent fired pixels. From any
cluster a single local reconstructed hit is generated, whose position is
evaluated by means of a charge sharing algorithm.

The track reconstruction is carried out on a single pot basis, using the set
of hits provided by the six planes. The procedure is iterative: it starts by
fitting all the hit combinations from the maximum number of planes available
in a region of interest (usually 6) and choosing the one that provides the
minimum $\chi^{2}$ over degree of freedom. If a track candidate with a
$\chi^{2}$ under a certain threshold is found, the track is saved, the used
hits removed and the iterative procedure restarted.

##  3 Alignment

The alignment procedure involves multiple steps. A special alignment
calibration LHC fill allows to determine the absolute position of the RPs with
respect to the beam. This calibration then serves as a reference for the
alignment of every physics fill with standard conditions. Once the tracking
RPs are aligned with respect the beam, the timing RPs are aligned with respect
to the tracking RPs. The alignment procedures depend on the LHC settings
(optics, collimators etc.). Therefore separate analyses are performed for
periods with different settings. The settings are often changed during
technical stops (TS), which are also used to adjust the RP configuration,
possibly influencing the alignment results.

The alignment fill is a calibration fill with standard optics but low beam
intensity which serves multiple purposes. Primarily, it allows to establish
the RP position with respect to the LHC collimators in a procedure analogous
to the LHC collimator alignment. This is a precondition for systematic RP
insertion close to the high-intensity LHC beams. Then, thanks to the low
intensity, one can insert both horizontal and vertical RPs very close to the
beam, at about 5$\sigma_{beam}$, where $\sigma_{beam}$ is the transverse size
of the beam. This yields an overlapping system of detectors that allows for
relative alignment among the RPs in each arm. The relative alignment, among
sensor planes in all RPs and among all RPs in one arm, is determined by
minimising residuals between hits and fitted tracks.

For each high-luminosity LHC fill (physics fill), the horizontal alignment is
obtained by matching observations from the fill to those from the reference
alignment fill. The chosen matching metrics is the slope of the profile of
$(y_{F}-y_{N})$ as a function of $x$, where $y_{N}$ and $y_{F}$ stand for the
vertical track positions in the near and far RP with respect to the IP. The
vertical alignment is obtained by extrapolating the observed vertical profile
of $x_{N}-x_{F}\textrm{ vs. }x_{N}$ to the horizontal beam position.

##  4 Optics calibration

Scattered protons are detected in the PPS RPs after having traversed a segment
of the LHC lattice containing 29 main and corrector magnets per beam. The
transverse proton position at the RP location with transverse position
($x^{*},y^{*}$) and angles ($\theta_{x}^{*},\theta_{y}^{*}$) at the IP is
described at first approximation as

$x=v_{x}(\xi)\cdot x^{*}+L_{x}(\xi)\cdot\theta_{x}^{*}+D_{x}(\xi)\cdot\xi$
$y=v_{y}(\xi)\cdot y^{*}+L_{y}(\xi)\cdot\theta_{y}^{*}+D_{y}(\xi)\cdot\xi.$

Therefore the reconstruction of the proton kinematics at the IP from the
measurements at the RPs needs a precise knowledge and calibration of the LHC
optics.

The first step of the optics calibration is to update the nominal transport
model using LHC databases. The proper version of the LHC magnet lattice
description for the data taking period under consideration is used, and the
nominal magnet strength file for a given beam optics is always updated using
measured data. For the evaluation of the proton kinematics, the horizontal
dispersion $D_{x}$ is the most important, because it allows to convert the
x-coordinate of the proton impact point into the fractional proton momentum
loss $\xi$. The determination of $D_{x}$, using proton trajectories measured
in the RPs, is described in details in [4]. The 2017 and 2018 optics
calibration procedure goes a step further and exploits (semi)-exclusive
$\mu\mu$ production in which the exclusivity of the process plays a key role.

In a third step the vertical dispersion $D_{y}$ is determined from minimum
bias RP data. The calibration of the dispersion functions is followed by the
calibration of the remaining optical functions in the transport equations,
namely the horizontal and vertical effective lengths $L_{x,y}(\xi)$, and the
corresponding magnification functions $v_{x,y}$; other optical functions are
less relevant from the point of view of the proton reconstruction.

##  5 Proton reconstruction

The proton reconstruction consists of back-propagating the protons from the
RPs, where they are measured, to the IP, where the kinematics is determined.
The propagation follows the LHC optics; inputs to it are the proton tracks
detected by the RPs and aligned with respect to the LHC beam. Two
complementary reconstruction strategies are used: single-RP and multi-RP.

The single-RP reconstruction is a simple approach that uses information from
single RPs only. Because of the reduced input information, only $x$ and
$\theta_{y}^{*}$ can be estimated via

$\xi\sim x/D_{x},\hskip 28.45274pt\theta_{y}^{*}=y/L_{y}(\xi).$

These equations reflect only the leading terms from the optics decomposition,
which implies a degraded resolution. On the other hand, notable advantages of
this approach are the applicability even when the proton track is not
available in the other RP(s) of the arm and a slightly smaller systematic
uncertainties with respect to the multi-RP method.

Figure 2: Comparison between single and multi-RP reconstruction using
simulated data: resolution (left) and systematic uncertainties (right).

The multi-RP reconstruction aims at exploiting the full potential of the
spectrometer: it searches for proton kinematics that matches best the
observations from all RPs and all projections by minimising the following
function:

$\chi^{2}=\sum_{i:RPs}\sum_{q:x,y}\left(\frac{d^{i}_{q}-(T^{i}q^{*})_{q}}{\sigma_{q}^{i}}\right)^{2}$

where the vector $d^{i}$ is the proton position at the i-th RP,
$\sigma_{q}^{i}$ its uncertainties, the vector $d^{*}$ denotes the proton
kinematics at the IP and the matrix $T^{i}$ is for the proton transport
between the IP and the i-th RP. The performance differences between the two
methods are summarized in Figure 2 where simulation studies are reported.

##  6 Future of PPS

PPS is preparing the installation of new detectors for LHC Run 3 that is
scheduled to start at the beginning of 2022. With respect to the setup used in
Run 2 an additional RP for timing measurements has been installed on both
sides. Moreover, another remarkable change is the installation of motors that
will allow the vertical movement of the pixel detectors, independently of the
RP position, in order to change the detector part facing the beam (that
undergoes most of the radiation damage), extending therefore the overall
duration of the sensor.

An Expression of Interest for a proton spectrometer to be built in the High
Luminosity phase of LHC (HL-LHC) has been published [5]. It would extend the
sensitivity of the spectrometer to CEP to lower masses and to lower production
cross-sections both in Standard Model measurements and in searches for new
physics.

##  7 Conclusion

The PPS detector of CMS has taken data during Run 2 and has provided the
kinematic parameters of the tagged proton by means of the reconstruction
strategy shortly described in this article. The experience achieved so far
will be exploited in Run 3 where PPS will take data with an updated version of
its detectors. An outlook on the longer term future of proton spectroscopy at
LHC has been given as well.

## References

  * [1] CMS and TOTEM Collaborations, CMS-TOTEM Precision Proton Spectrometer, CERN-LHCC-2014-021, TOTEM-TDR-003, CMS-TDR-13, http://cds.cern.ch/record/1753795
  * [2] CMS Collaboration, The CMS experiment at the CERN LHC, JINST 3 S08004 (2008), 10.1088/1748-0221/3/08/S08004
  * [3] TOTEM Collaboration, TOTEM TDR, CERN-LHCC-2004-002, TOTEM-TDR-001
  * [4] TOTEM Collaboration, LHC Optics Measurement with Proton Tracks Detected by the Roman Pots of the TOTEM Experiment, New J. Phys. 16 (2014) 103041, 10.1088/1367-2630/16/10/103041
  * [5] CMS Collaboration, Expression of Interest for a PPS at HL-LHC, CMS-NOTE-2020/008 https://cds.cern.ch/record/2750358