<figure><img src="image_1.png" /></figure>

HAL
open science

# Distributed Faulty Node Detection in Delay Tolerant Networks: Design and Analysis

Wenjie Li, Laura Galluccio, Francesca Bassi, Michel Kieffer

► To cite this version:

Wenjie Li, Laura Galluccio, Francesca Bassi, Michel Kieffer. Distributed Faulty Node Detection in Delay Tolerant Networks: Design and Analysis. IEEE Transactions on Mobile Computing, 2018, 17 (4), pp.831-844. 10.1109/TMC.2017.2743703. hal-01576587

HAL Id: hal-01576587

https://centralesupelec.hal.science/hal-01576587

Submitted on 23 Aug 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution - ShareAlike 4.0 International License
# Distributed Faulty Node Detection in Delay Tolerant Networks: Design and Analysis

Wenjie Li, Student Member, IEEE, Laura Galluccio, Member, IEEE, Francesca Bassi, Member, IEEE, and Michel Kieffer, Senior Member, IEEE

**Abstract**—Propagation of faulty data is a critical issue. In case of Delay Tolerant Networks (DTN) in particular, the rare meeting events require that nodes are efficient in propagating only correct information. For that purpose, mechanisms to rapidly identify possible faulty nodes should be developed. Distributed faulty node detection has been addressed in the literature in the context of sensor and vehicular networks, but already proposed solutions suffer from long delays in identifying and isolating nodes producing faulty data. This is unsuitable to DTNs where nodes meet only rarely. This paper proposes a fully distributed and easily implementable approach to allow each DTN node to rapidly identify whether its sensors are producing faulty data. The dynamical behavior of the proposed algorithm is approximated by some continuous-time state equations, whose equilibrium is characterized. The presence of misbehaving nodes, trying to perturb the faulty node detection process, is also taken into account. Detection and false alarm rates are estimated by comparing both theoretical and simulation results. Numerical results assess the effectiveness of the proposed solution and can be used to give guidelines for the algorithm design.

**Index Terms**—Delay Tolerant Networks; Fault detection; Iterative algorithms; Distributed estimation; Equilibrium analysis.

<figure><img src="image_3.png" /></figure>

## 1 INTRODUCTION

Delay Tolerant Networks (DTN) are challenging networks characterized by dynamic topology with frequent disconnections [1]. Examples of DTNs include Vehicular DTNs (VDTNs) [2] where two nodes can communicate with each other only when they are closely located. This connection is intermittent as the nodes are moving vehicles. Due to this sparse and intermittent connectivity, inference and learning over DTNs is much more complicated than in traditional networks, see, e.g., [3]–[8].

This paper considers the problem of distributed faulty node detection (DFD) in DTNs. A node is considered as faulty when one of its sensors frequently reports erroneous measurements. The identification of such faulty nodes is very important to save communication resources and to prevent erroneous measurements polluting estimates provided by the DTN. This identification problem is quite complicated in DTNs when interactions are mainly between pairs of encountering nodes. Most of the classical DFD algorithms are using measurements of spatially-correlated physical quantities collected by many nodes to determine the presence of outliers and identify the nodes producing these outliers. In case of pairwise interactions, mismatch between measurements provided by two different nodes

can still be detected, but identifying directly which node produces erroneous measurements is not possible.

This paper presents a fully distributed and easily implementable algorithm to allow each node of a DTN to determine whether its own sensors are defective. We assume as in [9] that nodes are not aware of the status (good or defective) of their sensors, while their computation and communication capabilities remain fine, even if some of their sensors are defective. Most of the nodes of the DTN are assumed to behave in a rational way and are willing to know the status of their sensors. Some nodes, however, may be misbehaving, trying to perturb the detection process.

As in [9]–[13], a Local Outlier Detection Test (LODT) is assumed to be able to detect the presence of outliers in a set of measurements, without necessarily being able to determine which are the outliers. This is a typical situation when only pairwise interactions are considered, where measurements from sensors of only two nodes are compared. The generic LODT is characterized by its probabilities of detection and false alarm. When two nodes meet, they exchange their local measurements and use them to perform the same LODT. The LODT results help both nodes to update their estimate of the status of their own sensors. When, for a given node, the proportion of meetings during which the LODT suggests the presence of outliers is larger than some threshold, this node decides its sensors may be defective. In this case, it becomes silent. Accordingly, it does not transmit any more its measurements to its neighbors, but keeps collecting measurements from nodes met and updates the estimate of the status of its sensors. It may then have the opportunity to change its estimate and communicate again. Although the LODT considered here are those of [9], this work differs significantly from [9] due to the communication conditions of DTNs, which require a totally

* W. Li, F. Bassi, and M. Kieffer are with Laboratoire des Signaux et Systemes (L2S, UMR CNRS 8506) CNRS-CentraleSupelec-Universitï¿½æ Paris-Sud 3, rue Joliot Curie 91192 Gif-sur-Yvette, France.
E-mail: firstname.name@lss.supelec.fr

* L. Galluccio is with Dipartimento di Ingegneria Elettrica, Elettronica e Informatica, University of Catania, Catania, Italy.
E-mail: laura.galluccio@dieei.unict.it

* M. Kieffer is also with LTCI Telecom ParisTech, 75013 Paris, France and Institut Universitaire de France, 75005 Paris, France

* F. Bassi is also with ESME-Sudria, 94200 Ivry-sur-Seine, France
different DFD algorithm. The analysis of the properties of the algorithm is also totally different. This paper shows that the behavior of the proposed DFD algorithm can be described using Markov models and tools borrowed from control theory and population dynamics.

More in depth, the belief of each node about the status of its sensors is quantized. The evolution of these quantized beliefs are then shown to follow two Markov chains. A continuous-time approximation of the evolution of the proportion of nodes with similar beliefs is then derived. Sufficient conditions on the decision parameters to ensure the existence and uniqueness of an equilibrium of the DFD algorithm are then provided. Given the characteristics of the LODT, upper and lower bounds of the detection rate, i.e., proportion of nodes which have effectively identified their sensors as defective, and of the false alarm rate, i.e., proportion of nodes which believe that their good sensors are in fact defective, are also obtained. The impact of misbehaving nodes, trying to perturb the results of the DFD algorithm, is also taken into account. These results provide guidelines to properly choose the parameters of the DFD algorithm.

The rest of the paper is organized as follows. Section 2 discusses some related work. Section 3 presents the system model and basic assumptions. Section 4 details the DFD algorithm for DTNs. Section 5 introduces the Markov model describing the behavior of the DFD algorithm and describes the transition probabilities between the node states. Section 6 develops the theoretical analysis of the macroscopic evolution of the proportion of nodes in different states. Section 7 analyzes the properties of the equilibrium obtained from the continuous-time state equations by approximating the stochastic evolution of the proportions of nodes with similar beliefs. Section 8 discusses the effect of having misbehaving nodes in the system. Section 9 provides some numerical results as well as a comparison with an alternative DFD algorithm and Section 10 concludes this paper. Notations are presented in Table 1. Proofs of propositions and lemmas are available in the Appendix.

## 2 RELATED WORK

DFD is a well-investigated topic when considering Wireless Sensor Networks (WSNs) (see [14]–[16] and references therein). The WSNs considered in most of the papers are dense and have a static topology. DFD in DTNs is much less investigated. Classical DFD algorithms usually consist of two phases. First, an LODT is performed using data collected from neighboring nodes. LODTs (based on majority voting [10], the median [11], or the mean [12] of the measurements, the modified three-sigma edit test [13], etc.) aim to decide which data is erroneous. Second, the outcomes of the LODTs are disseminated to improve the decision accuracy.

Nevertheless, when LODTs have to process measurements from two or three nodes only, it becomes difficult to identify the defective nodes. It may, however, still be possible to detect inconsistencies among measurements due to the presence of a node producing outliers. This is a typical situation in DTNs when there are mainly pairwise interactions: two nodes meet, take measurements, and share

these measurements. Applying directly classical DFD algorithms in DTNs may thus be quite ineffective. Moreover, usually the performance of DFD algorithms is characterized experimentally. A theoretical analysis of the equilibrium and convergence properties of these algorithms is seldom performed.

In the context of distributed estimation via consensus in a WSN, [17]–[20] have considered the simultaneous estimation of a common quantity from measurements corrupted by different levels of bias or of variance. A distributed ranking among nodes is performed according to the performance of their sensor. The proposed solution allows an identification of defective nodes with sensors producing measurements of high bias or high variance. Nevertheless, the proposed solution highly relies on the measurement models and on the communication conditions.

A problem related to DFD in DTNs has been considered in [21] in the context of VDTN. A large number of sensor nodes are fixed and some vehicles, called mobile carriers (MC) collect data from these sensors. The sensor nodes can only communicate with the MCs in their vicinity. A MC needs to collect enough measurements to perform a test to decide which have been produced by defective sensors. Once a node is deemed defective by a MC, it is added to its blacklist. The MC provides information to sensors about their status. MCs also exchange their blacklists to accelerate the faulty node detection.

In [22], a related problem of distributed malware detection in DTNs is addressed. Each node evaluates after a meeting with another node whether the latter has performed suspicious actions (malware transmission trial). When after several meetings with Node *j*, Node *i* detects for a given number of times suspicious activities, a cut-off decision is performed against Node *j*, which is ignored in next meetings. The main drawback of this approach is the long time required to identify and isolate misbehaving nodes. Misbehavior detection in DTNs is also considered in [6], [23], where the DTN is perturbed by routing misbehavior caused by selfish or malicious nodes. The identification approach in [6] is not distributed, since a trusted authority periodically checks the forwarding history of nodes to identify possible misbehavior. A collaborative approach is proposed in [23], where each node can detect whether the encountered node is selfish using a local watchdog. The detection result is disseminated over the network to increase the detection precision and to reduce the delay. Trust/Reputation management is another important aspect to help DTNs to resist various potential threats. For example, [24] provides an iterative trust management mechanism to fight against Byzantine attacks in which several nodes are totally controlled by the adversary. In [25], a defense against Sybil attacks [26] is introduced, which is based on the physical feature of the wireless propagation channel. A trust model in the scenario of underwater acoustic sensor networks is presented in [27] to take into account several trust metrics such as link trust, data trust, and node trust.

In this paper, differently from previous works in the field, we consider that in a distributed way each node performs a self-determination on whether its sensors are producing outliers in the context of DTNs. In this case, new issues arise, mainly related to the limited proximity time
<table><caption>TABLE 1
Symbols used in this paper</caption><thead><tr><th></th><th></th></tr></thead><tbody><tr><td>S<sub>0</sub>, S<sub>1</sub>, S<sub>2</sub></td><td>sets of good, defect, and misbehav. nodes</td></tr><tr><td>n<sub>S</sub></td><td>number of nodes</td></tr><tr><td>θ<sub>i</sub></td><td>status of node i</td></tr><tr><td>θ̂<sub>i</sub></td><td>estimate of θ<sub>i</sub></td></tr><tr><td>n<sub>θ</sub></td><td>number of nodes with status θ</td></tr><tr><td>p<sub>θ</sub></td><td>proportion of nodes with status θ</td></tr><tr><td>p<sup>θθ</sup></td><td>proportion of nodes with status θ and estimating their status as θ̂</td></tr><tr><td>p̂<sup>θθ</sup></td><td>value of p<sup>θθ</sup> at equilibrium</td></tr><tr><td>p̂<sup>θθ</sup></td><td>approximate value of p<sup>θθ</sup> at equilibrium</td></tr><tr><td>λ</td><td>inter-contact rate</td></tr><tr><td>ν</td><td>decision threshold</td></tr><tr><td>t</td><td>time</td></tr><tr><td>y<sub>i</sub></td><td>outcome of a LODT performed by node i</td></tr><tr><td>q<sub>D</sub></td><td>detection probability of a LODT</td></tr><tr><td>q<sub>FA</sub></td><td>false alarm probability of a LODT</td></tr><tr><td>c<sub>m,i</sub></td><td>number of LODTs performed by node i</td></tr><tr><td>c<sub>d,i</sub></td><td>number of LODTs by node i resulting in a detection of outliers</td></tr><tr><td>M</td><td>number of previous LODT results considered for the decision</td></tr><tr><td>x<sub>i</sub></td><td>state of node i, containing θ<sub>i</sub>, c<sub>m,i</sub>, and c<sub>d,i</sub></td></tr><tr><td>π<sup>δ<sub>m</sub>,δ<sub>d</sub></sup><sub>θ</sub>(t, c<sub>m</sub>, c<sub>d</sub>)</td><td>transition probability from state (θ, c<sub>m</sub>, c<sub>d</sub>) to state (θ, c<sub>m</sub> + δ<sub>m</sub>, c<sub>d</sub> + δ<sub>d</sub>)</td></tr><tr><td>X<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub></td><td>proportion of nodes of actual status θ with state x<sub>i</sub> = (θ, c<sub>m</sub>, c<sub>d</sub>)</td></tr><tr><td>X<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub></td><td>expected value of X<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub></td></tr><tr><td>X̃<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub></td><td>value of X<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub> at equilibrium</td></tr><tr><td>X̃<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub></td><td>approximate value of X<sup>c<sub>m</sub>,c<sub>d</sub></sup><sub>θ</sub> at equilibrium</td></tr></tbody></table>

of DTN nodes and the sporadic contacts which call for the consideration of the history of contacts in the identification process. Also, we provide a mathematical characterization of the problem and prove the convergence of the algorithm.

## 3 SYSTEM MODEL

Consider a set $S$ of $n_S$ moving nodes equipped with sensors. $S$ can be partitioned into three subsets, $S_0$, $S_1$, and $S_2$. $S_0$ contains all nodes with good sensors. $S_1$ is the subset of nodes with defective sensors producing outliers, i.e., measurements corrupted by a noise which has characteristics significantly different from those of the noise corrupting measurements provided by good sensors. Finally, $S_2$ represents the set of misbehaving nodes, trying to voluntarily perturb the behavior of the network.

The status of node $i$ is $θ_i(t) = 0$ (good node) if all its sensors are good, $θ_i(t) = 1$ (defective node) if at least one of them is defective, and by convention $θ_i(t) = 2$ (misbehaving node). The proportion and number of nodes with status $θ$ are respectively $p_θ$ and the number of nodes in status $θ$ is $n_θ = p_θ n_S$, with $p_0 + p_1 + p_2 = 1$. All nodes, except misbehaving nodes, are not aware of their own status. In what follows, we assume that over the time horizon of the experiment, the status of sensors does not change, i.e., $θ_i(t) = θ_i$.

Misbehaving nodes aim at disrupting network operations by causing congestion along paths, unreliable packet delivery, or erroneous data delivery [6], [24], [28]. Here, we assume that misbehaving nodes always transmit data to

their neighbors indicating that their sensors are good. Moreover, they try to provide measurements to the encountered nodes so that the LODTs performed by these nodes lead to the outcome of identifying outliers.

Our aim is (i) to design a distributed algorithm so that each node $i$ rapidly evaluates an accurate estimate $θ̂_i$ of its own status $θ_i$ despite the presence of misbehaving nodes, (ii) to provide a theoretical analysis of the behavior of this algorithm.

### 3.1 Communication model

Nodes can exchange information only during the limited time interval in which they are in vicinity. As in [7], [8], [23], [29], we assume that the time interval between two successive meetings follows an exponential distribution with an inter-contact rate $λ$. Moreover, we assume that each meeting involves only two nodes. When more than two nodes meet simultaneously, processing is performed pair-by-pair.

### 3.2 Local outlier detection test

As in [9], we consider a family of LODTs able to detect the presence of outliers in a set of $n$ data measurements $M = \{m_1, \dots, m_n\}$ but unable to identify which data is an outlier. Denote $y(M)$ the outcome of the LODT, i.e., $y(M) = 1$ if data corresponding to outliers are detected within $M$, otherwise, $y(M) = 0$.

