Machine Learning Approaches for Spectrum Management in Cognitive Radio Networks

Abstract

Cognitive radio (CR) provides a better way for utilization of spectrum resource by introducing an opportunistic usage of the frequency bands that are not heavily occupied by a licensed spectrum user or a primary user (PU). In cognitive radio, the detection and estimation of PU channel availability (unoccupied spectrum) are the key challenges that need to be overcome in order to prevent the interference with licensed spectrum user and improve spectrum resource utilization efficiency. This chapter focuses on developing new ways for detecting and estimating primary user channel availability based on machine-learning (ML) techniques.

1. Introduction

In this chapter, we study the problem of detection of unoccupied primary user spectrum (i.e., spectrum hole). We also introduce the methods for estimating the time when primary user channel state is available, so that the secondary spectrum user can adjust their transmission strategies accordingly.

The chapter is organized in two parts. The first part of the chapter focuses on the problem of detecting the unoccupied spectrum left by the primary user. In this part, we introduce the usage of machine-learning (ML) techniques as a fusion algorithm in cooperative spectrum sensing based on energy detector [1, 2]. In particular, we train a machine-learning classifier (i.e., K-nearest neighbor (KNN), support vector machine (SVM), Naive Bayes (NB), and Decision tree (DT)) over a set containing energy test statistics of PU channel frames along with their corresponding decisions about the presence or absence of PU transmission in the channel. Then, we use the trained classifier to predict the decisions for newly unseen PU channel frames [3]. The second part focuses on estimating the near future of PU channel state. In the literature, there are many proposals that have studied the problem of estimating PU channel state in cognitive radio (CR) [4, 5, 6]. However, most of these studies focused on predicting PU channel state in frequency domain by converting the received digital signals into frequency domain using fast Fourier transform (FFT). This increases the system complexity due to the FFT computations process. In the second part of the chapter, we introduce a new time-domain approach for PU channel state prediction based on time series prediction with some machine-learning prediction model. In particular, a time series is used to capture PU channel state detection sequence (PU channel “idle” or “occupied”) in time domain. Then, prediction models such as the hidden Markov model (HMM) and Markov switching model (MSM) are used to predict the behavior of the time series that used capture PU channel state [7].

2. Machine-learning fusion-based cooperative spectrum sensing

In this part, we, first, define the system model for energy detection-based spectrum sensing; then, we present the method of calculating the thresholds for energy detector with different fusion rules. Second, we formulate a machine-learning classification problem and present four machine-learning classifiers to solve it. Then, we evaluate the performance of these classifiers with simulation experiments.

2.1. Energy detection-based cooperative spectrum sensing

Figure 1 shows a block diagram of the system model used for energy detection cooperative spectrum sensing based on machine-learning fusion rule. In this model, we consider a cooperative CR network with $K$ cooperative nodes. Each node uses $N$ samples for energy detection, while $M\text{frames}$ are used for training the machine-learning (ML) classifier. The received signal of $i_{\mathit{th}}$ frame at the $j_{\mathit{th}}$ cooperative node $y_{\mathit{ij}}\left(n\right)$ ,$1\le n\le N$ ,$1\le i\le M$, $1\le j\le K$ is given by

