Mathuranathan, Author at GaussianWaves

Window function – figure of merits

Key focus: Window function smooths the observed signal over the edges. Analysis of some important parameters to help select the window for an application.

Spectral leakage

As we know, the DFT operation can be viewed as processing a signal through a set of filter banks with bandwidth Δf centered on the bin (frequency) of interest (Figure 1). To accurately resolve the frequency component in each bin, we desire the resolution bandwidth Δf to be as small as possible. In effect, each filter filter bank is viewed as filter having a brick-wall response with narrow bandwidth. If the signal under study contains several spectral components and broadband noise, the computed amplitude of the DFT output in each bin (filter bank) is significantly affected by the noise contained in the corresponding bandwidth of the filter bank.

*Figure 1: DFT (Discrete Fourier Transform) viewed as processing through a set of filter banks*

*Figure 1: DFT (Discrete Fourier Transform) viewed as processing through a set of filter banks*

In reality, signals are of time-limited nature and nothing can be known about the signal beyond the measured interval. When the time-limited signal slice is subjected to analysis using FFT algorithm (used for computing DFT), the FFT implicitly assumes that the signal essentially repeats itself after the observed interval. This may lead to discontinuities at the edges of each slice, that causes the energy contained in each frequency bin to spill into other bins. This phenomenon is called spectral leakage.

Equivalent Noise Bandwidth:^[2]

To suppress the spectral leakage, the signals are multiplied with a window function so as to smooth the discontinuity at the edges of the FFT slices. As a result, the choice of window function affects the amount of signal and noise that goes inside each filter bank. Hence the amount of noise that gets accumulated into each filter bank is affected by the bandwidth of the chosen window function. In effect, the amplitude estimate in each frequency bin is influenced by the amount of accumulated noise contributed by the window function of choice. In other words, the noise floor displayed at the FFT output, varies according to the window of choice.

Equivalent noise bandwidth (ENBW), often specified in terms of FFT bins, is the bandwidth of a fictitious brick-wall filter such that the amount of noise power accumulated inside the brick-wall filter is same as the noise power accumulated when processing white noise through the selected window.

*Figure 2: Illustrating equivalent noise bandwidth (ENBW) (figure not to scale*)

In discrete domain, the equivalent noise bandwidth B_enbw can be computed using time domain samples as

\[B_{enbw} = L \frac{\sum_{n}|w[n]|^2}{ \left| \sum_{n} w[n] \right|^2}\]

where L is the length of the window.

ENBW is a commonly used metric to characterize the variation in the displayed noise floor if a non rectangular window is used. If noise floor measurements are desired, the ENBW correction factor need to be applied to account for the variation in the noise floor [3]. Table 1 below, lists the ENBW (bins) and ENBW correction factor for some well known window functions.

Following figure illustrates the equivalent noise bandwidth illustrated for Hann window.

*Figure 3: Equivalent noise bandwidth (ENBW) illustrated for Hann window*

Coherent power gain:

In time-domain, windowing a signal reduces the amplitude of the signal at its boundaries (Figure 4). This causes a reduction in the amplitude of the spectral component when FFT used to visualize the signal in frequency domain. This reduction in amplitude in the spectral components is characterized as coherent power gain (which actually is a loss). So, when FFT is computed on a windowed signal, in order to compensate for the loss in the amplitude, we simply add the coherent power gain to the FFT output.

*Figure 4: Windowing a 10 Hz sinusoidal signal with Hann window of length L = 151 (sampling frequency F_s=160 Hz)*

Coherent power gain of a window w[n] of length L is given by

\[\text{Coherent power gain (dB)} = 20 \; log_{10} \left( \frac{\sum_n w[n]}{L} \right)\]

Following plot depicts the coherent power gain (i.e, the reduction in the FFT magnitude when the input signal is processed with a window) of a windowed sinusoidal signal of frequency 10 Hz. The length of the window is L = 151 points and the simulation assumes an oversampling factor of 16 (i.e, F_s=160 Hz). The FFT length is N_FFT =2048.

Table 1 below, lists the coherent power gain for some well known window functions.

*Figure 5: Coherent power gain – illustrated for various window functions*

Scalloping loss

As discussed above, the FFT (equivalently the DFT) operation can be seen as processing a signal through a set of filter banks with bandwidth Δf centered on the successive bin frequencies. The frequency resolution (Δf) of FFT depends on the size of FFT (N_FFT) and the sampling frequency F_s.

\[\Delta f = \frac{F_s}{N_{FFT}}\]

Therefore, the FFT process can measure amplitude of tones that are multiples of the frequency resolution Δf. That is, if the frequency of the observed signal matches a particular bin frequency, the signal amplitude will be maximum at that bin.

However, if the tone of the signal falls between two bins, the signal power is spread between the adjacent bins and hence the displayed signal amplitude/power is reduced.

Scalloping loss quantifies the worst case reduction in the signal level, when the tone of the observed signal falls exactly mid way between two bin frequencies. The worst case scalloping loss is calculated as

\[\text{Scalloping loss} = \frac{\sum_n w[n] e^{\left( -j \frac{\pi}{N_{FFT}}n \right)}}{\sum_n w[n]}\]

In reality, scalloping loss depends on the type of window and the relationship between the signal tones to the bin frequencies. Table 1 below, lists the scalloping loss for some well known window functions.

*Figure 6: Scalloping loss illustrated for processing a signal with Boxcar and Bartlett window*

Figure of merits for window functions

*Table 1: Figure of merits for various window functions*

List of window functions

Boxcar (Rectangular)

\[w[n] = \begin{cases} 1, \quad 0 \leq n \leq L-1 \\ 0, \quad \text{otherwise} \end{cases}\]

Figure 1: Boxcarr (rectangular) window *(L=51)*

Barthann

\[ w[n] = 0.62 – 0.48 \left| \frac{n}{L – 1} – 0.5 \right| + 0.38 cos \left(2 \pi \left| \frac{n}{L – 1} – 0.5 \right| \right), \quad 0 \leq n \leq L-1\]

Bartlett

\[w[n] = \frac{2}{L-1} \left(\frac{L-1}{2} – \left|n – \frac{L-1}{2}\right| \right),\quad 0 \leq n \leq L-1\]

Blackman

\[w[n] = 0.42 – 0.5 \cos(2\pi n/L) + 0.08 \cos(4\pi n/L), \quad 0 \leq n \leq L-1\]

Blackman-Harris

\[w[n]=a_0 – a_1 \cos \left ( \frac{2 \pi n}{L} \right)+ a_2 \cos \left ( \frac{4 \pi n}{L} \right)- a_3 \cos \left ( \frac{6 \pi n}{L} \right), \quad 0 \leq n \leq L-1\] \[a_0=0.35875;\quad a_1=0.48829;\quad a_2=0.14128;\quad a_3=0.01168\]

Bohman

\[w[n] = \left( 1 – |x|\right) cos \left( \pi |x| \right) + \frac{1}{\pi} sin \left( \pi |x| \right) , \quad -1 \leq x \leq 1 \]

Cosine

\[w[n] = sin \left( \frac{\pi}{ L } \left(n + 0.5 \right) \right),\quad 0 \leq n \leq L-1\]

Flattop

\[w[n]=a_0 – a_1 \cos \left ( \frac{2 \pi n}{L} \right)+ a_2 \cos \left ( \frac{4 \pi n}{L} \right)- a_3 \cos \left ( \frac{6 \pi n}{L} \right), \quad 0 \leq n \leq L-1 \\ a_0=0.21557895;\; a_1=0.41663158 \\ a_2=0.277263158;\; a_3=0.006947368 \]