LODTs can take various forms (see [9] and Example 1 below). The LODT is characterized by a false alarm probability $q_{FA}$ (the probability of having $y(M) = 1$ under the condition that none of the data in $M$ are produced by defective sensors) and by a detection probability $q_D$ (the probability of having $y(M) = 1$ under the condition that some data in $M$ are really produced by defective sensors). In $M$, let $n_0$ be the number of data produced by good sensors and $n_1$ be the number of data coming from defective sensors. We further assume that both $q_D$ and $q_{FA}$ in average depend only on the number of data involved in the LODT. As a consequence, we can denote $q_{FA}$ as $q_{FA}(n_0)$ and $q_D$ as $q_D(n_0, n_1)$. Each node performing a LODT on a set of data has not to know $n_0$ and $n_1$, but the performance of the LODT will depend on the actual values of $n_0$ and $n_1$, which are used in the analysis of the DFD algorithm.

**Example 1.** This example introduces a LODT inspired from bounded-error parameter estimation problem (see, e.g., [30]–[32]). It assumes that only pairwise interactions occur. Consider some sensor nodes taking temperature measurements in the same room, e.g., with actual value $t = 25°C$. For a non-defective sensor, suppose that its measurement error is bounded, e.g., $\pm 1°C$. Assume that two sensors provide $t_1 = 25.6°C$ and $t_2 = 23.5°C$ respectively. Assuming that none of the sensors is defective, and considering the bounded measurement noise, one deduces that $t \in t_1 = [t_1 - 1, t_1 + 1] = [24.6, 26.6]$ and $t \in t_2 = [22.5, 24.5]$. Since $t_1 \cap t_2 = \emptyset$, there exists at least one outlier, but one is unable to determine which sensor has produced an outlier.

### 3.3 Detection scenario

We assume that during each meeting of a pair of nodes $(i, j) \in S$, the nodes collect data with their sensors. Each
node may or may not transmit its data to the other node, as detailed in the DFD algorithm description. If a node has received data from its neighbor, it may run a LODT involving its own data and those received from its neighbor. We assume that the spatial and temporal correlation between data is such that only data collected during the meeting of two nodes can be exploited by a LODT. Therefore, previously collected data are not exploited. As a consequence, contrary to [9], where $n_0$ and $n_1$ may be large, in the DTN scenario, a LODT will involve $n_0, n_1 \in \{0, 1, 2\}$, with $n_0 + n_1 = 2$. In this paper, one furthermore assumes that $q_{FA}(2) < q_D(1, 1) \le q_D(0, 2)$, which is reasonable, since the outcome of a LODT is more likely to be 1 as the number of outliers involved increases.

## 4 DFD ALGORITHM

In the proposed DFD algorithm, each (good or defective) node $i$ manages two counters $c_{m,i}(t)$ and $c_{d,i}(t)$ initialized at 0 at $t = 0$. Using $c_{m,i}(t)$, node $i$ counts the number of meetings during which it has received data from its neighbor, and has been able to perform a LODT. Using $c_{d,i}(t)$, it counts the number of LODT resulting in the detection of outliers. When $c_{d,i}(t)/c_{m,i}(t) \ge \nu$, where $\nu$ is some constant decision threshold, node $i$ considers itself as carrying defective sensors, i.e., it sets its own estimate $\hat{\theta}_i(t) = 1$. Otherwise, it considers that its sensors are good, i.e., $\hat{\theta}_i(t) = 0$.

When a node with $\hat{\theta}_i(t) = 1$ encounters another node, it still takes measurements, but it does not send these data to the other node to avoid infecting the network with outliers. Any node, upon receiving data from another node, performs a LODT and updates $c_{m,i}(t)$ and $c_{d,i}(t)$. As a consequence, a node which meets another node considering itself as defective, transmits its data, but since it does not receive any data, it does not update $c_{m,i}(t)$ and $c_{d,i}(t)$ at the end of the meeting. Algorithm 1 summarizes the proposed DFD technique for an arbitrary reference node $i$.

The vector $\mathbf{x}_i(t) = (\theta_i, c_{m,i}(t), c_{d,i}(t))$ represents the (microscopic) state of each node $i$. As $t \to \infty$, one has $c_{m,i}(t) \to \infty$, which leads to an infinite number of possible values for $\mathbf{x}_i(t)$ and the global (macroscopic) behavior of the algorithm is difficult to analyze. To limit the number of possible states, we have chosen to consider the evolution of $c_{m,i}(t)$ and $c_{d,i}(t)$ over a sliding time window containing the time instants of the last $M$ meetings during which node $i$ has performed a LODT. Algorithm 2 is a modified version of Algorithm 1 accounting for this limited horizon. It involves an additional counter $\mu$ for the number of LODT performed by node $i$. As long as $\mu < M$, (5) is equivalent to (3).

Algorithm 2 is analyzed in the next sections.

## 5 EVOLUTION OF THE STATE OF A NODE

In this section, to simplify the presentation, the presence of misbehaving nodes is not taken into account. The impact of such nodes on the behavior of Algorithm 2 will be detailed in Section 8.

Consider the state $\mathbf{x}_i(t) = (\theta_i, c_{m,i}(t), c_{d,i}(t))$ of node $i$. Since $c_{m,i}(t) \in \{0, \dots, M\}$ and $c_{d,i}(t) \in \{0, \dots, c_{m,i}(t)\}$, the number of values that may be taken by the state of a

### Algorithm 1 DFD algorithm for node $i$

1) Initialize at $t_i^0 = 0, \hat{\theta}_i(t_i^0) = 0, c_{m,i}(t_i^0) = c_{d,i}(t_i^0) = 0, \kappa = 1$.

2) Do

$$ \left\{ \begin{array}{l} \hat{\theta}_i(t) = \hat{\theta}_i(t_i^{\kappa-1}), \\ c_{m,i}(t) = c_{m,i}(t_i^{\kappa-1}), \\ c_{d,i}(t) = c_{d,i}(t_i^{\kappa-1}), \end{array} \right. \qquad (1) $$

$$ t = t + \delta t \qquad (2) $$

until the $\kappa$-th meeting occurs at time $t = t_i^\kappa$ with Node $j^\kappa \in S \setminus \{i\}$.

3) Perform local measurement of data $m_i(t_i^\kappa)$.

4) If $\hat{\theta}_i(t_i^\kappa) = 0$, then transmit $m_i(t_i^\kappa)$ to node $j^\kappa$.

5) If data $m_{j^\kappa}$ have been received from node $j^\kappa$, then

a) Perform a LODT with outcome $y_i(t_i^\kappa)$.

b) Update $c_{m,i}$ and $c_{d,i}$ according to

$$ \left\{ \begin{array}{l} c_{m,i}(t_i^\kappa) = c_{m,i}(t_i^{\kappa-1}) + 1 \\ c_{d,i}(t_i^\kappa) = c_{d,i}(t_i^{\kappa-1}) + y_i(t_i^\kappa) \end{array} \right. \qquad (3) $$

c) Update $\hat{\theta}_i$ as follows

$$ \hat{\theta}_i(t_i^\kappa) = \begin{cases} 1 & \text{if } c_{d,i}(t_i^\kappa)/c_{m,i}(t_i^\kappa) \ge \nu, \\ 0 & \text{else.} \end{cases} \qquad (4) $$

6) $\kappa = \kappa + 1$.

7) Go to 2.

### Algorithm 2 Sliding-Window DFD algorithm for node $i$

1) Initialize $t_i^0 = 0, \hat{\theta}_i(t_i^0) = 0, c_{m,i}(t_i^0) = c_{d,i}(t_i^0) = 0, \kappa = 1$, and $\mu = 0$.

2) Do (1)-(2) until the $\kappa$-th meeting occurs at time $t_i^\kappa$ with Node $j^\kappa \in S \setminus \{i\}$.

3) Perform local measurement of data $m_i(t_i^\kappa)$.

4) If $\hat{\theta}_i(t_i^\kappa) = 0$, then transmit $m_i(t_i^\kappa)$ to node $j^\kappa$.

5) If data $m_{j^\kappa}$ have been received from node $j^\kappa$, then

a) $\mu = \mu + 1$. Perform a LODT with outcome $y_i^\mu$.

b) Update $c_{m,i}$ and $c_{d,i}$ as

$$ \left\{ \begin{array}{l} c_{m,i}(t_i^\kappa) = \min\{\mu, M\}, \\ c_{d,i}(t_i^\kappa) = \sum_{m=\max\{1, \mu-M+1\}}^{\mu} y_i^m. \end{array} \right. \qquad (5) $$

c) Update $\hat{\theta}_i$ using (4).

6) $\kappa = \kappa + 1$.

7) Go to 2.

node is $(M+1)(M+2)/2$. The evolution of $\mathbf{x}_i(t)$, conditioned by its status $\theta_i$, follows a Markov model with state transition diagram of the kind shown in Figure 1 for $M = 4$.

In particular, there are two chains, one conditioned by $\theta_i = 0$ and the other conditioned by $\theta_i = 1$. Both are characterized by a transient phase for state values with $c_{m,i}(t) < M$; then, a permanent regime starts when $c_{m,i}(t) = M$. With $c_{m,i}(t) = c_m$ and $c_{d,i}(t) = c_d$, the transitions from State $(\theta, c_m, c_d)$ to State $(\theta, c'_m, c'_d)$ are analyzed in the following.
<figure><img src="image_4.png" /><figcaption>Fig. 1. Example of Markov model for the evolution of the state $(\theta, c_{m,i}(t), c_{d,i}(t))$ of a node when $M = 4$.</figcaption></figure>

## 5.1 Analysis of some random events

In order to estimate the transition probability from a generic state of the Markov chain to another, one first calculates the probability a given node meets an other node believing its status good or its status is bad. One also evaluates the probability to perform a LODT outcome concluding in the absence or in the presence of outliers.

### 5.1.1 Probability of meeting a node believing its status is good/defective

Assume that at time $t$, some reference node $i$ meets an other node which index is represented by the random variable $J$ and define the random event $\mathcal{E}_1(t) = \{\hat{\theta}_J(t) = 0\}$, representing the event that the node met believes its status is good. As illustrated in (4), among the nodes with status $\theta$, the proportion of nodes that believe themselves as good is<fnref n="1" />

$$p^{\theta 0}(t) = \mathcal{X}_{\theta}^{c_{m}, 0}(t) + \sum_{c_{m} > 0, c_{d}: c_{d}/c_{m} < \nu} \mathcal{X}_{\theta}^{c_{m}, c_{d}}(t), \quad (6)$$

where $p^{10}(t)$ is in fact the non-detection rate (NDR) of the nodes with defective sensors at time $t$ and $\mathcal{X}_{\theta}^{c_{m}, c_{d}}(t)$ represents the proportion of nodes in state $(\theta, c_{m}, c_{d})$ defined as

$$\mathcal{X}_{\theta}^{c_{m}, c_{d}}(t) = |\mathcal{I}_{\theta}^{c_{m}, c_{d}}(t)|/n_{\theta}, \quad (7)$$

where

$$\mathcal{I}_{\theta}^{c_{m}, c_{d}}(t) = \{i : i \in S_{\theta}, \theta_{i} = \theta, c_{m,i}(t) = c_{m}, c_{d,i}(t) = c_{d}\},$$

and $|\cdot|$ denotes the cardinality of a set.

Assuming that the nodes are randomly spread, the probability that node $J$ believes it has only good sensors conditioned to its true status is

$$p^{\theta 0}(t) = \mathbb{P}(\hat{\theta}_J(t) = 0 | \theta_J(t) = \theta), \quad (8)$$

and then

$$\mathbb{P}\{\mathcal{E}_1(t)\} = p_0 p^{00}(t) + p_1 p^{10}(t). \quad (9)$$

Similarly, introduce $\mathcal{E}_1^*(t) = \{\hat{\theta}_J(t) = 1\}$, representing the event that the node met believes its status is defective. Among the nodes with sensors in status $\theta$, the proportion of nodes with $\hat{\theta}_J(t) = 1$ is

$$p^{\theta 1}(t) = \sum_{c_{m} > 0, c_{d}: c_{d}/c_{m} \ge \nu} \mathcal{X}_{\theta}^{c_{m}, c_{d}}(t), \quad (10)$$

where $p^{01}(t)$ and $p^{11}(t)$ represent the false alarm rate (FAR) and the detection rate (DR) respectively. From (10), one gets

$$\mathbb{P}\{\mathcal{E}_1^*(t)\} = p_0 p^{01}(t) + p_1 p^{11}(t). \quad (11)$$

### 5.1.2 Probability of detecting the presence of an outlier

Since node $i$ performs an LODT only when it meets a node $J$ with $\hat{\theta}_J(t) = 0$, one introduces the random event $\mathcal{E}_2^{\theta}(t) = \{Y_i(t) = 1 | \theta_i = \theta, \hat{\theta}_J(t) = 0\}$, for the reference node with actual status $\theta$. As discussed in Section 3.2, the statistical properties of the outcome $Y_i(t)$ of the LODT depend only on $\theta_i$ and $\theta_j$. For example, when node $i$ has good sensors, one has

$$\begin{align}
\mathbb{P}\{\mathcal{E}_2^{\theta}(t)\} &= \sum_{\varphi=0}^{1} \mathbb{P}\{Y_i(t) = 1, \theta_J = \varphi | \theta_i = 0, \hat{\theta}_J(t) = 0\} \nonumber \\
&\stackrel{(a)}{=} \sum_{\varphi=0}^{1} \mathbb{P}\{Y_i(t) = 1 | \theta_i = 0, \theta_J = \varphi\} \mathbb{P}\{\theta_J = \varphi | \hat{\theta}_J(t) = 0\} \nonumber \\
&\stackrel{(b)}{=} \frac{p_0 q_{\text{FA}}(2) p^{00}(t) + p_1 q_{\text{D}}(1, 1) p^{10}(t)}{p_0 p^{00}(t) + p_1 p^{10}(t)}. \tag{12}
\end{align}$$

In (12-a), one uses the fact that the LODT outcome is not influenced by the estimate of the status of a node and that in $\mathbb{P}\{\theta_J = \varphi | \theta_i = 0, \hat{\theta}_J(t) = 0\}$, the status of node $J$, does not depend on $\theta_i$. In (12-b),

