∎

11institutetext: Alberto De Santis, Stefano Lucidi, Tommaso Giovannelli,
Massimo Roma 22institutetext: Dipartimento di Ingegneria Informatica
Automatica e Gestionale “A. Ruberti”, SAPIENZA Università di Roma, via
Ariosto, 25 – 00185 Roma, Italy.
22email<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>33institutetext: Mauro
Messedaglia 44institutetext: ACTOR Start up of SAPIENZA Università di Roma,
via Nizza 45, 00198 Roma, Italy. 44email<EMAIL_ADDRESS>
# A simulation-based optimization approach for the calibration of a discrete
event simulation model of an emergency department

Alberto De Santis https://orcid.org/0000-0001-5175-4951 Tommaso Giovannelli
https://orcid.org/0000-0002-1436-5348 Stefano Lucidi
https://orcid.org/0000-0003-4356-7958 Mauro Messedaglia  Massimo Roma
https://orcid.org/0000-0002-9858-3616

###### Abstract

Accurate modeling of the patient flow within an Emergency Department (ED) is
required by all studies dealing with the increasing and well-known problem of
overcrowding. Since Discrete Event Simulation (DES) models are often adopted
with the aim of assessing solutions for reducing the impact of this worldwide
phenomenon, an accurate estimation of the service time of the ED processes is
necessary to guarantee the reliability of the results. However, simulation
models concerning EDs are frequently affected by data quality problems, thus
requiring a proper estimation of the missing parameters.

In this paper, a simulation-based optimization approach is used to estimate
the incomplete data in the patient flow within an ED by adopting a model
calibration procedure. The objective function of the resulting minimization
problem represents the deviation between simulation output and real data,
while the constraints ensure that the response of the simulation is
sufficiently accurate according to the precision required. Data from a real
case study related to a big ED in Italy is used to test the effectiveness of
the proposed approach. The experimental results show that the model
calibration allows recovering the missing parameters, thus leading to an
accurate DES model.

###### Keywords:

Simulation-based Optimization Emergency Department Discrete Event Simulation
Model Calibration Derivative-Free Optimization

††journal:

##  1 Introduction

Over the last few years, Emergency Departments (EDs) have been raising an
increasing attention in the Operations Research and Management Science
communities due to the international phenomenon of the overcrowding
hoot.aronsky:08 ; Paul.2010 ; Wiler:2011 ; disomma.2014 ; Nahhas.2017 , which
leads to longer waiting times, higher mortality rates, and lower patient
satisfaction pines.2008 . One of the most popular tools adopted to study this
crucial problem is Discrete Event Simulation (DES) Zeinali.2015 ; Joshi.2016 ;
Wong.2016 ; Bedoya-Valencia:2016 ; kramer.2020 ; fava.2021 , which is used to
represent the complex and stochastic patient flows in the ED in place of
analytical models (for a recent literature review of simulation modeling in
EDs, see Salmon:2018 ). Frequently, DES models are combined with an
optimization algorithm to make optimal decisions in relation to some Key
Performance Indicators (KPIs). The resulting Simulation-Based Optimization
(SBO) approach fu.2015 ; gosavi.2014 ; amaran.2014 , whose effective
application relies on the accuracy of the simulation model, has been
frequently applied to hospital EDs (for a recent review, see, e.g.,
yousefi.2020 ).

Many papers in the literature deal with the quality of the input data used in
ED simulation models. This issue, which strongly affects the reliability of
the results, has been carefully analyzed. Indeed, the cost, time, and
challenges required by collecting ED empirical data represents a serious
limitation for every simulation model Paul.2010 . In order to accurately
replicate and predict the patient flows within the ED, duma.2018  develops a
process mining approach which handles the noise factors in the dataset after
introducing assumptions on how to interpret the data. In particular, these
factors include the following cases: starting and ending time of each activity
performed in the ED may be unknown; information about urgent patients may be
registered after the activity is completed and, in general, the timestamps of
each activity may not be promptly recorded at the right time; for each
patient, the same activity may be recorded multiple times for technical
reasons, giving rise to misleading information. A framework to categorize all
the ED data quality issues is proposed in vanbrabant:2019 , which also
provides assessment techniques for each data quality problem category.
Moreover, this paper highlights that most of the works in the literature focus
on the problem of missing data, which is the problem addressed in this paper.
It consists in the lack of information on some or all the key timestamps that
define the activities performed in the ED, i.e., starting and ending time.
This well-known issue prevents gaining knowledge about the duration of each
activity, which is required to estimate the corresponding probability
distributions to use in the simulation model. Several SBO approaches are
proposed in the literature to tackle this crucial problem. The idea behind
each approach is to leverage the known information in order to estimate the
parameters of the probability distributions underlying the missing data. This
is accomplished by comparing the Key Performance Indicators (KPIs) computed
through the simulation model with the corresponding values derived by the data
collected in the real system. The resulting procedure is known as _model
calibration_.

The mathematical formulation of the optimization problem and the specific
algorithm used for determining the optimal solution are the two features that
distinguish the papers dealing with missing data. For instance, both Kuo.2016
and Guo.2016  adopt similar approaches that leverage the time differences
between the known timestamps. The former proposes an unconstrained
optimization problem where the objective function is a consistency measure
that compares the average, the standard deviation, and the proportions of the
time differences for each triage tag. The latter uses a constrained
optimization problem where the objective function is based on a modified
chi–square goodness of fit and the constraints are introduced to guarantee
each time difference the same level of accuracy. As regards the approaches
used to solve the problem, both papers use metaheuristic procedures. In
particular, Kuo.2016  considers both a descent method and a simulated
annealing algorithm, while Guo.2016  adopts an approach that combines a
genetic algorithm with simulated annealing and optimal computing budget
allocation. Moreover, in Guo.2016  the authors point out that the approach in
Kuo.2016  could be improved in several ways: including in the objective
function all the time differences defined by the available timestamps, since
delays in one activity may impact on downstream activities in the patient
flow; modeling the possibility for each patient to have more than one medical
visit, which may affect the time differences; using a more formal objective
function and solving the resulting optimization problem through a more
efficient algorithm.

Another paper in the literature on missing data that is worth being mentioned
is Liu.2017 , which assumes an agent–based simulation model, as opposed to the
DES model considered in Kuo.2016  and Guo.2016 . Although the framework
underlying this work is different, a simulation-based optimization problem is
used for the same goal of minimizing the deviation between real data and
simulation output in order to estimate the missing parameters of the
simulation model. The Length Of Stay (LOS), i.e., the difference between the
discharge time and the arrival time to the ED, is the time difference
considered in the objective function, which adopts the Jensen–Shannon
divergence. A systematic method based on the pattern search method APPSPACK
gray.2006  is used to find the optimal configuration of parameters.

The approach proposed in this paper aims to handle the problem of missing
timestamps, which affects many simulation models dealing with EDs, as
evidenced by the papers in the specific literature. Among the sources of noise
affecting the quality of the input data necessary to build a reliable
simulation model, the problem of missing timestamps has a strong impact on the
overall accuracy of the simulation, since it prevents the knowledge of the
values used to derive appropriate probability distributions. Therefore, while
other issues may be resolved by either carefully cleaning the dataset or
introducing assumptions on how to interpret the data, the unavailability of
timestamps requires more sophisticated procedures.

Compared to the other papers dealing with missing data in ED datasets, the
approach discussed in this paper aims both to propose a new formulation of the
resulting optimization problem and to improve on the optimization strategies
typically used in the literature. Since building a simulation model is a
process that requires a considerable amount of effort and thus is not expected
to be completed in a short time, an exact algorithm providing optimal
solutions may be preferable to metaheuristic procedures, whose final solutions
are returned faster but without optimality guarantees. However, although
global convergent algorithms appear to be a reasonable choice, metaheuristics
are the methods mainly adopted in the literature dedicated to missing data
(see, e.g., Kuo.2016 ; Guo.2016 ). To fill this gap, in this paper a SBO
approach is developed both to propose an alternative version of the
optimization problems used in Kuo.2016 ; Guo.2016  and to use a Derivative-
Free Optimization (DFO) strategy based on global convergence, which allows the
algorithm to find an optimal solution with optimality guarantees. Similarly to
Kuo.2016  and differently from Guo.2016 , Weibull distributions are adopted to
generate the values of the activity service times associated with missing
data. Indeed, this probability distribution is considered suitable when data
is unavailable (see, e.g., law:15 ). Moreover, the most critical patients, who
are excluded from the two approaches developed by Kuo.2016  and Guo.2016 , are
included to prevent the simulation model from returning inaccurate results for
this important class of patients. Finally, in the proposed constrained
optimization problem, the objective function is based on the comparison
between real and simulated Empirical Cumulative Distribution Functions
(ECDFs), which do not rely on intervals, unlike the functions adopted in
Kuo.2016 ; Guo.2016 ; Liu.2017 .

This paper is organized as follows. Section 2 reports the typical structure of
the datasets related to the patient flow within an ED, highlighting the
available information and the KPIs of interest. Section 3 defines the SBO
problem used to calibrate a general DES model concerning an ED. Section 4
describes the real ED used as a case study and the DES model implemented.
Section 5 reports the results obtained by applying the SBO approach to the
case study. Finally, in Section 6 we draw some concluding remarks and we
discuss directions for future work.

##  2 Data collection in ED

The timestamps defining the activities in the patient flow are represented in
Figure 2.1 and described in Table 2.1.

Figure 2.1: Timestamps collected throughout the patient flow. Table 2.1: Description of the timestamps collected throughout the patient flow. _$t_{0}$_ | starting time of triage.
---|---
_$t_{1}$_ | ending time of triage.
_$t_{2}$_ | starting time of medical visit.
_$t_{3}$_ | ending time of medical visit and starting time of examinations.
_$t_{4}$_ | ending time of the last exam or reassessment.
_$t_{5}$_ | patient receives the last medical report.
_$t_{6}$_ | patient is discharged and leaves the ED.

Although these timestamps are commonly recorded by the electronic systems used
to collect data, the actual patient flow may include extra times omitted from
this scheme. For example, since the arrival time to the ED is usually not
registered, before triage there may be an additional waiting time not included
in any record. Moreover, extra waiting times and observation periods may be
present both in the examination part and in the discharge phase. Since usually
there is no available data for these times, in Figure 2.1 both waiting and
service times for exams are included in a single box between $t_{3}$ and
$t_{4}$, which are associated with the end of the medical visit and the end of
the last exam or reassessment after a visit, respectively. After receiving the
medical report, patients may be either subject to additional medical visits or
discharged after a final waiting time, which is due to a further observation
or to the time required by the ED staff to prepare the discharge. It is
important to point out that the disposition, which refers to the decision to
admit a patient to a hospital ward, is associated with $t_{6}$, like the
discharge. This implies that the boarding time of the admitted patients (i.e.,
the waiting time before being transferred to the hospital ward) is assumed to
start at $t_{6}$, thus being excluded from the ED scope. This choice is in
accordance with the case study considered in Section 4. Contrarily, in some
EDs the discharge of patients admitted may be considered as the time of
admission to the hospital ward (see, e.g., Baumlin.2010 ).

In most cases, the available timestamps are $t_{0}$, $t_{2}$, and $t_{6}$,
which are marked in the scheme as known. Since this is a typical setting in
practice, throughout this paper they are supposed to be the only timestamps
recorded by the ED along with $t_{5}$ (which may not be always available). As
a consequence, for each patient it is possible to compute

  * •

the Door-to-Doctor Time (DOT), which is the time difference between the start
of the triage and the start of the medical visit, namely $t_{2}-t_{0}$;

  * •

the Doctor-to-Discharge Time (DIT), which is the time difference between the
start of the medical visit and the discharge, namely $t_{5}-t_{2}$ if $t_{5}$
is available, $t_{6}-t_{2}$ otherwise.

It is important to remark that when patients require additional medical visits
after the exams, $DIT$ can be interpreted as the time difference between the
timestamp of the last medical report (or the discharge if $t_{5}$ is not
available) and the starting time of the first medical visit. Note that $t_{4}$
and $t_{5}$ are equal if a patient does not need additional medical visits.