Hamming

\[w[n] = 0.54 – 0.46 \cos\left(\frac{2\pi{n}}{L-1}\right), \quad 0 \leq n \leq L-1\]

Hann

\[w[n] = 0.5 – 0.5 \cos\left(\frac{2\pi{n}}{L-1}\right), \quad 0 \leq n \leq L-1\]

Nuttall

\[ w[n]=a_0 – a_1 \cos \left ( \frac{2 \pi n}{L} \right)+ a_2 \cos \left ( \frac{4 \pi n}{L} \right)- a_3 \cos \left ( \frac{6 \pi n}{L} \right), \quad 0 \leq n \leq L-1 \\ a_0=0.3635819;\; a_1=0.4891775;\\ a_2=0.1365995;\; a_3=0.0106411\]

Parzen

\[w[n] = \begin{cases} 1- 6 \left( \frac{|n|}{L/2} \right)^2 + 6 \left( \frac{|n|}{L/2} \right)^3, & 0 \leq |n| \leq (L-1)/4 \\ 2 \left( 1 – \frac{|n|}{L/2} \right)^3, & (L-1)/4 < |n| \leq (L-1)/2 \end{cases}\]

Triangular

\[w[n] = 1 – \frac{|n|}{L/2}, \quad n = -L/2, \cdots, -1, 0, 1, \cdots, L/2\]

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

References:

[1] Harris, “On the Use of Windows for Harmonic Analysis with the Discrete Fourier Transform”, proceedings of IEEE, vol. 66, no. 1, January 1978 .↗
[2] Stefan Scholl, “Exact Signal Measurements using FFT Analysis”,Microelectronic Systems Design Research Group, TU Kaiserslautern, Germany.↗

Additional Resources:

[1] If you are particularly interested in calculating the Equivalent Noise Bandwidth of Analog Filters refer – Thomas J. Karras,”Equivalent Noise Bandwidth Analysis from Transfer Functions”,NASA Technical Note,NASA TN D-2842.↗
[2] Strategies for Choosing Windows : Michael Cerna and Audrey F. Harvey,”The Fundamentals of FFT-Based Signal Analysis and Measurement”,National Instruments Application Notes 041.↗

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Equivalent noise bandwidth (ENBW) of window functions

Key focus: Equivalent noise bandwidth (ENBW), is the bandwidth of a fictitious brick-wall filter that allows same amount of noise as a window function. Learn how to calculate ENBW in applications involving window functions and FFT operation.

FFT and spectral leakage

Here is an article with illustrations on – FFT slices and the phenomenon of spectral leakage,

Hence, to suppress the spectral leakage, the signals are multiplied with a window function so as to smooth the discontinuity at the edges of the FFT slices. As a result, the choice of window function affects the amount of signal and noise that goes inside each filter bank. Hence the amount of noise that gets accumulated into each filter bank is affected by the bandwidth of the chosen window function. In effect, the amplitude estimate in each frequency bin is influenced by the amount of accumulated noise contributed by the window function of choice. In other words, the noise floor displayed at the FFT output, varies according to the window of choice.

Equivalent noise bandwidth

In other words, in frequency domain, the power accumulated in the chosen window is given by

\[P_{window} = \int_{-f}^{f} |W(f)|^2 df\]

Then we wish to find the equivalent noise bandwidth B_enbw, such that

\[P_{brick-wall} = B_{enbw} |W(f_0)|^2 = P_{window}\]

where f₀ is the frequency at which the power peaks and P_brick-wall is the power accumulated in the fictitious rectangular window of bandwidth B_enbw.

Therefore, the equivalent noise bandwidth B_enbw is given by

\[B_{enbw} = \frac{\int_{-f}^{f} |W(f)|^2 df}{|W(f_0)|^2}\]

Translating to discrete domain, the equivalent noise bandwidth can be computed using DFT samples as

\[B_{enbw} =\frac{\sum_{k=0}^{N-1}|W[k]|^2}{|W[k_0]|^2}\]

where, k₀ is the index at which the magnitude of FFT output is maximum and N is the window length. Applying Parseval’s theorem and with a window of length L, B_enbw can also be computed using time domain samples as

\[B_{enbw} = L \frac{\sum_{n}|w[n]|^2}{ \left| \sum_{n} w[n] \right|^2}\]

Following figure illustrates the equivalent noise bandwidth illustrated for Hann window.

Python script

ENBW can be calculated from the time domain samples in a straightforward manner. Following Python 3.0 script calculates the ENBW for some well-known window functions provided as part of Scipy module↗.

#Author : Mathuranathan Viswanathan for gaussianwaves.com
#Date: 13 Sept 2020
#Script to calculate equivalent noise bandwidth (ENBW) for some well known window functions

import numpy as np
import pandas as pd
from scipy import signal

def equivalent_noise_bandwidth(window):
    #Returns the Equivalent Noise BandWidth (ENBW)
    return len(window) * np.sum(window**2) / np.sum(window)**2

def get_enbw_windows():
    #Return ENBW for all the following windows as a dataframe
    window_names = ['boxcar','barthann','bartlett','blackman','blackmanharris','bohman','cosine','exponential','flattop','hamming','hann','nuttall','parzen','triang']
    
    df = pd.DataFrame(columns=['Window','ENBW (bins)','ENBW correction (dB)'])
    for window_name in window_names:
        method_name = window_name
        func_to_run = getattr(signal, method_name) #map window names to window functions in scipy package
        L = 16384 #Number of points in the output window
        window = func_to_run(L) #call the functions
        
        enbw = equivalent_noise_bandwidth(window) #compute ENBW
        
        df = df.append({'Window': window_name.title(),'ENBW (bins)':round(enbw,3),'ENBW correction (dB)': round(10*np.log10(enbw),3)},ignore_index=True)
                       
    return df

df = get_enbw_windows() #call the function
df.head(14) #display dataframe

The resulting output is displayed in the following table (dataframe)

*Table 1: Table of equivalent noise bandwidth for some well known window functions*

As an example, the following plot depicts the difference in the noise floor of FFT output of a noise added sine wave that is processed through Boxcar and Flattop window.

*Figure 4: Noise floor and ENBW for flattop & boxcar window (FFT output) for noise added 10 Hz sinewave (oversampling factor = 16, Fs = 160 Hz, window length L =2048, SNR = 30 dB)*

Continue reading on Window function and figure of merits…

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

Reference

[1] Stefan Scholl, “Exact Signal Measurements using FFT Analysis”,Microelectronic Systems Design Research Group, TU Kaiserslautern, Germany.↗

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Parseval’s theorem – derivation

The Parseval’s theorem (a.k.a Plancherel theorem) expresses the energy of a signal in time-domain in terms of the average energy in its frequency components.

Suppose if the x[n] is a discrete-time sequence of complex numbers of length N : x_n={x₀,x₁,…,x_N-1}, its N-point discrete Fourier transform (DFT)[1] : X_k={X₀,X₁,…,X_N-1} is given by

\[X[k] = \sum_{n=0}^{N-1} x[n] e^{-j\frac{2 \pi}{N} k n} \]

The inverse discrete Fourier transform is given by

\[\tilde{x}[n] = \frac{1}{N} \sum_{k=0}^{N-1} X[k] e^{j \frac{2 \pi}{N} kn}\]

Suppose if x[n] and y[n] are two such sequences that follows the above definitions, the Parseval’s theorem is written as

\[\boxed{ \sum_{n=0}^{N-1} x[n] y^{\ast}[n] = \frac{1}{N} \sum_{k=0}^{N-1} X[k] Y^{\ast}[k]} \]