$$\left\{ \begin{aligned} \mathbb{P}\{Y_i(t) = 1 | \theta_i = 0, \theta_J = 0\} &= q_{\text{FA}}(2), \\ \mathbb{P}\{Y_i(t) = 1 | \theta_i = 0, \theta_J = 1\} &= q_{\text{D}}(1, 1). \end{aligned} \right. \quad (13)$$

Moreover,

$$\mathbb{P}\{\theta_J = \varphi | \hat{\theta}_J(t) = 0\} = \frac{\mathbb{P}\{\hat{\theta}_J(t) = 0 | \theta_J = \varphi\} \mathbb{P}\{\theta_J = \varphi\}}{\sum_{\phi=0}^{1} \mathbb{P}\{\hat{\theta}_J(t) = 0 | \theta_J = \phi\} \mathbb{P}\{\theta_J = \phi\}} = \frac{p_{\varphi} p^{\theta 0}(t)}{p_0 p^{00}(t) + p_1 p^{10}(t)}.$$

If node $i$ has defective sensors, one has

$$\mathbb{P}\{\mathcal{E}_2^*(t)\} = \frac{p_0 q_D(1, 1) p^{00}(t) + p_1 q_D(0, 2) p^{10}(t)}{p_0 p^{00}(t) + p_1 p^{10}(t)}. \quad (14)$$

Similarly, one may introduce the random event $\mathcal{E}_3^{\theta}(t) = \{Y_i(t) = 0 | \theta_i = \theta, \hat{\theta}_J(t) = 0\}$, and show that

$$\mathbb{P}\{\mathcal{E}_3^{\theta}(t)\} = \begin{cases} \frac{p_0(1-q_{\text{FA}}(2))p^{00}(t)+p_1(1-q_{\text{D}}(1,1))p^{10}(t)}{p_0p^{00}(t)+p_1p^{10}(t)}, & \text{if } \theta = 0, \\ \frac{p_0(1-q_{\text{D}}(1,1))p^{00}(t)+p_1(1-q_{\text{D}}(0,2))p^{10}(t)}{p_0p^{00}(t)+p_1p^{10}(t)}, & \text{if } \theta = 1. \end{cases} \quad (15)$$

## 5.2 Transition probabilities

One evaluates now the transition probabilities for the state of a node. More specifically, define as $\pi_{\theta}^{\delta_m, \delta_d}$ the transition probability from State $(\theta, c_m, c_d)$ to State $(\theta, c_m + \delta_m, c_d + \delta_d)$, where $\theta \in \{0, 1\}$. One has $\delta_m \in \{0, 1\}$ since $c_m$ may either increase ($\delta_m = 1$) in the transient regime or remain constant ($\delta_m = 0$) in the permanent regime. One has $\delta_d \in \{-1, 0, 1\}$, depending on the value of the last LODT outcome and on the value of the $M + 1$-th last LODT outcome, which is no more considered in the permanent regime.

<footnote n="1">1. For the sake of simplicity, the dependency of $p^{\theta 0}(t)$ in $\nu$ is omitted, as $\nu$ is constant during the DFD algorithm.</footnote>
Thus, $(\delta_m, \delta_d) \in \{(0,0), (0,1), (0,-1), (1,0), (1,1), (1,-1)\}$ are independently distributed over the time horizon corresponding to $m = \mu - M, \dots, \mu - 1$. One obtains then
Note that $\pi_{\theta}^{\delta_m, \delta_d}$ depends on the current state of the reference node, but also on the current proportion of active (good and defective) nodes. Therefore, the transition probabilities are denoted as $\pi_{\theta}^{\delta_m, \delta_d}(t, c_m, c_d)$, where $t$ is the time instant, $c_{m,i}(t) = c_m$, and $c_{d,i}(t) = c_d$. Depending on the value of $c_m$, two different cases are considered in Section 5.2.1 and in Section 5.2.2, respectively corresponding to the transient and permanent regimes.

### 5.2.1 Case I, $c_{m,i}(t) < M$

In the transient regime, when $c_{m,i}(t) < M$, $c_{m,i}(t)$ and $c_{d,i}(t)$ are updated according to (3) whenever node $J$ with $\hat{\theta}_J(t) = 0$ is met. The only possibility that leads to $\delta_m = 0$ is the event $\mathcal{E}_1^*(t)$, i.e., node $i$ meets node $J$ with $\hat{\theta}_J(t) = 1$. As a consequence, no LODT is performed by node $i$. Therefore, for any $\theta \in \{0,1\}$,

$$ \pi_{\theta}^{0,0}(t, c_m, c_d) = \mathbb{P}\{\mathcal{E}_1^*(t)\} = p_0 p^{01}(t) + p_1 p^{11}(t), \quad (16) $$

where $p^{01}(t)$ is defined by (10).

A state transition occurs with $(\delta_m, \delta_d) = (1,1)$ when node $i$ with status $\theta_i = \theta$ meets node $J$ with $\hat{\theta}_J(t) = 0$ and when the LODT yields $y_i(t) = 1$. Since the two events are independent, one has

$$ \begin{aligned} \pi_{\theta}^{1,1}(t, c_m, c_d) &= \mathbb{P}\{Y_i(t) = 1, \hat{\theta}_J(t) = 0 | \theta_i = \theta\} \\ &= \mathbb{P}\{\mathcal{E}_1(t)\} \mathbb{P}\{\mathcal{E}_2^{\theta}(t)\}. \end{aligned} \quad (17) $$

Depending on the value of $\theta_i$, using (9), (12), and (14), one may rewrite (17) as

$$ \pi_{\theta}^{1,1}(t, c_m, c_d) = \begin{cases} p_0 q_{\text{FA}}(2) p^{00}(t) + p_1 q_{\text{D}}(1,1) p^{10}(t), & \text{if } \theta = 0, \\ p_0 q_{\text{D}}(1,1) p^{00}(t) + p_1 q_{\text{D}}(0,2) p^{10}(t), & \text{if } \theta = 1. \end{cases} \quad (18) $$

Finally, $\pi_{\theta}^{1,0}(t, c_m, c_d) = \mathbb{P}\{Y_i(t) = 0, \hat{\theta}_J(t) = 0 | \theta_i = \theta\}$ is obtained similarly from (15)

$$ \pi_{\theta}^{1,0}(t, c_m, c_d) = \begin{cases} p_0(1-q_{\text{FA}}(2))p^{00}(t) + p_1(1-q_{\text{D}}(1,1))p^{10}(t), & \text{if } \theta = 0, \\ p_0(1-q_{\text{D}}(1,1))p^{00}(t) + p_1(1-q_{\text{D}}(0,2))p^{10}(t), & \text{if } \theta = 1. \end{cases} \quad (19) $$

### 5.2.2 Case II, $c_{m,i}(t) = M$

In the permanent regime, $c_{m,i}(t) = M$ and does not increase any more, thus $\delta_m = 0$. In Algorithm 2, $\mu$ is the number of LODTs performed by node $i$ up to time $t$. When $\mu \ge M$, only the last $M$ LODT outcomes are considered: LODT outcomes $y_i^m$ with $m \le \mu - M$ are no more considered.

To determine the value taken by $\delta_d \in \{-1,0,1\}$ after the $\mu$-th LODT, consider the random event

$$ \mathcal{E}_4^1(t) = \left\{ Y_i^{\mu-M} = 1 \mid \sum_{m=\mu-M}^{\mu-1} Y_i^m = c_d \right\}, \quad (20) $$

which corresponds to a situation where one knows that $c_d$ LODTs where positive among the last $M$ tests and the LODT that will be ignored, once the new LODT outcome is available, also concluded in the presence of defective sensors. $\mathbb{P}\{\mathcal{E}_1^1(t)\}$ is relatively complex to evaluate, since $\mathbb{P}\{Y_i^n = 1\}$ is time-varying according to (12-14). In what follows, we assume that LODT outcomes with $Y_i^m = 1$

$$ \mathbb{P}\{\mathcal{E}_4^1(t)\} = \frac{c_d}{M}. \quad (21) $$

This approximation is exact in steady-state, when the $\mathcal{X}_{\theta}^{c_m, c_d}$s do not vary any more.

Similarly, define

$$ \mathcal{E}_4^0(t) = \left\{ Y_i^{\mu-M} = 0 \mid \sum_{m=\mu-M}^{\mu-1} Y_i^m = c_d \right\}. \quad (22) $$

Considering the same assumption used to get (21), one has

$$ \mathbb{P}\{\mathcal{E}_4^0(t)\} = 1 - \mathbb{P}\{\mathcal{E}_4^1(t)\} \approx \frac{M-c_d}{M}. \quad (23) $$

Assume that the $(\mu - M)$-th LODT performed by node $i$ occurred at time $\tilde{t}$, then $y_i^{\mu-M}$ can also be denoted as $y_i(\tilde{t})$ and the transition related to $c_{d,i}$ is such that $\delta_d = y_i(t) - y_i(\tilde{t}) \in \{-1,0,1\}$.

To have $(\delta_m, \delta_d) = (0,1)$, three independent events have to occur: 1) the encountered node $J$ believes it is good at time $t$, i.e., $\mathcal{E}_1(t)$; 2) $y_i(t) = 1$, i.e., $\mathcal{E}_3^{\theta}(t)(t)$; 3) $y_i(t) = 0$, i.e., $\mathcal{E}_4^0(t)$. Thus the transition probability may be expressed as

$$ \pi_{\theta}^{0,1}(t, M, c_d) = \mathbb{P}\{\mathcal{E}_1(t)\} \mathbb{P}\{\mathcal{E}_3^{\theta}(t)\} \mathbb{P}\{\mathcal{E}_4^0(t)\}. \quad (24) $$

Using (9), (12), (14), and (21) in (24), one gets

$$ \begin{aligned} \pi_{\theta}^{0,1}(t, M, c_d) &= \left\{ \begin{aligned} & (p_0 q_{\text{FA}}(2) p^{00}(t) + p_1 q_{\text{D}}(1,1) p^{10}(t)) \frac{M-c_d}{M}, & \text{if } \theta = 0, \\ & (p_0 q_{\text{D}}(1,1) p^{00}(t) + p_1 q_{\text{D}}(0,2) p^{10}(t)) \frac{M-c_d}{M}, & \text{if } \theta = 1. \end{aligned} \right. \end{aligned} \quad (25) $$

Consider now $(\delta_m, \delta_d) = (0,-1)$. To have such transition, the three following independent events should occur: 1) $\mathcal{E}_1(t)$; 2) $y_i(t) = 0$, i.e., $\mathcal{E}_3^{\theta}(t)(t)$; 3) $y_i(t) = 1$, i.e., $\mathcal{E}_4^1(t)$. Thus, the transition probability is

$$ \begin{aligned} \pi_{\theta}^{0,-1}(t, M, c_d) &= \mathbb{P}\{\mathcal{E}_1(t)\} \mathbb{P}\{\mathcal{E}_3^{\theta}(t)\} \mathbb{P}\{\mathcal{E}_4^1(t)\} \\ &= \left\{ \begin{aligned} & (p_0(1-q_{\text{FA}}(2)) p^{00}(t) + p_1(1-q_{\text{D}}(1,1)) p^{10}(t)) \frac{c_d}{M}, & \text{if } \theta = 0, \\ & (p_0(1-q_{\text{D}}(1,1)) p^{00}(t) + p_1(1-q_{\text{D}}(0,2)) p^{10}(t)) \frac{c_d}{M}, & \text{if } \theta = 1. \end{aligned} \right. \end{aligned} \quad (26) $$

Finally, by substituting eqs. (25-26) it is possible to calculate $\pi_{\theta}^{0,0}(t, M, c_d)$ which is given by

$$ \pi_{\theta}^{0,0}(t, M, c_d) = 1 - \pi_{\theta}^{0,1}(t, M, c_d) - \pi_{\theta}^{0,-1}(t, M, c_d). \quad (27) $$

In this section, we have so far completely characterized the transition probabilities between any possible pair of states in the Markov chain. Accordingly, we are now able to completely describe the evolution of the DTN state components and, thus, the expected proportion of nodes in a specific state.

## 6 MACROSCOPIC EVOLUTION OF THE DTN STATE

All node state transition probabilities evaluated in Section 5 are now used to determine the evolution of the proportion of nodes in state $\theta$, i.e.

$$ \mathcal{X}_{\theta}(t) = (\mathcal{X}_{\theta}^{0,0}(t), \mathcal{X}_{\theta}^{1,0}(t), \mathcal{X}_{\theta}^{1,1}(t), \dots, \mathcal{X}_{\theta}^{M,0}(t), \dots, \mathcal{X}_{\theta}^{M,M}(t)) $$

and the corresponding expected values

$$ X_{\theta}(t) = (X_{\theta}^{0,0}(t), X_{\theta}^{1,0}(t), X_{\theta}^{1,1}(t), \dots, X_{\theta}^{M,0}(t), \dots, X_{\theta}^{M,M}(t)). $$
**Proposition 2.** The evolution of the DTN state components, i.e., the expected proportion of nodes $X_{\theta}^{c_m, c_d}(t)$ in the states $(\theta, c_m, c_d)$, with $\theta \in \{0, 1\}$, $c_m = 0, \dots, M$, and $c_d \le c_m$ is described by

$$ \left\{
\begin{aligned}
\frac{dX_{\theta}^{0,0}}{dt} & \stackrel{(a)}{=} -\lambda X_{\theta}^{0,0} (\pi_{\theta}^{1,0} (0,0) + \pi_{\theta}^{1,1} (0,0)), \\
\frac{dX_{\theta}^{c_m,0}}{dt} & \stackrel{(b)}{=} \lambda \left( -X_{\theta}^{c_m,0} (\pi_{\theta}^{1,0} (c_m,0) + \pi_{\theta}^{1,1} (c_m,0)) + X_{\theta}^{c_m-1,0} \pi_{\theta}^{1,0} (c_m-1,0) \right), \\
\frac{dX_{\theta}^{c_m,c_m}}{dt} & \stackrel{(c)}{=} \lambda \left( -X_{\theta}^{c_m,c_m} (\pi_{\theta}^{1,0} (c_m,c_m) + \pi_{\theta}^{1,1} (c_m,c_m)) + X_{\theta}^{c_m-1,c_m-1} \pi_{\theta}^{1,1} (c_m-1,c_m-1) \right), \\
\frac{dX_{\theta}^{M,0}}{dt} & \stackrel{(d)}{=} \lambda \left( -X_{\theta}^{M,0} \pi_{\theta}^{0,1} (M,0) + X_{\theta}^{M-1,0} \pi_{\theta}^{0,1} (M-1,0) + X_{\theta}^{M,1} \pi_{\theta}^{0,-1} (M,1) \right), \\
\frac{dX_{\theta}^{M,M}}{dt} & \stackrel{(e)}{=} \lambda \left( -X_{\theta}^{M,M} \pi_{\theta}^{0,-1} (M,M) + X_{\theta}^{M,M-1} \pi_{\theta}^{0,1} (M,M-1) + X_{\theta}^{M-1,M-1} \pi_{\theta}^{1,1} (M-1,M-1) \right),
\end{aligned}
\right.
\quad (28) $$

for any $c_m = 1, \dots, M-1$, with the initial conditions $X_{\theta}^{0,0}(0) = 1$ and $X_{\theta}^{c_m, c_d}(0) = 0, \forall c_m, c_d \neq 0$.

*Proof:* See Appendix A □

Kurtz's theorem [33], [34] can then be used to show that for all $\varepsilon > 0$, there exists $\alpha_1 \ge 0$ and $\alpha_2 (\varepsilon) > 0$ such that

$$ \mathrm{Pr} \left( \max_{t \in [0,T]} \| \mathcal{X}_{\theta}(t) - X_{\theta}(t) \| > \varepsilon \right) \le \alpha_1 \exp(-\alpha_2(\varepsilon) n_{\theta}). $$

As a consequence, $\mathcal{X}_{\theta}(t)$ converges in probability to $X_{\theta}(t)$ as $n_{\theta}$ goes to infinity. This is typically the approximation performed in the seminal work [35] where the SIR model was proposed. This model is the one used to characterize most widely studied classes of epidemic models. Accordingly, analogously to what was presented for example in [7], [35]–[40], the proposed system consists of ordinary differential equations approximating jump Markov processes.

The state equations in (28) are nonlinear, since each $\pi_{\theta}^{\delta_m, \delta_d}$ depends on $X_{\theta}^{c_m, c_d}$, see (6) and (10).

# 7 ANALYSIS OF THE DTN STATE EQUATIONS

In what follows, the asymptotic behavior of the DTN state equations (28) is characterized. Algorithm 2 may drive $X_{\theta}^{c_m, c_d}$ to an equilibrium $\bar{X}_{\theta}^{c_m, c_d}$ at which the proportions of nodes in different states $X_{\theta}^{c_m, c_d}(t)$ do not vary any more. As a consequence, $p^{\theta 0}(t)$ defined in (6) also tends to an equilibrium $\bar{p}^{\theta 0}$.

## 7.1 Equilibrium of $X_{\theta}^{c_m, c_d}$

One investigates first the evolution of $X_{\theta}^{c_m, c_d}(t)$ when $c_m < M$. As shown in the following proposition, the DTN state always reaches the permanent regime.

**Proposition 3.** For any $c_m < M$ and $c_d \le c_m$, $\lim_{t \to \infty} X_{\theta}^{c_m, c_d}(t) = 0$.

*Proof:* See Appendix B. □

From Proposition 3, the only possible value at equilibrium of $X_{\theta}^{c_m, c_d}(t)$ when $c_m < M$ is 0. Thus $\bar{p}^{\theta 0}$ may be written as

$$ \bar{p}^{\theta 0} = \sum_{c_d:c_d/M<\nu} \bar{X}_{\theta}^{M,c_d}. \quad (29) $$

Denote $\bar{p} = (\bar{p}^{00}, \bar{p}^{10}) \in \mathcal{P}_0$ with

$$ \mathcal{P}_0 = \{(x, y) \in [0, 1] \times [0, 1] \text{ and } (x, y) \neq (0, 0)\} \quad (30) $$

and consider the functions

$$ h_0(\bar{p}) = \frac{p_0 q_{\mathrm{FA}} (2) \bar{p}^{00} + p_1 q_{\mathrm{D}} (1, 1) \bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10}}, \quad (31) $$

$$ h_1(\bar{p}) = \frac{p_0 q_{\mathrm{D}} (1, 1) \bar{p}^{00} + p_1 q_{\mathrm{D}} (0, 2) \bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10}}, \quad (32) $$

$$ F_{\theta}(\bar{p}) = \sum_{c_d=0}^{[M\nu]-1} \binom{M}{c_d} (h_{\theta}(\bar{p}))^{c_d} (1 - h_{\theta}(\bar{p}))^{M-c_d}, \quad (33) $$