$$y_{\mathit{ij}}\left(n\right)=\left\{\begin{array}{c}{w}_{\mathit{ij}}\left(n\right)H0\\ \sqrt{{\gamma }_{\mathit{ij}}}{s}_{\mathit{ij}}\left(n\right)+w_{\mathit{ij}}\left(n\right)H1\end{array}\right.$$E1

where ${\mathrm{s}}_{\mathrm{ij}}\left(\mathrm{n}\right)$ is the PU signal which is assumed to follow Gaussian i.i.d random process (i.e., zero mean and ${{\mathrm{\sigma }}_{\mathrm{s}}}^2$ variance), $w_{\mathit{ij}}\left(n\right)$ is the noise which is also assumed to follow Gaussian i.i.d random process (zero mean and ${{\sigma }_u}^2$ variance) because $s_{\mathit{ij}}\left(n\right)$ and $w_{\mathit{ij}}\left(n\right)$ are independent. Due to the fact that all K nodes are sensing the same frame at a given time, the global decision about PU channel availability will be made at the fusion center only. Thus, the energy statistic for the $i_{\mathit{th}}$ frame at the $j_{\mathit{th}}$ cooperative node $Y_{\mathit{ij}}$ can be represented by the energy test statistic of the $i_{\mathit{th}}$ frame at the fusion center which is given by

$$Y_i=\frac{1}{N}\sum _{n=1}^N{\left|\left(y_{\mathit{ij}}\left(n\right)\right)\right|}^2,1\le i\le M$$E2

Figure 1.

Block diagram of machine-learning-based fusion rule spectrum sensing.

$Y_i$ is a random variable that has chi-square distribution probability density function ($2N$ degrees of freedom for complex value and with $N$ degrees of freedom for real value case). If we assume that the channel remains unchanged during the observation interval and there are enough number of samples observed $\left(N\ge 200\right)$ [8], then we can approximate $Y_i$ using Gaussian distribution as follows:

$$Y_i=\left\{\begin{array}{c}\left({{\sigma }_{\mathit{ij}}}^2,2{{\sigma }_{\mathit{ij}}}^4/N\right)H0\\ ({{\sigma }_{\mathit{ij}}}^2\left(1+{\gamma }_{\mathit{ij}}\right),2{{\sigma }_{\mathit{ij}}}^4{\left(1+{\gamma }_{\mathit{ij}}\right)}^2/NH1\end{array}\right.$$E3

where ${{\sigma }_{\mathit{ij}}}^2$, is the standard deviation of noise samples $w_{\mathit{ij}}\left(n\right)$, and ${\gamma }_{\mathit{ij}}$ is the observed signal-to-noise ratio (SNR) of the $i_{\mathit{th}}$ frame sensed at the $\mathrm{j}$ th cooperative node. Assuming that the noise variance and the SNR at the node remain unchanged for all $M$ frames, then ${\gamma }_{\mathit{ij}}={\gamma }_j$ and ${{\sigma }_{\mathit{ij}}}^2={{\sigma }_j}^2$ . For a chosen threshold ${\lambda }_j$ for each frame in the probability of the false alarm, $P_f$ as given in [9] can be written as

$$\begin{array}{c}{P}_f\left({\lambda }_j\right)=\mathit{Pr}\left(Y_i>{\lambda }_j|H0\right)\\ =\frac{1}{\sqrt{2\pi {\sigma }_j}}\underset{{\lambda }_j}{\overset{\infty }{\int }}{e}^{-{\left({\lambda }_j-{\sigma }_j\right)}^2/\sqrt{2}{{\sigma }_j}^2}\\ =\mathrm{Q}\left(\frac{{\mathrm{\lambda }}_{\mathrm{j}}}{{{\mathrm{\sigma }}_{\mathrm{j}}}^2}-1\right)\end{array}$$E4

and the probability of detection $P_d$ is given by

$$\begin{array}{c}{P}_d\left({\lambda }_j\right)=\mathit{Pr}\left(Y_i>{\lambda }_j|H1\right)\\ =Q\left(\left(\frac{{\lambda }_j}{{{\sigma }_j}^2\left(1+{\gamma }_j\right)}-1\right)\sqrt{\frac{N}{2}}\right)\end{array}$$E5

where $\mathrm{Q}\left(.\right)$ is the complementary distribution function of Gaussian distribution with zero mean and unit variance. To obtain the optimal threshold $\mathrm{\lambda }$ for $\mathrm{K}$ cooperative sensing nodes, data fusion rules are used. The calculation of the thresholds for single user and other fusion rules is presented in subsections 2.1.1 and 2.1.2.

2.1.1. The detection threshold for single-user-based sensing

For single user, sensing the number of the cooperative nodes is one (i.e., $\mathrm{K}$ = 1, ${{\mathrm{\sigma }}_{\mathrm{j}}}^2={{\mathrm{\sigma }}_{\mathrm{u}}}^2,{\mathrm{\gamma }}_{\mathrm{j}}={\mathrm{\gamma }}_{\mathrm{u}}. \text{From}$ Eq. (4) and for a given probability of false alarm ${\mathrm{P}}_{\mathrm{f}}$, the single-user threshold can be written as

$${\lambda }_{\mathit{\text{single}}}=\left(\sqrt{\frac{2}{N}}{Q}^{-1}\left(P_f\right)+1\right){{\sigma }_u}^2$$E6

where ${\mathrm{Q}}^{-1}\left(.\right)$ is the inverse of the $\mathrm{Q}\left(.\right)$ function, and the probability of the detection $P_{\mathit{\text{dsingle}}}$ can be written as

$$P_{\mathit{\text{dsingle}}}=Q\left(\left(\frac{{\lambda }_{\mathit{\text{single}}}}{{{\sigma }_u}^2\left(1+{\gamma }_u\right)}-1\right)\sqrt{\frac{N}{2}}\right)$$E7

2.1.2. The detection threshold for data fusion-based sensing

In a data fusion spectrum sensing scheme,$\mathrm{K}$ nodes cooperate in calculating the threshold that is used to make the global sensing decision. There are many fusion rules used to calculate the global sensing decision threshold, which are divided into: hard fusion rules including AND, OR, and majority rule and soft fusion rules including maximum ratio combining (MRC), equal gain combining (EGC), and square law selection (SLS).

2.1.2.1. AND fusion rule

The AND rule decides that the signal is present if all users have detected the signal. For a system with K cooperative nodes with the same false alarm probability $P_f$ cooperating using AND rule, the fusion center threshold can be expressed as

$${\lambda }_{\mathit{\text{AND}}}=\left(\sqrt{\frac{2}{N}}{Q}^{-1}\left({P_f}^{\frac{1}{K}}\right)+1\right){{\sigma }_u}^2$$E8

And the detection probability $P_{\mathit{\text{dAND}}}$ can be written as

$$P_{\mathit{\text{dAND}}}={\left(Q\left(\left(\frac{{\lambda }_{\mathit{\text{AND}}}}{{{\sigma }_u}^2\left(1+{\gamma }_u\right)}-1\right)\sqrt{\frac{N}{2}}\right)\right)}^K$$E9

2.1.2.2. OR fusion rule

The OR rule decides that a signal is present if any of the users detect a signal. The fusion center threshold for K cooperative nodes cooperate using OR fusion rule which can be expressed as

$${\lambda }_{\mathit{\text{OR}}}=\left(\sqrt{\frac{2}{N}}{Q}^{-1}\left((1-{\left(1-P_f\right)}^{\frac{1}{K}}\right)+1\right){{\sigma }_u}^2$$E10

And the detection probability $P_{\mathit{dOR}}$ is

$$P_{\mathit{dOR}}={\left(1-(1-Q\left(\left(\frac{{\lambda }_{\mathit{\text{OR}}}}{{{\sigma }_u}^2\left(1+{\gamma }_u\right)}1\right)\sqrt{\frac{N}{2}}\right)\right)}^K$$E11

2.1.2.3. Maximum ratio combination (optimal MRC) fusion rule

In soft combination fusion $K,$ cooperative nodes with noise variances $\left\{{{\sigma }_{11}}^2,{{\sigma }_{22}}^2,,\dots ,{{\sigma }_{\mathit{MK}}}^2\right\}$ and instantaneous SNRs {${\gamma }_{11},{\gamma }_{22},\dots ,{\gamma }_{\mathit{MK}}$} send their $i_{\mathit{th}}$ frame energy test statistics $Y_{\mathit{ij}}=\frac{1}{N}\sum _{n=1}^N{\left|\left(y_{\mathit{ij}}\left(n\right)\right)\right|}^2,1\le j\le K$ to the fusion center. The fusion center, weighs and adds them together after receiving these energy statistics as follows:

$${\mathit{Ys}}_i=\sum _{j=1}^K{w}_j{Y}_{\mathit{ij}},1\le i\le M$$E12

An assumption is made that SNRs and noise variances at the sensing node will remain unchanged for all the frames during the training process (i.e., ${\gamma }_{\mathit{ij}}={\gamma }_j$, ${{\sigma }_{\mathit{ij}}}^2={{\sigma }_j}^2$). For soft optimal linear combination, we need to find the optimum weight vector $w_j$ that maximizes the detection probability. For additive white Gaussian noise (AWGN) channels, the fusion threshold for MRC fusion rule is written as

$${\lambda }_{\mathit{MRC}}=\left(\sum _{j=1}^K{w}_j{{\sigma }_j}^2\right)Q^{-1}\left[P_f\right]+\sum _{j=1}^K{w}_j{{\sigma }_j}^2)$$E13

And the detection probability $P_{\mathit{\text{dMRC}}}$ is given by

$$P_{\mathit{\text{dMRC}}}=Q\left(\left(\frac{{\lambda }_{\mathit{MRC}}}{{\sum }_{j=1}^K\left(1+{\gamma }_j\right)w_j{{\sigma }_j}^2}-1\right)\sqrt{\frac{N}{2}}\right)$$E14

where the weighting coefficient vector ${\mathrm{w}}_{\mathrm{j}}\left\{{\mathrm{w}}_1,{\mathrm{w}}_2\dots {\mathrm{w}}_{\mathrm{K}}\right\}$ can be obtained by:

$${\mathrm{w}}_{\mathrm{j}}=\text{sign}\left({\mathrm{g}}^{\mathrm{T}}{\mathrm{w}}_0\right){\mathrm{w}}_0$$

where

$$w_0=\frac{{L_{H1}}^{-1/2}{\left[{L_{H1}}^{-1/T}\right]}^{T}g}{‖{L_{H1}}^{-1/2}{\left[{L_{H1}}^{-1/2}\right]}^{T}g‖}$$

where

$$L_{H1}=2\mathit{diag}\left({{\sigma }_1}^4{\left(1+{\gamma }_1\right)}^2,\dots ..{{\sigma }_k}^4{\left(1+{\gamma }_K\right)}^2\right)/N$$

$$g={\left[{{\sigma }_1}^2{\gamma }_1,{{\sigma }_2}^2{\gamma }_2,{{\sigma }_3}^2{\gamma }_3,{{\sigma }_4}^2{\gamma }_4,\dots .,{{\sigma }_K}^2{\gamma }_K\right]}^T$$

2.1.2.4. Equal gain combination (EGC) fusion rule

Equal weight linear combining employs straightforward averaging of the received soft decision statistics. In the equal gain combination, the received energies are equally weighted and then added together. The calculation of the threshold ${\mathrm{\lambda }}_{\mathrm{EGC}}$ and the detection probability ${\mathrm{P}}_{\text{dEGC}}$ follow Eqs. (13) and (14), respectively; the weighting vector is $\left\{w_j=w_1,\dots w_K\right\}$ where $w_1=w_2=w_3\dots =w_K=1/\sqrt{K}$ [10].

2.1.2.5. Square law selection (SLS) fusion rule

Here, the fusion center selects the node with the highest SNR ${\gamma }_{\mathit{SLS}}=\mathit{MAX}\left({\gamma }_1,{\gamma }_2,..{\gamma }_k\right)$ and considers the noise variance ${{\sigma }_{\mathit{SLS}}}^2$ associated with that node. Then the fusion center threshold is calculated as follows:

$${\lambda }_{\mathit{SLS}}=\left(\sqrt{\frac{2}{N}}{Q}^{-1}\left(1-{\left(1-P_f\right)}^{\frac{1}{K}}\right)+1\right){{\sigma }_{\mathit{SLS}}}^2$$E15

And the detection probability $P_{\mathit{\text{dSLS}}}$ is

$$P_{\mathit{\text{dSLS}}}=1-{\left((1-Q\left(\left(\frac{{\lambda }_{\mathit{SLS}}}{{{\sigma }_{\mathit{SLS}}}^2\left(1+{\gamma }_{\mathit{SLS}}\right)}-1\right)\sqrt{\frac{N}{2}}\right)\right)}^K$$E16

2.2. Machine-learning classification problem formulation

The $i_{\mathit{th}}$ frame energy test statistic ($Y_i$ for hard fusion or ${\mathit{Ys}}_i$ for soft fusion rule) given in Eq. (2) or (12) is compared to the sensing threshold to calculate the decision ${\mathrm{d}}_{\mathrm{i}}$ associated with $i_{\mathit{th}}$ frame in the training data set as follows:

$$d_i=\left\{\begin{array}{cc}1& Y_F\ge \lambda \\ -1& Y_F<\lambda \end{array}1\le i\le M\right.$$E17

where $\lambda \in \left({\lambda }_{\mathit{\text{single}}},{\lambda }_{\mathit{\text{and}}},{\lambda }_{\mathit{\text{OR}}},{\lambda }_{\mathit{MRC}},{\lambda }_{\mathit{EGC}},{\lambda }_{\mathit{SLS}}\right),Y_F\in \left\{Y_i,{\mathit{Ys}}_i\right\},$ M is the number of frames in the training set and “$-1$ ” represents the absence of primary user on the channel, and “1” represents the presence of the primary user transmission on the channel. The output of Eq. (17) gives a set of pairs $\left(Y_i,d_i\right),i=1,2\dots M,d_i\in \left(-1,1\right)$ that represent frame energy test statistics and their corresponding decisions. If we want to detect the decision (i.e., the class label) ${\mathrm{d}}_{\mathrm{x}}$ associated with a new frame energy test statistic $Y_x$, we can use one of the following machine-learning classifiers to solve this classification problem.

2.2.1. K-nearest neighbors (KNN) classifier

For K-nearest neighbors classifier,$K$ nearest points to $Y_x$ are used to predict the class label ${\mathrm{d}}_{\mathrm{x}}$ which corresponds to ${\mathrm{Y}}_{\mathrm{x}}$ [11]. for $\mathrm{K}=1$, the Euclidian distance $d_{\mathit{st}}$ between $Y_x$ and the training data points can be computed as

$$d_{\mathit{st}}\left(i\right)=\sqrt{{\left(Y_x-Y_i\right)}^2}=\left|Y_x-Y_i\right|i=1,2\dots M$$E18

and, the new $Y_x$ is classified with the label $d_x$ = $d_{\mathit{\text{in}}}$, where $d_{\mathit{\text{in}}}$ is the point that achieves the minimum Euclidian distance between $d_{\mathit{st}}$ and $Y_x.$

2.2.2. Naïve Bayes classifier

Under the assumption that $d_i=-1$ and $d_i=1$ are independent, the prior probabilities for $d_i=-1$ and $d_i=1$ given training example $\left(Y_i,d_i\right),i=1,2,\dots ,M$ can be calculated, and the class-conditional densities (likelihood probabilities) can also be estimated from the set $\left[Y_1,Y_2,\dots ,Y_k\right].\left[{\mathrm{Y}}_1,{\mathrm{Y}}_2,\dots ,{\mathrm{Y}}_{\mathrm{k}}\right]\text{in}$ which the new ${\mathrm{Y}}_{\mathrm{x}}$ is expected to fall in. And, the probability that the new ${\mathrm{Y}}_{\mathrm{x}}$ to be a member of either $d_i=-1$ or $d_i=1$ class is calculated using Naïve Bayes assumption and Bayes rule [12] as follows:

$$\mathit{\text{class}}\left(Y_x\right)=\underset{d_i}{\mathit{\text{argmax}}}\mathit{Pr}\left(d_i\right){\prod }_{j=1}^k\mathit{Pr}\left(Y_j/d_i\right)$$E19

where the prior probabilities are given to

$$\mathit{Pr}\left(d_i=-1\right)=\frac{\mathit{\text{number of}}{Y}_i\mathit{\text{with}}{\mathit{\text{class label}}}^{"}{1}^{"}}{\mathit{\text{total number of class labels}}}$$

$$\mathit{Pr}\left(d_i=1\right)=\frac{\mathit{\text{number of}}{Y}_i\mathit{\text{with}}a{\mathit{\text{class label}}}^{"}{0}^{"}}{\mathit{\text{total number of class labels}}}$$

Whereas the class-conditional densities “likelihood probabilities” can be estimated using Gaussian density function by:

$$\mathit{Pr}\left(Y_j/d_i\right)=\frac{1}{{\sigma }_j\sqrt{2\pi }}{e}^{\frac{-\left(Y-{\mu }_j\right)}{2{\sigma }_j}},Y_1<Y<Y_k,{\sigma }_j>0,$$

where ${\mu }_j,{\sigma }_j$ are mean and variance of the set $\left[Y_1,Y_2,\dots ,Y_k\right].$ Eq. (19) means that Naïve Bayes classifier will label the new $Y_x\text{with}$ the class label $d_i$ that achieves the highest posterior probability.

2.2.3. Support vector machine (SVM) classifier

For a given training set of pairs $\left(Y_i,d_i\right),i=1,2\dots M$, where $Y_i\in R$, and $d_i\in \left(+1,-1\right)$, the minimum weight $\mathrm{w}$ and a constant $\mathrm{b}$ that maximize the margin between the positive and negative class (i.e., $w{Y}_i+b=\pm 1$ ) with respect to the hyper-plane equation $w{Y}_i+b=0$ can be estimated using support vector machine classifier by performing the following optimization [13].

$$\underset{w,b}{\mathit{min}}\left(\frac{{‖w‖}^2}{2}\right),\mathit{\text{where}}{‖w‖}^2=w^{T}w$$E20

subject to $d_i\left(w{Y}_i+b\right)\ge 1i=1,2,\dots ,M.$

The solution of this quadratic optimization problem can be expressed using Lagrangian function as

$$L\left(w,b,\alpha \right)=\frac{{‖w‖}^2}{2}-\sum _{i=1}^M{\alpha }_i\left(d_i\left(w{Y}_i+b\right)-1\right),{\alpha }_i\ge 0$$E21

where $\mathrm{\alpha }=\left({\mathrm{\alpha }}_1,{\mathrm{\alpha }}_2,\dots ,{\mathrm{\alpha }}_{\mathrm{M}}\right)$ is the Lagrangian multipliers. IF we let $L\left(w,b,\alpha \right)=0$, we can get $w={\sum }_{i=1}^M{\alpha }_i{d}_i{Y}_i$ and ${\sum }_{i=1}^M{\alpha }_i{d}_i=0$, and by substituting them into Eq. (21), the dual optimization problem that describes the hyper-plane can be written as

$$\underset{\alpha }{\mathit{min}}\left(\frac{1}{2}\sum _{i=1}^M\sum _{j=1}^M{d}_i{d}_j\left(Y_i{Y}_j\right){\alpha }_i{\alpha }_j-\sum _{i=1}^M{\alpha }_j\right),{\alpha }_j\ge 0$$E22

From expression (22), we can assess $\mathrm{\alpha }$ and compute $\mathrm{w}$ using $w={\sum }_{i=1}^M{\alpha }_i{d}_i{Y}_i$. Then by choosing ${\alpha }_i>0,$ from the vector of $\mathrm{\alpha }=\left({\mathrm{\alpha }}_1,{\mathrm{\alpha }}_2,\dots ,{\mathrm{\alpha }}_{\mathrm{M}}\right)$ and calculating $b$ from b =$d_j-{\sum }_{i=1}^M{\alpha }_i{d}_i\left(Y_i{Y}_j\right)$, we classify the new instance $Y_x$ using the following classification function

$$\mathit{\text{class}}\left(Y_x\right)=\mathit{sign}\left(\sum _{i=1}^M{\alpha }_j{d}_j\left(Y_i{Y}_x\right)+b\right)$$E23

which means that the classification of new ${\mathrm{Y}}_{\mathrm{x}}$ can be expressed as dot product of $Y_x$ and the support vectors.

2.2.4. Decision tree (DT) classifier

For the training of a set of pairs of sensing decision $\left(Y_i,d_i\right),i=1,2,\dots ,M$ ,$d_i\in \left(-1,1\right)$, the decision tree classifier creates a binary tree based on either impurity or node error splitting rule in order to split the training set into separate subset. Then, it repeats the splitting rule recursively for each subset until the leaf of the subset becomes pure. After that, it minimizes the error in each leaf by taking the majority vote of the training set in that leaf [14]. For classifying a new example $Y_x$, DT classifier selects the leaf where the new $Y_x$ falls in and classifies the new $Y_x$ with the class label that occurs most frequently among that leaf.

2.3. Performance discussion

Figure 2 shows the receiver operating characteristic (ROC) curves for single-user soft and hard fusion rules under Additive White Gaussian Noise (AWGN) channel. In order to generate this figure, we assume a cognitive radio system with 7 cooperative nodes (i.e., $\mathrm{K}$ = 7) operate at SNR ${\mathrm{\gamma }}_{\mathrm{u}}$ = −22 dB. The local node decisions are made after observing1000 samples (i.e., energy detection samples $\mathrm{N}$ = 100). For soft fusion rules, the SNRs ${\mathrm{\gamma }}_{\mathrm{j}}$ for the nodes are equal to {−24.3, −21.8, −20.6, −21.6, −20.4, −22.2, −21.3} and the noise variances ${{\mathrm{\sigma }}_{\mathrm{j}}}^2$ are $\left\{1,1,1,1,1,1,1\right\}.$ We use a false alarm probability ${\mathrm{P}}_{\mathrm{f}}$ varied from 0 to 1 increasing by 0.025. The simulation results show that soft EGC and optimal MRC fusion rules perform better than other soft and hard fusion rules even though that soft EGC fusion rule does not need any channel state information from the nodes.

Figure 2.

ROC curves for the soft and hard fusion rules under the case of AWGN receiver noise $,{{\sigma }_u}^2=1$,${\gamma }_u$= −22 dB, $K$= 7 users and energy detection over $\mathrm{N}$=1000 samples.

Figure 3 shows the ROC curve depicting the performance of SVM classifier in classifying $1000$ new frames after training it over a set containing M = 1000 frames. The thresholds used for training SVM classifier (i.e., single-user threshold, AND, OR, MRC, SLS, and EGC fusion rule threshold) are obtained numerically by considering the cognitive system used to generate Figure 2; however, here, we set the false alarm probability to ${\mathrm{P}}_{\mathrm{f}}=0.1$ .

Figure 3.

ROC curves shows the performance of SVM classifier in predicting the decisions for 1000 new frames after training it over a set containing1000 frames when single user, AND, OR, MRC, SLS, and EGC thresholds are used for training process.

From Figure 3 and Table 1, we can notice that when training SVM classifier with anyone of the following thresholds: single user, OR, MRC, SLS, or EGC, it can detect 100% positive classes. We can also notice that training with EGC threshold can provide 90% precession in classifying the positive classes with 10% harmful interference, whereas training SVM with AND threshold can precisely classify the positive classes by 97.8%. Table 1 shows the classification accuracy of SVM classifier (i.e., the proportion of all true classifications over all testing examples) and the precession of classification (i.e., proportion of true positive classes over all positive classes) as well as the recall of classification (i.e., the effectiveness of the classifier in identifying positive classes).

Threshold	Single user (%)	AND rule (%)	OR rule (%)	MRC rule (%)	SLS rule (%)	EGC rule (%)
SVM
Accuracy	$96.1$	$98.3$	$98.1$	$97.6$	$98.9$	98.0
Precession	77.7	100	53.7	89.4	74.4	90.1
Recall	100	97.8	100	100	100	100

Table 1.

The accuracy, precession and the recall of SVM classifier.

Figure 4 shows ROC curves showing the comparison of four machine-learning classifiers: K-nearest neighbor (KNN), support vector machine (SVM), Naive Bayes and Decision tree when used to classify 1000 frames after training them over a set containing 1000 frames with single-user threshold (Note: the same system used to generate the simulation of Figure 3. is considered for computing the single-user threshold). We can notice from both Figure 4. and Table 2 that KNN and decision tree classifier perform better than Naïve Bayes and SVM classifier in terms of the accuracy of classifying the new frames.

Figure 4.

ROC curves shows a comparison of four machine learning classifiers: KNN, SVM, naive Bayes, and decision tree in classifying 1000 frames after training them over a set with 1000 frames using single user scheme threshold.

Classifier	Accuracy (%)	Precession (%)	Recall (%)
KNN	100	100	100
Decision Tree	100	100	100
Naïve Bayes	98.9	100	91.2
SVM	97.6	83.9	100

Table 2.

The accuracy, precession and recall of KNN, SVM, NB, and DT classifiers used in classifying 1000 new frames after being trained with 1000 frames.

Table 3 shows the accuracy, precession, and the recall for decision tree classifier when used to classify 3000 frames after training it over a set containing 1000 frames for the same cognitive system used to generate Figure 3. The single-user threshold is used for training the classifier. The simulation was run with different number of samples for energy detection process. It is clear from the table that decision tree can classify all of the 3000 frames correctly or achieve 100% detection rate using only 200 samples for the energy detection process. And, due to the fact that the sensing time is proportional to the number of samples taken by energy detector, a less number of samples used for energy detection leads to less sensing time. Thus, when we use machine-learning-based fusion, such as decision tree or KNN, we can reduce the sensing time from 200 to 40 μs for 5 MHz bandwidth channel as an example, while we still achieve 100% detection rate of the spectrum hole.

Number of samples	Accuracy (%)	Precession (%)	Recall (%)
200	100	100	100
400	100	100	100
600	100	100	100
800	100	100	100
1000	100	100	100

Table 3.

The accuracy, precession, and recall for decision tree classifier used in classifying 3000 frames for different number of samples.

3. Prediction of PU channel state based on hidden Markov and Markov switching model

In this part, the system model for forecasting the near future of PU channel state is divided into three models: (1) the model detecting the PU channel state (i.e., PU signal present or PU signal) which follows the conventional single-user energy detection (i.e., fusion techniques mentioned in Section 2.1 can also be considered here); (2) the model that generates a time series to capture PU channel state based on the detection sequence; and (3) the model for predicting the generated time series used to capture PU channel state based on hidden Markov model (HMM) and Markov switching model (MSM). The block diagram in Figure 5. illustrates these three models.

Figure 5.

Block diagram of PU channel state prediction model.

The PU channel state detection model can be written using Eq. (4); by giving probability of false alarm ${\mathrm{P}}_{\mathrm{f}}$, the detection threshold for single-user energy detector can be written as:

$$\mathrm{\lambda }=\left(\sqrt{\frac{2}{\mathrm{N}}}{\mathrm{Q}}^{-1}\left({\mathrm{P}}_{\mathrm{f}}\right)+1\right){{\mathrm{\sigma }}_{\mathrm{u}}}^2$$E24

where ${\mathrm{Q}}^{-1}\left(.\right)$ is the inverse of the $\mathrm{Q}\left(.\right)$ function.

And the decision of the sensing (i.e., PU detection sequence) over the time can be written as follows:

$${\mathrm{D}}_{\mathrm{t}}=\left\{\begin{array}{c}{}^{"}{0}^{"}\mathrm{PU}\text{absent}{\mathrm{Y}}_{\mathrm{t}}<\mathrm{\lambda }\\ {}^{"}{1}^{"}\mathrm{PU}\text{present}{\mathrm{Y}}_{\mathrm{t}}\ge \mathrm{\lambda }\end{array}\right.1\le \mathrm{t}\le \mathrm{T}$$E25

3.1. Time series generation model

Given PU channel state detection sequence over the time (i.e., PU absent, PU present), if we denote the period that the PU is inactive as “idle state,” and the period that PU is active as “occupied state,” our goal now is to predict when the detection sequence ${\mathrm{D}}_{\mathrm{t}}$ will change from one state to another (i.e., “idle” to “occupied “or vice versa) before that happens so that the secondary user can avoid interfering with primary user transmission. For this reason, we generate a time series ${\mathrm{z}}_{\mathrm{t}}\text{to}$ map each state of the detection sequence ${\mathrm{D}}_{\mathrm{t}}$ (i.e., “PU present” and “PU absent”) into another observation space using two different random variable distributions for each state (i.e., ${\mathrm{z}}_{\mathrm{t}}\in \left\{{\mathrm{v}}_1,{\mathrm{v}}_2\dots {\mathrm{v}}_{\mathrm{L}}\right\}$ represents PU absent or idle state and ${\mathrm{z}}_{\mathrm{t}}\in \left\{{\mathrm{v}}_{\mathrm{L}+1}\dots .{\mathrm{v}}_{\mathrm{M}}\right\}$ represents PU occupied or present), the time series ${\mathrm{z}}_{\mathrm{t}}$ can be written as

$${\mathrm{z}}_{\mathrm{t}}\in \left\{\begin{array}{cc}\left\{{\mathrm{v}}_1,{\mathrm{v}}_2,\dots ,{\mathrm{v}}_{\mathrm{L}}\right\}& {\mathrm{Y}}_{\mathrm{t}}<\mathrm{\lambda }\\ \left\{{\mathrm{v}}_{\mathrm{L}+1},..{\mathrm{v}}_{\mathrm{M}}\right\}& {\mathrm{Y}}_{\mathrm{t}}\ge \mathrm{\lambda }\end{array}1\le \mathrm{t}\le \mathrm{T}\right.$$E26

Now, supposing that we have given observations value $\mathrm{O}=\left\{{\mathrm{O}}_1,{\mathrm{O}}_2,{\mathrm{O}}_{\mathrm{t}},\dots {\mathrm{O}}_{\mathrm{T}}\right\}$ ,${\mathrm{O}}_{\mathrm{t}}\in \left\{{\mathrm{v}}_1,{\mathrm{v}}_2\dots {\mathrm{v}}_{\mathrm{M}}\right\}$ and a PU channel state at time step $\mathrm{t},{\mathrm{X}}_{\mathrm{t}}\in {\mathrm{s}}_{\mathrm{i}},\mathrm{i}=1,2\dots .\mathrm{K}$, ${\mathrm{s}}_{\mathrm{i}}\in \left\{0,1\right\}$ (i.e., 0 for PU idle and 1 PU occupied state), and we want to estimate the channel state at one time step ahead of the current state ${\mathrm{X}}_{\mathrm{t}+1}$. We can solve this problem using hidden Markov model Viterbi algorithm [15].

3.2. Primary users channel state estimation based on hidden Markov model

The generic HMM model can be illustrated by Figure 6.—in this figure, $\mathrm{X}=\left\{{\mathrm{X}}_1,{\mathrm{X}}_{\mathrm{t}},\dots {\mathrm{X}}_{\mathrm{T}}\right\}$ represents the hidden state sequence, where ${\mathrm{X}}_{\mathrm{t}}\in \left\{{\mathrm{s}}_1,{\mathrm{s}}_2,\dots ,{\mathrm{s}}_{\mathrm{K}}\right\}$, $\mathrm{K}$ represents the number of hidden states or Markov chain and $\mathrm{O}=\left\{{\mathrm{O}}_1,{\mathrm{O}}_{\mathrm{t}},\dots ,{\mathrm{O}}_{\mathrm{T}}\right\}$ represents the observation sequence where ${\mathrm{O}}_{\mathrm{t}}\in \left\{{\mathrm{v}}_1,{\mathrm{v}}_2,\dots ,{\mathrm{v}}_{\mathrm{M}}\right\}$ and $\mathrm{M}$ is the number of the observations in the observation space. $\mathrm{A}$ and $\mathrm{B}$ represent the transition probabilities matrix and the emission probabilities matrix, respectively, while $\mathrm{\pi }$ denotes the initial state probability vector. HMM can be defined by $\mathrm{\theta }=\left(\mathrm{\pi },\mathrm{A},\mathrm{B}\right)$ (i.e., the initial state probabilities, the transition probabilities, and emission probabilities) [15].

Figure 6.

Hidden Markov model.

Initial state probabilities for HMM can be written as

$$\mathrm{\pi }=\left({\mathrm{\pi }}_1,{\mathrm{\pi }}_2,{\mathrm{\pi }}_{\mathrm{i}}\dots .{\mathrm{\pi }}_{\mathrm{K}}\right)$$

$${\mathrm{\pi }}_{\mathrm{i}}=\mathrm{P}\left({\mathrm{X}}_1={\mathrm{s}}_{\mathrm{i}}\right),\mathrm{i}=1,2,\dots ,\mathrm{K}$$E27

For a HMM model with two hidden states $\mathrm{i}=2$ ,

$$\mathrm{\pi }=\left({\mathrm{\pi }}_1{\mathrm{\pi }}_2\right)$$

And the transition probabilities can be written as,

$$\mathrm{A}={\left({\mathrm{a}}_{\mathrm{ij}}\right)}_{\mathrm{K}\times \mathrm{K}}$$

$${\mathrm{a}}_{\mathrm{ij}}=\mathrm{P}\left({\mathrm{X}}_{\mathrm{t}+1}={\mathrm{s}}_{\mathrm{j}}|{\mathrm{X}}_{\mathrm{t}}={\mathrm{s}}_{\mathrm{i}}\right),\mathrm{i},\mathrm{j}=1,\dots ..,\mathrm{K}$$E28

where a_ij is the probability that next state equal ${\mathrm{s}}_{\mathrm{j}}$ when current state is equal to ${\mathrm{s}}_{\mathrm{i}}$. For HMM model with two states, the matrix $\mathrm{A}$ can be written as

$$\mathrm{A}=\left(\begin{array}{cc}{\mathrm{a}}_{00}& {\mathrm{a}}_{01}\\ {\mathrm{a}}_{10}& {\mathrm{a}}_{11}\end{array}\right)$$

The emission probabilities matrix for HMM model is written as

$$\mathrm{B}={\left({\mathrm{b}}_{\mathrm{jm}}\right)}_{\mathrm{K}\times \mathrm{M}}$$

$${\mathrm{b}}_{\mathrm{j}}\left(\mathrm{m}\right)=\mathrm{P}\left({\mathrm{O}}_{\mathrm{t}}={\mathrm{v}}_{\mathrm{m}}\right|{\mathrm{X}}_{\mathrm{t}}={\mathrm{s}}_{\mathrm{j}})\triangleq {\mathrm{b}}_{\mathrm{j}}\left({\mathrm{O}}_{\mathrm{t}}\right),\mathrm{j}=1\dots \mathrm{K},\mathrm{m}=1,\dots ..\mathrm{M}$$E29

$\mathrm{B}$ and b_j represent the probability that current observation is ${\mathrm{v}}_{\mathrm{m}}$ when current state is ${\mathrm{s}}_{\mathrm{j}}$. For example, in an HMM model with $\mathrm{M}=6$ and $\mathrm{K}=2$, $\mathrm{B}$ is written as

$$\mathrm{B}=\left(\begin{array}{cc}{\mathrm{b}}_{11}& {\mathrm{b}}_{12}{\mathrm{b}}_{13}\\ {\mathrm{b}}_{21}& {\mathrm{b}}_{22}{\mathrm{b}}_{23}\end{array}\begin{array}{cc}{\mathrm{b}}_{14}& {\mathrm{b}}_{15}{\mathrm{b}}_{16}\\ {\mathrm{b}}_{24}& {\mathrm{b}}_{25}{\mathrm{b}}_{26}\end{array}\right)$$

Now, for the problem we describe in subsection (3.1), if we assume that HMM parameters $\mathrm{\theta }=\left(\mathrm{\pi },\mathrm{A},\mathrm{B}\right)$ and the observations value $\mathrm{O}=\left\{{\mathrm{O}}_1,{\mathrm{O}}_2{\mathrm{O}}_{\mathrm{t}},\dots {\mathrm{O}}_{\mathrm{T}}\right\}$ are given. If we assume that the maximum probability of state sequence $\mathrm{t}$ that end at state $\mathrm{i}$ to be equal to ${\mathrm{\delta }}_{\mathrm{t}}\left(\mathrm{i}\right)$ where

$${\mathrm{\delta }}_{\mathrm{t}}\left(\mathrm{i}\right)=\underset{{\mathrm{X}}_1,\dots ,{\mathrm{X}}_{\mathrm{t}-1}}{\text{max}}\left\{\mathrm{P}\left({\mathrm{X}}_1,\dots ,{\mathrm{X}}_{\mathrm{t}}={\mathrm{s}}_{\mathrm{i}};{\mathrm{O}}_1,\dots ,{\mathrm{O}}_{\mathrm{t}}|\mathrm{\theta }\right)\right\}$$E30

And we let ${\mathrm{\psi }}_{\mathrm{t}}\left(\mathrm{i}\right)$ to be a vector that stores the arguments that maximize Eq. (30), we can write Viterbi algorithm to solve the problem mentioned in subsection (3.1) as follows:

1) step 1 initializes ${\mathrm{\delta }}_{\mathrm{t}}\left(\mathrm{i}\right)$ and ${\mathrm{\psi }}_{\mathrm{t}}\left(\mathrm{i}\right)$.

$${\mathrm{\delta }}_{\mathrm{t}}\left(\mathrm{i}\right)={\mathrm{\pi }}_{\mathrm{i}}{\mathrm{b}}_{\mathrm{i}}\left({\mathrm{O}}_1\right)$$

$${\mathrm{\psi }}_{\mathrm{t}}\left(\mathrm{i}\right)=0,\mathrm{i}=1,\dots ,\mathrm{K}$$E31

2) step 2 iterates to update ${\mathrm{\delta }}_{\mathrm{t}}\left(\mathrm{i}\right)$ and ${\mathrm{\psi }}_{\mathrm{t}}\left(\mathrm{i}\right)$.

$${\mathrm{\delta }}_{\mathrm{t}}\left(\mathrm{j}\right)=\underset{1\le \mathrm{i}\le \mathrm{K}}{\text{max}}\left[{\mathrm{\delta }}_{\mathrm{t}-1}\left(\mathrm{i}\right){\mathrm{a}}_{\mathrm{ij}}\right]{\mathrm{b}}_{\mathrm{j}}\left({\mathrm{O}}_{\mathrm{t}}\right),\mathrm{t}=2,\dots ,\mathrm{T},\mathrm{j}=1,\dots ,\mathrm{K}$$E32

$${\mathrm{\psi }}_{\mathrm{t}}\left(\mathrm{j}\right)=\underset{1\le \mathrm{i}\le \mathrm{K}}{\text{argmax}}\left[{\mathrm{\delta }}_{\mathrm{t}-1}\left(\mathrm{i}\right){\mathrm{a}}_{\mathrm{ij}}\right],\mathrm{t}=2,\dots ,\mathrm{T},\mathrm{j}=1,\dots ,\mathrm{K}$$E33

3) step 3 terminates the update and calculates the likelihood probability ${\mathrm{P}}^{\ast }$ and the estimated state ${{\mathrm{q}}_{\mathrm{T}}}^{\ast }$ at time $\mathrm{T}$ as

$${\mathrm{P}}^{\ast }=\underset{1\le \mathrm{i}\le \mathrm{K}}{\text{max}}\left[{\mathrm{\delta }}_{\mathrm{T}}\left(\mathrm{i}\right)\right]$$E34

$${{\mathrm{q}}_{\mathrm{T}}}^{\ast }=\underset{1\le \mathrm{i}\le \mathrm{K}}{\text{argmax}}\left[{\mathrm{\delta }}_{\mathrm{T}}\left(\mathrm{i}\right)\right]$$E35

In the above case, HMM parameters $\mathrm{\theta }=\left(\mathrm{\pi },\mathrm{A},\mathrm{B}\right)$ are unknown and need to be estimated. We estimate these parameters statistically using Baum-Welch algorithm [16].

3.2.1. Hidden Markov model parameters estimation using Baum-Welch algorithm

If we assume that we have given some training observations with length L $\left\{{\mathrm{O}}_1,{\mathrm{O}}_2{\mathrm{O}}_{\mathrm{t}},\dots {\mathrm{O}}_{\mathrm{L}}\right\}$ and want to approximate HMM parameters $\mathrm{\theta }=\left(\mathrm{\pi },\mathrm{A},\mathrm{B}\right)$ from them, we can use maximum likelihood estimation. In order to do that, we define ${\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)$ to be the probability of being in state ${\mathrm{s}}_{\mathrm{i}}$ at time $\mathrm{t},$ given t ${\mathrm{O}}_{\mathrm{t}}$ ,$\mathrm{t}=1,2\dots \mathrm{L}$ . ${\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)$ is written as

$${\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)=\mathrm{P}\left({\mathrm{X}}_{\mathrm{t}}={\mathrm{s}}_{\mathrm{i}}|{\mathrm{O}}_1,\dots ,{\mathrm{O}}_{\mathrm{L}},\mathrm{\theta }\right)$$E36

We also define ${\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)$ to be the probability of being in state ${\mathrm{s}}_{\mathrm{i}}$ at time $\mathrm{t}$ and transiting to state ${\mathrm{s}}_{\mathrm{j}}$ at time $\mathrm{t}+1,$ given ${\mathrm{O}}_{\mathrm{t}}$ ,$\mathrm{t}=1,2\dots \mathrm{L}$ . ${\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)$ is written as

$${\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)=\mathrm{P}\left(\left({\mathrm{X}}_{\mathrm{t}}={\mathrm{s}}_{\mathrm{i}};{\mathrm{X}}_{\mathrm{t}+1}={\mathrm{s}}_{\mathrm{j}}\right|\mathrm{O}{\mathrm{O}}_1,\dots ,{\mathrm{O}}_{\mathrm{L}},\mathrm{\theta }\right)$$E37

Given ${\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)$ and ${\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)$, the anticipated number of transitions from state ${\mathrm{s}}_{\mathrm{i}}$ during the path is written as

$$\mathrm{E}\left({\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)\right)=\sum _{\mathrm{t}=1}^{\mathrm{L}-1}{\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)$$E38

and the anticipated number of transitions from state ${\mathrm{s}}_{\mathrm{i}}$ to state ${\mathrm{s}}_{\mathrm{j}}$ during the path is written as

$$\mathrm{E}\left({\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)\right)=\sum _{\mathrm{j}=1}^{\mathrm{L}-1}{\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)$$E39

Given $\mathrm{E}\left({\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)\right)$ and $\mathrm{E}\left({\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)\right),$ we can extract the model parameters $\mathrm{\theta }=\left(\mathrm{\pi },\mathrm{A},\mathrm{B}\right)$ from the training sequence as given in [16] using the step listed below

1- for $\mathrm{i}=1,2,3\dots \mathrm{K}$, let ${\widehat{\mathrm{\pi }}}_{\mathrm{i}}=\text{expected frequency in state}{\mathrm{s}}_{\mathrm{i}}\mathrm{at}\text{time}\left(\mathrm{t}=1\right)$

$${\mathrm{\pi }}_{\mathrm{i}}={\mathrm{\gamma }}_1\left(\mathrm{i}\right)$$E40

2- for $\mathrm{i}=1,2,3\dots \mathrm{K}$ and $\mathrm{j}=1,2,3\dots \mathrm{K}$, compute

$${\widehat{\mathrm{a}}}_{\mathrm{ij}}=\frac{\text{Expected number of transitions fromstate}{\mathrm{s}}_{\mathrm{i}}\text{to state}{\mathrm{s}}_{\mathrm{j}}}{\text{Expected number of transitions from state}{\mathrm{s}}_{\mathrm{i}}}$$

$$=\frac{\mathrm{E}\left({\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)\right)}{\mathrm{E}\left({\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)\right)}=\frac{\sum _{\mathrm{j}=1}^{\mathrm{L}-1}{\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)}{\sum _{\mathrm{t}=1}^{\mathrm{L}-1}{\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)}$$E41

3- for $\mathrm{i}=1,2,3\dots \mathrm{K}$ and $\mathrm{j}=1,2,3\dots \mathrm{K}$, compute

$${\widehat{\mathrm{b}}}_{\mathrm{i}}\left(\mathrm{m}\right)=\frac{\text{Expected number of times in state}{\mathrm{s}}_{\mathrm{j}}\text{and observing}{\mathrm{v}}_{\mathrm{m}}}{\text{Expected number of times in state}{\mathrm{s}}_{\mathrm{j}}}$$

$$=\frac{{\sum }_{\begin{array}{c}\mathrm{t}=1\\ {\mathrm{O}}_{\mathrm{t}}={\mathrm{v}}_{\mathrm{m}}\end{array}}^{\mathrm{L}}{\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)}{{\sum }_{\mathrm{t}=1}^{\mathrm{L}}{\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)}$$E42

The estimation algorithm can be summarized in the following steps:

1. Get your observations ${\mathrm{O}}_1{\mathrm{O}}_2,\dots {\mathrm{O}}_{\mathrm{L}}$,

2. Set a guess of your first $\mathrm{\theta }$ estimate $\mathrm{\theta }\left(1\right),\mathrm{k}=1$

   $$\text{Update} \mathrm{k}=\mathrm{k}+1$$

3. Compute $\mathrm{\theta }\left(\mathrm{k}\right)$ based on ${\mathrm{O}}_1{\mathrm{O}}_2,\dots {\mathrm{O}}_{\mathrm{L}}$ and

   $${\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right),{\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)\forall 1\le \mathrm{t}\le \mathrm{L},\forall 1\le \mathrm{i}\le \mathrm{K},\forall 1\le \mathrm{j}\le \mathrm{K}$$

4. Compute $\mathrm{E}\left({\mathrm{\gamma }}_{\mathrm{t}}\left(\mathrm{i}\right)\right)$ and $\mathrm{E}\left({\mathrm{\zeta }}_{\mathrm{t}}\left(\mathrm{i},\mathrm{j}\right)\right)$ from Eqs. (38) and (39)

5. Compute according to 5 the new estimate of ${\mathrm{a}}_{\mathrm{ij}}$, ${\mathrm{b}}_{\mathrm{i}}\left(\mathrm{k}\right)$, ${\mathrm{\pi }}_{\mathrm{i}}$, and call them $\mathrm{\theta }$ (k + 1)

6. Go to 3 if not converged.

The prediction for a one-step ahead PU channel state can be done based on the trained parameters {$\mathrm{\pi },\mathrm{A},\mathrm{B}$ } with the help of Eqs. (31), (34), and (35) by setting $\mathrm{T}=1.$

3.3. Primary users channel state estimation based on Markov switching model

An alternative way to estimate PU channel state is to use Markov switching model (MSM). For the time series in Eq. (26), we assume that ${\mathrm{z}}_{\mathrm{t}}$ obeys two different Gaussian distributions $\mathrm{N}\sim \left({\mathrm{\mu }}_{\mathrm{z}0},{{\mathrm{\sigma }}_{\mathrm{z}0}}^2\right)$ or $\mathrm{N}\sim \left({\mathrm{\mu }}_{\mathrm{z}1},{{\mathrm{\sigma }}_{\mathrm{z}1}}^2\right)$ based on the sensed PU channel state “PU channel idle” or “PU occupied.” We can rewrite Eq. (26) as follows:

$${\mathrm{z}}_{\mathrm{t}}\sim \left\{\begin{array}{c}\left({\mathrm{\mu }}_{\mathrm{z}0},{{\mathrm{\sigma }}_{\mathrm{z}0}}^2\right){\mathrm{Y}}_{\mathrm{t}}<\mathrm{\lambda }\\ \left({\mathrm{\mu }}_{\mathrm{z}1},{{\mathrm{\sigma }}_{\mathrm{z}1}}^2\right){\mathrm{Y}}_{\mathrm{t}}\ge \mathrm{\lambda }\end{array}1\le \mathrm{t}\le \mathrm{T}\right.$$E43

It is obvious that Eq. (43) represents a two-state Gaussian regime switching time series which can be modeled using MSM [17]. In order to estimate the switching time of one state ahead of the current state for this time series, we need to derive MSM regression model for the time series and estimate its parameters.

3.3.1. Derivation of Markov switching model for Gaussian regime switching time series

A simple Markov switching regression model to describe the two-state Gaussian regime switching time series is given in Eq. (43). This model can be written by following Ref [17] as

$${\mathrm{z}}_{\mathrm{t}}={\mathrm{\mu }}_{{\mathrm{s}}_{\mathrm{t}}}+{\mathrm{ϵ}}_{\mathrm{t}}{\mathrm{ϵ}}_{\mathrm{t}}\sim \left(0,{{\mathrm{\sigma }}_{{\mathrm{s}}_{\mathrm{t}}}}^2\right)$$E44

where ${\mathrm{\mu }}_{{\mathrm{s}}_{\mathrm{t}}}$ is an array of predetermined variables measured at time t, which may include the lagged values of ${\mathrm{z}}_{\mathrm{t}}$, ${\mathrm{ϵ}}_{\mathrm{t}}$ is the white noise process,${\mathrm{s}}_{\mathrm{t}}=\left\{0,1\right\}$ is a hidden Markov chain which has a mean and standard deviation over the time equal to ${\mathrm{\mu }}_{\mathrm{st}}={\mathrm{\mu }}_0\left(1-{\mathrm{s}}_{\mathrm{t}}\right)+{\mathrm{\mu }}_1{\mathrm{s}}_{\mathrm{t}}$ and ${\mathrm{\sigma }}_{\mathrm{st}}={\mathrm{\sigma }}_0\left(1-{\mathrm{s}}_{\mathrm{t}}\right)+{\mathrm{\sigma }}_1{\mathrm{s}}_{\mathrm{t}}$, respectively (the state variable ${\mathrm{s}}_{\mathrm{t}}$ follows first order Markov chain (i.e., two-state Markov chain as in [18])). Given the past history of ${\mathrm{s}}_{\mathrm{t}}$, the probability of ${\mathrm{s}}_{\mathrm{t}}$ taking a certain value depends only on ${\mathrm{s}}_{\mathrm{t}-1}$, which takes the following Markov property:

$$\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{j}\left|{\mathrm{s}}_{\mathrm{t}-1}=\mathrm{i}\right)={\mathrm{P}}_{\mathrm{ij}}$$E45

where ${\mathrm{P}}_{\mathrm{ij}}\left(\mathrm{i};\mathrm{j}=0;1\right)$ denotes the transition probabilities of ${\mathrm{s}}_{\mathrm{t}}=\mathrm{j}$, given that ${\mathrm{s}}_{\mathrm{t}-1}=\mathrm{i}$. Clearly, the transition probabilities satisfy ${\mathrm{P}}_{\mathrm{i}0}+{\mathrm{P}}_{\mathrm{i}1}=1$ . We can gather the transition probabilities ${\mathrm{P}}_{\mathrm{ij}}$ into a transition matrix as follows:

$$\begin{array}{c}\text{P}=\left(\begin{array}{cc}\text{P}({\text{s}}_{\text{t}}=0\text{|}{\text{s}}_{\text{t}-1}=0\text{)}& \text{P}({\text{s}}_{\text{t}}=0\text{|}{\text{s}}_{\text{t}-1}=1\text{)}\\ \text{P}({\text{s}}_{\text{t}}=1\text{|}{\text{s}}_{\text{t}-1}=0\text{)}& \text{P}({\text{s}}_{\text{t}}=1\text{|}{\text{s}}_{\text{t}-1}=1\text{)}\end{array}\right)\\ =\left(\begin{array}{cc}{\text{P}}_{00}& {\text{P}}_{01}\\ {\text{P}}_{10}& {\text{P}}_{11}\end{array}\right)\end{array}$$E46

The transition matrix $\mathrm{P}$ is used to govern the behavior of the state variable ${\mathrm{s}}_{\mathrm{t}}$, and it holds only two parameters (${\mathrm{P}}_{00}$ and ${\mathrm{P}}_{11}$ ). Assuming that we do not observe ${\mathrm{s}}_{\mathrm{t}}$ directly, we only deduce its operation from the observed behavior of ${\mathrm{z}}_{\mathrm{t}}$. The parameters that need to be estimated to fully describe the probability law governing ${\mathrm{z}}_{\mathrm{t}}$ are the variance of the Gaussian innovation ${\mathrm{\sigma }}_0,{\mathrm{\sigma }}_1$, the expectation of the dependent variable ${\mathrm{\mu }}_0$, ${\mathrm{\mu }}_{1\prime }$ and the two-state transition probabilities ${\mathrm{P}}_{00}$ and ${\mathrm{P}}_{11}$.

3.3.2. Markov switching model parameters estimation via maximum likelihood estimation

There are many ways to estimate the parameters for the Markov switching model. Among these ways are Quasi-maximum likelihood estimation (QMLE) and Gibbs sampling. In this section, we focus on maximum likelihood estimation (MLE).

If we denote ${\mathrm{\psi }}_{\mathrm{t}-1}$ = {${\mathrm{z}}_{\mathrm{t}-1},{\mathrm{z}}_{\mathrm{t}},{\mathrm{z}}_{\mathrm{t}+1}\dots {\mathrm{z}}_1$} to be a vector of the training data until time $\mathrm{t}-1$ and denote $\mathrm{\theta }$ ={${\mathrm{\sigma }}_0,{\mathrm{\sigma }}_1$, ${\mathrm{\mu }}_0$, ${\mathrm{\mu }}_{1\prime }{\mathrm{P}}_{00},{\mathrm{P}}_{11}$} to be the vector of MSM parameters, then ${\mathrm{\psi }}_{\mathrm{L}}$ = {${\mathrm{z}}_{\mathrm{t}-1}, {\mathrm{z}}_{\mathrm{t}},$ …, ${\mathrm{z}}_{\mathrm{L}}$} to a vector of the available information with the length $\mathrm{L}$ sample (see Figure 7.). In order to evaluate the likelihood of the state variable ${\mathrm{s}}_{\mathrm{t}}$ based on the current trend of ${\mathrm{z}}_{\mathrm{t}},$ we need to assess its conditional expectations ${\mathrm{s}}_{\mathrm{t}}=\mathrm{i}$, $\left(\mathrm{i}=0,1\right)$ based on $\mathrm{\psi }$ and $\mathrm{\theta }$ . These conditional expectations include prediction probabilities $\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)$, which are based on the information prior to time $\mathrm{t}$, the filtering probabilities $\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\mid {\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }$) which are based on the past and current information, and finally the smoothing probabilities $\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)$ which are based on the full-sample information L. After getting these probabilities, we can obtain the log-likelihood function as a byproduct, and then we can compute the maximum likelihood estimates.

Figure 7.

Markov switching model.

Normally, the density of ${\mathrm{z}}_{\mathrm{t}}$ conditional on ${\mathrm{\psi }}_{\mathrm{t}-1}$ and ${\mathrm{s}}_{\mathrm{t}}=\mathrm{i}$, $\left(\mathrm{i}=0,1\right)$ is written as

$$\mathcal{F}({\mathrm{z}}_{\mathrm{t}}\left|{\mathrm{s}}_{\mathrm{t}}=\mathrm{i},{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)=\frac{1}{\sqrt{2\mathrm{\pi }{\mathrm{\sigma }}_{{\mathrm{s}}_{\mathrm{t}}}}}{\mathrm{e}}^{-\frac{{\left({\mathrm{z}}_{\mathrm{t}}-{\mathrm{\mu }}_{{\mathrm{s}}_{\mathrm{t}}}\right)}^2}{2{{\mathrm{\sigma }}_{\mathrm{st}}}^2}}$$E47

where $\mathcal{F}$ represent the probability density function. Given the prediction probabilities $\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right),$ the density of ${\mathrm{z}}_{\mathrm{t}}$ conditional on ${\mathrm{\psi }}_{\mathrm{t}-1}$ can be obtained from

$$\mathcal{F}({\mathrm{z}}_{\mathrm{t}}\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)=$$

$$=\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=0\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)\mathcal{F}\left({\mathrm{z}}_{\mathrm{t}}|{\mathrm{s}}_{\mathrm{t}}=0{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)$$

$$+\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=1\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)\mathcal{F}({\mathrm{z}}_{\mathrm{t}}\left|{\mathrm{s}}_{\mathrm{t}}=1{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)$$E48

For $\mathrm{i}=0;1$, the filtering probabilities of ${\mathrm{s}}_{\mathrm{t}}$ are given by:

$$\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)$$

$$=\frac{\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)\mathcal{F}({\mathrm{z}}_{\mathrm{t}}\left|{\mathrm{s}}_{\mathrm{t}}=\mathrm{i}{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)}{\mathcal{F}({\mathrm{z}}_{\mathrm{t}}\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)}$$E49

The prediction probabilities are:

$$\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)$$

$$=\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=0,{\mathrm{s}}_{\mathrm{t}+1}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)+\mathrm{P}\left({\mathrm{s}}_{\mathrm{t}}=1,{\mathrm{s}}_{\mathrm{t}+1}=\mathrm{i}|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)$$

$$={\mathrm{P}}_{0\mathrm{i}}\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=0\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)+{\mathrm{P}}_{1\mathrm{i}}\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=1\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)$$E50

