# Dynamic Network Service Selection in Intelligent Reflecting Surface-Enabled
Wireless Systems: Game Theory Approaches

Nguyen Thi Thanh Van, Nguyen Cong Luong, Feng Shaohan, Huy T. Nguyen, Kun Zhu,
Thien Van Luong, and Dusit Niyato

###### Abstract

In this paper, we address dynamic network selection problems of mobile users
in an Intelligent Reflecting Surface (IRS)-enabled wireless network. In
particular, the users dynamically select different Service Providers (SPs) and
network services over time. The network services are composed of IRS resources
and transmit power resources. To formulate the SP and network service
selection, we adopt an evolutionary game in which the users are able to adapt
their network selections depending on the utilities that they achieve. For
this, the replicator dynamics is used to model the service selection
adaptation of the users. To allow the users to take their past service
experiences into account their decisions, we further adopt an enhanced version
of the evolutionary game, namely fractional evolutionary game, to study the SP
and network service selection. The fractional evolutionary game incorporates
the memory effect that captures the users’ memory on their decisions. We
theoretically prove that both the game approaches have a unique equilibrium.
Finally, we provide numerical results to demonstrate the effectiveness of our
proposed game approaches. In particular, we have reveal some important
finding, for instance, with the memory effect, the users can achieve the
utility higher than that without the memory effect.

###### Index Terms:

Intelligent reflecting surface, next-generation wireless network, evolutionary
game, fractional game, dynamic network service selection.

##  I Introduction

Intelligent Reflecting Surface (IRS) is an emerging technology for the
development of the next-generation wireless networks [1], [2]. IRS consists of
passive elements that can reflect incident signals by intelligently adjusting
their phase-shifts corresponding to wireless channels. The signals reflected
by the IRS are added constructively with non-reflected signals, i.e., the
Line-of-Sight (LoS) signals, at the user receivers to boost the received
signal power and enhance the data rate at the users. As a result, IRS has
recently been proposed to be integrated with next-generation wireless
technologies such as terahertz (THz) communications. The reason is that THz
communication is able to provide data transmission rate up to terabit per
second (Tbps), but this technology is limited in distance due to the fact that
the THz waves are vulnerability to blockage and have severe path attenuation.
For this, network service providers (SPs) deploy multiple IRSs in the THz
networks to extend the network coverage and enhance the Quality of Service
(QoS) of the mobile users. This results in a high density of the IRSs and base
stations (BSs) in the THz networks, and the mobile users will more frequently
handover among IRSs, BSs, and even SPs to achieve their desired QoS with low
cost. In this case, the dynamic network selection of the mobile users in the
THz networks becomes critical.

Although there are some works, i.e., [3], [4], [5], [6], and [7], that have
recently investigated the IRS-enabled THz networks, they do not focus on the
network selection of the mobile users. In particular, the work in [5] is
proposed to determine phase shifts of IRS to maximize the data rate. Extending
the work in [5], the work in [3] aims to jointly optimize the IRS phase
shifts, beamforming at the BS, and spectrum allocation to maximize the data
rate. Similar to [3], the work in [6] aims to maximize the overall network
throughput by jointly optimizing the phase-shits of the IRS and beamforming at
the BS.

Also, there are some works that have recently investigated the network service
selection. However, these works are not considered under the THz networks, and
they do not account fir the dynamics of the network service selection. In
particular, the authors in [8] and [9] consider an IRS-enabled network in
which the BS is owned by an SP and the IRS belongs to a different SP. The
Stackelberg game is then adopted to maximize the individual utilities of the
SPs. Accordingly, the SP of IRS as the leader offers reflection modules as
network resources and decides their prices. Note that the authors consider the
allocation of reflection modules, i.e. instead of all the reflection elements,
to the users since triggering all the reflection elements frequently results
in an increased latency of adjusting phase-shift as well as the implementation
complexity. Given the price, the SP of BS as the follower selects the best
trigger reflection modules and determines their phase shifts and the transmit
beamforming at the BS. Although the proposed game approach is demonstrated by
the simulation results to be effective, the dynamics of the network service
selections are not modeled in the work. Therefore, a more effective approach
needs to be adopted to study the dynamic network service selection.

Evolutionary game [10] as an effective tool can be adopted to study the
dynamic selection and adaptation decision of a population of agents or
players. This game has significant advantages [11],[12] compared with the
traditional games. In particular, the traditional games allow players, e.g.
the mobile users in this work, to choose the desired solution immediately,
while in the evolutionary game, the players are able to gradually adjust their
strategies until they achieve a refined equilibrium solution. Especially, the
evolutionary game is able to track and capture the strategy dynamics of the
players as well as the strategy trends and behaviors of the players over time.
Therefore, the evolutionary game has been widely used to study the dynamics of
selection behaviors of users. The authors in [13] investigate the mining pool
selection of miners in a blockchain system. Accordingly, the miners compete to
solve a crypto-puzzle to win the reward of mining new blocks. Due to the
difficulty of the crypto-puzzle, the miners are willing to select mining pools
for their secure stable profits. To study the dynamic selection of mining pool
of the miners, the evolutionary game with the replicator dynamics is adopted
for modeling the strategy evolution of the miners. It is demonstrated by both
theory analysis and simulation results that in the case of two mining pools,
the evolutionary game approach exists a unique Nash equilibrium at which no
miner has an incentive to switch its pool selection since this will undermine
some other miner’s utility.

Next-generation wireless networks are expected to deploy different wireless
access technologies, and the wireless network technology selection of the
mobile users is crucial that impacts their QoS. For this, the evolutionary
game is adopted to effectively study the dynamic network selection of the
mobile users as proposed in [14]. The evolutionary game is also adopted to
model the network service selection of secondary transmitters in a
backscatter-based cognitive network [15]. In particular, the system model
includes multiple access points serving multiple secondary transmitters. Each
access point provides three network services, namely harvest-then-transmit
(HTT), backscatter, and HTT-backscatter, to the mobile users. The secondary
transmitters receive different utilities when choosing network services from
different access points. To model the access point and service adaptation of
the secondary transmitters, a series of ordinary differential equations is
used to formulate the replicator dynamic process. Both the theory and
numerical results show that the evolutionary game approach exists a unique
equilibrium at which the secondary transmitters achieve the same utility even
if they select different access points and network services. This demonstrates
that with the evolutionary game, the users can adapt their selections
gradually to reach the equilibrium. Especially, the complexity of algorithm to
implement the evolutionary game is low. In particular, as analyzed in [15],
the complexity of strategy adaptation at each secondary transmitter is $O(1)$,
that is suitable for the dynamic strategy selections of the users.

Given the aforementioned advantages, in this paper, we adopt the evolutionary
game theory to study the dynamic service selection strategies of mobile users
in the IRS-enabled terahertz network. The considered network consists of
multiple SPs that deploy BSs along with multiple IRSs to provide network
services to multiple mobile users. In particular, the SPs offer combinations
of IRS and transmit power resources as network services that the users can
select for their data transmissions. The network is thus considered to be a
user-centric network. To satisfy different QoS requirements of the mobile
users, similar to [16], [8], we assume that the SP divides its IRSs into
reflection modules. Furthermore, the SPs have different transmit power levels
that the mobile users can select. To model the SP and service adaptation of
the users in the network, we leverage the replicator dynamic process that is
expressed as a series of ordinary differential equations. Note that with the
classical evolutionary game, the users only consider the instantaneous utility
for their decision-making, i.e., the SP and service adaptation. This is not
natural and practical due to the fact that the user is typically aware of its
past network service experience when making the network selection. In other
words, the awareness of the users’ memory needs to be accounted. To address
the limitation of the classical evolutionary game, we further adopt the
fractional evolutionary game as a memory-aware economic process. The
fractional evolutionary game enables the users to incorporate the
instantaneous and past experiences of the users for their decisions. Both
theoretical analysis and simulation results show the effectiveness of the
proposed game approaches.

The main contributions of the paper include the followings:

  * •

We consider the IRS-enabled terahertz network in which multiple SPs deploy
IRSs to serve the mobile users. The SP offers IRS and transmit power resources
as network resources to the mobile users. Different combinations of the
network resources constitute different network services provided by the SPs.
The network is a user-centric network in which the users can select and adapt
the network services provided by the SPs over time to achieve their desired
utility. The IRS-enabled terahertz network introduces new transmission scheme
that makes the utility function of the users more complicated, and a new
solution is required to model the network service adaptation of the users. To
model the network service adaptation, we adopt the evolutionary game.

  * •

We consider the scenario in which the users use the delayed information for
their decisions. Such a delay can cause instability in the decision making
process. In this scenario, the delayed replicator dynamics is adopted to model
the SP and network service adaptation. We analyze the equilibrium region of
the delayed replicator dynamics and show in the simulation results that the
evolutionary game approach still reach an equilibrium with a small delay.

  * •

To capture the users’ memory on their decision-making, we incorporate the
users’ memory effect to reformulate the SP and network service selection
problem into a fractional evolutionary game. We then compare the network
selection strategies of the users between the classical game and fractional
game.

  * •

We theoretically prove that the fractional evolutionary game processes a
unique equilibrium. Then, the simulation results with the direction field of
the replicator dynamics are provided to verify the stability of the
equilibrium.

  * •

We provide performance evaluation to demonstrate the consistency with the
analytical results and to validate both the proposed game approaches. The
performance comparison between the two game approaches are also discussed and
analyzed.

The rest of the paper is organized as follows. In Section II, we present the
IRS-enabled terahertz system and utility functions of the mobile users. In
Section III, we formulate the dynamic SP and network service selection problem
as the classical evolutionary game and analytically derive the stability
region of the delayed replicator dynamics. In Section IV, we reformulate the
dynamic SP and network service selection problem as a fractional evolutionary
game, followed by the proofs of the existence, uniqueness, the stability of
the equilibrium of the game. The simulation results and discussions are
presented in Section V, and the conclusions are given in Section VI.

##  II System Model

This section presents the system model, channel models, and utility functions
of the users in the network. Typical notations used in this paper are
summarized in Table I.

###  II-A Network Model

TABLE I: List of frequency symbols used in this paper. Notation | Description
---|---
$M,N$ | Number of SPs, number of users
$I_{m},J_{m,j}$ |  Number of IRS of SP $m$, power level $j$ offered by SP $m$
$K_{m}$ |  Number of reflection elements of each IRS of SP $m$
$Q_{m}$ |  Number of modules of each IRS of SP $m$
$\theta_{m,l,k,e}$ |  Phase shift of element $e$ of subset $k$ in IRS $l$ of SP $m$
$N_{m,l,k,j}$ |  Number of users selecting power level $J_{m,j}$ and subset $k$ in IRS $l$ of SP $m$
$p_{m,l,k,j,i}$ |  Probability that user $i$ selects power level $J_{m,j}$ and subset $k$ in IRS $l$ of SP $m$
$\mathbf{h}_{m,i}$ |  Channel from BS $m$ to user $i$
$\mathbf{h}_{m,l,k,i}^{\rm{IU}}$ |  Channel from subset $k$ in IRS $l$ of SP $m$ to user $i$
$\mathbf{G}_{m,l,k}$ |  Channel from BS $m$ to subset $k$ in IRS $l$ of SP $m$
$\mathbf{w}_{m,j,i}$ |  Beamforming vector associated with user $i$ that selects power level $J_{m,j}$ of SP $m$
$f,\mu,\beta$ | Carrier frequency, learning rate, order of the Caputo fractional derivative
$\gamma_{m}^{I},\gamma_{m}^{P}$ | Price per IRS element, price per power unit
Figure 1: A system model with multiple SPs with multiple IRSs serving multiple
users.