and $\mathbf{F}(\bar{p}) = (F_0(\bar{p}), F_1(\bar{p}))$. The following proposition provides a non-linear equation that has to be satisfied by $\bar{p}$. The various $\bar{X}_{\theta}^{M, c_d}$ at equilibrium are easily deduced from the solutions of the mentioned equation.

**Proposition 4.** Assume that the dynamic system described by (28) admits some equilibrium $\bar{X}_{\theta}^{c_m, c_d}$, then $\bar{p} \in \mathcal{P}_0$ is the solution of

$$ \bar{p} = \mathbf{F}(\bar{p}), \quad (34) $$

and for any $\theta \in \{0, 1\}$ and $c_d \le c_m$,

$$ \bar{X}_{\theta}^{c_m, c_d} = \begin{cases} 0, & \forall c_m < M, \\ \binom{M}{c_d} (h_{\theta}(\bar{p}))^{c_d} (1 - h_{\theta}(\bar{p}))^{M-c_d}, & c_m = M. \end{cases} \quad (35) $$

*Proof:* See Appendix C. □

## 7.2 Existence and uniqueness of the equilibrium point

Now we investigate the existence and the uniqueness of the solution of (34), which is rewritten in detail in (36) at the top of the next page.

For that purpose, using fixed-point theorems, one may alternatively show that for all $\mathbf{p}(0) = (p^{00}(0), p^{10}(0)) \in \mathcal{P}_0$, the discrete-time system

$$ \left\{ \begin{aligned} p^{00}(n+1) &= F_0 (p^{00}(n), p^{10}(n)), \\ p^{10}(n+1) &= F_1 (p^{00}(n), p^{10}(n)). \end{aligned} \right. \quad (39) $$

converges to a unique equilibrium point $(\bar{p}^{00}, \bar{p}^{10})$, which is then solution of (36).

One first shows the existence of an equilibrium using Brouwer's fixed-point theorem [41] in the following proposition.

**Proposition 5.** For any $\nu \in [0, 1]$, (36) always admits a solution, which is an equilibrium point of the dynamical system (28).

Before proving Proposition 5, one first shows that $p^{00}(n)$ and $p^{10}(n)$ are contained in intervals with lower (and upper) bounds increasing (resp. decreasing) with $n$.

**Lemma 6.** For any $n \in \mathbb{N}^*$ and $\theta \in \{0, 1\}$, one has

$$ p_{\min}^{\theta 0}(n) \le p^{\theta 0}(n) \le p_{\max}^{\theta 0}(n), $$

with $p_{\min}^{\theta 0}(0) = 0, p_{\max}^{\theta 0}(0) = 1$, and

$$ \left\{ \begin{aligned} p_{\min}^{\theta 0}(n+1) &= F_{\theta} (p_{\min}^{00}(n), p_{\max}^{10}(n)), & \forall n \in \mathbb{N}^{+}, \\ p_{\max}^{\theta 0}(n+1) &= F_{\theta} (p_{\max}^{00}(n), p_{\min}^{10}(n)), & \forall n \in \mathbb{N}^{+}. \end{aligned} \right. \quad (40) $$

Moreover,

$$ p_{\min}^{00}(n+1) > p_{\min}^{00}(n), \quad p_{\max}^{00}(n+1) < p_{\max}^{00}(n). \quad (41) $$

*Proof:* See Appendix D. □

Using Lemma 6, one can now prove Proposition 5.

*Proof:* $F_0$ and $F_1$ are both continuous functions. For some $n > 0$, consider the set $\mathcal{P}_n = [p_{\min}^{00}(n), p_{\max}^{00}(n)] \times [p_{\min}^{10}(n), p_{\max}^{10}(n)]$, where $p_{\min}^{00}(n)$ and $p_{\max}^{00}(n)$ are defined
$$ \left\{ \begin{array}{l} \bar{p}^{00} = F_0 (\bar{p}^{00}, \bar{p}^{10}) = \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} \left( \frac{p_0 q_{\text{FA}}(2) \bar{p}^{00} + p_1 q_D(1,1) \bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10}} \right)^{c_d} \left( \frac{p_0(1-q_{\text{FA}}(2))\bar{p}^{00} + p_1(1-q_D(1,1))\bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10}} \right)^{M-c_d}, \\ \\ \bar{p}^{10} = F_1 (\bar{p}^{00}, \bar{p}^{10}) = \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} \left( \frac{p_0 q_D(1,1) \bar{p}^{00} + p_1 q_D(0,2) \bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10}} \right)^{c_d} \left( \frac{p_0(1-q_D(1,1))\bar{p}^{00} + p_1(1-q_D(0,2))\bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10}} \right)^{M-c_d}. \end{array} \right. \quad (36) $$

$$ c_0(q_{\text{FA}}(2), q_D(0,2), q_D(1,1), p_1, M, \nu, n) = \frac{M (q_D(1,1) - q_{\text{FA}}(2)) p_0 p_1 p_{\text{max}}^{00} (n) p_{\text{max}}^{10} (n)}{(p_0 p_{\text{min}}^{00} (n) + p_1 p_{\text{min}}^{10} (n)) ((1 - q_{\text{FA}}(2)) p_0 p_{\text{min}}^{00} (n) + (1 - q_D(1,1)) p_1 p_{\text{min}}^{10} (n))}, \quad (37) $$

$$ c_1(q_{\text{FA}}(2), q_D(0,2), q_D(1,1), p_1, M, \nu, n) = \frac{M (q_D(0,2) - q_D(1,1)) p_0 p_1 p_{\text{max}}^{00} (n) p_{\text{max}}^{10} (n)}{(p_0 p_{\text{min}}^{00} (n) + p_1 p_{\text{min}}^{10} (n)) ((1 - q_D(1,1)) p_0 p_{\text{min}}^{00} (n) + (1 - q_D(0,2)) p_1 p_{\text{min}}^{10} (n))}, \quad (38) $$

<figure><img src="image_5.png" /><figcaption>Fig. 2. Upper bounds of $\nu$ to satisfy (42), with $q_{\text{FA}}(2) = 0.05$, $q_D(0,2) = 0.9$, $q_D(1,1) \in \{0.5, 0.8\}$, $M \in \{4, 10\}$, and $p_1 \in [0.05, 0.5]$.</figcaption></figure>

in (40). For any $\mathbf{p} = (p^{00}, p^{10}) \in \mathcal{P}_n$, one can prove using Lemma 6 that $\mathbf{F}(\mathbf{p}) \in \mathcal{P}_n$. Thus $\mathbf{F}$ maps $\mathcal{P}_n$ to $\mathcal{P}_n$. Applying Brouwer's fixed-point theorem, $\mathbf{F}$ admits a fixed point and Proposition 5 is proved. $\square$

Sufficient conditions on $p_0, p_1, q_D, q_{\text{FA}}, M$ and $\nu$ are then provided to ensure the uniqueness of this equilibrium by applying Banach's fixed-point theorem [42].

**Proposition 7.** If there exists some $N'$, such that $\forall \theta \in \{0, 1\}$ and $\forall n > N'$, one has

$$ c_{\theta}(q_{\text{FA}}(2), q_D(0,2), q_D(1,1), p_1, M, \nu, n) < 1, \quad (42) $$

where $c_0$ and $c_1$ are defined in (37-38), then the discrete-time system (39) converges to a unique equilibrium point and the solution of (36) is unique.

*Proof:* See Appendix E. $\square$

Due to the monotonicity of $p_{\text{min}}^{00}(n)$ and $p_{\text{max}}^{00}(n)$ shown in Lemma 6, $c_{\theta}$ decreases with $n$. Hence, if a given $\nu$ satisfies (42) for some $N'$, then $\nu$ will satisfy (42) for all $n \ge N'$ and the equilibrium is unique. If the values of $p_1, q_D, q_{\text{FA}},$ and $M$ are fixed, then one may deduce sufficient conditions on the value of $\nu$ to have a unique equilibrium point. See Example 8.

**Example 8.** Consider $q_{\text{FA}}(2) = 0.05$, $q_D(0,2) = 0.9$, $q_D(1,1) \in \{0.5, 0.8\}$, $M \in \{4, 10\}$, and $p_1 \in [0.05, 0.5]$. One verifies whether (42) is satisfied considering $n = 10$ for different values of $\nu$. One obtains that (42) holds if $0 < \nu \le \nu_{\text{max}}$, where $\nu_{\text{max}}$ depends on the values of $p_1, q_D, q_{\text{FA}},$ and $M$. See Figure 2 for the numerical values of $\nu_{\text{max}}$ in each case.

### 7.3 Equilibrium point as $M \to \infty$

Both $\bar{p}^{00}$ and $\bar{p}^{10}$ can be seen as functions of $M$. As $M \to \infty$, Algorithm 2 turns into Algorithm 1. In this situation, if $\nu$ is

properly chosen, the probabilities of false alarm and non-detection tend to zero, as shown in Proposition 9.

**Proposition 9.** If $q_{\text{FA}}(2) < \nu < q_D(1,1)$, then (36) has a unique solution and

$$ \lim_{M \to \infty} \bar{p}^{00} = 1, \quad \lim_{M \to \infty} \bar{p}^{10} = 0. \quad (43) $$

*Proof:* See Appendix F. $\square$

### 7.4 Approximations of the Equilibrium

Closed-form expressions for $\bar{p}^{00}$ and $\bar{p}^{10}$ are difficult to obtain from (36). This section introduces an approximation of (36) from which some insights may be obtained on the way $\nu$ should be chosen.

Since $\bar{p}^{10}$ represents the expected proportion of nodes with defective sensors that have not detected their status, the value of $\bar{p}^{10}$ should be small. From (31-32) one sees that $\lim_{\bar{p}^{10} \to 0} h_0 = q_{\text{FA}}(2)$ and $\lim_{\bar{p}^{10} \to 0} h_1 = q_D(1,1)$, thus one may consider the following approximations

$$ h_0 \approx \tilde{h}_0 = q_{\text{FA}}(2), \quad h_1 \approx \tilde{h}_1 = q_D(1,1). \quad (44) $$

Therefore, (36) may be rewritten as

$$ \left\{ \begin{array}{l} \tilde{p}^{00} = \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} (q_{\text{FA}}(2))^{c_d} (1 - q_{\text{FA}}(2))^{M-c_d}, \\ \tilde{p}^{10} = \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} (q_D(1,1))^{c_d} (1 - q_D(1,1))^{M-c_d}. \end{array} \right. \quad (45) $$

from which one deduces approximate values $\tilde{X}_0^{M,c_d}$ of $X_0^{M,c_d}$ at equilibrium from eq. (35)

$$ \left\{ \begin{array}{l} \tilde{X}_0^{M,c_d} = \binom{M}{c_d} (q_{\text{FA}}(2))^{c_d} (1 - q_{\text{FA}}(2))^{M-c_d}, \\ \tilde{X}_1^{M,c_d} = \binom{M}{c_d} (q_D(1,1))^{c_d} (1 - q_D(1,1))^{M-c_d}. \end{array} \right. \quad (46) $$

For any fixed value of $M$, $q_{\text{FA}}(2)$, and $q_D(1,1)$, the values of detection rate ($\bar{p}^{11}$) and false alarm rate ($\bar{p}^{01}$) at equilibrium can be predicted using (45), since $\bar{p}^{01} = 1 - \bar{p}^{00}$ and $\bar{p}^{11} = 1 - \bar{p}^{10}$. Consider for example $M = 10$, $q_{\text{FA}}(2) = 0.05$, and $q_D(1,1) = 0.8$. Figure 3 presents $\bar{p}^{11}$ as a function of $\bar{p}^{01}$ for different values of $\nu$. This figure is helpful to choose the value of $\nu$ to meet different performance requirements. The actual values of $\bar{p}^{11}$ and $\bar{p}^{01}$ are also shown in Figure 3, which are very close to $\tilde{p}^{11}$ and $\tilde{p}^{01}$, in the region where $\bar{p}^{11}$ is close to 1.

## 8 INFLUENCE OF MISBEHAVING NODES

A LODT involving data coming from a misbehaving node will always result in the detection of an outlier. Thus, when a node $i$ with state $\mathbf{x}_i(t) = (\theta, c_{m,i}(t), c_{d,i}(t))$ meets a misbehaving node, the possible transitions are such that

$$ (\delta_m, \delta_d) = (1,0) \text{ or } (\delta_m, \delta_d) = (1,1) \text{ if } c_{m,i}(t) < M $$
<figure><img src="image_6.png" /><figcaption>Fig. 3. Approximate $p^{11}$ as a function of approximate $p^{01}$, for various $\nu$ and fixed $M = 10$.</figcaption></figure>

• $(\delta_m, \delta_d) = (0, 0)$ or $(\delta_m, \delta_d) = (0, 1)$ if $c_{m,i}(t) = M$ and $0 < c_{d,i}(t) < M$

• $(\delta_m, \delta_d) = (0, 0)$ if $c_{m,i}(t) = c_{d,i}(t) = M$.

Then, in the evaluation of the probability of the events $\mathcal{E}_1^*(t)$, $\mathcal{E}_2^\theta(t)$, and $\mathcal{E}_3^\theta(t)$ introduced in Section 5, one has to account for the probability of meeting a misbehaving node. For example, (11) can be rewritten as

$$ \mathbb{P}\{\mathcal{E}_1^*(t)\} = p_0 p^{01}(t) + p_1 p^{11}(t) + p_2. \quad (47) $$

The transition probabilities introduced in Sections 5.2.1 and 5.2.2 have to be updated accordingly. The form of the DTN state equations (28) remains the same.

Finally, the effect misbehaving nodes can be taken into account in (45) for the computation of the approximate expressions of $\tilde{p}^{00}$ and $\tilde{p}^{10}$. More specifically,

$$ \left\{ \begin{array}{l} \tilde{p}^{00} = \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} \left( \frac{p_0 q_{FA}(2) + p_2}{p_0 + p_2} \right)^{c_d} \left( 1 - \frac{p_0 q_{FA}(2) + p_2}{p_0 + p_2} \right)^{M - c_d}, \\ \tilde{p}^{10} = \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} \left( \frac{p_0 q_{D}(1,1) + p_2}{p_0 + p_2} \right)^{c_d} \left( 1 - \frac{p_0 q_{D}(1,1) + p_2}{p_0 + p_2} \right)^{M - c_d} \end{array} \right. \quad (49) $$

<figure><img src="image_7.png" /><figcaption>Fig. 4. Evolution of $X_0^{c_m, c_d}(t)$ (left) and $X_1^{c_m, c_d}(t)$ (right) obtained from (28), when $q_{FA}(0, 2) = 0.05$, $q_D(1, 1) = 0.8$, $q_D(0, 2) = 0.9$, $M = 4$, $\nu = 0.4$, and $\lambda = 1$.</figcaption></figure>

# 9 NUMERICAL RESULTS

In this Section we provide results aimed at assessing the convergence of the theoretical framework (Section 10.1), the appropriateness and accuracy of the framework also in case of specific mobility models such as the Brownian motion (Section 10.2) or other more realistic mobility models derived from user traces (Section 10.3), as well as to compare the proposed DFD methodology to other state-of-the-art solutions (Section 10.4). Finally, in Section 10.5 we investigate on the stability and accuracy of the Algorithm upon varying some key parameters.

## 9.1 Numerical verification of theoretical results

This section presents first the solution of the state equation (28) describing the evolution of the proportion of nodes in various states. Algorithm 2 is simulated considering a random displacement of nodes without any constraint on their speed. This allows to verify the correctness of the theoretical results presented in this paper.

Consider a LODT where $q_{FA}(0, 2) = 0.05$, $q_D(1, 1) = 0.8$, and $q_D(0, 2) = 0.9$. Besides $p_0 = 0.9$, $p_1 = 0.1$, $p_2 = 0$, $M = 4$, $\nu = 0.4$, and $\lambda = 1$. Figure 4 presents the evolution of the proportion of nodes with good sensors (left part) and defective sensors (right part) in different states, obtained

