Analysis of relay selection schemes in underlay cognitive radio non-orthogonal multiple access networks

Abstract

This article investigates the impacts of relay selection schemes on cooperative underlay cognitive radio non-orthogonal multiple access networks, where the partial relay selection scheme, the max–min relay selection scheme and the two-stage relay selection scheme are applied in the network. Moreover, decode-and-forward protocol is used at the transmission relays. What’s more, in order to show the effect of the schemes on the considered network, the closed-form expressions and asymptotic expressions for the outage probability of the system are derived. Furthermore, the outage performance under the effect of perfect and imperfect successive interference cancellation is analysed. Numerical results are given to illustrate the impacts of the relay selection schemes, the number of relays, the residual interference factor and the power allocation factor on the outage performance. Finally, Monte Carlo simulations are presented to validate the accuracy of the numerical results.

Introduction

Non-orthogonal multiple access (NOMA) has been proven as a promising multiple access technique to improve spectral efficiency (SE), enhance cell-edge throughput and decrease transmission latency, where multiple users with different power allocation (PA) levels are allowed to share the same frequency/time/code resources. ¹ The transmitter utilizes Superposition Coding (SC) to distribute the power of user signals, and the receiver employs successive interference cancellation (SIC) to distinguish different signals. ² Saito et al. ³ investigated the system-level performance of NOMA taking into account practical aspects of the cellular system. The results showed that the overall cell throughput, cell-edge user throughput and the degree of proportional fairness of NOMA are all superior to that for orthogonal multiple access (OMA).

Cooperative relaying technology can combat multipath fading and improve the reliability of wireless networks. Using relay nodes in NOMA can significantly improve SE and outage performance. Ding et al. ⁴ proposed a cooperative NOMA scheme, where the users with good channel conditions operated as relays. The results showed that the cooperative NOMA scheme outperforms both the non-cooperative NOMA and cooperative OMA. Abbasi et al. ⁵ considered a cooperative NOMA network in which the relay works in amplify-and-forward (AF) mode. Also, the approximated ergodic rate and asymptotic outage probability (OP) were derived. Pei et al. ⁶ investigated cooperative NOMA with a dedicated decode-and-forward (DF) relay, in which both full-duplex (FD) and half-duplex (HD) protocols are considered for the relay.

Relay selection (RS) techniques also remain an important issue due to the fact that the best relay can achieve performance improvements in reliability and throughput. Ding et al. ⁷ proposed a two-stage RS scheme for NOMA relay networks with fixed PA and derived closed-form expressions for the OP, besides, the authors demonstrated that the two-stage RS scheme can achieve minimal OP than other RS schemes. Xu et al. ⁸ further proposed two-stage weighted-max–min and two-stage max-weighted-harmonic-mean RS schemes in the case of fixed and adaptive PA at the relays, respectively. The outage probabilities of the two proposed RS schemes were also derived. In the research work, Yang et al. ⁹ proposed a new two-stage RS scheme and combined DF and AF relaying with this two-stage scheme, respectively. Besides, PA factors are related to the source–relay and relay–destination channels. Some works focus on the partial RS scheme on NOMA networks.^10,11 Lee et al. ¹⁰ considered the AF protocol and derived the closed-form expressions for the OP of terminal users. Asymptotic analysis at high signal-to-noise ratio (SNR) was also carried out. The simulation results showed that there is almost no gain in outage performance when the number of relays exceeds two for partial RS scheme. Hoang et al. ¹¹ combined NOMA with radio frequency (RF) energy harvesting (EH), besides, the closed-form expressions for the ergodic capacity and the OP were derived in the case of perfect and imperfect SIC. The results showed that the performance of the system is significantly influenced by the level of residual interference. The max–min RS scheme in NOMA system was studied in Kim ¹² where the candidate relay forwards the information to destination user in the presence of direct path. The results showed that the outage performance with max–min RS outperforms the random RS. Lee et al. ¹³ applied NOMA to an underlay cognitive radio (CR) scenario, where a user with strong channel gain was selected as a relay for assisting another user with poor channel gains. Besides, the impacts of the multiple antennas and the number of cooperative NOMA users on the system outage performance were investigated. Zhao et al. ¹⁴ proposed an NOMA-based joint relay-and-antenna selection scheme for Hybrid Satellite-Terrestrial Relay Network. The considered scheme can achieve the maximum communication rate of the secondary user when the primary user maintains the optimal outage performance.

With the rapid growth of data traffic, the shortage of spectrum resources will become the bottleneck of the development of wireless communication. The emergence of CR can improve spectrum utilization. ¹⁵ Particularly, underlay CR enables secondary users to access the spectrum of the primary network if the quality of service (QoS) of the primary users is guaranteed. ¹⁶ Hence, it can be foreseen that incorporating NOMA into underlay CR networks has the potential to increase SE and system capacity. The outage performance of the cooperative underlay cognitive radio non-orthogonal multiple access (CR-NOMA) networks was studied in some schemes,^17,18,19,20 where the interference temperature constraint (ITC) at the primary network was considered. In the research work, Chu and Zepernick ¹⁷ investigated the OP and ergodic capacity for secondary users and the whole system, what’s more, the impacts of the ITC, channel power gains and PA on the system performance were analysed. Im and Lee ¹⁸ considered the imperfect SIC in cooperative underlay CR-NOMA networks. Arzykulov et al. ¹⁹ showed that NOMA achieves better OP results compared to OMA, in which the OP of secondary users with imperfect channel state information (CSI) was investigated. In the research work, Nauryzbayev et al. ²⁰ considered the Nakagami-m fading channels and the closed-form expressions for the OP of user messages were derived.

The application of RS in underlay CR-NOMA networks has also been studied. Do and Le ²¹ discussed the AF-based opportunistic RS scheme and partial RS scheme. In addition, the influence of hardware imperfections on the outage performance is considered. The DF-based RS schemes have also been investigated.^22,23,24 Sultan ²² considered two NOMA transmission scenarios according to whether the near user assists transmission and proposed three RS schemes, besides, the closed-form expressions for the OP and bit error rate of the system are derived. In the research work, Sultan ²³ proposed a reactive RS scheme and derived the closed-form expressions for OP. The proposed RS scheme aims to maximize the second-hop SNR for the far user. Simulation results reveal the impacts of PA factor and interference threshold on outage performance. Do et al. ²⁴ considered an uplink–downlink CR-NOMA network, in which the partial RS scheme is exploited for both uplink and downlink communications. In the proposed model, the relays can exploit at FD or HD mode. The expressions of OP and ergodic capacity were also provided.

This article is the extension of Ding et al. ⁷ and Do et al. ²⁴ To the best of our knowledge, the max–min RS selection in CR-NOMA networks has not been discussed, and this is one of the main focuses of this study. In this article, we derive the closed-form expressions and asymptotic expressions for the OP of the system based on the proposed RS selection schemes. Besides, the imperfect SIC case is also investigated. The major contributions of this article are summarized as follows:

• First, three RS schemes in the underlay CR-NOMA network, namely, the partial RS scheme, the max–min RS scheme and the two-stage RS scheme are investigated. Moreover, the interference from the primary user to the secondary network is considered.

• Second, the closed-form expressions for the OP of the system are derived. Based on the numerical results, the outage performance of the three RS schemes is investigated. It is worth pointing out that an outage floor which depends on the ITC can be observed in the OP.

• Third, the asymptotic expressions for the OP are derived (1) when the transmit SNR goes to infinity and (2) when the ITC goes to infinity.

• Fourth, the outage performance under the effect of perfect SIC and imperfect SIC is analysed.

• Finally, the effect of the PA factor on system OP is analysed. Furthermore, Monte Carlo (MC) results are provided to validate the accuracy of the numerical results.

Performance analysis is the basis of engineering practice and is of great significance to guide the design of practical communication systems. From the theoretical and simulation results, it can be concluded that the two-stage RS scheme can achieve the best outage performance. Besides, the outage performance can be improved by increasing the number of relay nodes, while the performance gain can be neglected when the number of relays exceeds 2 for the partial RS scheme. In addition, the outage performance is severely affected by residual interference factor and PA factor.