The system model is an IRS-enabled terahertz MIMO system as shown in Fig. 1.
The system model consists of a set $\mathcal{M}$ of $M$ SPs and a set
$\mathcal{N}$ of $N$ single-antenna users. Without loss of generality, each SP
$m\in\mathcal{M}$ deploys a BS, i.e., BS $m$, that is equipped with $L_{m}$
antennas. Denote $B_{m}$ as the bandwidth allocated to SP $m$, i.e., BS $m$.
To provide flexible services to the users, BS $m$ has a set $\mathcal{P}_{m}$
of $P_{m}$ power levels, denoted by $\\{J_{m,1},\ldots,J_{m,P_{m}}\\}$, that
the users can select for their transmissions. Note that the assumption of the
discrete power levels is reasonable since in real networks, transmit power
control algorithms choose steps of power increment/decrement
[wu2001distributed]. Also, the number of power levels can be increased to
match with the real system implementation without causing much more
complexity. We also assume that $J_{m,1}<J_{m,2}<\dots<J_{m,P_{m}}$, where
$J_{m,P_{m}}$ is the maximum power of BS $m$. To improve the QoS for the
users, SP $m$ deploys a set $\mathcal{I}_{m}$ of $I_{m}$ IRSs. Let
$l\in\mathbb{Z}$ denote the index of IRS, and $1\leq
l\leq\max_{m}\\{I_{m}\\}$. We assume that IRSs belonging to the same SP have
the same size, and IRS $l$ of SP $m$ has $K_{m,l}=K_{m}$ reflection elements.
IRS $l$ of SP $m$ is divided into $Q_{m,l}=Q_{m}$ modules that are controlled
by parallel switches. Each module in IRS $l$ consists of $E_{m,l}=E_{m}$
elements. Note that during a time slot, one BS-IRS pair of the corresponding
SP can serve multiple users, but the user is associated with one BS-IRS pair.
Moreover, the user can select one or multiple modules, i.e., a subset of
modules, of the selected IRS. A network service is defined as a combination of
a power level and a subset of modules. In general, the data throughput
achieved by the user, say user $i$, depends on 1) the power level that the
user selects, 2) the bandwidth allocated to the BS that the user selects, 3)
the location of the selected IRS, 4) the number of modules of the selected
IRS, and 5) the interference caused by other users selecting the same BS with
user $i$. Note that the data throughput does not depend on indexes of the
modules of the selected IRS. As such, each IRS $l$ of SP $m$ has a set
$\mathcal{Q}_{m,l}$ including $Q_{m,l}$ of potential subsets of modules that
the user can select, and subset $k,1\leq k\leq Q_{m,l},$ of the IRS has $k$
modules. Denote $\boldsymbol{\Theta}_{m,l,k}$ as the phase-shift matrix
corresponding to the subset that the user selects, i.e., subset $k$ of IRS $l$
of SP $m$. Then, $\boldsymbol{\Theta}_{m,l,k}$ is a diagonal matrix in which
its main diagonal consists of phase-shifts of $kE_{m,l}$ reflection elements
of IRS $l$ of SP $m$. In particular, we have
$\boldsymbol{\Theta}_{m,l,k}=\text{diag}(\theta_{m,l,k,1},\ldots,\theta_{m,l,k,kE_{m,l}})$,
where $\theta_{m,l,k,e}$ is the phase-shift of reflection element $e$ of
subset $k$ in IRS $l$ of SP $m$,
$\theta_{m,l,k,e}=e^{j\varphi_{m,l,k,e}},\varphi_{m,l,k,e}\in[0,2\pi),e=\\{1,\ldots,kE_{m,l}\\}$.
With the assistance of subset $k$ of IRS $l$ of SP $m$, the signal received at
each user $i$ is the sum of 1) the received signal via the direct link, 2) the
received signal via the IRS-assisted link, and 3) the intra-interference
caused by other users that select the same BS with user $i$. Let
$\mathcal{N}_{m}$ denote the set of $N_{m}$ users selecting SP $m$, i.e., and
also BS $m$. Then, the received signal at user $i$ when selecting subset $k$
of IRS $l$ and power level $J_{m,j}$ offered by SP $m$ is determined as
follows:

$\displaystyle y_{i}=$
$\displaystyle\big{(}\mathbf{h}^{\text{H}}_{m,i}+(\mathbf{h}^{\rm{IU}}_{m,l,k,i})^{\text{H}}\boldsymbol{\Theta}^{\text{H}}_{m,l,k}\mathbf{G}_{m,l,k}\big{)}\mathbf{w}_{m,j,i}s_{i}+$
$\displaystyle+\sum_{n\in\mathcal{N}_{m},n\neq
i}\big{(}\mathbf{h}^{\text{H}}_{m,i}+(\mathbf{h}^{\rm{IU}}_{m,l,k,i})^{\text{H}}\boldsymbol{\Theta}^{\text{H}}_{m,l,k}\mathbf{G}_{m,l,k}\big{)}\mathbf{w}_{m,j,n}s_{n}+\omega_{i},$
(1)

where $s_{i}$ is the data symbol intended to user $i$,
$\mathbf{w}_{m,j,i}\in\mathbb{C}^{L_{m}\times 1}$ is the beamforming vector
associated with $s_{i}$ containing power level $J_{m,j}$ that the user
selects, $\mathbf{h}_{m,i}\in\mathbb{C}^{L_{m}\times 1}$ is the channel from
BS $m$ to user $i$,
$\mathbf{h}^{\rm{IU}}_{m,l,k,i}\in\mathbb{C}^{K_{m,l}\times 1}$ is the channel
from subset $k$ of IRS $l$ of SP $m$ to user $i$,
$\mathbf{G}_{m,l,k}\in\mathbb{C}^{K_{m,l}\times L_{m}}$ is the vector of
channels from BS $m$ to subset $k$ of IRS $l$, and $\omega_{i}$ is the
Gaussian noise at user $i$, $\omega_{i}\sim\mathcal{CN}(0,\sigma_{0}^{2})$,
where $\sigma_{0}^{2}$ is the variance. Similar to [18] and [19], we assume
that each BS $m$ has a perfect knowledge of channel state information (CSI) of
the channels. Note that we aim to model network selection strategies of the
users which is not influenced by this assumption. In fact, CSI estimation
algorithms such as [20] can be used to obtain the full CSI at all the BSs with
low training overhead. Since IRSs are typically deployed in static
environments, we can also assume that the quasi-static flat-fading model or
even static flat-fading model is applied for all channels [18]. The signal-to-
interference-plus-noise ratio (SINR) over bandwidth $B_{m}$ of user $i$ is
defined as follows:

$\displaystyle\eta_{m,l,k,j,i}=\frac{\left|(\mathbf{h}^{\text{H}}_{m,i}+(\mathbf{h}^{\rm{IU}}_{m,l,k,i})^{\text{H}}\boldsymbol{\Theta}^{\text{H}}_{m,l,k}\mathbf{G}_{m,l,k})\mathbf{w}_{m,j,i}\right|^{2}}{\left|\sum_{n\in\mathcal{N}_{m},n\neq
i}\big{(}\mathbf{h}^{\text{H}}_{m,i}+(\mathbf{h}^{\rm{IU}}_{m,l,k,i})^{\text{H}}\boldsymbol{\Theta}^{\text{H}}_{m,l,k}\mathbf{G}_{m,l,k}\big{)}\mathbf{w}_{m,j,n}\right|^{2}+B_{m}\sigma_{0}^{2}}.$
(2)

To remove the intra-interference among the users, the BSs can use time-
division multiple access for their users. In this case, $\eta_{m,l,k,i}$ can
be expressed by

$\eta_{m,l,k,j,i}=\frac{\left|(\mathbf{h}^{\text{H}}_{m,i}+(\mathbf{h}^{\rm{IU}}_{m,l,k,i})^{\text{H}}\boldsymbol{\Theta}^{\text{H}}_{m,l,k}\mathbf{G}_{m,l,k})\mathbf{w}_{m,j,i}\right|^{2}}{B_{m}\sigma_{0}^{2}}.$
(3)

###  II-B Channel Model

In this section, we describe the channel models for the IRS-enabled THz
network. The channel models for each user include the channel between the BS
that the user selects and the user, and the cascaded channel of the IRS-aided
link. Here, the cascaded channel of the IRS-aided link includes (1) the
channel between the BS and the subset of IRS modules that the user selects and
(2) the channel between the subset of IRS modules and the user. To model the
channels in the THz network, we adopt the Saleh-Valenzuela channel model [21].
Without loss of generality, we model the channels between user $i$ when it
selects BS $m$, subset $k$ of IRS $l$ of SP $m$. In particular, we determine
models of channels $\mathbf{h}_{m,i}$, $\mathbf{h}^{\rm{IU}}_{m,l,k,i}$, and
$\mathbf{G}_{m,l,k}$ that are given in (2). For an ease of presentation, we
remove the indices $m,l,k,i$ from the channels, and thus $\mathbf{h}_{m,i}$,
$\mathbf{h}^{\rm{IU}}_{m,l,k,i}$, and $\mathbf{G}_{m,l,k}$ can be expressed by
$\mathbf{h}$, $\mathbf{h}^{\rm{IU}}$, and $\mathbf{G}$, respectively. Also, we
assume that BS $m$ has $L$ antennas, subset $k$ of IRS $l$ that the user
selects has $K$ reflection elements.

####  II-B1 BS-user channel

The channel between the BS and the user is expressed by

$\mathbf{h}=\kappa^{(0)}\mathbf{a}(\Phi^{(0)})+\sum_{l=1}^{\mathcal{L}}\kappa^{(l)}\mathbf{a}(\Phi^{(l)}),$
(4)

where $\kappa^{(0)}\mathbf{a}(\Phi^{(0)})$ is the LoS element in which
$\kappa^{(0)}$ is the gain and $\mathbf{a}(\Phi^{(0)})$ is the spatial
direction, and $\kappa^{(l)}\mathbf{a}(\Phi^{(l)}),1\leq l\leq\mathcal{L}$, is
one of $\mathcal{L}$ non-LoS (NLoS) elements. $\mathbf{a}(\Phi^{(l)})$ is the
$L\times 1$ array steering vector corresponding to th-$l$ element that is
defined as follows:

$\mathbf{a}(\Phi^{(l)})=\frac{1}{\sqrt{L}}\big{[}e^{-j2\pi\Phi^{(l)}(-\frac{N-1}{2})},\ldots,e^{-j2\pi\Phi^{(l)}(\frac{N-1}{2})}\big{]},$
(5)

where $\Phi^{(l)}$ is the spatial direction of the signal corresponding to
component $l$, that is defined as $\Phi^{(l)}=\frac{d}{\lambda}\sin\xi^{(l)}$,
where $\xi^{(l)}\in[-\pi/2,\pi/2]$ is the angle-of-departure (AoD) of path $l$
corresponding to the BS and user, $\lambda$ is the signal wavelength, and $d$
is the distance between adjacent antennas of the BS or the distance between
adjacent IRS elements of the IRS that is typically defined as $d=\lambda/2$.

####  II-B2 BS-IRS-user channel

The channel between the BS and the subset of modules of IRS that the user
selects can be modeled as

$\mathbf{G}=\sqrt{\frac{LK}{\mathcal{L}_{\text{BS-I}}}}\sum_{l=1}^{\mathcal{L}_{\text{BS-I}}}\kappa_{\text{BS-I}}^{(l)}\mathbf{a}_{\text{I}}(\Phi^{(l)}_{\text{AOA}})\mathbf{a}_{\text{BS}}^{\text{H}}(\Phi^{(l)}_{\text{AOD}}),$
(6)

where $\mathcal{L}_{\text{BS-I}}$ denotes the scattering paths between the BS
and the subset of IRS that the user selects, $\kappa^{(l)}_{\text{BS-I}}$ is
the complex gain of path $l$, and $\Phi^{(l)}_{\text{AOD}}$ and
$\Phi^{(l)}_{\text{AOA}}$ are the spatial directions of path $l$ corresponding
to the BS and the subset of IRS, respectively. We consider the BS’s antennas
and the IRS’s reflection elements as uniform linear arrays (ULAs), and thus we
can determine
$\mathbf{a}_{\text{BS}}(\Phi^{(l)}_{\text{AoD}})\in\mathcal{C}^{L}$ and
$\mathbf{a}_{\text{I}}(\Phi^{(l)}_{\text{AOA}})\in\mathcal{C}^{K}$ as follows:
$\mathbf{a}_{\text{BS}}(\Phi^{(l)}_{\text{AOD}})=\frac{1}{\sqrt{L}}\big{[}e^{-j2\pi\Phi^{(l)}_{\text{BS}}(-\frac{L-1}{2})},\ldots,e^{-j2\pi\Phi^{(l)}_{\text{BS}}(\frac{L-1}{2})}\big{]}$,
and
$\mathbf{a}_{\text{I}}(\Phi^{(l)}_{\text{AOA}})=\frac{1}{\sqrt{K}}\big{[}e^{-j2\pi\Phi^{(l)}_{\text{I}}(-\frac{K-1}{2})},\ldots,e^{-j2\pi\Phi^{(l)}_{\text{I}}(\frac{K-1}{2})}\big{]}$.
Here, $\Phi^{(l)}_{\text{I}}=\frac{d}{\lambda}\sin\xi^{(l)}_{\text{I}}$ and
$\Phi^{(l)}_{\text{BS}}=\frac{d}{\lambda}\sin\xi^{(l)}_{\text{BS}}$, where
$\xi^{(l)}_{\text{BS}}\in[-\pi/2,\pi/2]$ and
$\xi^{(l)}_{\text{I}}\in[-\pi/2,\pi/2]$ are the AoD and the angle-of-arrival
(AoA) of path $l$ corresponding to the BS and the subset of IRS, respectively.

Similarly, we can determine the channel between the subset of IRS and the user
as follows:

$\mathbf{h}^{\text{IU}}=\sum_{l=0}^{\mathcal{L}_{\text{I-U}}}\kappa^{(l)}_{\text{I-U}}\mathbf{a}_{\text{I-U}}(\Phi^{(l)}_{\text{I-U}}),$
(7)

where $\mathcal{L}_{\text{I-U}}$ denotes the scattering paths between the
subset of IRS and the user, $\kappa^{(l)}_{\text{I-U}}$ is the complex gain of
path $l$, and
$\mathbf{a}_{\text{I}}(\Phi^{(l)}_{\text{I-U}})=\frac{1}{\sqrt{K}}\big{[}e^{-j2\pi\Phi^{(l)}_{\text{I-U}}(-\frac{K-1}{2})},\ldots,e^{-j2\pi\Phi^{(l)}_{\text{I-U}}(\frac{K-1}{2})}\big{]}$,
where $\Phi^{(l)}_{\text{I-U}}=\frac{d}{\lambda}\sin\xi^{(l)}_{\text{I-U}}$
with $\xi^{(l)}_{\text{I-U}}\in[-\pi/2,\pi/2]$ being the AoD of path $l$
corresponding to the subset of IRS.

####  II-B3 Path loss