by solving (28). Note that $\Delta t$ represents the duration of a unit time slot used in the simulation. One observes that the proportion of nodes in each state becomes almost constant as $t/\Delta t \ge 15$. For the nodes with $\theta = 0$, only $X_0^{4,0}$ and $X_0^{4,1}$ are larger than 0.05, while the others are very close to 0. For the nodes with $\theta = 1$, only $X_1^{4,4}$, $X_1^{4,3}$, and $X_1^{4,2}$ are relatively large as compared to the other states. Consider the two sets of $c_d$s for which $X_0^{M, c_d}$ and $X_1^{M, c_d}$ are large. These sets have no common $c_d$ and it is thus easy to choose a decision threshold to distinguish both sets. The accuracy of the algorithm is then very good. With $\nu = 0.4$, one has $\bar{p}^{00} = 0.985$ and $\bar{p}^{10} = 0.027$. Only 1.5% of the good nodes believe they are carrying defective sensors. Less than 3% of the nodes with defective sensors have not been detected.

Consider now a set $S$ of $n_S = 1000$ moving nodes uniformly distributed over a square of unit area. In the first displacement model (jump motion model): node $i$ randomly chooses its location at time instant $(k+1)\Delta t$, independently from its previous location at time $k\Delta t$. Two nodes communicate only at discrete time instants $k\Delta t$ when their distance is less than $r_0$. Node $i$ has its neighbors in the set $N_i = \{j \in S: 0 < r_{i,j} \le r_0\}$, where $r_{i,j}$ is the distance between Nodes $i$ and $j$. Furthermore, if $|N_i| > 1$, we assume that node $i$ communicates only with its closest neighbor. Denote $\rho = \pi r_0^2 n_S$ as the average value of $|N_i|$. The cardinality of $N_i$ approximately follows a Poisson distribution as $n_S$ is large enough, the inter-contact probability is thus

$$ \lambda \Delta t = \mathbb{P}\{|N_i| = 1\} = \rho \exp(-\rho). $$

In the Monte-Carlo simulations, we set $r_0 = 0.014$, so that $\rho \approx 0.6$ and $\lambda \Delta t \approx 0.33$. Using the same values of $p_0$, $M$, $\nu$, $q_D$, and $q_{FA}$ as in Figure 4, the simulation results for this jump motion model are shown in Figure 5. Comparing Figure 4 and Figure 5, one remarks that the state evolution in the transient phase has similar shape but with different convergence speed, which depends mainly on $\lambda$. Figure 6 shows a good match between theory and simulation for the proportions of states at equilibrium. The approximation of $\bar{X}_\theta^{4, c_d}$ using (46) is also presented in Figure 6, which is very close to its actual value. Note that the difference between the approximated value and those obtained through theory and simulation is less than 0.1%.

## 9.2 Simulations with Brownian motion model

Consider now a Brownian motion model where each node is moving with a random speed. Each node changes its
$$ \left\{ \begin{aligned} \bar{p}^{00} &= \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} \left( \frac{p_0 q_{FA}(2) \bar{p}^{00} + p_1 q_D(1,1) \bar{p}^{10} + p_2}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10} + p_2} \right)^{c_d} \left( \frac{p_0 (1-q_{FA}(2)) \bar{p}^{00} + p_1 (1-q_D(1,1)) \bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10} + p_2} \right)^{M-c_d}, \\ \bar{p}^{10} &= \sum_{c_d: c_d/M < \nu} \binom{M}{c_d} \left( \frac{p_0 q_D(1,1) \bar{p}^{00} + p_1 q_D(0,2) \bar{p}^{10} + p_2}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10} + p_2} \right)^{c_d} \left( \frac{p_0 (1-q_D(1,1)) \bar{p}^{00} + p_1 (1-q_D(0,2)) \bar{p}^{10}}{p_0 \bar{p}^{00} + p_1 \bar{p}^{10} + p_2} \right)^{M-c_d}. \end{aligned} \right. \quad (48) $$

<figure><img src="image_8.png" /><figcaption>Fig. 5. Evolution of $X_{c_d}^{c_m, c_d}(t)$ (left) and $X_{c_d}^{c_m, c_d}(t)$ (right) by simulations with the jump model, when $q_{FA}(0,2) = 0.05$, $q_D(1,1) = 0.8$, $q_D(0,2) = 0.9$, $M = 4$, $\nu = 0.4$, and $\lambda\Delta t = 0.33$.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 7. Evolution of $p^{01}$ and $p^{10}$ for the two moving models, with $\sigma \in \{0.1, 1\}$, $q_{FA}(2) = 0.05$, $q_D(1,1) = 0.8$, $q_D(0,2) = 0.9$, $M = 10$ and $\nu = 0.4$.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 6. Comparison of $\bar{X}_{\theta}^{4, c_d}$ at the equilibrium.</figcaption></figure>

orientation when it reaches the boundary of the unit square. Let $\chi_i = (\chi_i^1, \chi_i^2)$ be the location of Agent $i$. Consider a first-order Brownian mobility model

$$ \begin{aligned} \chi_1^i((k+1)\Delta t) &= \chi_1^i(k\Delta t) + \eta_1^i(k\Delta t) \\ \chi_2^i((k+1)\Delta t) &= \chi_2^i(k\Delta t) + \eta_2^i(k\Delta t) \end{aligned} $$

where $\eta_1^i(k\Delta t) \sim \mathcal{N}(0, (\sigma r_0)^2)$ and $\eta_2^i(k\Delta t) \sim \mathcal{N}(0, (\sigma r_0)^2)$.

Consider $\sigma \in \{0.1, 1\}$, $q_{FA}(2) = 0.05$, $q_D(1,1) = 0.8$, $q_D(0,2) = 0.9$, $M = 10$, and $\nu = 0.4$. Figure 7 compares the evolution of $p^{01}$ and $p^{10}$ as functions of time for the jump motion model and the Brownian motion model, with fixed $\rho \approx 0.6$. At equilibrium, the performance obtained for both models is quite close. However, the convergence speed depends on the inter-contact rate $\lambda$. When $\sigma = 0.1$, the algorithm converges slowly in the Brownian motion model. When $\sigma = 1$, which results in a larger value of $\lambda$, the evolution of $p^{01}$ and $p^{10}$ with the Brownian motion model are close to the jump motion model.

At the beginning of the algorithm, each node believes that its sensors are good, thus $p^{01}(0) = 0$ and $p^{10}(0) = 1$. During the algorithm, $p^{10}(t)$ decreases in the transient phase until it reaches the equilibrium. Whereas, $p^{01}(t)$ increases at first and then decreases to the equilibrium. This comes from the fact that $p^{10}(t)$ is large at the beginning and the LODT performed on a good node often detects outliers.

## 9.3 Simulations with real databases

In this section, Algorithm 2 is executed using some experimental databases instead of motion models. These databases, provided by the MIT Reality Mining Project [43] and the Haggle Project [44], have been used in several previous works, e.g., [4]. In this work, we use the following databases:

* Reality, where $n_S = 97$, lasts more than 200 days with about 111 inter-contacts per day.

* Infocom05, where $n_S = 41$, lasts 3 days with approximately 312 inter-contacts every hour.

More specifically, one is interested in the inter-contact trace, i.e., which pair of nodes have a meeting at which time. The traces were taken from [45], which are converted from the original databases [43], [44].

Consider again the following parameters: $q_{FA}(2) = 0.05$, $q_D(1,1) = 0.8$, $q_D(0,2) = 0.9$, $M = 10$, and $\nu = 0.4$. Monte-Carlo simulations are performed 500 times for each database. In each test, $n_1$ nodes with random index are chosen to be defective. One sets $n_1 = 10$ in Infocom05 and $n_1 = 20$ in Reality. Two cases are considered. First, no misbehaving node is introduced, i.e., $n_2 = 0$. In a second case, $n_2 = 1$ in Infocom05 and $n_2 = 2$ in Reality. At the top of Figure 8, the index of the active nodes (which have contact with the others) are presented at each time to show the frequency of the inter-contacts at different epochs. The evolution of $p^{10}$ and $p^{01}$ is plotted at the bottom of Figure 8. Interestingly, in absence of misbehaving nodes, both $p^{10}$ and $p^{01}$, obtained with both databases, decrease to $10^{-3}$ after a sufficiently long time. One also observes that the convergence speed of $p^{10}$ and $p^{01}$ is highly related to the inter-contact rate (reflected by the density of points in the sub-figures at the top). Considering Infocom05, Table 2 further shows the influence of misbehaving nodes for various $n_2$. In presence of misbehaving nodes, the performance of the DFD algorithms worsens, but remains satisfying if the value of $\nu$ is properly chosen.

When $n_2 = 0$, Figure 9 represents the states at equilibrium $\bar{X}_{\theta}^{M, c_d}$ obtained with the databases Reality and Infocom05.
<table><caption>TABLE 2
Values of $p^{01}$ and $p^{10}$ at the end of simulations for different number of misbehaving nodes</caption><thead><tr><th>$n_2$</th><th>1</th><th>2</th><th>3</th><th>4</th><th>5</th><th>6</th></tr></thead><tbody><tr><th>$\nu$</th><td>0.5</td><td>0.5</td><td>0.5</td><td>0.5</td><td>0.6</td><td>0.6</td></tr><tr><th>$p^{01}(\%)$</th><td>0.3</td><td>1.4</td><td>2.7</td><td>8.0</td><td>3.4</td><td>7.2</td></tr><tr><th>$p^{10}(\%)$</th><td>0.7</td><td>0.5</td><td>0.3</td><td>0.3</td><td>1.3</td><td>1.4</td></tr></tbody></table>

<figure><img src="image_11.png" /><figcaption>Fig. 9. Comparison of $\bar{X}_{\theta}^{10,c_d}$ at the equilibrium obtained using the Reality database, the Infocom05 database, and predicted by the approximation (46).</figcaption></figure>

com05, and also considering the approximation (46). There is an excellent match between the values at equilibrium predicted by theory and those obtained in practice. Note that the difference between the approximated value and those obtained through theory and simulation is less than 0.05%.

## 9.4 Comparison with state-of-the-art solutions

This section compares the proposed DFD algorithm to some closely related scheme in the literature. As mentioned in Section 2, classical DFD algorithms are difficult to apply in the context of DTN and no solutions have been presented so far in the literature for this specific scenario. Accordingly, in order to perform a meaningful comparison between our algorithm and a state-of-the-art approach, we have considered the gossip algorithm discussed in [18] which represents the most robust and efficient methodology in the context of classification and distributed estimation in dynamic scenarios like DTNs.

In [18], $n_S$ nodes are assumed to get a measurement

$$m_i = c + \theta_i + v_i, \forall i \in S, \quad (50)$$

of a common quantity $c$, where $v_i$ are realizations of independent zero-mean Gaussian random variables with variance $\sigma^2$ and $\theta_i \in \{0, 1\}$ denotes the bias of each node. Each node is interested in the joint estimation of $c$ and $\theta_i$. Since the measurements produced by the sensors with non-zero bias are more likely to have larger values, [18] proposes an estimator of $\theta_i$ based on a distributed ranking of the nodes according to their measurement $m_i$. Nodes with a large rank get an estimate $\hat{\theta}_i = 1$, while nodes with a small rank have $\theta_i = 0$.

In order to apply the proposed DFD algorithm to this problem, consider the following LODT

$$y_{i,j} = y_{j,i} = \begin{cases} 1 & \text{if } |m_i - m_j| > \delta, \\ 0 & \text{else} \end{cases},$$

<figure><img src="image_12.png" /><figcaption>Fig. 10. Comparison of the DFD part of the estimation algorithm proposed in [18] with the proposed DFD algorithm, when $\sigma = 0.2$ (top) and $\sigma = 0.3$ (bottom)</figcaption></figure>

where $\delta$ is a threshold that results in different values of probabilities $q_{FA}(2)$, $q_D(1, 1)$, and $q_D(0, 2)$; then the proper value of $\nu$ can be set accordingly.

Consider again Infocom05 for the simulation with $n_1 = 10$ nodes chosen randomly with $\theta = 1$ and without misbehaving node. Two scenarios are considered. In the first case, all nodes take a single measurement of $c$ at initialization. In the second case, nodes take measurements at each meeting. Results are obtained as the average of 200 independent Monte-Carlo simulations. Figure 10 compares the results when $\sigma = 0.2$ and $\sigma = 0.3$. The classification error and the estimation error are defined as $E_c = \sum |\theta_i - \hat{\theta}_i|/n_S$ and $E_e = \sum |c - \hat{c}_i|/n_S$. If nodes take a single measurement, the performance of the proposed algorithm is close to the reference method in terms of $E_c$ and $E_e$. When nodes take new measurements at each meeting, the proposed DFD algorithm performs better than the reference method: the value of $E_c$ decreases faster and turns to be much smaller. This is mainly due to the node ranking algorithm used in [18], which becomes less efficient when nodes update at each meeting the quantity according to which they are ranked.

## 9.5 Influence of the parameters

This section characterizes the influence of the parameters, such as $p_1$, $q_D(1, 1)$, and $M$, on the performance of Algorithm 2. The jump motion model is used throughout this section to describe the displacement of the nodes.

Consider fixed $q_{FA}(2) = 0.05$, $q_D(1, 1) = 0.8$, $q_D(0, 2) = 0.9$, and $p_2 = 0$. The evolution of $p^{10}$ and $p^{01}$ for various $p_1 \in \{0.1, 0.5\}$ and $M \in \{4, 10, 20\}$ is shown in Figure 11.
<figure><img src="image_13.png" /><figcaption>Fig. 8. Indexes of active nodes (having met another node) at different time (top) and evolution of $p^{01}$ and $p^{10}$ (bottom) obtained by using the Reality database (left) and the Infocom05 database (right), with $q_{FA}(2) = 0.05$, $q_D(1, 1) = 0.8$, $q_D(0, 2) = 0.9$, $M = 10$, and $\nu = 0.4$.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Fig. 11. Evolution of $p^{10}$ and $p^{01}$ for various $M \in \{4, 10, 20\}$ and $p_1 \in \{0.1, 0.5\}$, with $q_{FA}(2) = 0.05$, $q_D(1, 1) = 0.8$, $q_D(0, 2) = 0.9$.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Fig. 12. Achievable $\bar{p}^{10}$ or $\bar{p}^{01}$ for different values of the pair $(q_D(1, 1), q_{FA}(2))$ when $p_1 = 0.1$ (left) and for $p_1 = 0.5$ (right) and $p_2 = 0$.</figcaption></figure>

For each case, the value of $\nu$ is chosen such that it minimizes $\tilde{p}^{01} + \tilde{p}^{10}$. One observes that a large $M$ leads to a better performance at equilibrium. The price to be paid is a longer time required to reach equilibrium. When $M = 10$, both $\bar{p}^{10}$ and $\bar{p}^{01}$ are around $10^{-3}$. The proportion of the nodes with defective sensors has also an impact on the convergence speed of the algorithm. For example, when $p_1$ is large, more time is needed to achieve a given level of performance (in terms of $\bar{p}^{10}$ and $\bar{p}^{01}$).

To show the effectiveness of the proposed DFD algorithm, consider now $q_D(0, 2) = 0.9$ and $M = 10$. For $p_1 = 0.1$ and $p_1 = 0.5$, one is interested in the achievable $\bar{p}^{10}$ and $\bar{p}^{01}$ for $0 \le q_{FA}(2) < q_D(0, 2)$ and $q_{FA}(2) < q_D(1, 1) \le q_D(0, 2)$. Four areas are considered:

* Area 3: both $\bar{p}^{10}$ and $\bar{p}^{01}$ are less than $10^{-3}$;

* Area 2: both $\bar{p}^{10}$ and $\bar{p}^{01}$ are less than $10^{-2}$;

* Area 1: both $\bar{p}^{10}$ and $\bar{p}^{01}$ are less than $10^{-1}$;

* Area 0: either $\bar{p}^{10}$ or $\bar{p}^{01}$ cannot be less than $10^{-1}$.

Figure 12 shows partition of the $(q_D(1, 1), q_{FA}(2))$ triangle in four areas, represented in different colors. The ratio of defective nodes in the network has not a significant impact

on the performance at the equilibrium, even when 50% of nodes are defective.

This assesses the robustness of the approach.

# 10 CONCLUSION