Since in many cases only a subset of the timestamps related to the patient
flow are known vanbrabant:2019 , the service time cannot be computed for all
the activities. Apart from the final waiting time, which can be recovered for
each patient through the difference $t_{6}-t_{5}$, note that the durations of
all the other activities are not completely defined. In fact, in case of
triage and medical visit, the ending time is missing, while for exams both
timestamps are unknown.

##  3 Simulation-based optimization problem for ED model calibration

The goal of the approach proposed is to recover the information needed to
build an accurate simulation model by leveraging the known information through
the minimization of the deviation between real data and simulation output. In
order to define the simulation-based optimization problem, let us introduce
the sets below.

  * •

Let $C$ be the set of the triage tags.

  * •

Let $U(c)$ be the set of units where patients with tag $c\in C$ can be visited
and treated.

  * •

Let $\mathcal{T}=\\{DOT,DIT\\}$ be the set of the time differences considered.

In the absence of data, suitable probability distributions are Weibull and
lognormal (see law:15 ). Let us arbitrarily consider a Weibull distribution,
whose probability density function is

$f(w)=\begin{cases}\alpha\beta^{-\alpha}w^{\alpha-1}e^{-(w/\beta)^{\alpha}}\text{
if }w>0,\\\ 0\text{ otherwise,}\end{cases}$

where $\alpha>0$ and $\beta>0$ are the shape and the scale parameters,
respectively. Therefore, a Weibull distribution is assumed for the triage,
medical visit, and exams service times and, for each of them, different pairs
of shape and scale parameters are considered based on the triage tag $c$ and
on the unit $u$. This choice leads to a number of different pairs of
parameters equal to $\sum_{c\in C}|U(c)|$, where $|U(c)|$ refers to the
cardinality of the set $U(c)$. Let us denote as $x\in\mathbb{R}^{n_{1}}$,
$y\in\mathbb{R}^{n_{2}}$, and $z\in\mathbb{R}^{n_{3}}$ the corresponding
probability distribution parameters for the service times considered, where
$n_{v}\leq 2\,\sum_{c\in C}|U(c)|$ with $v\in\\{1,2,3\\}$. In particular, for
all $c\in C$ and for all $u\in U(c)$, the shape and scale parameters of the
Weibull distributions are

  * •

$x_{1}^{c\,u}$ and $x_{2}^{c\,u}$ for the triage probability distribution,

  * •

$y_{1}^{c\,u}$ and $y_{2}^{c\,u}$ for the medical visit probability
distribution,

  * •

$z_{1}^{c\,u}$ and $z_{2}^{c\,u}$ for the exams probability distribution.

Let $F^{sim}_{c\,u\,i}$ and $F^{real}_{c\,u\,i}$ be the ECDFs of the values of
the simulated and real time difference $i\in\mathcal{T}$ for patients with tag
$c$ visited in unit $u$. Moreover, let $k^{sim}_{c\,u\,i}$ and
$k^{real}_{c\,u\,i}$ be the number of such patients from the simulation and
from the real dataset. Hence, for all $j\in\\{sim,real\\}$, we have that

$F^{j}_{c\,u\,i}(t)=\frac{1}{k^{j}_{c\,u\,i}}\sum_{h=1}^{k^{j}_{c\,u\,i}}\textbf{1}_{\\{\tau_{{c\,u\,i}\,h}\leq
t\\}}\quad\text{ with }t\geq 0,$

where $\tau_{{c\,u\,i}\,h}$ is the value of the time difference $i$ recorded
for the $h$–patient, with $h\in\\{1,\ldots,k^{j}_{c\,u\,i}\\}$, considered
from $j\in\\{sim,real\\}$. It is important to point out that the values
$\tau_{{c\,u\,i}\,h}$ of the time differences computed from the simulation
depend on the service times of triage, medical visit, and exams drawn from the
Weibull distributions described above. Therefore, this dependence can be
written explicitly as $F^{sim}_{c\,u\,i}(t;x,y,z)$, where $x$, $y$, and $z$
are the vectors containing all the associated shape and scale parameters.

The mathematical problem formulation is reported as follows

$\begin{split}\min_{x,\,y,\,z}\text{ }&\sum_{c\in C}\sum_{u\in
U(c)}\sum_{i\in\mathcal{T}}\Big{(}\int_{0}^{\infty}(F^{sim}_{c\,u\,i}(t;x,y,z)-F^{real}_{c\,u\,i}(t))^{2}\,dt\Big{)}\\\
s.t.\text{ }\,&x\in{\cal P},\end{split}$ (3.1)

where the feasible set

$\begin{array}[]{l}{\cal
P}=\Bigl{\\{}(x,y,z)\in{\mathbb{R}}^{n_{1}}\times{\mathbb{R}}^{n_{2}}\times{\mathbb{R}}^{n_{3}}\leavevmode\nobreak\
|\vspace{0.2truecm}\nobreak\leavevmode\\\ \vspace{0.2truecm}\qquad\qquad\qquad
g_{c\,u\,i}(x,y,z)\leq 0,\\\ \vspace{0.2truecm}\qquad\qquad\qquad
h_{c\,u\,i}(x,y,z)\leq 0,\\\ \vspace{0.2truecm}\qquad\qquad\qquad l_{x}\leq
x\leq u_{x}\\\ \vspace{0.2truecm}\qquad\qquad\qquad l_{y}\leq y\leq u_{y}\\\
\vspace{0.2truecm}\qquad\qquad\qquad l_{z}\leq z\leq
u_{z}\bigr{.}\cr\cr\bigl{.}\hskip 19.91684pt\text{ for all }c\in C,\,u\in
U(c),\,i\in\mathcal{T}\Bigr{\\}}\end{array}$

is defined by the following functions

  * •

$g_{c\,u\,i}(x,y,z)=\big{|}\frac{\mu^{sim}_{c\,u\,i}(x,y,z)-\mu^{real}_{c\,u\,i}}{\mu^{real}_{c\,u\,i}}\big{|}-tol_{\mu}^{c\,u\,i}$,
which compares the sample means $\mu^{sim}_{c\,u\,i}(x,y,z)$ and
$\mu^{real}_{c\,u\,i}$ of the time difference $i$ computed through the
simulated (by averaging over the independent replications) and real data,
respectively,

  * •

$h_{c\,u\,i}(x,y,z)=\big{|}\frac{\sigma^{sim}_{c\,u\,i}(x,y,z)-\sigma^{real}_{c\,u\,i}}{\sigma^{real}_{c\,u\,i}}\big{|}-tol_{\sigma}^{c\,u\,i}$,
which compares the sample standard deviations $\sigma^{sim}_{c\,u\,i}(x,y,z)$
and $\sigma^{real}_{c\,u\,i}$ of the time difference $i$ computed through the
simulated (by averaging over the independent replications) and real data,
respectively,

and $l_{x},l_{y},l_{z},u_{x},u_{y}$, and $u_{z}$ are vectors defining bound
constraints. Positive thresholds, namely $tol_{\mu}^{c\,u\,i}$ and
$tol_{\sigma}^{c\,u\,i}$, are used to state the degree of accuracy required
for the simulation model. The objective function is the sum of the integrals
of the squared difference between $F^{sim}_{c\,u\,i}$ and $F^{real}_{c\,u\,i}$
over the sets $C$, $U(c)$, and $\mathcal{T}$. The choice of using ECDFs
instead of histograms, which are commonly adopted in the literature Kuo.2016 ;
Guo.2016 , is due to the different reliability of the information provided. In
particular, compared to ECDFs, the description of the data provided by
histograms is strongly affected by the choice of the width of each bin, which
is not straightforward and may lead to distributions with different shapes
based on how data is grouped. The decision variables are the Weibull
distribution parameters contained in the vectors $x,y$, and $z$. For the sake
of simplicity, the previous formulation does not include variables that are
not parameters of probability distributions, although they may be present. It
is important to remark that in the general framework described, the dependence
of each pair of parameters on both $c$ and $u$ implies that the triage service
time is affected by the triage tag and the ED unit. This assumption, which
turns out to be reasonable when applied to the duration of medical visit and
exams, may lead to an excessive number of variables for the triage service
time, if it does not hold. Indeed, while the color tag significantly affects
the triage duration since urgent patients undergo a faster triage than less
critical patients, the impact of the unit may be negligible. However, as
concerns the case study described in Section 4, interviews with the ED staff
have shown that slightly different procedures are adopted by the nurses in
charge of triage based on the unit where a patient is assigned. This means
that final conclusions should be drawn after analyzing the specific system
considered. Since the probability distribution parameters affect the dimension
of the optimization problem, avoiding unnecessary variables allows the
algorithm to benefit from a lower computational cost.

##  4 Case study: ED of Policlinico Umberto I in Rome

The case study concerns the ED of Policlinico Umberto I, which is a very large
hospital in Rome, Italy. It is the biggest ED in the region of Lazio in terms
of number of patient arrivals per year (about 140,000 on average). By using
the data collected from the patient flow through the ED for the whole year
2018, this case study is adopted to test the effectiveness of the approach
proposed in this paper.

###  4.1 Description of the ED

The ED of Policlinico Umberto I is divided into several areas, each one
associated with a medical specialty. The backbone of the ED is the _central
area_ , which is devoted to treating diseases and disorders related to
internal medicine and general surgery, which affect the majority of patients.
Separated from this main area, there are other parts of the ED that deal with
the following medical specialties: ophthalmology, obstetrics, pediatrics,
hematology, and dentistry. The focus of this paper is on the central area,
whose number of arrivals in 2018 amounted to more than 50,000.

In order to gain a complete understanding of the ED processes, the description
of the ED units and staff and the summary of the patient flow are reported in
the sequel. Other than a _triage area_ , where each incoming patient is
assigned a color tag by a nurse in charge of this task, the ED is composed by

  * $\bullet$

a _Medical Unit (MU)_ , devoted to patients needing specialized medications
and treatments, with areas dedicated to the most critical patients;

  * $\bullet$

a _Surgical Unit (SU)_ , devoted to patients needing either to receive a
surgical operation or to recover from it, with areas dedicated to the most
critical patients;

  * $\bullet$

a _Resuscitation Area (RA)_ , for the most acutely ill and injured patients,
who need timely treatments;

  * $\bullet$

a _Minor Injuries Unit (MIU)_ , for the least urgent patients, whose treatment
can be delayed or deferred;

  * $\bullet$

an _Orthopedic Unit (OU)_ , for patients suffering from orthopedic disorders.

Moreover, all of these units have rooms where patients can either wait for
exams or stay for observation. Red-tagged patients can be visited and treated
in RA or in dedicated areas within MU and SU, which are open 24 hours a day
and are provided with equipment and staff specialized for dealing with life-
threatening illnesses and injuries. In particular, 1 and 2 seats are available
in the dedicated areas of MU and SU, respectively, and further 2 seats are
available in RA. As concerns the medical treatment of the other patients, MU
and SU can host up to three and two patients during the day (8.00 a.m.–8.00
p.m.), respectively. At night, MU can host two patients, while in SU one seat
is available. Moreover, MIU has two seats, which are available from 8.00 a.m.
to 8.00 p.m., Monday through Saturday. When patients experience excessive
waiting times, two additional seats may be used to visit up to four patients
simultaneously. All this information is summarized in Tables 4.1–4.2.

Table 4.1: Number of seats available for medical visit and treatment in MU, SU, and MIU, the latter being open from Monday to Saturday. | MU | SU | MIU
---|---|---|---
Day (8.00 a.m.–8.00 p.m.) | 3 | 2 | 2
Night (8.00 p.m.–8.00 a.m.) | 2 | 1 | 0
Table 4.2: Number of seats available for medical visit and treatment of red-tagged patients in RA and in the dedicated areas of MU and SU. | MU | SU | RA
---|---|---|---
Day (8.00 a.m.–8.00 p.m.) | 1 | 2 | 2
Night (8.00 p.m.–8.00 a.m.) | 1 | 2 | 2
Table 4.3: Feasible assignments of patients to the ED units according to the color tag. A cross at the entry $(i,j)$ indicates that a patient with color tag $i$ can be assigned to the unit $j$.  | MU | SU | RA | MIU | OU
---|---|---|---|---|---
White | - | - | - | X | X
Green | X | X | - | X | X
Yellow | X | X | - | - | X
Red | X | X | X | - | -