where ${\mathrm{P}}_{0\mathrm{i}}=\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=\mathrm{i}\left|{\mathrm{s}}_{\mathrm{t}}=0\right)$ and ${\mathrm{P}}_{1\mathrm{i}}=\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=\mathrm{i}\left|{\mathrm{s}}_{\mathrm{t}}=1\right)$ are the transition probabilities. By setting the initial values as given in [19] assuming the Markov chain is presumed to be ergodic:

$$\mathrm{P}({\mathrm{s}}_0=\mathrm{i}\left|{\mathrm{\psi }}_0\right)=\frac{1-{\mathrm{P}}_{\mathrm{jj}}}{2-{\mathrm{P}}_{\mathrm{ii}}-{\mathrm{P}}_{\mathrm{jj}}}$$

we can iterate the Eqs. (49) and (50) to obtain the filtering probabilities $\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)$ and the conditional densities $\mathcal{F}\left({\mathrm{z}}_{\mathrm{t}}|{\mathrm{s}}_{\mathrm{t}}=0{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)$ for $\mathrm{t}=1,2,\dots ..\mathrm{T}$. Then we can compute the logarithmic likelihood function using

$$\text{log}\left(\mathrm{L}\left(\widehat{\mathrm{\theta }}\right)\right)=\sum _{\mathrm{t}=1}^{\mathrm{T}}\sum _{\mathrm{i}=1}^2\text{log}(\mathcal{F}\left({\mathrm{Z}}_{\mathrm{t}}\left|{\mathrm{S}}_{\mathrm{t}}=\mathrm{i},{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)\times \mathrm{P}\left({\mathrm{S}}_{\mathrm{t}}=\mathrm{i}|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)\right)$$E51

where $\mathrm{L}\left(\widehat{\mathrm{\theta }}\right)$ is the maximized value of the likelihood function. The model estimation can finally be obtained by finding the set of parameters $\widehat{\mathrm{\theta }}$ that maximize the Eq. (51) using numerical-search algorithm. The estimated filtering and prediction probabilities can then be easily calculated by plugging $\widehat{\mathrm{\theta }}$ into the equation formulae of these probabilities. We adopt the approximation in Ref [20] for computing the smoothing probabilities $\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)$

$$\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{s}}_{\mathrm{t}+1}=\mathrm{j},{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)\approx \mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{s}}_{\mathrm{t}+1},{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)$$