where, * indicates conjugate operation.

Deriving Parseval’s theorem

\[\begin{aligned} \sum_{n=0}^{N-1} x[n] y^{\ast}[n] &= \sum_{n=0}^{N-1} x[n] \left(\frac{1}{N} \sum_{k=0}^{N-1} Y[k] e^{j\frac{2 \pi}{N} k n} \right )^\ast \\ &= \frac{1}{N}\sum_{n=0}^{N-1} x[n] \sum_{k=0}^{N-1} Y^\ast[k] e^{-j\frac{2 \pi}{N} k n} \\ &= \frac{1}{N} \sum_{k=0}^{N-1} Y^\ast[k] \cdot \sum_{n=0}^{N-1} x[n] e^{-j\frac{2 \pi}{N} k n} \\ &= \frac{1}{N} \sum_{k=0}^{N-1} X[k] Y^\ast[k] \end{aligned}\]

If x[n] = y[n], then

\[ \begin{aligned} \sum_{n=0}^{N-1} x[n] x^{\ast}[n] &= \frac{1}{N} \sum_{k=0}^{N-1} X[k] X^\ast[k] \\ \Rightarrow \sum_{n=0}^{N-1} |x[n]|^2 &= \frac{1}{N} \sum_{k=0}^{N-1} |X[k]|^2\end{aligned} \]

Time-domain and frequency domain representations are equivalent manifestations of the same signal. Therefore, the energy of the signal computed from time domain samples must be equal to the total energy computed from frequency domain representation.

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

Reference

[1] Bertrand Delgutte and Julie Greenberg , “The discrete Fourier Transform”, Biomedical Signal and Image Processing Spring 2005, MIT web

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Rician flat-fading channel – simulation

In wireless environments, transmitted signal may be subjected to multiple scatterings before arriving at the receiver. This gives rise to random fluctuations in the received signal and this phenomenon is called fading. The scattered version of the signal is designated as non line of sight (NLOS) component. If the number of NLOS components are sufficiently large, the fading process is approximated as the sum of large number of complex Gaussian process whose probability-density-function follows Rayleigh distribution.

Rayleigh distribution is well suited for the absence of a dominant line of sight (LOS) path between the transmitter and the receiver. If a line of sight path do exist, the envelope distribution is no longer Rayleigh, but Rician (or Ricean). If there exists a dominant LOS component, the fading process can be represented as the sum of complex exponential and a narrowband complex Gaussian process g(t). If the LOS component arrive at the receiver at an angle of arrival (AoA) θ, phase ɸ and with the maximum Doppler frequency f_D, the fading process in baseband can be represented as (refer [1])

\[h(t)= \underbrace{\sqrt{\frac{K \Omega}{K +1}}}_\text{A:=} e^{\left( j2 \pi f_D cos(\theta)t+\phi \right)} + \underbrace{\sqrt{\frac{\Omega}{K+1}}}_\text{S:=}g(t)\]

where, K represents the Rician K factor given as the ratio of power of the LOS component A² to the power of the scattered components (S²) marked in the equation above.

\[K=\frac{A^2}{S^2}\]

The received signal power Ω is the sum of power in LOS component and the power in scattered components, given as Ω=A²+S². The above mentioned fading process is called Rician fading process. The best and worst-case Rician fading channels are associated with K=∞ and K=0 respectively. A Ricean fading channel with K=∞ is a Gaussian channel with a strong LOS path. Ricean channel with K=0 represents a Rayleigh channel with no LOS path.

The statistical model for generating flat-fading Rician samples is discussed in detail in chapter 11 section 11.3.1 in the book Wireless communication systems in Matlab (see the related article here). With respect to the simulation model shown in Figure 1(b), given a K factor, the samples for the Rician flat-fading samples are drawn from the following random variable

\[h= | X + jY |\]

where X,Y ~ N(μ,σ²) are Gaussian random variables with non-zero mean μ and standard deviation σ as given in references [2] and [3].

\[\mu = g_1 =\sqrt{\frac{K}{2\left(K+1\right)}} \quad \quad \sigma = g_2 = \sqrt{\frac{1}{2\left(K+1\right)}}\]

Kindly refer the book Wireless communication systems in Matlab for the script on generating channel samples for Ricean flat-fading.

*Figure 1: Simulation model for modulation and detection over flat fading channel*

Simulation and performance results

In chapter 5 of the book Wireless communication systems in Matlab, the code implementation for complex baseband models for various digital modulators and demodulator are given. The computation and generation of AWGN noise is also given in the book. Using these models, we can create a unified simulation for code for simulating the performance of various modulation techniques over Rician flat-fading channel the simulation model shown in Figure 1(b).

An unified approach is employed to simulate the performance of any of the given modulation technique – MPSK, MQAM or MPAM. The simulation code (given in the book) will automatically choose the selected modulation type, performs Monte Carlo simulation, computes symbol error rates and plots them against the theoretical symbol error rate curves. The simulated performance results obtained for various modulations are shown in the Figure 2.

*Figure 2: Performance of various modulations over Ricean flat fading channel*

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

References

[1] C. Tepedelenlioglu, A. Abdi, and G. B. Giannakis, The Ricean K factor: Estimation and performance analysis, IEEE Trans. Wireless Communication ,vol. 2, no. 4, pp. 799–810, Jul. 2003.↗
[2] R. F. Lopes, I. Glover, M. P. Sousa, W. T. A. Lopes, and M. S. de Alencar, A simulation framework for spectrum sensing, 13th International Symposium on Wireless Personal Multimedia Communications (WPMC 2010), Out. 2010.
[3] M. C. Jeruchim, P. Balaban, and K. S. Shanmugan, Simulation of Communication Systems, Methodology, Modeling, and Techniques, second edition Kluwer Academic Publishers, 2000.↗

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Statistical measures for stochastic signals

Key focus: Discuss statistical measures for stochastic signals : mean, variance, skewness, kurtosis, histogram, scatterplot, cross-correlation and auto-correlation.

Deterministic and stochastic signals

A deterministic signal is exactly predictable for the given time span of interest. It could be expressed using analytic form (example: x(t) = sin (2 π f_c t) ).

Many of the signals that are encountered in real world signal processing applications cannot be expressed or predicted using analytic equations, because they have an element of uncertainty associated with them. As a consequence, such signals are characterized using statistical concepts. Therefore, these signals are outside the realm of deterministic signals and are classified as stochastic signals.

For example, we look at an electrical circuit and monitor the voltage across a resistor for a period of time. Under an applied electric field, atomic particles inside resister tend to randomly move and it manifests as heat. This random thermal motion causes random fluctuation in the voltage measured across the resistor. Therefore, the measured voltage can be characterized by a probability distribution and can be analyzed using statistical methods, but it cannot be predicted with precision. This is an example of signal that is stochastic function of time. Such functions usually evolve according to certain probabilistic laws and are assumed to be generated by an underlying stochastic process (thermal motion in the example above).

*Figure 1: Examples for deterministic and stochastic signals*

Given the amplitude of the stochastic voltage signal at time , now we know, we cannot predict the value at t’. However, if we observed the signal for a sufficient amount of time, we can empirically determine its probability distribution based on which we should be able to answer questions like

Given the amplitude of the voltage at time t, what is the average (expected or mean) of the voltage at time t’?
How much can we expect the voltage at time t’, to fluctuate from the mean ? In other words, we are interested in the variance of the voltage at time t’.
What is the probability that the voltage at t’ exceeds or falls below a certain threshold value ?
How the voltages at times t and t’ are related ? In other words, we are interested in correlation.

Summary of descriptive statistical measures