As regards the patient flow, after arriving autonomously or by emergency
medical vehicles, all the incoming patients are admitted to the triage area,
where a nurse assigns the color tag. After the triage, the patients are
visited and treated in one of the units previously described, according to the
color tag assigned and the severity of the illness/injury. Table 4.3
represents a scheme showing the units where patients may be assigned based on
the color tags. In case of red-tagged patients, the medical visit is timely
performed in RA or in the dedicated areas of MU and SU. As concerns the other
color tags, the yellow and green tagged patients share the same resources in
MU and SU. However, while the former patients can be assigned only to MU and
SU, the latter may be sent to MIU during its opening hours if their health
conditions are deemed as not likely to worsen. If MIU is closed, all the
green-tagged patients are visited and treated in MU and SU. This diversion
allows the ED to reduce the occurrence of work overload in SU and MU, which
may give rise to overcrowded units. Moreover, it is important to point out
that the white triage tag is assigned only if MIU is open, otherwise the green
tag is used.

In many cases, after the medical visit, additional exams may be required.
Other than performing reassessments of the patients and requiring additional
exams, physicians may also require further observation periods. Finally, at
the end of the pathway, patients are discharged from the ED. This last stage
includes a final waiting time whose duration depends on the type of outcome.
Indeed, a longer wait is expected for patients that need to either be
hospitalized at a hospital ward or transferred to another hospital, while
patients discharged home can usually leave the ED in shorter time.

###  4.2 The Discrete Event Simulation model

This section describes the DES model of the ED of Policlinico Umberto I, which
has been implemented by using Simmer simmer.2019 , a process-oriented and
trajectory-based DES package for R. In the DES model, patients are represented
by the model entities, which are created according to the statistical model
used for modeling the arrival process. After arriving at the ED, the simulated
patient flow followed by each entity is based on trajectories determined by
the logical rules used to build the model. Each trajectory is associated with
the pathway followed by a patient according both to the color tag received at
triage and to the ED unit assigned. Figures 4.1–4.5 show the patient flow from
the arrival to the discharge according to the color tag assigned at triage. In
particular, Figure 4.1 focuses on the logic underlying the assignment of color
tag at triage. The other figures deal with the segments of the patient flow
following the triage, providing in brackets the information about the
resources seized in case the corresponding process is not a simple delay. Such
resources represent the seats at each ED unit, whose capacity can be either
fixed, as is the case with RA and the areas of MU and SU dedicated to red
patients, or based on schedule, as for MU, SU, and MIU.

The patient flow through the model can be described as follows. After being
created at the beginning of the simulation, the entities corresponding to
deceased patients leave the model, while the remaining entities enter the
triage area for the assignment of both color tag and ED unit, which are stored
as entity attributes. The discrete uniform probability distributions used for
assigning the color tag and the ED unit vary according to the time at which
the entity starts the triage. If this event happens in the daytime, the
corresponding probability distributions include also the white color tag and
MIU among the alternatives to sample. The triage phase is represented by as
many delay processes as the number of color tags and units. Each of these
delays is associated with a probability distribution that returns the value of
the triage duration for each patient. In particular, the 8 different processes
considered correspond to the possible pairs of color tags and units reported
in Table 4.3, with the exception of the pair given by the green color tag and
MIU. Indeed, the triage duration for a green-tagged patient eligible for MIU
depends on the specialty required, whether medical or surgical. Since this
information is unknown, in the DES model half of the MIU patients are subject
to the triage duration of MU patients, while the other half are subject to the
triage duration of SU patients. It is important to remark that, since the
processes are simple delay, no queue is generated before the triage. The
reason underlying this choice is twofold: on the one hand, this allows for
alignment between data e simulation model, since the timestamps denoted as
$t_{0}$ in Figure 2.1 are used to derive the arrival probability distribution,
which consequently returns the starting time of triage (and not the starting
time of the waiting time for triage); on the other hand, the waiting time
before triage can be considered negligible, according to the ED staff.

After the triage process, some entities leave the model, thus representing the
patients who leave without being seen, while the other ones start waiting for
the medical visit. The entity selected for the visit depends both on the
priority class, i.e., the color tag, and on the FIFO criterion. Then, when the
visit begins, one seat in the corresponding unit is seized and the duration of
the visit is returned by a suitable probability distribution based on the
entity color tag.

After the medical visit, the following phase includes exams and reassessments,
whose duration is generated by means of a suitable probability distribution.
Like the triage, a single delay process is used for this phase as well, due to
the lack of knowledge of the resources required. It is important to point out
that all possible further treatments are considered included in the service
time of this process. After the exams, the entities corresponding to patients
refusing hospitalization, leaving during exams, and being transferred to
another hospital are removed from the model, while the other entities proceed
to the next stage, i.e., the final waiting time before discharge.

Since in the dataset some timestamps defining the starting and ending time of
triage, medical visit, and exams are not available, the service times
corresponding to these activities cannot be recovered directly. Moreover,
additional delay processes are considered to reproduce all the service times
that cannot be directly computed through the data, such as the setup times
that are sometimes needed for sanitizing the ED areas and the idle times
caused by unexpected requests for personnel from other ED units, which give
rise to sudden activity disruption. In order to effectively model these times,
which impact on all the patients without predictable patterns, uniform
probability distributions are considered before the queue for the medical
visit, according to the color tag.

Figure 4.1: Plot of the patient flow in the ED from the arrival to the
assignment of the color tag. Figure 4.2: Plot of the red-tagged patient flow.

Figure 4.3: Plot of the yellow tagged patient flow.

Figure 4.4: Plot of the white tagged patient flow.

Figure 4.5: Plot of the green-tagged patient flow.

As regards the KPIs of interest, the focus is on the time differences $DOT$
and $DIT$ described in Figure 2.1, which can be computed through the data of
the real system and then compared with the corresponding values returned by
the simulation. Therefore, the KPIs of interest are

  * •

$DOT$, which is the time difference between the starting time of the triage
and the starting time of the visit, namely $t_{2}-t_{0}$ in Figure 2.1;

  * •

$DIT$, which is the time difference between the starting time of the visit and
the time of the discharge, namely $t_{6}-t_{2}$ in Figure 2.1.

###  4.3 Input analysis

To analyze the arrival process, we focus on the data collected from the 1st of
January to the 31st of March. In Figure 4.6, the weekly average hourly arrival
rate obtained by averaging the number of arrivals occurring in the same hourly
time slot over the 13 weeks considered is reported.

Figure 4.6: Plot of the weekly average hourly arrival rate for the first 13
weeks of the year.

It is worth observing that, in accordance with the literature (see, i.e.,
Kim.2014 ), the average arrival rate among the days of the week is
significantly different. Therefore, since averaging over these days would lead
to inaccurate results, the different days of the week must be considered
separately.

The stochastic process defined by the arrivals to the ED is assumed to be well
modeled by a Nonhomogeneous Poisson Process (NHPP), which is a standard
assumption in the literature (see, e.g.,Whitt.2017 ; Kuo.2016 ; Zeinali.2015 ;
Ahmed.2009 ; Ahalt.2018 ; Guo.2017 ). The validity of this hypothesis for the
case study considered is supported by the extensive analysis performed in
DGLMR:tr ; desantis.2021 , where statistical hypothesis testing is adopted. To
achieve an accurate representation of the arrival rate, 24 time slots of
unitary length are considered for each day of the week, thus obtaining a
piecewise constant approximation. While this approach allows taking into
account the within-day variation, the day-to-day variation is considered by
estimating the hourly arrival rate function separately for each day of the
week.

Although the collected data concerns the whole year 2018, the DES model of the
ED has been built by using the data related to January. The choice of focusing
on this month stems from the will of the ED management to reduce the
overcrowding level observed in the winter season, which exhibits longer
waiting times, as emerged by interviewing the ED staff. Among the winter
months, January is observed to suffer from the heaviest workload, which puts a
strain on the ED processes, thus requiring a careful analysis. However, the
simulation model could be also easily adapted to include input parameters
estimated from data related to different months.

From 00:00 of January 1 to 23:59 of January 31, 2018, the total number of
patients arrived to the ED is 4192. The timestamps recorded are the ones
marked as known in Figure 2.1. In Table 4.4, the number and the percentage of
color tags assigned at triage are reported along with the number of patients
who leave without being seen (LWBS). Since from Monday to Saturday MIU is open
from 8.00 a.m. to 8.00 p.m., the daytime and the night are analyzed separately
to reflect the change in the ED setting. Although, in some cases, re-
evaluations can lead to different patient tags at discharge with respect to
those assigned at triage, in the dataset considered this information is
unavailable.

Table 4.4: Number and percentage of color tags assigned at triage in the daytime (columns 2-3) and at night (column 4-5). | TAGS ASSIGNED AT TRIAGE |
---|---|---
| Day (8.00 a.m.–8.00 p.m.) | Night (8.00 p.m.–8.00 a.m.) | LWBS
White | 65 | 2.17% | - | - | 1
Green | 1306 | 44.32% | 498 | 40.42 % | 16
Yellow | 1420 | 48.18% | 638 | 51.79 % | -
Red | 157 | 5.33% | 95 | 7.71 % | -
| 2948 | 100% | 1231 | 100% | 17

Finally, Tables 4.5–4.6 show the number and proportion of the color tags among
the ED units. Although OU is out of the scope of this analysis, Table 4.5
includes this unit since OU patients share the triage station with the other
patients, thus affecting the counts reported in Table 4.4.

Table 4.5: Number (and percentage) of patients assigned to the ED units for each color tag. | MU | SU | RA | MIU | OU
---|---|---|---|---|---
White | - | - | - | 47 (73.44 %) | 17 (26.56 %)
Green | 248 (13.75 %) | 628 (34.81 %) | - | 260 (14.41 %) | 668 (37.03 %)
Yellow | 1316 (63.95 %) | 693 (33.67 %) | - | - | 49 (2.38 %)
Red | 191 (75.79 %) | 45 (17.86 %) | 16 (6.35 %) | - | -
| 1755 | 1366 | 16 | 307 | 734
Table 4.6: Number (and percentage) of green-tagged patients assigned to MU, SU, and MIU in the daytime (8.00 a.m. – 8.00 p.m.) and at night (8.00 p.m. – 8.00 a.m.). | MU | SU | MIU
---|---|---|---
Daytime | 132 (16.60 %) | 403 (50.69 %) | 260 (32.70 %)
Night | 116 (34.02 %) | 225 (65.98 %) | -
| 248 | 628 | 260

The percentage of patients requiring more than one medical visit is reported
in Table 4.7 for each color tag and ED unit. Since the number of patients
needing more than visit is small, in the DES model all the entities are
assumed to undergo one medical visit and, accordingly, each patient flow is
represented sequentially. The probability distributions used to generate
values for the final waiting time before the discharge are reported in Table
4.7 as well (note that the base unit for time adopted in the DES model is the
hour).

Table 4.7: Percentage of patients requiring more than one visit and probability distribution of the final waiting time before discharge for each color tag and ED unit. |  | MORE THAN ONE VISIT | FINAL WAITING TIME
---|---|---|---
White | MIU | 0% | Weib(1.12, 0.41)
| MIU | 1.61% | Weib(0.83, 0.57)
Green | SU | 2.44% | Weib(0.61, 1.07)
| MU | 4.8% | Weib(0.57, 4.91)
Yellow | SU | 2.04% | Weib(0.63, 2.33)
MU | 3.86% | Weib(0.62, 7.22)
| RA | 6.25% | Weib(0.67, 9.26)
Red | SU | 6.52% | Weib(0.71, 4.47)
| MU | 7.9% | Weib(0.69, 6.65)

###  4.4 Model verification/validation and design of experiments