The rest of this article is organized as follows. Section ‘System model’ illustrates the system model and the three RS schemes. The OP results achieved by the three RS schemes are shown in section ‘Formulation analysis in the case of perfect SIC’ and section ‘Formulation Analysis in the case of imperfect SIC’, respectively. Section ‘Performance analysis’ provides MC simulation results and necessary discussions. Finally, the conclusions are given in section ‘Conclusion’.

System model

As illustrated in Figure 1, consider a downlink DF underlay CR-NOMA network which includes a primary user $\mathrm{PU}$ , a secondary source $S$ , two destination users $(D_1\mathrm{and}{D}_2)$ and $N$ DF relays $\{R_n{\}}_{n=1}^N$ . It is assumed that the direct link is not available due to deep shadow fading. In addition, all of the channels experience independent Rayleigh fading. It is also assumed that all nodes are equipped with a single antenna and operate in HD mode. It is assumed that short distances between the relays compared to the distances $S-\{R_n\}$ and $\{R_n\}-D_i$ , $i\in \{1,2\}$ , respectively. Thus, we denote $d_{\mathrm{SR}}=d_{S{R}_n}$ , $d_{\mathrm{RP}}=d_{R_{n}P}$ and $d_i=d_{R_n{D}_i}$ for all $n$ , where $d$ represents the distance of the link. The complex channel coefficients of the link $S\to \mathrm{PU}$ , $S\to R_n$ and $R_n\to \mathrm{PU}$ are denoted by $h_{\mathrm{SP}}~\mathrm{CN}(0,{\lambda }_{\mathrm{SP}}=d_{\mathrm{SP}}^{-\alpha })$ , $h_{S{R}_n}~\mathrm{CN}(0,{\lambda }_0=d_0^{-\alpha })$ and $h_{R_{n}P}~\mathrm{CN}(0,{\lambda }_{\mathrm{RP}}=d_{\mathrm{RP}}^{-\alpha })$ , respectively. $g_{R_n{D}_i}~\mathrm{CN}(0,{\lambda }_i=d_i^{-\alpha })$ is the complex channel coefficient between $R_n$ and $D_i$ , $i\in \{1,2\}$ . Far user and near user can be distinguished by a suitable channel feedback mechanism. Suppose that $D_1$ is the far user while $D_2$ is the near user, thus the path loss and shadowing effect of the link $R\to D_1$ are more severe than those of $R\to D_2$ , that is, ${\lambda }_1<{\lambda }_2$ . Specifically, $\alpha$ represents the path loss factor.

Figure 1.

Illustration of the system model.

In the underlay CR network, the secondary users are allowed to use the frequency band of $\mathrm{PU}$ if the interference caused by secondary users is tolerable. Thus, the transmit power at $S$ and $R_n$ is restricted as ²⁵

$$P_S=min\{\frac{I}{{\left|h_{\mathrm{SP}}\right|}^2},P\}$$

(1)

$$P_R=min\{\frac{I}{{\left|h_{R_{n}P}\right|}^2},\beta P\}$$

(2)

where $I$ denotes the ITC at $\mathrm{PU}$ while $P$ indicates the maximum average transmit power of $S$ . Positive constant $\beta$ indicates the relationship between the maximum average transmit power of $S$ and the maximum average transmit power of $R_n$ . In addition, $P_I$ represents the interference from $\mathrm{PU}$ to the secondary network, which can be seen as additive white Gaussian noise (AWGN) with $\mathrm{CN}(0,\eta {\sigma }^2)$ . ¹⁹ What’s more, it is assumed that all secondary nodes obtain the same $P_I$ for simplicity.

RS schemes

1. Partial RS scheme: the best relay is selected based on the channels of the first hop ¹⁰

   $$R_b=\arg \underset{n=1,2,\dots ,N}{max}{\left|h_{S{R}_n}\right|}^2$$

   (3)

2. Max–min RS scheme: this scheme considers the channels of the first hop and the second hop which can be formulated as follows ⁷

   $$\begin{array}{l}{R}_b=\arg \underset{n=1,2,\dots ,N}{\max }\\ \left\{\min \right\{{\left|h_{S{R}_n}\right|}^2,{\left|g_{R_n{D}_1}\right|}^2,{\left|g_{R_n{D}_2}\right|}^2\}\}\end{array}$$

   (4)

3. Two-stage RS scheme: the first stage is to select the relays that can guarantee user 1’s targeted data rate ⁷

$$\begin{array}{c}{S}_r=\{n:1\le n\le N,\frac{1}{2}\log (1+\frac{a_1{\rho }_S{\left|h_{S{R}_n}\right|}^2}{a_2{\rho }_S{\left|h_{S{R}_n}\right|}^2+\eta +1})\ge u_1\hfill \\ \frac{1}{2}\log (1+\frac{a_1{\rho }_R{\left|g_{R_n{D}_1}\right|}^2}{a_2{\rho }_R{\left|g_{R_n{D}_1}\right|}^2+\eta +1})\ge u_1\hfill \\ \frac{1}{2}\log (1+\frac{a_1{\rho }_R{\left|g_{R_n{D}_2}\right|}^2}{a_2{\rho }_R{\left|g_{R_n{D}_2}\right|}^2+\eta +1})\ge u_1\}\hfill \end{array}$$

(5)

In other words, the QoS requirements of $D_1$ can always be met in this phase.

On the basis of satisfying the above condition, the second stage is to choose the relay which can maximize user 2’s data rate

$$\begin{array}{c}{R}_b=\arg \underset{n}{max}\{min\{\log (1+\frac{a_2{\rho }_S{\left|h_{S{R}_n}\right|}^2}{\eta +1}),\hfill \\ \log (1+\frac{a_2{\rho }_R{\left|g_{R_n{D}_2}\right|}^2}{\eta +1})\},n\in S_r\}\hfill \end{array}$$

(6)

Based on the principle of NOMA, the communication from $S$ to $D_i$ consists of two time periods. During the first time slot, $S$ transmits the superimposed signal $\sum _{i=1}^2\sqrt{P_S{a}_i}{x}_i$ to the selected relay, where $x_i$ denotes the signal for $D_i$ and $a_i$ represents the PA factor of $x_i$ with $a_1+a_2=1$ . Considering the QoS requirements and throughput fairness among individual users, it is required that $a_1>a_2$ . Therefore, the received signal at $R_n$ can be written as

$$y_{R_n}=h_{S{R}_n}\sqrt{P_S}(\sqrt{a_1}{x}_1+\sqrt{a_2}{x}_2)+P_I+w_R$$

(7)

where $w_{(·)}$ denotes the AWGN at each receive node with zero mean and variance ${\sigma }^2$ . According to the principle of SIC, the selected relay $R_b$ first decodes and removes $x_1$ by treating $x_2$ as noise, and then will decode $x_2$ .

The transmit SNR at $S$ and $R_n$ are defined as ${\rho }_S=P_S/{\sigma }^2$ and ${\rho }_R=P_R/{\sigma }^2$ , respectively. Random variables $X={\rho }_S{\left|h_{S{R}_b}\right|}^2$ and $Y_i={\rho }_R{\left|g_{R_n{D}_i}\right|}^2$ denote accordingly instantaneous SNRs at $R_n$ and $D_i$ . Thus, the signal-to-interference-plus-noise ratio (SINR) for $R_b$ to decode $x_1$ and $x_2$ can be given by

$${\gamma }_{R,1}=\frac{a_1{\rho }_S{\left|h_{S{R}_b}\right|}^2}{a_2{\rho }_S{\left|h_{S{R}_b}\right|}^2+\eta +1}=\frac{a_{1}X}{a_{2}X+\eta +1}$$

(8)

$${\gamma }_{R,2}=\frac{a_2{\rho }_S{\left|h_{S{R}_b}\right|}^2}{\eta +1}=\frac{a_{2}X}{\eta +1}$$

(9)

During the second time slot, $R_b$ forwards the superimposed signal $\sum _{i=1}^2\sqrt{P_R{a}_i}{\tilde{x}}_i$ to $D_i$ . Therefore, $D_i$ observes

$$y_i=g_{R_b{D}_i}\sqrt{P_R}(\sqrt{a_1}{\tilde{x}}_1+\sqrt{a_2}{\tilde{x}}_2)+P_I+w_{D_i}$$

