Decision Fusion for Large-Scale Sensor Networks with Nonideal Channels

Abstract

Since there has been an increasing interest in the areas of Internet of Things (IoT) and artificial intelligence that often deals with a large number of sensors, this chapter investigates the decision fusion problem for large-scale sensor networks. Due to unavoidable transmission channel interference, we consider sensor networks with nonideal channels that are prone to errors. When the fusion rule is fixed, we present the necessary condition for the optimal sensor rules that minimize the Monte Carlo cost function. For the K-out-of-L fusion rule chosen very often in practice, we analytically derive the optimal sensor rules. For general fusion rules, a Monte Carlo Gauss-Seidel optimization algorithm is developed to search for the optimal sensor rules. The complexity of the new algorithm is of the order of OLN compared with OLNL of the previous algorithm that was based on Riemann sum approximation, where L is the number of sensors and N is the number of samples. Thus, the proposed method allows us to design the decision fusion rule for large-scale sensor networks. Moreover, the algorithm is generalized to simultaneously search for the optimal sensor rules and the optimal fusion rule. Finally, numerical examples show the effectiveness of the new algorithms for large-scale sensor networks with nonideal channels.

1. Introduction

Distributed detection has been an active research area in the past decades [1, 2, 3, 4, 5, 6, 7]. It involves the design of decision rules for the sensors¹ and the fusion rule [8]. Early work on distributed detection mainly focused on conditionally independent sensor observations, such as [2, 4, 9, 10], and the resulting optimal sensor decision rules, as well as the fusion rule, were likelihood ratio tests (LRTs). Details on distributed detection with conditionally independent sensor observations can be seen in [1, 6, 7] and references therein.

In this chapter, we focus on conditionally dependent observations in sensor networks. In [5], the computational difficulty of obtaining the optimal sensor rules was shown by a rigorous mathematical approach. Some early progress was made on the derivation of sensor rules for the dependent observation case such as in [11, 12, 13, 14, 15]. More recently, a hierarchical conditional independence model was provided that was applicable to some specific classes of multisensor detection problems with dependent observations [16]. Copula-based distributed decision fusion methods have been proposed to deal with dependent observations in sensor networks, such as [17, 18, 19] and references therein. Given a fusion rule, Monte Carlo methods were proposed to reduce the computational complexity of deriving sensor decision rules with ideal channels in [20, 21], and the optimal sensor rules were obtained analytically for the $K$-out-of-$L$ fusion rule in [20].

Some works on the derivation of optimal fusion rules can be seen in [15, 22, 23, 24]. For some specific parallel network decision systems, a unified fusion rule was presented in [15]. Some further results on the problem are available in [25, 26]. In [27], the authors provided methods that search for the sensor rules and the fusion rule simultaneously by combining the methods of [2] and [15] in order to attain near-optimal system performance.

The works discussed thus far assumed the availability of ideal channels in sensor networks. However, channel errors between the sensors and the fusion center are omnipresent in practical multisensor detection networks, and, therefore, studies on multisensor detection in the presence of nonideal channels have attracted some recent interest, such as in [8, 28, 29, 30, 31, 32, 33]. Under the Neyman-Pearson criterion, the design of sensor rules in the presence of nonideal channels was addressed in [32]. The parallel fusion structure was extended by incorporating the fading channel layer and two alternative fusion schemes were presented based on fixed sensor rules in [28]. It was shown that the optimal sensor decision rule that minimizes the error probability at the fusion center is equivalent to a local LRT for independent sensor observations in [29]. Under Neyman-Pearson and Bayesian criteria, the work was generalized to dependent and noisy channels, respectively, in [8]. In [31], the authors considered the optimal sensor rules with channel errors via Riemann sum approximation under a given fusion rule for general dependent sensor observations. Although the method based on the Riemann sum approximation has been developed for dependent observations with channel errors, it is too computationally expensive to be of practical use in large-scale sensor networks.

In this chapter, a Monte Carlo importance sampling method is provided to reduce the computational complexity of multisensor detection fusion with channel errors. Based on the strong law of large numbers, the Bayesian cost function is approximated by its empirical average through the Monte Carlo importance sampling method. The main contributions of this chapter are listed below:

1. When the fusion rule is fixed, we derive a necessary condition for the optimal sensor rules that minimize the approximated Bayesian cost function. A Monte Carlo Gauss-Seidel optimization algorithm is developed and it is shown to be finitely convergent. The complexity of the new algorithm is shown to be of the order of $O\left(\mathit{LN}\right)$ compared with $O\left({\mathit{LN}}^L\right)$ of the previous algorithm based on the Riemann sum approximation.

2. When the fusion rule is the $K$-out-of-$L$ rule, we prove that there exists an analytical form for the optimal sensor rules in the presence of nonideal channels. Thus, the proposed method allows us to design decision rules for large-scale sensor networks.

3. The Monte Carlo Gauss-Seidel optimization algorithm is extended to simultaneously search for the optimal sensor rules and the optimal fusion rule.

Numerical examples show the effectiveness of the new algorithms for large-scale sensor networks with dependent observations and channel errors.

The rest of this chapter is organized as follows: In Section 2, the parallel binary Bayesian detection network with channel errors is formulated and the Monte Carlo cost function is introduced. In Section 3, the necessary condition for the optimal sensor rules is presented. For the $K$-out-of-$L$ fusion rule, the analytical form for the optimal sensor rules is provided. In Section 4, the Monte Carlo Gauss-Seidel iterative algorithm and its convergence analysis are presented. The extension to search for the optimal sensor rules and the optimal fusion rule are simultaneously described in Section 5. Simulation results are provided in Section 6. Conclusions are contained in Section 7.

2. Preliminaries

2.1 Problem formulation

The $L$-sensor parallel Bayesian detection network structure with two hypotheses H₀ and H₁ in the presence of nonideal channels is considered (see Figure 1). Assume that $y_1$, $y_2$, …, $y_L$ are sensor observations and the $j$th sensor compresses the $n_j$-dimension vector observation $y_j$ to one bit: $I_j\left(y_j\right):{\mathbb{R}}^{n_j}\to \left\{0,1\right\},j=1,\dots ,L$. For notational convenience, $n_j=1$ in the following description. The $L$ sensors transmit the compressed data to the fusion center and the fusion center makes the decision between $H_0$ and $H_1$. Since external interference and internal errors may occur, the channels are not reliable and the fusion center may not correctly receive the symbol $I_j$ sent by the $j$th sensor. Let $I_j^0$ denote the received bit by the fusion center for $j=1,2,\dots ,L$. Generally speaking, $I_j^0$ may not be equal to $I_j$. The definition and assumptions on channel errors (see e.g., [29, 31]) are summarized below:

Figure 1.

The $L$-sensor parallel binary Bayesian detection network structure in the presence of nonideal channels.

Definition 1: The channel errors between the $j$th sensor and the fusion center are described as $P_j^{\mathit{ce}1}=P\left(I_j^0=0|I_j=1\right)$ and $P_j^{\mathit{ce}0}=P\left(I_j^0=1|I_j=0\right)$ for $j=1,2,\dots ,L$, where $P_j^{\mathit{ce}1}$ is the probability of channel error when the $j$th sensor sends 1 but the fusion center receives 0, and $P_j^{\mathit{ce}0}$ is the probability of channel error when the $j$th sensor sends 0 but the fusion center receives 1.

Assumption 1: The probabilities of channel error are statistically independent of the hypotheses, namely $P\left(I_j^0|I_j,H_{\nu }\right)=P\left(I_j^0|I_j\right),\nu =0,1$.

Remark 1: Assumption 1 is due to the hierarchical structure based on the Markov property (see [29]).