$$=\frac{\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i},{\mathrm{s}}_{\mathrm{t}+1}\left|{\mathrm{\psi }}_{\mathrm{t}-1};\mathrm{\theta }\right)}{\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=\mathrm{j}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)}$$

$$=\frac{{\mathrm{P}}_{0\mathrm{i}}\mathrm{P}\left({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)}{\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=\mathrm{j}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)}$$

And, for $\mathrm{i};\mathrm{j}=0;1,$ smoothing probabilities is expressed as:

$$\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)$$

$$=\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=0\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=1\left|{\mathrm{s}}_{\mathrm{t}+1}=0,{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)$$

$$+\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=1\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{s}}_{\mathrm{t}+1}=1,{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)$$

$$=\mathrm{P}({\mathrm{s}}_{\mathrm{t}}=\mathrm{i}\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)\times \left(\frac{{\mathrm{P}}_{\mathrm{i}0}\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=0\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)}{\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=0\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)}+\frac{{\mathrm{P}}_{\mathrm{i}1}\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=1\left|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)}{\mathrm{P}({\mathrm{s}}_{\mathrm{t}+1}=1\left|{\mathrm{\psi }}_{\mathrm{t}};\mathrm{\theta }\right)}\right)$$E52

Using $\mathrm{P}\left({\mathrm{S}}_{\mathrm{L}}=\mathrm{i}|{\mathrm{\psi }}_{\mathrm{L}};\mathrm{\theta }\right)$ as the initial value, we can iterate the equations regressively for filtering and prediction probabilities along with the equation above to get the smoothing probabilities for $\mathrm{t}=\mathrm{L}-1,··,\mathrm{k}+1$.