In THz communication systems, the non-LoS elements is proved to be much weaker
than the LoS element, i.e., lower than $20$ dB [22]. Therefore, in the IRS-
enabled THz network, we consider the LoS elements of the involved channels.
Without loss of generality, we calculate the path loss of the LoS element
between the BS and the user. This can be applied to calculating the path loss
of the LoS elements between the BS and the subset of IRS as well as that
between the subset of IRS and the user. The path loss of the LoS element,
denoted by $\chi_{\text{LoS}}$, is a function of spreading loss and molecular
absorption loss, denoted by $\chi_{\text{abs}}$. Then, the path loss is
determined as [22]

$\chi_{\text{LoS}}(f)=\chi_{\text{spr}}(f)\chi_{\text{abs}}(f)e^{-j2\pi
f\tau_{\text{LoS}}},$ (8)

where $f$ is the carrier frequency, and $\tau_{\text{LoS}}=r/c$ is the LoS
propagation time, $r$ is the distance between the BS and the user, and $c$ is
the speed of light. $\chi_{\text{spr}}$ is the spreading loss that is
determined by

$\chi_{\text{spr}}(f)=\frac{c}{4\pi fr}.$ (9)

While, $\chi_{\text{abs}}(f)$ is the molecular absorption loss that is
determined as follows:

$\chi_{\text{abs}}(f)=e^{-0.5\zeta(f)r},$ (10)

where $\zeta(f)$ is the medium absorption coefficient that depends on carrier
frequency and the composition of the transmission medium at a molecular level.
For example, given $f=2$ THz, the molecular absorption coefficient
$\zeta=2.3\times 10^{-5}$ m-1 for oxygen (O2) [23].

###  II-C Utility Functions

This section presents the utility functions of the users when they select
different SP and network services. There are totally $N$ users, $M$ BSs,
$P_{m}$ power levels, and $\sum_{m=1}^{M}I_{m}$ IRSs in the network. Users
selecting the same SP, the same subset of the IRS, and the same power level
are grouped into a group. Thus, there are totally $G$ groups in the network,
where $G=\sum_{m=1}^{M}\sum_{l=1}^{I_{m}}P_{m}I_{m}Q_{m,l}$. Without loss of
generality, we can assume that group $g,1\leq g\leq G$, consists of users that
select SP $m$, subset $k$ of IRS $l$, i.e., the corresponding phase-shift
matrix $\boldsymbol{\Theta}_{m,l,k}$, and power level $J_{m,j}$. For this, we
can denote $g$ as the combination of indexes $(m,l,k,j)$ for the expression
simplification. Let $\mathcal{N}_{m,l,k,j}$, i.e., $\mathcal{N}_{g}$, be a set
of $N_{m,l,k,j}$, i.e., $N_{g}$, of users in group $g$. We have
$\sum_{m=1}^{M}\sum_{l=1}^{I_{m}}\sum_{k=1}^{Q_{m,l}}\sum_{j=1}^{P_{m}}N_{m,l,k,j}=N$,
and user $i\in\mathcal{N}_{g}$ selects SP $m$, IRS $l$,
$\boldsymbol{\Theta}_{m,l,k}$, and $J_{m,j}$ at a probability of
$p_{m,l,k,j,i}=N_{m,l,k,j}/N$. As the expected number of the users selecting
SP $m$, IRS $l$, $\boldsymbol{\Theta}_{m,l,k}$, and $J_{m,j}$ is
$p_{m,l,k,j,i}N$ and BS $m$ adopts the time-division multiple access, each
user in the group will access the channel with a probability of
$\frac{1}{p_{m,l,k,j,i}N}$ for every time slot. Therefore, the expected data
rate that the user in the group can achieve is

$\overline{R}_{m,l,k,j,i}=\frac{B_{m}}{p_{m,l,k,j,i}N}\log_{2}\big{(}1+\eta_{m,l,k,j,i}\big{)},$
(11)

where $\eta_{m,l,k,j,i}$ is given in (3). Let $v_{m,l,k,j,i}$ denote the value
of unit data to user $i$ in group $g$ when selecting SP $m$, subset $k$ of IRS
$l$, and power level $J_{m,j}$. Denote $\gamma_{m}^{I}$ as the price per
element in IRSs of SP $m$ and $\gamma_{m}^{P}$ as the price per unit power.
The prices, i.e., $\gamma_{m}^{I}$ and $\gamma_{m}^{P}$, are set by SP $m$
that are constant. Since users in each group share the same resources, they
should share the resource cost. Then, the utility of the user is given by

$u_{m,l,k,j,i}=v_{m,l,k,j,i}\overline{R}_{m,l,k,j,i}-\frac{\gamma_{m}^{I}\lVert\boldsymbol{\Theta}_{m,l,k}\rVert_{0}-\gamma_{m}^{P}J_{m,j}}{p_{m,l,k,j,i}N},$
(12)

where the $l_{0}$-norm is used to count the number of non-zero elements of a
diagonal matrix that here refers to the number of active reflection elements
of IRS $l$ of SP $m$ that the user selects.

##  III Evolutionary Game Formulation

In this section, we leverage the evolutionary game to model the dynamic SP and
network service selection of the users. We prove that the game can achieve the
evolutionary equilibrium at which no user has an incentive to change their
network service strategy.

###  III-A Game Formulation

Each user in the network is able to adapt their network selection over time,
and it can achieve different utility at different time points. Thus, by taking
the SP and network service selection strategies, the expected or average
utility of user $i$ at time $t$ is

$\overline{u}_{i}=\sum_{m=1}^{M}\sum_{l=1}^{I_{m}}\sum_{k=1}^{Q_{m,l}}\sum_{j=1}^{P_{m}}p_{m,l,k,j,i}u_{m,l,k,j,i}.$
(13)

To model the SP and service adaptation of the users, we leverage the
replicator dynamic process that is expressed as a series of ordinary
differential equations as follows [15], [17]:

$\displaystyle\dot{p}_{m,l,k,j,i}(t)=$ $\displaystyle\mu
p_{m,l,k,j,i}(t)\left[u_{m,l,k,j,i}(t)-\overline{u}_{i}(t)\right],$
$\displaystyle
m\in\mathcal{M},l\in\mathcal{I}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g},\forall
t,$ (14)

where $\dot{p}_{m,l,k,j,i}(t)$ represents the first derivative of
$p_{m,l,k,j,i}$ with respect to $t$, and
$p_{m,l,k,j,i}(t_{0})=p^{0}_{m,l,k,j,i}$ is the initial strategy of the user
in group $g$ at $t_{0}$. The factor $\mu$ is the learning rate of the users
that evaluates the strategy adaptation rate. The replicator dynamics process
given in (14) represents the population strategy evolution of the users in the
network. That is, the population of users evolves over time, and the game
converges to the evolutionary equilibrium. This means that the users select an
SP and its service with higher utility over time, and the evolutionary
equilibrium can be defined as the set of stable fixed points of the replicator
dynamics.

To show the existence of the evolutionary equilibrium of the game as defined
in (14), we use the following theorem. First, we let $f_{g,i}(t,p_{g,i})=\mu
p_{g,i}(t)\left[u_{g,i}(t)-\overline{u}_{i}(t)\right]$, where $g$ denotes as
the combination of indexes $(m,l,k,j)$, and then we rewrite equation (14) as
follows:

$\dot{p}_{g,i}(t)=f_{g,i}(t,p_{g,i}),\text{ }p_{g,i}(t_{0})=p^{0}_{g,i},\text{
}g\in\\{1,\ldots,G\\},i\in\mathcal{N}_{g}.$ (15)

######  Theorem 1.

[24] Suppose that functions $f_{g,i}(t,p_{g,i})$ $\frac{\partial
f_{g,i}}{\partial p_{g,i}}(t,p_{g,i})$ are continuous in some open rectangle
$\left\\{(t,p_{g,i}):0\leq t\leq\tau,0<p_{g,i}\leq 1\right\\}$ that contains
point $(t_{0},p^{0}_{g,i})$. Then, the problem in (15) has a unique solution
in the interval of $I=\left[t_{0}-\varsigma,t_{0}+\varsigma\right]$, where
$\varsigma>0$. Moreover, the Picard iteration defined by

$p^{<z+1>}_{g,i}(t)=p^{0}_{g,i}+\int_{t_{0}}^{t}f_{g,i}(t,p^{<z>}_{g,i}(t))\;\mathrm{d}t$
(16)

produces a sequence of functions ${p^{<z>}_{g,i}(t)}$ that converges to the
solution uniformly on $I$.

###### Proof.

Theorem 1 means that problem in (15) converges to a unique solution, i.e., the
equilibrium of the game defined in (15), given that function
$f_{g,i}(t,p_{g,i})$ and its derivative are continuous with respect to time
$t$. Therefore, we first show that $f_{g,i}(t,p_{g,i})$ and $\frac{\partial
f_{g,i}}{\partial p_{g,i}}(t,p_{g,i})$ are continuous functions in the
rectange $\left\\{(t,p_{g,i}):0\leq t\leq\tau,0\leq p_{g,i}\leq 1\right\\}$.
Indeed, it is clear that function $p_{g,i}(t)=\frac{N_{g,i}(t)}{N}$ is
continuous at every $t_{0}\in[0,\tau]$. Moreover, due to the static flat-
fading channel model, variables $\mathbf{h}_{m,i}(t)$,
$\mathbf{h}_{g,i}^{\rm{IU}}(t)$, and $\mathbf{G}_{m,l,k}(t)$ are constant and
thereby continuous at every $t_{0}\in[0,\tau]$. Correspondingly,
$\eta_{g,i}(t)$ is continuous at every $t_{0}\in[0,\tau]$, and functions
$\overline{R}_{g,i}(t),u_{g,i}(t)$, and $\overline{u}_{i}(t)$ are also
continuous at every $t_{0}\in[0,\tau]$ if $p_{g,i}(t_{0})\neq 0$. Since
$f_{g,i}(t,p_{g,i})=\mu p_{g,i}(t)[u_{g,i}(t)-\overline{u}_{i}(t)]$ and
$\frac{\partial f_{g,i}}{\partial
p_{g,i}}(t,p_{g,i})=\mu\left[u_{g,i}(t)-\overline{u}_{i}(t)\right]$, then
$f_{g,i}(t,p_{g,i})$ and $\frac{\partial f_{g,i}}{\partial
p_{g,i}}(t,p_{g,i})$ are continuous functions in the open rectangle
$\left\\{(t,p_{g,i}):0\leq t\leq\tau,0<p_{g,i}\leq 1\right\\}$.

Given that the continuity of $f_{g,i}(t,p_{g,i})$ and its derivative, there
are many ways to prove Theorem 1. One of them is leveraging the Banach Fixed
Point Theorem (BFPT) [25] to approximate a solution, i.e., a fixed point, to
(15) by constructing a sequence of functions that converges to a unique
solution. The proof of Theorem 1 using the BFPT is well explained and
presented in [24]. The unique solution refers to the game equilibrium at which
1) all the users achieve the same utility and 2) no user has an incentive to
change its network service selection. ∎

The overall process of the network service selection is summarized as follows.
Initially, each user randomly selects an SP and a service of the SP. Given the
user selection, the SP determines the optimal phase-shift and beamforming for
its associated users according to Algorithm 1. The user computes its utility
according to (12) and transmits the utility information to the SP. The user
compares its utility and the expected utility determined by (13) and can
change its network service selection to achieve a higher utility value. After
all the users achieve the same utility by choosing any strategies, then no
user has an incentive to change its network service selection and the game
converges to the evolutionary equilibrium.

Algorithm 1 Optimizing phase-shift and beamforming for user $i$ when selecting
SP $m$, subset $k$ of IRS $l$, and power level $J_{m,j}$ [26] (in this
algorithm, function $\rm{unt}(\mathbf{a})$ is defined as
$\rm{unt}(\mathbf{a})=[a_{1}/|a_{1}|,\ldots,a_{n}/|a_{n}|]$)

1:Output: $\mathbf{w}_{m,j}$ and $\mathbf{\Theta}_{m,l,k}$ for each user;

2:Initialize: $t=0$, $\epsilon_{1}=0$, $\mathbf{\Theta}_{m,l,k}^{0}$;

3:Calculate
$\mathbf{R}_{11}=\text{diag}(({\mathbf{h}^{\rm{IU}}_{m,l,k,i})^{\rm{H}})\mathbf{G}_{m,l,k}\mathbf{G}_{m,l,k}^{\rm{H}}\text{diag}(\mathbf{h}^{\rm{IU}}_{m,l,k,i}})$;

4:Construct $\mathbf{R}$ =
$\begin{bmatrix}\mathbf{R}_{11}&\text{diag}(({\mathbf{h}^{\rm{IU}}_{m,l,k,i}})^{\rm{H}})\mathbf{G}_{m,l,k}\mathbf{h}_{m,l,k,i}\\\
\mathbf{h}_{m,i}^{\rm{H}}\mathbf{G}_{m,l,k}^{\rm{H}}\text{diag}(\mathbf{h}^{\rm{IU}}_{m,l,k,i})&0\end{bmatrix}$

5:repeat

6: $\mathbf{v}^{(t)}=[\rm{diag}(\mathbf{\Theta}_{m,l,k}),t]^{\top}$;

7: Calculate $\mathbf{v}^{(t+1)}$ according by
$\mathbf{v}^{(t+1)}$=unt$(\mathbf{Rv}^{(t)})$;

8: $t\leftarrow t+1$;