Assumption 2: The channels that connect the sensors to the fusion center are independent, i.e., $P\left(I_1^0,I_2^0,\dots ,I_L^0|I_1,I_2,\dots ,I_L\right)={\prod }_{j=1}^{L}P\left(I_j^0|I_j\right)$.

We consider the parallel binary Bayesian detection network with nonideal channels that is built on the above definition and assumptions. The final decision is made by the fusion center based on the received binary bits $\left(I_1^0,I_2^0,\dots ,I_L^0\right)$ from the $L$ sensors. From the definition of a general Bayesian cost function given in [25], the $L$-sensor binary Bayesian cost function with channel errors at the fusion center can be written as follows:

$$\begin{array}{l}C\left(I_1^0\left(y_1\right),\dots ,I_L^0\left(y_L\right);F^0\right)=C_{00}{P}_{0}P\left(F^0=0|H_0\right)+C_{01}{P}_{1}P\left(F^0=0|H_1\right)\\ +C_{10}{P}_{0}P\left(F^0=1|H_0\right)+C_{11}{P}_{1}P\left(F^0=1|H_1\right)\end{array}$$E1

$$=c+\mathit{aP}\left(F^0=0|H_1\right)-\mathit{bP}\left(F^0=0|H_0\right),$$E2

where $C_{\mathit{\alpha \beta }},\alpha ,\beta =0,1$ are suitable cost coefficients, $P_0$ and $P_1$ are the prior probabilities for the hypotheses $H_0$ and $H_1$, respectively, $F^0$ is the fusion rule, and $P\left(F^0=\mu |H_{\nu }\right),\mu ,\nu =0,1$ denotes the conditional probability of the event that the fusion center decides in favor of hypothesis $\mu$ when the real hypothesis is $H_{\nu }$. The cost function (1) is simplified to (2) by defining $\mathit{c}={\mathit{C}}_{10}{\mathit{P}}_0+{\mathit{C}}_{11}{\mathit{P}}_1,\mathit{a}={\mathit{P}}_1\left({\mathit{C}}_{01}-{\mathit{C}}_{11}\right),$$\mathit{b}={\mathit{P}}_0\left({\mathit{C}}_{10}-{\mathit{C}}_{00}\right)$. $F^0$ is actually a function of the disjoint set of all possible binary messages $\left(I_1^0,I_2^0,\dots ,I_L^0\right)$. The received decisions are divided into two sets denoted as $H_0^0$ and $H_1^0$ which are given by

$$\begin{array}{l}{H}_0^0=\left\{\left(u_1^0,u_2^0,\dots ,u_L^0\right):F^0\left(\left(I_1^0,I_2^0,\dots ,I_L^0\right)\right)=0,I_j^0=u_j^0,u_j^0=0/1,j=1,\dots ,L\right\};\\ H_1^0=\left\{\left(u_1^0,u_2^0,\dots ,u_L^0\right):F^0\left(\left(I_1^0,I_2^0,\dots ,I_L^0\right)\right)=1,I_j^0=u_j^0,u_j^0=0/1,j=1,\dots ,L\right\}.\end{array}$$

Obviously, $H^0=\left\{\left(u_1^0,u_2^0,\dots ,u_L^0\right):I_j^0=u_j^0,u_j^0=0/1,j=1,\dots ,L\right\}=H_0^0\cup H_1^0$. For any binary decisions $\left(I_1^0,I_2^0,\dots ,I_L^0\right)$ received by the fusion center, the original sensor decision bits before transmission are $\left(I_1,I_2,\dots ,I_L\right)$ and they consist of the set $H=\left\{\left(u_1,u_2,\dots ,u_L\right):I_j=u_j,u_j=0/1,j=1,\dots ,L\right\}.$ Therefore, based on the law of total probability, the conditional probability formula, and Assumption 1:

$$P\left(F^0=0|H_{\nu }\right)=\sum _{s^0\in H_0^0}P\left(D^0|H_{\nu }\right)=\sum _{s^0\in H_0^0}\sum _{s\in H}P\left(D^0|D\right)P\left(D|H_{\nu }\right),$$E3

where $D^0=\left(I_1^0,I_2^0,\dots ,I_L^0\right)$, $s^0=\left(s^0\left(1\right),\dots ,s^0\left(L\right)\right)$, $I_j^0=s^0\left(j\right)$, and $s^0\left(j\right)=0/1$ is a specific value of $I_j^0$; in the same way, $D=\left(I_1,I_2,\dots ,I_L\right)$, $s=\left(s\left(1\right),\dots ,s\left(L\right)\right)$, $I_j=s\left(j\right)$, and $s\left(j\right)=0/1$ is a specific value of $I_j$. Strictly speaking, we should use $P\left(D^0=s^0|H_{\nu }\right)$ to represent $P\left(D^0|H_{\nu }\right)$ and we use the latter for notational simplicity. It is similar to $P\left(D|H_{\nu }\right)$. Based on Assumption 2:

$$P\left(D^0|D\right)=\prod _{j=1}^{L}P\left(I_j^0|I_j\right),$$E4

where for any $1\le j\le L$

$$P\left(I_j^0|I_j\right)=\left(1-P_j^{\mathit{ce}0}\right)\left(1-I_j^0\right)\left(1-I_j\right)+P_j^{\mathit{ce}0}{I}_j^0\left(1-I_j\right)+\left(1-P_j^{\mathit{ce}1}\right)I_j^0{I}_j+P_j^{\mathit{ce}1}\left(1-I_j^0\right)I_j.$$E5

Thus, the cost function (2) becomes

$$\begin{array}{l}C\left(I_1^0\left(y_1\right),\dots ,I_L^0\left(y_L\right);F^0\right)=c+\sum _{s^0\in H_0^0}\sum _{s\in H}P\left(D^0|D\right)\left[\mathit{aP}\left(D|H_1\right)-\mathit{bP}\left(D|H_0\right)\right]\end{array}$$E6

$$\triangleq C\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right),$$E7

where $P^{\mathit{ce}0}=\left(P_1^{\mathit{ce}0},\dots ,P_L^{\mathit{ce}0}\right),P^{\mathit{ce}1}=\left(P_1^{\mathit{ce}1},\dots ,P_L^{\mathit{ce}1}\right)$. Hence, the cost function now becomes a function of the sensor rules $\left(I_1,\dots ,I_L\right)$, the probabilities of channel errors $P^{\mathit{ce}0},P^{\mathit{ce}1}$, and the fusion rule $F^0$. The goal of this chapter is to optimize the sensor rules and the fusion rule so as to minimize the cost function with known probabilities of channel errors.

We rewrite $\mathit{aP}\left(D|H_1\right)-\mathit{bP}\left(D|H_0\right)$ as follows:

$$\begin{array}{l}\mathit{aP}\left(D|H_1\right)-\mathit{bP}\left(D|H_0\right)={\int }_{{\mathrm{\Omega }}_s}\mathit{ap}\left(y_1,\dots ,y_L|H_1\right)-\mathit{bp}\left(y_1,\dots ,y_L|H_0\right){\mathit{dy}}_1\cdots {\mathit{dy}}_L\\ =\int I_{{\mathrm{\Omega }}_s}\left[\mathit{ap}\left(y_1,\dots ,y_L|H_1\right)-\mathit{bp}\left(y_1,\dots ,y_L|H_0\right)\right]{\mathit{dy}}_1\cdots {\mathit{dy}}_L,\end{array}$$E8