The aim of the DES model of the ED of Policlinico Umberto I is to provide a
reliable tool for assessing the impact of changes to the current settings on
the overcrowding level. To this end, verification and validation are important
steps to be carefully considered. While standard techniques have been adopted
to verify the model, such as debugging and model trace, validation has
required a more sophisticated analysis. Indeed, a model calibration approach
has been used to achieve an accurate simulation output that well reproduces
the real system data.

Important KPIs to describe the ED status are the two time differences
contained in ${\cal T}$ and the number of entities generated by the simulation
and associated with each color tag. While the former KPIs are part of the
objective function used in the calibration procedure and they are calculated
across the iterations of the optimization algorithm, the latter do not require
a continuous monitoring, since their values are not affected by the
calibration. Indeed, the only variations are due to the different number of
entities discharged before the simulation ends, which leads to different
results in terms of KPI computation. Table 4.8 reports the values of the
patient counts obtained from the simulation of the starting point of the
calibration procedure along with their confidence intervals. By comparing
these values with the counts in Table 4.5, it can be observed that the
simulation model provides an accurate output in terms of number of patients.

Table 4.8: Output values (with confidence interval) for the number of patients returned by the simulation of the starting point of the calibration procedure. | MU | SU | MIU
---|---|---|---
White | - | - | $41\pm 2.82$
Green | $259\pm 5.33$ | $612\pm 9.73$ | $239\pm 6.98$
Yellow | $1335\pm 20.55$ | $716\pm 9.39$ | -
Red | $213\pm 9.39$ | $57\pm 5.35$ | -

Since the objective and constraint functions compare the real system data and
the simulation output, achieving adequate accuracy is imperative for a fair
comparison. To this end, the simulation output is estimated through 30
independent simulation replications, each of them 38 days long, with a warm up
period of 7 days, thus matching the 31 days of January.

##  5 Experimental results

The sets defined in Section 3 can be adapted to the case study as follows.

  * •

Let $C=\\{\text{W,G,Y,R}\\}$ be the set of the triage tags.

  * •

Let $U(c)\subseteq\\{\text{MU},\text{SU},\text{MIU}\\}$ be the set of the ED
units where patients with tag $c\in C$ can be visited and treated. In
particular,

    * –

$U(W)=\\{\text{MIU}\\}$,

    * –

$U(G)\ =\\{\text{MU},\text{SU},\text{MIU}\\}$,

    * –

$U(Y)\ =\\{\text{MU},\text{SU}\\}$,

    * –

$U(R)\ =\\{\text{MU},\text{SU}\\}$.

For the specific instance represented by this case study, 8 different pairs of
parameters are considered for the Weibull distributions representing the
service times of medical visit and exams. Each pair is associated with an
element of the sets $U(c)$. Instead, 7 pairs are used for the probability
distributions of triage for the reasons stated in Section 4.2, thus resulting
in 23 overall pairs of parameters, i.e., 46 unknown parameters to be
determined through the calibration procedure.

A crucial and preliminary step is the choice of the starting point of the
optimization, given by $x^{0}$, $y^{0}$, and $z^{0}$ (where $x$, $y$, and $z$
are the vectors containing all the corresponding shape and scale parameters,
as defined in Section 3), which appears not to be straightforward since the
missing data prevents the computation of good initial values for the
parameters. Instead of starting from randomly generated values, a better
strategy is to leverage the known information in a reasonable manner. While
the parameters $x_{c\,u}$ of the triage probability distributions are
arbitrarily fixed so that the generated service times are in accordance with
the ED staff suggestions, for the parameters of the medical visit and exams
distributions a different approach is followed. Indeed, in addition to the
timestamps shown in Figure 2.1, which are the most commonly known timestamps
in every ED, for this specific case study further available information is
represented by the time at which the physician requires exams. In general,
this happens during the medical visit, although sometimes exams are required
throughout the patient flow when periodic check-ups are performed. To provide
the parameters $y_{c\,u}$ of the medical visit probability distributions with
initial values, a good starting point can be obtained by computing for each
patient the duration between the start of the visit and the latest time at
which an exam is required. The latter time is considered within a reasonable
period from the start of the visit that should reflect the medical visit’s
maximum service time. Once this timestamp is identified and the durations are
consequently available, their values can be used to initialize the parameters
by fitting the corresponding Weibull probability distributions. Moreover, the
time between the presumptive end of the visit and the discharge can be used in
the same manner to obtain good initial values for the parameters $z_{c\,u}$ of
the exams distributions.

Although Problem 3.1 is a continuous optimization problem, its variables are
considered as discrete to efficiently solve the problem. In particular, at the
time of solving the problem, $x^{c\,u}_{p}$, $y^{c\,u}_{p}$, and
$z^{c\,u}_{p}$, with $p\in\\{1,2\\}$, are treated as granular variables, i.e.,
variables with a controlled number of decimals (see, e.g., audet.2019 ). This
choice is motivated by the fact that the output of the simulation model is
insensitive to small changes in the values of the decision variables, thus
making discrete variables preferable. To this end, a specific granularity is
considered based on the role of the parameter associated with the variable in
the corresponding probability distribution. All the variables are pairs of
shape and scale parameters, which have a different impact on the simulation
output. Indeed, while the former parameters determine the shape of the
probability distribution, the latter affect the scale of the values generated.
With respect to the Weibull distribution adopted, in practice the behavior
observed is that the larger the value of the shape parameter, the larger the
dispersion of the values of the random variables associated. Moreover, such
random variables take on values with a larger scale as the value of the scale
parameter increases. Since some preliminary analyses have shown that the
impact on the simulation output, measured in terms of the two KPIs $DOT$ and
$DIT$, is more noticeable for variations in the scale parameter, for each pair
of variables the granularity $\delta^{min}_{p}$ of the values is fixed to
$10^{-3}$ if $p=1$ (i.e., shape parameter) and to $10^{-4}$ if $p=2$ (i.e.,
scale parameter). Note that treating the variables as granular requires the
new variables to be equal to $x^{c\,u}_{p}/\delta^{min}_{p}\in\mathbb{Z}$,
thus meaning that $x^{c\,u}_{p}$ have to take on real values that are multiple
of $\delta^{min}_{p}$. The same type of constraint applies to the variables of
the visit and exams probability distributions.

To solve the SBO problem for calibrating the simulation model of the ED in
hand, the approach of the Sample Average Approximation (SAA) pagnoncelli.2009
; kim.2015-b  is adopted. As a consequence, the resulting optimization problem
is deterministic and it can be solved by applying an algorithm from the class
of DFO conn.2009 ; audet.2017 . In particular, the optimization algorithm
proposed in liuzzi.2020  is used for solving the problem by adopting its
default parameters. This algorithm has been successfully applied to the same
case study in desantis.2021  and it is suited for efficiently solving integer
black-box constrained problems, like Problem 3.1, thanks to unconventional
search directions, an effective penalty approach, and strong global
convergence properties. The maximum number of function evaluations, which
represents the stopping condition, is set to 3000. Note that by using the SAA
approach, the empirical cumulative distribution functions used in the
optimization problem are estimated through the corresponding sample means over
the 30 independent simulation replications.

As concerns the optimization problem to solve, the values of the lower and
upper bounds introduced in Problem (3.1) are reported in Table 5.1. Moreover,
the tolerance $\tau^{c\,u\,i}_{\mu}$ and $\tau^{c\,u\,i}_{\sigma}$ are both
fixed to $0.35$ for yellow and green-tagged patients visited and treated in MU
and SU, $0.2$ otherwise.

Table 5.1: Values of the lower and upper bounds for each pair of variables $(x_{1}^{c\,u},x_{2}^{c\,u})$, $(y_{1}^{c\,u},y_{2}^{c\,u})$, and $(z_{1}^{c\,u},z_{2}^{c\,u})$, for all $c\in C$ and $u\in U(c)$.  $l_{x^{c\,u}_{1}}$ | $l_{x^{c\,u}_{2}}$ | $l_{y^{c\,u}_{1}}$ | $l_{y^{c\,u}_{2}}$ | $l_{z^{c\,u}_{1}}$ | $l_{z^{c\,u}_{2}}$
---|---|---|---|---|---
0.01 | 0.01 | 0.01 | 0.01 | 0.01 | 0.01
$u_{x^{c\,u}_{1}}$ | $u_{x^{c\,u}_{2}}$ | $u_{y^{c\,u}_{1}}$ | $u_{y^{c\,u}_{2}}$ | $u_{z^{c\,u}_{1}}$ | $u_{z^{c\,u}_{2}}$
1000 | 0.5 | 1000 | 4 | 1000 | 40

Given these settings, the experimental results are obtained by using a PC with
Intel Core i7-4790K quad-core 4.00 GHz Processor and 32 GB RAM. Table 5.2
reports the service time distributions of triage, medical visit, and exams and
reassessments at the final solution determined by the optimization algorithm
(we recall that the base unit for time adopted in the DES model is the hour).

Table 5.2: Service time distributions of triage, medical visit, and exams and reassessments at the final solution determined by the optimization algorithm. Note that the triage service time of green-tagged patients in MIU is missing since in the DES model these patients are subject to the Weibull distributions used for SU and MU. |  | TRIAGE | MEDICAL VISIT | EXAMS
---|---|---|---|---
White | MIU | Weib(0.61, 0.42) | Weib(0.99, 0.57) | Weib(0.62, 0.23)
| MIU | - | Weib(1.06, 0.41) | Weib(1.29, 1.61)
Green | SU | Weib(1000, 0.50) | Weib(0.64, 0.31) | Weib(0.68, 2.11)
| MU | Weib(0.60, 0.49) | Weib(1000, 1.63) | Weib(0.69, 12.26)
Yellow | SU | Weib(1000, 0.28) | Weib(0.51, 0.22) | Weib(0.72, 6.32)
MU | Weib(0.55, 0.50) | Weib(0.62, 7.22) | Weib(0.78, 25.46)
Red | SU | Weib(0.92, 0.10) | Weib(0.88, 0.75) | Weib(1.41, 19.34)
MU | Weib(0.80, 0.19) | Weib(2.24, 0.36) | Weib(0.78, 30.92)

Two types of plots are shown in Appendices A and B to assess the accuracy of
the calibrated simulation model. The former appendix focuses on the KPI given
by the average hourly number of patients inside the ED after the start of the
medical visit, namely, patients satisfying two conditions: $(i)$ they have
already started the medical visit; $(ii)$ they have not been discharged yet.
In particular, the plots report the comparison between the value of the KPI
computed through the real data and the corresponding value derived from the
simulation output along with its $95\%$ confidence interval. The latter
appendix reports the comparison between the Empirical Cumulative Distribution
Functions (ECDFs) of the values of the time differences $DOT$ and $DIT$
computed for each patient through the real data and the simulation output
resulting from the model calibration procedure. In both appendices, the
comparisons are performed for all the color tags $c\in C$ and for all the ED
units $u\in U(c)$ where a patient with tag $c$ can be assigned. Moreover, both
the average hourly number of patients inside the ED after the start of the
medical visit and the ECDFs are obtained as an average over the independent
simulation replications. Both types of plots show that, in general, the values
of the KPIs corresponding to the system data are well reproduced by the
simulation output. However, some dissimilarities may be still present in spite
of the calibration, especially for green and yellow-tagged patients. This is
due to the difficulty in reproducing the patient flow when sharing of the
resources is involved between patients with different triage tag, as is the
case for the medical visit of green and yellow-tagged patients at both $MU$
and $SU$. Moreover, the experimental results show that, in general, the
comparison between real data and simulation output is more accurate for the
time differences DIT. Indeed, larger errors are observed in both types of
plots when time differences DOT are considered due to the difficulty in
reproducing the actual waiting times of patients, which depend both on the
service time of the medical visit and on the number of patients that are in
queue. Hence, it is the result of a seize-delay-release process used within
the simulation, while DIT is given by the mere sum of service times.