9:until
$\parallel{\mathbf{Rv}^{(t+1)}\parallel_{1}-\parallel\mathbf{Rv}^{(t)}\parallel_{1}}\leq\epsilon$;

10:Take first $kE_{m,l}$ elements of ${(\mathbf{v}^{t+1})}^{*}$ as the main
diagonal of $\mathbf{\Theta}_{m,l,k}$;

11:Compute
$\mathbf{w}_{m,j,i}=\sqrt{J_{m,j}}\frac{\mathbf{G}_{m,l,k}^{\rm{H}}\rm{diag}(\mathbf{h}^{\rm{IU}}_{m,l,k,i})\mathbf{\Theta}^{\rm{H}}_{m,l,k}+\mathbf{h}_{m,l,k,i}}{\parallel\mathbf{G}_{m,l,k}^{\rm{H}}\rm{diag}(\mathbf{h}^{\rm{IU}}_{m,l,k,i})\mathbf{\Theta}^{\rm{H}}_{m,l,k}+\mathbf{h}_{m,i}\parallel}$.

The computational complexity of the algorithm is mainly caused from 1) the
phase-shift and beamforming optimization implemented at the SPs (BSs) side and
2) the utility computation implemented at the user side. When the users select
a network service of SP $m$, the SP optimizes the phase-shift matrix and
beamforming for each user using Algorithm 1 that requires
$5K_{m}^{2}L_{m}+K_{m}^{2}(5L_{m}+2)+K_{m}(L_{m}+2)+2L_{m}+1$ multiplications
and additions. When the size increases to infinite, the complexity of the
algorithm for each user is $\mathcal{O}(n^{3})$. Each user calculates its own
utility based on the prices and network services that it selects. Thus, the
computational complexity of the user does not increase with the total numbers
of users and the SPs. Thus, the complexity of the algorithm implemented at
each user is $\mathcal{O}(1)$. This implies that the game approach is
computationally efficient and highly scalable.

###  III-B Delay in Replicator Dynamics

In the game model discussed in the previous section, to make the decision on
SP and service selections, the users need information about the average
utility, i.e., $\overline{u}_{i}$, and the proportion of users choosing
different strategies, i.e., $p_{m,l,k,j,i}$, from the BSs. However, the up-to-
date information may not be available at the users due to the communication
latency. Thus, at time instance $t$, the users may need to use the information
at time $t-\delta$, i.e., delay for $\delta$ time units, to make the SP and
service selections. Thus, the delayed replicator dynamic process is expressed
as

$\displaystyle\dot{p}_{m,l,k,j,i}(t)=$ $\displaystyle\mu
p_{m,l,k,j,i}(t-\delta)\left[u_{m,l,k,j,i}(t-\delta)-\overline{u}_{i}(t-\delta)\right],$
$\displaystyle
m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g},\forall
t.$ (17)

Note that as delay $\delta$ is large, the decisions of the users based on the
outdated information tend to be inaccurate. In this case, the SP and service
selections may not converge. How to determine $\delta^{*}$ such that the
selections converge is challenging. As an example, consider a simple scenario
with $M=2$, and SP $m$ offers one service including subset
$\boldsymbol{\Theta}_{m}$ and power level $J_{m}$:

######  Theorem 2.

The evolutionary game can converge to a stable equilibrium if the value of
$\delta$ is satisfied the following condition:

$\delta^{*}<\frac{\pi}{2\mu\sum_{m\in\mathcal{M}}\frac{v_{m}B_{m}\log_{2}(1+\eta_{m})-\gamma_{m}^{I}||\boldsymbol{\Theta}_{m}||_{0}-\gamma_{m}^{P}J_{m}}{N}}.$
(18)

###### Proof.

The delayed replicator dynamics in (III-B) can be rewritten as

$\mathbf{\dot{p}}(t)=\mathbf{A}\mathbf{p}(t-\delta)+\mathbf{c},$ (19)

where
$\mathbf{\dot{p}}(t)={[{\dot{p}_{1}}(t),\ldots,{\dot{p}_{M}}(t)]}^{\top}$,
$\mathbf{p}(t)={[{p}_{1}(t-\delta),\ldots,{p}_{M}(t-\delta)]}^{\top}$,
$\mathbf{c}={[\frac{\mu a_{1}}{N},\ldots,\frac{\mu a_{M}}{N}]}^{\top}$ with
$a_{m}=\frac{v_{m}B_{m}\log_{2}(1+\eta_{m})-\gamma_{m}^{I}||\boldsymbol{\Theta}_{m}||_{0}-\gamma_{m}^{P}J_{m}}{N}$,
and $\mathbf{A}=-\kappa\mathbf{I}$. Here, $I$ is the identity matrix of size
$M$, and $\kappa$ is defined as

$\displaystyle\kappa=\mu\sum_{m\in\mathcal{M}}\frac{v_{m}B_{m}\log_{2}(1+\eta_{m})-\gamma_{m}^{I}||\boldsymbol{\Theta}_{m}||_{0}-\gamma_{m}^{P}J_{m}}{N}.$

Otherwise, the evolutionary game with the delayed replicator dynamics can
converge to a stable equilibrium if the real parts of all the roots are
negative [27]. This is equivalently the condition $\kappa\delta<\pi$, and thus
we have

$\displaystyle\delta^{*}<\frac{\pi}{2\kappa}=\frac{\pi}{2\mu\sum_{m\in\mathcal{M}}\frac{v_{m}B_{m}\log_{2}(1+\eta_{m})-\gamma_{m}^{I}||\boldsymbol{\Theta}_{m}||_{0}-\gamma_{m}^{P}J_{m}}{N}}.$

∎

Theorem 18 means that the evolutionary game is guaranteed to converge to the
equilibrium as the users use information at $t<\delta^{*}$ for their
decisions.

##  IV Fractional Evolutionary Game Formulation

In this section, we discuss the use of the fractional evolutionary game to
model the SP and network service selection of the users with memory effect in
the IRS-enabled terahertz system. In particular, we first present the concept
of memory-aware economic process. Then, we present how to cast the
evolutionary game that describes the SP and network service selection into a
fractional evolutionary game by using the memory-aware economic processes.
Finally, we analyze the equilibrium of the game.

###  IV-A Memory-aware Economic Process

With the classical evolutionary game as presented in Section III, each user,
say user $i$, decides on the SP and network service selection according to its
instantaneously achievable utility functions, i.e., functions
$u_{m,l,k,j,i}(t)$ and $\overline{u}_{i}(t)$ at time instant $t$. In reality,
the users take into account their memory, i.e., of service experience, on
their strategy decisions. Specifically, the selection decision of the users at
time $t$ is based not only on the information about the state of the process
at time $t$, but also on the information about the process states at previous
time instants $\tau\in[0,t]$. This is considered to be a memory-aware economic
process [28], [29]. To describe the memory-aware economic process, we consider
a typical economic model with two variables, namely exogenous variable and
endogenous variable. The exogenous variable and endogenous variable are the
input and output of the economic model, respectively. This means that the
endogenous variable depends on the exogenous variable, and they are similar to
the independent and dependent variables, respectively. Denote $X(t)$ as the
exogenous variable and $Y(t)$ as the endogenous variable variable, in which
the exogenous variable changes depends on the changes of the endogenous
variable. Then, the economic process is typically expressed by
$Y(t)=F^{t}_{0}(X(\tau))+Y_{0}$, where $\tau\in[0,t]$, $Y_{0}$ is the initial
state of the output of the process, and $F^{t}_{0}$ is an operator. To enable
the memory awareness of the economic process, the operator is defined as
$F^{t}_{0}(X(\tau)):=\int_{0}^{t}M_{\beta}(t-\tau)X(\tau)\rm{d}t$, where
$M_{\beta}(t-\tau)$ is the weight function that represents how the input
$X(\tau)$ at time $\tau$ impacts on the output $Y(t)$ at time $t$. In general,
function $M_{\beta}(t-\tau)$ changes with respect to $\tau$ so as to capture
the dynamic characteristic of the memory. Furthermore, by taking the time
derivative of $Y(t)$, we have
$\frac{\rm{d}}{\rm{d}t}Y(t)=M_{\beta}(t)X(0)+\int_{0}^{t}M_{\beta}(t-\tau)[\frac{\rm{d}}{\rm{d}t}Y(t)]X(\tau)\rm{d}\tau$
that depends on both $X(t)$ and $X(\tau)$ with $\tau\in[0,t)$. The formulation
of $M_{\beta}(t-\tau)$ is
$M_{\beta}(t-\tau)=\frac{1}{\Gamma(\beta)(t-\tau)^{1-\beta}}$, where
$\Gamma(\cdot)$ is the gamma function that is defined by
$\Gamma(\beta)=\int_{0}^{+\infty}x^{\beta-1}e^{-x}\rm{d}x$.

Since $Y(t)$ depends on both $X(t)$ and $X(\tau)$, the economic process is
namely memory-aware economic process that can be expressed in the fractional
equation by taking the derivation of $Y(t)$ at the order of $\beta$ through
the left-sided Caputo fractional derivative as follows:

${}^{C}_{0}D^{\beta}_{t}Y(t)=X(t),$ (20)

where $Y(0)=Y_{0}$ is the initial state, and ${}^{C}_{0}D^{\beta}_{t}Y(t)$ is
the left-sided Caputo fractional derivative [30] of $Y(t)$ at the order of
$\beta$ that is given by:

${}^{C}_{0}D^{\beta}_{t}Y(t)=\frac{1}{\Gamma(\lceil\beta\rceil-\beta)}\int^{t}_{0}\frac{Y^{(\lceil\beta\rceil)}(\tau)}{(t-\tau)^{\beta+1-\lceil\beta\rceil}}\rm{d}\tau,$
(21)

where $\lceil\cdot\rceil$ is the ceiling function.

The memory-aware economic process given in (21) has two key properties. First,
the past experiences of the user at different time instances have different
impacts on its decision-making so as to capture dynamically the memory of the
user. Second, the user is affected by the experience within the memory rather
than that at the current time, and consequently the memory-aware users can
make network selection decisions differently from the memory-unaware users.
Given the properties, we incorporate the memory awareness of the economic
process when modeling the SP and network service selection of the users. The
memory-aware economic process can be modeled as the fractional evolutionary
game that is presented in the next section.

###  IV-B Fractional Game Formulation

For convenience, we rewrite the replicator dynamic process of the users in the
classical evolutionary game as expressed in (14) as follows:

$\displaystyle\dot{p}_{m,l,k,j,i}(t)=$ $\displaystyle\mu
p_{m,l,k,j,i}(t)\left[u_{m,l,k,j,i}(t)-\overline{u}_{i}(t)\right],$
$\displaystyle
m\in\mathcal{M},l\in\mathcal{I}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g},\forall
t.$ (22)

Then, given the utility functions and the average utility of the users defined
in (12) and (13), respectively, and by incorporating the memory characteristic
of the users, we can formulate the fractional evolutionary game as follows:

$\displaystyle{}^{C}_{0}D^{\beta}_{t}p_{m,l,k,j,i}(t)=$ $\displaystyle\mu
p_{m,l,k,j,i}(t)[u_{m,l,k,j,i}(t)-\overline{u}_{i}(t)],$ $\displaystyle
m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g},\forall
t,$ (23)

where $\beta\in(0,2),\beta\neq 1$ is the order of the Caputo fractional
derivative, and it is called memory effect coefficient. The physical meaning
of the left-sided Caputo is further explained and discussed in Section V. The
equilibrium analysis of the fractional evolutionary game is presented in the
next section.

###  IV-C Equilibrium Analysis

In this section, we theoretically discuss the existence and the uniqueness of
the equilibrium, and the unique and stable equilibrium is admitted as the
solution of the fractional evolutionary game defined in (IV-B). The specific
steps are as follows. First, we transfer the fractional game defined in (IV-B)
into an equivalent problem, i.e, (25), and the equivalence between which is
verified in Theorem 3. Then, to prove the existence and uniqueness of
equilibrium of the game defined in (IV-B), we provide the proof of the
uniqueness of the solution to the equivalent problem defined in (25).

For the ease of presentation, we let
$\mathbf{P}(t)=[p_{m,l,k,j,i}\left(t\right)]_{m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}}$
and $\mathbf{E}(\mathbf{P}(t))=\big{[}\mu
p_{m,l,k,j,i}(t)[u_{m,l,k,j,i}(t)-\overline{u}_{i}(t)]\big{]}_{m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}}$,
and reorganize the fractional evolutionary game defined in (IV-B) as follows:

${}^{C}_{0}D^{\beta}_{t}\mathbf{P}(t)=\mathbf{E}(\mathbf{P}(t)),$ (24)

with the initial strategy
$\mathbf{P}(0)=\mathbf{P}^{0}=[p^{0}_{m,l,k,j,i}]_{m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}}$
and the time horizon $\mathcal{T}=[0,T]$.

######  Theorem 3.

If all the elements of vector $\mathbf{E}$ in (24), i.e., $e_{n}$ (element $n$
of vector $\mathbf{E}$) for all
$n\in\left\\{\left.\left(m,l,k,j,i\right)\right|m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}\right\\}$,
can satisfy the following two conditions:

  * •

$e_{n}\in{\cal{C}}^{2}$ with ${\cal{C}}^{2}$ being the set of the twice
differentiable functions;

  * •

$\frac{\partial}{\partial p_{m,l,k,j,i}}e_{n}$ exists and is bounded for all
$m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}$.