This paper presents a fully distributed algorithm allowing each node of a DTN to estimate the status of its own sensors using LODT performed during the meeting of nodes. The DFD algorithm is analyzed considering a Markov model of the evolution of the proportion of nodes with a given belief in their status. This model is then used to derive a system of ordinary differential equations approximating the evolution of the proportions of the nodes in different states. The existence and uniqueness of an equilibrium is discussed. Interestingly, the proportions at the equilibrium follow a binomial distribution. The approximations of these proportions of nodes at equilibrium provide insight to properly choose the decision parameter of the DFD algorithm. In the simulations, a jump motion model, a Brownian motion model, as well as databases containing traces of inter-contact time instants are considered. The results show a good match
with theory. The convergence speed of the DFD algorithm depends on the inter-contact rate and on the proportion of nodes with defective sensors $p_1$. Nevertheless, $p_1$ has not a significant impact on the non-detection and false alarm rates at equilibrium, showing the robustness of the approach also in case of a large number of defective nodes. The impact of the presence of misbehaving nodes has also been considered, showing the robustness of the proposed DFD algorithm.

## REFERENCES

[1] M. J. Khabbaz, C. M. Assi, and W. F. Fawaz, "Disruption-tolerant networking: A comprehensive survey on recent developments and persisting challenges," *IEEE Com. Surv. & Tut.*, vol. 14, no. 2, pp. 607-640, 2012.

[2] P. R. Pereira, A. Casaca, J. J. Rodrigues, V. N. Soares, J. Triay, and C. Cervelló-Pastor, "From delay-tolerant networks to vehicular delay-tolerant networks," *IEEE Com. Surv. & Tut.*, vol. 14, no. 4, pp. 1166-1182, 2012.

[3] K. Wei, M. Dong, J. Weng, G. Shi, K. Ota, and K. Xu, "Congestion-aware message forwarding in delay tolerant networks: a community perspective," *Concurrency and Computation: Practice and Experience*, vol. 27, no. 18, pp. 5722-5734, 2015.

[4] P. Hui, J. Crowcroft, and E. Yoneki, "Bubble rap: Social-based forwarding in delay-tolerant networks," *IEEE Trans. Mob. Comp.*, vol. 10, no. 11, pp. 1576-1589, 2011.

[5] V. N. Soares, J. J. Rodrigues, and F. Farahmand, "Geospray: A geographic routing protocol for vehicular delay-tolerant networks," *Information Fusion*, vol. 15, pp. 102-113, 2014.

[6] H. Zhu, S. Du, Z. Gao, M. Dong, and Z. Cao, "A probabilistic misbehavior detection scheme toward efficient trust establishment in delay-tolerant networks," *IEEE Trans. Paral. and Dist. Syst.*, vol. 25, no. 1, pp. 22-32, 2014.

[7] L. Galluccio, B. Lorenzo, and S. Glisic, "Sociality-aided new adaptive infection recovery schemes for multicast dtns," *IEEE Trans. Veh. Tech.*, vol. 65, no. 5, pp. 3360-3376, 2016.

[8] M. Panda, A. Ali, T. Chahed, and E. Altman, "Tracking message spread in mobile delay tolerant networks," *IEEE Trans. Mob. Comp.*, vol. 14, no. 8, pp. 1737-1750, 2015.

[9] W. Li, F. Bassi, D. Dardari, M. Kieffer, and G. Pasolini, "Defective sensor identification for WSNs involving generic local outlier detection tests," *IEEE Trans. Sig. and Inf. Proc. over Netw.*, vol. 2, no. 1, pp. 29-48, 2016.

[10] J. Chen, S. Kher, and A. Somani, "Distributed fault detection of wireless sensor networks," in *Proc. DIWANS*, New York, NY, 2006, pp. 65-72.

[11] J.-L. Gao, Y.-J. Xu, and X.-W. Li, "Weighted-median based distributed fault detection for wireless sensor networks," *Jnl of Softw.*, vol. 18, no. 5, pp. 1208-1217, 2007.

[12] S. Ji, S.-F. Yuan, T.-H. Ma, and C. Tan, "Distributed fault detection for wireless sensor based on weighted average," in *Proc. NSWCTC*, Wuhan, China, 2010, pp. 57-60.

[13] M. Panda and P. Khilar, "Distributed self fault diagnosis algorithm for large scale wireless sensor networks using modified three sigma edit test," *Ad Hoc Networks*, vol. 25, pp. 170-184, 2015.

[14] Y. Zhang, N. Meratnia, and P. Havinga, "Outlier detection techniques for wireless sensor networks: A survey," *IEEE Com. Surv. & Tut.*, vol. 12, no. 2, pp. 159-170, 2010.

[15] A. Mahapatro and P. M. Khilar, "Fault diagnosis in wireless sensor networks: A survey," *IEEE Com. Surv. & Tut.*, vol. 15, no. 4, pp. 2000-2026, 2013.

[16] H. Dong, Z. Wang, S. X. Ding, and H. Gao, "A survey on distributed filtering and fault detection for sensor networks," *Math. Prob. in Eng.*, 2014.

[17] E. F. Nakamura, A. A. Loureiro, and A. C. Frery, "Information fusion for wireless sensor networks: Methods, models, and classifications," *ACM Comp. Surv.*, vol. 39, no. 3, p. 9, 2007.

[18] A. Chiuso, F. Fagnani, L. Schenato, and S. Zampieri, "Gossip algorithms for simultaneous distributed estimation and classification in sensor networks," *IEEE Jnl Sel. Top. Sig. Proc.*, vol. 5, no. 4, pp. 691-706, 2011.

[19] F. Fagnani, S. M. Fosson, and C. Ravazzi, "A distributed classification/estimation algorithm for sensor networks," *SIAM Jnl Contr. and Opt.*, vol. 52, no. 1, pp. 189-218, 2014.

[20] ——, "Consensus-like algorithms for estimation of gaussian mixtures over large scale networks," *Math. Mod. and Meth. in App. Sciences*, vol. 24, no. 02, pp. 381-404, 2014.

[21] B. Zhu, W. Zhang, W. Feng, and L. Zhang, "Distributed faulty node detection and isolation in delay-tolerant vehicular sensor networks," in *Proc. PIMRC*, Sept 2012, pp. 1497-1502.

[22] W. Peng, F. Li, X. Zou, and J. Wu, "Behavioral malware detection in delay tolerant networks," *IEEE Trans. on Paral. and Dist. Syst.*, vol. 25, no. 1, pp. 53-63, 2014.

[23] E. Hernandez-Orallo, M. D. Serrat Olmos, J.-C. Cano, C. T. Calafate, and P. Manzoni, "Cocowa: A collaborative contact-based watchdog for detecting selfish nodes," *IEEE Trans. Mob. Comp.*, vol. 14, no. 6, pp. 1162-1175, 2015.

[24] E. Ayday and F. Fekri, "An iterative algorithm for trust management and adversary detection for delay-tolerant networks," *IEEE Trans. Mob. Comp.*, vol. 11, no. 9, pp. 1514-1531, 2012.

[25] Y. Liu, D. R. Bild, R. P. Dick, Z. M. Mao, and D. S. Wallach, "The mason test: A defense against sybil attacks in wireless networks without trusted authorities," *IEEE Trans. Mob. Comp.*, vol. 14, no. 11, pp. 2376-2391, 2015.

[26] J. R. Douceur, "The sybil attack," in *Peer-to-peer Systems*. Springer, 2002, pp. 251-260.

[27] G. Han, J. Jiang, L. Shu, and M. Guizani, "An attack-resistant trust model based on multidimensional trust metrics in underwater acoustic sensor network," *IEEE Trans. Mob. Comp.*, vol. 14, no. 12, pp. 2447-2459, 2015.

[28] M. Abdelhakim, L. E. Lightfoot, J. Ren, and T. Li, "Distributed detection in mobile access wireless sensor networks under byzantine attacks," *IEEE Trans. on Paral. and Dist. Syst.*, vol. 25, no. 4, pp. 950-959, 2014.

[29] H. Zhu, L. Fu, G. Xue, Y. Zhu, M. Li, and L. Ni, "Recognizing exponential inter-contact time in vanets," in *Proc. INFOCOM*, March 2010, pp. 1-5.

[30] J. P. Norton, Ed., Special Issue on Bounded-Error Estimation: Issue 1, 1994, *Int. Jnl Adapt. Cont. Sig. Proc.* 8(1):1-118.

[31] ——, Special Issue on Bounded-Error Estimation: Issue 2, 1995, *Int. Jnl Adapt. Cont. Sig. Proc.* 9(1):1-132.

[32] M. Milanese, J. Norton, H. Piet-Lahanier, and E. Walter, Eds., *Bounding Approaches to System Identification*. New York, NY: Plenum Press, 1996.

[33] A. Shwartz and A. Weiss, "Large deviations for performance analysis," 1995.

[34] T. G. Kurtz, "The central limit theorem for markov chains," *The Annals of Prob.*, pp. 557-560, 1981.

[35] W. O. Kermack and A. G. McKendrick, "A contribution to the mathematical theory of epidemics," in *Proc. Royal Soc. London, Ser. A*, vol. 115, no. 772, 1927, pp. 700-721.

[36] Z. J. Haas and T. Small, "A new networking model for biological applications of ad hoc sensor networks," *IEEE/ACM Trans. Netw.*, vol. 14, no. 1, pp. 27-40, 2006.

[37] E. Renshaw, *Modelling biological populations in space and time*. Cambridge University Press, 1993, vol. 11.

[38] M. E. Newman, "Spread of epidemic disease on networks," *Phys. Rev. E*, vol. 66, no. 1, July 2002.

[39] T. Spyropoulos, T. Turletti, and K. Obraczka, "Routing in delay-tolerant networks comprising heterogeneous node populations," *IEEE Trans. Mob. Comp.*, vol. 8, no. 8, pp. 1132-1147, 2009.

[40] Y. Lin, B. Li, and B. Liang, "Stochastic analysis of network coding in epidemic routing," *IEEE Jnl Sel. Ar. Com.*, vol. 26, no. 5, 2008.

[41] A. Granas and J. Dugundji, *Fixed point theory*. Springer, 2013.

[42] S. Banach, "Sur les opérations dans les ensembles abstraits et leur application aux équations intégrales," *Fund. Math.*, vol. 3, no. 1, pp. 133-181, 1922.

[43] N. Eagle and A. Pentland, "Reality mining: sensing complex social systems," *Pers. Ubiq. Comp.*, vol. 10, no. 4, pp. 255-268, 2006.

[44] J. Scott, R. Gass, J. Crowcroft, P. Hui, C. Diot, and A. Chaintreau, "CRAWDAD dataset cambridge/haggle (v. 2009-05-29)," Downloaded from http://crawdad.org/cambridge/haggle/20090529, May 2009.

[45] M. Orlinski, "Encounter traces for the one simulator," Downloaded from http://www.shigs.co.uk/index.php?page=traces.

[46] T. M. Cover and J. A. Thomas, *Elements of information theory*. John Wiley & Sons, 2012.
<figure><img src="image_16.png" /></figure>

information theory.

**Wenjie Li** received the B.S. degree from the Huazhong University of Science and Technology, China, in 2011, and the M.S. degree in Advanced Systems of Wireless Communication from the University of Paris-Sud, France, in 2013. He is currently pursuing the Ph.D. degree at the Laboratoire des Signaux et Systemes, France. He has been visiting Ph.D. student at the University of Bologna, Italy, in 2014 and 2015. His research interests lie in the areas of distributed signal processing over networks and

<figure><img src="image_17.png" /></figure>

**Laura Galluccio** received her laurea degree in Electrical Engineering from University of Catania, Catania, Italy, in 2001. In March 2005 she got her Ph.D. in Electrical, Computer and Telecommunications Engineering at the same university. Since 2002 she is also at the Italian National Consortium of Telecommunications (CNIT), where she worked as a Research Fellow within the VICOM (Virtual Immersive Communications) and the SATNEX Projects. Since November 2010 she is Assistant Professor at University of Catania. Her research interests include ad hoc and sensor networks, protocols and algorithms for wireless networks, and network performance analysis. From May to July 2005 she has been a Visiting Scholar at the COMET Group, Columbia University, NY. In September 2015 she has been Visiting Professor at CentraleSupelec, Gif-sur-Yvette, France. She is member of the IEEE, Comsoc and ACM N2Women Group.

<figure><img src="image_18.png" /></figure>

**Francesca Bassi** received the Master of Science degree in telecommunications engineering from the University of Bologna, Bologna, Italy, in 2006, and the Ph.D. degree from University Paris-Sud, Orsay, France, and from the University of Bologna, in 2010. In 2007 and 2009, she has been a visiting Ph.D. student at the University of Bologna. From 2006 to 2010, she was with Laboratoire des Signaux et des Systemes, Gif-sur-Yvette, France. Between 2011 and 2012, she was a Postdoctoral Researcher at the University of Michigan, Ann Arbor, MI, USA. Since 2012, she has been an Assistant Professor at ESME-Sudria, Ivry-sur-Seine, France, and an Associate Researcher with the Laboratoire des Signaux et des Systemes.

<figure><img src="image_19.png" /></figure>

**Michel Kieffer** is a full professor in signal processing for communications at the Paris-Sud University and a researcher at the Laboratoire des Signaux et Systemes, Gif-sur-Yvette. Since 2009, he is also part-time invited professor at the Laboratoire Traitement et Communication de l'Information, Telecom ParisTech, Paris. His research interests are in signal processing for multimedia, communications, and networking, distributed source coding, network coding, joint source-channel coding and decoding, joint source-network coding. Michel Kieffer is co-author of more than 150 contributions in journals, conference proceedings, or books. He is one of the co-authors of the book *Applied Interval Analysis* published by Springer-Verlag in 2001, and of the book *Joint source-channel decoding: A crosslayer perspective with applications in video broadcasting* published by Academic Press in 2009. He serves as associate editor of Signal Processing since 2008 and of the IEEE Transactions on Communications from 2012 to 2016. From 2011 to 2016, Michel Kieffer was junior member of the Institut Universitaire de France.
<figure><img src="image_20.png" /><figcaption>Fig. 13. Transient regime: Possible state transitions from and to state $(\theta, c_m, c_d)$ when $0 < c_m < M$ and $0 < c_d < c_m$</figcaption></figure>

# APPENDIX A
## PROOF OF PROPOSITION 2

At time $t$, remind that $\mathcal{X}_{\theta}^{c_m, c_d}(t)$ is the proportion of nodes in state $(\theta, c_m, c_d)$. For $\theta \in \{0, 1\}$, the process

$$ \mathcal{X}_{\theta}(t) = (\mathcal{X}_{\theta}^{0,0}(t), \mathcal{X}_{\theta}^{1,0}(t), \mathcal{X}_{\theta}^{1,1}(t), \dots, \mathcal{X}_{\theta}^{M,0}(t), \dots, \mathcal{X}_{\theta}^{M,M}(t)) $$

is a jump Markov process with random jump time instants and with a jump size $1/n_{\theta}$. In order to get the expected proportions of nodes $X_{\theta}^{c_m, c_d}(t)$, one will consider an inter-contact rate $\lambda$ and a well-mixed population of nodes. During a short time interval $[t, t + \delta t]$ the number of nodes with state $(\theta, c_m, c_d)$ that will meet another node is $\lambda p_{\theta} n_S \mathcal{X}_{\theta}^{c_m, c_d}(t) \delta t$.

When $0 \le c_m < M$, and thus also $c_d \le c_m < M$, nodes with state $(\theta, c_m, c_d)$ will switch to the states $(\theta, c_m + \delta_m, c_d + \delta_d)$, with $(\delta_m, \delta_d) \in \{(0,0), (1,0), (1,1)\}$ with a probability $\pi_{\theta}^{\delta_m, \delta_d}(t, c_m, c_d)$, see Figure 13.

As a consequence, at time $t + \delta t$, the number of nodes in State $(\theta, c_m, c_d)$ may be expressed as follows

$$ \begin{aligned} p_{\theta} n_S \mathcal{X}_{\theta}^{c_m, c_d}(t + \delta t) = & p_{\theta} n_S \mathcal{X}_{\theta}^{c_m, c_d}(t) \\ & + \lambda \delta t p_{\theta} n_S (-\mathcal{X}_{\theta}^{c_m, c_d}(t) (\pi_{\theta}^{1,0}(t, c_m, c_d) + \pi_{\theta}^{1,1}(t, c_m, c_d)) \\ & + \mathcal{X}_{\theta}^{c_m-1, c_d-1}(t) \pi_{\theta}^{1,1}(t, c_m-1, c_d-1) + \mathcal{X}_{\theta}^{c_m-1, c_d}(t) \pi_{\theta}^{1,0}(t, c_m-1, c_d)). \end{aligned} \quad (51) $$