By observing the ECDF plots in Appendix B, it is possible to gain useful
insight that histograms, which are commonly adopted in the literature, fail to
provide. For instance, the ECDF related to yellow-tagged patients in SU shows
that the time differences DOT collected from the simulation model are
systematically larger than the corresponding real values, meaning that either
the estimated parameters of the probability distribution of the medical visit
service time or the modeling of the interaction with green-tagged patients in
SU may be improved with further research. To avoid ambiguities, it is
important to remark that since the simulated ECDFs are obtained as an average
over independent simulation replications, the jumps (or steps) in these
simulated functions correspond to patients observed across all the
replications. This is a consequence of averaging different step functions,
such as the ECDFs, which have jumps at different points in their domain.

The same conclusions on the accuracy of the simulation model can be drawn also
by observing Figures 5.1–5.2, which report the comparison between the average
current and simulated values of DOT and DIT for all the triage tags and all
the ED units. As already mentioned, the simulation output related to the time
differences DIT is a good approximation of the real system values since for
each color tag the current values are either within the corresponding
confidence interval or close to it. Instead, even though for the time
differences DOT the accuracy of the simulation responses as estimates of the
current values is lower, the relative error is within the tolerances
$\tau^{c\,u\,i}_{\mu}$ chosen. Indeed, all the constraints of Problem (3.1)
are (of course) satisfied at the optimal solution obtained by applying the
model calibration approach. Note that reducing the values of the tolerances
potentially allows achieving more accurate solutions, however, the
optimization algorithm is not guaranteed to determine a feasible solution.

Figure 5.1: Plot of current values and simulation output of DOT with the
confidence interval. Figure 5.2: Plot of current values and simulation output
of DIT with the confidence interval.

The experimental results discussed above show that the model calibration
procedure enables estimating the missing parameters in order for the
simulation model to satisfactorily reproduce the ED operations despite the
unavailable timestamps, which represent the main hurdle for achieving a high
level of accuracy. The constraints on the sample means used in Problem (3.1)
ensure that the responses of the simulation model are on average close to the
corresponding real values. Moreover, the constraints on the sample variances
of the KPIs help the optimization algorithm to avoid points associated with
inaccurate results. However, when looking more closely into the details of the
numerical results, some dissimilarities between real data and simulation
output may be observed as a result of the impact of the problem of missing
data, which necessarily undermines the overall accuracy. Notwithstanding, the
simulation model may be deemed sufficiently reliable with respect to the
specific objectives of the analysis.

##  6 Conclusions

The SBO approach proposed in this paper addresses one of the data quality
problems that frequently affect the DES models concerning EDs. In particular,
the focus is on the problem of missing data, which consists in the
unavailability of data related to some of the starting and ending times of the
activities performed in the ED. This well-known issue is responsible for the
lack of knowledge of the service time of some processes, which is required to
estimate the corresponding probability distributions to use in the simulation
model. For this purpose, after assuming suitable probability distributions for
the processes associated with the unknown service times, a model calibration
procedure is proposed to estimate the missing parameters of such
distributions.

The preliminary experimental results on the ED of a big hospital in Italy show
that the model calibration procedure enables estimating the missing parameters
in order for the simulation model to satisfactorily reproduce the ED
operations despite the unavailable timestamps, which represent the main hurdle
for achieving a high level of accuracy. Although proper constraints ensure
that the responses of the simulation model are on average close to the
corresponding real values, when looking more closely into the details of the
experimental results, some dissimilarities between real data and simulation
output may be observed as a result of the impact of the problem of missing
data, which necessarily undermines the overall accuracy. Notwithstanding, the
simulation model may be deemed sufficiently reliable with respect to the
specific objectives of the analysis. Since there is still a margin for
improvement, several ideas are considered for further research. For instance,
as the worst results are observed for the time differences DOT, a significant
improvement may be obtained by using different weights for the terms of the
objective function of Problem (3.1) and by assigning larger values to the
terms associated with DOT. Moreover, different objective functions and
starting points may be taken into account also to assess the robustness of the
approach. Finally, increasing the number of function evaluations used as
stopping condition may lead to better solutions allowing a more thorough
exploration of the feasible region.

## Acknowledgment

The authors are grateful to Prof. F. Romano (Chief Medical Officer) and Dr. L.
De Vito (Medical Director of ED) of Policlinico Umberto I of Rome for their
availability to carry out this study.

##  Appendix A Model calibration - Plots I

This appendix focuses on the KPI given by the average hourly number of
patients satisfying two conditions: $(i)$ they have already started the
medical visit; $(ii)$ they have not been discharged yet. In particular, the
plots report the comparison between the value of the KPI computed through the
real data and the corresponding value derived from the simulation output along
with its $95\%$ confidence interval. The comparisons are performed for all
color tags $c\in C$ and for all the ED units $u\in U(c)$ where a patient with
color tag $c$ can be assigned.

##  Appendix B Model calibration - Plots II

This appendix reports the comparison between the Empirical Cumulative
Distribution Functions (ECDFs) of the time differences $DOT$ and $DIT$
computed through the real data and the simulation output resulting from the
calibration procedure. The comparisons are performed for all color tags $c\in
C$ and for all ED units $u\in U(c)$. The colored curves correspond to the
simulation output.

## Conflict of interest

The authors declare that they have no conflict of interest.

## References

  * (1) Ahalt, V., Argon, N., Strickler, J., Mehrotra, A.: Comparison of emergency department crowding scores: a discrete-event simulation approach.  Health Care Management Science 21, 144–155 (2018)
  * (2) Ahmed, M.A., Alkhamis, T.M.: Simulation optimization for an emergency department healthcare unit in Kuwait.  European Journal of Operational Research 198(3), 936 – 942 (2009)
  * (3) Amaran, S., Sahinidis, N., Sharda, B., Bury, S.: Simulation optimization: a review of algorithms and applications.  4OR 12, 301–333 (2014)
  * (4) Audet, C., Hare, W.: Derivative-Free and Blackbox Optimization.  Springer International Publishing (2017).  DOI 10.1007/978-3-319-68913-5
  * (5) Audet, C., Le Digabel, S., Tribes, C.: The mesh adaptive direct search algorithm for granular and discrete variables.  SIAM Journal on Optimization 29, 1164–1189 (2019).  DOI 10.1137/18M1175872
  * (6) Baumlin, K., Shapiro, J., Weiner, C., Gottlieb, B., Chawla, N., Richardson, L.: Clinical information system and process redesign improves emergency department efficiency.  The Joint Commission Journal on Quality and Patient Safety 36(4), 179–185 (2010)
  * (7) Bedoya-Valencia, L., Kirac, E.: Evaluating alternative resource allocation in an emergency department using discrete event simulation.  Simulation 92, 1041–1051 (2016)
  * (8) Conn, A., Scheinber, K., Vicente, L.: Introduction to Derivative-Free Optimization.  Society for Industrial and Applied Mathematics, Philadelphia (2009)
  * (9) De Santis, A., Giovannelli, T., Lucidi, S., Messedaglia, M., Roma, M.: An optimal non–uniform piecewise constant approximation for the patient arrival rate for a more efficient representation of the emergency departments arrival process.  Technical Report 1–2020, Dipartimento di Ingegneria Informatica Automatica e Gestionale “A. Ruberti”, SAPIENZA Università di Roma (2020)
  * (10) De Santis, A., Giovannelli, T., Lucidi, S., Messedaglia, M., Roma, M.: Determining the optimal piecewise constant approximation for the Nonhomogeneous Poisson Process rate of Emergency Department patient arrivals.  arXiv e-prints arXiv:2101.11138 (2021)
  * (11) Di Somma, S., Paladino, L., Vaughan, L., Lalle, I., Magrini, L., Magnanti, M.: Overcrowding in emergency department: an international issue.  Internal and Emergency Medicine 10 (2014).  DOI 10.1007/s11739-014-1154-8
  * (12) Duma, D., Aringhieri, R.: An ad hoc process mining approach to discover patient paths of an emergency department.  Flexible Services and Manufacturing Journal 32 (2018).  DOI 10.1007/s10696-018-9330-1
  * (13) Fava, G., Giovannelli, T., Messedaglia, M., Roma, M.: Effect of different patient peak arrivals on an Emergency Department via discrete event simulation.  arXiv e-prints arXiv:2101.12432 (2021)
  * (14) Fu, M.C. (ed.): Handbook of Simulation Optimization.  Springer, New York (2015)
  * (15) Gosavi, A.: Simulation-Based Optimization: Parametric Optimization Techniques and Reinforcement Learning, 2nd edn.  Springer Publishing Company, Incorporated (2014)
  * (16) Gray, G.A., Kolda, T.G.: Algorithm 856: Appspack 4.0: Asynchronous parallel pattern search for derivative-free optimization.  ACM Trans. Math. Softw. 32(3), 485–507 (2006).  DOI 10.1145/1163641.1163647.  URL https://doi.org/10.1145/1163641.1163647
  * (17) Guo, H., Gao, S., Tsui, K., Niu, T.: Simulation optimization for medical staff configuration at emergency department in Hong Kong.  IEEE Transactions on Automation Science and Engineering 14(4), 1655–1665 (2017)
  * (18) Guo, H., Goldsman, D., Tsui, K.L., Zhou, Y., Wong, S.Y.: Using simulation and optimisation to characterise durations of emergency department service times with incomplete data.  International Journal of Production Research 54(21), 6494–6511 (2016)
  * (19) Hoot, N., Aronsky, D.: Systematic review of emergency department crowding: causes, effects, and solutions.  Annals of Emergency Medicine 52(2), 126–136 (2008)
  * (20) Joshi, V., Lim, C., Teng, S.G.: Simulation study: Improvement for non-urgent patient processes in the emergency department.  Engineering Management Journal 28:3, 145–157 (2016)
  * (21) Kim, S., Pasupathy, R., Henderson, S.G.: A Guide to Sample Average Approximation, pp. 207–243.  Springer New York, New York, NY (2015).  DOI 10.1007/978-1-4939-1384-8_8.  URL https://doi.org/10.1007/978-1-4939-1384-8_8
  * (22) Kim, S.H., Whitt, W.: Are call center and hospital arrivals well modeled by nonhomogeneous Poisson process ?  Manufactory & Service Operations Management 16, 464–480 (2014)
  * (23) Kramer, A., Dosi, C., Iori, M., Vignoli, M.: Successful implementation of discrete event simulation: the case of an italian emergency department (2020)
  * (24) Kuo, Y.H., Rado, O., Lupia, B., Leung, J.M.Y., Graham, C.A.: Improving the efficiency of a hospital emergency department: a simulation study with indirectly imputed service-time distributions.  Flexible Services and Manufacturing Journal 28(1), 120–147 (2016)
  * (25) Law, A.: Simulation Modeling and Analysis, fifth edn.  McGraw-Hill (2015)
  * (26) Liu, Z., Rexachs, D., Epelde, F., Luque, E.: A simulation and optimization based method for calibrating agent-based emergency department models under data scarcity.  Comput. Ind. Eng. 103(C), 300–309 (2017).  DOI 10.1016/j.cie.2016.11.036.  URL https://doi.org/10.1016/j.cie.2016.11.036
  * (27) Liuzzi, G., Lucidi, S., Rinaldi, F.: An algorithmic framework based on primitive directions and nonmonotone line searches for black-box optimization problems with integer variables.  Mathematical Programming Computation 12(4), 673–702 (2020).  DOI 10.1007/s12532-020-00182-7.  URL https://doi.org/10.1007/s12532-020-00182-7
  * (28) Nahhas, A., Alwadi, A., Reggelin, T.: Simulation and the emergency department overcrowding problem.  Procedia Engineering 178, 368–376 (2017)
  * (29) Pagnoncelli, B.K., Ahmed, S., Shapiro, A.: Sample Average Approximation Method for Chance Constrained Programming: Theory and Applications.  Journal of Optimization Theory and Applications 142(2), 399–416 (2009).  DOI 10.1007/s10957-009-9523-6
  * (30) Paul, S.A., Reddy, M.C., De Flitch, C.J.: A systematic review of simulation studies investigating emergency department overcrowding.  Simulation 86(8-9), 559–571 (2010)
  * (31) Pines, J.M., Iyer, S., Disbot, M., Hollander, J.E., Shofer, F.S., Datner, E.M.: The effect of emergency department crowding on patient satisfaction for admitted patients.  Academic Emergency Medicine 15(9), 825–831 (2008).  DOI https://doi.org/10.1111/j.1553-2712.2008.00200.x
  * (32) Salmon, A., Rachuba, S., Briscoe, S., Pitt, M.: A structured literature review of simulation modelling applied to emergency departments: Current patterns and emerging trends.  Operations Research for Health Care 19, 1–13 (2018)
  * (33) Ucar, I., Smeets, B., Azcorra, A.: simmer: Discrete-event simulation for R.  Journal of Statistical Software 90(2), 1–30 (2019).  DOI 10.18637/jss.v090.i02
  * (34) Vanbrabant, L., Martin, N., Ramaekers, K., Braekers, K.: Quality of input data in emergency department simulations: Framework and assessment techniques.  Simulation Modelling Practice and Theory 91, 83 – 101 (2019)
  * (35) Whitt, W., Zhang, X.: A data-driven model of an emergency department.  Operations Research for Health Care 12, 1 – 15 (2017)
  * (36) Wiler, J., Griffey, R., Olsen, T.: Review of modeling approaches for emergency department patient flow and crowding research.  Academic Emergency Medicine 18, 1371–1379 (2011)
  * (37) Wong, Z.S.Y., Lit, A.C.H., Leung, S.Y., Tsui, K.L., Chin, K.S.: A discrete-event simulation study for emergency room capacity management in a hong kong hospital.  2016 Winter Simulation Conference (WSC) pp. 1970–1981 (2016)
  * (38) Yousefi, M., Yousefi, M., Fogliatto, F.S.: Simulation-based optimization methods applied in hospital emergency departments: A systematic review.  SIMULATION 96(10), 791–806 (2020).  DOI 10.1177/0037549720944483
  * (39) Zeinali, F., Mahootchi, M., Sepehri, M.: Resource planning in the emergency departments: a simulation-base metamodeling approach.  Simulation Modelling Practice and Theory 53, 123–138 (2015)