Different statistical measures are available to gather, review, analyze and draw conclusions from stochastic signals. Some of the statistical measures are summarized below:

Quantitative measures of shape:

In many statistical analysis, the fundamental task is to get a general idea about the shape of a distribution and it is done by using moments.

Measure of central tendency – mean – the first moment:

Measures of central tendency attempt to identify the central position in the distribution of samples that make up the stochastic signal. Mean, mode and median are different measures of central tendency. In signal processing, we are mostly interested in computing mean which is the average of values of the given samples. Given a discrete random variable X with probability mass function p_X(x) the mean is denoted by

\[\mu = E \left[ X\right] = \sum_{x: p_X(x) > 0} x p_X(x) \]

Measure of dispersion – variance – the second moment:

Measures of dispersion describe how the values in the signal samples are spread out around a central location. Variance and standard deviation are some of the measures of dispersion. If the values are widely dispersed, the central location is said to be less representative of the values as a whole. If the values are tightly dispersed, the central location is considered more reliable.

For a discrete random variable X with probability mass function p_X(x) , the variance (σ_X²) is given by the second central moment. The term central moment implies that this is measured relative to mean. For an electrical signal, the second moment is proportional to the average power.

\[\sigma_X^2 = E \left[ \left(X – \mu \right)^2\right] = \sum_{x: p_X(x) > 0} \left(x – \mu \right)^2 p_X(x) \]

The square root of variance is standard deviation

\[\sigma_X = \sqrt{E \left[ \left(X – \mu \right)^2\right]}\]

Figure 2, demonstrates the use of histogram for getting a visual perception of shape of the distribution of samples from two different stochastic signals. Though the signals look similar in nature in time domain, their histogram reveals different picture altogether. The central location (mean) of the first signal is around zero (no DC shift in the time domain view) and for the second signal the mean is at 0.75. The average power of first signal varies widely (histogram is widely spread out) compared to that of the second signal.

*Figure 2: Histogram is a visual method that provides a general idea about the shape of a distribution.*

Higher order moments – skewness and kurtosis:

Further characterization of the shape of a distribution includes higher order moments : skewness and kurtosis. They identify anomalies and outliers in many signal processing applications.

Skewness provides a measure to quantify the presence of asymmetry in the shape of the distribution. Actually, skewness measures the relative size of the tails in a distribution. Presence of asymmetry manifests as non-zero value for the standardized third central moment. The term “standardized” stands for the normalization of the third central moment by σ³.

\[Skewness = \frac{E \left[ \left( X – \mu \right)^3 \right]}{\sigma^3}\]

*Figure 3: A random sequence showing positive skewness (asymmetry in shape) in its distribution*

Kurtosis measures the amount of probability in the two tails of a distribution, relative to a normal distribution of same variance. Kurtosis is 3 for normal distribution (perfect bell shaped curve). If the kurtosis of a distribution is greater than 3, it implies that the tails are heavier compared to that of normal distribution. A value less than 3 for kurtosis implies lighter tails compared to that of the normal distribution. Essentially, kurtosis is a measure of outliers. Kurtosis is calculated as the standardized fourth central moment.

\[Kurtosis = \frac{E \left[ \left( X – \mu \right)^4 \right]}{\sigma^4}\]

*Figure 4: Histogram showing kurtosis of a random sequence vs. kurtosis of normal distribution*

Measures of association

Statistics, such as measures of central tendency, dispersion and higher order moments, describe about a single distribution are called univariate statistics. If we are interested in the relationship between two or more variables, we have to move to at least the realm of bivariate statistics.

Measures of association, summarize the size of association between two variables. Most measures of associates are scaled to a range of values. For example, a measure of association can be construed to have a range 0 to 1, where the value 0 implies no relationship and a value of 1 implies perfect relationship between the two variables under study. In another case, the measure can range from -1 to 1, which can help us determine if the two variables have negative or positive relationship between each other.

Scatter plot

Scatter plot is a great way to visually assess the nature of relationship between two variables. Following figure contains the scatter plots between different stochastic signals, exhibiting different strengths of relationships between the signals.

Correlation

Correlation functions are commonly used in signal processing, for measuring the similarity of two signals. They are especially used in signal detection and pattern recognition.

Cross-correlation

Cross-correlation is commonly used for measuring the similarity between two different signals. In signal processing, cross-correlation measures the similarity of two waveforms as a function of time lags applied to one of the waveform.

For discrete-time waveforms – x[n] and y[n], cross correlation is defined as

\[Corr_{xy}[l] = \sum_{m=-\infty}^{\infty} x[n]^{\ast} y[n+l] = \sum_{m=-\infty}^{\infty} x[n-l]^{\ast} y[n]\]

where, * denotes complex conjugate operation and l is the discrete time lags applied to one of the waveforms. That is, the cross-correlation can be calculated as the dot product of a sequence with each shifted version of another sequence.

Auto-correlation

Auto-correlation is the cross-correlation of a waveform/sequence with itself.

For discrete-time waveform – x[n], auto-correlation is defined as

\[Corr_{xx}[l] = \sum_{m=-\infty}^{\infty} x[n]^{\ast} x[n+l] = \sum_{m=-\infty}^{\infty} x[n-l]^{\ast} x[n] \]

where, * denotes complex conjugate operation and l is the discrete time lags applied to the copy of the same waveform.

Correlation properties are useful for identifying/distinguishing a known bit sequence from a set of other possible known sequences. For example, in GPS technology, each satellite is assigned a unique 10-bit Gold code sequence (2¹⁰ = 1023 possible combinations). Cross-correlation between different Gold code sequence is very low, since the Gold codes are orthogonal to each other. However, the auto-correlation of a Gold code is maximum at zero lag. The satellites are identified using these correlation properties of Gold codes. (I have described the hardware implementation of Gold codes, it can be accessed here).

Following plots illustrate the application of auto-correlation for audio analysis. The auto-correlation of an audio signal will have a peak at zero lag (i.e, where there is no time shifting when computing the correlation) as shown in Figure 6. Figure 7 contains the same audio file that is synthetically processed to produce reverberation characteristics. By looking at the time series plot in Figure 7, we cannot infer anything. However, the auto-correlation plot reveals the reverberation characteristics embedded in the audio.

*Figure 6: Normalized auto-correlation of an original sound*

*Figure 7: Normalized auto-correlation of sound with simulated reverberation characteristics (synthetically processed through an IIR comb filter)*

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

For further reading

[1] Steven M. Kay, “Fundamentals of Statistical Signal Processing, Volume I: Estimation Theory”, ISBN: 978-0133457117, Prentice Hall, Edition 1, 1993.↗

BPSK bit error rate simulation in Python & Matlab

Key focus: Simulate bit error rate performance of Binary Phase Shift Keying (BPSK) modulation over AWGN channel using complex baseband equivalent model in Python & Matlab.

Why complex baseband equivalent model

The passband model and equivalent baseband model are fundamental models for simulating a communication system. In the passband model, also called as waveform simulation model, the transmitted signal, channel noise and the received signal are all represented by samples of waveforms. Since every detail of the RF carrier gets simulated, it consumes more memory and time.

In the case of discrete-time equivalent baseband model, only the value of a symbol at the symbol-sampling time instant is considered. Therefore, it consumes less memory and yields results in a very short span of time when compared to the passband models. Such models operate near zero frequency, suppressing the RF carrier and hence the number of samples required for simulation is greatly reduced. They are more suitable for performance analysis simulations. If the behavior of the system is well understood, the model can be simplified further.

Passband model and converting it to equivalent complex baseband model is discussed in this article.