where ${\mathrm{\Omega }}_s=\left\{\left(y_1,\dots ,y_L\right):I_1\left(y_1\right)=s\left(1\right),\dots ,I_L\left(y_L\right)=s\left(L\right)\right\}$, $I_{{\mathrm{\Omega }}_s}$ is an indicator function on ${\mathrm{\Omega }}_s$, and the region of integration in (8) is the full space. Assume that $p\left(y_1,y_2,\dots ,y_L|H_{\nu }\right),\nu =0,1$ (or $p\left(y|H_{\nu }\right)$) are the known conditional joint probability density functions. If not, we can learn the joint probability density functions from training data using copula functions (see, e.g., [17]). Note that $I_1\left(y_1\right),\dots ,I_L\left(y_L\right)$ are indicator functions and $s\left(j\right)=0/1,j=1,\dots ,L$,

$$\begin{array}{c}{I}_{{\mathrm{\Omega }}_s}=I_{\left\{\left(y_1,\dots ,y_L\right):I_1\left(y_1\right)=s\left(1\right),\dots ,I_L\left(y_L\right)=s\left(L\right)\right\}}\\ =I_{\left\{y_1:I_1\left(y_1\right)=s\left(1\right)\right\}}\cdots I_{\left\{y_L:I_L\left(y_L\right)=s\left(L\right)\right\}}\\ =\left[\left(1-I_1\right)\left(1-s\left(1\right)\right)+I_{1}s\left(1\right)\right]\cdots \left[\left(1-I_L\right)\left(1-s\left(L\right)\right)+I_{L}s\left(L\right)\right].\end{array}$$E9

For simplicity, denote $Q_j\left(I_j\right)=\left(1-I_j\right)\left(1-s\left(j\right)\right)+I_{j}s\left(j\right)$. Substituting (8) into (6),

$$\begin{array}{l}C\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\\ =c+\sum _{s^0\in H_0^0}\sum _{s\in H}P(\left.D^0|D\right)\cdot \int Q_1\left(I_1\right)\cdots Q_L\left(I_L\right)[\left.\mathit{ap}\left(y|H_1\right)-\mathit{bp}\left(y|H_0\right)\right]\mathit{dy}\\ =c+\int P_{H_0^0}\widehat{L}\left(y\right)\mathit{dy},\end{array}$$E10

where $P_{H_0^0}={\sum }_{s^0\in H_0^0}{\sum }_{s\in H}P\left(D^0|D\right)Q_1\left(I_1\right)\cdots Q_L\left(I_L\right)$ and $\widehat{L}\left(y\right)=\mathit{ap}\left(y|H_1\right)-\mathit{bp}\left(y|H_0\right)$. Note that from the definition of $H_0^0,H_1^0,H^0$, and $H$, we have

$$\begin{array}{l}{P}_{H_0^0}=\sum _{s^0\in H^0}\left[1-F^0\left(D^0\right)\right]\sum _{s\in H}P\left(D^0|D\right)Q_1\left(I_1\right)\cdots Q_L\left(I_L\right)\\ =\sum _{k^{\prime }=1}^{2^L}\sum _{k=1}^{2^L}\left[1-F^0\left(s_{k^{\prime }}\right)\right]P\left(s_{k^{\prime }}|s_k\right)\cdot \prod _{j=1}^L\left\{\left[1-I_j\left(y_j\right)\right]\left[1-s_k\left(j\right)\right]+I_j\left(y_j\right)s_k\left(j\right)\right\},\end{array}$$E11

where $s_{k^{\prime }}$ is the element of $H^0$ and $s_k$ is the element of $H$. $F^0\left(D^0\right)=0/1$ is used in the first equality. The second equality holds since there are $2^L$ elements in both $H$ and $H^0$.

2.2 Monte Carlo cost function

An essential difficulty of the Bayesian cost function (10) is the required high dimensional integration when dealing with large-scale sensor networks. Monte Carlo importance sampling is an attractive method to deal with this problem. In this subsection, we approximate the Bayesian cost function (10) by the Monte Carlo importance sampling method (see, e.g., [34, 35]). According to (10),

$$\begin{array}{l}C\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\\ =c+\int P_{H_0^0}\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\frac{\widehat{L}\left(y\right)}{g\left(y\right)}\cdot g\left(y\right)\mathit{dy}\end{array}$$E12

$$={\mathbb{E}}_g\frac{P_{H_0^0}\left(Y\right)\widehat{L}\left(Y\right)}{g\left(Y\right)}+c,$$E13

where $y=\left(y_1,y_2,\dots ,y_L\right)$, and $g\left(y\right)$ is a given importance sampling density such that (12) is well-defined (i.e., $g\left(y\right)>0$). In (13), the expectation is taken with respect to the importance sampling density $g$. Consequently, assume that $N$ samples $Y_1,\dots ,Y_N$ are generated from the density $g$, that is, $Y\sim g\left(y\right)$, where $Y_i=\left[Y_{i1},Y_{i2},\dots ,Y_{\mathit{iL}}\right]$. Then

$$\begin{array}{l}C\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\approx \frac{1}{N}\sum _{i=1}^N\frac{P_{H_0^0}\left(Y_{i1},Y_{i2},\dots ,Y_{\mathit{iL}}\right)\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}+c\end{array}$$E14

$$\triangleq C_{\mathit{MC}}\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right).$$E15

Based on the strong law of large numbers, the expectation (13) can be approximated by the empirical average (14). Denote (14), namely the Monte Carlo cost function, as $C_{\mathit{MC}}\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)$. The optimal importance sampling density is $g\left(y_1,y_2,\dots ,y_L\right)\propto \mid P_{H_0^0}\widehat{L}(y_1,y_2,$$\dots ,y_L)\mid$ (see, e.g., [34, 35]).

The initial goal is to minimize the Bayesian cost function (10). Instead, we can minimize the Monte Carlo cost function (15) by selecting a set of optimal sensor rules $I_1\left(y_1\right),I_2\left(y_2\right),\dots ,I_L\left(y_L\right)$ and an optimal fusion rule $F^0$. In this manner, the high-dimensional integration problem is converted to a problem where we need to deal with the single summation objective function for large-scale sensor networks. Thus, for dependent observations with channel errors, the computational complexity is reduced significantly by the Monte Carlo importance sampling method. In the following sections, we assume that the samples drawn from the importance sampling density are fixed so that $C_{\mathit{MC}}(I_1,\dots ,I_L;F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},$$N)$ does not have any randomness, since only deterministic decision rules are considered in this chapter.

3. A necessary condition for the optimal sensor rules

In this section, when the fusion rule is fixed, we derive a necessary condition for the optimal sensor rules that minimize the Monte Carlo cost function. First, we need some equivalent transformations for $P_{H_0^0}$. Then based on the transformations, the necessary condition can be obtained. At the same time, an analytical result is obtained when the fusion rule is the $K$-out-of-$L$ rule.

3.1 Necessary condition

First, we need some equivalent transformations for $P_{H_0^0}$.

Lemma 1 $P_{H_0^0}$ can be rewritten as follows:

$$\begin{array}{l}{P}_{H_0^0}\triangleq \left[1-I_j\left(y_j\right)\right]P_{j1}\left(I_1\left(y_1\right),\dots ,I_{j-1}\left(y_{j-1}\right),I_{j+1}\left(y_{j+1}\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\\ +P_{j2}\left(I_1\left(y_1\right),\dots ,I_{j-1}\left(y_{j-1}\right),I_{j+1}\left(y_{j+1}\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right),\end{array}$$E16

where for $j=1,2,\dots ,L$,

$$P_{j1}\left(\cdot \right)\triangleq \sum _{k^{\prime }=1}^{2^L}\left[1-F\left(s_{k^{\prime }}\right)\right]\left(1-P_j^{\mathit{ce}0}-P_j^{\mathit{ce}1}\right)\left(1-2s_{k^{\prime }}\left(j\right)\right)P_{m\ne j},$$E17

$$P_{j2}\left(\cdot \right)\triangleq \sum _{k^{\prime }=1}^{2^L}\left[1-F\left(s_{k^{\prime }}\right)\right]\left[s_{k^{\prime }}\left(j\right)+P_j^{\mathit{ce}1}\left(1-2s_{k^{\prime }}\left(j\right)\right)\right]P_{m\ne j},$$E18

$$\begin{array}{l}{P}_{m\ne j}\triangleq \prod _{m,m\ne j}^L\left\{\left(1-P_m^{\mathit{ce}0}\right)\left(1-s_{k^{\prime }}\left(m\right)\right)\left(1-I_m\right)+P_m^{\mathit{ce}0}{s}_{k^{\prime }}\left(m\right)(1-I_m\right.)\\ \left.+\left(1-P_m^{\mathit{ce}1}\right)s_{k^{\prime }}\left(m\right)I_m+P_m^{\mathit{ce}1}\left(1-s_{k^{\prime }}\left(m\right)\right)I_m\right\}.\end{array}$$E19

Proof: If $s_k\left(m\right)=I_m\left(y_m\right)$ for all $m=1,\dots ,L$, then the continued product term ${\prod }_{m=1}^L\left\{\left[1-I_m\left(y_m\right)\right]\left[1-s_k\left(m\right)\right]+I_m\left(y_m\right)s_k\left(m\right)\right\}=1$ in $P_{H_0^0}$. Otherwise, it is 0. Thus, $P_{H_0^0}$ can be rewritten as $P_{H_0^0}={\sum }_{k^{\prime }=1}^{2^L}\left[1-F\left(s_{k^{\prime }}\right)\right]P\left(s_{k^{\prime }}|\left(I_1,I_2,\dots ,I_L\right)\right),$ where the terms that equal zero are omitted and $P\left(s_{k^{\prime }}|\left(I_1,I_2,\dots ,I_L\right)\right)={\prod }_{j=1}^{L}P\left(s_{k^{\prime }}\left(j\right)|I_j\right).$ Recalling the conditional probability formula (5), we rewrite $P\left(s_{k^{\prime }}\left(j\right)|I_j\right)$ as $P\left(s_{k^{\prime }}\left(j\right)|I_j\right)=\left[1-I_j\right]\left(1-P_j^{\mathit{ce}0}-P_j^{\mathit{ce}1}\right)\left(1-2s_{k^{\prime }}\left(j\right)\right)+s_{k^{\prime }}\left(j\right)+P_j^{\mathit{ce}1}\left(1-2s_{k^{\prime }}\left(j\right)\right).$ Based on these transformations, $P_{H_0^0}$ can be decomposed as (16).

Remark 2: Note that $P_{j1}\left(\cdot \right)$ and $P_{j2}\left(\cdot \right)$ are both independent of $I_j\left(y_j\right)$ for $j=1,\dots ,L$. In addition, they can also be applied in the Riemann sum approximation (see, e.g., [31]). Compared with [36], the sum of $2^L$ terms about $s_k$ is eliminated and it greatly reduces the computational time. In addition, the expression for $P_{j1}\left(\cdot \right)$ given in (17) is also a key equation in the following results:

Substituting the transformations (16) into (15), we obtain

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ =c+\frac{1}{N}\sum _{i=1}^N\left\{\left[1-I_j\left(Y_{\mathit{ij}}\right)\right]P_{j1}(I_1\left(Y_{i1}\right),\dots ,I_{j-1}\left(Y_{i\left(j-1\right)}\right),I_{j+1}\left(Y_{i\left(j+1\right)}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1})\right.\\ \left.+P_{j2}(I_1\left(Y_{i1}\right),\dots ,I_{j-1}\left(Y_{i\left(j-1\right)}\right),I_{j+1}\left(Y_{i\left(j+1\right)}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1})\right\}\cdot \frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)},\end{array}$$E20

where $Y_i=\left(Y_{i1},Y_{i2},\dots ,Y_{\mathit{iL}}\right)$. According to (20), the necessary condition for the optimal sensor rules that minimize $C_{\mathit{MC}}(I_1\left(y_1\right),$$\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N)$ is stated in the following lemma:

Lemma 2: Let $\left\{I_1\left(y_1\right),\dots ,I_L\left(y_L\right)\right\}$ be the set of optimal sensor rules, i.e., they minimize $C_{\mathit{MC}}(I_1\left(y_1\right),$$\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N)$ in the parallel Bayesian detection network, then they must satisfy the following equations:

$$I_1\left(Y_{i1}\right)=I\left[P_{11}\left(I_2\left(Y_{i2}\right),I_3\left(Y_{i3}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\cdot \widehat{L}\left(Y_i\right)\right],$$E21

$$I_2\left(Y_{i2}\right)=I\left[P_{21}\left(I_1\left(Y_{i1}\right),I_3\left(Y_{i3}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\cdot \widehat{L}\left(Y_i\right)\right],$$E22

$$\cdots \cdots$$

$$I_L\left(Y_{\mathit{iL}}\right)=I\left[P_{L1}\left(I_1\left(Y_{i1}\right),I_2\left(Y_{i2}\right),\dots ,I_{L-1}\left(Y_{i\left(L-1\right)}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\cdot \widehat{L}\left(Y_i\right)\right],$$E23

where $P_{j1}\left(\cdot \right)$ are defined by (17) and $I\left[\cdot \right]$ is an indicator function defined as follows:

$$I\left[x\right]=\left\{\begin{array}{ll}1,& \text{if}x\ge 0;\\ 0,& \text{if}x<0.\end{array}\right.$$E24

Proof: Note that both $P_{j1}\left(\cdot \right)$ and $P_{j2}\left(\cdot \right)$ are independent of $I_j\left(y_j\right)$ for $j=1,\dots ,L$. If $I_1\left(Y_{i1}\right)$ minimizes the Monte Carlo cost function under the given $I_2\left(Y_{i2}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right)$, we only need to minimize the first term of the summation in (20), that is, $\left[1-I_1\left(Y_{i1}\right)\right]P_{11}\left(I_2\left(Y_{i2}\right),I_3\left(Y_{i3}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}.$ Note that the value of $I_1\left(Y_{i1}\right)$ is $0$ or $1$ and $g\left(y\right)$ is well defined, that is, $g\left(Y_i\right)>0$, $I_1\left(Y_{i1}\right)$ should be equal to 1 when $P_{11}(I_2\left(Y_{i2}\right),$$I_3\left(Y_{i3}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right);F^0;$$P^{\mathit{ce}0},P^{\mathit{ce}1})\widehat{L}\left(Y_i\right)\ge 0$ for $i=1,\dots ,N$, otherwise it should be equal to 0. Therefore, we obtain (21) by the definition of $I\left[x\right]$ in (24). Similarly, we obtain (22) and (23) by minimizing (20).

3.2 An analytical result for the $K$-out-of-$L$ rule

When the fusion rule is a $K$-out-of-$L$ rule, we would obtain an analytical result in the presence of nonideal channels. It is described as follows:

Theorem 1.1: If the fusion rule is a $K$-out-of-$L$ rule and the probabilities of channel errors are less than 0.5 (i.e., $0<P_j^{\mathit{ce}0}<\mathrm{0.5,0}<P_j^{\mathit{ce}1}<0.5$) for each channel, the optimal sensor rules are $I_j\left(Y_{\mathit{ij}}\right)=I\left[\widehat{L}\left(Y_i\right)\right]$ for $i=1,\dots ,N$ and $j=1,\dots ,L$.

Proof: From Lemma 1, we know

$$\begin{array}{l}{P}_{j1}\left(\cdot \right)=\sum _{k^{\prime }=1}^{2^L}\left[1-F\left(s_{k^{\prime }}\right)\right]P_{m\ne j}\cdot \left(1-P_j^{\mathit{ce}0}-P_j^{\mathit{ce}1}\right)\left[1-2s_{k^{\prime }}\left(j\right)\right]\\ =\left(1-P_j^{\mathit{ce}0}-P_j^{\mathit{ce}1}\right)\sum _{k^{\prime }=1}^{2^{L-1}}\left[\left(1-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)\right)-\left(1-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)\right)\right]P_{m\ne j}\\ =\left(1-P_j^{\mathit{ce}0}-P_j^{\mathit{ce}1}\right)\sum _{k^{\prime }=1}^{2^{L-1}}\left[F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)\right]P_{m\ne j}.\end{array}$$E25

Since $0<P_j^{\mathit{ce}0}<\mathrm{0.5,0}<P_j^{\mathit{ce}1}<0.5$, we have $1-P_j^{\mathit{ce}0}-P_j^{\mathit{ce}1}>0$. Obviously, $P_{m\ne j}>0$ holds from its definition. If $\left[F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)\right]\ge 0$, $P_{j1}\left(\cdot \right)\ge 0$ can be derived. When the fusion rule is a $K$-out-of-$L$ rule, $F\left(s_{k^{\prime }}\right)=I\left[{\sum }_{j=1}^L{s}_{k^{\prime }}\left(j\right)-K\right]$. Thus,

$$\begin{array}{l}F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)=I\left[\sum _{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+1-K\right],\\ F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)=I\left[\sum _{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+0-K\right].\end{array}$$

If ${\sum }_{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+0-K\ge 0$, then ${\sum }_{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+1-K\ge 0$, and we can get that $F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)=0$. If ${\sum }_{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+0-K<0$, then $F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)\ge 0$. In a word, $F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)\ge 0$ is derived, thus $P_{j1}\ge 0$. It is easy to find a $s_{k^{\prime }}\left(m\right),m\ne j$ so that ${\sum }_{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+0=K-1$ and ${\sum }_{m=1,m\ne j}^L{s}_{k^{\prime }}\left(m\right)+1=K$. Thus, there must exist $F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=1\right)-F\left(s_{k^{\prime }}|s_{k^{\prime }}\left(j\right)=0\right)>0$. Therefore, the $P_{j1}>0$ is derived. Recalling the necessary condition for the optimal sensor rules, that is, $I\left(Y_{\mathit{ij}}\right)=I\left[P_{j1}\left(\cdot \right)\cdot \widehat{L}\left(Y_i\right)\right]$, the proof is completed.

Remark 3: The $K$-out-of-$L$ rule counts the number of sensors that vote in favor of $H_1$ and compares it with a given threshold $K$ [37]. It is also referred to as the counting rule or voting rule and is widely used in the practical decision fusion area [38, 39]. It encompasses a general class of fusion rules such as AND, OR, and Majority Boolean fusion rules [40]. The reason we assume that the probabilities of channel errors are less than 0.5 is based on practical considerations. If the probabilities of channel errors are greater than or equal to 0.5, the channel is totally unreliable and the performance is not better than a random decision. Obviously, the analytical solution is very efficient to tackle large-scale sensor networks with dependent observations and channel errors.

4. Monte Carlo Gauss-Seidel iterative algorithm and its convergence

For general fusion rules that do not have the form of a $K$-out-of-$L$ rule, an efficient algorithm can be obtained that is inspired by Lemma 2. Next, we present a Monte Carlo Gauss-Seidel iterative algorithm and derive its convergence, when the fusion rule is fixed.

4.1 Monte Carlo Gauss-Seidel iterative algorithm

Based on the fixed-point type necessary condition given in Lemma 2, the Monte Carlo Gauss-Seidel iterative algorithm is presented in Algorithm 1.

Algorithm 1: Optimization of the sensor rules.

Given the fusion rule $F^0$:

• Step 1: Generate $N$ samples: $Y_1,\dots ,Y_N\sim g\left(y\right)$, where $g\left(y\right)$ is an importance sampling density and $Y_i=\left[Y_{i1},Y_{i2},\dots ,Y_{\mathit{iL}}\right]$.

• Step 2: Initialize the $L$ sensor rules, for $j=1,2,\dots ,L$ and $i=1,\dots ,N$,

  $$I_j^{\left(0\right)}\left(Y_{\mathit{ij}}\right)=0/1.$$E26

• Step 3: Iteratively search for the $L$ sensor rules until a termination criterion in Step 4 is satisfied. The $n+1$th iteration is given as follows: for $i=1,\dots ,N$,

  $$I_1^{\left(n+1\right)}\left(Y_{i1}\right)=I\left[P_{11}\left(I_2^{\left(n\right)}\left(Y_{i2}\right),I_3^{\left(n\right)}\left(Y_{i3}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\widehat{L}\right(Y_i],$$E27
  $$I_2^{\left(n+1\right)}\left(Y_{i2}\right)=I\left[P_{21}\left(I_1^{\left(n+1\right)}\left(Y_{i1}\right),I_3^{\left(n\right)}\left(Y_{i3}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\widehat{L}\left(Y_i\right)\right],$$E28
  $$I_L^{\left(n+1\right)}\left(Y_{\mathit{iL}}\right)=I\left[P_{L1}\left(I_1^{\left(n+1\right)}\left(Y_{i1}\right),\dots ,I_{L-1}^{\left(n+1\right)}\left(Y_{i\left(L-1\right)}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\widehat{L}\left(Y_i\right)\right].$$E29

• Step 4: For $i=1,\dots ,N$, the termination criterion of the iteration process is

  $$\begin{array}{c}{I}_1^{\left(n+1\right)}\left(Y_{i1}\right)=I_1^{\left(n\right)}\left(Y_{i1}\right),\\ I_2^{\left(n+1\right)}\left(Y_{i2}\right)=I_2^{\left(n\right)}\left(Y_{i2}\right),\\ \cdots \cdots \\ I_L^{\left(n+1\right)}\left(Y_{\mathit{iL}}\right)=I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right).\end{array}$$E30

Remark 4: When we obtain $I_1\left(Y_{i1}\right)$ for $i=1,\dots ,N$, we can compress $y_1$ by defining $I_1\left(y_1\right)=I_1\left(Y_{i1}\right)$ if the distance $\Vert y_1-Y_{i1}\Vert \le \Vert y_1-Y_{i^{\prime }1}\Vert$ for all $i^{\prime }\ne i$. In the same way, we can compress $y_j$ for $j=2,\dots ,L$. In fact, the method is to find one nearest neighbor of $y_j$ for $j=1,\dots ,L$ and use the corresponding compression rule. Moreover, we can utilize the k-nearest neighbor (knn) to compress $y_j$ (see more in [41]).

Remark 5: The main computation burden of Algorithm 1 is included in (27)–(29). If we let the number of discretized points $N_1=N_2=\dots =N_L=N$ in [31], then $P_{j1}\left(\cdot \right)\widehat{L}\left(Y_i\right)$, $j=1,\dots ,L$, and $i=1,\dots ,N$ are computed ${\mathit{LN}}^L$ times for each iteration, as in [31]. But they only need to be computed $\mathit{LN}$ times in Algorithm 1. Thus, the computational complexity of Algorithm 1, i.e., $O\left(\mathit{LN}\right)$ is much less than that in [31], that is, $O\left({\mathit{LN}}^L\right)$. It is more efficient to tackle large-scale sensor networks with dependent observations and channel errors.

4.2 Convergence of the iterative algorithm

Now, we show that Algorithm 1 must converge to a stationary point and the algorithm cannot oscillate infinitely, that is, it terminates after a finite number of iterations.

Lemma 3: Given the fusion rule $F^0$, for any initial values $\left(I_1^{\left(0\right)},\dots ,I_L^{\left(0\right)}\right)$ in (26), $C_{\mathit{MC}}\left(I_1^{\left(n\right)},\dots ,I_j^{\left(n\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)$ must converge to a stationary point after a finite number of iterations.

Proof: For $j=1,\dots ,L$, we denote $C_{\mathit{MC}}$(20) in the $n+1$th iteration process by

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_j^{\left(n+1\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ =\frac{1}{N}\sum _{i=1}^N\left\{\left[1-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]P_{j1}(I_1^{\left(n+1\right)}\left(Y_{i1}\right),\dots ,I_{j-1}^{\left(n+1\right)}\left(Y_{i\left(j-1\right)}\right),I_{j+1}^{\left(n\right)}\left(Y_{i\left(j+1\right)}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^0;\right.\\ \left.P^{\mathit{ce}0},P^{\mathit{ce}1},N)+P_{j2}(I_1^{\left(n+1\right)}\left(Y_{i1}\right),\dots ,I_{j-1}^{\left(n+1\right)}\left(Y_{i\left(j-1\right)}\right),I_{j+1}^{\left(n\right)}\left(Y_{i\left(j+1\right)}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N)\right\}\\ \cdot \frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}+c.\end{array}$$E31

Similarly, we denote the $\left(n+1\right)$th iteration process of the iterative items $P_{j1}\left(\cdot \right)\widehat{L}\left(\cdot \right)$ in (27)–(29) by

$$G_j^i=P_{j1}\left(I_1^{\left(n+1\right)}\left(Y_{i1}\right),\dots ,I_{j-1}^{\left(n+1\right)}\left(Y_{i\left(j-1\right)}\right),I_{j+1}^{\left(n\right)}\left(Y_{i\left(j+1\right)}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\widehat{L}\left(Y_i\right),$$E32

for $i=1,\dots ,N$ and $j=1,\dots ,L.$ Plugging $G_j^i$ into (31), we know

$$C_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_j^{\left(n+1\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)=\frac{1}{N}\sum _{i=1}^N\frac{\left[1-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]}{g\left(Y_i\right)}{G}_j^i+C_j^i,$$E33

where $C_j^i=c+\frac{1}{N}{\sum }_{i=1}^N{P}_{j2}(I_1^{\left(n+1\right)}\left(Y_{i1}\right),\dots ,I_{j-1}^{\left(n+1\right)}\left(Y_{i\left(j-1\right)}\right),I_{j+1}^{\left(n\right)}\left(Y_{i\left(j+1\right)}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^0;$$P^{\mathit{ce}0},P^{\mathit{ce}1},N)\frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}$ is independent of $I_j^{\left(n\right)}$ and $I_j^{\left(n+1\right)}$. Splitting $1-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)$ into two terms, we obtain

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_j^{\left(n+1\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ =\frac{1}{N}\sum _{i=1}^N\frac{\left[1-I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)\right]+\left[I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]}{g\left(Y_i\right)}{G}_j^i+C_j^i\\ =\frac{1}{N}\sum _{i=1}^N\frac{\left[1-I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)\right]}{g\left(Y_i\right)}{G}_j^i+C_j^i+\frac{1}{N}\sum _{i=1}^N\frac{\left[I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]}{g\left(Y_i\right)}{G}_j^i\\ =C_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_{j-1}^{\left(n+1\right)},I_j^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)+D_j^{\left(n+1\right)},\end{array}$$E34

where

$$D_j^{\left(n+1\right)}=\frac{1}{N}\sum _{i=1}^N\frac{\left[I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]}{g\left(Y_i\right)}{G}_j^i.$$E35

Note that (27)–(29) imply that $I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)=0$ if and only if $G_j^i<0$ and $I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)=1$ if and only if $G_j^i\ge 0$ for $i=1,\dots ,N,j=1,\dots ,L$. It means

$$\left[I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]G_j^i\le 0.$$E36

Thus, for $\forall i,j$

$$\left[I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]G_j^i/g\left(Y_i\right)\le 0,$$E37

where the inequality (35) holds since $g\left(\cdot \right)$ is well-defined (i.e., $g\left(\cdot \right)>0$). Substituting (35) into (33) yields $D_j^{\left(n+1\right)}\le 0$. Thus, for $\forall j\le L$,

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_j^{\left(n+1\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ \le C_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_{j-1}^{\left(n+1\right)},I_j^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right).\end{array}$$E38

Furthermore,

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1^{\left(n+1\right)},I_2^{\left(n+1\right)},\dots ,I_L^{\left(n+1\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ \le C_{\mathit{MC}}\left(I_1^{\left(n\right)},I_2^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right).\end{array}$$E39

It means $C_{\mathit{MC}}$ is nonincreasing. Note that $C_{\mathit{MC}}(I_1^{\left(n\right)},I_2^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N)$ is a finite value. We conclude that it must converge to a stationary point after a finite number of iterations.

Theorem 1.2: Given the fusion rule $F^0$, the sensor rules $I_1^{\left(n\right)},I_2^{\left(n\right)},\dots ,I_L^{\left(n\right)}$ are finitely convergent, i.e., Algorithm 1 converges after a finite number of iterations.

Proof: By Lemma 3, $C_{\mathit{MC}}$ must attain a stationary point after a finite number of iterations. It means that the value of $C_{\mathit{MC}}$ cannot change after $n$th iteration, that is,

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_j^{\left(n+1\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ =C_{\mathit{MC}}\left(I_1^{\left(n+1\right)},\dots ,I_{j-1}^{\left(n+1\right)},I_j^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right).\end{array}$$E40

Using (32) and (37), we derive that $D_j^{\left(n+1\right)}=0$. Combining (33)–(35), we know

$$\left[I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)\right]G_j^i=0,\mathit{\text{for}}i=1,\dots ,N,$$E41

which implies either $I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)-I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)=0,i.e.,I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)=I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)$ or $G_j^i=0,i.e.,I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)=1,I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)=0.$ It follows that when $C_{\mathit{MC}}$ converges to a stationary point, either $I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)$ is invariant or $I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)=1,I_j^{\left(n\right)}\left(Y_{\mathit{ij}}\right)=0$. Namely, $I_j^{\left(n+1\right)}\left(Y_{\mathit{ij}}\right)$ can only change from 0 to 1 at most a finite number of times. Therefore, the $I_1^{\left(n\right)},I_2^{\left(n\right)},\dots ,I_L^{\left(n\right)}$ are finitely convergent.

5. Extension for simultaneous search for the optimal sensor rules and fusion rule

In this section, we extend the Monte Carlo method to search for the optimal sensor rules and the optimal fusion rule simultaneously. Firstly, the necessary condition is generalized to search for the optimal sensor rules and the optimal fusion rule simultaneously. Secondly, we describe a generalized Monte Carlo Gauss-Seidel iterative algorithm. We also give the convergence of the iterative algorithm.

5.1 A necessary condition for the optimal sensor rules and the optimal fusion rule

Note that (15) can be rewritten as follows:

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ =c+\frac{1}{N}\sum _{i=1}^N\sum _{k^{\prime }=1}^{2^L}\sum _{k=1}^{2^L}\left[1-F^0\left(s_{k^{\prime }}\right)\right]P\left(s_{k^{\prime }}|s_k\right)\cdot P_{s_k}(\mathbf{I}\left(Y_i\right)\frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}\\ =c+\frac{1}{N}\sum _{k^{\prime }=1}^{2^L}\left[1-F^0\left(s_{k^{\prime }}\right)\right]\sum _{i=1}^N\sum _{k=1}^{2^L}P\left(s_{k^{\prime }}|s_k\right)\cdot P_{s_k}(\mathbf{I}\left(Y_i\right)\frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)},\end{array}$$E42

where $P_{s_k}\left(\mathbf{I}\left(Y_i\right)\right)\triangleq {\prod }_{j=1}^L\left[s_k\left(j\right)I_j\left(Y_{\mathit{ij}}\right)+\left(1-s_k\left(j\right)\right)\left(1-I_j\left(Y_{\mathit{ij}}\right)\right)\right]$ and $\mathbf{I}\left(Y_i\right)=\left(I_1\left(Y_{i1}\right),I_2\left(Y_{i2}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right)\right)$. Since $P_{s_k}(\mathbf{I}\left(Y_i\right)=1$ if and only if $I_j=s_k\left(j\right)$ for all $j=1,\dots ,L$, (39) can be simplified as follows:

$$\begin{array}{l}{C}_{\mathit{MC}}\left(I_1\left(y_1\right),\dots ,I_L\left(y_L\right);F^0;P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)\\ =c+\frac{1}{N}\sum _{k^{\prime }=1}^{2^L}\left[1-F^0\left(s_{k^{\prime }}\right)\right]\cdot \sum _{i=1}^{N}P\left(s_{k^{\prime }}|\left(I_1\left(Y_{i1}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right)\right)\right)\frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)},\end{array}$$E43

where the terms $P_{s_k}\left(\mathbf{I}\left(Y_i\right)\right)=0$ are eliminated.

Remark 6: According to (20) and (40), the necessary condition for the optimal sensor rules is similar to Lemma 2 and the necessary condition for the optimal fusion rule is given by

$$F^0\left(s_{k^{\prime }}\right)=I\left[\sum _{i=1}^{N}P\left(s_{k^{\prime }}|\left(I_1\left(Y_{i1}\right),\dots ,I_L\left(Y_{\mathit{iL}}\right)\right)\right)\cdot \frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}\right]$$E44

for $k^{\prime }=1,\dots ,2^L$. The proofs are similar to Lemma 2.

5.2 Generalized Monte Carlo Gauss-Seidel iterative algorithm

Based on the fixed-point type necessary condition, the generalized Monte Carlo Gauss-Seidel iterative algorithm is presented in Algorithm 2.

Remark 7: For any initial values $\left(I_1^{\left(0\right)},\dots ,I_L^{\left(0\right)};F^{0\left(0\right)}\right)$, the Monte Carlo cost function $C_{\mathit{MC}}\left(I_1^{\left(n\right)},\dots ,I_j^{\left(n\right)},I_{j+1}^{\left(n\right)},\dots ,I_L^{\left(n\right)};F^{0\left(n\right)};P^{\mathit{ce}0},P^{\mathit{ce}1},N\right)$ must converge to a stationary point and Algorithm 2 terminates after a finite number of iterations. The proofs are similar to those of Lemma 3 and Theorem 1.2.

Algorithm 2: Simultaneous optimization of the sensor rules and the fusion rule.

• Step 1: Generate $N$ samples: $Y_1,\dots ,Y_N\sim g\left(y\right)$, where $g\left(y\right)$ is an importance sampling density and $Y_i=\left[Y_{i1},Y_{i2},\dots ,Y_{\mathit{iL}}\right]$.

• Step 2: Initialize the $L$ sensor rules and the fusion rule, respectively, for $j=1,2,\dots ,L$, $i=1,\dots ,N$, and $k^{\prime }=1,\dots ,2^L$,

  $$I_j^{\left(0\right)}\left(Y_{\mathit{ij}}\right)=0/1,F^{0\left(0\right)}\left(s_{k^{\prime }}\right)=0/1.$$

• Step 3: Iteratively search for the $L$ sensor rules and the fusion rule until a termination criterion in Step 4 is satisfied. The $n+1$th iteration is given as follows: for $i=1,\dots ,N$ and $k^{\prime }=1,\dots ,2^L$

  $$\begin{array}{l}{I}_1^{\left(n+1\right)}\left(Y_{i1}\right)=I\left[P_{11}\left(I_2^{\left(n\right)}\left(Y_{i2}\right),I_3^{\left(n\right)}\left(Y_{i3}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^{0\left(n\right)};P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\cdot \widehat{L}\left(Y_i\right)\right],\\ I_2^{\left(n+1\right)}\left(Y_{i2}\right)=I\left[P_{21}\left(I_1^{\left(n+1\right)}\left(Y_{i1}\right),I_3^{\left(n\right)}\left(Y_{i3}\right),\dots ,I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);F^{0\left(n\right)};P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\cdot \widehat{L}\left(Y_i\right)\right],\\ \cdots \cdots \\ I_L^{\left(n+1\right)}\left(Y_{\mathit{iL}}\right)=I\left[P_{L1}\left(I_1^{\left(n+1\right)}\left(Y_{i1}\right),I_2^{\left(n+1\right)}\left(Y_{i2}\right),\dots ,I_{L-1}^{\left(n+1\right)}\left(Y_{i\left(L-1\right)}\right);F^{0\left(n\right)};P^{\mathit{ce}0},P^{\mathit{ce}1}\right)\cdot \widehat{L}\left(Y_i\right)\right],\\ F^{0\left(n+1\right)}\left(s_{k^{\prime }}\right)=I\left[\sum _{i=1}^{N}P\left(s_{k^{\prime }}|\left(I_1^{\left(n+1\right)}\left(Y_{i1}\right),I_2^{\left(n+1\right)}\left(Y_{i2}\right),\dots ,I_L^{\left(n+1\right)}\left(Y_{\mathit{iL}}\right)\right)\right)\frac{\widehat{L}\left(Y_i\right)}{g\left(Y_i\right)}\right].\end{array}$$

• Step 4: For $i=1,\dots ,N$ and $k^{\prime }=1,2,\dots ,2^L$, the termination criterion of the iteration process is

  $$\begin{array}{c}{I}_1^{\left(n+1\right)}\left(Y_{i1}\right)=I_1^{\left(n\right)}\left(Y_{i1}\right),\\ I_2^{\left(n+1\right)}\left(Y_{i2}\right)=I_2^{\left(n\right)}\left(Y_{i2}\right),\\ \cdots \cdots \\ I_L^{\left(n+1\right)}\left(Y_{\mathit{iL}}\right)=I_L^{\left(n\right)}\left(Y_{\mathit{iL}}\right);\\ F^{0\left(n+1\right)}\left(s_{k^{\prime }}\right)=F^{0\left(n\right)}\left(s_{k^{\prime }}\right).\end{array}$$

6. Numerical examples

In this section, in order to evaluate the performance of Algorithms 1 and 2, we present some examples with a Gaussian signal $s$ observed in the presence of Gaussian sensor noises.

The random signal $s$ and observation noises $v_1$, $v_2$, …, $v_L$ are as follows:

$$H_0:y_j=v_j;H_1:y_j=s+v_j,\mathit{\text{for}}j=1,\dots ,L,$$E45

where $v_1$, $v_2$ …, $v_L$, $s$ are all mutually independent and

$$v_j\sim N\left(\mathrm{0,0.6}\right),s\sim N\left(\mathrm{1,0.4}\right),\mathit{\text{for}}j=1,\dots ,L.$$

Thus, given $H_0$ and $H_1$, the two conditional probability density functions are

$$p\left(y_1,y_2,\dots ,y_L|H_0\right)\sim N\left({\mu }_0,{\mathrm{\Sigma }}_0\right),p\left(y_1,y_2,\dots ,y_L|H_1\right)\sim N\left({\mu }_1,{\mathrm{\Sigma }}_1\right),$$

where ${\mu }_0,{\mu }_1,{\mathrm{\Sigma }}_0,{\mathrm{\Sigma }}_1$ are easily obtained from the relationship of $s,v_1$, $v_2$, …, $v_L$.

Assume that each sensor is required to transmit a bit through a channel with probabilities of $P_j^{\mathit{ceo}}=P_j^{\mathit{ce}1}=p$, where $p=\mathrm{0.05,0.15,0.3}$, for $j=1,2,\dots ,L$. In the cost function (2), let the cost coefficients $C_{00}=C_{11}=0$ and $C_{10}=C_{01}=1$. The receiver operating characteristics (ROC) curves are used to evaluate the performance of the algorithms. $\mathit{Pf}$ and $\mathit{Pd}$ denote the probability of false alarm and the probability of detection, respectively.

6.1 Two-sensor network

We compare the Monte Carlo Gauss-Seidel iterative algorithm with the centralized algorithm and the iterative algorithm based on the Riemann sum approximation in [31] by using the receiver operating characteristics (ROC) curves.

In this case, we know ${\mu }_0={\left[0,0\right]}^T,{\mu }_1=\left[1,1\right]){}^T$ and ${\mathrm{\Sigma }}_0=\left[\mathrm{0.6,0};\mathrm{0,0.6}\right]$, ${\mathrm{\Sigma }}_1=\left[\mathrm{1,0.4};\mathrm{0.4,1}\right]$. Some discrete values of $a$ and $b$ are used to plot ROC curves. We refer to the optimal importance sampling density $g\left(y\right)\propto \mid P_{H_0^0}\left(y\right)\widehat{L}\left(y\right)\mid$ in Section 2.2 and $\mid \widehat{L}\left(y\right)\mid =\mid \mathit{ap}\left(y|H_1\right)-\mathit{bp}\left(y|H_0\right)\mid$. The form is similar to the mixture-Gaussian distribution. Therefore, the importance sampling density $g\left(y\right)$ is chosen to be the mixture-Gaussian distribution. The effects of choosing different $g\left(y\right)$ in terms of the performance of the Monte Carlo method were shown in [21] via numerical examples. For Algorithm 1, we take $N=200$ samples from the density $g\left(y\right)$. For the Riemann sum approximation iterative algorithm in [31], we take a discretized step-size $\Delta =0.09$, $y_i\in \left[-8,10\right]$, i.e., $N_1=N_2=N=200$. The ROC curves for three important fusion rules: AND, OR, and XOR rules with $p=\mathrm{0.05,0.15,0.3}$ are plotted in Figure 2. We compare the computational time of the two algorithms with $p=0.15$ in Figure 3. Note that the analytical solution is used for the AND rule and the OR rule. Since the XOR rule is not a $K$-out-of-$L$ rule, we use Algorithm 1 to search for the sensor rules.

Figure 2.

Two-sensor ROC curves with the probability of channel errors $p=\mathrm{0.05,0.15,0.3}$.

Figure 3.

Two-sensor computational time as $N$ increases with the probability of channel errors $p=0.15$.

Some observations in Figures 2 and 3 are presented as follows:

• Given the fusion rule, the two points $\left(0,0\right)$ and $\left(1,1\right)$ may not be the beginning or ending points of the ROC curves, which is different from the case in the ideal channel cases. In addition, the larger the probability of channel errors is, the farther away from $\left(0,0\right)$ or $\left(1,1\right)$ the ROC curves are. A possible reason is that the detection probability is not equal to 0 or 1, even when the false alarm probability is 0 or 1 in the presence of channel errors.

• From Figure 2, when the probability of channel transmission errors increases, the decision fusion performance of different methods using the optimal sensor rules decreases.

• It can be seen in Figure 2 that the ROC curves of the new Monte Carlo approach are very close to those of the previous algorithm based on the Riemann sum approximation. However, from Figure 3, the computational time of the Monte Carlo importance sampling approximation is much less than that of the Riemann sum approximation for the three different fusion rules. It also implies that the new method can be used to deal with large-scale sensor networks.

• Note that the computational time of the AND rule and the OR rule is less than that of the XOR rule for the Monte Carlo importance sampling approximation from Figure 3. The reason is that the AND rule and the OR rule belong to the $K$-out-of-$L$ rules. The analytical form is used for the AND rule and the OR rule, therefore, the corresponding computation time is relatively lower.

6.2 Ten-sensor network

We consider a larger sensor network with 10 sensors, which cannot be dealt with by the previous decision fusion algorithm based on the Riemann sum approximation due to its heavy computation requirements. For different probabilities of channel errors, the ROC curves of the AND rule, the OR rule, the 4-out-of-10 rule, the 6-out-of-10 rule, and Algorithm 2 are plotted in Figure 4.

Figure 4.

Ten-sensor ROC curves with the probability of channel errors $p=\mathrm{0.05,0.15,0.3}$.

Some observations in Figure 4 are presented as follows:

• The ROC curves for the ten-sensor network exhibit similar behavior as those for the two-sensor network.

• Given the fusion rule and the probability of channel errors, the decision fusion performance of the AND rule is better than the other rules for a small false alarm probability and the decision fusion performance of the OR rule is better than the other rules for a large false alarm probability. The reason may be that both of them are extreme cases of the fusion rules. For other cases, the 4-out-of-10 rule and the 6-out-of-10 rule perform better than the AND rule and the OR rule, respectively.

• Regardless of the centralized detection algorithm, Figure 4 shows that the ROC curves generated by Algorithm 2 obtain almost the best performance for different probabilities of channel errors. It implies that a simultaneous search for the sensor rules and the fusion rule would provide better decision fusion performance.

6.3 One-hundred-sensor network

We consider a large-scale network with one hundred sensors. The parameter settings are similar to Section 6.2. The ROC curves of the 20-out-of-100 rule, the 40-out-of-100, the 50-out-of-100, the 60-out-of-100 rule, and the 80-out-of-100 rule are plotted in Figure 5.

Figure 5.

One-hundred-sensor ROC curves with the probability of channel errors $p=\mathrm{0.05,0.15,0.3}$.

From Figure 5, it can be seen that the dramatically lower computational requirement of our method enables us to handle a large sensor network consisting of one hundred sensors. This is due to the fact that we have shown that there exist analytical forms of the optimal sensor rules for the $K$-out-of-$L$ rule. In addition, the decision fusion performance of different methods is improved, as the number of sensors becomes large.

7. Conclusion

By employing the Monte Carlo importance sampling technique, decision fusion algorithms have been provided for large-scale sensor networks with dependent observations and channel errors. The Bayesian cost function is approximated by the Monte Carlo cost function. The necessary conditions for the optimal sensor rules and the optimal fusion rule that minimize the Monte Carlo cost function have been obtained. Computationally efficient Monte Carlo Gauss-Seidel iterative algorithms have been proposed to search for the optimal sensor rules and the optimal fusion rule. These algorithms have been shown to converge after a finite number of iterations. The computational complexity of the new algorithm (i.e., $O\left(\mathit{LN}\right)$) is much less than that of the previous algorithm based on Riemann sum approximation (i.e., $O\left({\mathit{LN}}^L\right)$). For the $K$-out-of-$L$ rule, an analytical solution has been presented for the optimal sensor rules. Simulations have demonstrated the effectiveness of Algorithms 1 and 2. Future work will include the decision fusion algorithms under the Monte Carlo framework for other networks such as tandem networks, tree networks, and sensor networks under Byzantine attack.

Acknowledgments

This work was supported in part by the Sichuan Youth Science and Technology Innovation Team under Grant 2022JDTD0014, and Grant 2021JDJQ0036.