Then, (24) can be equivalently transformed into the following problem

$\mathbf{P}(t)=\mathbf{P}^{0}+{}_{0}\rm{I}^{\beta}_{t}\mathbf{E}(\mathbf{P}(t)),\forall
t\in\mathcal{T}.$ (25)

The second condition means that for all
$n\in\left\\{\left.\left(m,l,k,j,i\right)\right|m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}\right\\}$,
there exists $L\in\mathbb{R}^{+}$ such that
$|e_{n}({\hat{\mathbf{P}}}(t))-e_{n}({\tilde{\mathbf{P}}}(t))|<L||({\hat{\mathbf{P}}}(t))-{\tilde{\mathbf{P}}}(t)||_{\mathcal{L}_{1}}$,
which implies the satisfaction of the Lipschitz condition.

###### Proof.

According to (25), the $\lceil\beta\rceil$-th derivative of $\mathbf{P}(t)$
with respect to $t$ is as follows:

$\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)=\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\big{[}{}_{0}\rm{I}^{\beta}_{t}\mathbf{E}(\mathbf{P}(t))\big{]}={}^{RL}_{0}D^{\beta-\lceil\beta\rceil}_{t}{\mathbf{E}}\left(\mathbf{P}(t)\right),$
(26)

with ${}^{RL}_{0}D^{\beta}_{t}{\mathbf{E}}\left(\mathbf{P}(t)\right)$ being
defined as the left-sided Riemann-Liouville fractional derivative with respect
to $t$, and the following derivation is satisfied

$\displaystyle{}^{RL}_{0}D^{\lceil\beta\rceil-\beta}_{t}{\mathbf{E}}\left(\mathbf{P}(t)\right)$
$\displaystyle=\frac{\rm{d}}{\rm{d}t}\Big{[}\frac{1}{\Gamma(1-\lceil\beta\rceil+\beta)}\int^{t}_{0}\frac{\mathbf{E}(\mathbf{P}(\tau))}{(t-\tau)^{\lceil\beta\rceil-\beta}}\rm{d}\tau\Big{]}$
(27)
$\displaystyle=\frac{1}{\Gamma(1-\lceil\beta\rceil+\beta)}\frac{\rm{d}}{\rm{d}t}\int^{t}_{0}\theta^{\beta-\lceil\beta\rceil}\mathbf{E}(\mathbf{P}(t-\theta))\rm{d}\theta$
$\displaystyle=\frac{1}{\Gamma(1-\lceil\beta\rceil+\beta)}\Big{[}t^{\beta-\lceil\beta\rceil}\mathbf{E}(\mathbf{P}^{0})+\int^{t}_{0}\theta^{\beta-\lceil\beta\rceil}\frac{\rm{d}}{\rm{d}t}\mathbf{E}(\mathbf{P}(t-\theta))\rm{d}\theta\Big{]}$
$\displaystyle=\frac{1}{\Gamma(1-\lceil\beta\rceil+\beta)}\Big{[}t^{\beta-\lceil\beta\rceil}\mathbf{E}(\mathbf{P}^{0})+\int^{t}_{0}(t-\tau)^{\beta-\lceil\beta\rceil}\frac{\rm{d}}{\rm{d}\tau}\mathbf{E}(\mathbf{P}(\tau))\rm{d}\tau\Big{]}$
$\displaystyle=\frac{t^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})+{}_{0}\rm{I}^{\beta}_{t}\Big{[}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{E}(\mathbf{P}(t))\Big{]}.$

Let $\sigma\in(0,t)$, the $\mathcal{L}^{1}$ norm of
$\frac{t^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})$
exists an upper bound that is derived as follows:

$\|\frac{t^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})\|_{\mathcal{L}^{1}}\leq\|\frac{\sigma^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})\|_{\mathcal{L}^{1}}.$
(28)

Using the condition in Theorem 3 that $\frac{\delta}{\delta
p_{m,l,k,j,i}}e_{n}$ exists and is bounded for all
$m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}$,
and (26) and (27) as well as (28), we have

$\displaystyle\Big{\|}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}_{\mathcal{T}}$
$\displaystyle<\|\frac{\sigma^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})\|_{\mathcal{L}^{1}}+\Big{\|}{}_{0}\rm{I}_{t}^{\beta}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{E}(\mathbf{P}(t))\Big{\|}_{\mathcal{T}}$
(29)
$\displaystyle<\|\frac{\sigma^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})\|_{\mathcal{L}^{1}}+\Big{\|}{}_{0}\rm{I}_{t}^{\beta}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}_{\mathcal{T}}AL,$

where $\|z\|_{\mathcal{T}}=\int_{\mathcal{T}}\exp\left(-\mu
t\right)\|z\|_{\mathcal{L}^{1}}\rm{d}t$ and $A$ is the cardinality of
$\left\\{\left.\left(m,l,k,j,i\right)\right|m\in\mathcal{M},l\in\mathcal{L}_{m},k\in\mathcal{Q}_{m,l},\right.$
$\left.j\in\mathcal{P}_{m},i\in\mathcal{N}_{g}\right\\}$. For the last term in
(29), we have

$\displaystyle\Big{\|}{}_{0}\rm{I}_{t}^{\beta}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}_{\mathcal{T}}$
$\displaystyle=\int^{T}_{0}\exp\left(-\mu
t\right)\Big{\|}{}_{0}\rm{I}_{t}^{\beta}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}\rm{d}t\leq\int_{0}^{T}\exp\left(-\mu
t\right)\int^{t}_{0}\frac{1}{\Gamma(\beta)}\frac{\|\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}s^{\lceil\beta\rceil}}\mathbf{P}(s)\|}{(t-s)^{1-\beta}}\rm{d}s\rm{d}t$
(30) $\displaystyle=\int_{0}^{T}\frac{\exp\left(-\mu
s\right)}{\Gamma(\beta)}\|\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}s^{\lceil\beta\rceil}}\mathbf{P}(s)\|\int^{T}_{s}\exp\left(-\mu(t-s)\right)(t-s)^{\beta-1}\rm{d}t\rm{d}s$
$\displaystyle=\int_{0}^{T}\frac{\exp\left(-\mu
s\right)}{\Gamma(\beta)}\|\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}s^{\lceil\beta\rceil}}\mathbf{P}(s)\|\int^{\mu(T-s)}_{0}\exp\left(-\Psi\right)(\frac{\Psi}{\mu})^{\beta-1}\rm{d}(\frac{\Psi}{\mu})\rm{d}s$
$\displaystyle<\frac{1}{\mu^{\beta}\Gamma(\beta)}\int_{0}^{T}\exp\left(-\mu
s\right)\|\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}s^{\lceil\beta\rceil}}\mathbf{P}(s)\|\rm{d}s\int^{+\infty}_{0}\exp\left(-\Psi\right)\Psi^{\beta-1}\rm{d}\Psi$
$\displaystyle=\frac{1}{\mu^{\beta}}\|\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}s^{\lceil\beta\rceil}}\mathbf{P}(t)\|_{\cal{T}}.$

Then, we substitute (30) into (29) as follows

$\displaystyle\Big{\|}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}_{\mathcal{T}}<\|\frac{\sigma^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})\|_{\mathcal{L}^{1}}+\frac{AL}{\mu^{\beta}}\|\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}s^{\lceil\beta\rceil}}\mathbf{P}(t)\|_{\cal{T}}$
(31) $\displaystyle\Leftrightarrow$
$\displaystyle\Big{\|}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}_{\mathcal{T}}<\frac{1}{1-\frac{AL}{\mu^{\beta}}}\|\frac{\sigma^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})\|_{\mathcal{L}^{1}},$

which implies that if $\mu$ is sufficiently large such that
$\frac{AL}{\mu^{\beta}}<1$,
$\Big{\|}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)\Big{\|}_{\mathcal{T}}$
has an upper bound. In this case, the fractional derivative of $\mathbf{P}(t)$
with the order of $\beta\in\left(0,1\right)\cup\left(1,2\right)$ exists and
can be obtained in the following:

${}^{C}_{0}D^{\beta}_{t}\mathbf{P}(t)={}_{0}\rm{I}^{{\lceil\beta\rceil}-\beta}_{t}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{P}(t)={}_{0}\rm{I}^{{\lceil\beta\rceil}-\beta}_{t}\left\\{\frac{t^{\beta-\lceil\beta\rceil}}{\Gamma(1-\lceil\beta\rceil+\beta)}\mathbf{E}(\mathbf{P}^{0})+{}_{0}\rm{I}^{\beta}_{t}\Big{[}\frac{\rm{d}^{\lceil\beta\rceil}}{\rm{d}t^{\lceil\beta\rceil}}\mathbf{E}(\mathbf{P}(t))\Big{]}\right\\}=\mathbf{E}(\mathbf{P}(t)),$
(32)

∎

and hence the equivalence between (24) and (25) has been verified, which
completes the proof.

######  Theorem 4.

Given the conditions in Theorem 3, the uniqueness of the solution to the
problem defined in (25) can be guaranteed.

###### Proof.

First, by defining an operator
$\Lambda:\mathcal{D}_{\mathbf{P}}\mapsto\mathcal{D}_{\mathbf{P}}$, where
$\mathcal{D}_{\mathbf{P}}$ is the feasible domain of ${\mathbf{P}}$, there
exists an inequality expression as follows:

$\|\Lambda{\hat{\mathbf{P}}}(t)-\Lambda{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}}<\frac{AL}{\mu^{\beta}}\|{\hat{\mathbf{P}}}(t)-{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}},$
(33)

and the specific derivation of which has been shown as follows

$\displaystyle\|\Lambda{\hat{\mathbf{P}}}(t)-\Lambda{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}}=\int^{T}_{0}\exp\left(-\mu
t\right)\|{}_{0}\mathbf{I}^{\beta}_{t}\mathbf{E}({\hat{\mathbf{P}}}(t))-{}_{0}\mathbf{I}^{\beta}_{t}\mathbf{E}({\tilde{\mathbf{P}}}(t))\|\rm{dt}$
(34) $\displaystyle<AL\Big{[}\int^{T}_{0}\exp\left(-\mu
t\right)\|{}_{0}\mathbf{I}^{\beta}_{t}{\hat{\mathbf{P}}}(t)-{}_{0}\mathbf{I}^{\beta}_{t}{\tilde{\mathbf{P}}}(t)\|\rm{dt}\Big{]}$
$\displaystyle\leq AL\Big{[}\int^{T}_{0}\exp\left(-\mu
t\right)\int^{t}_{0}\frac{\|{\hat{\mathbf{P}}}(s)-{\tilde{\mathbf{P}}}(s)\|}{\Gamma(\beta)(t-s)^{1-\beta}}{\rm{d}}s{\rm{d}}t\Big{]}$
$\displaystyle=\frac{AL}{\Gamma(\beta)}\Big{[}\int^{T}_{0}\int^{T}_{s}\exp\left(-\mu
t\right)\frac{\|{\hat{\mathbf{P}}}(s)-{\tilde{\mathbf{P}}}(s)\|}{(t-s)^{1-\beta}}\rm{dtds}\Big{]}$
$\displaystyle=\frac{AL}{\Gamma(\beta)}\Big{[}\int^{T}_{0}\exp\left(-\mu
s\right)\|{\hat{\mathbf{P}}}(s)-{\tilde{\mathbf{P}}}(s)\|\int^{T}_{s}\frac{\exp\left(-\mu(t-s)\right)}{(t-s)^{1-\beta}}{\rm{d}}s{\rm{d}}t\Big{]}$
$\displaystyle=\frac{AL}{\mu^{\beta}\Gamma(\beta)}\Big{[}\int^{T}_{0}\exp\left(-\mu
s\right)\|{\hat{\mathbf{P}}}(s)-{\tilde{\mathbf{P}}}(s)\|\int^{\mu(T-s)}_{0}\exp\left(-\psi\right)\psi^{\beta-1}{\rm{d}}\psi{\rm{d}}s\Big{]}$
$\displaystyle<\frac{AL}{\mu^{\beta}\Gamma(\beta)}\|{\hat{\mathbf{P}}}(t)-{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}}\int^{+\infty}_{0}\exp\left(-\sigma\right)\sigma^{\beta-1}\rm{d}\sigma=\frac{AL}{\mu^{\beta}}\|{\hat{\mathbf{P}}}(t)-{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}}.$

Based on (34), we can conclude that
$\|\Lambda{\hat{\mathbf{P}}}(t)-\Lambda{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}}<\|{\hat{\mathbf{P}}}(t)-{\tilde{\mathbf{P}}}(t)\|_{\mathcal{T}}$
if $\mu^{\beta}\geq AL$. In this case, the operator $\Lambda$ satisfies the
fixed point theorem, which indicates the uniqueness of the solution to (25).
By following this, there exists a unique solution to the fractional
evolutionary game defined in (24), which completes this proof.

∎

##  V Performance Evaluation