Simulation of bit error rate performance of BPSK using passband simulation model is discussed in this article.

BPSK constellation

In binary phase shift keying, all the information gets encoded in the phase of the carrier signal. The BPSK modulator accepts a series of information symbols drawn from the set m ∈ {0,1}, modulates them and transmits the modulated symbols over a channel.

The general expression for generating a M-PSK signal set is given by

Here, M denotes the modulation order and it defines the number of constellation points in the reference constellation. The value of M depends on the parameter k – the number of bits we wish to squeeze in a single M-PSK symbol. For example if we wish to squeeze in 3 bits (k=3) in one transmit symbol, then M = 2^k = 2³ = 8 and this results in 8-PSK configuration. M=2 gives BPSK (Binary Phase Shift Keying) configuration. The parameter A is the amplitude scaling factor, f_cis the carrier frequency and g(t) is the pulse shape that satisfies orthonormal properties of basis functions.

Using trigonometric identity, equation (1) can be separated into cosine and sine basis functions as follows

Therefore, the signaling set {s_i,s_q} or the constellation points for M-PSK modulation is given by,

For BPSK (M=2), the constellation points on the I-Q plane (Figure 1) are given by

Simulation methodology

Note: If you are interested in knowing more about BPSK modulation and demodulation, kindly visit this article.

In this simulation methodology, there is no need to simulate each and every sample of the BPSK waveform as per equation (1). Only the value of a symbol at the symbol-sampling time instant is considered. The steps for simulation of performance of BPSK over AWGN channel is as follows (Figure 2)

Generate a sequence of random bits of ones and zeros of certain length (N_sym typically set in the order of 10000)
Using the constellation points, map the bits to modulated symbols (For example, bit ‘0’ is mapped to amplitude value A, and bit ‘1’ is mapped to amplitude value -A)
Compute the total power in the sequence of modulated symbols and add noise for the given E_bN₀ (SNR) value (read this article on how to do this). The noise added symbols are the received symbols at the receiver.
Use thresholding technique, to detect the bits in the receiver. Based on the constellation diagram above, the detector at the receiver has to decide whether the receiver bit is above or below the threshold 0.
Compare the detected bits against the transmitted bits and compute the bit error rate (BER).
Plot the simulated BER against the SNR values and compare it with the theoretical BER curve for BPSK over AWGN (expressions for theoretical BER is available in this article)

*Figure 2: Simulation methodology for performance of BPSK modulation over AWGN channel*

Let’s simulate the performance of BPSK over AWGN channel in Python & Matlab.

Simulation using Python

Following standalone code simulates the bit error rate performance of BPSK modulation over AWGN using Python version 3. The results are plotted in Figure 3.

For more such examples refer the book (available as PDF and paperback) – Digital Modulations using Python

#Eb/N0 Vs BER for BPSK over AWGN (complex baseband model)
# © Author: Mathuranathan Viswanathan (gaussianwaves.com)
import numpy as np #for numerical computing
import matplotlib.pyplot as plt #for plotting functions
from scipy.special import erfc #erfc/Q function

#---------Input Fields------------------------
nSym = 10**5 # Number of symbols to transmit
EbN0dBs = np.arange(start=-4,stop = 13, step = 2) # Eb/N0 range in dB for simulation
BER_sim = np.zeros(len(EbN0dBs)) # simulated Bit error rates

M=2 #Number of points in BPSK constellation
m = np.arange(0,M) #all possible input symbols
A = 1; #amplitude
constellation = A*np.cos(m/M*2*np.pi)  #reference constellation for BPSK

#------------ Transmitter---------------
inputSyms = np.random.randint(low=0, high = M, size=nSym) #Random 1's and 0's as input to BPSK modulator
s = constellation[inputSyms] #modulated symbols

fig, ax1 = plt.subplots(nrows=1,ncols = 1)
ax1.plot(np.real(constellation),np.imag(constellation),'*')

#----------- Channel --------------
#Compute power in modulatedSyms and add AWGN noise for given SNRs
for j,EbN0dB in enumerate(EbN0dBs):
    gamma = 10**(EbN0dB/10) #SNRs to linear scale
    P=sum(abs(s)**2)/len(s) #Actual power in the vector
    N0=P/gamma # Find the noise spectral density
    n = np.sqrt(N0/2)*np.random.standard_normal(s.shape) # computed noise vector
    r = s + n # received signal
    
    #-------------- Receiver ------------
    detectedSyms = (r <= 0).astype(int) #thresolding at value 0
    BER_sim[j] = np.sum(detectedSyms != inputSyms)/nSym #calculate BER

BER_theory = 0.5*erfc(np.sqrt(10**(EbN0dBs/10)))

fig, ax = plt.subplots(nrows=1,ncols = 1)
ax.semilogy(EbN0dBs,BER_sim,color='r',marker='o',linestyle='',label='BPSK Sim')
ax.semilogy(EbN0dBs,BER_theory,marker='',linestyle='-',label='BPSK Theory')
ax.set_xlabel('$E_b/N_0(dB)$');ax.set_ylabel('BER ($P_b$)')
ax.set_title('Probability of Bit Error for BPSK over AWGN channel')
ax.set_xlim(-5,13);ax.grid(True);
ax.legend();plt.show()

*Figure 3: Bit error rate performance of BPSK modulation over AWGN using Python*

Simulation using Matlab

Following code simulates the bit error rate performance of BPSK modulation over AWGN using basic installation of Matlab. You will need the add_awgn_noise function that was discussed in this article. The results will be same as Figure 3.

For more such examples refer the book (available as PDF and paperback) – Digital Modulations using Matlab: build simulation models from scratch

%Eb/N0 Vs BER for BPSK over AWGN (complex baseband model)
% © Author: Mathuranathan Viswanathan (gaussianwaves.com)
clearvars; clc;
%---------Input Fields------------------------
nSym=10^6;%Number of symbols to transmit
EbN0dB = -4:2:14; % Eb/N0 range in dB for simulation

BER_sim = zeros(1,length(EbN0dB));%simulated Symbol error rates
    
M=2; %number of constellation points in BPSK
m = [0,1];%all possible input bits
A = 1; %amplitude
constellation = A*cos(m/M*2*pi);%constellation points

d=floor(M.*rand(1,nSym));%uniform random symbols from 1:M
s=constellation(d+1);%BPSK modulated symbols
    
for i=1:length(EbN0dB)
    r  = add_awgn_noise(s,EbN0dB(i));%add AWGN noise
    dCap = (r<=0);%threshold detector
    BER_sim(i) = sum((d~=dCap))/nSym;%SER computation
end

semilogy(EbN0dB,BER_sim,'-*');
xlabel('Eb/N0(dB)');ylabel('BER (Pb)');
title(['Probability of Bit Error for BPSK over AWGN']);

Rate this article: Note: There is a rating embedded within this post, please visit this post to rate it.

Reference

[1] Andrea Goldsmith, “Wireless Communications”, ISBN: 978-0521837163, Cambridge University Press; 1 edition, August 8, 2005.↗

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Constellation diagram – investigate phase transitions

The phase transition properties of the different variants of QPSK schemes and MSK, are easily investigated using constellation diagram. Let’s demonstrate how to plot the signal space constellations, for the various modulations used in the transmitter.

Typically, in practical applications, the baseband modulated waveforms are passed through a pulse shaping filter for combating the phenomenon of intersymbol interference (ISI). The goal is to plot the constellation plots of various pulse-shaped baseband waveforms of the QPSK, O-QPSK and π/4-DQPSK schemes. A variety of pulse shaping filters are available and raised cosine filter is specifically chosen for this demo. The raised cosine (RC) pulse comes with an adjustable transition band roll-off parameter α, using which the decay of the transition band can be controlled.