The evolution of the expected value $X_{\theta}^{c_m, c_d}(t)$ of $\mathcal{X}_{\theta}^{c_m, c_d}(t)$ is then described by the following differential equation²

$$ \frac{dX_{\theta}^{c_m, c_d}}{dt} = -\lambda X_{\theta}^{c_m, c_d} (\pi_{\theta}^{1,0} (c_m, c_d) + \pi_{\theta}^{1,1} (c_m, c_d)) \\ + \lambda X_{\theta}^{c_m-1, c_d-1} \pi_{\theta}^{1,1} (c_m-1, c_d-1) + \lambda X_{\theta}^{c_m-1, c_d} \pi_{\theta}^{1,0} (c_m-1, c_d). \quad (52) $$

When $c_m = M$ and $0 < c_d < M$, nodes in state $(\theta, M, c_d)$ will switch to the states $(\theta, M, c_d + \delta_d)$, $\delta_d \in \{-1, 0, 1\}$ with a probability $\pi_{\theta}^{0, \delta_d}(t, M, c_d)$. Nodes in the states $(\theta, M-1, c_d-1)$ and $(\theta, M-1, c_d)$ that have met an other node in the time interval $[t, t + \delta t]$ may reach state $(\theta, M, c_d)$, respectively with a probability $\pi_{\theta}^{1,1}(t, M-1, c_d-1)$ and $\pi_{\theta}^{1,0}(t, M-1, c_d)$, see Figure 14. As a consequence, the evolution of the expected value $X_{\theta}^{M, c_d}(t)$ of $\mathcal{X}_{\theta}^{M, c_d}(t)$ can be described by

$$ \frac{dX_{\theta}^{M, c_d}}{dt} = -\lambda X_{\theta}^{M, c_d} (\pi_{\theta}^{0,1}(M, c_d) + \pi_{\theta}^{0,-1}(M, c_d)) \\ + \lambda X_{\theta}^{M-1, c_d-1} \pi_{\theta}^{1,1}(M-1, c_d-1) + \lambda X_{\theta}^{M-1, c_d} \pi_{\theta}^{1,0}(M-1, c_d) \\ + \lambda X_{\theta}^{M, c_d-1} \pi_{\theta}^{0,1}(M, c_d-1) + \lambda X_{\theta}^{M, c_d+1} \pi_{\theta}^{0,-1}(M, c_d+1). \quad (53) $$

2. Notice that to lighten notations, time dependency is omitted. Moreover, $\pi_{\theta}^{\delta_m, \delta_d}(c_m, c_d)$ and $p^{\theta\theta}$ in the rest of the paper represent the expected values, as they can be represented as functions of $\mathcal{X}_{\theta}^{c_m, c_d}(t)$ or $X_{\theta}^{c_m, c_d}(t)$.

<figure><img src="image_21.png" /><figcaption>Fig. 14. Permanent regime: Possible state transitions from and to State $(\theta, M, c_d)$ when $0 < c_d < M$</figcaption></figure>

Similar derivations can be made for the remaining DTN state components to obtain (28).

# APPENDIX B
## PROOF OF PROPOSITION 3

For the proof, one considers first the following lemmas.

**Lemma 10.** If

$$ \lim_{t \to \infty} \int_0^t (p_0 p^{00}(\tau) + p_1 p^{10}(\tau)) d\tau = \infty \quad (54) $$

then $p_0 p^{00}(t) + p_1 p^{10}(t) > 0$ for all $t \in \mathbb{R}^+$.

*Proof:* Since $p_0 > 0$, $p_1 > 0$, $p^{00} \ge 0$, and $p^{10} \ge 0$, it suffices to prove that

$$ p^{00}(t) + p^{10}(t) \neq 0 \quad \forall t > 0. \quad (55) $$

Assume that there exists a time instant $t^* > 0$, such that $p^{00}(t^*)+p^{10}(t^*) = 0$. As a consequence, at time $t^*$, all nodes in the network believe themselves as carrying defective sensors. As a consequence, no node will transmit its data to its neighbors. No LODTs will be performed after time $t^*$ and the state of nodes will remain constant. Hence, if $p^{00}(t^*) + p^{10}(t^*) = 0$ for some $t^*$, then $p^{00}(t) + p^{10}(t) = 0$ for all $t > t^*$. Consequently,

$$ \lim_{t \to \infty} \int_0^t (p_0 p^{00}(\tau) + p_1 p^{10}(\tau)) d\tau = \int_0^{t^*} (p_0 p^{00}(\tau) + p_1 p^{10}(\tau)) d\tau, $$

which contradicts (54). □

**Lemma 11.** The property (54) is always satisfied.

*Proof:* From (28-a), one has

$$ X_{\theta}^{0,0}(t) = \exp \left( -\lambda \int_{0}^{t} (p_{0}p^{00}(\tau) + p_{1}p^{10}(\tau)) d\tau \right). \quad (56) $$

Assume that there exists $C^* > 0$ such that

$$ \lim_{t \to \infty} \int_0^t (p_0 p^{00}(\tau) + p_1 p^{10}(\tau)) d\tau \le C^* \quad (57) $$

then $\forall t \ge 0$, one has

$$ \int_0^t (p_0 p^{00}(\tau) + p_1 p^{10}(\tau)) d\tau \le C^*. \quad (58) $$

Combining (56) and (58), one gets

$$ X_{\theta}^{0,0}(t) \ge \exp(-\lambda C^*) > 0. \quad (59) $$

Moreover, from (6), one has $p^{00}(\tau) \ge X_{\theta}^{0,0}(\tau)$, leading to

$$ \int_0^t (p_0 p^{00}(\tau) + p_1 p^{10}(\tau)) d\tau \ge \int_0^t p_0 X_{\theta}^{0,0}(\tau) d\tau > p_0 \exp(-\lambda C^*) t. \quad (60) $$

Since $\exp(-\lambda C^*) t \to \infty$ as $t \to \infty$, (60) leads to a violation of the hypothesis (57). Hence, one always has (54). □
The proof of Proposition 3 is then by induction. Starting with (28-a), one has (56). Since (54) is satisfied according to Lemma 11, for any $\xi > 0$, there exists $t_{00} > 0$ such that $t > t_{00}$ implies $X_{\theta}^{0,0}(t) < \xi$ and $\lim_{t \to \infty} X_{\theta}^{0,0}(t) = 0$.

Then, assume that for any $c_m \le M-1$, and $\xi > 0$, there exists $t_{(c_m-1)0} > \dots > t_{00}$ such that $t > t_{(c_m-1)0}$ implies $X_{\theta}^{j,0}(t) < \xi$ for $j = 0, \dots, c_m-1$. One has to show now that there exists $t_{c_m0} > t_{(c_m-1)0}$ such that $X_{\theta}^{c_m0}(t) < \xi$ for all $t > t_{c_m0}$.

Define $Z_{\theta}^{c_m,0}(t) = \sum_{j=0}^{c_m} X_{\theta}^{j,0}(t)$. From (28a) and (28b), one has

$$ \frac{dZ_{\theta}^{c_m,0}}{dt} = -\lambda \left( v(t) Z_{\theta}^{c_m-1,0}(t) + (p_0 p^{00}(t) + p_1 p^{10}(t)) X_{\theta}^{c_m,0}(t) \right), $$

where $v(t) = \pi_{\theta}^{1,1}(t, c_m, c_d)$, since $\pi_{\theta}^{1,0}$ and $\pi_{\theta}^{1,1}$ do not depend on $c_m$ and $c_d$ when $c_m < M$. Using (55) one has $dZ_{\theta}^{c_m,0}/dt < 0$ for any $X_{\theta}^{c_m0} > 0$. As a consequence, $Z_{\theta}^{c_m,0}(t)$ decreases until $X_{\theta}^{c_m0}(t)$ reaches 0. Hence, for any $\xi > 0$, there exists $t_{c_m0} > t_{(c_m-1)0}$, such that $X_{\theta}^{c_m0} < \xi$ and $\lim_{t \to \infty} X_{\theta}^{c_m0}(t) = 0$.

In the same way, using (28c) and the previous results that $X_{\theta}^{c_m c_d}(t) \to 0$ with $c_d = 1, \dots, M-2$ and $c_m = c_d, \dots, M-2$, one can prove that for any $c_d = 1, \dots, M-1$, $X_{\theta}^{c'_m, (c_d+1)}(t)$ tends to zero as $t \to \infty$, with any $c'_m = c_d+1, \dots, M-1$.

# APPENDIX C

## PROOF OF PROPOSITION 4

According to Proposition 3, one has $\bar{X}_{\theta}^{c_m, c_d} = 0$, for all $c_m < M$ and $c_d \le c_m$. To evaluate $\bar{X}_{\theta}^{M, c_d}$, one thus considers the following simplified dynamics derived from (28) for $\theta \in \{0, 1\}$,

$$ \left\{ \begin{aligned} \frac{d\bar{X}_{\theta}^{M,0}}{dt} &= \lambda \left( -X_{\theta}^{M,0} \pi_{\theta}^{0,1}(M,0) + X_{\theta}^{M,1} \pi_{\theta}^{0,-1}(M,1) \right), \\ \frac{d\bar{X}_{\theta}^{M,M}}{dt} &= \lambda \left( -X_{\theta}^{M,M} \pi_{\theta}^{0,-1}(M,M) + X_{\theta}^{M,M-1} \pi_{\theta}^{0,1}(M,M-1) \right), \\ \frac{d\bar{X}_{\theta}^{M,c_d}}{dt} &= \lambda \left( -X_{\theta}^{M,c_d} \left( \pi_{\theta}^{0,-1}(M,c_d) + \pi_{\theta}^{0,1}(M,c_d) \right) + X_{\theta}^{M,c_d+1} \pi_{\theta}^{0,-1}(M,c_d+1) + X_{\theta}^{M,c_d-1} \pi_{\theta}^{0,1}(M,c_d-1) \right). \end{aligned} \right. \quad (61) $$

At equilibrium, one has $d\bar{X}_{\theta}^{M, c_d}(t)/dt = 0$ for all $c_d \le M$. Moreover, the transition probabilities do not vary any more.

Let $\bar{X}_{\theta}^{M} = (\bar{X}_{\theta}^{M,1}, \dots, \bar{X}_{\theta}^{M,M})^T$, $a_{\theta}(c_d) = \pi_{\theta}^{0,1}(M, c_d)$, and $b_{\theta}(c_d) = \pi_{\theta}^{0,-1}(M, c_d)$. From (61), one deduces that the vector $\bar{X}_{\theta}^{M}$ should satisfy $\Psi_{\theta} \cdot \bar{X}_{\theta}^{M} = 0$ where

$$ \Psi_{\theta} = \begin{pmatrix} -a_{\theta}(0) & b_{\theta}(1) & & \\ a_{\theta}(0) & -a_{\theta}(1) - b_{\theta}(1) & b_{\theta}(2) & \\ & \ddots & \ddots & \\ & & a_{\theta}(M-1) & -b_{\theta}(M) \end{pmatrix}. $$

Summing Lines 1 to $c_d+1$, for all $c_d = 0, \dots, M-1$, one obtains $a_{\theta}(c_d) \bar{X}_{\theta}^{M, c_d} = b_{\theta}(c_d+1) \bar{X}_{\theta}^{M, c_d+1}$, which leads to

$$ \bar{X}_{\theta}^{M, c_d} = \bar{X}_{\theta}^{M,0} \prod_{j=0}^{c_d-1} \frac{a_0(j)}{b_0(j+1)}. \quad (62) $$

One evaluates

$$ \frac{a_{\theta}(j)}{b_{\theta}(j+1)} = \frac{\pi_{\theta}^{0,1}(M, j)}{\pi_{\theta}^{0,-1}(M, j+1)} = \eta_{\theta} \frac{M-j}{j+1}, \quad (63) $$

where using (25) and (26), one has

$$ \left\{ \begin{aligned} \eta_0 &= \frac{p_0 q_{\mathrm{FA}}(2) \bar{p}^{00} + p_1 q_D(1,1) \bar{p}^{10}}{p_0(1-q_{\mathrm{FA}}(2)) \bar{p}^{00} + p_1(1-q_D(1,1)) \bar{p}^{10}}, \\ \eta_1 &= \frac{p_0 q_D(1,1) \bar{p}^{00} + p_1 q_D(0,2) \bar{p}^{10}}{p_0(1-q_D(1,1)) \bar{p}^{00} + p_1(1-q_D(0,2)) \bar{p}^{10}}. \end{aligned} \right. \quad (64) $$

with $\bar{p}^{00}$ and $\bar{p}^{10}$ defined in (29).

From (62) and (63), one deduces

$$ \begin{aligned} \bar{X}_{\theta}^{M, c_d} &= \bar{X}_{\theta}^{M,0} \prod_{j=0}^{c_d-1} \left( \eta_{\theta} \frac{M-j}{j+1} \right) \\ &= \bar{X}_{\theta}^{M,0} \eta_{\theta}^{c_d} \frac{M \cdot (M-1) \cdot (M-c_d+1)}{1 \cdot 2 \cdots c_d} = \binom{M}{c_d} \eta_{\theta}^{c_d} \bar{X}_{\theta}^{M,0}. \end{aligned} \quad (65) $$

Since $\sum_{c_d=0}^{M} \bar{X}_{\theta}^{M, c_d} = 1$, one has

$$ 1 = \sum_{c_d=0}^{M} \binom{M}{c_d} \eta_{\theta}^{c_d} \bar{X}_{\theta}^{M,0} = (\eta_{\theta} + 1)^M \bar{X}_{\theta}^{M,0}. \quad (66) $$

From (65) and (66), $\forall c_d = 0, \dots, M$,

$$ \bar{X}_{\theta}^{M, c_d} = \binom{M}{c_d} \left( \frac{\eta_{\theta}}{\eta_{\theta} + 1} \right)^{c_d} \left( \frac{1}{\eta_{\theta} + 1} \right)^{M-c_d} = \binom{M}{c_d} (h_{\theta})^{c_d} (1-h_{\theta})^{M-c_d} \quad (67) $$

with $h_{\theta} = \frac{\eta_{\theta}}{\eta_{\theta}+1}$. Introducing (67) into (29), one obtains (34) with $F_{\theta}$ defined in (33). Thus one needs to solve (34) to determine $\bar{p}$, which is then used to deduce $\bar{X}_{\theta}^{M,d}$ using (67).

# APPENDIX D

## PROOF OF LEMMA 6

To prove Lemma 6, one needs first to investigate the monotonicity of $F_{\theta}$. To lighten the notations, let $\alpha = q_{\mathrm{FA}}(2)$, $\beta = q_D(1,1)$ and $\gamma = q_D(0,2)$. Then $h_0$ and $h_1$ defined in (31-32) can be rewritten as

$$ h_0(x, y) = \frac{\alpha p_0 x + \beta p_1 y}{p_0 x + p_1 y}, \quad h_1(x, y) = \frac{\beta p_0 x + \gamma p_1 y}{p_0 x + p_1 y}, \quad (68) $$

with $(x, y) \in \mathcal{P}_0$. One starts showing some monotonicity properties.

**Lemma 12.** If $\alpha < \beta < \gamma$, then $h_0$ and $h_1$ are decreasing with $x$ and increasing with $y$, for all $(x, y) \in \mathcal{P}_0$. If $\beta = \gamma$, then $h_1 = \beta = \gamma$ is a constant.

*Proof:* Since $\alpha < \beta \le \gamma$, one has

$$ \frac{\partial h_0}{\partial x} = \frac{(\alpha - \beta) p_0 p_1 y}{(p_0 x + p_1 y)^2} \le 0, \quad \frac{\partial h_0}{\partial y} = \frac{(\beta - \alpha) p_0 p_1 x}{(p_0 x + p_1 y)^2} \ge 0, $$

$$ \frac{\partial h_1}{\partial x} = \frac{(\beta - \gamma) p_0 p_1 y}{(p_0 x + p_1 y)^2} \le 0, \quad \frac{\partial h_1}{\partial y} = \frac{(\gamma - \beta) p_0 p_1 x}{(p_0 x + p_1 y)^2} \ge 0. $$

then Lemma 12 can be proved. $\square$

**Lemma 13.** For $z \in [0, 1]$, the family of functions

$$ f_i(z) = z^i (1-z)^{M-i}, \quad i = 0, \dots, M. \quad (69) $$

are increasing over $[0, \frac{i}{M}]$ and decreasing over $[\frac{i}{M}, 1]$.

*Proof:* Consider three possible situations: 1) If $i=0$, $f_0(z) = (1-z)^M$ is decreasing over $[0, 1]$. 2) If $i=M$, $f_M(z) = z^M$ is increasing over $[0, 1]$. 3) If $1 \le i \le M-1$,