In this section, we present and discuss simulation results obtained by the
proposed evolutionary game approaches. To evaluate the game approaches, we
consider three cases, i.e., $\beta=1$ corresponding to the classical
evolutionary game, and $\beta=1.1$ and $0.8$ corresponding to the fractional
evolutionary games. For the comparison purpose, we consider a network that
consists of two SPs, namely SP 1 and SP 2, and $100$ users. Each SP deploys a
BS that is equipped with $4$ antennas. SP 1 deploys 2 IRSs, namely IRS $1_{1}$
and $1_{2}$, and SP 2 deploys 1 IRS, namely IRS $2$. SP 1 divides each IRS
into two modules and offers $1$ power level, i.e., $J_{1,1}=30$ dBm. SP 2 does
not divide its IRS and offers $2$ power levels, i.e., $J_{2,1}=25$ dBm and
$J_{2,2}=35$ dBm. As such, SP 1 offers 4 services, and SP 2 offers 2 services
that the users can select. Correspondingly, $100$ users are divided into $6$
groups. Note that the SPs can offer more services and our proposed game
approaches are scalable since the complexity of the algorithm implemented at
each user is $\mathcal{O}(1)$ as analyzed in Section III-A. The locations of
the BSs, IRSs, and users are shown in Fig. 2. The simulation parameters are
provided in Table II. In particular for the involved channels, the non-LoS
elements is proved to be much weaker than the LoS element, i.e., lower than
$20$ dB [22], and thus similar to [32], we mainly consider the channel with
only LoS element, i.e.,
$\mathcal{L}=\mathcal{L}_{\text{BS-I}}=\mathcal{L}_{\text{I-U}}=1$.

Figure 2: Coordinates of BSs, IRSs, and users. TABLE II: Simulation parameters. Parameters | Value | Parameters | Value | Parameters | Value | BSs, IRSs and user | Coordinate $(m)$
---|---|---|---|---|---|---|---
$M$ | $2$ | $B_{1}=B_{2}$ |  $2$ MHz |  $v_{1}$, $v_{2}$ | $10^{-6}$ |  BS $1$ | [100 20]
$L_{1}=L_{2}$ | $4$ | $f$ |  $0.3$ THz | $\gamma^{I}_{1}$ | $2.5\times 10^{-4}$ | BS 2 | [0 20]
$P_{1}$ | $1$ | $c$ |  $3\times 10^{8}$ m/s | $\gamma^{I}_{2}$ | $10^{-3}$ |  IRS $1_{1}$ | [80 20]
$P_{2}$ | $2$ | $N$ | $100$ | $\gamma^{P}_{1}$ | $0.05$ |  IRS $1_{2}$ | [40 20]
$J_{1,1}$ |  $30$ dBm | $Q_{1,1}=Q_{1,2}$ | $2$ | $\gamma^{P}_{2}$ | $0.03$ |  IRS $2_{1}$ | [30 20]
$J_{2,1}$ |  $25$ dBm | $Q_{2,1}$ | $1$ | $\gamma^{P}_{2}$ | $0.03$ | User group | [50 0]
$J_{2,2}$ |  $35$ dBm |  $E_{1,1}$, $E_{1,2}$, $E_{2,1}$ | $8$ | $(\mu,\zeta)$ | $(e^{-2},0)$ |  |

(a)

(b)

(c)

(d)

Figure 3: Proportion of users selecting different SPs and services.

First, we discuss strategies that the users choose different SPs and services
over evolutionary time. Figures 3(a), (b), and (c) illustrate the results for
the evolutionary games with $\beta=1,1.1$ and $0.8$, respectively. As seen,
the strategies of the users choosing different SPs and services eventually
converge to an equilibrium point over time. Moreover, during the initial
phase, the users’ strategies in the evolutionary games with $\beta=1$ and
$0.8$ fluctuate in a range smaller than those in the evolutionary game with
$\beta=1.1$. The results indicate that the adaptations of the users’
strategies in the game with $\beta=1.1$ are faster than those in the games
with $\beta=1$ and $0.8$. These results are further verified in Fig. 3(d) in
which the strategy adaptation frequency of the users in the game with
$\beta=1.1$ is higher than those in the games with $\beta=1$ and $0.8$. Note
that as the users’ strategies have a larger fluctuation, the convergence speed
can be slower. Thus, as we can observe from Fig. 3(d), the game with
$\beta=1.1$ converges to the equilibrium more slowly than those in the games
with $\beta=1$ and $0.8$.

(a)

(b)

Figure 4: The impact of $\beta$ on the evolutionary games.

Now, we discuss how the memory effect coefficient, i.e., $\beta$, impacts on
the evolutionary games. As shown in Fig. 4(a), with $\beta<1$ and $\beta>1$,
as $\beta$ increases, the convergence time of the user’s strategy is shorter.
This implies that as $\beta$ increases, the replicator dynamics converges
faster. Moreover, as shown in Fig. 4(b), as $\beta$ increases, the rate that
the corresponding games converge to the vicinity of the equilibrium is faster.
This means that the adaptation rate of the user’s strategy increases with the
increase of $\beta$. This results is also consistent with descriptions in
Figs. 3(a), (b), (c), and (d).

Next, it is important to show the utility that the users can achieve when
different games are used. For this, we vary the values of $\beta$, i.e., the
memory effect, and we evaluate the total utility the the users achieve. As
shown in Fig. 4(b), the utility for the users with $\beta=0.8<1$ is worse than
that for the users with $\beta=1$. Meanwhile, with $\beta=1.1>1$, the users
achieve higher utility values when the users have no memory effect, i.e.,
$\beta=1$. Since the users achieve higher utility values with $\beta=1.1$, we
can say that the memory effect with $\beta>1$ is a positive effect. This
further implies that to achieve a higher utility value, the users should
incorporate both the past and instantaneously achievable experiences for their
network selection. In addition, as seen from Fig. 4(b), the total utility of
users with IRSs in the games, i.e., $\beta=0.8$, $1$, and $1.1$ is much higher
than that of users without IRSs. This result demonstrates that deploying IRSs
increases the throughput of the users.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 5: Direction field of the replicator dynamics when $t=200$.

To show that the users’ strategies in the proposed game approaches can be
stabilized at the equilibrium, we present the direction field of the
replicator dynamics. As illustrated in Fig. 5, the strategies of the users
eventually reach the equilibrium strategy after a certain time, i.e., $t=200$,
that is represented by the black circles. We can take the results shown in Fig
5(a) as an example. In the figure, we show the replicator dynamics of
selection strategy of SP 1, service 1 and SP 2, service 1 $(p_{1,1}$ and
$p_{2,1}$). Assuming that the strategies the users select services provided by
SP 2 (i.e., $p_{2,1}$ and $p_{2,2}$) achieve the equilibrium, and
$p_{1,1}=1-p_{1,2}-p_{1,3}-p_{1,4}-p_{2,1}-p_{2,2}$. As seen, the users are
able to adapt their strategies by following the directions of the arrows.
Furthermore, any initial strategy eventually reach the equilibrium that
verifies the stability of our proposed game approaches. That is similar to
Figs. 5(b), (c), (d), (e), and (f).

(a)

(b)

(c)

Figure 6: Utility of user groups.

Now, we discuss how the utilities of the users obtained at the equilibrium.
Figures 6(a), (b), and (c) show the utilities that the users achieve by
selecting different SPs and network services over evolutionary time. As seen,
the utilities of the users vary until the equilibrium is reached. At the
equilibrium, the users have the same utility even if they select different SPs
and services. The reason is that the evolutionary equilibrium is reached only
when the utilities of the users choosing any SP and any service are equal to
their expected utility.

(a)

(b)

(c)

Figure 7: Proportion of users selecting different service as the size of IRS
provided by SP 2 varies.

Next, we investigate the impact of sizes of IRSs on the proportions of users
selecting different SPs and services. Figures 7 (a), (b), and (c) show the
results obtained from the evolutionary game with $\beta=1,1.1$ and $0.8$,
respectively. In particular, we vary the size $K_{2}$ of IRS 2 of SP 2. As
shown in Fig. 7(a), as $K_{2}$ increases, the proportions of users selecting
services provided by SP 2 increase since the throughput and utility obtained
by the users selecting services provided by SP 2 increase. However, as the
size of IRS 2 is large, the increasing rate tends to be slower. This is
because of that the users pay a very high resource cost if they select the
services provided by SP 2, and thus they tend to select the services provided
by SP 1. Figures 7(b) and (c) have the same pattern as Fig. 8(a), and the
results can be explained in the same way.

(a)

(b)

(c)

Figure 8: Impact of the number of users on the total utility and the time to
reach the equilibrium.

Note that the time to reach the equilibrium can be different depending on the
learning rate $\mu$ and the number of users. Figures 8(a), (b), and (c) show
the results obtained from the fractional evolutionary games with $\beta=1,1.1$
and $0.8$, respectively. As seen, the evolutionary equilibrium is reached
faster as $\mu$ is higher since the frequency of the strategy adaptation of
the users is higher. Moreover, the games need more time to converge to the
equilibrium as the number of users $N$ increases. Note that as $N$ increases,
the total utility of the users decreases. This can be explained based on (11),
more users share the fixed amount of bandwidth that results in reducing the
throughput of the users. Moreover, from Fig. 8(a), the total utility of users
with IRSs is much higher than that of users without IRSs. The reason is that
deploying IRSs increases the throughput of the users.

(a)

(b)

(c)

Figure 9: Proportions of users selecting different SPs vs. distance.

Next, we discuss how the mobility of the users impacts the selection
strategies of the users. In particular, we evaluate the proportions of users
selecting different SPs as the distance between the users and IRS 1 provided
by SP 1 varies. The results for the evolutionary games with $\beta=1,1.1$ and
$0.8$ are respectively shown in Figs. 9(a), (b), and (c). As observed from
Figs. 9(a), (b), and (c), the proportion of users selecting SP 2 increases as
the distance between the users and IRS 1 of SP 1 increases. This is because of
that the throughput obtained by the users selecting services of SP 1
decreases. Thus, the users are willing to select services of SP 2. In
addition, we consider the case as SP 1 deploys one IRS. As observed from Figs.
9(a), (b), and (c), as SP 1 deploys 1 IRS, the proportions that the users
select SP 1 is lower than those that SP 1 deploys 2 IRSs. Especially,
decreases more slowly than those that SP 1 deploys 2 IRSs. Interestingly, as
SP 1 deploys 1 IRS, the proportions that the users select SP 1 decrease faster
than those that SP 1 deploys 2 IRSs. This implies that by deploying more IRSs,
the SP can further improve the QoS of the users and can prevent the users to
select the network service of other SPs.

(a)

(b)

(c)

Figure 10: Proportion of the users selecting SP1 and Service 1 with different
delays.

Next, we discuss how the proportions of users selecting different SPs and
services depend on the information delay $\delta$. For the evaluation purpose,
we consider the proportion of the users selecting Service 1 provided by SP 1
as shown in Figs. 10(a) (b), and (c). First, we discuss the results obtained
by the classical evolutionary game shown in Figs. 10(a). This figure shows the
proportion of the users choosing Service 1 of SP 1 when the users use
information for their decisions at $\delta=0,40$ and $80$. As seen, when
$\delta>0$, there is a fluctuating dynamics of strategy adaptation. In
particular, as $\delta=40$, the game can sill converge to the equilibrium that
is the same as the case when $\delta=0$. However, as $\delta=80$, the game
cannot reach the equilibrium. These results mean that when the users use
information with a small delay for their decisions, the game is still
guaranteed for the convergence. In addition, we find in the figure that as
$\delta=40$, the proportions of the users have more fluctuations and the game
needs more time to reach the equilibrium, compared with the game where
$\delta=0$. This implies that the convergence speed is slower as the users use
information with larger delay. The results obtained by the fractional games
with $\beta=0.8$ and $1.1$ are the same as that obtained by the classical
evolutionary game. However, it seems to be that the service selection cannot
reach the evolutionary equilibrium even when $\delta$ is small. For example,
for fractional evolutionary game with $\beta=1.1$, as $\delta\leq 2$, the
service selection cannot reach the evolutionary equilibrium. This is also a
shortcoming of the fractional game in which the very outdated information may
not be used for the user decisions.

(a)

(b)

(c)

Figure 11: Proportion of the users selecting different SPs as the size of IRS
and the price of IRS element provided by SP2 varies.

We finally discuss how the sizes of IRSs and the price of IRS 2 elements
impact on the SP and service selection of the users. Figure 11(a) shows that,
for a given price, the proportion of users selecting SP 2 increases as the
number of elements of IRS 2 increases. This is because of that the throughput
obtained by the users selecting services of SP 2 increases. Thus, the users
are willing to select services of SP 2. This figure also shows the proportion
of users selecting SP 1 increases as the price of IRS 2 elements increases.
This is because of that the utility of users choose SP 2 will fall down as the
price of IRS 2 elements increases. Therefore, the users are willing choose
services of SP 1. The same explanations can be applied to the results obtained
as the values of $\beta$ are $0.8$ and $1.1$ (Figs. 11(b) and (c)).

##  VI Conclusions

We have proposed dynamic game frameworks for modeling the dynamic network
selection of mobile users in the IRS-enabled terahertz network. First, we have
adopted the classical evolutionary game in which the SP and service adaptation
of the users is modeled as replicator dynamics. We have further considered the
scenario in which the users use delayed information for their decision-making.
In this scenario, we have analyzed the stability region of the delayed
replicator dynamics. Furthermore, we have adopted the fractional evolutionary
game that incorporates the memory effect to model the SP and service
adaptation of the users. The proof of the existence and uniqueness of the game
equilibrium has been provided. We have finally provided simulation results
obtained by the proposed evolutionary game approaches. In addition, we have
further discussed the selection behaviors of the users and compared the
performance obtained by the proposed evolutionary games. For future work, we
will study the SP and network service selections in a heterogeneous network
that includes different types of relaying, i.e., IRS and active relay devices,
and different communication technologies, i.e., terahertz and millimeter wave
communications.