This article is part of the following books
● Digital Modulations using Matlab : Build Simulation Models from Scratch, ISBN: 978-1521493885
● Digital Modulations using Python ISBN: 978-1712321638
All books available in ebook (PDF) and Paperback formats

The RC pulse shaping function is expressed in frequency domain as

Equivalently, in time domain, the impulse response corresponds to

A simple evaluation of the equation (2) produces singularities (undefined points) at p(t = 0) and p(t = ±T_sym/(2α)). The value of the raised cosine pulse at these singularities can be obtained by applying L’Hospital’s rule [1] and the values are

Using the equations above, the raised cosine filter is implemented as a function (refer the books Digital Modulations using Python and Digital Modulations using Matlab for the code).

The function is then tested. It generates a raised cosine pulse for the given symbol duration T_sym = 1s and plots the time-domain view and the frequency response as shown in Figure 1. From the plot, it can be observed that the RC pulse falls off at the rate of 1/|t|³ as t→∞, which is a significant improvement when compared to the decay rate of a sinc pulse which is 1/|t|. It satisfies Nyquist criterion for zero ISI – the pulse hits zero crossings at desired sampling instants. The transition bands in the frequency domain can be made gradual (by controlling α) when compared to that of a sinc pulse.

***Figure 1: Raised-cosine pulse and its manifestation in frequency domain***

Plotting constellation diagram

Now that we have constructed a function for raised cosine pulse shaping filter, the next step is to generate modulated waveforms (using QPSK, O-QPSK and π/4-DQPSK schemes), pass them through a raised cosine filter having a roll-off factor, say α = 0.3 and finally plot the constellation. The constellation for MSK modulated waveform is also plotted.

**Figure 2: Constellations plots for: (a) a = 0.3 RC-filtered QPSK, (b) α = 0.3 RC-filtered O-QPSK, (c) α = 0.3 RC-filtered π/4-DQPSK and (d) MSK**

Conclusions

The resulting simulated plot is shown in the Figure 2. From the resulting constellation diagram, following conclusions can be reached.

Conventional QPSK has 180° phase transitions and hence it requires linear amplifiers with high Q factor
The phase transitions of Offset-QPSK are limited to 90° (the 180° phase transitions are eliminated)
The signaling points for π/4-DQPSK is toggled between two sets of QPSK constellations that are shifted by 45° with respect to each other. Both the 90° and 180° phase transitions are absent in this constellation. Therefore, this scheme produces the lower envelope variations than the rest of the two QPSK schemes.
MSK is a continuous phase modulation, therefore no abrupt phase transition occurs when a symbol changes. This is indicated by the smooth circle in the constellation plot. Hence, a band-limited MSK signal will not suffer any envelope variation, whereas, the rest of the QPSK schemes suffer varied levels of envelope variations, when they are band-limited.

References

[1] Clay S. Turner, Raised Cosine and Root Raised Cosine Formulae, Wireless Systems Engineering, Inc, (May 29, 2007) V1.2↗

In this chapter

Digital Modulators and Demodulators - Passband Simulation Models
● Introduction
● Binary Phase Shift Keying (BPSK)
□ BPSK transmitter
□ BPSK receiver
□ End-to-end simulation
● Coherent detection of Differentially Encoded BPSK (DEBPSK)
● Differential BPSK (D-BPSK)
□ Sub-optimum receiver for DBPSK
□ Optimum noncoherent receiver for DBPSK
● Quadrature Phase Shift Keying (QPSK)
□ QPSK transmitter
□ QPSK receiver
□ Performance simulation over AWGN
● Offset QPSK (O-QPSK)
● π/p=4-DQPSK
● Continuous Phase Modulation (CPM)
□ Motivation behind CPM
□ Continuous Phase Frequency Shift Keying (CPFSK) modulation
□ Minimum Shift Keying (MSK)
● Investigating phase transition properties
● Power Spectral Density (PSD) plots
● Gaussian Minimum Shift Keying (GMSK)
□ Pre-modulation Gaussian Low Pass Filter
□ Quadrature implementation of GMSK modulator
□ GMSK spectra
□ GMSK demodulator
□ Performance
● Frequency Shift Keying (FSK)
□ Binary-FSK (BFSK)
□ Orthogonality condition for non-coherent BFSK detection
□ Orthogonality condition for coherent BFSK
□ Modulator
□ Coherent Demodulator
□ Non-coherent Demodulator
□ Performance simulation
□ Power spectral density

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Matplotlib histogram and estimated PDF in Python

Key focus: Shown with examples: let’s estimate and plot the probability density function of a random variable using Python’s Matplotlib histogram function.

Generation of random variables with required probability distribution characteristic is of paramount importance in simulating a communication system. Let’s see how we can generate a simple random variable, estimate and plot the probability density function (PDF) from the generated data and then match it with the intended theoretical PDF. Normal random variable is considered here for illustration.

Step 1: Generate random samples

A survey of commonly used fundamental methods to generate a given random variable is given in [1]. For this demonstration, we will consider the normal random variable with the following parameters : μ – mean and σ – standard deviation. First generate a vector of randomly distributed random numbers of sufficient length (say 100000) with some valid values for μ and σ. There are more than one way to generate this. Two of them are given below.

● Method 1: Using the in-built numpy.random.normal() function (requires numpy package to be installed)

import numpy as np

mu=10;sigma=2.5 #mean=10,deviation=2.5
L=100000 #length of the random vector

#Random samples generated using numpy.random.normal()
samples_normal = np.random.normal(loc=mu,scale=sigma,size=(L,1)) #generate normally distributted samples

● Method 2: Box-Muller transformation [2] method produces a pair of normally distributed random numbers (Z₁, Z₂) by transforming a pair of uniformly distributed independent random samples (U₁,U₂). The algorithm for transformation is given by

#Samples generated using Box-Muller transformation

U1 = np.random.uniform(low=0,high=1,size=(L,1)) #uniformly distributed random numbers U(0,1)
U2 = np.random.uniform(low=0,high=1,size=(L,1)) #uniformly distributed random numbers U(0,1)

a = np.sqrt(-2*np.log(U1))
b = 2*np.pi*U2

Z = a*np.cos(b) #Standard Normal distributed numbers
samples_box_muller= Z*sigma+mu #Normal distribution with mean and sigma

Step 2: Plot the estimated histogram

Typically, if we have a vector of random numbers that is drawn from a distribution, we can estimate the PDF using the histogram tool. Matplotlib’s hist function can be used to compute and plot histograms. If the density argument is set to ‘True’, the hist function computes the normalized histogram such that the area under the histogram will sum to 1. Estimate and plot the normalized histogram using the hist function.

#For plotting
import matplotlib.pyplot as plt
%matplotlib inline
plt.style.use('ggplot')

fig, ax0 = plt.subplots(ncols=1, nrows=1) #creating plot axes
(values, bins, _) = ax0.hist(samples_normal,bins=100,density=True,label="Histogram of samples") #Compute and plot histogram, return the computed values and bins

Step 3: Theoretical PDF:

And for verification, overlay the theoretical PDF for the intended distribution. The theoretical PDF of normally distributed random samples is given by

Theoretical PDF for normal distribution is readily obtained from stats.norm.pdf() function in the SciPy package.

from scipy import stats
bin_centers = 0.5*(bins[1:] + bins[:-1])
pdf = stats.norm.pdf(x = bin_centers, loc=mu, scale=sigma) #Compute probability density function
ax0.plot(bin_centers, pdf, label="PDF",color='black') #Plot PDF
ax0.legend()#Legend entries
ax0.set_title('PDF of samples from numpy.random.normal()');

