In this post we take a first look at the spectral parameters of higher-order cyclostationarity (HOCS). In previous posts, I have introduced the topic of HOCS and have looked at the temporal parameters, such as cyclic cumulants and cyclic moments. Those temporal parameters have proven useful in modulation classification and parameter estimation settings, and will likely be an important part of my ultimate radio-frequency scene analyzer.
The spectral parameters of HOCS have not proven to be as useful as the temporal parameters unless you include the trivial case where the moment/cumulant order is equal to two. In that case, the spectral parameters reduce to the spectral correlation function, which is extremely useful in CSP (see the TDOA and signal-detection posts for examples).
When we started thinking about generalizing second-order cyclostationarity (SOCS) to HOCS in the time domain, it seemed natural to consider higher-order moment functions. This followed from thinking of the time-varying autocorrelation function as a second-order moment. So, what can be learned from third-order time-varying moments, or fourth-order, or th-order? That led us to cumulants, and in particular cyclic cumulants, which have the same highly desirable signal-selectivity properties as does the spectral correlation function. Higher-order temporal moments, on the other hand, do not always possess the properties we seek.
Now we wish to look at HOCS in the frequency domain. We can parallel our development in the time domain by recalling that the spectral correlation function is a second-order frequency-domain moment. Any guesses about which class of functions–spectral moments or spectral cumulants–will be most useful? Most well-behaved?
Second-Order Spectral Moments and Cumulants
Let’s start off by reviewing what we did for second order: the spectral correlation function. This is a kind of second-order moment involving a frequency-domain representation of two components of the signal. Recall we first defined the second-order product of interest, the cyclic periodogram, which is the product of two downconverted narrowband components of the signal under study ,
as the amount of processed data increases without bound, and then the spectral resolution () is allowed to decrease to zero,
In these equations, I define the narrowband spectral component in the following way
Higher-Order Spectral Moments
Let’s simply look at a higher-order product of spectral components,
where is the approximate bandwidth of the spectral components, and
denotes an infinite time-average operation. Let me try to justify the appearance of the optional conjugations
and corresponding optional negations
. Consider some arbitrary collections of time-series (signals)
. An
-th order spectral moment for these
signals is simply
Now if we make the particular choice , we end up with (4). So the particular conjugation configuration determines which factors in the moment have conjugations and corresponding negated frequencies.
It isn’t too hard to show that the spectral moment function (4) is nonzero only for a certain subset of the frequencies
. The frequencies must sum to an
th-order cycle frequency (corresponding to
optional conjugations) for
,
That is, the sum of the frequencies has to equal a cycle frequency for the set
.
What we name the actual -th order spectral moment function (SMF) is the function (4) as the bandwidth of the narrowband components approaches zero,
or
,
Restricting our attention to those subsets of -dimension frequency space for which the SMF is potentially non-zero, we see from (6) that the
frequencies are not independent. For a given cycle frequency
, the
frequencies can be written as
independently chosen frequencies followed by the
th frequency which is determined by the others and
,
For convenience, let’s define the vector as the first
elements of
,
Then the reduced-dimension spectral moment function (RD-SMF) is defined by
It turns out that this function is the -dimensional Fourier transform of the reduced-dimension cyclic temporal moment function (RD-CTMF),
although in many cases it contains impulses (Dirac delta functions) or even products of impulses, so it exists only in a formal mathematical sense. This happens for order as well. Consider a signal that contains an additive sine-wave component, such as
where the signal has a continuous spectrum (no additive sine-wave components) and
is AWGN. Then we know that the spectrum (the ‘second-order spectral moment’ using our terminology in this post) must contain an impulse due to the sine wave with frequency
. So we’re used to impulses in our second-order spectral moments. For the higher orders, though, there are generally more impulses and there can be products of impulses in the spectral moment function.
The reduced-dimension spectral moment and the spectral moment are related by an impulsive function,
So the spectral moment is zero unless its frequencies sum to an
th-order cycle frequency (with
conjugations). And in that case, it is impulsive.
Recalling our discussions of the higher-order temporal characterizations of cyclostationary signals, we see that the spectral moment (and its reduced-dimension version) don’t help us clearly answer questions about the value of higher-order measures relative to second order. We have to go to cumulants for that.
Higher-Order Spectral Cumulants
Since the spectral moment was obtained using the infinite-time average of the finite-resolution spectral components (4) and then letting the spectral resolution approach zero in (7), we start our spectral cumulant investigation by looking at the cumulant of the
narrowband components,
Then we let the spectral resolution , which means that the narrowband spectral components are estimated with an increasingly long data segment (length
), and obtain the spectral cumulant function
This spectral cumulant is the sum of products of lower-order (spectral) moments, just like all cumulants. But each lower-order spectral moment is either zero or impulsive. It is impulsive if the subset of corresponding to its partition itself sums to a lower-order cycle frequency, and is zero otherwise. Since the
th-order cumulant can always be written as the
th-order moment plus several weighted products of lower-order moments, this means that the spectral moment function will equal the spectral cumulant function, and will be nonzero, if the sum of the
frequencies is an
th-order cycle frequency and none of the partition-induced subsets of
themselves sum to the appropriate-order cycle frequencies.
I’m avoiding writing all that down because the mathematical notation gets really complex. You can see it all, though, in My Papers [5,6] and my downloadable dissertation (item 2 on the Downloads page).
The subsets of -dimensional space for which the sum-of-lower-order spectral moments is not zero are called the
-submanifolds. When we pick
th-order frequency vectors
that do not lie on the
-submanifolds, then we can estimate the spectral cumulant by estimating the spectral moment.
When we reduce the dimension of the spectral cumulant, as we did for the spectral moment, we obtain a similar characterization,
where is the
th-order cyclic polyspectrum. It can be shown that the cyclic polyspectrum is the
dimensional Fourier transform of the
th-order reduced-dimension cyclic temporal cumulant function,
If the reduced-dimension cyclic cumulant is integrable, then the Fourier transform exists and has no impulsive components. This, in turn, will be the case if the cumulant decays in all directions in -dimensional space. All that means is that the underlying signal
possesses a correlation property such that widely separated samples become uncorrelated as the separation gets large, a property that almost all signals of practical interest possess.
Putting it All Together
We’ve looked at the second-order temporal parameters of cyclostationary signals, the second-order spectral parameters, the higher-order temporal parameters, and finally in this post, we’ve looked at the higher-order spectral parameters. The relationships between the various parameters are summarized in the following figure:

The set of mathematical objects we’ve defined and evaluated encompass many familiar objects from various parts of signal processing (and beyond), including the autocorrelation, power spectrum, cyclic autocorrelation, spectral correlation function, conjugate cyclic autocorrelation, conjugate spectral correlation function, conventional polyspectrum, conventional temporal cumulant, cyclic cumulant, and cyclic polyspectrum.
There are other objects of great interest to us when we apply these basic tools, such as the spectral coherence function, which do not appear in the golden diamond figure above.
We’ll discuss some of the connections to more commonly known mathematical objects in a future post, as a means of potentially harnessing previously gained intuition from practical signal-processing experience and formal study in stochastic-process theory.
Estimators and the
-Submanifold Problem
There are at least three ways to estimate the cyclic polyspectrum (My Papers [6]).
The first is to estimate the reduced-dimension cyclic cumulant on a hypercube (many values of ), then perform an
-dimensional Fourier transform.
The second is to estimate the spectral moment function, being careful to avoid -submanifolds, with a generalization of the frequency-smoothing method of spectral correlation estimation. This also requires estimation on a hypercube of dimension
only this time the hypercube corresponds to
.
The third is a generalization of the time-smoothing method of spectral correlation estimation. This method simply averages, over time, products of narrowband spectral components, being careful to avoid -submanifolds. It can be used for arbitrary frequency vectors
. That is, to obtain an estimate of the cyclic polyspectrum on some one-dimensional subspace of
, one need not estimate anything on any hypercubes.
So the time-smoothing method is the only one of the three that has a low computational cost because it can be directly used to estimate the cyclic polyspectrum on a low-dimensional slice through its high-dimensional domain.
The Cyclic Polyspectrum for Digital QAM/PSK
The cyclic cumulant function and the cyclic polyspectrum are known (My Papers [5,6]) for the case of textbook digital QAM/PSK. Such signals have the model
The cycle frequencies for such signals take the form
So we can speak of the “ cyclic cumulant” or the “
cyclic polyspectrum.” The cyclic cumulant and cyclic polyspectrum for a particular set of
may be zero, however, depending on the nature of the constellation embodied by the sequence
, and on the nature of the pulse function
. The values of the symbol cumulant function are tabulated in this post for several modulation types.
The theoretical cyclic polyspectrum for is given by
The function is the cumulant for the set of variables
, and
is the Fourier transform of
.
The post on signal processing operations can be used to find the generalization of (19) for the case where is not zero.
Numerical Examples
We use the generalization of the time-smoothing method for our cyclic polyspectrum estimator and evaluate (19) for the theoretical function.
The signal is QPSK, with pulse function given by the square-root raised-cosine pulse with roll-off (excess bandwidth of
%). The signal uses
and
.
Let’s start with . We have an excellent comparison between measured and theoretical cyclic polyspectra:

That second-order result helps to validate the estimator. Notice that the features for
are shifted away from the carrier, which has value of zero. This is because the CP for order two is a frequency shifted version of the spectral correlation function.
Now let’s turn to the various fourth-order cyclic polyspectra. Each is estimated for a range of and for specific values, as stated below, for the remaining frequencies
and
. The fourth frequency is determined by
and
and the cycle frequency
.
The cyclic polyspectrum for
and
is shown here:

The zero near is due to a
submanifold. The cyclic polyspectrum is not estimated for this frequency. One could optionally interpolate between neighboring points to fill this in.
The cyclic polyspectrum for
and
and the
cyclic polyspectrum for
and
is shown next:

And finally the cyclic polyspectrum for
and
is shown here:

Final Thoughts
The cyclic polyspectrum is difficult to interpret and to compute due to the presence of the -submanifolds, which become more numerous as the number of lower-order cycle frequencies increases (think rectangular pulse
or DSSS) and as the order increases.
The application of the cyclic polyspectrum that is most appealing is channel estimation. The equation (19) is the key. The effect of a linear time-invariant propagation channel leads to an new effective pulse
Taking the logarithm of both sides transforms the product of pulse transforms involving into a sum of those transforms. Choosing to evaluate the resulting equation for multiple choices for
leads to a linear system of equations in the logarithm of the channel transfer function, which might be invertible. I’ve made some progress on that front and hope to reveal it in a future post.
As always, let me know if you find any errors in the post or leave a comment if you have experience with the cyclic polyspectrum. I’d like to hear from you!