## References

  * [1] Q. Wu and R. Zhang, “Towards smart and reconfigurable environment: Intelligent reflecting surface aided wireless network,” _IEEE Communications Magazine_ , vol. 58, no. 1, pp. 106–112, Jan. 2019.
  * [2] M. Di Renzo _et al._ , “Smart radio environments empowered by reconfigurable ai meta-surfaces: An idea whose time has come,” _EURASIP Journal on Wireless Communications and Networking_ , vol. 2019, no. 1, pp. 1–20, May 2019.
  * [3] W. Chen, X. Ma, Z. Li, and N. Kuang, “Sum-rate maximization for intelligent reflecting surface based terahertz communication systems,” in _2019 IEEE/CIC International Conference on Communications Workshops in China (ICCC Workshops)_. IEEE, 2019, pp. 153–157.
  * [4] X. Ma, Z. Chen, W. Chen, Y. Chi, Z. Li, C. Han, and Q. Wen, “Intelligent reflecting surface enhanced indoor terahertz communication systems,” _Nano Communication Networks_ , vol. 24, p. 100284, 2020.
  * [5] X. Ma, Z. Chen, W. Chen, Z. Li, Y. Chi, C. Han, and S. Li, “Joint channel estimation and data rate maximization for intelligent reflecting surface assisted terahertz mimo communication systems,” _IEEE Access_ , 2020.
  * [6] W. Hao, G. Sun, M. Zeng, Z. Zhu, Z. Chu, O. A. Dobre, and P. Xiao, “Robust design for intelligent reflecting surface assisted mimo-ofdma terahertz communications,” _arXiv preprint arXiv:2009.05893_ , 2020.
  * [7] X. Ma, Z. Chen, W. Chen, Y. Chi, L. Yan, C. Han, and S. Li, “Joint hardware design and capacity analysis for intelligent reflecting surface enabled terahertz mimo communications,” _arXiv preprint arXiv:2012.06993_ , 2020\.
  * [8] Y. Gao, C. Yong, Z. Xiong, D. Niyato, Y. Xiao, and J. Zhao, “A stackelberg game approach to resource allocation for intelligent reflecting surface aided communications,” _arXiv preprint arXiv:2003.06640_ , 2020.
  * [9] ——, “Resource allocation for intelligent reflecting surface aided cooperative communications,” _arXiv preprint arXiv:2012.10229_ , 2020.
  * [10] J. Hofbauer and K. Sigmund, “Evolutionary game dynamics,” _Bulletin of the American mathematical society_ , vol. 40, no. 4, pp. 479–519, 2003.
  * [11] Z. Han, D. Niyato, W. Saad, T. Başar, and A. Hjørungnes, _Game theory in wireless and communication networks: theory, models, and applications_. Cambridge university press, 2012.
  * [12] N. Quijano, C. Ocampo-Martinez, J. Barreiro-Gomez, G. Obando, A. Pantoja, and E. Mojica-Nava, “The role of population games and evolutionary dynamics in distributed control systems: The advantages of evolutionary game theory,” _IEEE Control Systems Magazine_ , vol. 37, no. 1, pp. 70–97, 2017.
  * [13] X. Liu, W. Wang, D. Niyato, N. Zhao, and P. Wang, “Evolutionary game for mining pool selection in blockchain networks,” _IEEE Wireless Communications Letters_ , vol. 7, no. 5, pp. 760–763, 2018.
  * [14] D. Niyato and E. Hossain, “Dynamics of network selection in heterogeneous wireless networks: An evolutionary game approach,” _IEEE transactions on vehicular technology_ , vol. 58, no. 4, 2008.
  * [15] X. Gao, S. Feng, D. Niyato, P. Wang, K. Yang, and Y.-C. Liang, “Dynamic access point and service selection in backscatter-assisted rf-powered cognitive networks,” _IEEE Internet of Things Journal_ , vol. 6, no. 5, pp. 8270–8283, 2019.
  * [16] Y. Gao, C. Yong, Z. Xiong, J. Zhao, Y. Xiao, and D. Niyato, “Reflection resource management for intelligent reflecting surface aided wireless networks,” _arXiv preprint arXiv:2002.00331_ , 2020.
  * [17] S. Feng, D. Niyato, X. Lu, P. Wang, and D. I. Kim, “Dynamic game and pricing for data sponsored 5g systems with memory effect,” _IEEE Journal on Selected Areas in Communications_ , vol. 38, no. 4, pp. 750–765, 2020.
  * [18] Q. Wu and R. Zhang, “Intelligent reflecting surface enhanced wireless network via joint active and passive beamforming,” _IEEE Transactions on Wireless Communications_ , vol. 18, no. 11, pp. 5394–5409, 2019.
  * [19] G. Zhou, C. Pan, H. Ren, K. Wang, and A. Nallanathan, “Intelligent reflecting surface aided multigroup multicast miso communication systems,” _IEEE Transactions on Signal Processing_ , to appear.
  * [20] C. Liaskos, A. Tsioliaridou, A. Pitilakis, G. Pirialakos, O. Tsilipakos, A. Tasolamprou, N. Kantartzis, S. Ioannidis, M. Kafesaki, A. Pitsillides _et al._ , “Joint compressed sensing and manipulation of wireless emissions with intelligent surfaces,” in _International Conference on Distributed Computing in Sensor Systems (DCOSS)_ , 2019, pp. 318–325.
  * [21] C. Lin and G. Y. Li, “Indoor terahertz communications: How many antenna arrays are needed?” _IEEE Transactions on Wireless Communications_ , vol. 14, no. 6, pp. 3097–3107, 2015.
  * [22] C. Han, A. O. Bicen, and I. F. Akyildiz, “Multi-ray channel modeling and wideband characterization for wireless communications in the terahertz band,” _IEEE Transactions on Wireless Communications_ , vol. 14, no. 5, pp. 2402–2412, 2014.
  * [23] J. M. Jornet and I. F. Akyildiz, “Channel modeling and capacity analysis for electromagnetic wireless nanonetworks in the terahertz band,” _IEEE Transactions on Wireless Communications_ , vol. 10, no. 10, pp. 3211–3221, 2011\.
  * [24] D. Gutermuth, “Picard’s existence and uniqueness theorem,” _notes of Fundamental of Differential equations. https://embedded. eecs. berkeley. edu/eecsx44/lectures7Spring2013/Picard. pdf_.
  * [25] K. Ciesielski _et al._ , “On stefan banach and some of his results,” _Banach Journal of Mathematical Analysis_ , vol. 1, no. 1, pp. 1–10, 2007\.
  * [26] X. Yu, D. Xu, and R. Schober, “Miso wireless communication systems via intelligent reflecting surfaces,” in _IEEE International Conference on Communications in China_ , 2019, pp. 735–740.
  * [27] K. Gopalsamy, _Stability and oscillations in delay differential equations of population dynamics_. Springer Science & Business Media, 2013.
  * [28] V. V. Tarasova and V. E. Tarasov, “Concept of dynamic memory in economics,” _Communications in Nonlinear Science and Numerical Simulation_ , vol. 55, pp. 127–145, 2018.
  * [29] V. Tarasova and V. Tarasov, “A generalization of the concepts of the accelerator and multiplier to take into account of memory effects in macroeconomics,” _J. Econ. Entrep_ , vol. 10, pp. 1121–1129, 2016.
  * [30] V. V. Tarasova and V. E. Tarasov, “Logistic map with memory from economic model,” _Chaos, Solitons & Fractals_, vol. 95, pp. 84–91, 2017.
  * [31] S. Kang and J. B. Cheek, _Numerical solution of differential equations_. Waterways Experiment Station, 1972.
  * [32] X. Gao, L. Dai, Y. Zhang, T. Xie, X. Dai, and Z. Wang, “Fast channel tracking for terahertz beamspace massive mimo systems,” _IEEE Transactions on Vehicular Technology_ , vol. 66, no. 7, pp. 5689–5696, 2016.


11institutetext: CNR-ISTI – National Research Council, Institute of
Information Science and Technologies “Alessandro Faedo”, 56124 Pisa, Italy
11email<EMAIL_ADDRESS>22institutetext: OpenAIRE AMKE, Athens,
Greece

# “Knock knock! Who’s there?”
A study on scholarly repositories’ availability

Andrea Mannocci 11 0000-0002-5193-7851 Miriam Baglioni 11 0000-0002-2273-9004
Paolo Manghi 1122 0000-0001-7291-3210

###### Abstract

Scholarly repositories are the cornerstone of modern open science, and their
availability is vital for enacting its practices. To this end, scholarly
registries such as FAIRsharing, re3data, OpenDOAR and ROAR give them presence
and visibility across different research communities, disciplines, and
applications by assigning an identifier and persisting their profiles with
summary metadata. Alas, like any other resource available on the Web,
scholarly repositories, be they tailored for literature, software or data, are
quite dynamic and can be frequently changed, moved, merged or discontinued.
Therefore, their references are prone to link rot over time, and their
availability often boils down to whether the homepage URLs indicated in
authoritative repository profiles within scholarly registries respond or not.

For this study, we harvested the content of four prominent scholarly
registries and resolved over 13 thousand unique repository URLs. By performing
a quantitative analysis on such an extensive collection of repositories, this
paper aims to provide a global snapshot of their availability, which
bewilderingly is far from granted.

###### Keywords:

Scholarly repositories Availability HTTP resolution Scholarly communication
Open science

##  1 Introduction

Scholarly repositories are a vital part of the scholarly infrastructure and
therefore are a cornerstone of modern open science practices. Scholarly
registries, such as FAIRsharing111FAIRsharing – https://fairsharing.org [12],
re3data222re3data – https://re3data.org [11], OpenDOAR333OpenDOAR –
https://v2.sherpa.ac.uk/opendoar, and ROAR444ROAR –
http://roar.eprints.org/information.html, facilitate the discovery and
referencing of scholarly repositories by assigning them an identifier (either
local or persistent, as for FAIRsharing and re3data) and maintaining a public
profile displaying summary metadata, as in
https://fairsharing.org/FAIRsharing.wy4egf.

However, like any other resource on the Web [1, 2], scholarly repositories as
well, be they tailored for literature, software or data, are dynamic and often
changed, moved, merged or discontinued. Indeed, as scholars soon got to
understand, referencing consistently scholarly resources is far from
straightforward [5, 8, 4, 9]. Therefore, references to repositories are
vulnerable to link rot over time, and their availability often boils down,
especially for research infrastructures and scholarly communication services
aggregating content from registries, to whether their homepages are live and
able to respond or not.

In this paper, we harvested the content of four prominent scholarly
registries, namely FAIRsharing, re3data, OpenDOAR, and ROAR, and distilled a
comprehensive, longitudinal, repository-type-agnostic collection of over 13
thousand unique URLs pointing to scholarly repositories worldwide. Then, each
URL has been requested via HTTP with multiple methods and the response codes,
redirection lists and ultimately resolved URLs have been tracked, thus
providing a global snapshot of the availability of scholarly repositories.

The analysis highlights that about the 25% of repository URLs and homepages
registered in scholarly registries are problematic. Even more so, the requests
terminating successfully from a syntactical standpoint can still be incorrect
semantically, as the content served has nothing to relate with the original
scope of the repository. In our opinion, such results are bewildering, as the
URLs considered in this analysis are those available in official repository
profiles maintained by scholarly registries and not “vanilla” repository
references parsed “in the wild” from the research literature.

##  2 Related work

Web reference and link rot inside and outside the scholarly domain have been
extensively studied over the years [1, 2, 5, 6, 8, 4, 9]. However, to the best
of our knowledge, no prior study focused on the availability of scholarly
repositories by examining the URLs contained in repository profiles registered
by repository managers into scholarly registries.