3.4. Results and discussions

Figure 8 shows the training detection sequence which we generate as a training observation using randomly distributed PU channel state “idle and occupied” over T = 250 ms simulation time. We use this training observations to train Baum-Welch algorithm in order to estimate HMM model parameters $\mathrm{\theta }=\left(\mathrm{\pi },\mathrm{A},\mathrm{B}\right)$, assuming that the first estimate of $\mathrm{\theta }\left(1\right)$ is:

$${\mathrm{\pi }}_1=\left(10\right)$$

$${\mathrm{A}}_1=\left(\begin{array}{cc}0.85& 0.15\\ 0.10& 0.90\end{array}\right)$$

$${\mathrm{B}}_1=\left(\begin{array}{ccccccc}0.17& 0.16& 0.17& 0.16& 0.17& 0.17& 0.17\\ 0.60& 0.08& 0.08& 0.08& 0.08& 0.08& 0.08\end{array}\right)$$

Figure 8.

The training detection sequence for HMM and MSM.

Figure 9a shows the performance of HMM algorithm in estimating the PU channel states (i.e., PU idle or PU occupied) of the time series that capture the detection sequence for a single-user cognitive radio network. Figure 9a contains three plots; the top plot shows the randomly distributed PU channel states over time $\mathrm{T}=500\mathrm{ms}$ . The middle plot shows the generated time series following the distribution ${\mathrm{z}}_{\mathrm{t}}\in \left\{1,2,3\right\}$ for idle states and ${\mathrm{z}}_{\mathrm{t}}\in \left\{4,5,6\right\}$ for occupied states (note: we can construct the observation space from these two distributions as ${\mathrm{O}}_{\mathrm{t}}\in \left\{1,2,3,4,5,6\right\},\mathrm{t}=1,2\dots 500\mathrm{ms}$). The bottom plot shows performance of HMM algorithm in forecasting the time series generated to capture PU detection sequence.