$$ \frac{df_i}{dz} = z^{i-1} (1-z)^{M-i-1} (i-Mz), \quad (70) $$
and $df_i/dz \ge 0$ when $z \in [0, \frac{i}{M}]$ and $df_i/dz \le 0$ when $z \in [\frac{i}{M}, 1]$. Therefore, Lemma 13 holds $\forall i = 0, \dots, M$. $\square$

**Lemma 14.** If $0 < \nu < 1$, the function

$$g(z) = \sum_{i:i/M < \nu} \binom{M}{i} f_i(z) = \sum_{i:i/M < \nu} \binom{M}{i} z^i (1-z)^{M-i}, \quad (71)$$

is decreasing for all $z \in [0, 1]$.

*Proof:* First, consider $z \in [\nu, 1]$. In (71), each $i$ in the sum is such that $\frac{i}{M} < \nu \le z$. From Lemma 13, $f_i(z)$ is a decreasing function for any $\frac{i}{M} < z$, thus $g(z)$ is also decreasing with $z$.

Now, consider $z \in [0, \nu]$, one rewrites (71) as

$$g(z) = 1 - \sum_{i:i/M \ge \nu} \binom{M}{i} f_i(z), \quad (72)$$

in which each $i$ in the sum is such that $z < \nu \le \frac{i}{M}$. Applying again Lemma 13, since $f_i(z)$ is an increasing function for any $z \le \frac{i}{M}$, the sum in (72) is also increasing with $z$ and $g(z)$ is decreasing. Thus $g(z)$ decreases over $[0, 1]$. $\square$

Considering the functions $h_\theta$ and $g$, then one may rewrite $F_\theta$ as $F_\theta(x, y) = g(h_\theta(x, y))$, $\forall \theta \in \{0, 1\}$. The monotonicity of $F_0$ and $F_1$ is shown in the following lemma.

**Lemma 15.** If $\alpha < \beta < \gamma$, then $F_0$ and $F_1$ are increasing functions of $x$ and decreasing functions of $y$, for all $(x, y) \in \mathcal{P}_0$. If $\beta = \gamma$, then $F_1 = g(\beta) = g(\gamma)$ is a constant.

*Proof:* The proof of obtained by combining Lemma 12 and Lemma 14. $\square$

The proof of Lemma 6 is by induction. At the beginning, one has $0 \le p^{\theta_0}(0) \le 1$, thus $p^{\theta_0}_{\min}(0) = 0$ and $p^{\theta_0}_{\max}(0) = 1$. Using Lemma 15, one has $F_\theta(0, 1) \le F_\theta(p^{\theta_0}(0), p^{\theta_0}(0)) \le F_\theta(1, 0)$, thus

$$\left\{ \begin{array}{l} p_{\min}^{\theta_0}(1) = F_0(0, 1) = g(\beta) > 0 = p_{\min}^{\theta_0}(0), \\ p_{\max}^{\theta_0}(1) = F_0(1, 0) = g(\alpha) < 1 = p_{\max}^{\theta_0}(0), \\ p_{\min}^{\theta_0}(1) = F_1(0, 1) = g(\gamma) > 0 = p_{\min}^{\theta_0}(0), \\ p_{\max}^{\theta_0}(1) = F_1(1, 0) = g(\beta) < 1 = p_{\max}^{\theta_0}(0), \end{array} \right. \quad (73)$$

thus (40) and (41) are true for $n = 1$.

Consider then an arbitrary $n \in \mathbb{N}^*$ and $n > 1$. Assume that (40) and (41) are satisfied for any $n' < n$ and $n' \in \mathbb{N}^*$, one needs to see whether (40) and (41) are still satisfied for $n$. Applying Lemma 15 again, one obtains

$$p_{\min}^{\theta_0}(n) = F_{\theta} \left( p_{\min}^{\theta_0}(n-1), p_{\max}^{\theta_0}(n-1) \right) \\ > F_{\theta} \left( p_{\min}^{\theta_0}(n-2), p_{\max}^{\theta_0}(n-2) \right) = p_{\min}^{\theta_0}(n-1),$$

and

$$p_{\max}^{\theta_0}(n) = F_{\theta} \left( p_{\max}^{\theta_0}(n-1), p_{\min}^{\theta_0}(n-1) \right) \\ < F_{\theta} \left( p_{\max}^{\theta_0}(n-2), p_{\min}^{\theta_0}(n-2) \right) = p_{\max}^{\theta_0}(n-1),$$

Similarly, one gets $p_{\min}^{\theta_0}(n) > p_{\min}^{\theta_0}(n-1)$ and $p_{\max}^{\theta_0}(n) < p_{\max}^{\theta_0}(n-1)$.

F is contracting, i.e., that for any pairs $\mathbf{p} = (x, y) \in \mathcal{P}_n$ and $\mathbf{p} + \mathbf{\delta} = (x + \delta_x, y + \delta_y) \in \mathcal{P}_n$, one has

$$|\mathbf{F}(\mathbf{p} + \mathbf{\delta}) - \mathbf{F}(\mathbf{p})| < |\mathbf{\delta}|. \quad (74)$$

A sufficient condition to have (74) is that the eigenvalues of the matrix

$$\mathbf{A} = \begin{pmatrix} \frac{\partial F_0(x, y)}{\partial x} & \frac{\partial F_0(x, y)}{\partial y} \\ \frac{\partial F_1(x, y)}{\partial x} & \frac{\partial F_1(x, y)}{\partial y} \end{pmatrix}$$

have module less than 1. The eigenvalues of $\mathbf{A}$ are the solutions of

$$z^2 - \left( \frac{\partial F_0}{\partial x} + \frac{\partial F_1}{\partial y} \right) z + \left( \frac{\partial F_0}{\partial x} \frac{\partial F_1}{\partial y} - \frac{\partial F_0}{\partial y} \frac{\partial F_1}{\partial x} \right) = 0. \quad (75)$$

As in Appendix D, denote $\alpha = q_{\mathrm{FA}}(2)$, $\beta = q_{\mathrm{D}}(1, 1)$ and $\gamma = q_{\mathrm{D}}(0, 2)$. First, one evaluates

$$\frac{\partial F_0}{\partial x} \frac{\partial F_1}{\partial y} - \frac{\partial F_0}{\partial y} \frac{\partial F_1}{\partial x} = \frac{\partial g}{\partial h_0} \frac{\partial g}{\partial h_1} \left( \frac{\partial h_0}{\partial x} \frac{\partial h_1}{\partial y} - \frac{\partial h_0}{\partial y} \frac{\partial h_1}{\partial x} \right) \stackrel{(a)}{=} 0,$$

where (a) comes from $\frac{\partial h_0}{\partial x} \frac{\partial h_1}{\partial y} = \frac{\partial h_0}{\partial y} \frac{\partial h_1}{\partial x}$, using the partial derivatives calculated in the proof of Lemma 12. Then, the solutions of (75) are $z_1 = \frac{\partial F_0}{\partial x} + \frac{\partial F_1}{\partial y}$ and $z_2 = 0$. Hence, it suffices to prove that $|z_1| < 1$.

We begin with the evaluation of an upper bound of the partial derivative of $F_0(x, y)$ with respect to $x$

$$\begin{align*} \frac{\partial F_0(x, y)}{\partial x} &= \frac{\partial g(h_0(x, y))}{\partial x} = \frac{\partial g}{\partial h_0} \cdot \frac{\partial h_0}{\partial x} \\ &\stackrel{(a)}{=} \frac{(\beta - \alpha) p_0 p_1 y}{(p_0 x + p_1 y)^2} \sum_{i:i/M < \nu} \binom{M}{i} h_0^i (1 - h_0)^{M-i} \frac{h_0 M - i}{h_0 (1 - h_0)} \\ &\stackrel{(b)}{\le} \frac{(\beta - \alpha) p_0 p_1 y}{(p_0 x + p_1 y)^2} F_0(x, y) \frac{M}{1 - h_0} \le c_0 (\alpha, \beta, \gamma, M, \nu, n), \quad (76) \end{align*}$$

where (a) is obtained using (70), (b) comes from $i \ge 0$, and $c_0$ is defined in (37). Meanwhile, from Lemma 15, one has $\partial F_0(x, y) / \partial x \ge 0$, as $F_0$ is an increasing function of $x$.

Similarly,

$$\begin{align*} \frac{\partial F_1(x, y)}{\partial y} &= \frac{\partial g(h_1(x, y))}{\partial y} = \frac{\partial g}{\partial h_1} \cdot \frac{\partial h_1}{\partial y} \\ &= \frac{(\gamma - \beta) p_0 p_1 x}{(p_0 x + p_1 y)^2} \sum_{i:i/M < \nu} \binom{M}{i} h_1^i (1 - h_1)^{M-i} \frac{i - h_1 M}{h_1 (1 - h_1)} \\ &\ge \frac{(\gamma - \beta) p_0 p_1 x}{(p_0 x + p_1 y)^2} F_1(x, y) \frac{-M}{1 - h_1} \ge -c_1 (\alpha, \beta, \gamma, M, \nu, n), \quad (77) \end{align*}$$

and $\partial F_1(x, y) / \partial y \le 0$ as $F_1$ is a non-decreasing function of $y$. One concludes that

$$-c_1 \le \frac{\partial F_0(x, y)}{\partial x} + \frac{\partial F_1(x, y)}{\partial y} \le c_0,$$

thus $c_0 < 1$ and $c_1 < 1$ lead to $|z_1| < 1$, which ensures the uniqueness of the equilibrium.

# APPENDIX E

## PROOF OF PROPOSITION 7

As seen in the proof of Proposition 5, $\forall n \in \mathbb{N}^*$, $\mathbf{F}(\mathbf{p})$ maps $\mathcal{P}_n$ to $\mathcal{P}_n$, with $\mathcal{P}_n = [p_{\min}^{\theta_0}(n), p_{\max}^{\theta_0}(n)] \times [p_{\min}^{\theta_0}(n), p_{\max}^{\theta_0}(n)]$. In order to apply Banach's fixed-point theorem [42] to prove Proposition 7, it suffices to show that

# APPENDIX F

## PROOF OF PROPOSITION 9

First, one shows that if $\nu < q_D(1, 1)$, then for any $\varepsilon > 0$, there exists $M > M'$, such that $\bar{p}^{10} < \varepsilon$.
From Lemma 6, $\bar{p}^{10}$ can be bounded as

$$ 
\begin{aligned}
\bar{p}^{10} &= F_1(\bar{p}^{00}, \bar{p}^{10}) \\
&< \sum_{c_d:c_d/M < \nu} \binom{M}{c_d} (q_D(1,1))^{c_d} (1-q_D(1,1))^{M-c_d}
\end{aligned}
\qquad (78)
$$

Consider $\Phi_1, \Phi_2, \dots$ an infinite sequence of i.i.d. binary random variables with $\mathbb{P}\{\Phi_m = 1\} = q_D(1,1)$. For any $\varrho \in [0,1]$ such that $\varrho M \in \mathbb{N}^+$, one has

$$ \mathbb{P}\left\{\frac{\sum_{m=1}^{M} \Phi_m}{M} = \varrho\right\} = \binom{M}{\varrho M} (q_D(1,1))^{\varrho M} (1-q_D(1,1))^{M(1-\varrho)}. $$

According to the weak law of large numbers [46], for $\varepsilon > 0$, there exists $M'$, such that for any $M > M'$, one has

$$ \mathbb{P}\left\{\left|\frac{\sum_{m=1}^{M} \Phi_m}{M} - q_D(1,1)\right| > q_D(1,1)\right\} < \varepsilon. \qquad (79) $$

From (79), one also has

$$ 
\begin{aligned}
& \sum_{c_d:c_d/M < (q_D(1,1)-\varepsilon)} (q_D(1,1))^{c_d} (1-q_D(1,1))^{M-c_d} \\
&= \mathbb{P}\left\{\frac{\sum_{m=1}^{M} \Phi_m}{M} - q_D(1,1) < -\varepsilon\right\} \le \mathbb{P}\left\{\left|\frac{\sum_{m=1}^{M} \Phi_m}{M} - q_D(1,1)\right| > \varepsilon\right\} \\
&< \varepsilon.
\end{aligned}
\qquad (80)
$$

If $\nu < q_D(1,1) - \varepsilon$, then using (80), the bound of $\bar{p}^{10}$ in (78) may be further written as

$$ 
\begin{aligned}
\bar{p}^{10} &< \sum_{c_d:c_d/M < \nu} \binom{M}{c_d} (q_D(1,1))^{c_d} (1-q_D(1,1))^{M-c_d} \\
&\le \sum_{c_d:c_d/M < (q_D(1,1)-\varepsilon)} (q_D(1,1))^{c_d} (1-q_D(1,1))^{M-c_d} < \varepsilon. \quad (81)
\end{aligned}
$$

From Lemma 12 and the fact that $q_{\text{FA}}(2) \le \bar{p}^{10} \le q_D(1,1)$ and $0 \le \bar{p}^{10} < \varepsilon$, one has $h_0(\bar{p}^{00}, \bar{p}^{10}) \in [q_{\text{FA}}(2), \chi(\varepsilon)]$, with

$$ \chi(\varepsilon) = \frac{p_0 (q_{\text{FA}}(2))^2 + p_1 q_D(1,1) \varepsilon}{p_0 q_{\text{FA}}(2) + p_1 \varepsilon}. \qquad (82) $$

Thus, according to Lemma 14,

$$ 
\begin{aligned}
\bar{p}^{00} &= F_0(\bar{p}^{00}, \bar{p}^{10}) = g(h_0(\bar{p}^{00}, \bar{p}^{10})) \\
&\ge g(\chi(\varepsilon)) = \sum_{c_d:c_d/M < \nu} \binom{M}{c_d} (\chi(\varepsilon))^{c_d} (1-\chi(\varepsilon))^{M-c_d}. \quad (83)
\end{aligned}
$$

Using derivations similar to those leading to (80), one gets

$$ \sum_{c_d:c_d/M > (\chi(\varepsilon)+\varepsilon)} \binom{M}{c_d} (\chi(\varepsilon))^{c_d} (1-\chi(\varepsilon))^{M-c_d} < \varepsilon, \qquad (84) $$

which leads to

$$ \sum_{c_d:c_d/M \le (\chi(\varepsilon)+\varepsilon)} \binom{M}{c_d} (\chi(\varepsilon))^{c_d} (1-\chi(\varepsilon))^{M-c_d} \ge 1-\varepsilon. \qquad (85) $$

If $\nu > \chi(\varepsilon) + \varepsilon$, then

$$ 
\begin{aligned}
\bar{p}^{00} &\ge \sum_{c_d:c_d/M < \nu} \binom{M}{c_d} (\chi(\varepsilon))^{c_d} (1-\chi(\varepsilon))^{M-c_d} \\
&\ge \sum_{c_d:c_d/M \le (\chi(\varepsilon)+\varepsilon)} (\chi(\varepsilon))^{c_d} (1-\chi(\varepsilon))^{M-c_d} \ge 1-\varepsilon. \quad (86)
\end{aligned}
$$

As a conclusion, for any $\varepsilon > 0$, if $\chi(\varepsilon) + \varepsilon < \nu < q_D(1,1) - \varepsilon$, then $\bar{p}^{00} \ge 1 - \varepsilon$ and $\bar{p}^{10} < \varepsilon$. Since $\lim_{\varepsilon \to 0} \chi(\varepsilon) = q_{\text{FA}}(2)$, one concludes that if $q_{\text{FA}}(2) < \nu < q_D(1,1)$, one obtains (43).