*Figure 1: Estimated PDF (histogram) and the theoretical PDF for samples generated using numpy.random.normal() function*

The histogram and theoretical PDF of random samples generated using Box-Muller transformation, can be plotted in a similar manner.

#Samples generated using Box-Muller transformation
from numpy.random import uniform
U1 = uniform(low=0,high=1,size=(L,1)) #uniformly distributed random numbers U(0,1)
U2 = uniform(low=0,high=1,size=(L,1)) #uniformly distributed random numbers U(0,1)
a = np.sqrt(-2*np.log(U1))
b = 2*np.pi*U2
Z = a*np.cos(b) #Standard Normal distribution
samples_box_muller = Z*sigma+mu #Normal distribution with mean and sigma

#Plotting
fig, ax1 = plt.subplots(ncols=1, nrows=1) #creating plot axes
(values, bins, _) = ax1.hist(samples_box_muller,bins=100,density=True,label="Histogram of samples") #Plot histogram
bin_centers = 0.5*(bins[1:] + bins[:-1])
pdf = stats.norm.pdf(x = bin_centers, loc=mu, scale=sigma) #Compute probability density function
ax1.plot(bin_centers, pdf, label="PDF",color='black') #Plot PDF
ax1.legend()#Legend entries
ax1.set_title('Box-Muller transformation');

*Figure 2: Estimated PDF (histogram) and theoretical PDF for samples generated using Box-Muller transformation for normal distribution*

References:

[1] John Mount, ‘Six Fundamental Methods to Generate a Random Variable’, January 20, 2012
[2] Thomas, D. B., Luk. W., Leong, P. H. W., and Villasenor, J. D. 2007. Gaussian random number generators. ACM Comput. Surv. 39, 4, Article 11 (October 2007), 38 pages DOI = 10.1145/1287620.1287622 http://doi.acm.org/10.1145/1287620.1287622

Rate this post: Note: There is a rating embedded within this post, please visit this post to rate it.

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Articles in this series
[1] Fibonacci series in python
[2] Central Limit Theorem – a demonstration
[3] Moving Average Filter in Python and Matlab
[4] How to plot FFT in Python – FFT of basic signals : Sine and Cosine waves
[5] How to plot audio files as time-series using Scipy python
[6] How to design a simple FIR filter to reject unwanted frequencies
[7] Analytic signal, Hilbert Transform and FFT
[8] Non-central Chi-squared Distribution
[9] Simulation of M-PSK modulation techniques in AWGN channel (in Matlab and Python)
[10] QPSK modulation and Demodulation (with Matlab and Python implementation)

Hidden Markov Models (HMM) – Simplified !!!

Markov chains are useful in computing the probability of events that are observable. However, in many real world applications, the events that we are interested in are usually hidden, that is we don’t observe them directly. These hidden events need to be inferred. For example, given a sentence in a natural language we only observe the words and characters directly. The parts-of-speech from the sentence are hidden, they have to be inferred. This bring us to the following topic– the hidden Markov models.

Hidden Markov models enables us to visualize both observations and the associated hidden events. Let’s consider an example for understanding the concept.

The cheating casino and the gullible gambler

Consider a dishonest casino that deceives it player by using two types of die : a fair die (F) and a loaded die (L). For a fair die, each of the faces has the same probability of landing facing up. For the loaded die, the probabilities of the faces are skewed as given next

When the gambler throws the die, numbers land facing up. These are our observations at a given time t (denoted as O_t = {1,2,3,4,5,6}). At any given time t, whether these number are rolled from a fair die (state S_t = F) or a loaded die (S_t = L), is unknown to an observer and therefore they are the hidden events.

Emission probabilities

The probabilities associated with the observations are the observation likelihoods, also called emission probabilities (B).

Initial probabilities

The initial probability of starting (at time = 0) with any of fair die or loaded die (hidden event) is 50%.

Transition probabilities

A gullible gambler switches from the fair die to loaded die with 10% probability. He switches back from loaded die to fair die with 5% probability.

The probabilities of transitioning from one hidden event to another is described by the transition probability matrix (A). The elements of the probability transition matrix, are the transition probabilities (p_ij) of moving from one hidden state i to another hidden state j.

The transition probabilities from time t-1 to t, for the hidden events are

Therefore, the transition probability matrix is

Based on the given information so far, a probability model is constructed. This is the Hidden Markov Model (HMM) for the given problem.

Figure 1: Hidden Markov Model for the cheating Casino problem

Assumptions

We saw, in previous article, that the Markov models come with assumptions. Similarly, HMMs models also have such assumptions.

1. Assumption on probability of hidden states

In the model given here, the probability of a given hidden state depends only on the previous hidden state. This is a typical first order Markov chain assumption.

2. Assumption on Output

The probability of any observation (output) depends on the hidden state that produce it and not on any other hidden state or output observations.

Problems and Algorithms

Let’s briefly discuss the different problems and the related algorithms for HMMs. The algorithms will be explained in detail in the future articles.

In the dishonest casino, the gambler rolls the following numbers:

1. Evaluation

Given the model of the dishonest casino, what is the probability of obtaining the above sequence ? This is a typical evaluation problem in HMMs. Forward algorithm is applied for such evaluation problems.

2. Decoding

What is the most likely sequence of die (hidden states) given the above sequence ? Such problems are addressed by Viterbi decoding.

What is the probability of fourth die being loaded, given the above sequence ? Forward-backward algorithm to our rescue.

3. Learning

Learning problems involve parametrization of the model. In learning problems, we attempt to find the various parameters (transition probabilities, emission probabilities) of the HMM, given the observation. Baum-Welch algorithm helps us to find the unknown parameters of a HMM.

Some real-life examples

Here are some real-life examples of HMM applications:

Speech recognition: HMMs are widely used in speech recognition systems to model the variability of speech sounds. In this application, the observable events are the acoustic features of the speech signal, while the hidden states represent the phonemes or words that generate the speech signal.
Handwriting recognition: HMMs can be used to recognize handwritten characters by modeling the temporal variability of the pen strokes. In this application, the observable events are the coordinates of the pen on the writing surface, while the hidden states represent the letters or symbols that generate the handwriting.
Stock price prediction: HMMs can be used to model the behavior of stock prices and predict future price movements. In this application, the observable events are the daily price movements, while the hidden states represent the different market conditions that generate the price movements.
Gene prediction: HMMs can be used to identify genes in DNA sequences. In this application, the observable events are the nucleotides in the DNA sequence, while the hidden states represent the different regions of the genome that generate the sequence.
Natural language processing: HMMs are used in many natural language processing tasks, such as part-of-speech tagging and named entity recognition. In these applications, the observable events are the words in the text, while the hidden states represent the grammatical structures or semantic categories that generate the text.
Image and video analysis: HMMs can be used to analyze images and videos, such as for object recognition and tracking. In this application, the observable events are the pixel values in the image or video, while the hidden states represent the object or motion that generates the pixel values.
Bio-signal analysis: HMMs can be used to analyze physiological signals, such as electroencephalograms (EEGs) and electrocardiograms (ECGs). In this application, the observable events are the signal measurements, while the hidden states represent the physiological states that generate the signal.
Radar signal processing: HMMs can be used to process radar signals and detect targets in noisy environments. In this application, the observable events are the radar measurements, while the hidden states represent the presence or absence of targets.

Rate this post: Note: There is a rating embedded within this post, please visit this post to rate it.

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Markov Chains – Simplified !!

