# Revisiting the effect of lens mass models in cosmological applications of
strong gravitational lensing

Christopher Harvey-Hawes, David L. Wiltshire
School of Physical & Chemical Sciences, University of Canterbury, Private Bag
4800, Christchurch 8041, New Zealand
E-mail<EMAIL_ADDRESS>0009-0002-6751-2695
0000-0003-1992-6682

(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Strong gravitational lens system catalogues are typically used to constrain a
combination of cosmological and empirical lens mass model parameters, even
though the simplest singular isothermal sphere (SIS) models yield a $\chi^{2}$
per degree of freedom $\simeq 2$. To date, this problem has been alleviated by
introducing additional empirical parameters in the extended power law (EPL)
models and constraints from high resolution imagery. The EPL parameters are
taken to vary from lens to lens, rather than defining universal lens profiles.
We investigate these lens models using Bayesian methods through a novel
alternative that treats spatial curvature via the
non–Friedmann–Lemaître–Robertson–Walker (non–FLRW) Timescape cosmology. We
apply Markov Chain Monte Carlo methods using the catalogue of 161 lens systems
of Chen et al. (2019) to simulate large mock catalogues for: (i) the standard
$\Lambda$CDM model with zero spatial curvature; and (ii) the Timescape model.
Furthermore, this methodology can be applied to any cosmological model. In
agreement with previous results we find that in combination with SIS
parameters, models with zero FLRW spatial curvature fit better as the free
parameter approaches an unphysical empty universe, $\Omega_{\mathrm{M}0}\to
0$. By contrast, the Timescape cosmology is found to prefer parameter values
in which its cosmological parameter, the present void fraction, is driven to
$f_{\mathrm{v}0}\to 0.73$ matches, close to values found to best fit
independent cosmological data sets: supernovae Ia distances and cosmic
microwave background. This conclusion holds for a large range of seed values
$f_{\mathrm{v}0}\to 0.73\in\\{0.1,0.9\\}$, and for Timescape fits to both
Timescape and FLRW mocks. Regardless of cosmology, unphysical estimates of the
distance ratios given from power-law lens models result in poor goodness of
fit. Nonetheless, the results are consistent with non-FLRW spatial curvature
evolution.

###### keywords:

gravitational lensing: strong – cosmology: cosmological parameters – galaxies:
haloes

††pubyear: 2023††pagerange: Revisiting the effect of lens mass models in
cosmological applications of strong gravitational lensing–A

##  1 Introduction

Strong gravitational lensing (SGL) occurs when the light path from a distant
source is warped around an intermediary lens body on its way to an observer,
producing several images, arcs or rings. Over the last couple of decades,
strong gravitational lenses have been observed with an increased number. The
current largest compilation of all strong gravitational lenses which can be
used viably for constraining lens properties and cosmological parameters
consists of 161 systems (Chen et al., 2019). The number of such strong
gravitational lenses is set to increase by orders of magnitude over the coming
decade: the Rubin Observatory Large Synoptic Survey Telescope (LSST) and
Euclid are projected to observe $1.2\times 10^{6}$ and $1.7\times 10^{6}$
galaxy-galaxy scale lenses respectively (Collett, 2015). Since the standard
cosmology faces increasing challenges (Buchert et al., 2016; Peebles, 2022;
Aluri et al., 2023) it is important to understand the interplay of lens models
and cosmological models. In this paper, we will test such assumptions using
the existing Chen et al. (2019) catalogue.

It is possible to constrain the properties of the lens galaxies given an
assumed background cosmology, typically using the standard spatially flat
$\Lambda$ Cold Dark Matter ($\Lambda$CDM ) model with Planck value parameters
(Aghanim et al., 2020). Alternatively, one can constrain the parameters of the
assumed cosmological model given a particular choice of lens mass density
profile. Attempts to simultaneously fit both the lens density profile and
cosmological parameters have been conducted (Cao et al., 2015; Chen et al.,
2019), inevitably running into degeneracies that require further observational
data to resolve.

There are several different methods to perform statistical cosmological
analysis with strong gravitational lenses. Time-delay cosmography as pioneered
by Refsdal (1964) is the most sensitive to cosmological parameters and deals
with each system individually rather than trying to deal with an ‘average’
ideal lens model. This method is particularly useful in determining the Hubble
constant, $H_{0}\equiv H(t_{0})$, as the measured time-delay between multiply
imaged events can be used to determine their different path lengths, and
thereby a value for the Hubble expansion (Wong et al., 2020; Millon et al.,
2020). Despite time-delay cosmography being a powerful technique, only a small
fraction of the observed lensing systems have well-defined time-delays, thus
limiting the potential application of this approach. The vast increase in
available lens systems from the next generation of telescopes is now set to
overcome this hurdle.

To date most observed galaxy-scale lensing systems involve distant quasars
sources. As compared to observations of smaller numbers of lensed SneIa, the
uncertainties in time-delay measurements of quasars are significantly larger.
Given the paucity of data available for accurate time-delay cosmography,
alternate approaches should be considered to make full use of the observed
lens systems.

This paper will focus in particular on the distance ratio

$\frac{d_{ls}}{d_{s}}=\sqrt{1-k\,d_{l}^{2}}-\frac{d_{l}}{d_{s}}\,\sqrt{1-k\,d_{s}^{2}}\quad,$
(1)

introduced in the distance sum rule test of Räsänen et al. (2015). Here
$d_{s}$, $d_{l}$ and $d_{ls}$ denote the dimensionless comoving distances,
$d(z_{l},z_{s})$, from the source to observer, lens to observer and source to
lens respectively. The comoving distance is related to the angular diameter
distances, $D_{A}$, and dimensionless luminosity distance, $D_{L}$ by

$\displaystyle d$
$\displaystyle=\frac{H_{0}}{c}D_{A}(1+z)=\frac{H_{0}}{c}\frac{D_{L}}{1+z}$ (2)
$\displaystyle=\frac{1}{\sqrt{\Omega_{\rm k0}}}\sinh{\sqrt{\Omega_{\rm
k0}}\int_{1}^{1+z}\frac{{\rm d}x}{\sqrt{\Omega_{\Lambda 0}+\Omega_{\rm
k0}\,x^{2}+\Omega_{\rm M0}\,x^{3}}}},$ (3)

where $\Omega_{\rm k0}=-kc^{2}/H_{0}^{2}=1-\Omega_{\rm M0}-\Omega_{\Lambda
0}$, and $\Omega_{\rm M0}$ and $\Omega_{\Lambda 0}$ are the present epoch
matter and cosmological constant density parameters.

The relation (1) provides an exact expression for the ratio of two observed
angular distances in the case of a FLRW metric with constant spatial curvature
$k$. It is then typically used as a self-consistency check of the standard
FLRW cosmology, or to constrain the global spatial curvature of some assumed
FLRW universe (Liao et al., 2017; Qi et al., 2019). However, any cosmological
model – whether FLRW or not – can be substituted on the r.h.s. of (1) provided
the model has a suitable definition of the distance ratio in terms of its
underlying parameters. This then allows a comparison of different cosmological
models through Bayesian inference (Cardone et al., 2016). Such Bayesian
analysis might involve two FLRW models with different priors on cosmological
parameters, or alternatively an FLRW model compared to a non-FLRW model such
as the Timescape (Wiltshire, 2009b).

The distance ratio is independent of the Hubble parameter, and therefore
sidesteps the issue of its current epoch value $H_{0}=H(t_{0})$. This is a
major advantages of this technique, as the Hubble tension continues to pose a
significant challenge to the $\Lambda$CDM model (Riess et al., 2023).
Furthermore, since SGL bypasses the cosmic distance ladder, it directly allows
us to infer cosmological parameters without relying on other astronomical
distance measures and the systematics they depend upon.

While the distance ratio has been used to test the self-consistency of FLRW
models, to date it has not been applied to non-FLRW metrics. This paper will
fill that gap by revisiting the analysis of Chen et al. (2019) and extending
it to include the Timescape model. We will perform a Bayesian comparison
between the spatially flat FLRW model and the Timescape model, and determine
the corresponding best fit cosmological parameters. In addition, we
investigate three different parametrisations of the lens galaxies mass profile
and their effect on cosmological fits. Strong gravitational lensing is heavily
dependent on the choice of lens model, as well as the background cosmology the
systems are embedded within. Therefore, it is vital that lens galaxies are
well modelled in advance of new observations by the next generation of
telescopes.

An overview of the paper is as follows. In Section 2, we will briefly review
the Timescape model. In Section 3 we will discuss the different lens models
used in this paper, as well as provide a description of the catalogue of
lensing systems, identifying the uncertainties involved and model biases.
Section 4 outlines the methodology for determining cosmological parameters
within a given lens model and gives the results of the fitting procedure.
Section 5 presents the results of simulations of mock data generated from the
catalogue data with randomly added noise. In section 6 we discuss our results
and the conclusions that can be drawn.

##  2 Timescape Model Overview

Standard $\Lambda$CDM cosmology, and indeed all FLRW models, are based upon
the cosmological principle that on average the universe is both spatially
homogeneous and isotropic. While some notion of spatial homogeneity and
isotropy certainly holds for whole sky averages of very distant objects, the
universe is inhomogeneous on scales $\lesssim$100$\,h^{-1}\hbox{Mpc}$
(Eisenstein et al., 2005; Sylos Labini et al., 2009; Scrimgeour et al., 2012)
where it is dominated by a cosmic web of over-dense filaments and under-dense
voids. The standard cosmology has been incredibly successful. Nonetheless,
challenges to the $\Lambda$CDM model remain unresolved despite of the growth
of available data and increased observational precision (Buchert et al.
(2016); Aluri et al. (2023); Peebles (2022)).

The Timescape Model assumes that cosmology should not invoke a single global
reference background with a unique split of space and time (Wiltshire, 2007a,
b, 2009a). To model the universe effectively, we have to account for
backreaction of inhomogeneities, namely deviations from average FLRW evolution
at the present epoch. This is due to structure formation of filaments and
voids and the resulting cosmic web, which grows in complexity in the late
universe. Timescape makes use of Buchert spatial averages (Buchert, 2001) to
provide the average evolution of the Einstein equations with backreaction.

The interpretation of Buchert’s formalism has been much debated; critics
pointed out that as a fraction of the averaged energy density, a term Buchert
denotes the “kinematical backreaction” $\bar{\Omega}_{\cal Q}$, should be
small (Ishibashi & Wald, 2006). A universal feature of all viable backreaction
proposals is that, unlike FLRW, average spatial curvature does not scale as a
simple power of the average cosmic scale factor $\bar{a}$,
$\bar{\Omega}_{K}\neq-k/(\bar{a}^{2}\bar{H}^{2})$, where $k$ is constant and
$\bar{H}$ an average expansion rate. In the Timescape model, in particular,
$\bar{\Omega}_{K}\propto f_{\rm v}^{1/3}/(\bar{a}^{2}\bar{H}^{2})$, where the
variable void volume fraction, $f_{\rm v}\to 0$ at early times but is
significant in the late epoch universe. This means that conclusions about
spatial curvature derived from CMB data in FLRW models are not relevant for
backreaction models: when fit to the angular diameter distance of the CMB
(with $\Lambda=0$) it is curvature, $\bar{\Omega}_{K}$, that dominates in the
late epoch universe, (Duley et al., 2013, Fig. 1). Furthermore, while the
average energy-density parameters for matter, radiation, curvature and
kinematical backreaction satisfy a simple sum rule

$\bar{\Omega}_{M}+\bar{\Omega}_{R}+\bar{\Omega}_{K}+\bar{\Omega}_{\cal
Q}=1\,,$ (4)

the quantities appearing in (4) are statistical volume average quantities, not
directly related to local observables. Consequently, conclusions about
apparent cosmic acceleration do not follow simply from determining the
magnitude of $\bar{\Omega}_{\cal Q}$. The Timescape ansatz fixes the
relationship of these quantities to local observables, resulting in a viable
phenomenology with $\bar{\Omega}_{K}\,{\mathop{\sim}\limits}\,0.86$,
$\bar{\Omega}_{M}\,{\mathop{\sim}\limits}\,0.17$, $\bar{\Omega}_{\cal
Q}\,{\mathop{\sim}\limits}\,-0.03$ at present, $z=0$.

Since the local geometry varies dramatically depending on the position of an
observer in the cosmic web, local geometry is generally very different to the
global volume average. Local structure not only affects the average evolution
of the Einstein equations, but also the calibration of distance and time
measures relative to the volume average111A volume average observer in the
Timescape model represents an observer whose local spatial curvature coincides
with that of the largest spatial averages. Such an observer is necessarily at
an average position by volume - in a void - and systematically different from
observers in bound structures. A typical location by volume is not a typical
location by mass. observer. Timescape allows for variation of the calibration
of regional clocks relative to one another due to the gravitational energy
cost of the gradients in spatial curvature between filaments and voids. Ideal
observers, who see a close to isotropic CMB, in different local environments
would age differently, hence the name Timescape.

The accelerated expansion of the universe is now an apparent effect derived
from fitting an FLRW model with constant spatial curvature and global cosmic
time to an inhomogeneous non-FLRW universe. An ideal void observer not
gravitationally bound to any structure will not infer a cosmic acceleration at
all, whereas an observer in a bound structure such as galaxy will use a
different time parameter and will infer an accelerated expansion at late
cosmic epochs. The Timescape model can eliminate the need for dark energy
entirely and does so without introducing new ad hoc scalar fields or
modifications to the gravitational action; rather it revisits averaging
procedures that build upon Einstein’s general relativity.

Since the $\Lambda$CDM has been so successful in diverse cosmological tests,
any phenomenologically viable model must yield similar predictions,
particularly with regard to large scale averages for distances beyond the
scale of inhomogeneity. For over a decade, timescape has consistently given
fits to type Ia supernovae which are essentially statistically
indistinguishable from $\Lambda$CDM by Bayesian comparison (Leith et al.,
2007; Smale & Wiltshire, 2011; Dam et al., 2017). The most recent analysis of
the Pantheon+ catalogue shows, furthermore, that Timescape predictions for
cosmic expansion below $\,{\mathop{\sim}\limits}\,100\,h^{-1}\hbox{Mpc}$
scales are consistent with observations in a manner which may provide a self-
consistent resolution of the Hubble tension (Lane et al., 2023). The
difference between the distance–redshift relations of the Timescape model and
those of some $\Lambda$CDM model with fixed $\Omega_{M0}$, $\Omega_{\Lambda
0}$, is only $1$–$3$% over a small range of redshifts, but Timescape
effectively interpolates between $\Lambda$CDM models with different parameters
over larger redshift ranges – which is why it can resolve the Hubble tension
in a natural fashion. Of course, this must be tested with completely
independent cosmological tests.

Among the many different definition of distance in the Timescape model, here
we focus on the matter dominated era, where a simple ‘tracking solution’ is
found to apply. Instead of the parameters $\Omega_{\rm M0}$ in a spatially
flat FLRW model, distances are instead parameterised in terms of the present
epoch void fraction $f_{\rm v0}$. The effective comoving distances directly
comparable to those in the distance sum rule (1) are222Care must be taken with
choices of units. With an explicit factor $c$, $b$ and $\mathcal{F}$ must have
dimensions of $t$ and $t^{1/3}$ respectively, which is different to a
convention often assumed for Timescape (Wiltshire, 2014).

$\begin{cases}d_{s}=H_{0}t_{s}^{2/3}[\mathcal{F}(t_{0})-\mathcal{F}(t_{s})](1+z_{s})\\\
d_{l}=H_{0}t_{l}^{2/3}[\mathcal{F}(t_{0})-\mathcal{F}(t_{l})](1+z_{l})\\\
d_{ls}=H_{0}t_{s}^{2/3}[\mathcal{F}(t_{l})-\mathcal{F}(t_{s})](1+z_{s})\end{cases}$
(5)

where

$\mathcal{F}(t)=2\,t^{1/3}+\frac{b^{1/3}}{6}\ln\left(\frac{(t^{1/3}+b^{1/3})^{2}}{t^{2/3}-b^{1/3}\,t^{1/3}+b^{2/3}}\right)\\\
+\frac{b^{1/3}}{\sqrt{3}}\arctan\left(\frac{2t^{1/3}-b^{1/3}}{\sqrt{3}\,b^{1/3}}\right)\,,$
(6) $b=\frac{2\,(1-f_{\rm v0})(2+f_{\rm v0})}{9\,f_{\rm v0}\,\bar{H}_{0}}\,,$
(7)

and $H_{0}$ and $\bar{H}_{0}$ are the dressed and bare Hubble constants
respectively.333The bare or volume average Hubble parameter is defined by the
average volume expansion - not a global scale factor - and uses a time
parameter relative to an idealised volume average observer. This time
parameter differs systematically from that of observers in bound systems. The
bare Hubble parameter is found to be related to the dressed Hubble parameter
by  $H(t)=\frac{[4f_{\rm v}^{2}(t)+f_{\rm v}(t)+4]\bar{H}(t)}{2(2+f_{\rm
v}(t))}\,.$ and $f_{\rm v0}\equiv f_{\rm v}(t_{0})$ etc denote present epoch
values. The Buchert volume averaged time parameter $t$ has been interpreted in
different ways in different backreaction models. For the Timescape tracker
solution observers in gravitationally bound systems measure a time parameter
related to volume-average time according to

$\tau=\frac{2}{3}t+\frac{2(1-f_{\rm v0})(2+f_{\rm v0})}{27f_{\rm
v0}\bar{H}_{0}}\,\ln{\left[1+\frac{9f_{\rm v0}\bar{H}_{0}}{2(1-f_{\rm
v0})(2+f_{\rm v0})}\right]}\,.$ (8)

The volume-average time parameter is related to our observed redshift by

$z+1=\frac{2^{4/3}\,t^{1/3}\,(t+b)}{f_{\rm
v0}^{1/3}\,\bar{H}_{0}\,t\,(2t+3b)^{4/3}}\,,$ (9)

(Larena et al., 2009). The distance ratio is then given as

$\mathcal{D}^{\text{Time}}\equiv\frac{d_{ls}}{d_{s}}=\frac{\mathcal{F}(t_{l})-\mathcal{F}(t_{s})}{\mathcal{F}(t_{0})-\mathcal{F}(t_{s})}\,.$
(10)

Further details on the motivation and derivation of these equations is given
by Wiltshire (2009b, 2014). It is now simple to make comparisons between
Timescape and spatially flat FLRW models, as they have the same number of free
parameters: $f_{v0}$ and $\Omega_{M0}$ are direct analogues due to their
effect on distance measurements.

##  3 Catalogue, Observables and Lens Models

###  3.1 Catalogue Data

The data set used in this investigation is taken from the catalogue of 161
strongly lensed systems from Chen et al. (2019). It is the largest catalogue
compiled to date, and is itself built upon a previous catalogue of 118 systems
from Cao et al. (2015). This catalogue is a compilation of multiple surveys
consisting of the LSD (Treu & Koopmans, 2002, 2004; Koopmans & Treu, 2002,
2003), SL2S (Ruff et al., 2011; Sonnenfeld et al., 2013a; Sonnenfeld et al.,
2013b, 2015), SLACS (Bolton et al., 2008; Auger et al., 2009, 2010), S4TM (Shu
et al., 2017, 2015), BELLS (Brownstein et al., 2012) and BELLS Gallery surveys
(Shu et al., 2016a, b). Each of these surveys has differing instrumental
parameters and thus has different uncertainties in the measurement of the
central velocity dispersions444This is obtained by converting the redshift
dispersion from spectroscopic measurements around the galaxy using the radial
special relativistic formula $z=\sqrt{(c+v)/(c-v)}-1$. of the lens galaxies,
which contributes the largest portion of the uncertainty in this method.

Velocity dispersions of lens galaxies are measured spectroscopically within a
given aperture radius $\theta_{\rm ap}$. However, the shape and size of the
aperture used in each survey differs and thus an aperture correction formula
is applied so that velocity dispersions are determined as if they were found
with a normalized typical circular aperture $\sigma_{0}$ (Jorgensen et al.,
1995). To transform rectangular apertures into an equivalent circular aperture
for comparison, the following formula is used,

$\theta_{\rm ap}\approx 1.025\sqrt{\frac{\theta_{x}\theta_{y}}{\pi}}\,,$ (11)

where $\theta_{x}$ and $\theta_{y}$ are the width and height of the
rectangular aperture respectively. Since aperture radii vary even within the
class of circular apertures, the following normalization is also applied along
the line of sight,

$\sigma_{\parallel}\equiv\sigma_{0}=\sigma_{\text{ap}}\left(\frac{\theta_{\text{eff}}}{2\theta_{\text{ap}}}\right)^{\eta}\,.$
(12)

The value of $\eta$ is taken from Cappellari et al. (2006) of $\eta=-0.066\pm
0.035$ is found empirically via fitting to individual galaxy profiles and is
in agreement with other values determined by Jorgensen et al. (1995) and
Mehlert et al. (2003) of $\eta=-0.04$ and $\eta=-0.06$ respectively.

The uncertainty in $\eta$ is the greatest source of uncertainty in the
distance ratio inferred from gravitational lensing. In addition, (12) contains
a statistical error in $\sigma_{\text{ap}}$ as quoted in the surveys, and a
further systematic error source. The latter is understood as the effect of
interfering matter along the line of sight and whilst the surveys are
restricted to isolated systems, the sytematic error is estimated at 3% in all
cases (Bernardi et al., 2003). The combined total error budget is then given
by,

$\Delta\sigma_{0}^{\text{tot}}=\sqrt{(\Delta\sigma_{0}^{\eta})^{2}+(\Delta\sigma_{0}^{\text{sys}})^{2}+(\Delta\sigma_{0}^{\text{stat}})^{2}}\,.$
(13)

Further cuts were also applied in compiling the catalogue. Not only did
systems have to be isolated and have no significant substructure, but lens
galaxies must also be early-type galaxies, either elliptical (E) or lenticular
(S0).

###  3.2 Lens Models

When modelling lens mass distributions individually for use in time-delay
cosmography, one of two descriptions are used: (i) power-law models that take
into account both baryonic and dark matter simultaneously; and (ii) a
combination of a luminous baryonic mass profile with a Navarro-Frank-White
dark matter halo (Navarro et al., 1996; Wagner, 2020). Constraints on the
value of $H_{0}$ from either class of model are very similar and statistically
consistent with one another (Millon et al., 2020). However, to model a
statistically average lens galaxy within large catalogues, more general power-
law models are used (Koopmans, 2006).

In this section, we review the three power-law models applied to the (Chen et
al., 2019) catalogue, commonly used to model the matter distribution of
elliptical lens galaxies: Extended Power-Law (EPL), Spherical Power-Law (SPL)
and Singular Isothermal Sphere (SIS).

All three lens models discussed in this paper rely upon the assumption that
while galaxies may differ in shape, when taking averages the dominant
component of the matter density is spherically symmetric. However, the angular
structure of a lens galaxies mass distribution is very significant to the
image separation produced. In fact, quadruply imaged sources are not possible
with spherically symmetric mass distributions. Since the luminous matter in
individual galaxies is not spherically symmetric, whether the angular
component of all matter in each lens galaxy can be averaged out in this
fashion is an assumption which can be tested.

Solving the radial Jeans equation, following the procedure given in Appendix
A, the distance ratio for the most general Extended Power Law (EPL) model is
found to be

$\mathcal{D}^{\text{obs}}\equiv\frac{d_{ls}}{d_{s}}=\frac{c^{2}\theta_{E}}{2\sqrt{\pi}\sigma_{0}^{2}}\left(\frac{\theta_{\text{eff}}}{2\theta_{E}}\right)^{2-\gamma}F(\gamma,\delta,\beta)\,,$
(14)

where $F(\gamma,\delta,\beta)$ depends on the total mass density slope,
$\gamma$, the luminous matter denisty slope, $\delta$, and the stellar orbital
anisotropy, $\beta$.

In deriving (14) several physical assumptions have been made, viz.:

  1. 1.

The creation, destruction and collisions of stars are neglected.

  2. 2.

The thin lens approximation applies as the distances from observer to lens and
lens to source are far larger than the width of the lens object itself. A
projected mass density can then be defined in the lens plane that is
responsible for the deflection of light.

  3. 3.

The system is stationary, i.e., $\partial_{t}=0$. The properties of the lens
galaxy vary insignificantly over the short periods in which lensed images are
observed.

  4. 4.

The weak field limit of general relativity applies (linearized gravity) with
an asymptotically flat (Minkowski) background. For the vast majority of
lensing scenarios one assumes $\Phi/c^{2}\ll 1$.

  5. 5.

Angles of deflection due to gravitational lensing are small, so that the small
angles approximation is invoked.

  6. 6.

The Born approximation applies as the deflection angle is small. The
gravitational potential along the deflected and undeflected light paths can be
considered to be approximately the same.

####  3.2.1 Extended Power-Law Model

For the EPL model (Koopmans, 2006), the most nuanced and complex model
considered in this paper, the function $F$ in (14) takes the form

$F(\gamma,\delta,\beta)=\frac{3-\delta}{(\xi-2\beta)(3-\xi)}\,\times\\\
\left[\frac{\Gamma\quantity(\frac{\xi-1}{2})}{\Gamma\quantity(\frac{\xi}{2})}-\beta\frac{\Gamma\quantity(\frac{\xi+1}{2})}{\Gamma\quantity(\frac{\xi+2}{2})}\right]\left[\frac{\Gamma\quantity(\frac{\gamma}{2})\Gamma\quantity(\frac{\delta}{2})}{\Gamma\quantity(\frac{\gamma-1}{2})\Gamma\quantity(\frac{\delta-1}{2})}\right]\,,$
(15)

where $\xi=\gamma+\delta-2$.

Through a specific choice of $\gamma$, $\delta$ and $\beta$, which
parameterise the total mass density profile $\rho$, the luminous mass density
profile $\nu$ and the stellar orbital anisotropy $\beta$,

$\displaystyle\rho(r)=\rho_{\lower 2.0pt\hbox{$\scriptstyle 0$}}(r/r_{\lower
2.0pt\hbox{$\scriptstyle 0$}})^{-\gamma}\,,$ (16)
$\displaystyle\nu(r)=\nu_{\lower 2.0pt\hbox{$\scriptstyle 0$}}(r/r_{\lower
2.0pt\hbox{$\scriptstyle 0$}})^{-\delta}\,,$ (17) $\displaystyle\beta\
=1-\sigma^{2}_{\theta}/\sigma^{2}_{r}\,,$ (18)

simpler models are recovered. Both $\gamma=\gamma_{0}$ and

$\gamma=\gamma_{0}+\gamma_{z}z_{l}+\gamma_{s}\log\tilde{\Sigma}$ (19)

parameterisations are considered for the total matter density profile where
the normalised surface mass density
$\tilde{\Sigma}\propto\sigma_{0}^{2}/R_{\rm eff}$.

####  3.2.2 Spherical Power-Law Model

If both $\delta=2$ and $\beta=0$ are fixed, the SPL model is obtained. In this
spherically symmetric model, the only free parameter is the total mass density
(dark matter halo) profile, $\gamma$, and

$F:=F(\gamma,2,0)=\frac{1}{\sqrt{\pi}}\frac{1}{\gamma(3-\gamma)}\,.$ (20)

The SPL model is often generalized further to allow for the variation of
$\gamma$ with redshift or surface brightness density of the lens galaxy. For
our investigation, we limit $\gamma$ to a constant value, $\gamma_{0}$.
However, Chen et al. (2019) explore different parametrizations of $\gamma$,
which are discussed in Section 4.1.

####  3.2.3 Singular Isothermal Sphere Model

The SIS model is the simplest discussed, where $\gamma=\delta=2$, $\beta=0$,
and

$F:=F(2,2,0)=\frac{1}{2\sqrt{\pi}}\,.$ (21)

It assumes that the halo mass of the lens galaxy is in isothermal equilibrium,
with gravitational attraction balanced entirely by the pressure associated
with the constituent stars’ interactions. The total mass density relation is
then

$\rho(r)=\frac{\sigma^{2}}{2\pi Gr^{2}}\,,$ (22)

which leads to the distance ratio,

$\frac{d_{ls}}{d_{s}}=\frac{c^{2}\theta_{E}}{4\pi\sigma_{0}^{2}}\,.$ (23)

##  4 Parameter Determination and Model Preference Results

Once the purely cosmology dependent and lens model dependent distance ratios
have been calculated, the defining parameters of each model can be determined
through a Markov Chain Monte Carlo (MCMC) sampling procedure. We initially
take wide uniform priors for all parameters to avoid biases in the values
determined. However, we later take the luminosity density profile $\delta$ to
have priors informed by profiles fitted to Hubble Space Telescope imagery,
giving $2.003<\delta<2.343$. The stellar orbital anisotropy $\beta$ of all
lens galaxies in the catalogue is not measured, so we assume a $2\sigma$
Gaussian prior $\beta=0.18\pm 0.26$, established from a study of nearby
elliptical galaxies and adopted by Chen et al. (2019).

To distinguish between the spatially flat FLRW and Timescape models, which
have the same number of free parameters, a straightforward Bayes factor, $B$,
can be calculated by integrating the likelihoods $\mathcal{L}\sim
e^{-\chi^{2}/2}$ over the $2\sigma$ priors $f_{\rm v0}\in(0.5,0.799)$ and
$\Omega_{\rm M0}\in(0.143,0.487)$ (Dam et al., 2017; Trotta, 2017). These
priors are determined using non-parametric fits of the Planck CMB data
(Aghamousa & Shafieloo, 2015) constrained by estimates of the angular BAO
scale from the BOSS survey data and Lyman alpha forest statistics (Alam et
al., 2017; Delubac et al., 2015). Following Dam et al. (2017) we use wide
priors for both models so as not to unfairly disadvantage $\Lambda$CDM . The
priors usually adopted for $\Lambda$CDM, including a precise estimate
$\Omega_{\rm M0}=0.315\pm 0.007$ (Aghanim et al., 2020), arise from
constraints on perturbation theory on an FLRW background. A Timescape
equivalent has yet to be developed to constrain the equivalent free parameter,
$f_{\rm v0}$.

We adopt the Jeffrey’s scale for the Bayes factor

$B=\frac{\int\mathcal{L}^{\rm timescape}\,\differential f_{\rm
v0}\,\differential\gamma\,\differential\delta\,\differential\beta}{\int\mathcal{L}^{\rm
FLRW}\,\differential\Omega_{\rm
M0}\,\differential\gamma\,\differential\delta\,\differential\beta}$ (24)

so that values of $B>1$ will indicate preference for Timescape and $B<1$ for
FLRW. By the standard interpretation, evidence with $|\ln B|<1$ is ‘not worth
more than a bare mention’ (Kass & Raftery, 1995) or ‘inconclusive’ (Trotta,
2007), while $1\leq|\ln B|<3$, $3\leq|\ln B|<5$ and $|\ln B|\geq 5$ indicate
‘positive’, ‘strong’ and ‘very strong’ evidences respectively (Kass & Raftery,
1995).

Bayes Factors (B)
---
Lens Model |  $B$ |  $\absolutevalue{\ln{B}}$
EPL |  0.4928 |  0.7077
SPL |  $9.38\times 10^{-7}$ |  13.9
SIS |  $3.28\times 10^{-3}$ |  5.72
Table 1: Bayes factors (24), where $B<1$ favours the spatially flat FLRW
model.

(a) $k=0$ FLRW

(b) Timescape

Figure 1: Parameter probability distributions of the EPL lens model for both
spatially flat FLRW and Timescape models, constrained using 161 lensing
systems. Dashed lines representing the $2\sigma$ bounds and median value of
MCMC samples.

(a) $k=0$ FLRW

(b) Timescape

Figure 2: Parameter probability distributions of the EPL lens model for both
spatially flat FLRW and Timescape models, constrained using 161 lensing
systems. Priors are unchanged from Fig.1 except that of the luminosity density
slope $\delta$, which are uniform priors from $2.003<\delta<2.343$. These
values are informed by fitting luminosity profiles to Hubble Space Telescope
imagery of local elliptical galaxies. Dashed lines representing the $2\sigma$
bounds and median value of MCMC samples

(a) $k=0$ FLRW

(b) Timescape

Figure 3: Parameter probability distributions of the fully extended EPL lens
model for both spatially flat FLRW and Timescape models, constrained using 161
lensing systems.  FLRW ($k=0$)
---
Lens Model |  $\beta$ |  $\delta$ |  $\gamma$ |  $\Omega_{\rm M0}$
EPL (full parametrisation)* |  $0.0493^{+5.93\times 10^{-7}}_{-5.34\times 10^{-7}}$ |  $2.3430^{+2.11\times 10^{-13}}_{-9.02\times 10^{-13}}$ |  $\gamma_{0}=2.5000^{+4.16\times 10^{-4}}_{-4.16\times 10^{-4}}$ |  $1.0000^{+6.00\times 10^{-3}}_{-6.00\times 10^{-3}}$
|  |  |  $\gamma_{z}=-0.1682^{+6.91\times 10^{-4}}_{-6.91\times 10^{-4}}$ |
|  |  |  $\gamma_{s}=0.1179^{+1.76\times 10^{-4}}_{-1.76\times 10^{-4}}$ |
EPL ($\gamma=\gamma_{0}$) |  $0.2578^{+1.71\times 10^{-2}}_{-2.42\times 10^{-2}}$ |  $2.4646^{+1.63\times 10^{-2}}_{-2.29\times 10^{-2}}$ |  $2.0588^{+1.02\times 10^{-2}}_{-1.02\times 10^{-2}}$ |  $0.0259^{+2.57\times 10^{-2}}_{-2.57\times 10^{-2}}$
EPL* |  $0.0544^{+9.09\times 10^{-3}}_{-2.16\times 10^{-2}}$ |  $2.2708^{+9.01\times 10^{-3}}_{-2.17\times 10^{-2}}$ |  $2.0588^{+3.19\times 10^{-6}}_{-3.19\times 10^{-6}}$ |  $0.0259^{+3.21\times 10^{-6}}_{-3.21\times 10^{-6}}$
SPL |  – |  – |  $1.9953^{+7.30\times 10^{-3}}_{-7.30\times 10^{-3}}$ |  $0.0208^{+3.17\times 10^{-2}}_{-3.17\times 10^{-2}}$
SIS |  – |  – |  – |  $0.0159^{+3.02\times 10^{-2}}_{-3.02\times 10^{-2}}$
Table 2: Ideal parameter values obtained from MCMC sampling of velocity
dispersions obtained using FLRW distances and one of the listed lens models.
Wide uniform priors are taken in general, except for EPL* where priors on
$\delta$ are informed by HST imagery.  Timescape
---
Lens Model |  $\beta$ |  $\delta$ |  $\gamma$ |  $f_{\rm v0}$
EPL (full parametrisation) |  $0.1886^{+8.68\times 10^{-3}}_{-3.17\times 10^{-2}}$ |  $2.3430^{+5.97\times 10^{-12}}_{-1.53\times 10^{-11}}$ |  $\gamma_{0}=1.8870^{+4.31\times 10^{-4}}_{-4.31\times 10^{-4}}$ |  $0.0000^{+7.67\times 10^{-4}}_{-7.67\times 10^{-4}}$
|  |  |  $\gamma_{z}=-0.3019^{+7.93\times 10^{-4}}_{-7.93\times 10^{-4}}$ |
|  |  |  $\gamma_{s}=0.2300^{+5.33\times 10^{-4}}_{-5.33\times 10^{-4}}$ |
EPL ($\gamma=\gamma_{0}$) |  $0.1576^{+8.68\times 10^{-3}}_{-3.17\times 10^{-2}}$ |  $2.4154^{+5.45\times 10^{-3}}_{-2.02\times 10^{-2}}$ |  $2.0073^{+5.37\times 10^{-3}}_{-5.37\times 10^{-3}}$ |  $0.7685^{+2.38\times 10^{-2}}_{-2.38\times 10^{-2}}$
EPL* |  $-0.0362^{+1.08\times 10^{-2}}_{-5.82\times 10^{-3}}$ |  $2.2836^{+7.92\times 10^{-3}}_{-4.31\times 10^{-3}}$ |  $2.0073^{+6.52\times 10^{-6}}_{-6.52\times 10^{-6}}$ |  $0.7685^{+1.54\times 10^{-5}}_{-1.54\times 10^{-5}}$
SPL |  – |  – |  $1.8837^{+6.20\times 10^{-3}}_{-6.20\times 10^{-3}}$ |  $0.6751^{+1.73\times 10^{-2}}_{-1.73\times 10^{-2}}$
SIS |  – |  – |  – |  $0.7037^{+1.66\times 10^{-2}}_{-1.66\times 10^{-2}}$
Table 3: Ideal parameter values obtained from MCMC sampling of velocity
dispersions obtained using Timescape distances and one of the listed lens
models. Wide uniform priors are taken in general, except for EPL* where priors
on $\delta$ are informed by HST imagery.

###  4.1 Discussion of Results

The Bayes factors favour a spatially flat FLRW model over the Timescape model
in all cases but with varying degrees of strength, as shown in Table 1. For
the SIS and SPL lens models, the preference is a strong. However, both FLRW
and Timescape have a minimum $\chi^{2}$ per degree of freedom $\sim 2$, which
shows both cases give a very poor fit to the data, indicating a potential
problem with the lens models rather than the cosmology. Consequently, the
Bayes factors of Table 1 should not be interpreted naïvely. It is important to
note that the lower $\chi^{2}$ for FLRW models comes at the expense of an
unphysical matter density at the extremes $\Omega_{\rm M0}\simeq 0$ or
$\Omega_{\rm M0}\simeq 1$. By contrast, the values of $f_{\rm v0}$ predicted
are within the $2\sigma$ priors $f_{\rm v0}\in(0.5,0.799)$ for Timescape and
remain physically plausible for the majority of lens model parametrisations.

Chen et al. (2019) already noted that $\Omega_{\rm M0}\simeq 0$ for specific
lens models. However, by considering

  1. 1.

an alternate parametrisation of $\gamma$ to include a dependency on redshift
and normalised surface mass density of each lens galaxy,
$\gamma=\gamma_{0}+\gamma_{z}z_{l}+\gamma_{s}\log{\tilde{\Sigma}}$; and

  2. 2.

$\delta$ as an observable for each lens galaxy

more realistic values of $\Omega_{\rm M0}$ are inferred. One can only fit
$\delta$ given high-resolution imagery of a lens galaxy; for the Chen et al.
(2019) catalogue this reduces the sample of lensing systems from 161 to 130.
Future surveys are predicted to observe several orders of magnitude more
lensing systems. The requirement of follow up high resolution imaging will
therefore face significantly greater challenges on account of the vast
increase of data volume.

Physically plausible values of $\Omega_{\rm M0}$ for FLRW have only been found
when steps (i) and (ii) are applied to the further parametrised EPL model. If
$\delta$ is not constrained for each system in the EPL model, then
$\Omega_{\rm M0}\simeq 1$, the other unphysical extreme. Timescape also
returns the physically implausible extreme value for $f_{\rm v0}\simeq 0$
which corresponds exactly to the Einstein de Sitter universe ($\Omega_{\rm
M0}=1$) for the full paramaterisation of $\gamma$. Our values of
$\gamma_{0}$,$\gamma_{z}$ and $\gamma_{s}$ differ from those of Chen et al.
(2019) due to different normalisations on $\tilde{\Sigma}$.

The values of $\delta$ determined are far higher than values expected from
observations of nearby elliptical galaxies for both Timescape and FLRW models,
when wide uniform priors, $1.2<\delta<2.8$, are applied. This is shown in Fig
2, where $\delta>2.4$ are obtained. When priors informed by observations are
applied, $2.003<\delta<2.343$, $\gamma$ and $\\{\Omega_{\rm M0},f_{\rm v0}\\}$
values are essentially unchanged (see Table 2 and 3).

##  5 Mock Catalogues and Parameter Fitting

###  5.1 Methodology

The lens galaxy mass distribution is fixed to be the simplest case of an SIS
model (22). We use the catalogue data and some initial seed parameters
($\Omega_{\rm M0}$ and $f_{\rm v0}$) to generate a mock catalogue. We then
determine whether the original parameters are recovered from the fitting
procedure. In fact, we find that regardless of the initial input values of the
cosmological parameters, this procedure always lowers $\Omega_{\rm M0}$
significantly for FLRW. For Timescape, a value of $f_{\rm v0}\simeq 0.73$ is
returned irrespective of the input value.

The cosmological distance ratios (i.e., $\mathcal{D}^{\text{FLRW}}$ and
$\mathcal{D}^{\text{Time}}$) are calculated using the catalogue of lens and
source redshifts with seed cosmological parameters in the range $\Omega_{\rm
M0}\in\\{0.1,0.9\\}$ for FLRW, and similarly $f_{\rm v0}\in\\{0.1,0.9\\}$ for
Timescape. From these distance ratios, the ‘ideal’ velocity dispersion for
each system is found using the Einstein radius in (23). To generate the
distance ratios of a mock catalogue, we then apply (23) for the new parameters
drawn from each distribution. Each of the input distributions is assumed to be
Gaussian with standard deviation (13), but with the systematic uncertainty
removed.

The $\chi^{2}$ value is then minimized to determine the corresponding ideal
cosmological parameters for each mock, using

$\chi^{2}=\sum^{161}_{i=1}\quantity(\frac{\mathcal{D}_{i}^{\rm
model}-\mathcal{D}_{i}^{\rm mock}}{\Delta\mathcal{D}_{i}})^{2}\,.$ (25)

To remove the effect of statistical outliers we created $10^{4}$ mock samples,
with the best fit parameters being binned into a histogram for each individual
mock. However, no significant convergence is observed between $10^{3}$ and
$10^{4}$ mocks. For this reason we anticipate generating $10^{3}$ mocks will
be sufficient with larger data sets in future.

###  5.2 Simulation Results and Discussion

Figure 4: Kernel density estimate plots of the best fit parameters for: (a)
Timescape fit to Timescape data, (b) FLRW fit to Timescape data, (c) FLRW fit
to FLRW data and (d) Timescape fit to FLRW data. The dashed vertical red line
indicates the initial parameter value used to generate the mock data, $f_{\rm
v0}=0.76$ and $\Omega_{\rm M0}=0.315$ for Timescape and FLRW models
respectively.

(a) FLRW

(b) Timescape

Figure 5: Kernel density estimate plots showing the best fit cosmological
parameters for both FLRW (a) and Timescape (b) models fit against data
generated with varying initial parameter values.

The histograms produced from the mock data provide interesting results, shown
in Fig 4. Naïvely, one would assume that if the cosmological model plays a
significant role in the pipeline then the generated cosmological parameters
would match the input values. However, this is only the case when the
uncertainty in the velocity dispersion is set artificially low.

The Timescape fit finds values of $f_{\rm v0}\approx 0.73$ which are
consistent with the model expectation based on constraints from SneIa
distances, CMB, etc (Duley et al., 2013). In fact, this is still true for
Timescape fits of the FLRW mocks. The FLRW mocks have a significantly wider
cosmological parameter spread, as shown in Fig. 5(a) as compared to Fig. 5(b).
Furthermore, the FLRW model yields an $\Omega_{\rm M0}$ more in-line with the
expected value of $\Omega_{\rm M0}$ from the Planck CMB data (Aghanim et al.,
2020) when it is fit against data generated with the Timescape model. For FLRW
the maximum likelihood of $\Omega_{\rm M0}$ is found to increase as the input
$\Omega_{\rm M0}$ is lowered, resulting in an overall maximum likelihood
$\Omega_{\rm M0}\to 0$, corresponding to an unphysical Milne universe. For
Timescape the maximum likelihood increases for seed values close to $f_{\rm
v0}\simeq 0.73$.

Both the data and mock catalogues produce some outlying unphysical distance
ratios $\mathcal{D}>1$ with closer source than lens to observer (see Fig. 6).
This is an artefact of choosing a global power-law lens model and fitting to a
data set with high variations in the observed velocity dispersions. Some
systems are not well described by the lens model choice, and thus produce
abnormally high distance ratios that skew the distribution. Fitting any curve
$D(\sigma)$ to the data of Fig. 6 will be skewed by the unphysical values.

The implications can by understood by considering Fig. 7(a) and Fig. 7(b)
where distance ratios $\mathcal{D}^{\text{FLRW}}$ and
$\mathcal{D}^{\text{Time}}$ are shown for particular fixed source and lens
redshifts, $z_{s}$ and $z_{l}$. Varying $z_{s}$ and $z_{l}$ we see that while
the distance ratio generally increases with increasing $z_{s}$, the maximum is
always found at $\Omega_{\rm M0}\simeq 0$ for FLRW and $f_{\rm v0}\simeq 0.73$
for Timescape. With a high enough uncertainty in the velocity dispersion, much
larger distance ratios are possible in the simulations. Therefore, the fitting
of cosmological parameters is forced towards the values that enable the
highest possible distance ratios.

Whilst our simulations were performed using a SIS lens model, the entire class
of lens models in our investigation return unphysical distance ratios,
$\mathcal{D}\geq 1$, when using the ideal parameters determined by MCMC
sampling. Even the EPL model with the full parametrisation (19), which enables
a slight variation of the power law model between systems, often yields
unphysical distance ratios.

Figure 6: Distance ratios against velocity dispersion ($\rm kms^{-1}$) for a
Single Isothermal Sphere lens model using catalogue data for 161 systems.

(a) FLRW

(b) Timescape

Figure 7: Distance ratios for fixed lens and source redshifts in both
spatially flat FLRW and Timescape models with varying $\Omega_{M0}$ and
$f_{v0}$ respectively.

##  6 Conclusions

Astrophysicists outside the SGL community might naturally question why lens
system catalogues are used for cosmology when the simplest SIS models, closest
to the best known dark matter profiles, produce fits with $\chi^{2}$ per
degree of freedom $\approx 2$. All choices of power-law model naïvely prefer a
spatially flat FLRW cosmology, but as a consequence delegate all the mass
along the line of sight to the lens galaxy, resulting in $\Omega_{\rm
M0}\simeq 0$. For the full parametrisation (19) of the EPL model the fits
return an Einstein de Sitter universe with $\Omega_{\rm M0}=1$, which is also
clearly unphysical.

The Timescape cosmology gives reasonable void fraction values of $0.68\leq
f_{\rm v0}\leq 0.77$ in all cases except that of (19) when it also returns the
same unphysical Einstein de Sitter universe ($f_{\rm v0}=0$). This highlights
the imperative of exploring alternatives to the standard $\Lambda$CDM
cosmology. Even if any particular non-FLRW cosmology is incorrect, such
comparisons may give insight into breaking the many degeneracies present in
SGL analyses to date.

It is clear that either

  1. 1.

power-law models cannot be applied globally to large surveys of lensing
systems as they make unphysical predictions of the distance ratios, or,

  2. 2.

our measurements of the velocity dispersions of lens galaxies are not
currently at a high enough accuracy and precision to enable this kind of
analysis.

Navarro-Frenk-White profiles are generally taken as the standard for modelling
the dark matter haloes of elliptical galaxies, rather than the power-law
models used here for large SGL catalogues. The use of a general lens model
across all systems in the catalogue, assuming that on average they will be
well-defined, is evidently not the case. Many lens galaxies have large
uncertainties in their respective velocity dispersion measurements, leading to
unphysical distance ratios for many lensing systems regardless of the choice
of power-law used.

Systems that deviate from these power-law models negatively skew the
distribution of distance ratios, leading to biasing in the determination of
cosmological parameters. In fact, even an individualised EPL model, with the
full parametrisation of $\gamma$, fails to produce physically plausible
distance ratios for all systems. What results is that the ‘best fit’
cosmological parameters found are those which can produce the highest possible
distance ratios to match the unphysical predictions of the lens models.

Without additional imaging data, the technique of using distance ratios as a
test statistic is weak at constraining cosmological parameters. However, it
has the advantage of a far larger dataset than that of using time-delay
distances. With the next generation of telescopes (JWST, Euclid, LSST etc)
predicted to observe even more lensing systems, it is of paramount importance
that appropriate lens models are determined in advance of the upcoming
catalogues if one aims to use strong lensing for robust statistical
constraints in cosmology.

Time-delay cosmography is able to constrain the lens mass distribution with
far greater precision than the distance sum rule test. This arises from the
extra constraints given by the magnification and time-delay between each of
the images. In addition, time-delay cosmography models lensing systems
individually, whereas the distance ratio technique fits global lens models to
large catalogues of highly variable systems. At present, despite the smaller
number of systems available, time-delay cosmography may present the only
viable way to use strong gravitational lensing to constrain both cosmological
models and lens matter models with the precision of measurements available in
ongoing surveys.

## Acknowledgements

We warmly thank Jenny Wagner for her perceptive physical insights in defining
this project, and for introducing us to the concepts, methodology, statistical
and systematic issues that dominate cosmological applications of strong
gravitational lensing. We also thank Zachary Lane, Marco Galoppo, Morag Hills
and Michael Williams for our many discussions and their helpful comments.

## Data Availability

The data used for this analysis was taken from the catalogue of lensing
surveys compiled in Chen et al. (2019).

## References

  * Aghamousa & Shafieloo (2015) Aghamousa A., Shafieloo A., 2015, Journal of Cosmology and Astroparticle Physics, 06, 003
  * Aghanim et al. (2020) Aghanim N., et al., 2020, Astronomy & Astrophysics, 641, A6
  * Alam et al. (2017) Alam S., et al., 2017, Monthly Notices of the Royal Astronomical Society, 470, 2617
  * Aluri et al. (2023) Aluri P. K., et al., 2023, Classical and Quantum Gravity, 40, 094001
  * Auger et al. (2009) Auger M. W., Treu T., Bolton A. S., Gavazzi R., Koopmans L. V. E., Marshall P. J., Bundy K., Moustakas L. A., 2009, Astrophysical Journal, 705, 1099
  * Auger et al. (2010) Auger M. W., Treu T., Bolton A. S., Gavazzi R., Koopmans L. V. E., Marshall P. J., Moustakas L. A., Burles S., 2010, Astrophysical Journal, 724, 511
  * Bernardi et al. (2003) Bernardi M., et al., 2003, The Astronomical Journal, 125, 1817
  * Bolton et al. (2008) Bolton A. S., Burles S., Koopmans L. V. E., Treu T., Gavazzi R., Moustakas L. A., Wayth R., Schlegel D. J., 2008, Astrophysical Journal, 682, 964
  * Brownstein et al. (2012) Brownstein J. R., et al., 2012, Astrophysical Journal, 744, 41
  * Buchert (2001) Buchert T., 2001, General Relativity and Gravitation, 33, 1381
  * Buchert et al. (2016) Buchert T., Coley A. A., Kleinert H., Roukema B. F., Wiltshire D. L., 2016, International Journal of Modern Physics D, 25, 1630007
  * Cao et al. (2015) Cao S., Biesiada M., Gavazzi R., Piórkowska A., Zhu Z.-H., 2015, Astrophysical Journal, 806, 185
  * Cappellari et al. (2006) Cappellari M., et al., 2006, Monthly Notices of the Royal Astronomical Society, 366, 1126
  * Cardone et al. (2016) Cardone V. F., Piedipalumbo E., Scudellaro P., 2016, Monthly Notices of the Royal Astronomical Society, 455, 831
  * Chen et al. (2019) Chen Y., Li R., Shu Y., Cao X., 2019, Monthly Notices of the Royal Astronomical Society, 488, 3745
  * Collett (2015) Collett T. E., 2015, Astrophysical Journal, 811, 20
  * Dam et al. (2017) Dam L. H., Heinesen A., Wiltshire D. L., 2017, Monthly Notices of the Royal Astronomical Society, 472, 835
  * Delubac et al. (2015) Delubac T., et al., 2015, Astronomy & Astrophysics, 574, A59
  * Duley et al. (2013) Duley J. A. G., Nazer M. A., Wiltshire D. L., 2013, Classical and Quantum Gravity, 30, 175006
  * Eisenstein et al. (2005) Eisenstein D. J., et al., 2005, Astrophysical Journal, 633, 560
  * Ishibashi & Wald (2006) Ishibashi A., Wald R. M., 2006, Class. Quant. Grav., 23, 235
  * Jorgensen et al. (1995) Jorgensen I., Franx M., Kjaergaard P., 1995, Monthly Notices of the Royal Astronomical Society, 276, 1341
  * Kass & Raftery (1995) Kass R. E., Raftery A. E., 1995, J. Am. Statist. Assoc., 90, 773
  * Koopmans (2006) Koopmans L. V. E., 2006, EAS Publications Series, 20, 161
  * Koopmans & Treu (2002) Koopmans L. V. E., Treu T., 2002, Astrophysical Journal, 568, L5
  * Koopmans & Treu (2003) Koopmans L. V. E., Treu T., 2003, Astrophysical Journal, 583, 606
  * Lane et al. (2023) Lane Z. G., Seifert A., Ridden-Harper R., Wiltshire D. L., 2023, Cosmological foundations revisited with Pantheon+ (arXiv:2311.01438)
  * Larena et al. (2009) Larena J., Alimi J.-M., Buchert T., Kunz M., Corasaniti P.-S., 2009, Physical Review D, 79, 083011
  * Leith et al. (2007) Leith B. M., Ng S. C. C., Wiltshire D. L., 2007, Astrophysical Journal, 672, L91
  * Liao et al. (2017) Liao K., Li Z., Wang G.-J., Fan X.-L., 2017, Astrophysical Journal, 839, 70
  * Mehlert et al. (2003) Mehlert D., Thomas D., Saglia R. P., Bender R., Wegner G., 2003, Astronomy & Astrophysics, 407, 423
  * Millon et al. (2020) Millon M., et al., 2020, Astronomy & Astrophysics, 639, A101
  * Navarro et al. (1996) Navarro J. F., Frenk C. S., White S. D. M., 1996, Astrophysical Journal, 462, 563
  * Peebles (2022) Peebles P. J. E., 2022, Annals of Physics, 447, 169159
  * Qi et al. (2019) Qi J., Cao S., Biesiada M., Zheng X., Ding X., Zhu Z.-H., 2019, Physical Review D, 100, 023530
  * Räsänen et al. (2015) Räsänen S., Bolejko K., Finoguenov A., 2015, Physical Review Letters, 115, 101301
  * Refsdal (1964) Refsdal S., 1964, Monthly Notices of the Royal Astronomical Society, 128, 307
  * Riess et al. (2023) Riess A. G., et al., 2023, Crowded No More: The Accuracy of the Hubble Constant Tested with High Resolution Observations of Cepheids by JWST (arxiv:2307.15806), doi:10.48550/arXiv.2307.15806
  * Ruff et al. (2011) Ruff A. J., Gavazzi R., Marshall P. J., Treu T., Auger M. W., Brault F., 2011, Astrophysical Journal, 727, 96
  * Scrimgeour et al. (2012) Scrimgeour M. I., et al., 2012, Monthly Notices of the Royal Astronomical Society, 425, 116
  * Shu et al. (2015) Shu Y., et al., 2015, Astrophysical Journal, 803, 71
  * Shu et al. (2016a) Shu Y., et al., 2016a, Astrophysical Journal, 824, 86
  * Shu et al. (2016b) Shu Y., et al., 2016b, Astrophysical Journal, 833, 264
  * Shu et al. (2017) Shu Y., et al., 2017, Astrophysical Journal, 851, 48
  * Smale & Wiltshire (2011) Smale P. R., Wiltshire D. L., 2011, Monthly Notices of the Royal Astronomical Society, 413, 367
  * Sonnenfeld et al. (2013a) Sonnenfeld A., Gavazzi R., Suyu S. H., Treu T., Marshall P. J., 2013a, Astrophysical Journal, 777, 97
  * Sonnenfeld et al. (2013b) Sonnenfeld A., Treu T., Gavazzi R., Suyu S. H., Marshall P. J., Auger M. W., Nipoti C., 2013b, Astrophysical Journal, 777, 98
  * Sonnenfeld et al. (2015) Sonnenfeld A., Treu T., Marshall P. J., Suyu S. H., Gavazzi R., Auger M. W., Nipoti C., 2015, Astrophysical Journal, 800, 94
  * Sylos Labini et al. (2009) Sylos Labini F., Vasilyev N. L., Pietronero L., Baryshev Y. V., 2009, EPL (Europhysics Letters), 86, 49001
  * Treu & Koopmans (2002) Treu T., Koopmans L., 2002, Astrophysical Journal, 575, 87
  * Treu & Koopmans (2004) Treu T., Koopmans L. V. E., 2004, Astrophysical Journal, 611, 739
  * Trotta (2007) Trotta R., 2007, MNRAS, 378, 72
  * Trotta (2017) Trotta R., 2017, Bayesian Methods in Cosmology (arxiv:1701.01467), doi:10.48550/arXiv.1701.01467
  * Wagner (2020) Wagner J., 2020, International Journal of Modern Physics D, 29, 2043017
  * Wiltshire (2007a) Wiltshire D. L., 2007a, New Journal of Physics, 9, 377
  * Wiltshire (2007b) Wiltshire D. L., 2007b, Physical Review Letters, 99, 251101
  * Wiltshire (2009a) Wiltshire D. L., 2009a, International Journal of Modern Physics D, 18, 2121
  * Wiltshire (2009b) Wiltshire D. L., 2009b, Physical Review D, 80, 123512
  * Wiltshire (2014) Wiltshire D. L., 2014, in Perez Bergliaffa S., Novello M., eds, , Proceedings of the XVth Brazilian School of Cosmology and Gravitation. Cambridge Scientific Publishers, pp 203–244 (arXiv:1311.3787)
  * Wong et al. (2020) Wong K. C., et al., 2020, Monthly Notices of the Royal Astronomical Society, 498, 1420

##  Appendix A Model derivation and explanation of parameters

In this section we will discuss the various assumptions that go into the
extended power-law model of Koopmans (2006) for elliptical galaxies as
gravitational lenses. The extended power-law model, as well as variants SPL
and SIS, are derived from the radial Jeans equation here given in spherical
coordinates $(r,\theta,\phi)$:

$\frac{\differential}{\differential
r}[\nu(r)\,\sigma^{2}_{r}]+\frac{2\beta}{r}\,\nu(r)\,\sigma^{2}_{r}=-\nu(r)\frac{\differential\Phi}{\differential
r}\,,$ (26)

where

$\frac{\differential\Phi}{\differential r}=\frac{GM(r)}{r^{2}}$ (27)

Here, the velocity dispersion $\sigma_{r}$ is defined as follows with $f$
being the distribution function of the stars within the galaxy.

$\sigma_{r}^{2}=\frac{1}{\nu}\int v_{r}f\differential^{3}v$ (28)

The Jeans equation looks at the average motion of the distribution of stars
defined by the luminous matter density distribution $\nu(r)$ within a
Newtonian gravitational potential created by the total matter density
$\rho(r)$, including dark matter. The parameter $\beta(r)$ denotes the
anisotropy of the stellar velocity dispersion and is also known as the stellar
orbital anisotropy. These are parametrised as such for the extended power-law,

$\displaystyle\rho(r)=\rho_{\lower 2.0pt\hbox{$\scriptstyle 0$}}(r/r_{\lower
2.0pt\hbox{$\scriptstyle 0$}})^{-\gamma}$ (29)
$\displaystyle\nu(r)=\nu_{\lower 2.0pt\hbox{$\scriptstyle 0$}}(r/r_{\lower
2.0pt\hbox{$\scriptstyle 0$}})^{-\delta}$ (30)
$\displaystyle\beta(r)=1-\sigma^{2}_{\theta}/\sigma^{2}_{r}\,,$ (31)

where $\sigma_{\theta}$ and $\sigma_{r}$ are the tangential and radial
velocity dispersions respectively. With $\beta(r)\neq 0$, we can account the
average orbits of stars deviating from circular paths, however the total mass
density will remain spherical - it does not disrupt the spherical dark matter
halo of the galaxy.

In order to solve the Jeans equation, numerous assumptions are made:

  * •

the system is static i.e., $\partial_{t}=0$. This feels like a reasonable
assumption to make as the structure of the elliptical galaxy should not be
expected to change over the course of lensing measurements.

  * •

the stress tensor $\sigma_{ij}$ is diagonal, such that
$\expectationvalue{v_{i}v_{j}}=0$ if $i\neq j$.

  * •

spherical symmetry applies such that
$\sigma_{r\phi}^{2}=\sigma_{r\theta}^{2}=0$ and $\rho$ and $\nu$ only depend
on $r$ with no angular component.

  * •

collisionless matter. The Jeans equation is originally derived from the
collisionless Boltzman equation and as such does not permit interactions
between the constituent stars within galaxies.

After making these assumptions, the radial Jeans equation can be solved for
the radial velocity dispersion in terms of the luminous matter distribution
and the total mass inside a sphere of radius $r$,

$\sigma^{2}_{r}(r)=\frac{G\int^{\infty}_{r}\differential
r^{\prime}{r^{\prime}}^{2\beta-2}\,\nu(r^{\prime})M(r^{\prime})}{r^{2\beta}\,\nu(r)}\,.$
(32)

We can then define the mass contained within a cylinder of radius $R_{E}$, the
Einstein radius, which quantifies the mass of the lens galaxy,

$M_{E}=\int^{R_{E}}_{0}\differential R\,2\pi R^{\prime}\Sigma(R^{\prime})\,,$
(33)

where $\Sigma(R)$ is the mass density projected into the lens plane,

$\Sigma(R)=\int^{\infty}_{-\infty}\differential
Z\,\rho(r)=\int^{\infty}_{-\infty}\differential Z\,\rho_{\lower
2.0pt\hbox{$\scriptstyle 0$}}\,r_{0}^{\gamma}(Z^{2}+R^{2})^{-\gamma/2}\,,$
(34)
$\Sigma(R)=\sqrt{\pi}\,R^{1-\gamma}\frac{\Gamma((\gamma-1)/2)}{\Gamma(\gamma/2)}\rho_{\lower
2.0pt\hbox{$\scriptstyle 0$}}\,r_{0}^{\gamma}\,.$ (35)

Here $R$ is the radius of the galaxy in the lens plane and $Z$ is the distance
along the line of sight perpendicular to the lens plane, $r$ is the spherical
radius as previous such that $r^{2}=R^{2}+Z^{2}$.

The mass contained inside the Einstein radius is therefore,

$M_{E}=2\pi^{3/2}\,\frac{R_{E}^{3-\gamma}}{3-\gamma}\frac{\Gamma((\gamma-1)/2)}{\Gamma(\gamma/2)}\rho_{\lower
2.0pt\hbox{$\scriptstyle 0$}}\,r_{0}^{\gamma}\,.$ (36)

The total mass within a sphere of radius $r$ can also be calculated,

$M(r)=\int^{r}_{0}\differential
r^{\prime}4\pi\,{r^{\prime}}^{2}\,\rho(r^{\prime})=4\pi\,\rho_{\lower
2.0pt\hbox{$\scriptstyle 0$}}\,r_{0}^{\gamma}\frac{r^{3-\gamma}}{3-\gamma}$
(37)

which can be written in terms of the mass $M_{E}$ contained within the
cyclinder as

$\sigma^{2}_{r}(r)=\frac{2}{\sqrt{\pi}}\frac{GM_{E}}{R_{E}}\frac{1}{\xi-2\beta}\frac{\Gamma((\gamma-1)/2)}{\Gamma(\gamma/2)}\left(\frac{r}{R_{E}}\right)^{2-\gamma}$
(38)

where $\xi=\gamma+\delta-2$. The velocity dispersion from observation is the
component of the luminosity weighted average along the line of sight and pver
the effective spectroscopic aperture $R_{A}$. This can be written as:

$\sigma^{2}_{\parallel}(\leq
R_{A})=\frac{\int^{R_{A}}_{0}\int^{\infty}_{-\infty}\differential
R\,\differential Z\,2\pi
R\,\sigma^{2}_{r}(r)\left(1-\beta\frac{R^{2}}{r^{2}}\right)\nu(r)}{\int^{R_{A}}_{0}\int^{\infty}_{-\infty}\differential
R\,\differential Z\,2\pi R\nu(r)}\,.$ (39)

By substituting in the expression for $\sigma^{2}_{r}$ and the power-law of
$\nu(r)$ we arrive at

$\sigma^{2}_{\parallel}(\leq
R_{A})=\frac{2}{\sqrt{\pi}}\frac{GM_{E}}{R_{E}}\frac{3-\delta}{(\xi-2\beta)\,(3-\xi)}\left[\frac{\Gamma\left(\frac{\xi-1}{2}\right)}{\Gamma\left(\frac{\xi}{2}\right)}-\beta\frac{\Gamma\left(\frac{\xi+1}{2}\right)}{\Gamma\left(\frac{\xi+2}{2}\right)}\right]\\\
\times\frac{\Gamma\left(\frac{\gamma}{2}\right)\Gamma\left(\frac{\delta}{2}\right)}{\Gamma\left(\frac{\gamma-1}{2}\right)\Gamma\left(\frac{\delta-1}{2}\right)}\left(\frac{R_{A}}{R_{E}}\right)^{2-\gamma}\,.$
(40)

Using the deflecting mass found from the lensing equation

$M_{E}=\frac{c^{2}\,\theta_{E}}{4G}\frac{D_{s}\,D_{l}}{D_{ls}}\,,$ (41)

$R_{A}=D_{l}\theta_{A}$ and $R_{E}=\theta_{E}D_{l}$ we get,

$\sigma^{2}_{\parallel}(\leq
R_{A})=\frac{2}{\sqrt{\pi}}\frac{D_{s}\,\theta_{E}\,c^{2}}{D_{ls}}\frac{3-\delta}{(\xi-2\beta)(3-\xi)}\left[\frac{\Gamma\left(\frac{\xi-1}{2}\right)}{\Gamma\left(\frac{\xi}{2}\right)}-\beta\frac{\Gamma\left(\frac{\xi+1}{2}\right)}{\Gamma\left(\frac{\xi+2}{2}\right)}\right]\\\
\times\frac{\Gamma\left(\frac{\gamma}{2}\right)\Gamma\left(\frac{\delta}{2}\right)}{\Gamma\left(\frac{\gamma-1}{2}\right)\Gamma\left(\frac{\delta-1}{2}\right)}\left(\frac{\theta_{A}}{\theta_{E}}\right)^{2-\gamma}\,,$
(42)

which is the expression we call the extended power-law. The use of the lens
equation also is built upon several assumptions:

  * •

the weak field limit of general relativity applies and that $\Phi/c^{2}\ll 1$,

  * •

small deflection angles,

  * •

the Born approximation applies.


# What You See Is (not) What You Get: A VR Framework For Correcting Robot
Errors

Maciej K. Wozniak<EMAIL_ADDRESS>KTH Royal Institute of TechnologyStockholm,
Sweden ,  Rebecca Stower<EMAIL_ADDRESS>KTH Royal Institute of
TechnologyStockholm, Sweden ,  Patric Jensfelt<EMAIL_ADDRESS>KTH Royal
Institute of TechnologyStockholm, Sweden and  Andre Pereira<EMAIL_ADDRESS>KTH
Royal Institute of TechnologyStockholm, Sweden

(2023)

###### Abstract.

Many solutions tailored for intuitive visualization or teleoperation of
virtual, augmented and mixed (VAM) reality systems are not robust to robot
failures, such as the inability to detect and recognize objects in the
environment or planning unsafe trajectories. In this paper, we present a novel
virtual reality (VR) framework where users can (i) recognize when the robot
has failed to detect a real-world object, (ii) correct the error in VR, (iii)
modify proposed object trajectories and, (iv) implement behaviors on a real-
world robot. Finally, we propose a user study aimed at testing the efficacy of
our framework. Project materials can be found in the OSF
repository111https://osf.io/c8wap/?view_only=e1c799b564b043d0aeaac289513ebff0.
The code will be added after user studies are completed..

robotics, human-robot interaction, VR, AR, perception

††journalyear: 2023††copyright: rightsretained††conference: Companion of the
2023 ACM/IEEE International Conference on Human-Robot Interaction; March
13–16, 2023; Stockholm, Sweden††booktitle: Companion of the 2023 ACM/IEEE
Int’l Conference on Human-Robot Interaction (HRI ’23 Companion), March 13–16,
2023, Stockholm, Sweden††doi: 10.1145/3568294.3580081††isbn:
978-1-4503-9970-8/23/03††ccs: Human-centered computing Virtual reality

##  1\. Introduction

Robots designed for human-robot-interaction (HRI), such as cobots, can
increasingly be seen in homes (Robinson et al., 2014), offices (Hawes et al.,
2017) or common areas like restaurants or bars (Foster et al., 2014).
Nevertheless, such robots are not perfect and can suffer from many potential
failures. Although deep learning models used in the robot perception modules
can perform well on the datasets they are trained and evaluated on, they often
fail when deployed in the real world (Honig and Oron-Gilad, 2018). Simple
failures, such as not detecting an object, can have severe consequences. First
of all, if the robot fails a simple task, the user may get angry and annoyed,
resulting in lower trust and motivation to collaborate (Lucas et al., 2018;
Khavas et al., 2020). Secondly, a robot may collide or crash with objects if
it fails to detect them.

Another potential issue relates to how the robot plans and executes the
trajectories. Users may have different preferences about how the robot should
complete a specific task. A robot reaching for an object in a collaborative
task (e.g., joint assembly) may be perceived as unsafe from the user’s
perspective, even if the movement is executed successfully from the robot
point of view (i.e., no collision occurs). Allowing users to see the robot’s
actions in advance, and potentially modify them, could then help improve the
overall user experience.

In our framework, we therefore consider two main sources of failure that can
occur in human-robot-interactions; failure to detect an object (potentially
causing a collision in the real world), or failure to plan an acceptable
trajectory for the user. Rather than preventing these failures entirely, we
aim to provide users with the opportunity to directly correct the robot
themselves when such failures occur. We provide two main functions for (i)
correcting the robot’s understanding about its surroundings by modifying the
virtual environment created from its perception module (ii) modifying the
robot’s trajectory to adjust to the users’ preferences or to avoid obstacles.
In doing so, we contribute a virtual reality framework for human-robot
collaboration that considers that robots are not perfect and makes it easy for
the user to correct their mistakes.

Figure 1. Example of the framework with Franka Panda robotic arm. From the
left: (1) the VR user interface; (2) user interacting with the robot through
Oculus Quest 2 headset; (3) real-world environment.

##  2\. Related Work

Several recent projects have focused on improving human-robot collaboration
using virtual, augmented and mixed reality (VAM) devices. Teleoperation
(remotely controlling the robot) is one of the most common forms of human-
robot interaction (Nenna and Gamberini, 2022; Naceri et al., 2019; Rosen et
al., 2018). Ostatin (Ostanin and Klimchik, 2018) and Togias (Togias et al.,
2021) both use VR to plan the robotic arm’s trajectory. Others, like Chandan
et al. (Chandan et al., 2021) use AR to visualize the states and intentions of
the robots. Xu et al. (Xu et al., 2022) used VR to steer the robot by moving
its end effector in the VR space. Kennel-Maushart et al. (Kennel-Maushart et
al., 2022) and Barentine et al. (Barentine et al., 2021) also focused on
steering the robot, however, the former did so by applying force in AR by
pressing virtual objects held by the virtual robot, making the real world
robot move. The latter uses VR for steering the robot from a first-person
view. Lastly, Wozniak and Jensfelt (Wozniak and Jensfelt, 2022) presented a VR
simulation framework for interacting with the virtual representation of a
robot and modifying its trajectory.

Many of these projects focus only on the execution of planned robot
trajectories and do not consider the perception module (and associated
potential failures). In addition, the visualization and manipulation of the
robot’s movements are limited to either only choosing the final position of
the robot’s end-effector or adding no-go zones. Finally, none of these
projects follow the whole pipeline of first visualizing, then correcting and
testing the planned trajectory in simulation, followed by deploying it to a
real-world robot.

Our framework could also be classified as a digital twin project, which is a
virtual representation of a robot’s actions and abilities (Liu et al., 2021).
Modern physics engines can imitate reality in great detail, allowing a digital
twin to be an accurate test-bed for the real-world (Malik and Brem, 2021).
There are various applications and benefits of using VR interfaces for digital
twin projects, such as improved factory safety or workers’ training (Wang et
al., 2021; Kaarlela et al., 2020), collaborative tasks (Ortenzi et al., 2022)
or assembly process (Guzzi et al., 2022). However, these methods so far focus
primarily on personnel training for the manufacturing industry.

In addition, what is still missing in these works is an explicit understanding
of the robot’s perception of the working space and planned course of action.
Our framework addresses this challenge by allowing the user to (i) assess how
the robot perceives the environment, (ii) modify and correct the robot’s
actions and understanding of the environment, and (iii) verify the proposed
actions and trajectories before deploying them and testing in the real world.

Sensor dataSensorsMotorsDeep Learning moduleDatatransfernodeVR
environmentUserRobot controller

Figure 2. Proposed system architecture. The elements within the orange box
correspond to nodes and parts connected with the robot’s motion and
perception, whereas the ones in the green box are explicitly corresponding to
its hardware. Nodes inside the blue box are connected to the user and VR
setting.

##  3\. Technical Implementation

Our proposed framework integrates several systems; the cameras and associated
perception module, the VR headset and user interface, and the real-world
robot.

###  3.1. Camera views and perception module

We use two cameras to create and examine the environment: one fixed, pointing
at the table (main camera) and one on the robot’s end effector (secondary
camera). The images detected by the cameras are fed into the deep learning
model (Wu et al., 2019) for scene segmentation and object detection to obtain
segmentation masks and bounding boxes of the detected objects222Any other
object detection architectures could be used as well.. These tell us what the
object is and where it is located. Our system then uses this information to
correctly position and add objects to the scene in VR.

The VR environment is then built from the output of the perception module. The
user can only ask the robot to interact with objects shown in the VR
environment. The user can switch between the VR environment view and the two
distinct camera views to verify whether the robot has successfully detected
all the objects.

In case an object is occluded, the robot’s arm with the secondary camera can
be moved (since from another angle the object may be more clearly visible) to
detect it and add it to the environment. Alternatively, if the object is still
not detectable or the perception module itself fails, the user can add the
missing object to the VR environment by bringing up a menu and adding the
object to the scene (as shown in the supplementary videoFootnote 1).

###  3.2. Robot and data transfer

We have so far tested our framework with the Franka Panda and Niryo One
robotic arms; however our framework can accommodate most robotic arms,
provided the URDF files. We use the MoveIt library333https://moveit.ros.org/
for planning and Robot Operating System (ROS Melodic)444https://www.ros.org/
for transferring data and running different tasks on separate nodes.

Since the user interface (UI) of our framework is in the VR headset, other
operations, such as planning and object detection, must be run on separate
machines. This creates the challenge of efficiently transferring the data
between different nodes of the system and the VR headset itself. In order to
facilitate the exchange of information, improve the user experience, and
minimize the necessary bandwidth, we only transfer the output of a detection
network between the robot and the virtual environment. When the robot detects
and classifies an object, it sends its class, location, and estimated size to
the VR UI application. In the virtual environment, we have multiple prefabs
(3D models of the objects we built into the project) corresponding to the
detected classes, from which we can quickly create a 3D representation of the
environment when the message from the perception module (containing detected
object classes and poses) is received.

###  3.3. VR Interface

The VR interface is built in Unity ($2020.3$ version) and used the Oculus
Quest 2 headset555Our framework can also be adapted to other headsets by
changing the target device while building the Unity project. In our UI design
(shown in Fig. 3(d) and the supplementary videoFootnote 1) the user can bring
up two different menus with primary and secondary buttons on the left
controller and choose an option by pointing at the specific button with the
right-hand controller. The first menu has options for creating and modifying
the environment. The second menu contains commands for the robot (such as
selecting the object or sending trajectories to the real robot).

##  4\. Interaction Scenario

In this section we show an example scenario of our framework with a 7 DoF
Franka Panda robotic arm. While we use a real-world robot, we want to ensure
that the framework is available to labs without access to a physical robot. In
such a case, communication with a real-world robot can be simulated by running
the same Unity project on another machineFootnote 1.

There are three main actions that users take to arrive at real-world object
manipulation with the robot. First, the VR environment depicting both the
robot and the objects detected and recognized by the perception module is
created (and corrected by the users, if necessary). Next, the users ask the
virtual robot to perform an action in VR. If the users are satisfied with how
the virtual robot acts, they can send the trajectory to the real–world robot.
If not, they can modify the trajectory beforehand. All the steps described in
this section can be found in the video in attached materialsFootnote 1.

###  4.1. Common environment understanding

The VR environment corresponds to what the robot understands about the real-
world environment. Objects detected by the robot’s perception module appear in
the UI. Users can examine whether the objects correspond to what is in the
real-world environment in two ways. First, as shown in Fig. 3(c), they can
activate passthrough where they will see the main camera image projected on
the table where the objects are located. Secondly, they can go to camera view,
where they can see and compare the VR environment with the output from both
cameras’ deep learning modules as shown in Fig. 3(a) and Fig. 3(b).

However, the perception module is imperfect and often partially fails, as
described in Section 3.1. There are two ways to recover from these failures.
The users can move the robotic arm around and try to detect the objects using
the camera on the robot’s arm. This can help with, e.g., partially occluded
objects or objects that cannot be recognized from the main camera perspective
(Fig. 3(a) vs. Fig. 3(b)). Another way is to insert virtual objects into the
environment. Users can choose the object from the UI and grab and place it in
the VR environment where it belongs, as shown in Fig. 3(d). The passthrough
view on the table allows the user to see where the objects are in the real
world and place a matching virtual object in the correct position. Objects in
the environment are then considered for collision avoidance while planning the
trajectory.

(a) Main camera view

(b) Secondary camera view

(c) Passthrough view

(d) Adding missing objects to the environment

Figure 3. Different ways to see the real world view in the VR. Users can
verify if the environment was correctly created and correct potential flaws of
the robot’s perception module by comparing these views with the state of the
VR environment.

###  4.2. Trajectory manipulation

When the users are satisfied that the VR environment’s current state is an
accurate representation of reality, they can proceed to ask the virtual robot
to move objects. Whilst the virtual robot is executing a movement, it
generates a trajectory and waypoints in real-time in the VR environment. After
the robot finishes its movement, the user then has the option to either send
the action to the real robot, or bring the virtual robot and the selected
object back to their initial positions and edit the trajectory. While the
planner considers objects in the scene and plans the trajectory so that the
robotic arm will not collide with them (except for the object we want to pick
up), this might not always align with user’s preferences. To accommodate this,
we first show the user the planned trajectory containing multiple waypoints in
VR. If the user is not satisfied with the proposed trajectory, they can edit
it by moving the waypoints. Now, the user can ask the robot to perform the
task again, following the modified trajectory (see video for
demonstrationFootnote 1). This additional functionality could be used to e.g.,
help home robots learn users’ preferences.

###  4.3. Real robot

When the users are satisfied with the state of the environment and the way the
robot performs the task, they can translate the movement to the real robot.
After the real robot finishes the task, they can ask it to start another
interaction in the same manner as before: by planning and verifying it first
in VR and then sending it to the real robot.

##  5\. Proposed User Study

In order to test the efficacy of our framework, we plan to conduct a user
study comparing our VR system with a screen interface using a keyboard and
mouse (hereon referred to as screen interface). In this study, we will focus
only on failures within the perception module, rather than editing robot
trajectories. Our goal is to investigate whether VR systems are preferred to
screen-based systems for perceiving and correcting robot perceptual errors.

As our independent variable, we will manipulate the type of interface
participants use to interact with the robot (VR, screen). We chose the screen
interface as our point of comparison because (i) screen interfaces are often
used to teleoperate robots (Liu et al., 2017; Wonsick and Padır, 2021; Szafir,
2019), and (ii) there are very few direct VR - screen comparisons in the HRI
literature to date (Mara et al., 2021; Li, 2015). We plan to use a within-
groups study design, where all participants will be exposed to both
conditions. Whether participants start with the VR system first or screen
first will be randomly assigned. Before beginning the study, participants will
be given a pre-questionnaire assessing their familiarity with robots and VAM
reality systems.

Participants will first be given a training phase where they can familiarise
themselves with the system setup and controls. In particular, participants
will be shown the two different cameras (main and secondary) with which they
can view the real-world scene, as well as the passthrough option. Participants
will then have 5 experimental trials where they interact with the robot. Each
trial will consist of the robot failing to recognize an object. That is, the
user can see that the object exists in the real world (via the external
cameras/passthrough), but it is not represented in the virtual world (either
in VR for the VR condition, or on the screen in the screen condition). The
user can then correct this error in one of two different ways. First, they can
move the secondary camera in an effort to detect the missing object.
Alternatively, they can manually add the object to the environment by
selecting the missing object from an array of different possible objects and
placing it (on the screen / in VR) as close to the real-life position of the
object as possible. Once the user is satisfied with their object positioning,
they can submit their feedback to the virtual robot, view the planned
trajectory, and finally, execute the movement in the real world (as described
in section 4.2) before moving on to the next trial. Although the overall
framework also allows users to view and edit the planned trajectory, for this
study we will disable this function, as we aim to focus first on how users
correct the robot when an object detection failure occurs.

As dependent variables, we will measure how often users choose each correction
method (moving the secondary camera, or adding virtual objects), the time
taken per trial, and the user’s subjective perceptions of each interface. We
will also measure the accuracy of the object placement by looking at where
users choose to manually add the object to the environment. For the subjective
perceptions, we will use the Technology Acceptance Model (TAM) (Venkatesh and
Davis, 2000) and the performance trust subscale from the Multi-Dimensional
Measure of Trust (MDMT) scale (Malle and Ullman, 2021). Participants will
complete the questionnaires directly after interacting with each system. After
participants have seen both conditions, we will also ask a binary forced-
choice question about which system they preferred to use. Finally, we will
interview participants about the perceived strengths and weaknesses of each
system.

We aim to recruit 50 participants, based on an a-priori power analysis with
$\alpha=0.05,~{}\beta=0.8$ and $d=0.4$. All hypotheses and planned analyses
will be pre-registered on the Open Science FrameworkFootnote 1 prior to data
collection.

###  5.1. Hypotheses

Based on previous studies which show that VR systems lead to better task
performance (Liu et al., 2017) and have higher perceived usability (Whitney et
al., 2020) than 2D interactions, we formulate the following hypotheses:

  * H1

The VR interface will be preferred over the screen interface for interacting
with the robot

  * H2

Participants will have better task performance with the VR system in terms of:

    * 2a

More accurate object placement when adding virtual objects to the environment

    * 2b

Shorter average time taken per trial

  * H3

Subjective ratings on the self–report questionnaires (technology acceptance
and trust) will be higher for the VR system compared to the screen interface.

##  6\. Conclusions and future work

This paper presents a virtual reality human-robot collaboration framework for
object detection and interaction accounting for robot perceptual failures.
Although we are interested in using our framework to enhance human robot
collaboration, it could also be used for other tasks.

As briefly mentioned in Section 4.2, trajectory manipulation can be adapted to
the preference learning task. While future home robots might be able to
successfully perform many different tasks, the users’ satisfaction will vary.
For example, different groups of users might have different preferences about
how to do certain things. People’s level of comfort being close to a moving
robot may also vary. A similar tool to what we present can be used in the
transition when a user buys a new robot and wants to collect data on how users
correct a robot’s trajectory so that the system learns to perform tasks in the
way that this particular user likes.

The framework can also be used to collect data and improve deep learning
models responsible for robot perception. When users realize a specific object
is not detected, they could, e.g., draw a bounding box around that object or
even just place it in the correct place in the VR, triggering bounding box
generation. Then, such images could be collected and the model retrained to
perform better on the tasks it failed before. Such a solution can help the
users improve the robot’s functionality, even if they do not have any
programming experience.

In sum, in this paper we presented a technical implementation and interaction
scenario of a VR framework for correcting robot perception and planning
errors, as well as proposed a future user study. Future work will continue to
develop and test this framework further.

##  7\. Acknowledgments

This work was partially supported by the Wallenberg AI, Autonomous Systems and
Software Program (WASP) funded by the Knut and Alice Wallenberg Foundation and
the Vinnova Competence Center for Trustworthy Edge Computing Systems and
Applications.

## References

  * (1)
  * Barentine et al. (2021) Christian Barentine, Andrew McNay, Ryan Pfaffenbichler, Addyson Smith, Eric Rosen, and Elizabeth Phillips. 2021.  A VR Teleoperation Suite with Manipulation Assist. In _Companion of the 2021 ACM/IEEE International Conference on Human-Robot Interaction_ (Boulder, CO, USA) _(HRI ’21 Companion)_. Association for Computing Machinery, New York, NY, USA, 442–446.  https://doi.org/10.1145/3434074.3447210
  * Chandan et al. (2021) Kishan Chandan, Vidisha Kudalkar, Xiang Li, and Shiqi Zhang. 2021\.  ARROCH: Augmented reality for robots collaborating with a human. In _2021 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 3787–3793.
  * Foster et al. (2014) Mary Ellen Foster, Simon Keizer, and Oliver Lemon. 2014\.  Towards action selection under uncertainty for a socially aware robot bartender. In _Proceedings of the 2014 ACM/IEEE international conference on Human-robot interaction_. 158–159.
  * Guzzi et al. (2022) Jérôme Guzzi, Gabriele Abbate, Antonio Paolillo, and Alessandro Giusti. 2022.  Interacting with a Conveyor Belt in Virtual Reality Using Pointing Gestures. In _Proceedings of the 2022 ACM/IEEE International Conference on Human-Robot Interaction_ (Sapporo, Hokkaido, Japan) _(HRI ’22)_. IEEE Press, 1194–1195.
  * Hawes et al. (2017) Nick Hawes, Christopher Burbridge, Ferdian Jovan, Lars Kunze, Bruno Lacerda, Lenka Mudrova, Jay Young, Jeremy Wyatt, Denise Hebesberger, Tobias Kortner, et al. 2017\.  The strands project: Long-term autonomy in everyday environments.  _IEEE Robotics & Automation Magazine_ 24, 3 (2017), 146–156.
  * Honig and Oron-Gilad (2018) Shanee Honig and Tal Oron-Gilad. 2018.  Understanding and resolving failures in human-robot interaction: Literature review and model development.  _Frontiers in psychology_ 9 (2018), 861.
  * Kaarlela et al. (2020) Tero Kaarlela, Sakari Piesk, and Tomi Pitkaho. 2020\.  Digital Twin and Virtual Reality for Safety Training. In _2020 11th IEEE International Conference on Cognitive Infocommunications (CogInfoCom)_. 000115–000120.  https://doi.org/10.1109/CogInfoCom50765.2020.9237812
  * Kennel-Maushart et al. (2022) Florian Kennel-Maushart, Roi Poranne, and Stelian Coros. 2022\.  Multi-Arm Payload Manipulation via Mixed Reality. In _2022 International Conference on Robotics and Automation (ICRA)_. 11251–11257.  https://doi.org/10.1109/ICRA46639.2022.9811580
  * Khavas et al. (2020) Zahra Rezaei Khavas, S Reza Ahmadzadeh, and Paul Robinette. 2020\.  Modeling trust in human-robot interaction: A survey. In _International Conference on Social Robotics_. Springer, 529–541.
  * Li (2015) Jamy Li. 2015.  The benefit of being physically present: A survey of experimental works comparing copresent robots, telepresent robots and virtual agents.  _International Journal of Human-Computer Studies_ 77 (2015), 23–37.
  * Liu et al. (2021) Mengnan Liu, Shuiliang Fang, Huiyue Dong, and Cunzhi Xu. 2021\.  Review of digital twin about concepts, technologies, and industrial applications.  _Journal of Manufacturing Systems_ 58 (2021), 346–361.
  * Liu et al. (2017) Oliver Liu, Daniel Rakita, Bilge Mutlu, and Michael Gleicher. 2017.  Understanding human-robot interaction in virtual reality. In _2017 26th IEEE international symposium on robot and human interactive communication (RO-MAN)_. IEEE, 751–757.
  * Lucas et al. (2018) Gale M Lucas, Jill Boberg, David Traum, Ron Artstein, Jonathan Gratch, Alesia Gainer, Emmanuel Johnson, Anton Leuski, and Mikio Nakano. 2018.  Getting to know each other: The role of social dialogue in recovery from errors in social robots. In _Proceedings of the 2018 acm/ieee international conference on human-robot interaction_. 344–351.
  * Malik and Brem (2021) Ali Ahmad Malik and Alexander Brem. 2021.  Digital twins for collaborative robots: A case study in human-robot interaction.  _Robotics and Computer-Integrated Manufacturing_ 68 (2021), 102092\.
  * Malle and Ullman (2021) Bertram F Malle and Daniel Ullman. 2021.  A multidimensional conception and measure of human-robot trust.  In _Trust in human-robot interaction_. Elsevier, 3–25.
  * Mara et al. (2021) Martina Mara, Jan-Philipp Stein, Marc Erich Latoschik, Birgit Lugrin, Constanze Schreiner, Rafael Hostettler, and Markus Appel. 2021\.  User responses to a humanoid robot observed in real life, virtual reality, 3D and 2D.  _Frontiers in Psychology_ 12 (2021), 633178.
  * Naceri et al. (2019) Abdeldjallil Naceri, Dario Mazzanti, Joao Bimbo, Domenico Prattichizzo, Darwin G. Caldwell, Leonardo S. Mattos, and Nikhil Deshpande. 2019.  Towards a Virtual Reality Interface for Remote Robotic Teleoperation. In _2019 19th International Conference on Advanced Robotics (ICAR)_. 284–289.  https://doi.org/10.1109/ICAR46387.2019.8981649
  * Nenna and Gamberini (2022) Federica Nenna and Luciano Gamberini. 2022.  The Influence of Gaming Experience, Gender and Other Individual Factors on Robot Teleoperations in VR. In _Proceedings of the 2022 ACM/IEEE International Conference on Human-Robot Interaction_. 945–949.
  * Ortenzi et al. (2022) Valerio Ortenzi, Maija Filipovica, Diar Abdlkarim, Tommaso Pardi, Chie Takahashi, Alan Wing, Massimiliano Di Luca, and Katherine J Kuchenbecker. 2022.  Robot, Pass Me the Tool: Handle Visibility Facilitates Task-oriented Handovers. In _Proceedings of the ACM/IEEE International Conference on Human-Robot Interaction (HRI)_. 1–9.
  * Ostanin and Klimchik (2018) Mikhail Ostanin and Alexandr Klimchik. 2018.  Interactive robot programing using mixed reality.  _IFAC-PapersOnLine_ 51, 22 (2018), 50–55.
  * Robinson et al. (2014) Hayley Robinson, Bruce MacDonald, and Elizabeth Broadbent. 2014\.  The role of healthcare robots for older people at home: A review.  _International Journal of Social Robotics_ 6, 4 (2014), 575–591.
  * Rosen et al. (2018) Eric Rosen, David Whitney, Elizabeth Phillips, Daniel Ullman, and Stefanie Tellex. 2018\.  Testing robot teleoperation using a virtual reality interface with ROS reality. In _Proceedings of the 1st International Workshop on Virtual, Augmented, and Mixed Reality for HRI (VAM-HRI)_. 1–4.
  * Szafir (2019) Daniel Szafir. 2019\.  Mediating human-robot interactions with virtual, augmented, and mixed reality. In _International Conference on Human-Computer Interaction_. Springer, 124–149.
  * Togias et al. (2021) Theodoros Togias, Christos Gkournelos, Panagiotis Angelakis, George Michalos, and Sotiris Makris. 2021\.  Virtual reality environment for industrial robot control and path design.  _Procedia CIRP_ 100 (2021), 133–138.
  * Venkatesh and Davis (2000) Viswanath Venkatesh and Fred D Davis. 2000.  A theoretical extension of the technology acceptance model: Four longitudinal field studies.  _Management science_ 46, 2 (2000), 186–204.
  * Wang et al. (2021) Xi Wang, Ci Jyun Liang, Carol C Menassa, and Vineet R Kamat. 2021\.  Interactive and immersive process-level digital twin for collaborative human–robot construction work.  _Journal of Computing in Civil Engineering_ 35, 6 (2021), 04021023\.
  * Whitney et al. (2020) David Whitney, Eric Rosen, Elizabeth Phillips, George Konidaris, and Stefanie Tellex. 2020\.  Comparing robot grasping teleoperation across desktop and virtual reality with ROS reality.  In _Robotics Research_. Springer, 335–350.
  * Wonsick and Padır (2021) Murphy Wonsick and Taşkın Padır. 2021.  Human-humanoid robot interaction through virtual reality interfaces. In _2021 IEEE Aerospace Conference (50100)_. IEEE, 1–7.
  * Wozniak and Jensfelt (2022) Maciej K Wozniak and Patric Jensfelt. 2022.  Virtual Reality Framework for Better Human-Robot Collaboration and Mutual Understanding.  (2022).
  * Wu et al. (2019) Yuxin Wu, Alexander Kirillov, Francisco Massa, Wan-Yen Lo, and Ross Girshick. 2019.  Detectron2.  https://github.com/facebookresearch/detectron2.
  * Xu et al. (2022) Shiyu Xu, Scott Moore, and Akansel Cosgun. 2022.  Shared-Control Robotic Manipulation in Virtual Reality.  _arXiv preprint arXiv:2205.10564_ (2022).