$Z^{-1}AZ=\begin{bmatrix}A_{1}&A_{21}&A_{22}\\\ 0&A_{41}&0\\\
0&0&A_{44}\end{bmatrix}.$

Then $\sigma(A_{44})\cap\Omega=\varnothing$ implies that

$\widehat{A}:=Z^{-1}AZ\sim\begin{bmatrix}A_{1}&A_{21}&0\\\ 0&A_{41}&0\\\
0&0&A_{44}\end{bmatrix}.$

Similarly, $B$ is similar to

$\widehat{B}:=\begin{bmatrix}B_{1}&0&0\\\ B_{21}&B_{41}&0\\\
0&0&B_{44}\end{bmatrix},$

where, as we noted above, $A_{44}$ is similar to $B_{44}$.

Now

  1. (i)

$\widehat{A}\oplus\widehat{A}$ is similar to $\widehat{B}\oplus\widehat{B}$;

  2. (ii)

$\sigma(\begin{bmatrix}A_{1}&A_{21}\\\
0&A_{41}\end{bmatrix})=\Omega=\sigma(\begin{bmatrix}B_{1}&0\\\
B_{21}&B_{41}\end{bmatrix});\mbox{\ \ \ \ \ and }$

  3. (iii)

$\Omega\cap\sigma(A_{44})=\varnothing$,
$\Omega\cap\sigma(B_{44})=\varnothing$.

It follows that $(\begin{bmatrix}A_{1}&A_{21}\\\
0&A_{41}\end{bmatrix})^{(2)}\sim(\begin{bmatrix}B_{1}&0\\\
B_{21}&B_{41}\end{bmatrix})^{(2)}$. But these are finite matrices, so by
considering their Jordan forms, we immediately conclude that

$\begin{bmatrix}A_{1}&A_{21}\\\
0&A_{41}\end{bmatrix}\sim\begin{bmatrix}B_{1}&0\\\
B_{21}&B_{41}\end{bmatrix},$

from which we deduce that $\widehat{A}\sim\widehat{B}$, and thus $A$ is
similar to $B$.

$\Box$

When $A$ and $B$ are not assumed to be compact, we may no longer conclude that
$S_{1}$ is Fredholm of index zero, and the above proof immediately breaks
down. Indeed, the question of whether or not $A\oplus A$ similar to $B\oplus
B$ implies that $A$ is similar to $B$ seems to be very difficult. There exist
some cases, however, where we can assert positive results. First, we propose a
possible “line of attack” to solving the problem, in the hope that the reader
may have more luck than we did.

### 4.11.

Proposition. Let $A,B\in\mathcal{B}(\mathcal{H})$.

  1. (a)

Every operator $Z\in\mathcal{B}(\mathcal{H})$ in the commutant of $J_{2}(A)$
is of the form

$\begin{bmatrix}Z_{1}&Z_{2}\\\ Z_{2}A&Z_{1}\\\ \end{bmatrix},$

where $Z_{1}$ and $Z_{2}$ commute with $A$.

  2. (b)

Suppose that $S\in\mathcal{B}(\mathcal{H}\oplus\mathcal{H})$ is invertible and
that $J_{2}(B)=SJ_{2}(A)S^{-1}$. By Lemma 4.9(a), we may write
$S=\begin{bmatrix}S_{1}&S_{2}\\\ S_{2}A&S_{1}\end{bmatrix}$.

Then $A$ is similar to $B$ if and only if there exists an invertible operator
$Z=\begin{bmatrix}Z_{1}&Z_{2}\\\ Z_{3}&Z_{4}\end{bmatrix}$ in the commutant of
$J_{2}(A)$ such that

$W_{1}:=S_{2}AZ_{2}+S_{1}Z_{1}$

is invertible.

Proof.

  1. (a)

This follows immediately from Lemma 4.9 (a), setting $B=A$ in that Lemma.

  2. (b)

First, suppose that there exists an invertible operator
$Z=\begin{bmatrix}Z_{1}&Z_{2}\\\ Z_{3}&Z_{4}\end{bmatrix}$ in the commutant of
$J_{2}(A)$ such that $W_{1}:=S_{2}AZ_{2}+S_{1}Z_{1}$ is invertible. Since $Z$
commutes with $J_{2}(A)$, we observe that

$J_{2}(B)=SJ_{2}(A)S^{-1}=SZJ_{2}(A)Z^{-1}S^{-1}.$

By part (a) above, we have that $Z_{3}=Z_{2}A$ and $Z_{4}=Z_{1}$, so that
$Z=\begin{bmatrix}Z_{1}&Z_{2}\\\ Z_{2}A&Z_{1}\end{bmatrix}$. Now, by Lemma
4.9(a),

$SZ=\begin{bmatrix}S_{1}&S_{2}\\\
S_{2}A&S_{1}\end{bmatrix}\begin{bmatrix}Z_{1}&Z_{2}\\\
Z_{2}A&Z_{1}\end{bmatrix}=\begin{bmatrix}W_{1}&W_{2}\\\
W_{2}A&W_{1}\end{bmatrix}$

where $W_{2}=S_{1}Z_{2}+S_{2}Z_{1}$. By considering the $(1,1)$-entry of
$(SZ)J_{2}(A)=J_{2}(B)(SZ)$, we see that $W_{1}A=BW_{1}$. Since $W_{1}$ is
invertible by hypothesis, $B$ is similar to $A$.

Conversely, suppose that $B=W_{1}AW_{1}^{-1}$, where
$W_{1}\in\mathcal{B}(\mathcal{H})$ is invertible. With $W:=W_{1}\oplus W_{1}$,
it follows that $J_{2}(B)=WJ_{2}(A)W^{-1}$.

Next, suppose that $S\in\mathcal{B}(\mathcal{H}\oplus\mathcal{H})$ and

$SJ_{2}(A)S^{-1}=J_{2}(B)=WJ_{2}(A)W^{-1}.$

Then $S^{-1}WJ_{2}(A)=J_{2}(A)S^{-1}W$, and so $Z:=S^{-1}W$ is invertible and
lies in the commutant of $J_{2}(A)$. From (a) above, we may write
$Z=\begin{bmatrix}Z_{1}&Z_{2}\\\ Z_{2}A&Z_{1}\end{bmatrix}$ where
$Z_{1},Z_{2}\in\mathcal{B}(\mathcal{H})$ commute with $A$. Using Lemma 4.9(a),
we see that we may also write $S=\begin{bmatrix}S_{1}&S_{2}\\\
S_{2}A&S_{1}\end{bmatrix}$, and therefore

$\begin{bmatrix}W_{1}&0\\\
0&W_{1}\end{bmatrix}=W=SZ=\begin{bmatrix}S_{1}Z_{1}+S_{2}Z_{2}A&S_{1}Z_{2}+S_{2}Z_{1}\\\
S_{2}AZ_{1}+S_{1}Z_{2}A&S_{2}AZ_{2}+S_{1}Z_{1}\end{bmatrix}.$

In particular, $W_{1}$ is invertible and

$W_{1}=S_{1}Z_{1}+S_{2}Z_{2}A=S_{2}AZ_{2}+S_{1}Z_{1}.$

$\Box$

If the spectrum of $A$ (and hence of $B$) is sufficiently nice, we can again
obtain positive results.

### 4.12.

Proposition. Suppose that $A,B\in\mathcal{B}(\mathcal{H})$ are invertible, and
that there exists a net $(p_{\alpha})_{\alpha}$ of polynomials such that
$X:=\textsc{wot}-\lim_{\alpha}p_{\alpha}(A)$ and
$Y:=\textsc{wot}-\lim_{\alpha}p_{\alpha}(B)$ both exist, $X^{2}=A$, and
$Y^{2}=B$. If

$J_{2}(A)\sim J_{2}(B),$

then $A\sim B$.

Proof. Note that $(\sigma(X))^{2}=\sigma(X^{2})=\sigma(A)$ implies that $X$ is
invertible. Now

$J_{2}(A)\sim\begin{bmatrix}I&0\\\ 0&X^{-1}\end{bmatrix}\begin{bmatrix}0&I\\\
A&0\end{bmatrix}\ \begin{bmatrix}I&0\\\ 0&X\end{bmatrix}=\begin{bmatrix}0&X\\\
X^{-1}A&0\end{bmatrix}=\begin{bmatrix}0&X\\\ X&0\end{bmatrix}\sim X\oplus-X.$

Similarly, $J_{2}(B)\sim Y\oplus-Y$. Thus

$X\oplus-X\sim J_{2}(A)\sim J_{2}(B)\sim Y\oplus-Y.$

Choose $R\in\mathcal{B}(\mathcal{H}\oplus\mathcal{H})$ invertible such that

$R^{-1}(X\oplus-X)R=Y\oplus-Y.$

Squaring both sides, we obtain

$R^{-1}(A\oplus A)R=B\oplus B.$

Let $(p_{\alpha})_{\alpha}$ be a net of polynomials such that
$\textsc{wot}-\lim_{\alpha}p_{\alpha}(A)=X$. Then

$\displaystyle Y\oplus Y$
$\displaystyle=\textsc{wot}-\lim_{\alpha}p_{\alpha}(B\oplus B)$
$\displaystyle=\textsc{wot}-\lim_{\alpha}p_{\alpha}(R^{-1}(A\oplus A)R)$
$\displaystyle=R^{-1}(\textsc{wot}-\lim_{\alpha}(p_{\alpha}(A\oplus A)))R$
$\displaystyle=R^{-1}(X\oplus X)R.$

Thus

$\displaystyle R^{-1}(2X\oplus 0)R$ $\displaystyle=R^{-1}((X\oplus-X)+(X\oplus
X))R$ $\displaystyle=R^{-1}(X\oplus-X)R+R^{-1}(X\oplus X)R$
$\displaystyle=(Y\oplus-Y)+(Y\oplus Y)$ $\displaystyle=2Y\oplus 0.$

Since $X$ and $Y$ are invertible, this implies (by Lemma 2.16) that $X$ is
similar to $Y$, and therefore that $A=X^{2}$ is similar to $Y^{2}=B$.