(10)

where ${\tilde{x}}_i$ is the detected message of $R_b$ . As discussed earlier, the SINR for $D_1$ to decode $x_1$ can be written as

$${\gamma }_{D1,1}=\frac{a_1{\rho }_R{\left|g_{R_b{D}_1}\right|}^2}{a_2{\rho }_R{\left|g_{R_b{D}_1}\right|}^2+\eta +1}=\frac{a_1{Y}_1}{a_2{Y}_1+\eta +1}$$

(11)

Correspondingly, the SINR for $D_2$ to decode $x_1$ and $x_2$ can be given by

$${\gamma }_{D2,1}=\frac{a_1{\rho }_R{\left|g_{R_b{D}_2}\right|}^2}{a_2{\rho }_R{\left|g_{R_b{D}_2}\right|}^2+\eta +1}=\frac{a_1{Y}_2}{a_2{Y}_2+\eta +1}$$

(12)

$${\gamma }_{D2,2}=\frac{a_2{\rho }_R{\left|g_{R_b{D}_2}\right|}^2}{\eta +1}=\frac{a_2{Y}_2}{\eta +1}$$

(13)

It should be noted that the required condition to decode $x_i$ at each receive node is $(1/2)\log (1+{\gamma }_{(·),i})\ge u_i$ , that is, ${\gamma }_{(·),i}\ge {\epsilon }_i=2^{2u_i-1}$ , where ${\epsilon }_i$ is the receive SNR threshold and $u_i$ denotes the predefined targeted data rate for $D_i$ .

Formulation analysis in the case of perfect SIC

In this part, we consider the scenario that SIC is perfect at $R_n$ and $D_2$ , and the closed-form expressions of the OP for the system under the three RS schemes are derived.

Formulation analysis of partial RS

This section focuses on the outage performance of the partial RS scheme. The cumulative distribution function (CDF) of the random variable $X$ is given by

$$\begin{array}{c}\multicolumn{1}{c}{F_X\left(\epsilon \right)=P(\rho {\left|h_{S{R}_b}\right|}^2<\epsilon ,\rho <{\theta }_1)}\\ +P({\theta }_1{\left|h_{S{R}_b}\right|}^2<\epsilon ,\rho >{\theta }_1)\hfill \\ \multicolumn{1}{c}{=P({\left|h_{S{R}_b}\right|}^2<\frac{\epsilon }{\rho },{\left|h_{\mathrm{SP}}\right|}^2<\frac{{\rho }_I}{\rho })}\\ +P(\frac{{\left|h_{S{R}_b}\right|}^2}{{\left|h_{\mathrm{SP}}\right|}^2}<\frac{\epsilon }{{\rho }_I},{\left|h_{\mathrm{SP}}\right|}^2>\frac{{\rho }_I}{\rho })\hfill \\ \multicolumn{1}{c}{=\underset{0}{\overset{\frac{\epsilon }{\rho }}{\int }}{f}_{{\left|h_{S{R}_b}\right|}^2}\left(x\right)\mathrm{dx}\underset{0}{\overset{\frac{{\rho }_I}{\rho }}{\int }}{f}_{{\left|h_{\mathrm{SP}}\right|}^2}\left(y\right)\mathrm{dy}}\\ +\underset{\frac{{\rho }_I}{\rho }}{\overset{\infty }{\int }}{f}_{{\left|h_{\mathrm{SP}}\right|}^2}\left(y\right)\mathrm{dy}\underset{0}{\overset{\frac{\epsilon y}{{\rho }_I}}{\int }}{f}_{{\left|h_{S{R}_b}\right|}^2}\left(x\right)\mathrm{dx}\hfill \end{array}$$

(14)

where $\rho =P/{\sigma }^2$ , ${\theta }_1={\rho }_I/{\left|h_{\mathrm{SP}}\right|}^2$ , ${\rho }_I=I/{\sigma }^2$ .

To derive the close-form expression in equation (14), it is required to get the probability density function (PDF) of ${\left|h_{S{R}_b}\right|}^2$ .

First, the PDF of the ordered variable ²⁶ ${\left|h_m\right|}^2$ is given by

$$\begin{array}{c}{f}_{{\left|h_m\right|}^2}\left(x\right)=\frac{N!}{(N-m)!(m-1)!}{f}_{{\left|{\tilde{h}}_i\right|}^2}\left(x\right){\left[F_{{\left|{\tilde{h}}_i\right|}^2}\left(x\right)\right]}^{m-1}\hfill \\ \multicolumn{1}{c}{\times {[1-F_{{\left|{\tilde{h}}_i\right|}^2}\left(x\right)]}^{N-m}}\end{array}$$

(15)

where ${\left|{\tilde{h}}_i\right|}^2$ is the unordered variable and its PDF is expressed as

$$f_{{\left|{\tilde{h}}_i\right|}^2}\left(x\right)=\frac{1}{{\lambda }_0}{e}^{-\frac{x}{{\lambda }_0}}$$

(16)

As $R_b$ is selected with the largest ${\left|h_{S{R}_n}\right|}^2,n\in \{1,2,...,N\}$ , the PDF of ${\left|h_{S{R}_b}\right|}^2$ can be derived by assigning $m$ in equation (14) to $N$

$$\begin{array}{c}\multicolumn{1}{c}{f_{{\left|h_{S{R}_b}\right|}^2}\left(x\right)=N{\left[F_{{\left|{\tilde{h}}_i\right|}^2}\left(x\right)\right]}^{N-1}{f}_{{\left|{\tilde{h}}_i\right|}^2}\left(x\right)}\\ \multicolumn{1}{c}{=\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}\frac{n}{{\lambda }_0}{e}^{-\frac{\mathrm{nx}}{{\lambda }_0}}}\end{array}$$

(17)

Next, the PDF of ${\left|h_{\mathrm{SP}}\right|}^2$ is given by

$$f_{{\left|h_{\mathrm{SP}}\right|}^2}\left(x\right)=\frac{1}{{\lambda }_{\mathrm{SP}}}{e}^{-\frac{x}{{\lambda }_{\mathrm{SP}}}}$$

(18)

Then, by substituting equations (17) and (18) into equation (14), the CDF of random variable $X$ can be rewritten as