Figure 9.

(a) Shows the performance of HMM algorithm in predicting the generated time series to capture PU channel state detection sequence. (b) Shows the performance of MSM algorithm in predicting the generated time series to capture the same PU detection sequence in Figure 8.

Figure 9b shows the performance of MSM algorithm in predicting the switching process between the two PU channel states for the same PU detection sequence given in Figure 9a $\mathrm{T}=500\mathrm{ms}$ . The top graph in Figure 9b shows the generated time series with the following distribution ${\mathrm{z}}_{\mathrm{t}}\sim \left(\mathrm{0.1,0.5}\right)$ for idle states and ${\mathrm{z}}_{\mathrm{t}}\sim \left(\mathrm{0.01,0.2}\right)$ for occupied states and the bottom graph shows the prediction performance using MSM. As it is clear from the figure, the prediction performance is smoother than HMM approach.

4. Conclusions

In this chapter, we have presented a per-frame decision-based cooperative spectrum sensing based on machine-learning classifier-based fusion approach. The simulation and numerical results have shown that the machine-learning classifier-based fusion algorithm performs same as conventional fusion rules in terms of sensing accuracy with less sensing time, overheads, and extra operations that limit achievable cooperative gain among cognitive radio users. In addition, we have also studied the problem of primary user channel state prediction in cognitive radio network and introduced Markov model and Markov Switching Model to solve this problem. We finally showed by the means of simulation that both hidden Markov model and Markov switching model perform very well in predicting the time series that capture the actual primary user channel state.

Acknowledgments

I gratefully acknowledge the funding received from Shanghai Jiao Tong University to undertake my PhD. I also thank Prof. Bin Guo for his encouragement and help on the topic.