$\Box$

### 4.13.

Corollary. Let $A,B\in\mathcal{B}(\mathcal{H})$ be invertible operators, and
suppose that there exists a simply connected domain
$\Omega\subseteq\mathbb{C}\setminus\\{0\\}$ such that
$\sigma(A)\subseteq\Omega$. If $J_{2}(A)\sim J_{2}(B)$, then $A$ is similar to
$B$.

Proof. First note that – as we have already seen – $J_{2}(A)\sim J_{2}(B)$
implies that $\sigma(B)=\sigma(A)$.

Since $0\not\in\Omega$, it follows that the identity function $\iota(z)=z$,
$z\in\Omega$ admits a holomorphic square root, say $q$. Thus $q(A)$ is a norm-
limit of polynomials in $A$, and $q(B)$ is a norm-limit of the same
polynomials in $B$. The result now follows from Proposition 4.12.

$\Box$

In the case where one of the operators $A,B$ is an isometry, we can dispense
with the stronger hypothesis that $J_{2}(A)$ is similar to $J_{2}(B)$ and we
can obtain a positive solution to Kaplansky’s original problem in the context
of similarity (see Theorem 4.17 below). Before doing so, however, we need to
recall a few definitions and results from the literature. The monograph [22]
is an excellent reference.

### 4.14.

Definition. Recall that an operator $T\in\mathcal{B}(\mathcal{H})$ is said to
be polynomially bounded if there exists a constant $\kappa>0$ such that for
all polynomials $p\in\mathbb{C}[x]$, we have

$\|p(T)\|\leq\kappa\|p\|_{\infty},$

where $\|p\|_{\infty}:=\sup\\{|p(z)|:|z|\leq 1\\}.$ We say that $T$ is
completely polynomially bounded if there exists a constant $\kappa>0$ such
that for all $q=[q_{ij}]\in\mathbb{M}_{n}(\mathbb{C}[x])$, we have

$\|[q_{ij}(T)]\|\leq\kappa\|q\|_{\infty},$

where $\|q\|_{\infty}:=\sup\\{\|q_{ij}(z)\|:|z|\leq 1\\}$.

### 4.15.

We note that von Neumann’s Inequality [22, Corollary 1.2] is the assertion
that every contraction $T\in\mathcal{B}(\mathcal{H})$ is polynomially bounded,
while Arveson [2, 3] observed that the Sz.-Nagy Dilation theorem implies that
every operator $X\in\mathcal{B}(\mathcal{H})$ which is similar to a
contraction is in fact completely polynomially bounded.

The converse of this, namely that every completely polynomially bounded
operator is similar to a contraction, is due to Paulsen [21].

We also remind the reader that an operator $T\in\mathcal{B}(\mathcal{H})$ is
said to be quasinormal if $T$ commutes with $T^{*}T$. It follows from the work
of Brown [7] that $T$ is quasinormal if and only if $T$ is unitarily
equivalent to an operator of the form $N\oplus(S\otimes P)$, where $S$ denotes
the unilateral forward shift, $N$ is normal and $P$ is a positive definite
operator.

If $W\in\mathcal{B}(\mathcal{H})$ is an isometry, then $W^{*}W=I$, and so
clearly $W$ is quasinormal. The Wold Decomposition for isometries acting on
infinite-dimensional, separable Hilbert spaces asserts that there exists a
unitary operator $U_{W}$ (perhaps acting on a space of dimensional zero) and a
unique cardinal number $\alpha\in\mathbb{N}\cup\\{0,\aleph_{0}\\}$ such that
$W\simeq U_{W}\oplus S^{(\alpha)}$. Observing that
$S^{(\alpha)}\simeq(S\otimes I_{\alpha})$, where $\alpha$ denotes the identity
operator acting on a Hilbert space of dimension $\alpha$, we see for
isometries, the corresponding positive definite operator in Brown’s structure
theorem is nothing more than $I_{\alpha}$.

We shall need one more result, due to Popescu [23].

### 4.16.

Theorem. _[Popescu]_ Let $A\in\mathcal{B}(\mathcal{H})$ be such that

$P_{A}:=\textsc{sot}-\lim_{m}(A^{*})^{m}A^{m}$

exists. Then $A$ is similar to an isometry if and only if $P_{A}$ is
invertible.

### 4.17.

Theorem. Let $T\in\mathcal{B}(\mathcal{H})$ and $n\geq 1$ be an integer and
let $W\in\mathcal{B}(\mathcal{H})$ be an isometry. If $T^{(n)}$ is similar to
$W^{(n)}$, then $T$ is similar to $W$. In particular, if $J_{n}(T)$ is similar
to $J_{n}(W)$, then $T$ is similar to $W$.

Proof. Note that $\|W\|=1=\|W^{(n)}\|$. It follows from Arveson’s results that
$T^{(n)}$ is completely polynomially bounded, which is easily seen to imply
that $T$ is completely polynomially bounded, and thus by Paulsen’s Theorem,
that $T$ is similar to a contraction $X$.

If $A\in\mathcal{B}(\mathcal{H})$ is a contraction, then the sequence
$((A^{*})^{m}A^{m})_{m}$ is a decreasing sequence of positive operators, and
so $P_{A}:=\textsc{sot}-\lim_{m}(A^{*})^{m}A^{m}$ exists. In our case, we
conclude that

$P_{X}:=\textsc{sot}-\lim_{m}(X^{*})^{m}X^{m}\mbox{\ \ \ \ and that \ \ \ \
}P_{X^{(n)}}:=\textsc{sot}-\lim_{m}((X^{(n)})^{*})^{m}(X^{(n)})^{m}$

both exist, and clearly we must have

$P_{X^{(n)}}=(P_{X})^{(n)}.$

But $X^{(n)}$ is similar to the isometry $W^{(n)}$, and so from Theorem 4.16,
$P_{X^{(n)}}$ is invertible. This in turn implies that $P_{X}$ is invertible.
Once again applying Theorem 4.16, we conclude that $X$ is similar to an
isometry $V$, and consequently that $T$ is similar to $V$. By the Wold
decomposition, $V$ is unitarily equivalent to $U_{V}\oplus S^{(\alpha)}\simeq
U_{V}\oplus(S\otimes I_{\alpha})$, where $U_{V}$ is unitary
$\alpha\in\mathbb{N}\cup\\{0,\aleph_{0}\\}$, and $W$ is unitarily equivalent
to $U_{W}\oplus S^{(\beta)}\simeq U_{W}\oplus(S\otimes I_{\beta})$, where
$U_{W}$ is unitary and $\beta\in\mathbb{N}\cup\\{0,\beta\\}$. In particular,

$U_{V}^{(n)}\oplus(S\otimes I_{n\cdot\alpha})\simeq V^{(n)}\sim T^{(n)}\sim
W^{(n)}\simeq U_{W}^{(n)}\oplus(S\otimes I_{n\cdot\beta}).$

Now $V^{(n)}$ and $W^{(n)}$ are quasinormal, and so by [9, Theorem 2.1],
$U_{V}^{(n)}\simeq U_{W}^{(n)}$ and for all $\lambda\in\mathbb{C}$,

$\mathrm{nul}(I_{n\cdot\alpha}-\lambda
I_{n\cdot\alpha})=\mathrm{nul}(I_{n\cdot\beta}-\lambda I_{n\cdot\beta}).$

By the results of Kadison and Singer [16] (see Theorem 4.2 above),
$U_{V}\simeq U_{W}$, and the condition on the nullity for $\lambda=1$ implies
that $\alpha=\beta$. Thus $V$ is unitarily equivalent to $W$, and so $T$ is
similar to $W$.

The last statement easily follows from the above argument and the fact that if
$J_{n}(T)\sim J_{n}(W)$, then $T^{(n)}=J_{n}(T)^{n}\sim J_{n}(W)^{n}=W^{(n)}$.

$\Box$

##  5\. Remarks and observations.

### 5.1.

Let $T\in\mathcal{B}(\mathcal{H})$ and $2\leq n\in\mathbb{N}$. As we have seen
multiple times, $J_{n}(T)^{n}=T^{(n)}$, and thus if $T$ is
$\textsc{j}_{n}^{[\textsc{u}]}$-stable (resp.
$\textsc{j}_{n}^{[\textsc{s}]}$-stable), then $T^{n}\simeq T^{(n)}$ (resp.
$T^{n}\sim T^{(n)}$). In their paper [11], the authors studied operators
satisfying this last relation in the case where $n=2$, referring to these
phenomena as condition U (i.e. $T^{2}\simeq T^{(2)}$) and as condition S (i.e.
$T^{2}\sim T^{(2)}$) respectively. Unfortunately, they mistakenly state that
“a scalar multiple of an operator satisfying condition S (or condition U)
satisfies the same condition only if the scalar is either $0$ or $1$.”

A counterexample to this claim comes in the form of the forward unilateral
shift operator $S$, as for any $\alpha\in\mathbb{T}$,

$(\alpha S)^{2}=\alpha^{2}S^{2}\simeq S^{2}\simeq S\oplus S\simeq(\alpha
S\oplus\alpha S).$

### 5.2.

Example. There exist (positive, invertible) operators
$H_{1},H_{2}\in\mathcal{B}(\mathcal{H})$ such that $H_{1}\sim H_{2}$ and
$H_{1}\simeq_{a}H_{2}$, but $H_{1}$ is not unitarily equivalent to $H_{2}$.

Choose two non-atomic Borel measures $\mu_{1},\mu_{2}$ on the interval
$[1,2]\subseteq\mathbb{R}$, with the supports of $\mu_{1}$ and $\mu_{2}$ both
equal to $[1,2]$, such that $\mu_{1}$ and $\mu_{2}$ are not mutually
absolutely continuous. Let $H_{k}$ denote the (hermitian) multiplication
operator $M_{x}$ acting on $L^{2}([1,2],\mu_{k}),\ k=1,2$.

Then $\sigma(H_{1})=[1,2]=\sigma(H_{2})$, and so $H_{1}$ and $H_{2}$ are
approximately unitarily equivalent, by the Weyl-von Neumann-Berg Theorem [12,
Theorem II.4.4], but they are not unitarily equivalent, since the underlying
measures $\mu_{1}$ and $\mu_{2}$ are not mutually absolutely continuous [12,
Corollary II.3.6].

Let $Z_{1}=S\otimes H_{1}$, $Z_{2}=S\otimes H_{2}$, where $S$ is the usual
unilateral forward shift operator. It is easy to check that $Z_{1}$ and
$Z_{2}$ are approximately unitarily equivalent since $H_{1}$ and $H_{2}$ are.
Furthermore, it follows from [9, Theorem 2.1] that $Z_{1}$ and $Z_{2}$ are
similar.

Note, however, that $Z_{1}$ and $Z_{2}$ are not unitarily equivalent. For
suppose otherwise. Then

$H_{1}^{2}\oplus
0=(Z_{1}^{*}Z_{1}-Z_{1}Z_{1}^{*})\simeq(Z_{2}^{*}Z_{2}-Z_{2}Z_{2}^{*})=H_{2}^{2}\oplus
0,$

and so by [16, Theorem 1], we have that $H_{1}^{2}\simeq H_{2}^{2}$. Since
$H_{1},H_{2}$ are invertible and positive, we conclude that $H_{1}\simeq
H_{2}$, a contradiction.

### 5.3.

Remark. Let $A,B\in\mathcal{B}(\mathcal{H})$. The problem of whether or not
$J_{2}(A)\sim J_{2}(B)$ implies that $A\sim B$ remains open. In the hope that
it may be of use to someone wishing to try their hand at this problem, we make
a minor observation, which may be of some use: if $J_{2}(A)\sim J_{2}(B)$ and
$\kappa\in\mathbb{C}$, then $J_{2}(\kappa A)\sim J_{2}(\kappa B)$.

The result is obvious if $\kappa=0$. If $\kappa\not=0$, then $J_{2}(A)\sim
J_{2}(B)$ implies that

$\kappa^{1/2}J_{2}(A)\sim\kappa^{1/2}J_{2}(B).$