$$\begin{array}{c}{F}_X\left(\epsilon \right)=1-\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}\hfill \\ \multicolumn{1}{c}{\times (e^{-\frac{n\epsilon }{{\lambda }_0\rho }}-\frac{n\epsilon {\lambda }_{\mathrm{SP}}}{n\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{e}^{-\frac{n\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{{\lambda }_0{\lambda }_{\mathrm{SP}}\rho }})}\end{array}$$

(19)

Similar to equation (14), the CDF of random variable $Y_i$ is given by

$$\begin{array}{c}{F}_{Y_i}\left(\epsilon \right)=P(\beta \rho {\left|g_{R_b{D}_i}\right|}^2<\epsilon ,\beta \rho <{\theta }_2)\hfill \\ \multicolumn{1}{c}{+P({\theta }_2{\left|g_{R_b{D}_i}\right|}^2<\epsilon ,\beta \rho >{\theta }_2)}\\ \multicolumn{1}{c}{=\underset{0}{\overset{\frac{\epsilon }{\beta \rho }}{\int }}{f}_{{\left|g_{R_b{D}_i}\right|}^2}\left(x\right)\mathrm{dx}\underset{0}{\overset{\frac{{\rho }_I}{\beta \rho }}{\int }}{f}_{{\left|h_{R_{b}P}\right|}^2}\left(y\right)\mathrm{dy}}\\ \multicolumn{1}{c}{+\underset{\frac{{\rho }_I}{\beta \rho }}{\overset{\infty }{\int }}{f}_{{\left|h_{R_{b}P}\right|}^2}\left(y\right)\mathrm{dy}\underset{0}{\overset{\frac{\epsilon y}{{\rho }_I}}{\int }}{f}_{{\left|g_{R_b{D}_i}\right|}^2}\left(x\right)\mathrm{dx}}\\ \multicolumn{1}{c}{=1-e^{-\frac{\epsilon }{{\lambda }_i\beta \rho }}+\frac{\epsilon {\lambda }_{\mathrm{RP}}}{\epsilon {\lambda }_{\mathrm{RP}}+{\lambda }_i{\rho }_I}{e}^{-\frac{\epsilon {\lambda }_{\mathrm{RP}}+{\lambda }_i{\rho }_I}{{\lambda }_i{\lambda }_{\mathrm{RP}}\beta \rho }}}\end{array}$$

(20)

where ${\theta }_2={\rho }_I/{\left|h_{R_{b}P}\right|}^2$ .

The OP of the system can be expressed as

$$\begin{array}{c}\mathrm{OP}=1-P({\gamma }_{R,1}>{\epsilon }_1,{\gamma }_{R,2}>{\epsilon }_2,{\gamma }_{D1,1}>{\epsilon }_1,{\gamma }_{D2,1}>{\epsilon }_1\hfill \\ ,{\gamma }_{D2,2}>{\epsilon }_2)\hfill \\ {=}_{\left(a\right)}1-P(X>{\varphi }_1,Y_1>{\xi }_1,Y_2>{\varphi }_1)\hfill \\ =1-[1-F_X\left({\varphi }_1\right)][1-F_{Y_1}\left({\xi }_1\right)][1-F_{Y_2}\left({\varphi }_1\right)]\hfill \end{array}$$

(21)

where ${\xi }_1={\epsilon }_1(1+\eta )/a_1-a_2{\epsilon }_1$ , ${\xi }_2={\epsilon }_2(1+\eta )/a_2$ and ${\varphi }_1=max\{{\xi }_1,{\xi }_2\}$ . In addition, step (a) requires that $a_1>a_2{\epsilon }_1$ , otherwise $\mathrm{OP}=1$ . Finally, by substituting equations (19) and (20) into equation (21), the OP of the system can be rewritten as

$$\begin{array}{c}\mathrm{OP}=1-\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}\hfill \\ \multicolumn{1}{c}{\times (e^{-\frac{n{\varphi }_1}{{\lambda }_0\rho }}-\frac{n{\varphi }_1{\lambda }_{\mathrm{SP}}}{n{\varphi }_1{\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{e}^{-\frac{n{\varphi }_1{\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{{\lambda }_0{\lambda }_{\mathrm{SP}}\rho }})}\\ \multicolumn{1}{c}{\times (e^{-\frac{{\xi }_1}{{\lambda }_1\beta \rho }}-\frac{{\xi }_1{\lambda }_{\mathrm{RP}}}{{\xi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_1{\rho }_I}{e}^{-\frac{{\xi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_1{\rho }_I}{{\lambda }_1{\lambda }_{\mathrm{RP}}\beta \rho }})}\\ \multicolumn{1}{c}{\times (e^{-\frac{{\varphi }_1}{{\lambda }_2\beta \rho }}-\frac{{\varphi }_1{\lambda }_{\mathrm{RP}}}{{\varphi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_2{\rho }_I}{e}^{-\frac{{\varphi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_2{\rho }_I}{{\lambda }_2{\lambda }_{\mathrm{RP}}\beta \rho }})}\end{array}$$

(22)

Asymptotic analysis with $\rho \to \infty$ for partial RS

The exponential function can be approximated as

$$\underset{x\to 0}{lim}{e}^{-x}\approx 1-x$$

(23)

Therefore, the following simplifications can be performed

$$e^{-\frac{n\epsilon }{{\lambda }_0\rho }}\approx 1-\frac{n\epsilon }{{\lambda }_0\rho }$$

(24)

$$e^{-\frac{n\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{{\lambda }_0{\lambda }_{\mathrm{SP}}\rho }}\approx 1-\frac{n\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{{\lambda }_0{\lambda }_{\mathrm{SP}}\rho }$$

(25)

The CDF of random variable X in equation (19) can be approximated as

$$F_X\left(\epsilon \right)\approx \sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}\frac{n\epsilon {\lambda }_{\mathrm{SP}}}{n\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}$$

(26)

The CDF of random variable $Y_i$ in equation (20) can be approximated as

$$F_{Y_i}\left(\epsilon \right)\approx \frac{\epsilon {\lambda }_{\mathrm{RP}}}{\epsilon {\lambda }_{\mathrm{RP}}+{\lambda }_i{\rho }_I}$$

(27)

Using the above approximations, the asymptotic expression for the OP of the system is given by

$$\begin{array}{c}\underset{\rho \to \infty }{lim}\mathrm{OP}=1-\hfill \\ [1-\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}\frac{n{\varphi }_1{\lambda }_{\mathrm{SP}}}{n{\varphi }_1{\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}]\hfill \\ \times (1-\frac{{\xi }_1{\lambda }_{\mathrm{RP}}}{{\xi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_1{\rho }_I})·(1-\frac{{\varphi }_1{\lambda }_{\mathrm{RP}}}{{\varphi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_2{\rho }_I})\hfill \end{array}$$

(28)

Asymptotic analysis with $I\to \infty$ for partial RS

When the ITC $I$ goes to infinity, the asymptotic expression for the OP of the system is as follows

$$\begin{array}{c}\underset{I\to \infty }{lim}\mathrm{OP}=1-\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}{e}^{-\frac{n{\varphi }_1}{{\lambda }_0\rho }}\hfill \\ \times e^{-\frac{{\xi }_1}{{\lambda }_1\beta \rho }}{e}^{-\frac{{\varphi }_1}{{\lambda }_2\beta \rho }}\hfill \end{array}.$$

(29)

Formulation analysis of max–min RS

The PDF of ${\left|h_{S{R}_b}\right|}^2$ based on the max–min RS scheme can be expressed as (see Appendix)

$$\begin{array}{c}{f}_{{\left|h_{S{R}_b}\right|}^2}^M\left(x\right)=\sum _{i=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{i}\end{array}\right){(-1)}^{i-1}i\hfill \\ \times [\frac{1}{{\lambda }_0}{e}^{-\frac{\mathrm{ix}}{{\tau }_Z}}+\frac{{\tau }_Z}{{\tau }_0(i{\lambda }_0-{\tau }_Z)}(e^{-\frac{x}{{\lambda }_0}}-e^{-\frac{\mathrm{ix}}{{\tau }_Z}})]\hfill \end{array},$$

(30)

where ${\tau }_Z={\lambda }_0{\lambda }_1{\lambda }_2/{\lambda }_0{\lambda }_1+{\lambda }_0{\lambda }_2+{\lambda }_1{\lambda }_2$ , ${\tau }_0={\lambda }_1{\lambda }_2/{\lambda }_1+{\lambda }_2$ . By replacing $f_{{\left|h_{S{R}_b}\right|}^2}(x)$ with $f_{{\left|h_{S{R}_b}\right|}^2}^M(x)$ in equation (14), the CDF of random variable $X$ can be calculated as in equation (31)

$$\begin{array}{l}{F}_X^M\left(x\right)=\sum _{i=1}^N\left(\begin{array}{c}N\\ i\end{array}\right){\left(-1\right)}^{i-1}\{\frac{i{\lambda }_0{\tau }_Z}{{\tau }_0\left(i{\lambda }_0-{\tau }_Z\right)}\left(1-e^{-\frac{x}{{\lambda }_0\rho }}\right)\left(1-e^{-\frac{{\rho }_I}{{\lambda }_{SP}\rho }}\right)\\ +\frac{{\tau }_0{\tau }_Z\left(i{\lambda }_0-{\tau }_Z\right)-{\lambda }_0{{\tau }_Z}^2}{{\lambda }_0{\tau }_0\left(i{\lambda }_0-{\tau }_Z\right)}\left(1-e^{-\frac{ix}{{\tau }_Z\rho }}\right)\left(1-e^{-\frac{{\rho }_I}{{\lambda }_{SP}\rho }}\right)\\ +\frac{i{\lambda }_0{\tau }_Z}{{\lambda }_{SP}{\tau }_0\left(i{\lambda }_0-{\tau }_Z\right)}\left({\lambda }_{SP}{e}^{-\frac{{\rho }_I}{{\lambda }_{SP}\rho }}-\frac{{\lambda }_{SP}{\lambda }_0{\rho }_I}{{\lambda }_0{\rho }_I+x{\lambda }_{SP}}{e}^{-\frac{{\lambda }_0{\rho }_I+x{\lambda }_{SP}}{{\lambda }_{SP}{\lambda }_0\rho }}\right)\\ +\frac{{\tau }_0{\tau }_Z\left(i{\lambda }_0-{\tau }_Z\right)-{\lambda }_0{{\tau }_Z}^2}{{\lambda }_0{\tau }_0\left(i{\lambda }_0-{\tau }_Z\right)}{e}^{-\frac{{\rho }_I}{{\lambda }_{SP}\rho }}\\ -\frac{{\tau }_0{{\tau }_Z}^2{\rho }_I\left(i{\lambda }_0-{\tau }_Z\right)-{\lambda }_0{{\tau }_Z}^3{\rho }_I}{{\lambda }_0{\tau }_0\left(i{\lambda }_0-{\tau }_Z\right)\left({\tau }_Z{\rho }_I+ix{\lambda }_{SP}\right)}{e}^{-\frac{{\tau }_Z{\rho }_I+ix{\lambda }_{SP}}{{\lambda }_{SP}{\tau }_Z\rho }}\}\end{array}$$

(31)

In the same way, the CDF of $Y_i$ is given by equation (32), where ${\tau }_1={\lambda }_0{\lambda }_2/{\lambda }_0+{\lambda }_2$ , ${\tau }_2={\lambda }_0{\lambda }_1/{\lambda }_0+{\lambda }_1$

$$\begin{array}{l}{F}_{Y_i}^M\left(x\right)=\sum _{j=1}^N\left(\begin{array}{l}N\\ j\end{array}\right){\left(-1\right)}^{j-1}\{\frac{j{\lambda }_i{\tau }_Z}{{\tau }_i\left(j{\lambda }_i-{\tau }_Z\right)}\left(1-e^{-\frac{x}{{\lambda }_i\beta \rho }}\right)\left(1-e^{-\frac{{\rho }_I}{{\lambda }_{RP}\beta \rho }}\right)\\ +\frac{{\tau }_i{\tau }_Z\left(j{\lambda }_i-{\tau }_Z\right)-{\lambda }_i{{\tau }_Z}^2}{{\lambda }_i{\tau }_i\left(j{\lambda }_i-{\tau }_Z\right)}\\ \left(1-e^{-\frac{jx}{{\tau }_Z\beta \rho }}\right)\left(1-e^{-\frac{{\rho }_I}{{\lambda }_{RP}\beta \rho }}\right)\\ +\frac{j{\lambda }_i{\tau }_Z}{{\lambda }_{RP}{\tau }_i\left(j{\lambda }_i-{\tau }_Z\right)}\\ \left({\lambda }_{RP}{e}^{-\frac{{\rho }_I}{{\lambda }_{RP}\beta \rho }}-\frac{{\lambda }_{RP}{\lambda }_i{\rho }_I}{{\lambda }_i{\rho }_I+x{\lambda }_{RP}}{e}^{-\frac{{\lambda }_i{\rho }_I+x{\lambda }_{RP}}{{\lambda }_{RP}{\lambda }_i\beta \rho }}\right)\\ +\frac{{\tau }_i{\tau }_Z\left(j{\lambda }_i-{\tau }_Z\right)-{\lambda }_i{{\tau }_Z}^2}{{\lambda }_i{\tau }_i\left(j{\lambda }_i-{\tau }_Z\right)}{e}^{-\frac{{\rho }_I}{{\lambda }_{RP}\beta \rho }}\\ -\frac{{\tau }_i{{\tau }_Z}^2{\rho }_I\left(j{\lambda }_i-{\tau }_Z\right)-{\lambda }_i{{\tau }_Z}^3{\rho }_I}{{\lambda }_i{\tau }_i\left(j{\lambda }_i-{\tau }_Z\right)\left({\tau }_Z{\rho }_I+jx{\lambda }_{RP}\right)}{e}^{-\frac{{\tau }_Z{\rho }_I+jx{\lambda }_{RP}}{{\lambda }_{RP}{\tau }_Z\beta \rho }}\}\end{array}$$

(32)

The OP of the system can be expressed as

$$\begin{array}{c}\mathrm{O}{\mathrm{P}}^M=1-P({\gamma }_{R,1}>{\epsilon }_1,{\gamma }_{R,2}>{\epsilon }_2,{\gamma }_{D1,1}>{\epsilon }_1\hfill \\ \multicolumn{1}{c}{,{\gamma }_{D2,1}>{\epsilon }_1,{\gamma }_{D2,2}>{\epsilon }_2)}\\ \multicolumn{1}{c}{=1-[1-F_X^M\left({\varphi }_1\right)][1-F_{Y_1}^M\left({\xi }_1\right)][1-F_{Y_2}^M\left({\varphi }_1\right)]}\end{array}$$

(33)

By substituting equations (31) and (32) into equation (33), the OP of the system can be obtained.

Asymptotic analysis with $\rho \to \infty$ for max–min RS

Based on equation (23), the asymptotic expressions for $F_X^M(x)$ and $F_{Y_i}^M(x)$ are given by equations (34) and (35), respectively, which is

$$F_X^{M,\rho }\left(x\right)=\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{i}\end{array}\right){(-1)}^{i-1}[\frac{i{\lambda }_0{\tau }_Z{\lambda }_{\mathrm{SP}}x}{{\tau }_0(i{\lambda }_0-{\tau }_Z)({\lambda }_0{\rho }_I+x{\lambda }_{\mathrm{SP}})}+\frac{{\tau }_0{\tau }_Z(i{\lambda }_0-{\tau }_Z)-{\lambda }_0{\tau }_Z^2}{{\lambda }_0{\tau }_0(i{\lambda }_0-{\tau }_Z)}-\frac{{\tau }_0{\tau }_Z^2{\rho }_I(i{\lambda }_0-{\tau }_Z)-{\lambda }_0{\tau }_Z^3{\rho }_I}{{\lambda }_0{\tau }_0(i{\lambda }_0-{\tau }_Z)({\tau }_Z{\rho }_I+\mathrm{ix}{\lambda }_{\mathrm{SP}})}]$$

(34)

$$F_{Y_i}^{M,\rho }\left(x\right)=\sum _{j=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{j}\end{array}\right){(-1)}^{j-1}[\frac{j{\lambda }_i{\tau }_Z{\lambda }_{\mathrm{RP}}x}{{\tau }_i(j{\lambda }_i-{\tau }_Z)({\lambda }_i{\rho }_I+x{\lambda }_{\mathrm{RP}})}+\frac{{\tau }_i{\tau }_Z(j{\lambda }_i-{\tau }_Z)-{\lambda }_i{\tau }_Z^2}{{\lambda }_i{\tau }_i(j{\lambda }_i-{\tau }_Z)}-\frac{{\tau }_i{\tau }_Z^2{\rho }_I(j{\lambda }_i-{\tau }_Z)-{\lambda }_i{\tau }_Z^3{\rho }_I}{{\lambda }_i{\tau }_i(j{\lambda }_i-{\tau }_Z)({\tau }_Z{\rho }_I+\mathrm{jx}{\lambda }_{\mathrm{RP}})}]$$

(35)

By substituting equations (34) and (35) into equation (33), the asymptotic expression for the OP of the system is given by

$$\begin{array}{c}\underset{\rho \to \infty }{lim}\mathrm{O}{\mathrm{P}}^M=1-[1-F_X^{M,\rho }\left({\varphi }_1\right)][1-F_{Y_1}^{M,\rho }\left({\xi }_1\right)]\hfill \\ \multicolumn{1}{c}{\times [1-F_{Y_2}^{M,\rho }\left({\varphi }_1\right)]}\end{array}$$

(36)

Asymptotic analysis with $I\to \infty$ for max–min RS

When the ITC $I$ goes to infinity, the asymptotic expressions for $F_X^M(x)$ and $F_{Y_i}^M(x)$ are given by equations (37) and (38), respectively, which are

$$\begin{array}{c}\multicolumn{1}{c}{F_X^{M,I}\left(x\right)=\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{i}\end{array}\right){(-1)}^{i-1}}\\ \multicolumn{1}{c}{[\frac{i{\lambda }_0{\tau }_Z}{{\tau }_0(i{\lambda }_0-{\tau }_Z)}(1-e^{-\frac{x}{{\lambda }_0\rho }})+\frac{{\tau }_0{\tau }_Z(i{\lambda }_0-{\tau }_Z)-{\lambda }_0{\tau }_Z^2}{{\lambda }_0{\tau }_0(i{\lambda }_0-{\tau }_Z)}(1-e^{-\frac{\mathrm{ix}}{{\tau }_Z\rho }})]}\end{array}$$

(37)

$$\begin{array}{c}\multicolumn{1}{c}{F_{Y_i}^{M,I}\left(x\right)=\sum _{j=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{j}\end{array}\right){(-1)}^{j-1}}\\ \multicolumn{1}{c}{[\frac{j{\lambda }_i{\tau }_Z}{{\tau }_i(j{\lambda }_i-{\tau }_Z)}(1-e^{-\frac{x}{{\lambda }_i\beta \rho }})+\frac{{\tau }_i{\tau }_Z(j{\lambda }_i-{\tau }_Z)-{\lambda }_i{\tau }_Z^2}{{\lambda }_i{\tau }_i(j{\lambda }_i-{\tau }_Z)}(1-e^{-\frac{\mathrm{jx}}{{\tau }_Z\beta \rho }})]}\end{array}$$

(38)

Therefore, the asymptotic expression for the OP of the system is given by

$$\begin{array}{c}\underset{I\to \infty }{lim}\mathrm{O}{\mathrm{P}}^M=1-[1-F_X^{M,I}\left({\varphi }_1\right)][1-F_{Y_1}^{M,I}\left({\xi }_1\right)]\hfill \\ \multicolumn{1}{c}{\times [1-F_{Y_2}^{M,I}\left({\varphi }_1\right)]}\end{array}$$

(39)

Formulation analysis of two-stage RS

The CDF of ${\left|h_{S{R}_b}\right|}^2$ based on the two-stage RS scheme is given by

$$F_{{\left|h_{S{R}_b}\right|}^2}^T\left(x\right)=1-e^{-\frac{x}{{\lambda }_0}}$$

(40)

Accordingly, the CDF of random variable $X$ in equation (14) can be rewritten as

$$F_X^T\left(\epsilon \right)=1-e^{-\frac{\epsilon }{{\lambda }_0\rho }}+\frac{\epsilon {\lambda }_{\mathrm{SP}}}{\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{e}^{-\frac{\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{{\lambda }_0{\lambda }_{\mathrm{SP}}\rho }}$$

(41)

Defining $O_1$ as the event that the goal of the first stage cannot be achieved, and $O_2$ as the event that the goal of the second stage cannot be achieved. Thus, the OP of the system can be expressed as

$$\mathrm{O}{\mathrm{P}}^T=P\left(O_1\right)+P\left(O_2\right)$$

(42)

Term $P(O_1)$ is given by

$$\begin{array}{c}P\left(O_1\right)=P(\left|S_r\right|=0)\hfill \\ =\underset{n=1}{\overset{N}{\Pi }}[1-P({\gamma }_{R,1}\ge {\epsilon }_1,{\gamma }_{D1,1}\ge {\epsilon }_1,{\gamma }_{D2,1}\ge {\epsilon }_1)]\hfill \\ =\underset{n=1}{\overset{N}{\Pi }}[1-(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))]\hfill \end{array}$$

(43)

where $\left|S_r\right|$ denotes the size of $S_r$ .

Term $P(O_2)$ can be expressed as

$$\begin{array}{c}P\left(O_2\right)=P(max\{min\{\log (1+\frac{a_2{\rho }_S{\left|h_{S{R}_n}\right|}^2}{\eta +1})\hfill \\ \multicolumn{1}{c}{,\log (1+\frac{a_2{\rho }_R{\left|g_{R_n{D}_2}\right|}^2}{\eta +1})\}\}<2u_2,\left|S_r\right|>0)}\\ \multicolumn{1}{c}{=P(max\{min\{\frac{a_2{\rho }_S{\left|h_{S{R}_n}\right|}^2}{\eta +1}}\\ \multicolumn{1}{c}{,\frac{a_2{\rho }_R{\left|g_{R_n{D}_2}\right|}^2}{\eta +1}\}\}<{\epsilon }_2,\left|S_r\right|>0)}\end{array}$$

(44)

Defining

$$\begin{array}{c}{t}_n=min\{\frac{a_2{\rho }_S{\left|h_{S{R}_n}\right|}^2}{\eta +1},\frac{a_2{\rho }_R{\left|g_{R_n{D}_2}\right|}^2}{\eta +1}\}\hfill \\ =min\{{\gamma }_{R,2},{\gamma }_{D2,2}\}\hfill \end{array}$$

(45)

and

$$t_b=max\{t_i,i\in S_r\}$$

(46)

The CDF of $t_n$ can be calculated as follows

$$\begin{array}{c}\multicolumn{1}{c}{F_{t_n}\left({\epsilon }_2\right)=P(min\{{\gamma }_{R,2},{\gamma }_{D2,2}\}<{\epsilon }_2|\left|S_r\right|>0)}\\ \multicolumn{1}{c}{=P(min\{{\gamma }_{R,2},{\gamma }_{D2,2}\}<{\epsilon }_2|X>{\xi }_1,Y_1>{\xi }_1,Y_2>{\xi }_1)}\\ \multicolumn{1}{c}{=1-\frac{P(min\{{\gamma }_{R,2},{\gamma }_{D2,2}\}>{\epsilon }_2,X>{\xi }_1,Y_1>{\xi }_1,Y_2>{\xi }_1)}{P(X>{\xi }_1,Y_1>{\xi }_1,Y_2>{\xi }_1)}}\\ \multicolumn{1}{c}{=1-\frac{P(X>{\varphi }_1,Y_1>{\xi }_1,Y_2>{\varphi }_1)}{P(X>{\xi }_1,Y_1>{\xi }_1,Y_2>{\xi }_1)}}\\ \multicolumn{1}{c}{=1-\frac{[1-F_X^T\left({\varphi }_1\right)][1-F_{Y_2}\left({\varphi }_1\right)]}{[1-F_X^T\left({\xi }_1\right)][1-F_{Y_2}\left({\xi }_1\right)]}}\end{array}$$

(47)

Thus, $P(O_2)$ can be calculated as

$$\begin{array}{c}P\left(O_2\right)=\sum _{l=1}^{N}P(t_b<{\epsilon }_2,\left|S_r\right|=l)\hfill \\ \multicolumn{1}{c}{=\sum _{l=1}^{N}P(t_b<{\epsilon }_2|\left|S_r\right|=l)P(\left|S_r\right|=l)}\\ \multicolumn{1}{c}{=\sum _{l=1}^N{\left(F_{t_n}\left({\epsilon }_2\right)\right)}^{l}P(\left|S_r\right|=l)}\end{array}$$

(48)

In the following, the probability that $\left|S_r\right|$ equals to $l$ is given by

$$\begin{array}{c}P(\left|S_r\right|=l)=\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{l}\end{array}\right)\underset{n=1}{\overset{N-l}{\Pi }}[1-P({\gamma }_{R,1}>{\epsilon }_1)P({\gamma }_{D1,1}>{\epsilon }_1)\hfill \\ \multicolumn{1}{c}{\times P({\gamma }_{D2,1}>{\epsilon }_1)]}\\ \multicolumn{1}{c}{\times \underset{n=N-l+1}{\overset{N}{\Pi }}[P({\gamma }_{R,1}>{\epsilon }_1)P({\gamma }_{D1,1}>{\epsilon }_1)}\\ \multicolumn{1}{c}{\times P({\gamma }_{D2,1}>{\epsilon }_1)]}\\ \multicolumn{1}{c}{=\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{l}\end{array}\right)\underset{n=1}{\overset{N-l}{\Pi }}[1-(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))}\\ \multicolumn{1}{c}{\times (1-F_{Y_2}\left({\xi }_1\right))]}\\ \multicolumn{1}{c}{\times \underset{n=N-l+1}{\overset{N}{\Pi }}[(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))}\\ \multicolumn{1}{c}{\times (1-F_{Y_2}\left({\xi }_1\right))]}\end{array}$$

(49)

Finally, from equations (42)–(49), the OP of the system can be derived as in equation (50), which is

$$\begin{array}{c}\mathrm{O}{\mathrm{P}}^T=\underset{n=1}{\overset{N}{\Pi }}[1-(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))]\hfill \\ \multicolumn{1}{c}{+\sum _{l=1}^N{[1-\frac{(1-F_X^T\left({\varphi }_1\right))(1-F_{Y_2}\left({\varphi }_1\right))}{(1-F_X^T\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))}]}^l\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{l}\end{array}\right)}\\ \multicolumn{1}{c}{{[1-(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))]}^{N-l}}\\ \multicolumn{1}{c}{\times {\left[(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))\right]}^l}\end{array}$$

(50)

Asymptotic analysis with $\rho \to \infty$ for two-stage RS

Based on equation (23), the asymptotic approximation for $F_X^T(\epsilon )$ is given by

$$F_X^T\left(\epsilon \right)\approx \frac{\epsilon {\lambda }_{\mathrm{SP}}}{\epsilon {\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}$$

(51)

The asymptotic approximation for $F_{Y_i}(\epsilon )$ is given by

$$F_{Y_i}\left(\epsilon \right)\approx \frac{\epsilon {\lambda }_{\mathrm{RP}}}{\epsilon {\lambda }_{\mathrm{RP}}+{\lambda }_i{\rho }_I}$$

(52)

Then, by substituting equations (51) and (52) into equation (50), the asymptotic expression for the OP of the system is obtained.

Asymptotic analysis with $I\to \infty$ for two-stage RS

The asymptotic approximation for $F_X^T(\epsilon )$ is given by

$$F_X^T\left(\epsilon \right)\approx 1-e^{-\frac{\epsilon }{{\lambda }_0\rho }}$$

(53)

The asymptotic approximation for $F_{Y_i}(\epsilon )$ is given by

$$F_{Y_i}\left(\epsilon \right)\approx 1-e^{-\frac{\epsilon }{{\lambda }_i\beta \rho }}$$

(54)

By substituting equations (53) and (54) into equation (50), the asymptotic expression for $O{P}^T$ is obtained.

Formulation analysis in the case of imperfect SIC

In this part, we consider the scenario that symbol $x_1$ is incompletely removed at $R_b$ and $D_2$ . The SINRs of $x_2$ at $R_b$ and $D_2$ are, respectively, given by ¹¹

$${\gamma }_{R,2}^{\mathrm{ISIC}}=\frac{a_2{\rho }_S{\left|h_{S{R}_b}\right|}^2}{k_1{a}_1{\rho }_S{\left|h_{S{R}_b}\right|}^2+\eta +1}=\frac{a_{2}X}{k_1{a}_{1}X+\eta +1}$$

(55)

$${\gamma }_{D2,2}^{\mathrm{ISIC}}=\frac{a_2{\rho }_R{\left|g_{R_b{D}_2}\right|}^2}{k_2{a}_1{\rho }_R{\left|g_{R_b{D}_2}\right|}^2+\eta +1}=\frac{a_2{Y}_2}{k_2{a}_1{Y}_2+\eta +1}$$

(56)

where $0\le k_i\le 1(i\in 1,2)$ represents the level of residual interference due to imperfect SIC.

Formulation analysis of partial RS

From equations (55) and (56), the OP of the system in the case of imperfect SIC is given by

$$\begin{array}{c}\mathrm{O}{\mathrm{P}}^{\mathrm{ISIC}}=1-P({\gamma }_{R,1}>{\epsilon }_1,{\gamma }_{R,2}^{\mathrm{ISIC}}>{\epsilon }_2,{\gamma }_{D1,1}>{\epsilon }_1\hfill \\ ,{\gamma }_{D2,1}>{\epsilon }_1,{\gamma }_{D2,2}^{\mathrm{ISIC}}>{\epsilon }_2)\hfill \\ {=}_{\left(b\right)}1-P(X>{\varphi }_{k1},Y_1>{\xi }_1,Y_2>{\varphi }_{k2})\hfill \\ =1-[1-F_X\left({\varphi }_{k1}\right)][1-F_{Y_1}\left({\xi }_1\right)][1-F_{Y_2}\left({\varphi }_{k2}\right)]\hfill \end{array}$$

(57)

where ${\xi }_{k1}={\epsilon }_2(1+\eta )/a_2-k_1{a}_1{\epsilon }_2$ , ${\xi }_{k2}={\epsilon }_2(1+\eta )/a_2-k_2{a}_1{\epsilon }_2$ , ${\varphi }_{k1}=max\{{\xi }_1,{\xi }_{k1}\}$ and ${\varphi }_{k2}=max\{{\xi }_1,{\xi }_{k2}\}$ . In addition, step (b) requires that $a_2-k_1{a}_1{\epsilon }_2>0$ , $a_2-k_2{a}_1{\epsilon }_2>0$ and $a_1-a_2{\epsilon }_1>0$ , otherwise $O{P}^{\mathrm{ISIC}}=1$ .

Finally, by substituting equations (19) and (20) into equation (57), $O{P}^{\mathrm{ISIC}}$ can be expressed as

$$\begin{array}{c}\mathrm{O}{\mathrm{P}}^{\mathrm{ISIC}}=1-\sum _{n=1}^N\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{n}\end{array}\right){(-1)}^{n-1}\hfill \\ \multicolumn{1}{c}{\times (e^{-\frac{n{\varphi }_{k1}}{{\lambda }_0\rho }}-\frac{n{\varphi }_{k1}{\lambda }_{\mathrm{SP}}}{n{\varphi }_{k1}{\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{e}^{-\frac{n{\varphi }_{k1}{\lambda }_{\mathrm{SP}}+{\lambda }_0{\rho }_I}{{\lambda }_0{\lambda }_{\mathrm{SP}}\rho }})}\\ \multicolumn{1}{c}{\times (e^{-\frac{{\xi }_1}{{\lambda }_1\beta \rho }}-\frac{{\xi }_1{\lambda }_{\mathrm{RP}}}{{\xi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_1{\rho }_I}{e}^{-\frac{{\xi }_1{\lambda }_{\mathrm{RP}}+{\lambda }_1{\rho }_I}{{\lambda }_1{\lambda }_{\mathrm{RP}}\beta \rho }})}\\ \multicolumn{1}{c}{\times (e^{-\frac{{\varphi }_{k2}}{{\lambda }_2\beta \rho }}-\frac{{\varphi }_{k2}{\lambda }_{\mathrm{RP}}}{{\varphi }_{k2}{\lambda }_{\mathrm{RP}}+{\lambda }_2{\rho }_I}{e}^{-\frac{{\varphi }_{k2}{\lambda }_{\mathrm{RP}}+{\lambda }_2{\rho }_I}{{\lambda }_2{\lambda }_{\mathrm{RP}}\beta \rho }})}\end{array}$$

(58)

Formulation analysis of max–min RS

From equations (33), (55) and (56), the OP of the system in the case of imperfect SIC is given by

$$\begin{array}{c}\mathrm{O}{\mathrm{P}}^{M,\mathrm{ISIC}}=1-P({\gamma }_{R,1}>{\epsilon }_1,{\gamma }_{R,2}^{\mathrm{ISIC}}>{\epsilon }_2,{\gamma }_{D1,1}>{\epsilon }_1\hfill \\ \multicolumn{1}{c}{,{\gamma }_{D2,1}>{\epsilon }_1,{\gamma }_{D2,2}^{\mathrm{ISIC}}>{\epsilon }_2)}\\ \multicolumn{1}{c}{=1-[1-F_X^M\left({\varphi }_{k1}\right)][1-F_{Y_1}^M\left({\xi }_1\right)][1-F_{Y_2}^M\left({\varphi }_{k2}\right)]}\end{array}$$

(59)

By substituting equations (31) and (32) into equation (59), the OP of the system can be obtained.

Formulation analysis of two-stage RS

For the two-stage RS scheme in the case of imperfect SIC, the expression (6) is given as

$$\begin{array}{c}{R}_b=\arg \underset{n}{max}\{min\{\log (1+\frac{a_2{\rho }_S{\left|h_{S{R}_n}\right|}^2}{k_1{a}_1{\rho }_S{\left|h_{S{R}_b}\right|}^2+\eta +1}),\log (1+\frac{a_2{\rho }_R{\left|g_{R_n{D}_2}\right|}^2}{k_2{a}_1{\rho }_R{\left|g_{R_b{D}_2}\right|}^2+\eta +1})\},n\in S_r\}\hfill \end{array}$$

(60)

$P(O_2)$ in equation (44) should be reformulated. Defining

$$t_n^{\mathrm{ISIC}}=min\{{\gamma }_{R,2}^{\mathrm{ISIC}},{\gamma }_{D2,2}^{\mathrm{ISIC}}\}$$

(61)

and

$$t_b^{\mathrm{ISIC}}=max\{t_i^{\mathrm{ISIC}},i\in S_r\}$$

(62)

The CDF of $t_n^{\mathrm{ISIC}}$ is given by

$$F_{t_n^{\mathrm{ISIC}}}\left({\epsilon }_2\right)=1-\frac{[1-F_X^T\left({\varphi }_{k1}\right)][1-F_{Y_2}\left({\varphi }_{k2}\right)]}{[1-F_X^T\left({\xi }_1\right)][1-F_{Y_2}\left({\xi }_1\right)]}$$

(63)

Therefore, we have

$$\begin{array}{c}{P}^{\mathrm{ISIC}}\left(O_2\right)=\sum _{l=1}^{N}P(t_b^{\mathrm{ISIC}}<{\epsilon }_2,\left|S_r\right|=l)\hfill \\ \multicolumn{1}{c}{=\sum _{l=1}^N{\left(F_{t_n^{\mathrm{ISIC}}}\left({\epsilon }_2\right)\right)}^{l}P(\left|S_r\right|=l)}\end{array}$$

(64)

From equations (42) and (64), the OP of system is derived as in equation (65), which is

$$\begin{array}{c}\mathrm{O}{\mathrm{P}}^{T,\mathrm{ISIC}}=\underset{n=1}{\overset{N}{\Pi }}[1-(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))]\hfill \\ \multicolumn{1}{c}{+\sum _{l=1}^N{[1-\frac{(1-F_X^T\left({\varphi }_{k1}\right))(1-F_{Y_2}\left({\varphi }_{k2}\right))}{(1-F_X^T\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))}]}^l\left(\begin{array}{c}N\hfill \\ \multicolumn{1}{c}{l}\end{array}\right)}\\ \multicolumn{1}{c}{{[1-(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))]}^{N-l}}\\ \multicolumn{1}{c}{\times {\left[(1-F_X^T\left({\xi }_1\right))(1-F_{Y_1}\left({\xi }_1\right))(1-F_{Y_2}\left({\xi }_1\right))\right]}^l}\end{array}$$

(65)

Performance analysis

This section presents MATLAB simulations to investigate the impacts of system parameters on outage performance. All the simulations are obtained by performing ${10}^6$ channel realizations. The settings of the system parameters are as follows: $d_{\mathrm{SP}}=d_{\mathrm{RP}}=d_0=d_2=d$ , where $d$ is assumed to be unity for simplicity; without loss of generality, $d_1=2d_2$ , ${\sigma }^2=1$ , $\alpha =3$ , $\beta =1$ , $\eta =0.5$ , $a_1=0.8$ , $a_2=0.2$ and ${\epsilon }_1={\epsilon }_2=3\mathrm{dB}$ , extracted from. ¹⁹

Figure 2 shows the OP results versus the transmit SNR for the three RS schemes in the case of perfect SIC, where the number of relays $N=3$ and $I=30\mathrm{dB}$ . In addition, the asymptotic cases without ITC $(I\to \infty )$ and high SNR $(\rho \to \infty )$ are also plotted. It is observed that outage floors occur in the high region of the transmit SNR, which indicates that the transmit power of $S$ and $R_n$ is limited by ITC. Furthermore, the two-stage RS scheme achieves the best outage performance, while the max–min RS scheme outperforms the partial RS scheme regardless of the presence or absence of ITC. This is because the two-stage RS scheme considers the QoS requirements of both users, the max–min RS scheme is based on the channels of the entire links, and the partial RS scheme only considers the channels of the first hop. It is noticed that the simulation results are in excellent agreement with the analytical results.

Figure 2.

OP versus transmit SNR with different RS schemes, where $N=3$ and $I=30\mathrm{dB}$ .

To compare the outage performance in the case of perfect SIC and imperfect SIC, in Figure 3, the OP results are depicted based on the partial RS scheme, assuming $N=3$ and $I=30\mathrm{dB}$ . Moreover, the OP of OMA based on the partial RS scheme is illustrated as a benchmark to be compared with the NOMA. For OMA, four time slots are needed. Thus, the data rate requirement for OMA is set as twice higher as for NOMA for a fair comparison. Besides, the power allocated for each data transmission is equal to $1/2P_S$ . In perfect SIC case, NOMA obtains better outage performance compared to that of OMA. It is noticed that the residual interference factor $(k_i)$ significantly degrades the OP of the system, when $k_i=0.12$ , OMA is even better than NOMA.

Figure 3.

Comparison of OP for the perfect and imperfect SIC, where $N=3$ and $I=30\mathrm{dB}$ .

In Figure 4, the OP versus the number of relays in the case of perfect SIC is illustrated, where the SNR is $\rho =30\mathrm{dB}$ and the ITC is $I=20\mathrm{dB}$ . First, for the partial RS scheme, it is noticed that the improvement of outage performance is no longer obvious when the number of relays exceeds 2, while the performance gain is significant for the other two RS schemes. It can be concluded that there is no need to employ more than two relays for NOMA system with partial RS. Furthermore, the gaps between the curves achieved by the three RS schemes become large when increasing the number of relays.

Figure 4.

OP versus the number of relays with different RS schemes, where $\rho =30\mathrm{dB}$ and $I=20\mathrm{dB}$ .

The impacts of the PA on the outage performance with different RS schemes are investigated in Figure 5 in the case of perfect SIC, where $N=3$ , $\rho =30\mathrm{dB}$ and the ITC is $I=20\mathrm{dB}$ . As shown in the figure, the plots first decrease and then increase with increasing $a_1$ when ${\epsilon }_1<a_1/a_2$ . This can be explained by the fact that the larger $a_1$ is, the easier it is to decode $x_1$ , and the more difficult it is to decode $x_2$ . Hence, the combined effect of decoding $x_1$ and $x_2$ leads to the changes in the curves. In addition, it can be observed that the optimal PA factors are different for the three RS schemes. The optimal PA factors for different RS schemes are shown in Table 1.

Figure 5.

OP versus $a_1$ with different RS schemes, where $N=3$ , $\rho =30\mathrm{dB}$ and $I=20\mathrm{dB}$ .

Table 1.

Optimal PA factors for different RS schemes, where $N=3$ , $\rho =30\mathrm{dB}$ and $I=20\mathrm{dB}$ .

RS schemes	Partial RS	Max–min RS	Two-stage RS
$a_1$	0.823	0.774	0.798

RS: relay selection.

Conclusion

This article investigated the impacts of RS schemes on DF cooperative underlay CR-NOMA networks. Particularly, the partial RS scheme, the max–min RS scheme and the two-stage RS scheme were studied, where the closed-form expressions and asymptotic expressions for the OP were obtained. Moreover, theoretical and simulation results revealed the impacts of ITC, PA factor, residual interference factor and the number of relay nodes on the outage performance. Some important conclusions can be summarized from the results. (1) The two-stage RS scheme outperforms the other two schemes in outage performance. (2) The outage performance of the system can be improved by increasing the number of relay nodes, but for the partial RS scheme, the performance gain can be neglected when the number of relays exceeds 2; therefore, the large number of relay nodes is unnecessary for the partial RS scheme. (3) Selecting the appropriate PA factors for different RS schemes can significantly improve outage performance. (4) The outage performance is severely affected by residual interference. When the residual interference factor is large, the outage performance of NOMA is even worse than that of OMA. The above conclusions can guide the design of practical communication systems.