Despite registries are aware of such issues, and they best-effort notify their
users (e.g., https://fairsharing.org/1724 or
https://www.re3data.org/repository/r3d100011299), in several cases, repository
profiles could be not yet flagged as problematic. Furthermore, the
authoritative nature of such URLs can open to implications regarding
repository management (and not just referencing) that are worth investigating.

In this work, we took inspiration from the methodology introduced in [5] for
studying the persistence of DOIs, which we adapted to the peculiarities of the
case study at hand.

##  3 Data and methods

For this analysis, we selected four major scholarly repository registries,
namely FAIRsharing, re3data, OpenDOAR and ROAR, whose details are briefly
summarised in Table 1. The content of the four registries was dumped by
various means in February 2022. Each repository profile in the four registries
was processed, and its homepage accrued in a list of 13,356 unique URLs.

Registry | Dump date | Dump method | Registry Licence | # repos
---|---|---|---|---
FAIRsharing | Feb 2022 | JSON (rest API) | CC-BY-SA | 1,853
re3data | Feb 2022 | OpenAIRE | CC-BY | 2,793
OpenDOAR | Feb 2022 | OpenAIRE | CC-BY-NC-ND | 5,811
ROAR | Feb 2022 | CSV (website) | CC-BY | 5,444
Table 1: Overview of the four registries considered in this analysis.

Each URL was requested via HTTP transactions, as documented in RFC 7231 [3],
comprehending a client issuing an HTTP request consisting of a request method
and request headers, and a server replying to such request with response
headers and (optionally) a response body. In this study, the HTTP requests
adhere to the two most common HTTP request methods, i.e., GET and HEAD. As
documented in the RFC, the main difference between HEAD and GET is that the
server returns a resource representation in the response body for GET
requests, which is instead omitted for HEAD requests. As per RFC 7231, we
expect identical results for the two request methods, except for negligible
differences in latency. HTTP requests were issued via Python, were allowed to
follow a maximum of 30 redirects, and had a timeout set to 30 seconds. For
each HTTP request successfully concluded, a number of response headers was
tracked, namely  (i) original URL;  (ii) final URL (same as the original or
redirected);  (iii) final status code;  (iv) redirection chain (in case of
redirects);  (v) redirection status codes (in case of redirects);  (vi)
latency.  If the server instead could not respond within this interval, the
request was marked as unsuccessful, and the error message (e.g., timeout,
connection refused) was noted. It is worth mentioning that we did not alter
the protocol in any way, be it http:// or https://, nor did we alter any other
part of the original URLs such as ports, as doing so could alter the perceived
error rate, despite reducing further the URLs space to process.

For the sake of open science and reproducibility, the content of the
registries, the URL list, the code issuing the requests, the collected data,
and the Jupyter notebooks for the analysis of the results are available on
Zenodo [10].

##  4 Results

We start our analysis by examining the last HTTP status code returned by the
last accessible link in the redirection chain. If no redirection takes place,
only one status code is returned. We aggregated the returned status codes in
Table 2, where rows indicate HTTP request methods, while columns indicate the
number of status code occurrences aggregated per response classes (i.e., 2xx
for successful requests, 4xx for client errors, and 5xx for server errors).

Method | 2xx | 4xx | 5xx | Err
---|---|---|---|---
HEAD | 74.83% (9,995) | 6.11% (816) | 0.86% (115) | 18.2% (2,431)
GET | 76.05% (10,158) | 5.19% (694) | 0.76% (102) | 18% (2,403)
Table 2: Final HTTP status codes.

The largest part of the performed HTTP requests returned 2xx-class status
codes (74.83% and 76.05% for HEAD and GET, respectively), while the remaining
requests exhibited issues of various kinds. More specifically, about 18% of
HEAD and GET requests failed due to time out, excessive retries or redirects,
connection reset by peer, malformed URLs, and so on. Furthermore, a
considerable amount of 4xx-class status codes was returned (i.e., 6.11% for
HEAD, 5.19% for GET), while a much smaller yet non-negligible quantity of 5xx-
class status codes was collected (i.e., 0.86% for HEAD, 0.76% for GET). Please
notice that the 3xx-class is relative to redirection, and therefore no final
HTTP status code is expected to belong to this class as all redirections were
followed, eventually resolving into a non-3xx status code or an error (e.g.,
for maxing out the number of redirects allowed, or for a timeout). A
significant divergence across the two request methods cannot be observed, even
though the slight difference suggests that requests were served differently
based on the method used, which infringes the RFC recommendations. Indeed,
https://idr.openmicroscopy.org triggers an error due to an exceeding number of
redirects for HEAD while returning a 200 OK for GET requests.

In total, it is safe to say that about one-quarter of the performed HTTP
requests failed to hit their target and resolve successfully, a result that we
found bewildering. This is particularly true, especially if we consider that
HTTP requests resolving into 2xx-class status codes, despite being
syntactically correct from the protocol perspective, can still be semantically
wrong in relation to the content of the web page served. Indeed, on manual
inspection of randomly sampled URLs, we verified that, despite the syntactic
correctness, the page resolved had nothing to do with academia or with the
repository/project anymore. For example, http://www.dlese.org/library/
(https://v2.sherpa.ac.uk/id/repository/425) resolves, without further
redirects, to a web page suggesting that the domain has not been renewed for
its original purposes and that someone else put it back onto the market.
Similarly, http://ejournal.windeng.net (http://roar.eprints.org/1530/),
together with five more URLs, redirects to http://survey-smiles.com via HEAD,
which ostensibly has nothing to do with the initial purpose of the repository.

This event led us to investigate redirections further, which per se do not
signal a sure malfunction (e.g., redirect from http:// to https://) while can
still contribute to semantic incorrectness as they could serve an unrelated
page (e.g., the repository has been moved somewhere else and the original URL
redirects to an error page). Redirection happens in about 32% of the requests,
i.e., 4,330 and 4,359 for HEAD and GET requests, respectively. We also noticed
that 378 redirection chains for HEAD requests ended up with a status code
other than 2xx-class status codes (338 in the 4xx-class, 40 in the 5xx-class),
while, for GET requests, 344 redirections chains shared a similar epilogue
(307 in the 4xx-class, 37 in the 5xx-class).

No redirection chain longer than five was observed, except for those maxing
out the number of redirections allowed. This was observed sporadically, eight
times for HEAD and once for GET, leaving little room for further speculations.
Most redirections resolved in one hop, and the number of N-long redirection
chains quickly decreases as N increases. We experimentally verified that the
most frequent reason for one-hop-long redirection chains was to switch from
http:// to https://, prepend www, append a trailing slash, or combinations of
such and similar alterations.

Incidentally, we noticed that the URL space dramatically shrinks, i.e., from
the over 13 thousand unique URLs we started with, excluding the errors, 10,926
of them resolved to 9,331 unique URLs via HEAD (-14.59%), while 10,954 of them
resolved to 9,353 unique URLs via GET (-14.61%). More specifically, for HEAD
requests, a pool of 3,051 unique repository URLs redirected at least twice
onto the same URL, for a total of 1,456 unique URLs after redirection. For GET
requests, we registered similar numbers: 3,060 unique repository URLs
redirected at least twice onto one of 1,459 unique URLs after redirection. As
it can be noticed, there is no substantial difference between the two methods,
yet the behaviour is not identical as it should be expected as per RFC
recommendations. This is also confirmed by 146 unique URLs for which HEAD and
GET returned two different locations.

By aggregating the final URLs after redirections and counting the number of
conflated original URLs, we discovered something interesting. In several
cases, for which unfortunately do not exist systematic identification
criteria, redirections hid issues worth commenting on. For example, 16
distinct URLs555http://brage.bibsys.no/hia/ {ks, hig, hiak, politihs, hsf,
hive, misjon, hinesna, hvo, hibo, histm, dhh, hint, hibu, bdh} redirected to
https://www.unit.no/ugyldig-lenke-til-dokument-i-vitenarkiv. The latter alerts
about the incorrectness of the (original) URLs (i.e., “You are trying to reach
a document via a link that is not valid”) and leaves the user with no further
indication of the repository’s actual location.

Similarly, 6 distinct URLs666https://toxnet.nlm.nih.gov/cgi-
bin/sis/htmlgen?CCRIS, https://toxnet.nlm.nih.gov/cgi-bin/sis/htmlgen?HSDB,
https://www.nlm.nih.gov/toxnet/index.html, https://toxmap.nlm.nih.gov/toxmap/,
https://toxnet.nlm.nih.gov/newtoxnet/tri.htm,
https://toxnet.nlm.nih.gov/newtoxnet/cpdb.htm conflated onto
https://www.nlm.nih.gov/toxnet/index.html, where a note claims that “most of
NLM’s toxicology information services have been integrated into other NLM
products and services”, thus leaving no trace of the previous databases and
relative information about the content. Similarly,
http://csdb.glycoscience.ru/help/migrate.html indicates that two original
databases redirecting here have been merged into a new one.

On a slightly different note, 7 distinct and seemingly unrelated,
URLs777http://content.wsulibs.wsu.edu/cdm/,
http://idahohistory.cdmhost.com/cdm/,
http://ccdl.libraries.claremont.edu/collection.php?alias=irw,
http://gettysburg.cdmhost.com, http://cdm16378.contentdm.oclc.org,
http://content.wsulibs.wsu.edu/, http://trinity.cdmhost.com/index.php all
redirected to https://www.oclc.org/url/notfound, whose domain is held by a
company positioned in the digital libraries & services market segment.

Alas, such problems can manifest for any number of conflated URLs, and making
a comprehensive list of such cases would equate to examining all the URLs by
hand, which quickly grows unfeasible.

##  5 Discussion

The current study suffers from a few limitations. Firstly, a one-shot
resolution provides a snapshot that institution-wide or nationwide
infrastructure outages could temporarily distort. As of June 2022, we tested
this hypothesis by checking Ukrainian repositories (i.e., .ua domain) and
verified that 36 out of 146 (24.6%) were unreachable, of which several are in
Kharkiv, east Ukraine. Similarly, other repositories can have downtime for
technical reasons. Therefore, repeating the experiment over an extended
timeframe would help recover from spurious transients, provide a better
overview of the actual situation, and put scholarly repositories’ availability
in a long-term perspective.

Furthermore, we resolved the URLs in this analysis via HEAD and GET only.
However, we empirically noticed that in some cases, URLs behave differently
when accessed via browser, a behaviour already observed in [7]. Therefore, we
could extend the study by including HTTP parameters (e.g., user agent,
accepted cookies) and by simulating human-like browser interactions via
Selenium888Selenium WebDriver – https://www.selenium.dev.

Nonetheless, the results show that about one out of four URLs from repository
profiles in scholarly registries is problematic. Moreover, this result is a
lower bound, as problems could be related to the content served after a
successfully served request. In fact, as we came to realise, many HTTP
requests returned, both with and without redirection, a web page unrelated to
the repository pertaining to the original URL. The methodology described in
[4] could help to assess the occurrence of such events; however, this would
apply only to the cases for which pre- and post- repository registration
snapshots exist. Unfortunately, an exhaustive estimation of content drift
would entail a fully manual inspection of repository URLs and the relevant
entries on the registries.

As a final remark, repository and IT infrastructure managers have a role in
such dysfunctions and should be aware that changes in the infrastructure or,
by any means, in the repository lifecycle should be notified and reflected
onto registries to strive for accountability and that an incorrect metadata
description of a repository can have an impact on reproducibility and open
science practices.

## Acknowledgements

This work was partially funded by the EC H2020 OpenAIRE-Nexus (Grant agreement
101017452).

## References

  * [1] Bar-Yossef, Z., Broder, A.Z., Kumar, R., Tomkins, A.: Sic transit gloria telae: Towards an understanding of the web’s decay. In: Proceedings of the 13th Conference on World Wide Web - WWW ’04. p. 328. ACM Press, New York, NY, USA (2004). https://doi.org/10.1145/988672.988716
  * [2] Cho, J., Garcia-Molina, H.: Estimating frequency of change. ACM Transactions on Internet Technology 3(3), 256–290 (Aug 2003). https://doi.org/10.1145/857166.857170
  * [3] Fielding, R., Reschke, J.: Hypertext transfer protocol (HTTP/1.1): Semantics and content (2014)
  * [4] Jones, S.M., de Sompel, H.V., Shankar, H., Klein, M., Tobin, R., Grover, C.: Scholarly Context Adrift: Three out of Four URI References Lead to Changed Content. PLOS ONE 11(12), e0167475 (Dec 2016). https://doi.org/10.1371/journal.pone.0167475
  * [5] Klein, M., Balakireva, L.: On the Persistence of Persistent Identifiers of the Scholarly Web. In: Hall, M., Merčun, T., Risse, T., Duchateau, F. (eds.) Digital Libraries for Open Knowledge. pp. 102–115. Lecture Notes in Computer Science, Springer International Publishing, Cham (2020). https://doi.org/10.1007/978-3-030-54956-5_8
  * [6] Klein, M., Balakireva, L.: An extended analysis of the persistence of persistent identifiers of the scholarly web. International Journal on Digital Libraries 23(1), 5–17 (Mar 2022). https://doi.org/10.1007/s00799-021-00315-w
  * [7] Klein, M., Balakireva, L., Shankar, H.: Who is asking? Humans and Machines experience a different Scholarly Web. In: 16th International Conference on Digital Preservation. Open Science Framework, Amsterdam (2019). https://doi.org/10.17605/OSF.IO/SMCY2
  * [8] Klein, M., de Sompel, H.V., Sanderson, R., Shankar, H., Balakireva, L., Zhou, K., Tobin, R.: Scholarly Context Not Found: One in Five Articles Suffers from Reference Rot. PLOS ONE 9(12), e115253 (Dec 2014). https://doi.org/10.1371/journal.pone.0115253
  * [9] Lawrence, S., Pennock, D., Flake, G., Krovetz, R., Coetzee, F., Glover, E., Nielsen, F., Kruger, A., Giles, C.: Persistence of Web references in scientific research. Computer 34(2), 26–31 (Feb 2001). https://doi.org/10.1109/2.901164
  * [10] Mannocci, A.: Analysis of scholarly repositories’ availability. Data and notebooks. (Jul 2022). https://doi.org/10.5281/zenodo.6906885, https://doi.org/10.5281/zenodo.6906885
  * [11] Pampel, H., Vierkant, P., Scholze, F., Bertelmann, R., Kindling, M., Klump, J., Goebelbecker, H.J., Gundlach, J., Schirmbacher, P., Dierolf, U.: Making Research Data Repositories Visible: The re3data.org Registry. PLOS ONE 8(11), e78080 (Nov 2013). https://doi.org/10.1371/journal.pone.0078080
  * [12] Sansone, S.A., McQuilton, P., Rocca-Serra, P., Gonzalez-Beltran, A., Izzo, M., Lister, A.L., Thurston, M.: FAIRsharing as a community approach to standards, repositories and policies. Nature Biotechnology 37(4), 358–367 (Apr 2019). https://doi.org/10.1038/s41587-019-0080-8