But then

$J_{2}(\kappa A)\sim\kappa^{1/2}J_{2}(A)\sim\kappa^{1/2}J_{2}(B)\sim
J_{2}(\kappa B).$

### 5.4.

Example. (due to J. Bell) It is natural to ask whether Kaplansky’s second
problem can be positively answered in the general setting of rings, which
would obviously imply a positive answer in the setting of operators on a
Hilbert space. More specifically, is it true that if $R$ is a unital ring, and
$x,y\in R$ are such that $x\oplus x$ is _conjugate_ (the ring-theoretic
nomenclature for similarity) to $y\oplus y$ in $\mathbb{M}_{2}(R)$, then $x$
and $y$ must be conjugate? In general, the answer is no. In fact, in the
general setting of unital rings, the analogue of Theorem 4.10 fails. The
following counterexample is due to our colleague Jason Bell, and we thank him
for allowing us to reproduce it here.

There exists a unital ring $R$ and elements $x,y\in R$ such that if
$J_{2}(x):=\begin{bmatrix}0&1\\\ x&0\end{bmatrix}$ and
$J_{2}(y):=\begin{bmatrix}0&1\\\ y&0\end{bmatrix}$ then $J_{2}(x)$ is
conjugate to $J_{2}(y)$ in $\mathbb{M}_{2}(R)$, but $x$ is not conjugate to
$y$ in $R$.

Let $S=\mathbb{Q}[a,b,c]/(a^{2}-bc=1)$. Notice that $S$ is an integral domain
since $a^{2}-bc$ is irreducible. We claim that the group of units of $S$ is
$\mathbb{Q}^{*}$. To see this, notice that $S$ admits a $\mathbb{Z}$-grading,
by giving $a$ degree $0$, $b$ degree $1$ and $c$ degree $-1$. Since $S$ is a
$\mathbb{Z}$-graded domain, the only units of $S$ are contained in $S_{0}$,
the homogeneous part of $S$ of degree $0$. Indeed, since $S$ is a domain the
set of units is homogeneous. Consider $S$ as a module over $S_{0}$, and let
$S_{n}$ denote the set of elements of $S$ of degree $n$, $n\in\mathbb{Z}$.
Then for $n\geq 1$, $S_{n}=S_{0}b^{n}$, and $S_{-n}=S_{0}c^{n}$. Note further
that $b$ and $c$ are not units since $bc=a^{2}-1$ is not a unit in $S_{0}$.
Hence the group of units of $S$ lies in $S_{0}$.

Observe that $S_{0}$ is the subring of $S$ generated by $a$ and $bc$. Thus

$S_{0}\cong\mathbb{Q}[a,bc]/(a^{2}-bc-1)\cong\mathbb{Q}[a].$

Since the group of units of $\mathbb{Q}[a]$ is $\mathbb{Q}^{*}$, the claim
follows.

Next, let $T$ denote the free product (in the category of associative
$\mathbb{Q}$-algebras) of $S$ and $\mathbb{Q}\\{x,y\\}$, the free associative
algebra on $x$ and $y$. Consider the ideal $I$ of $T$ generated by $ax-ya$,
$bx-yb$, and $cx-yc$. If we let $R=T/I$, then in $\mathbb{M}_{2}(R)$, the
matrix $\begin{bmatrix}a&b\\\ c&a\end{bmatrix}$ is invertible (since it lies
in $\mathbb{M}_{2}(S)$, $S$ is commutative, and the determinant is
$a^{2}-bc=1$). Furthermore,

$\begin{bmatrix}a&b\\\ c&a\end{bmatrix}(x\oplus x)=(y\oplus
y)\begin{bmatrix}a&b\\\ c&a\end{bmatrix}\in\mathbb{M}_{2}(R),$

by our definition of $I$. Thus $x\oplus x$ and $y\oplus y$ are conjugate in
$\mathbb{M}_{2}(R)$.

Notice that $R$ admits an $\mathbb{N}$-grading (not compatible with the one on
$S$) obtained by assigning degree $0$ to $a,b$ and $c$, and degree $1$ to $x$
and $y$, since the generators of $I$ and the relation $a^{2}-bc=1$ are all
homogeneous. Now $R/(x,y)\cong S$ and so $S$ is the degree $0$ component of
$R$.

We now argue by contradiction. Suppose that $x$ and $y$ are conjugate in $R$.
Then there is some unit $u$ of $R$ such that $ux=yu$. Write
$u=u_{0}+\cdots+u_{d}$, where $u_{i}\in R$ is homogeneous of degree $i$. Let
$v=v_{0}+\cdots+v_{e}$ be the inverse of $u$. By considering the degree $0$
terms of both sides of $1=uv$, we see that $u_{0}v_{0}=1$ and since
$u_{0},v_{0}$ are in the degree $0$ component of $R\cong S$, we see that
$u_{0},v_{0}$ must be in $\mathbb{Q}^{*}$.

Now we have that $(u_{0}+\dots+u_{d})x=y(u_{0}+\dots+u_{d})$ and so looking at
the degree one part of both sides we see that $u_{0}x=yu_{0}$. But $u_{0}$ is
a nonzero scalar, from which it follows that $x=y$ in $R$. Thus to show that
$x$ and $y$ are not conjugate in $R$, it suffices to show that $x\neq y$ in
$R$. To do this, we shall define a homomorphism
$\Phi:R\to\mathbb{M}_{2}(\mathbb{Q})$ such that $\Phi(x)\neq\Phi(y)$.

We consider the map $\Phi$ induced by setting

  * •

$\Phi(1)=\begin{bmatrix}1&0\\\ 0&1\end{bmatrix}$;

  * •

$\Phi(a)=\begin{bmatrix}0&1\\\ 1&0\end{bmatrix}$;
$\Phi(b)=\Phi(c)=\begin{bmatrix}0&0\\\ 0&0\end{bmatrix}$;

  * •

$\Phi(x)=\begin{bmatrix}1&0\\\ 0&-1\end{bmatrix}$; and
$\Phi(y)=\begin{bmatrix}-1&0\\\ 0&1\end{bmatrix}$.

By definition of $R$ as the quotient of the free product $T$ of $S$ and
$\mathbb{Q}\\{x,y\\}$ by the ideal $I$, to see that $\Phi$ extends to a ring
homomorphism on $R$, it suffices to show that

$\displaystyle 0$ $\displaystyle=\Phi(a)\Phi(x)-\Phi(y)\Phi(a)$
$\displaystyle=\Phi(b)\Phi(x)-\Phi(y)\Phi(b)$
$\displaystyle=\Phi(c)\Phi(x)-\Phi(y)\Phi(c),\mbox{ and }$
$\displaystyle\Phi(1)$ $\displaystyle=\Phi(a)^{2}-\Phi(b)\Phi(c),$

all of which are routine calculations.

This implies that $\Phi$ indeed induces a homomorphism from
$R\to\mathbb{M}_{2}(\mathbb{Q})$. But $\Phi(x)\neq\Phi(y)$, and thus $x\neq y$
in $R$, completing the argument.

As a corollary to this example, we note that $x\oplus x=J_{2}(x)^{2}$ is
conjugate to $y\oplus y=J_{2}(y)^{2}$, and yet $x$ is not conjugate to $y$,
and so this example also shows that this formulation of Kaplansky’s second
problem admits a negative answer in the setting of unital rings.

Statements and Declarations

Data Availability There is no dataset associated to this paper, as it is a
study in pure mathematics.

Competing Interests There are no competing interests for the authors to
disclose. The authors contribute equally to this piece of theoretical work.

## References

  * [1] C. Apostol, C. Foiaş, and D. Voiculescu.  Some results on nonquasitriangular operators. IV.  Rev. Roumaine Math. Pures Appl., 18:487–514, 1973.
  * [2] W.B. Arveson.  Subalgebras of ${C}^{*}$-algebras.  Acta Math., 123:141–224, 1969.
  * [3] W.B. Arveson.  Subalgebras of ${C}^{*}$-algebras, II.  Acta Math., 128:271–308, 1972.
  * [4] E.A. Azoff.  Test problems for operator algebras.  Trans. Amer. Math. Soc., 347:2989–3001, 1995.
  * [5] H. Bauer.  Measure and integration theory.  Berlin, New York: De Gruyter, 2001.
  * [6] V.I. Bogachev.  Weak convergence of measures, volume 234 of Mathematical Surveys and Monographs.  American Mathematical Society, Providence, R.I., 2018.
  * [7] A. Brown.  On a class of operators.  Proc. Amer. Math. Soc., 4:723–728, 1953.
  * [8] L.G. Brown.  The universal coefficient theorem for Ext and quasidiagonality.  In Operator algebras and group representations, (Neptun, 1980), volume 1 of Monogr. Stud. Math., pages 60–64, Boston, MA, 1984. Pitman.
  * [9] K.Y. Chen, D.A. Herrero, and P.Y. Wu.  Similarity and quasisimilarity of quasinormal operators.  J. Operator Th., 27:385–412, 1992.
  * [10] J.B. Conway.  A course in functional analysis, volume 96.  Springer Verlag, New York, 1985.
  * [11] J.B. Conway, G. Prǎjiturǎ, and A. Rodríguez-Martínez.  Powers and direct sums.  J. Math. Anal. Appl., 413:880–889, 2014.
  * [12] K.R. Davidson.  ${C}^{*}$-algebras by example, volume 6 of Fields Institute Monographs.  Amer. Math. Soc., Providence, RI, 1996.
  * [13] K. Falconer.  Sub-self-similar sets.  Trans. Amer. Math. Soc., 347:3121–3129, 1995.
  * [14] K. Guo, Y. Ji, and S. Zhu.  A $C^{*}$-algebra approach to complex symmetric operators.  Trans. Amer. Math. Soc., 367:6903–6942, 2015.
  * [15] D.A. Herrero.  Approximation of Hilbert space operators I, volume 224 of Pitman Research Notes in Math. Longman Scientific and Technical, Harlow, New York, second edition, 1989\.
  * [16] R.V. Kadison and I.M. Singer.  Three test problems in operator theory.  Pacific J. Math., 7:1101–1106, 1957.
  * [17] I. Kaplansky.  Infinite abelian groups.  University of Michigan Press, Ann Arbor, MI, 1954.
  * [18] L.W. Marcoux.  On the distance between unitary orbits of weighted shifts.  Trans. Amer. Math. Soc., 326:585–612, 1991.
  * [19] L.W. Marcoux, H. Radjavi, and Y. Zhang.  Around the closure of the set of commutators of idempotent elements of $\mathcal{B}(\mathcal{H}$): biquasitriangularity and factorisation.  J. Funct. Anal. 284: Paper No. 109854, 2023.
  * [20] D.P. O’Donovan.  Weighted shifts and covariance algebras.  Trans. Amer. Math. Soc., 208:1–25, 1975.
  * [21] V. Paulsen.  Every completely polynomially bounded operator is similar to a contraction.  J. Funct. Anal., 55:1–17, 1984.
  * [22] V. Paulsen.  Completely Bounded Maps and Operator Algebras, volume 78 of Cambridge Studies in Advanced Mathematics.  Cambridge University Press, Cambridge, 2002.
  * [23] G. Popescu.  On similarity of operators to isometries.  Michigan Math. J., 39:385–393, 1992.
  * [24] D. Sherman.  Divisible operators in von Neumann algebras.  Illinois J. Math., 54:567–600, 2010.
  * [25] A.L. Shields.  Weighted shift operators and analytic function theory, volume 13 of Math. Surveys, pages 49–128.  Amer. Math. Soc., Providence, RI, 1974.
  * [26] W. Specht.  Zur Theorie der Matrizen, II.  Jahresber. Dtsch. Math-Ver., 50:19–23, 1940.
  * [27] Cédric Villani.  Topics in optimal transportation, volume 58 of Graduate Studies in Mathematics.  American Mathematical Society, Providence, RI, first edition, 2003.
  * [28] L.N. Wasserstein.  Markov processes over denumerable products of spaces describing large systems of automata.  Probl. Inform. Transmission, 5:47–52, 1969.