Key focus: Markov chains are a probabilistic models that describe a sequence of observations whose occurrence are statistically dependent only on the previous ones.

● Time-series data like speech, stock price movements.
● Words in a sentence.
● Base pairs on the rung of a DNA ladder.

States and transitions

Assume that we want to model the behavior of a driver behind the wheel. The possible behaviors are

● accelerate (state 1)
● constant speed (state 2)
● idling (engine running slowly but the vehicle is not moving – (state 3))
● brake (state 4)

Let’s refer each of these behaviors as a state. In the given example, there are N=4 states, refer them as Q = {q₁,q₂,q₃,q₄}.

We observe the following pattern in the driver’s behavior (Figure 1). That is, the driver operates the vehicle through a certain sequence of states. In the graph shown in Figure 1, the states are represented as nodes and the transitions as edges.

Figure 1: Driver’s behavior – operating the vehicle through a sequence of states

We see that, sometimes, the driver changes the state of the vehicle from one state to another and sometimes stays in the same state (as indicated by the arrows).

We also note that either the vehicle stays in the same state or changes to the next state. Therefore, from this model, if we want to predict the future state, all that matters is the current state of the vehicle. The past states has no bearing on the future state except through the current state. Take note of this important assumption for now.

Probabilistic model

We know that we cannot be certain about the driver’s behavior at any given point in time. Therefore, to model this uncertainty, the model is turned into a probabilistic model. A probabilistic model allows us to account for the likelihood of the behaviors or change of states.

An example for a probabilistic model for the given problem is given in Figure 2.

Figure 2: Driver’s behavior – a probabilistic model (transition matrix shown)

In this probabilistic model, we have assigned probability values to the transitions.These probabilities are collectively called transition probabilities. For example, considering the state named “idling”, the probability of the car to transition from this state to the next state (accelerate) is 0.8. In probability mathematics this is expressed as a conditional probability conditioned on the previous state.

p(state = “accelerate” | previous state = “idling” ) = 0.8

Usually, the transition probabilities are formulated in the form of matrix called transition probability matrix. The transition probability matrix is shown in Figure 2. In a transition matrix, denoted as A, each element a_ij represent the probability of transitioning from state i to state j. The elements of the transition matrix satisfy the following property.

That is, the sum of transition probabilities leaving any given state is 1.

As we know, in this example, the driver cannot start car in any state (example, it is impossible to start the car in “constant speed” state). He can only start the car from at rest (i.e, brake state). To model this uncertainty, we introduce π_i – the probability that the Markov chain starts in a given state i. The set of starting probabilities for all the N states are called initial probability distribution (π = π₁, π₂, …, π_N). In Figure 3, the starting probabilities are denoted by green arrows.

Figure 3: Markov Chain model for driver’s behavior

Markov Assumption

As noted in the definition, the Markov chain in this example, assumes that the occurrence of each event/observation is statistically dependent only on the previous one. This is a first order Markov chain (or termed as bigram language model in natural language processing application). For the states Q = {q₁, …, q_n}, predicting the probability of a future state depends only on the current observation, all other previous observations do not matter. In probabilistic terms, this first order Markov chain assumption is denoted as

Extending the assumption for m^th order Markov chain, predicting the probability of a future observation depends only on the previous m observations. This is an m-gram model.

Given a set of n arbitrary random variables/observations Q = {q₁, …, q_n}, their joint probability distribution is usually computed by applying the following chain rule.

However, if the random observations in Q are of sequential in nature and follows the generic m^th order Markov chain model, then the computation of joint probability gets simplified.

The Markov assumptions for first and second order of Markov models are summarized in Figure 4.Figure 4: Assumptions for 1st order and 2nd order Markov chains

Hidden Markov Model (HMM)

Markov chains are useful in computing the probability of events that are observable. However, in many real world applications, the events that we are interested in are usually hidden, that is we don’t observe them directly. These hidden events need to be inferred. For example, given a sentence in a natural language we only observe the words and characters directly. The parts-of-speech from the sentence are hidden, they have to be inferred. This brings us to the next topic of discussion – the hidden Markov models.

Rate this post: Note: There is a rating embedded within this post, please visit this post to rate it.

Books by the author

Wireless Communication Systems in Matlab Second Edition(PDF) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Python (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart	Digital Modulations using Matlab (PDF ebook) Note: There is a rating embedded within this post, please visit this post to rate it. Checkout Added to cart
Hand-picked Best books on Communication Engineering Best books on Signal Processing

Cookie	Duration	Description
cookielawinfo-checbox-analytics	11 months	This cookie is set by GDPR Cookie Consent plugin. The cookie is used to store the user consent for the cookies in the category "Analytics".
cookielawinfo-checbox-analytics	11 months	This cookie is set by GDPR Cookie Consent plugin. The cookie is used to store the user consent for the cookies in the category "Analytics".
cookielawinfo-checbox-functional	11 months	The cookie is set by GDPR cookie consent to record the user consent for the cookies in the category "Functional".
cookielawinfo-checbox-functional	11 months	The cookie is set by GDPR cookie consent to record the user consent for the cookies in the category "Functional".
cookielawinfo-checbox-others	11 months	This cookie is set by GDPR Cookie Consent plugin. The cookie is used to store the user consent for the cookies in the category "Other.
cookielawinfo-checbox-others	11 months	This cookie is set by GDPR Cookie Consent plugin. The cookie is used to store the user consent for the cookies in the category "Other.
cookielawinfo-checkbox-necessary	11 months	This cookie is set by GDPR Cookie Consent plugin. The cookies is used to store the user consent for the cookies in the category "Necessary".
cookielawinfo-checkbox-performance	11 months	This cookie is set by GDPR Cookie Consent plugin. The cookie is used to store the user consent for the cookies in the category "Performance".
viewed_cookie_policy	11 months	The cookie is set by the GDPR Cookie Consent plugin and is used to store whether or not user has consented to the use of cookies. It does not store any personal data.

Spectral leakage

Equivalent Noise Bandwidth:[2]

Coherent power gain:

Scalloping loss

Figure of merits for window functions

List of window functions

Boxcar (Rectangular)

Barthann

Bartlett

Blackman

Blackman-Harris

Bohman

Cosine

Flattop

Hamming

Hann

Nuttall

Parzen

Triangular

References:

Similar articles

Additional Resources:

Books by the author

FFT and spectral leakage

Equivalent noise bandwidth

Python script

Reference

Similar articles

Books by the author

Deriving Parseval’s theorem

Reference

Books by the author

Simulation and performance results

References

Books by the author

Deterministic and stochastic signals

Summary of descriptive statistical measures

Quantitative measures of shape:

Measure of central tendency – mean – the first moment:

Measure of dispersion – variance – the second moment:

Higher order moments – skewness and kurtosis:

Measures of association

Scatter plot

Correlation

Cross-correlation

Auto-correlation

Books by the author

For further reading

Why complex baseband equivalent model

BPSK constellation

Simulation methodology

Simulation using Python

Simulation using Matlab

Reference

Books by the author

Plotting constellation diagram

Conclusions

References

In this chapter

Books by the author

Step 1: Generate random samples

Step 2: Plot the estimated histogram

Step 3: Theoretical PDF:

References:

Books by the author

Similar topics

The cheating casino and the gullible gambler

Emission probabilities

Initial probabilities

Transition probabilities

Assumptions

1. Assumption on probability of hidden states

2. Assumption on Output

Problems and Algorithms

1. Evaluation

2. Decoding

3. Learning

Some real-life examples

Books by the author

Similar topics

Equivalent Noise Bandwidth:^[2]