Technical Report: Compressive Temporal Higher Order Cyclostationary Statistics

T ec hnical Rep ort: Compressiv e T emp oral Higher Order Cyclostationary Statistics Chia W ei Lim and Mic hael B. W akin ∗ June 2, 2021 Abstract The application of nonlinear transformations to a cyclostationary signal for the purp ose of rev ealing hidden perio dicities has pro ven to be useful for applications requiring signal selectiv- it y and noise tolerance. The fact that the hidden p eriodicities, referred to as cyclic momen ts, are often compressible in the F ourier domain motiv ates the use of compressive sensing (CS) as an eﬃcien t acquisition proto col for capturing such signals. In this w ork, we consider the class of T emp oral Higher Order Cyclostationary Statistics (THOCS) estimators when CS is used to acquire the cyclostationary signal assuming compressible cyclic momen ts in the F ourier domain. W e dev elop a theoretical framework for estimating THOCS using the low-rate non uniform sam- pling proto col from CS and illustrate the p erformance of this framework using simulated data. 1 In tro duction 1.1 T emp oral Higher Order (Cyclostationary) Statistics A cyclostationary signal (as deﬁned in [8]) is one whic h, after undergoing nonlinear transformations, will con tain certain perio dicities. Comm unication signals are not, in general, p erio dic. Ho wev er, man y comm unication signals are cyclostationary . F or example, it is kno wn [17] that b y taking a Binary Phase-Shift-Keying (BPSK) signal and multiplying it by itself, one can generate p eriodicities in the signal at twice its carrier frequency . Figure 1a shows the sp ectrum of a BPSK signal while Figure 1b sho ws the sp ectrum of the squared signal. Likewise, if one raises a BPSK signal to the fourth p o w er, p eriodicities at four time s its carrier frequency can be generated. These “hidden” p eriodicities, referred to as higher or der statistics (HOS) [17], are generated when suc h nonlinear transformations (namely , raising the signal to some p o wer) are applied to the signal. The pioneering w ork of Gardner et al. [10, 20], in formulating the theory of temp or al higher or der cyclostationary statistics (THOCS), refers to the hidden perio dicities as the cyclic moments of the signal. Giv en the cyclic moments of the signal, one can compute its corresp onding cyclic cumulants whic h ha ve b een shown to b e signal selectiv e and robust to Gaussian noise. The theory of THOCS subsumes that of HOS and deﬁnes additional statistical functions such as cyclic cumulan ts. 1.2 Compressiv e Sensing (CS) The theory of c ompr essive sensing (CS) has in recent years oﬀered the great promise of eﬃciently capturing essen tial signal information with low-rate sampling proto cols, often below the minimum ∗ Departmen t of Electrical Engineering and Computer Science, Colorado School of Mines. E-mail: clim,m wakin@mines.edu. This w ork w as partially supported b y DSO National Laboratories of Singapore. A prelim- inary v ersion of this w ork w as ﬁrst presen ted at the 2012 SPIE Defense, Security , and Sensing Conference [13]. 1 (a) BPSK signal (b) squared BPSK signal Figure 1: Periodicities (p eaks) are hidden in the sp ectrum of a BPSK signal but rev ealed in the spectrum of the squared signal. rate required by the Nyquist sampling theorem. The use of CS can dramatically reduce b oth the complexit y of a s ignal receiv er and the amount of data that must be stored and/or transmitted do wnstream [1, 7, 24]. CS exploits the fact that some high-bandwidth signals are c ompr essible in that they can b e well describ ed by a small num b er of coeﬃcients when expressed in an appropriate domain. Early work in CS [1] established the feasibility of collecting a random low-rate stream of nonuniform signal samples (as opp osed to a deterministic high-rate stream of uniform signal samples) to acquire a signal that has relativ ely few large coeﬃcients in its F ourier-domain spectrum. Giv en a suﬃcient n umber of nonuniform samples, subsequent reconstruction of the signal is p ossible by solving a con vex optimization problem. What is remark able ab out CS is that the requisite num b er of samples dep ends on the compressibilit y of the signal’s sp ectrum instead of the width of its s pectrum (the latter is what is required b y the Nyquist sampling theorem). 1.3 Challenges of THOCS and the Promise of CS While THOCS hav e prov en to b e useful, there are key c hallenges that one faces when attempting to estimate THOCS from traditional uniform samples of a cyclostationary signal. The higher the order of THOCS (that is, the higher the order of the nonlinear transformations applied to the signal), the more one m ust ov ersample the signal relative to its Nyquist rate to preven t aliasing. T o accurately estimate THOCS, one also needs to observ e the signal for an adequate (often long) time in terv al. These restrictions translate to hea vy sampling burdens and data storage requirements. As noted in Section 1.1, a typical cyclostationary signal is not compressible in the frequency domain. Ho wev er, the signal ma y b e compressible after undergoing a nonlinear transformation, and the dominant p eaks in the F ourier sp ectrum (see, for example, Figure 1b) are in fact the THOCS that one migh t wish to estimate. This motiv ates the question of whether CS in general, and low-rate non uniform sampling in particular, can somehow b e used to simplify the data acquisition pro cess in problems in volving THOCS. The immediate con undrum that one faces is the follo wing: since the signal coming in to the receiver is not compressible, one cannot simply apply traditional CS sampling and reconstruction proto cols to that signal (and then, say , estimate THOCS using conv en tional metho ds). While one could consider applying nonlinear transformations to the incoming signal and then collecting nonuniform samples of the transformed signal, this could signiﬁcan tly increase the 2 complexit y of a receiver. In this pap er we consider the problem of estimating THOCS directly from low-rate nonuniform samples of the original cyclostationary signal as well as low-rate nonuniform samples of delay ed v ersions of the signal. In short, rather than collecting nonuniform samples of a nonlinearly trans- formed signal, one can simply apply these same nonlinear transformations to nonuniform samples of the original cyclostationary signal and non uniform samples of the dela yed versions of it. With these transformed non uniform samples one could consider applying traditional CS reconstruction proto cols and then extracting THOCS features. How ev er, these reconstruction algorithms could result in additional pro cessing latency and require signiﬁcan t computational resources. Therefore, w e b orro w ideas from the ﬁeld of compressive signal processing (CSP) [4] and prop ose a very sim- ple alternativ e metho d for estimating THOCS directly from the transformed nonuniform samples without the need for full-scale signal reconstruction. In summary , w e prop ose a traditional CS receiver fron t end that collects random lo w-rate streams of non uniform samples of the incoming signal and delay ed v ersions of it. While the com- plexit y of such a receiver (see, e.g., [24] for a simpliﬁed version) is increased due to the additional signal dela y paths, it has the potential to dramatically reduce the data acquisition, storage, and transmission burdens compared to conv en tional THOCS estimators that are based on high-rate uniform sampling. Our back end process ing, in whic h w e estimate THOCS from the non uniform samples, diﬀers from traditional CS and is computationally very simple. W e refer to the resulting estimated THOCS as c ompr essive temp or al higher or der cyclostationary statistics (CTHOCS). This pap er describ es the acquisition/estimation framew ork underlying CTHOCS and analyzes its eﬀec- tiv eness theoretically and empirically in comparison with conv en tional THOCS based on high-rate uniform sampling. Our paper is organized as follo ws. In Section 2, we ﬁrst provide brief reviews of the theories of THOCS and lo w-rate nonuniform sampling. In Section 3, w e discuss how THOCS can b e estimated from lo w-rate non uniform samples and subsequently ev aluate the p erformance of such a CTHOCS estimator. W e demonstrate the signal selectivity of CTHOCS in Section 4. W e conclude with a ﬁnal discussion in Section 5. 2 Preliminaries While the conten t in this section is intended to mak e this pap er as self-contained as p ossible, it is not complete. F or completeness, the reader may refer to [1, 7, 10, 12, 19, 20, 24]. 2.1 THOCS Preliminaries 2.1.1 F raction-of-Time Probability (FOT) Mathematical analysis of HOCS is based on the fr action-of-time (F OT) probability framew ork [9, 10, 12], where statistical parameters or functions are deﬁned based on a single realization of an inﬁnite duration signal. This is in con trast to the stochastic pro cess framew ork where statistical parameters or functions are deﬁned based on the ensemble a verage of m ultiple realizations of the sto c hastic signal. (F or ergo dic signals, statistical functions coincide b et ween the t wo framew orks.) The F OT probability framework is useful in applications where only a single realization of the signal is a v ailable. As we will discuss, it facilitates the deﬁnitions of exp ectation, momen ts and cum ulants which are important for the study of cyclostationary time series. 3 F or a given time-series f ( u ), the m ultiple sine-wa v e extraction op erator 1 is deﬁned [10, 20] as ˆ E [ f ( u )] ( t ) , X α ∈ A ( f ( u )) h f ( u ) e − j 2 π αu i u e j 2 π αt where h·i u denotes the time av eraging op eration with resp ect to the v ariable u and the set A ( f ( u )) con tains all frequencies α ∈ R for which h f ( u ) e − j 2 π αu i u 6 = 0; this set is assumed to b e countable for the time-series w e will deal with. (Equiv alently , A ( f ( u )) con tains the set of frequencies α at whic h the F ourier transform of f ( u ) contains a Dirac impulse.) The op erator ˆ E [ · ] is said to extract all additive sine-wa ve comp onen ts from f . It can b e viewed as ﬁrst analyzing a signal (in the F ourier domain) for all p erio dic comp onen ts and then rebuilding the signal only from these same perio dic comp onen ts. 2 When f ( u ) is a time series deﬁned in the F OT probabilit y framew ork, the ˆ E [ · ] operator extracts the deterministic component of its argumen t and hence is analogous to the exp ectation operator E [ · ] in the sto c hastic pro cess framew ork. 2.1.2 THOCS Statistical F unctions W e no w provide the deﬁnitions of v arious THOCS statistical functions based on the F OT analysis framew ork. The THOCS statistical functions, collectively referred to as moments and cumulants , are analogous to the momen ts and cum ulants deﬁned in the sto chastic analysis framework. W e note here that while there exists another set of HOCS in the sp ectral domain [10, 20], THOCS shall b e the only set of HOCS considered in this framew ork. 1. Time-v arying n th order moments The time-varying n th or der moments of a signal x ( t ) are deﬁned [10, 20] as R x ( t, τ ) n,q , ˆ E " n Y i =1 x ( ∗ ) i ( u − τ i ) # ( t ) , (1) where the i th factor x ( u − τ i ) 3 could b e (optionally) conjugated and this conjugation is denoted b y ( ∗ ) i . The v ariable q is the total num b er of conjugations and τ , [ τ 1 . . . τ n ] T is a v ector of lags. 4 The argumen t to ˆ E [ · ] in (1) is termed the lag pr o duct of the signal x ( t ). In other w ords, giv en a signal x ( t ), its time-v arying n th order moments can b e obtained by ﬁrst analyzing its lag pro duct for all p erio dic comp onen ts and then rebuilding its lag pro duct from only these comp onen ts. 2. n th order cyclic moments The time-v arying n th order moments of a signal x ( t ) can be further expressed in terms of its lag-dep endent F ourier co eﬃcien ts (also kno wn as its n th or der cyclic moments ) as R x ( t, τ ) n,q = P α R α x ( τ ) n,q e j 2 π αt , where the n th order cyclic moments satisfy R α x ( τ ) n,q = h R x ( t, τ ) n,q e − j 2 π αt i t . (2) As mentioned in [10], due to R α x ( τ ) n,q b eing sinusoidal join tly in the n translation v ariables τ , setting τ 1 = 0 retains all the information present in R α x ( τ ) n,q . In general, R α x ( τ ) n,q is not in tegrable with resp ect to τ . 1 This op erator is denoted by ˆ E { α } in [10, 20]. 2 In general, the signal ma y contain non p eriodic comp onen ts as well. 3 The use of a minus sign here is to ensure consistency with regards to delaying the signal. 4 F ollowing [10, 20], (1) do es not indicate which q terms get conjugated, but throughout our pap er, w e will alwa ys assume the ﬁrst q terms in the pro duct are conjugated. 4 D n { 1 }{ 2 }{ 3 } { 1 2 }{ 3 } { 1 3 }{ 2 } { 2 3 }{ 1 } { 1 2 3 } d 3 2 2 2 1 v i v 1 = { 1 } v 1 = { 1 , 2 } v 1 = { 1 , 3 } v 1 = { 2 , 3 } v 1 = { 1 , 2 , 3 } v 2 = { 2 } v 2 = { 3 } v 2 = { 2 } v 2 = { 1 } v 3 = { 3 } T able 1: Partitions for n = 3 . F or a giv en n , all distinct partitions D n can be generated using a recursive form ula [19]. 3. Time-v arying n th order cumulan ts W e ﬁrst deﬁne some key parameters that the time-v arying n th order cumulan ts are dep enden t on. Distinct partitions of the index set { 1 , 2 , . . . , n } are referred to as D n , d is the num b er of elements in a partition, and the set of indices belonging to a partition is denoted by v i , i ∈ { 1 , 2 , . . . , d } . T able 1 sho ws the key parameters with their resp ectiv e v alues when n = 3. Using the moment-to-cum ulan t (M-C) con version formula 5 [10], the time-varying n th or der cumulants of a signal x ( t ) are given by C x ( t, τ ) n,q = X D n " ( − 1) d − 1 ( d − 1)! d Y i =1 R x ( t, τ v i ) n i ,q i # , (3) where n i = | v i | , q i is the num b er of conjugations used to compute the time-v arying momen ts deﬁned b elo w, R x ( t, τ v i ) n i ,q i , ˆ E   Y k ∈ v i x ( ∗ ) k ( u − τ k )   ( t ) , (4) τ v i is the vector of lags in τ with indices in v i , and ( ∗ ) k again denotes optional conjugation of the k th factor x ( u − τ k ). 4. n th order cyclic cumulan ts The n th or der cyclic cumulants are the F ourier co eﬃcien ts of the time-v arying n th order cum ulants, C β x ( τ ) n,q , h C x ( t, τ ) n,q e − j 2 π β t i t . It can b e further sho wn [19] that C β x ( τ ) n,q = X D n   ( − 1) d − 1 ( d − 1)! X α T 1 = β d Y i =1 R α i x ( τ v i ) n i ,q i   , (5) The inner sum in (5) is ov er all p ossible vectors α = [ α 1 . . . α d ] T whose entries sum to β and satisfy R α i x ( τ v i ) n i ,q i 6 = 0 for i = 1 , 2 , . . . , d . (The v alue of d changes depending on the partition of D n in the outer sum.) Thus, given the cyclic momen ts of a signal, one can use (5) to compute its corresp onding cyclic cumulan ts. In practice, for a giv en n i , q i , and τ v i , one can determine the candidate frequencies for α i b y iden tifying peaks in the sp ectrum of the lag pro duct app earing in (4). 5 There also exists a cumulan t-to-moment conv ersion formula [10]. 5 Giv en that C β x ( τ ) n,q is also sinusoidal join tly in the n translation v ariables τ , setting τ 1 = 0 again retains all the information present in C β x ( τ ) n,q . In general, C α x ( τ ) n,q is in tegrable with resp ect to τ for a large class of cyclostationary signals of practical interest [20]. 2.2 CS Preliminaries 2.2.1 Sparsit y and Compressive Measurements CS builds on the framework of sparse (and compressible) signal represen tations in order to simplify signal acquisition. Supp ose w e ha ve a ﬁnite length- L vector y whose en tries are samples of an analog signal acquired at or ab ov e the Nyquist rate. If w e represent y in terms of an L × L dictionary (or basis) Φ, where y = Φ η , y is said to hav e an s -sparse representation in Φ if the length L -co eﬃcien t v ector η has just s  L nonzero entries. If η has s signiﬁcan t co eﬃcien ts (but the rest are v ery small), y is said to b e compressible in the dictionary Φ. Cand ` es et al. [1] show ed that if one collects just P < L compressiv e measurements of a sparse signal y via a P × L sensing matrix R suc h that w = Ry = R Φ η =: Aη , exact recov ery of η (and therefore, y ) from w ma y b e p ossible. F or example, it is p ossible to solve the under-determined system of equations w = Aη (ha ving P equations and L unkno wns), if the P × L matrix A satisﬁes a condition known as the r estricte d isometry pr op erty (RIP). A matrix A satisﬁes the RIP of order s with isometry constant δ s ∈ (0 , 1) if the follo wing holds for all s -sparse v ectors η : (1 − δ s ) k η k 2 ≤ k Aη k 2 ≤ (1 + δ s ) k η k 2 . Exact recov ery of η from w can b e ac hieved b y solving a con vex optimization problem (namely , ` 1 -minimization [2]) if A satisﬁes the RIP of order 2 s and δ 2 s is small. A v ariety of greedy iterative algorithms (e.g., CoSaMP [16]) ha ve also b een prop osed. Moreo ver, these recov ery algorithms are robust to noise and provide an appro ximate solution if η is not exactly sparse but compressible. All of the algorithms mentioned ab o ve attempt to ﬁnd, among all η that satisfy the constraint w = Aη , the sparsest p ossible solution. This reco very problem w ould b e straightforw ard if A w ere a square matrix and its columns were orthogonal—one could uniquely identify η simply by computing A H w , where ( · ) H denotes the conjugate transp ose op erator. In a CS problem, A is not square (it is underdetermined) and so its columns cannot b e orthogonal. How ever, the RIP prop ert y requires that small sets of columns of A b e roughly orthogonal, and so A H w will provide a rough appro ximation to the true v ector η . Most CS reconstruction algorithms—particularly greedy algorithms—can b e in terpreted as correcting the estimate A H w to accoun t for the minor correlations among the columns of A . Ho wev er, it is also w orth noting that CSP theory [4] sup p orts the fact t hat the uncorrected estimate A H w can itself provide a reasonable approximation of η when η is sparse or compressible, and using this e stimate can b e computationally attractive b ecause it is m uch simpler to compute than a full-scale (but more accurate) signal reconstruction. 2.2.2 Lo w-rate Non uniform Sampling F or a signal that is sparse or compressible in the frequency domain—that is, one whose F ourier sp ectrum contains relatively few large p eaks—Cand ` es et al. [1] prop ose a CS acquisition strategy based on lo w-rate nonuniform sampling. This strategy in volv es collecting a random lo w-rate stream of nonuniform samples of an incoming analog signal, as opp osed to the deterministic high-rate stream of uniform signal samples acquired b y a conv entional Nyquist-rate receiver [24]. In some implemen tations of lo w-rate non uniform sampling, the random time in terv als are integer multiples of a ﬁner sampling in terv al (at least as ﬁne as Nyquist rate of the signal); we assume this is the case in our subsequent discussions. In practice, one can generate these random time interv als via a pseudorandom n umber generator. 6 Lo w-rate nonuniform sampling is compatible with the form ulation and notation deﬁned ab o ve. Recall the length- L vector y whose entries are uniform samples of an analog signal acquired at or ab o ve the Nyquist rate (denote the sampling interv al as T s ). A collection of P non uniform samples of the same analog signal can b e represen ted as w = Ry , where R is a P × L binary matrix con taining a single 1 in a random position in eac h row. F or example, if P = 4 and L = 10, one p ossible R matrix would b e R =    1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1    . In such an instance, w e w ould hav e w = [ w (0) w (1) w (2) w (3)] T = [ y (0) y (3 T s ) y (4 T s ) y (9 T s )] T . It w as sho wn in [18] that if Φ is the L × L DFT matrix and R is a random selection matrix as shown ab o v e, then A = R Φ will satisfy the RIP of order P / log 4 ( L ) with high probabilit y . (T ec hnically , the RIP holds for the appropriately scaled matrix A = p L/P · R Φ.) Hence, one can recov er an s -sparse co eﬃcient v ector η and consequently the Nyquist-rate sample v ector y from only a small m ultiple times s nonuniform samples. 3 CTHOCS F ramew ork In this section, w e explain our CTHOCS framework, in whic h THOCS are estimated directly from lo w-rate nonuniform samples of the original cyclostationary signal as w ell as low-rate nonuniform samples of dela yed versions of the signal. 3.1 Signal Acquisition In a t ypical scenario one wishes to estimate THOCS for v arious orders n , optional conjugations q , and delay v ectors τ ; recall from (1) that each combination of n , q , and τ deﬁnes a unique lag pro duct of the incoming signal x ( t ). The question is how to estimate these statistics without excessiv e sampling or sp ecialized hardw are. The idea of using low-rate non uniform sampling for estimating THOCS is inspired b y the compressibility of the lag pro duct. F or the sake of discussion, Figure 2a sho ws one hypothetical wa y that low-rate nonuniform sampling could b e used assuming sp ecialized analog hardware was av ailable to generate v arious (compressible) signal lag pro ducts. This system consists of n t c hannels, where n t is the total n umber of unique lag pro ducts of in terest. In each channel j ∈ { 1 , 2 , . . . , n t } , an analog hardware device H/W j computes a lag pro duct of x ( t ) corresp onding to some desired v alues of n , q , and τ , all of which ma y v ary from c hannel to c hannel. This lag pro duct is then passed through a lo w-rate non uniform sample-and-hold (s/h) device (to b e formally deﬁned in Section 3.4), and the output w lag ,j is a non uniform sample stream of the j th lag pro duct. Unfortunately , such an acquisition setup w ould require a signiﬁcant amoun t of sp ecialized hardware and its complexity w ould scale linearly with the desired n umber n t of unique lag pro ducts. As an alternative, we prop ose the signal acquisition setup sho wn in Figure 2b. F or the sake of simplicit y , this ﬁgure presents the acquisition sc heme for a single lag pro duct as app earing in (1), i.e., for a single com bination of n , q , and τ . In eac h channel i ∈ { 1 , 2 , . . . , n } , an analog delay pro duces x ( t − τ i ), where τ i is the i th en try of the dela y vector τ . This dela yed signal is then passed through a low-rate non uniform s/h device, whose output w ch i is a nonuniform sample stream of the i th delay ed input signal. In this sc heme w e assume that all s/h devices are driv en b y a common 7 x ( t ) w lag , 1 w lag , 2 w lag ,n t H/W 1 H/W 2 H/W n t s/h s/h s/h (a) Sp ecialized hardware sampling sc heme. x ( t ) w ch 1 w ch 2 w ch 3 w ch n NU CLK ( · ) ∗ ( · ) ∗ ( · ) ∗ ( · ) ∗ × w lag τ 1 τ 2 τ 3 s/h s/h s/h s/h τ n (b) Prop osed CS m ulti-channel sampling scheme. Figure 2: Comparison of t wo sampling sc hemes. (a) Lag products computed with sp ecialized analog hard- w are. (b) Prop osed scheme using analog signal dela ys and low-rate non uniform samplers. non uniform clo c k (NU CLK) which dictates the sample times. The output w lag of the acquisition system is a digital pro duct of the non uniform sample streams from each channel: w lag [ k ] = n Y i =1 w ( ∗ ) i ch i [ k ] , (6) where the optional conjugations matc h those in (1). This multi-c hannel acquisition setup can b e easily implemen ted using only analog signal dela ys and low-rate nonuniform samplers. At ﬁrst sigh t, the prop osed sampling sc heme may seem at o dds with our ob jectiv e of reducing the o verall sampling rate necessary for THOCS estimation due to the use of n channels, eac h with its o wn s/h device. How ev er, we note that if the delay vector τ con tains an y rep eated entries (for example, if τ 2 = τ 3 ), the duplicate channel can be remov ed. Moreo ver, this system can easily be extended to allo w for m ultiple lag pro ducts (e.g., n t lag pro ducts as in Figure 2a). The only c hannels that would need to b e added would b e those corresp onding to dela ys τ not already represented in the system. In fact, in most practical cases of interest known to the authors, THOCS estimation is very often restricted to a small set of unique delays τ of in terest. Examples include the Delay-R e duc e d Classiﬁer as prop osed in [23], [22] and the set of classiﬁcation features prop osed in [21] where all τ i = 0; in such a case, only one sampling c hannel w ould be needed. Nonetheless, we note here that one cannot exp ect signiﬁcant av erage sampling rate reductions in the case when a wide range of delays τ are required since the trade-oﬀ b et w een a verage sampling rate reduction and increased n umber of analog signal dela ys required no longer justiﬁes suc h a multi-c hannel sampling scheme. 3.2 The Lag Pro duct Measuremen ts W e now describ e how the lag pro duct measurements w lag (coming from the system in Figure 2b) can b e represented mathematically . Although the following hypothetical vector is not physically measured, for eac h channel i ∈ { 1 , 2 , . . . , n } , let us deﬁne 6 y ch i ,  x (0 − τ i ) x ( T s − τ i ) · · · x (( L − 1) T s − τ i )  T (7) 6 Classical windowing techniques could b e used in our framework to concentrate main lob e energy and reduce side lob e energy . 8 to be a collection of L uniform time samples, ha ving sampling in terv al T s (where 1 T s is at least the Nyquist rate), of the analog cyclostationary signal x ( t − τ i ) for the i th channel. A t the i th channel output, the nonuniform sample v ector w ch i can be though t of as b eing constructed b y ph ysically measuring only certain, randomly selected en tries of y ch i for each channel. F or m = 0 , 1 , . . . , L − 1, w e deﬁne a binary random v ariable a m indicating whether en try m of y ch i is selected for inclusion in w ch i ; we assume these a m are I ID Bernoulli random v ariables with P [ a m = 1] = γ and P [ a m = 0] = 1 − γ for some constan t γ ∈ (0 , 1]. As indicated in our diagram, we assume that all s/h devices are driv en by a common nonuniform clo c k; that is, we assume that a m tak es on the same v alue across diﬀeren t c hannels for the same entry m . Thus, on a verage, the parameter γ denotes the fraction of the total Nyquist rate samples retained by the nonuniform sampling proto col p er c hannel. W e let P denote the actual total num b er of nonuniform samples of x ( t ) collected p er channel; w e will ha ve P ≈ γ L , and we can write the sample vector as w ch i = Ry ch i , where R is a P × L selection matrix as explained in Section 2.2.2 and is identical across all n channels. As noted in Section 1.1, a typical cyclostationary signal is not compressible in the frequency domain. In particular, consider the uniform high-rate sample vector y ch i deﬁned in (7). Letting Φ denote the L × L DFT matrix, the DFT co eﬃcien ts of y ch i , given b y η ch i = Φ H y ch i , will not t ypically b e sparse or compressible. This raises the question of how we can hop e to apply CS reconstruction tec hniques given the non uniform samples w ch i = Ry ch i . The answer comes from the fact that w e do not attempt to reconstruct the full rate sample vec- tors y ch i in the individual c hannels. Rather, w e note that the ov erall system output w lag corresp onds to non uniform samples of the signal’s lag pro duct, which is compressible. Deﬁning x lag ( t ) , n Y i =1 x ( ∗ ) i ( t − τ i ) to b e the lag pro duct app earing in (1), we can deﬁne y lag to b e a h yp othetical collection of L uniform time samples of x lag ( t ): y lag [ k ] , x lag (( k − 1) T s ) = n Y i =1 y ( ∗ ) i ch i [ k ] . (8) F rom (6) and (8), it follows that w lag = Ry lag , where R is the same P × L selection matrix used in the individual channels. So, the system output w lag can b e viewed as a subset of the high-rate collection of L uniform time samples from the signal’s lag pro duct. W e note that such a computation is p ossible thanks to the commutativ e relationship of the lag product op eration with the non uniform sampling op erator: it is not ne c essary to c ompute the analo g lag pr o duct b efor e sampling the signal . F or suitable c hoices of n and q , the DFT co eﬃcien ts of y lag , given by η lag = Φ H y lag , may b e compressible. Indeed, when dealing with ﬁnite data in practice it is common to estimate THOCS directly from η lag : the lo cations and v alues of p eaks in the DFT sp ectrum provide estimates of the true α v alues and corresp onding cyclic moments R α x ( τ ) n,q , resp ectively . In situations where η lag is indeed compressible, we should b e able to apply CS techniques for reconstructing or estimating η lag from w lag . This is the approach w e adopt in the CTHOCS framework. 3.3 A Simple Estimator of the Lag Pro duct Sp ectrum W e can express the collected the non uniform samples of the lag pro duct as w lag = Ry lag = R Φ η lag = Aη lag , where A = R Φ. As w e explained in Section 2.2.1, when η lag is sparse or compressible, the 9 estimate b η lag = A H w lag = Φ H R H w lag , (9) can provide a suitable appro ximation to η lag that is v ery simple to compute. In particular, when R is a random selection matrix and Φ is a DFT matrix, the estimate b η lag = Φ H R H w lag can b e computed simply b y zero-padding w lag (creating a length- L vector containing the P en tries of w lag at the non uniform sample lo cations) and taking the FFT of this zero-padded v ector. This estimation pro cedure is not only simple (b oth conceptually and computationally), but as w e show in Section 3.4 it is also compatible with the THOCS analytical framew ork. In particular, w e propose to extract statistics from the lo cations and v alues of the p eaks of b η lag . Using the RIP , we can argue that the lo cations and v alues of the p eaks of b η lag should coincide roughly with those of η lag . Lemma 3.1. L et s denote the numb er of sp e ctr al p e aks in η lag and supp ose for simplicity that L = ms for some p ositive inte ger m . L et η lag ,s denote a length- L ve ctor c ontaining only the s p e aks fr om η lag (and zer os elsewher e), and let η lag ,t = η lag − η lag ,s denote the r emaining L − s entries of η lag . If A satisﬁes the RIP of or der s + 1 with isometry c onstant δ , the DFT estimate b η lag = A H w lag appr oximates the true DFT c o eﬃcients η lag up to the fol lowing ac cur acy:    b η lag [ k ] − η lag [ k ]    ≤ δ k η lag ,s k 2 + δ ( m − 1) √ s k η lag ,t k ∞ . (10) Pro of : See App endix A. As noted in Section 2.2.2, when Φ is an L × L DFT matrix and R is a random selection matrix with P ro ws, then A = R Φ (appropriately rescaled) will satisfy the RIP of order P / log 4 ( L ) with high probabilit y . Thus, if η lag is highly concentrated at its s peak v alues, one can expect η lag ,t to b e small, leading to a rough guaran tee in (10) that b η lag [ k ] ≈ η lag [ k ] for eac h co eﬃcien t k . While suc h a simple estimator may not b e the most accurate possible—the peak v alue of η lag ,t ma y not b e small unless one inv okes Lemma 3.1 with a large v alue of s to account for man y large entries in η lag , and most CS reconstruction algorithms are more sophisticated to accoun t for correlations among the columns of A —we opt for this FFT-based estimation technique b ecause of its simplicity and its compatibilit y with the THOCS analytical framework (see Section 3.4). As an illustration of this tec hnique, the solid blue curv e in Figure 3 sho ws (a zoomed p ortion of ) the magnitude of the DFT η lag of the lag pro duct of a 2PSK communication signal for τ = 0 . This curv e w as generated from uniform samples acquired at the Nyquist rate for this signal. The dashed red curv e in this plot sho ws b η lag , the zero-padded DFT of the lag pro duct computed from nonuniform samples of the signal with γ = 0 . 2 (so that P ≈ 0 . 2 L ). As we ha ve discussed, the lo cation and v alue of the p eak in the blue curve can be used to infer one of the cyclic momen ts of the signal. (Figure 3 only shows one of a few p eaks.) Con venien tly , we see a p eak in the red curve that— although w eaker—coincides with the p osition of the p eak in the blue curv e. This gives some hop e that information ab out the cyclic moments can b e recov ered from the low-rate non uniform samples of the signal. W e note here that while the w eaker p eak seems at o dds with the approximation error bound giv en in Lemma 3.1, the display ed red curv e has not b een appropriately normalized (rescaled). This shall b e discussed in the next section. W e shall also explain formally what statistics w e extract from b η lag and discuss how these relate to THOCS of the original cyclostationary signal in the next section. 10 Figure 3: Lag product DFT (solid blue curv e) and estimate obtained from lo w-rate nonuniform samples (dashed red curve). 3.4 Deﬁnition of CTHOCS Recall that η lag corresp onds to the DFT of uniform samples of the lag product of the original signal. Con ven tionally , one w ould estimate THOCS directly from η lag : the lo cations and v alues of p eaks in the DFT sp ectrum pro vide estimates of the true α v alues and corresp onding cyclic moments R α x ( τ ) n,q , resp ectiv ely . Although THOCS are formally deﬁned in terms of analog signals and con tinuous-time F ourier transforms, it is p ossible to sho w [20] that cyclic moments of the analog signal and cyclic momen ts of its uniformly sampled counterpart are equiv alent up to a p erio dic rep etition in the frequency domain, as o ccurs with any uniform sampling of an analog signal. When w e compute b η lag b y taking a zero-padded DFT of w lag , one w ay to view b η lag is simply as an estimate of η lag . Another wa y to view b η lag , how ever, is that it is the DFT of zero ed-out uniform samples of the lag pro duct of the original signal (with zeros in lo cations not sampled in the non uniform sampling proto col). W e adopt this latter viewpoint for the purpose of formally deﬁning CTHOCS: in short, we refer to THOCS of this zero ed-out signal as CTHOCS, w e can relate these CTHOCS to THOCS of the original signal, and w e can estimate CTHOCS from the lo cations and v alues of the p eaks in b η lag . 3.4.1 Preliminaries As explained ab o ve, the cyclic momen ts of the analog cyclostationary signal and cyclic momen ts of its uniformly sampled counterpart are eﬀectively equiv alent [20]. In our discussion, to av oid the complexities of dealing with Dirac delta functions (esp ecially when raising these functions to higher p ow ers), we frame our discussion around a s/h version of the delay ed analog cyclostationary signal for eac h c hannel of the prop osed multi-c hannel sampling scheme as follo ws. W e ﬁrst assume uniform sampling at each channel of the multi-c hannel sampler and deﬁne the s/h v ersion of the dela yed signal x ( t − τ i ) for the i th c hannel to b e x i s/h ( t, τ i ) , x ( t − τ i ) " ∞ X m = −∞ p s/h ( t − mT s ) # , (11) where m is the sample time index, T s is the sampling interv al (b ounded by the Nyquist sampling in terv al), and p s/h ( t ) is the s/h pulse shap e used in the prop osed m ulti-channel sampling scheme. F or conv enience, w e will use p s/h ( t ) = rect  t υ  , where υ is the pulse duration and υ  T s . Using (11), the s/h version of the lag product of signal x ( t ) for a ﬁxed τ (corresponding to the τ i ’s in 11 Figure 2b), x lag s/h ( t, τ ) n,q is deﬁned as x lag s/h ( t, τ ) n,q , n Y i =1 x ( ∗ ) i i s/h ( t, τ i ) . The F ourier co eﬃcien ts of x lag s/h ( t, τ ) n,q at the cycle frequencies will approximate the cyclic mo- men ts of x ( t ), and one can make this approximation arbitrarily accurate b y letting υ → 0 and normalizing appropriately . W e next consider nonuniform sampling at eac h c hannel of the multi-c hannel sampler and simi- larly deﬁne the s/h v ersion of the zero ed-out delay ed signal for the i th channel as follo ws: x i ncs ( t, τ i ) , x ( t − τ i ) " ∞ X m = −∞ a m p s/h ( t − mT s ) # , (12) where a m are the nonuniform sample indicators as deﬁned in Section 3.1, and the subscript “ncs” refers to nonuniform compressiv e samples. F or use in our analysis, w e note that the a m ’s are indep enden t, identically distributed (I ID) Bernoulli random v ariables where E [ a n m ] = γ , ∀ n ∈ Z + . In the F OT sense, a m can be though t of as taking on the v alue 1 with F OT probabilit y γ and taking on the v alue 0 with FOT probabilit y 1 − γ . Again using (12), the s/h v ersion of the zero ed-out lag pro duct of signal x ( t ) for a ﬁxed τ (corresp onding to the τ i ’s in Figure 2b), x lag ncs ( t, τ ) n,q can b e easily computed using the v arious channel outputs as x lag ncs ( t, τ ) n,q = n Y i =1 x ( ∗ ) i i ncs ( t, τ i ) . 3.4.2 CTHOCS Statistical F unctions W e deﬁne CTHOCS statistical functions, for a giv en τ , n and q , in the context of x lag ncs ( t, τ ) n,q , the s/h v ersion of the zero ed-out lag product of the incoming cyclostationary signal, and we compare these to THOCS of x lag s/h ( t, τ ) n,q , the s/h version of the lag pro duct of the same signal. Let us remind the reader that the main assumption underlying the deﬁnition of CTHOCS is that the lag pro duct of the cyclostationary signal is compressible in the frequency domain, i.e., the signal’s lag pro duct only contains a small n umber of unique dominant cycle frequencies (p erio dicities). 1. Time-v arying n th order compressive moments The time-v arying n th order moments of a cyclostationary signal are deﬁned in Section 2.1.2. W e denote the time-v arying n th order momen ts computed using the s/h version of the lag pro ducts x lag s/h ( t, τ ) n,q and x lag ncs ( t, τ ) n,q b y R x s/h ( t, τ ) n,q and R x ncs ( t, τ ) n,q , resp ectiv ely . Sp eciﬁcally , R x s/h ( t, τ ) n,q , ˆ E  x lag s/h ( t, τ ) n,q  and R x ncs ( t, τ ) n,q , ˆ E [ x lag ncs ( t, τ ) n,q ] . Our CTHOCS framew ork relies heavily on the fact that these statistics are strongly related. Lemma 3.2. The time-varying n th or der moments obtaine d using the nonuniformly sample d lag pr o duct signal x lag ncs ( t, τ ) n,q for a ﬁxe d τ , n and q satisfy R x ncs ( t, τ ) n,q = γ R x s/h ( t, τ ) n,q . (13) 12 Pro of : See App endix B. In other w ords, the time-v arying n th order moments computed using the collection of non uni- formly sampled signals { x i ncs ( t, τ i ) } are simply scale d b y a factor of γ compared to the time-v arying n th order momen ts computed using the collection of uniformly sampled signals { x i s/h ( t, τ i ) } for a giv en τ , n and q . In light of (13), w e deﬁne our ﬁrst CTHOCS quan tity: the time-varying n th or der c ompr essive moments of our signal are deﬁned to b e R c x ncs ( t, τ ) n,q , 1 γ R x ncs ( t, τ ) n,q . (14) (All CTHOCS quan tities will b e denoted with a sup erscript “ c ”.) Thanks to (13), it follo ws that the “compressiv e momen ts” R c x ncs ( t, τ ) n,q w e deﬁne are actually equal to the traditional momen ts R x s/h ( t, τ ) n,q of the original uniformly sampled lag pro duct of the signal for a given τ , n and q . In practice, one must approximate R c x ncs ( t, τ ) n,q based on a ﬁnite n umber of nonuniform samples of the signal at the output of each channel in Figure 2b. W e prop ose to do this b y computing the zero-padded DFT b η lag of the v ector w lag , preserving just the p eaks in b η lag , rescaling those p eaks by 1 γ , and then inv erting the DFT. The resulting estimate of R c x ncs ( t, τ ) n,q doubles as an appro ximation to the true moments R x s/h ( t, τ ) n,q . It is imp ortan t to note that the s trength of the p eaks of the zero-padded DFT reduces prop ortionally with γ (recall Figure 3), and so when dealing with a ﬁnite set of data, a peak detection scheme will b e necessary to iden tify the lo cations of the p eaks. When dealing with noisy data and high compression ratios ( γ  1) some p eaks may actually fall b elow the noise ﬂo or and will thus b e missed. The implications of dealing with ﬁnite data are discussed more thoroughly in Section 3.5. 2. n th order compressive cyclic moments W e deﬁne the n th or der c ompr essive cyclic moments of our signal by applying the cyclic momen t formula (2) to the time-v arying n th order compressive moments w e deﬁned in (14): R α,c x ncs ( τ ) n,q , h R c x ncs ( t, τ ) n,q e − j 2 π αt i t . (15) Again, thanks to (13) and (14), these “compressiv e cyclic momen ts” R α,c x ncs ( τ ) n,q are equal to the cyclic momen ts R α x s/h ( τ ) n,q of the uniformly sampled lag pro duct of the signal for a given τ , n and q . In practice, when dealing with ﬁnite data one can estimate these compressive cyclic moments simply b y using the lo cations and rescaled v alues of the p eaks identiﬁed in the estimation of the compressiv e moments R c x ncs ( t, τ ) n,q ab o v e. By setting τ 1 = 0 here again, all the information present in R α,c x ncs ( τ ) n,q is retained. 3. Time-v arying n th order compressive cumulan ts W e deﬁne the time-varying n th or der c ompr essive cumulants of our signal b y applying the momen t-to-cumulan t conv ersion form ula (3) to the time-v arying n th order compressiv e mo- men ts we deﬁned in (14): C c x ncs ( t, τ ) n,q , X D n " ( − 1) d − 1 ( d − 1)! d Y i =1 R c x ncs ( t, τ ) n i ,q i # , (16) 13 where each term R c x ncs ( t, τ ) n i ,q i can b e computed from the channel outputs in Figure 2b. Again, b ecause of (13) and (14), these “compressive cum ulants” C c x ncs ( t, τ ) n,q are equal to the traditional cumulan ts C x s/h ( t, τ ) n,q of the uniformly sampled lag pro duct of the signal for a giv en τ , n and q . 7 When dealing with ﬁnite data one can approximate the compressive moments R c x ncs ( t, τ ) n,q as discussed ab o ve and then plug these estimates in to (16) to obtain an estimate of the compressiv e cumulan ts since cumulan ts are not directly measurable quantities. If some p eaks w ere missed in the zero-padded DFT spectrum, the estimated compressive cumulan ts could diﬀer substan tially from the original signal cumulan ts C x s/h ( t, τ ) n,q . Again, b y setting τ 1 = 0, all the information presen t in C c x ncs ( t, τ ) n,q is retained. 4. n th order compressive cyclic cumulan ts W e deﬁne the n th or der c ompr essive cyclic cumulants of our signal b y applying the momen t- to-cum ulant con version form ula (5) to the n th order compressive cyclic momen ts w e deﬁned in (15): C β ,c x ncs ( τ ) n,q , X D n   ( − 1) d − 1 ( d − 1)! X α T 1 = β d Y i =1 R α i ,c x ncs ( τ ) n i ,q i   . (17) Again, thanks to (13) and (14), these “compressive cyclic cum ulants” are equal to the tra- ditional cyclic cumulan ts C β x s/h ( n, q ) of the uniformly sampled signal for a given τ , n and q . In practice, when dealing with ﬁnite data one can appro ximate the compressiv e cyclic mo- men ts R α,c x ncs ( n, q ) as discussed abov e and then plug these estimates into (17) to obtain an es- timate of the compressive cyclic cum ulants. Again, if some p eaks were missed, the estimated compressiv e cyclic cumulan ts could diﬀer substantially from the original signal’s cyclic cumu- lan ts C β x s/h ( n, q ). W e will demonstrate in Section 4, ho wev er, that the estimated compressiv e cyclic cum ulants still retain their signal selective prop ert y as compared to their conv en tional cyclic cum ulants. 3.5 Implications of Finite Observ ations W e conclude this section with a few additional remarks concerning the estimates of THOCS and CTHOCS from ﬁnite, rather than inﬁnite, data streams. Using ﬁnite observ ations of the uniformly sampled channel outputs in Figure 2b, an estimator of the cyclic momen t R α x s/h ( τ ) n,q referred to as b R α x s/h ( t, τ ) n,q can b e deﬁned as b R α x s/h ( t, τ ) n,q , 1 b T Z t + b T 2 − τ ∗ 0 2 t − b T 2 + τ ∗ 0 2 x lag s/h ( u, τ ) n,q e − j 2 π αu du, (18) where b T is the duration of a ﬁnite time windo w centered at time t applied to the incoming analog signal x ( u ) b efore its b eing fed into the multi-c hannel sampler of Figure 2b, τ l , max { τ 1 , · · · , τ n } , 7 Because the pro duct in (16) contains a v ariable num b er of terms d , the “compressive cumulan ts” C c x ncs ( t, τ ) n,q will not hav e a simple linear relationship to the traditional cumulan ts C x ncs ( t, τ ) n,q of the zero ed-out lag pro duct signal. 14 τ r , min { τ 1 , · · · , τ n } and τ ∗ 0 , τ r − τ l . The time-v arying mean of this estimator is equiv alent to that prop osed in [19] and (from [19]) is giv en by 1 b T R α x s/h ( τ ) n,q ( τ ∗ 0 + b T ) + 1 b T X b 6 = α R b x s/h ( τ ) n,q e − iπ ( α − b )( τ l + τ r ) ( τ ∗ 0 + b T )sinc h ( α − b )( τ ∗ 0 + b T ) i e − i 2 π ( α − b ) t . (19) This equation reveals a bias in the estimator due to contributions by other (cycle) frequencies. This bias v anishes as b T → ∞ . By extending this analysis, one can show that when estimating a compressiv e cyclic moment R α,c x ncs ( τ ) n,q from observ ations of the nonuniformly sampled channel outputs in Figure 2b with ﬁnite duration b T cen tered at time t , e.g., by computing the estimator b R α,c x ncs ( t, τ ) n,q , 1 γ b T Z t + b T 2 − τ ∗ 0 2 t − b T 2 + τ ∗ 0 2 x lag ncs ( u, τ ) n,q e − j 2 π αu du, (20) the mean of the estimator is again given by (19). The time-a veraged v ariance of a cyclostationary time series x ( t ) is deﬁned [20] as V ar 0 { x ( t ) } , h| x ( t ) | 2 i t − h| ˆ E { x ( t ) }| 2 i t . (21) Hence, the time-av eraged v ariance of the compressiv e cyclic moment estimator b R α,c x ncs ( t, τ ) n,q is giv en b y V ar 0 n b R α,c x ncs ( t, τ ) n,q o =     b R α,c x ncs ( t, τ ) n,q    2  t −     ˆ E n b R α,c x ncs ( t, τ ) n,q o    2  t , (22) where     b R α,c x ncs ( t, τ ) n,q    2  t = * 1 γ 2 b T 2 Z t + b T 2 − τ ∗ 0 2 t − b T 2 + τ ∗ 0 2 Z t + b T 2 − τ ∗ 0 2 t − b T 2 + τ ∗ 0 2 x lag ncs ( r , τ ) n,q  x lag ncs ( s, τ ) n,q  ∗ e − j 2 π α ( r − s ) dr ds + t (23) is the mean square absolute v alue of the estimator and ˆ E n b R α,c x ncs ( t, τ ) n,q o is giv en by (19). Lemma 3.3. The me an squar e absolute value of the estimator     b R α,c x ncs ( t, τ ) n,q    2  t for a ﬁxe d τ c an b e shown to satisfy     b R α,c x ncs ( t, τ ) n,q    2  t = 1 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t ×   1 γ * ∞ X m = −∞ p n s/h ( r + t − mT s ) p n s/h ( s + t − mT s ) + t + * X u 6 = v p n s/h ( r + t − uT s ) p n s/h ( s + t − v T s ) + t   × e − j 2 π α ( r − s ) dr ds. (24) Pro of : See App endix C. W e note that when γ = 1,     b R α,c x ncs ( t, τ ) n,q    2  t equals     b R α,c x s/h ( t, τ ) n,q    2  t whic h results in iden tical v ariance when comparing b oth estimators. Since (24) can b e split in to tw o integrals, 15 whic h m ust add to a real, nonnegativ e n umber for all γ ∈ (0 , 1], it follo ws that     b R α,c x ncs ( t, τ ) n,q    2  t m ust increase as γ decreases. This, in turn, causes V ar 0 n b R α,c x ncs ( t, τ ) n,q o to increase if b T is ﬁxed and γ decreases; ho wev er, this increase in v ariance can b e mitgated by increasing b T . While our main goal in this section is to establish the v ariance of our compressiv e cyclic moment estimator in relation to its full rate uniform sampling cyclic moment estimator coun terpart, the interested reader is referred to [3] and [19] for more detailed discussions on the v ariance of the full rate uniform sampling cyclic momen t estimator as w ell as the cyclic cum ulant estimator. In the following section, w e fo cus our discussions on the signal selectivit y asp ect of CTHOCS and v alidate our theorems with sim ulations. 4 CTHOCS Signal Selectivit y As highligh ted in our introduction and in the seminal pap ers of Gardner et al. [10, 20], t wo highly desirable prop erties of cyclic cum ulants are signal selectivit y and tolerance to Gaussian noise (when n > 2). Therefore, w e devote this section to empirically v alidating the signal selectivity of our prop osed compressive cyclic cum ulants. Our prop osed CTHOCS only requires the lag pro duct of the cyclostationary signal to b e com- pressible in the frequency domain. One example of suc h a cyclostationary signal is a baseband digital Quadrature Amplitude Mo dulation (Digital QAM) comm unication signal, as we explain in the next section. 4.1 Baseband Digital QAM Signal Mo del W e assume the cyclostationary signal to b e acquired is a baseband digital QAM comm unication signal of the form s ( t ) = ∞ X k = −∞ as k p ( t − k T − t 0 ) e j (2 π ∆ f c t + θ c ) , (25) where a is the signal amplitude, s k is the k th transmitted sym b ol, p ( t ) is the signaling pulse shap e, T is the sym b ol p erio d, ∆ f c is the residual carrier frequency oﬀset, and θ c is the carrier phase. W e mak e the common assumption that ∆ f c  1 T , and we also assume that symbols s k are I ID, uniformly c hosen from some ﬁnite dictionary depending on the modulation t yp e of the signal. In the sections that follow, the notation s ( t ) is reserved sp eciﬁcally for an incoming cyclostationary signal that ob eys the mo del ab o ve. Henceforth, in the subsequent sections, w e shall restrict ourselves to signals b elonging to the signal mo del of (25). 4.2 Signal Selectivity of Cyclic Cumulan ts W e now pro vide a brief literature review of previous w orks fo cused on the signal selectivity of cyclic cum ulants for the Digital QAM signal mo del (25). W e note here that while other pap ers may focus on the signal selectivit y of other signal statistics suc h as non-cyclic cumulan ts or spectral HOCS for our signal mo del (25), our aim here is to establish the v alidity of compressive cyclic cumulan ts using cyclic cum ulants as a b enc hmark. Sp ooner considered the classiﬁcation of b oth single and multiple carrier (co-c hannel) signals b elonging to (25) using cyclic cum ulants in [21]. The signal selective features proposed were F s = n max τ    C β s ( τ ) n,q    o , (26) 16 for n = 2 , 4 , . . . , N , q = 0 , 1 , . . . , n , τ = [ τ 1 . . . τ n ], β = ( n − 2 q )∆ f c + k T , k = 0 , . . . , K , and C β s ( τ ) n,q = a n T C s n,q e j θ c ( n − 2 q ) e j 2 π ∆ f c P n m =1 ( − ) m τ m e − j 2 π k T t 0 " Z ∞ −∞ n Y ` =1 p ( t + τ ` ) ! e − j 2 π k T t dt # . (27) Here, C s n,q = P D n ( − 1) d − 1 ( d − 1)! Q d i =1 R s,v i is known as the sym b ol cumulan ts parameter, R s,v i are known as the symbol moments, and ( − ) m denotes an optional min us sign dep ending on the optional conjugation. The symbol dep enden t parameters R s,v i can b e computed from the symbols s k in (25) using R s,v i = E h s | v i |− q i k ( s ∗ k ) q i i , where q i is the n umber of conjugated indices in the set v i . The prop osed features from the incoming signal were subsequently estimated assuming kno wn cycle frequencies (v alues of β ), normalized, group ed in vector form and classiﬁed using a minimum Euclidean distance approac h with theoretical feature vectors. Classiﬁcation results were obtained for scenarios in volving six (single carrier) and ten classes (four co-c hannel carriers). W e note here that for each n , only the maximum v alue of (27) w as used and this would o ccur for the case when τ = 0 (the length- n v ector of zeros) based on insp ection of (27) and giv en the deca ying nature of typical comm unication signal pulse shap es p ( t ) [6]. As a revision to [21], Spo oner et al. demonstrated the v alidity of cyclic cumulan ts, for b oth single carrier and multiple co-channel carriers, when estimated with no prior information about the signals in [23] using warped cyclic cumulan ts (a v ariant of the cyclic cumulan t features prop osed in [21]). A costly maximum likelihoo d (ML) classiﬁer w as derived and subsequently approximated with tw o realizable alternativ es, namely , the order-reduced classiﬁer (ORC) and the delay-reduced classiﬁer (DR C). Classiﬁcation results were obtained for scenarios in volving ﬁve classes (single carrier) and up to three co-c hannel carriers. In [22], Sp o oner extended the use of the fourth order cyclic cumulan ts (used in [21, 23]) to the sixth order and demonstrated that classiﬁcation p erformance c an b e further improv ed with its use. Classiﬁcation results w ere obtained for up to four class scenarios in volving tw o co-channel signals. Marc hand et al. considered the use of a linear com bination of cyclic cum ulants containing second and fourth order cyclic cum ulants for classiﬁcation of digital QAM signals in [14, 15]. Finally , Dobre et al. inv estigated the use of cyclic cumulan t features for classiﬁcation of digital QAM signals up to the eighth order in the single carrier scenario in [6] assuming all parameters of the signal are kno wn. 4.3 Signal Selectivity of Compressiv e Cyclic Cumulan ts As mentioned in the previous section, our goal here is to establish the v alidity of compressive cyclic cum ulants using cyclic cumulan ts as a b enc hmark. Therefore, we shall b e using the signal selective cyclic cum ulant features as prop osed in [21] as well as the minimum distance classiﬁcation scheme prop osed in [21]. W e note here that in subsequent sim ulations, a four class scenario comprising of 2PSK, 4PSK, 8PSK and 16QAM single carrier classiﬁcation will b e considered. T herefore, a single feature minim um classiﬁcation scheme, namely using F s =    C β ,c s ncs ( 0 ) 4 , 0    , (28) where β = 4∆ f c suﬃces. In this instance, using (27), the theoretical features can b e computed as F s = a 4 T | C s 4 , 0 |     Z ∞ −∞ p ( t ) 4 dt     , (29) 17 where | C s 4 , 0 | takes the v alues 2, 1, 0, and 0 . 68 for 2PSK, 4PSK, 8PSK and 16QAM respectively . In our simulations, we assume the signal amplitude a is known and the interested reader is referred to [21] for a detailed discussion on feature normalization when v arious signal parameters of (25) are unkno wn. 4.4 W alkthrough CTHOCS Estimation Example W e now presen t an example of how b C β ,c s ncs ( τ ) 4 , 2 , an estimate of C β ,c s ncs ( τ ) 4 , 2 can b e obtained from ﬁnite sets of non uniform samples of a cyclostationary signal belonging to (25) obtained from our prop osed multi-c hannel sampling scheme. 1. Using the collection of w ch i lo w-rate nonuniform sample v ectors, compute its lag pro duct w lag using (6) for n = 4 and q ∈ { 0 , 1 , 2 } . F or each q ∈ { 0 , 1 , 2 } , compute b η lag , the zero-padded DFT of w lag , using (9). 2. F or each q ∈ { 0 , 1 , 2 } , detect the p eaks in b η lag , rescale these p eaks b y 1 γ , and denote these rescaled v alues as the estimated compressiv e cyclic momen ts b R α,c s ncs ( τ ) n,q . F or eac h peak detected, α denotes the cycle frequency at which the p eak was detected. 3. F orm an estimate of the compressiv e cyclic cumu lant C β ,c s ncs ( τ ) 4 , 2 (denote it as b C β ,c s ncs ( τ ) 4 , 2 ) b y plugging the estimated compressiv e cyclic moments b R α,c s ncs ( τ ) n,q in to (17). F or the sake of clarit y , we pro vide the formula for computing b C β ,c s ncs ( τ ) 4 , 2 : b C β ,c s ncs ( τ ) 4 , 2 = X D 4   ( − 1) d − 1 ( d − 1)! X α T 1 = β d Y i =1 b R α i ,c s ncs ( τ ) n i ,q i   . Before proceeding further, it may b e helpful to examine the set D 4 , whic h is the set con taining all p ossible partitions of the index set { 1 2 3 4 } , as well as the P term inside the square brac kets. Since w e consider the case where n = 4 and q = 2, w e shall assume the ﬁrst t wo indices are the conjugated terms and denote the index set as { 1 ∗ 2 ∗ 3 4 } . Then D 4 con tains the follo wing partitions of { 1 ∗ 2 ∗ 3 4 } : (a) { 1 ∗ 2 ∗ 3 4 } : F or this partition, d = 1 (since there is only one factor) and the only compressiv e cyclic moment estimate w e require is b R β ,c s ncs ( τ ) 4 , 2 . (b) { 1 ∗ 2 ∗ }{ 3 4 } : F or this partition, d = 2 (since there are tw o factors) and the P term has the form P b R α i ,c s ncs ( τ ) 2 , 2 b R α j ,c s ncs ( τ ) 2 , 0 o ver all α i and α j suc h that α i + α j = β . W e require the compressiv e cyclic moment estimates b R α i ,c s ncs ( τ ) 2 , 2 and b R α j ,c s ncs ( τ ) 2 , 0 due to the factors { 1 ∗ 2 ∗ } and { 3 4 } , resp ectiv ely . (c) { 1 ∗ 3 }{ 2 ∗ 4 } : F or this partition, d = 2 (since there are tw o factors) and the P term has the form P b R α i ,c s ncs ( τ ) 2 , 1 b R α j ,c s ncs ( τ ) 2 , 1 o ver all α i and α j suc h that α i + α j = β . (d) { 1 ∗ 4 }{ 2 ∗ 3 } : F or this partition, d = 2 and the P term has the form P b R α i ,c s ncs ( τ ) 2 , 1 b R α j ,c s ncs ( τ ) 2 , 1 o ver all α i and α j suc h that α i + α j = β . Recall that THOCS (and thus CTHOCS) are nonzero only for even v alues of n and this holds for signals b elonging to the mo del of (25), hence we hav e eliminated partitions containing factors with o dd num b er of indices. 18 4.5 CTHOCS V alidation 4.5.1 Sim ulation T o v alidate the inference quality and signal selectivity quality of CTHOCS, Monte Carlo (MC) sim ulations were p erformed for 2PSK, 4PSK, 8PSK and 16QAM signal mo dulation types. The follo wing signal parameters w ere k ept constan t for all sim ulations: signal amplitude ( a = 1), residual carrier frequency oﬀset (∆ f c = 23 . 0625Hz), symbol rate ( 1 T = 12999 . 5625 Hz, 6499 . 5625 Hz, 3249 . 5625 Hz, 1624 . 5625 Hz, 799 . 5625 Hz and 399 . 5625 Hz for 13000, 6500, 3250, 1625, 800, and 400 pro cessing symbols resp ectiv ely), R C pulse shap e with roll-oﬀ factor of 0 . 3, and signal duration ( b T = 1 second), and all of these parameters were assumed to b e known. While prior kno wledge of the parameters simpliﬁed the estimation of the CTHOCS features, a p eak detection algorithm was necessary (and used in the simulations) since the actual n um b er of p eaks (in the lag sp ectrum of the signal) to be detected for a sp eciﬁc signal parameter conﬁguration was not known b eforehand. Sp eciﬁcally , eac h lag pro duct sp ectrum w as tested for the presence of a p eak up to the 6th harmonic frequency (a total of 12 including the negative frequencies). The absolute v alue at each harmonic frequency was then compared against those of its neigh b oring frequencies with a p eak detected if the absolute v alue was greater than an empirically pre-determined threshold. Due to the v arious symbol rates considered, sampling frequencies 1 T s = 131072 Hz, 65536 Hz, 32768 Hz, 16384 Hz, 8192 Hz, and 4096 Hz w ere chosen for symbol rates 12999 . 5625 Hz, 6499 . 5625 Hz, 3249 . 5625 Hz, 1624 . 5625 Hz, 799 . 5625 Hz, and 399 . 5625 Hz for all simulations. F or each signal t yp e, 50 trials w ere p erformed. 4.5.2 Compressiv e Cyclic Cumulan ts Inference Quality Before w e pro ceed to establish the signal selectivit y of compressiv e cyclic cum ulants, in this section, w e v alidate their inference qualit y . Sp eciﬁcally , we compare estimated compressiv e cyclic cum ulants using v arious pro cessing lengths (in terms of num b er of symbols) for select v alues of n , q across a broad range of γ against cyclic cumulan ts (which corresp ond to γ = 1 in the preceding plots) using the normalized mean square error (NMSE) NMSE = 1 N tr N tr X i =1 h b C β ,c s ncs ( 0 ) n,q − C β ,c s ncs ( 0 ) n,q i 2 h C β ,c s ncs ( 0 ) n,q i 2 , (30) where N tr is the num b er of trials, β = ( n − 2 q )∆ f c and C β ,c s ncs ( 0 ) n,q is computed using (27). W e note that since there is only one unique v alue of τ i = 0, the prop osed multi-c hannel sampling sc heme of Figure 2b collapses in to a single c hannel and subsequent compression rates obtained represent the b est compression rates achiev able. The plots in Figures 4, 5, 6, and 7 sho w the NMSE for the situations when no noise is added to the signal and also when noise is added suc h that resulting carrier-to-noise ratio (CNR) is 9dB, 6dB and 3dB resp ectiv ely . The NMSE for each 8PSK signal has b een omitted from the plots in the ﬁrst column in Figures 4–7 since an 8PSK cyclic cumulan t for n = 4 and q = 0 is zero. W e make the following observ ations with regards to the inference quality of compressive cyclic cum ulants. First, as n increases, the NMSE of compressiv e cyclic cumulan t estimates increases across the entire range of γ v alues considered. This also o ccurs in the case of cyclic cumulan t estimates. Second, the NMSE of compressiv e cyclic cum ulant estimates increases as the n umber of pro cessing sym b ols decreases. This also o ccurs in the case of cyclic cum ulant estimates. Third, as the compression factor increases (i.e., as γ decreases), the NMSE of the compressiv e cyclic cumulan t 19 estimates increases at a minimal rate of 1 γ as explained in the sequel. Equation (24) gives the time- a veraged v ariance of the compressive cyclic moment estimator whic h can b e used to predict the rate of increase in the NMSE of compressiv e cyclic cumulant estimates (Figures 4–7), for 4PSK and 16QAM signals (when n = 4 and q = 0) since for these signals, C α,c s ncs ( τ ) 4 , 0 = R α,c s ncs ( τ ) 4 , 0 . As the n umber of pro cessing symbols decreases (i.e., as b T decreases), for a ﬁxed γ , the increase in NMSE can b e explained b y the increase in time-a veraged v ariance of the compressiv e cyclic cumulan t estimates due to the 1 b T 2 scaling in both terms of (22) and the increase in bias due to (19). On the other hand, for a ﬁxed num b er of pro cessing symbols (i.e., ﬁxed b T ) the NMSE increases at rate 1 γ due to the 1 γ scaling in (24). F or the other compressive cyclic cumulan ts such as C β ,c s ncs ( τ ) 4 , 2 and C β ,c s ncs ( τ ) 6 , 3 , their NMSE increases at a rate greater than 1 γ whic h is exp ected since C β ,c s ncs ( τ ) 4 , 2 and C β ,c s ncs ( τ ) 6 , 3 are b oth functions of higher p o w ers of low er order compressive cyclic moments. As we shall discuss in the next section, in some regimes, compressive cyclic cum ulants can giv e reliable appro ximations to their cyclic cumulan t counterparts. 4.5.3 Compressiv e Cyclic Cumulan t Signal Selectivity Qualit y In this section, w e establish the signal selectivity quality of b F s , an estimate of the compressive cyclic cum ulant feature F s (cf. (28)). Figures 8–11 show the confusion matrix when the minimum distance classiﬁcation sc heme was used to classify a four class scenario in the noiseless setting and with CNR of 9dB, 6dB and 3dB resp ectiv ely . In each of the confusion matrices, × , +, # and 2 represent 2PSK, 4PSK, 8PSK and 16QAM resp ectiv ely . Each row in these ﬁgures shows the confusion matrices for γ = 0 . 1, 0.5 and 1. W e make the following observ ations with regards to the classiﬁcation p erformance using com- pressiv e cyclic cumulan t estimates. First, as the compression factor increases (i.e., as γ decreases), the classiﬁcation p erformance degrades for a given CNR and n umber of processed symbols. The degradation of classiﬁcation p erformance due to the decreasing γ v alue is expected since the strengths of the sp ectral p eaks (compressive cyclic momen ts) decrease as predicted by Lemma 3.2, and this results in increased estimation errors of b F s . Second, the classiﬁcation p erformance de- grades with decreasing CNR as well as decreasing num b er of pro cessed symbols. This is due to the increased estimation errors asso ciated with increased noise p o wer and a reduced observ ation in terv al. This trend is also observed in the classiﬁcation p erformance of cyclic cum ulant based classiﬁers. The observed trends are eviden t in Figure 12, which sho ws the a verage correct classiﬁcation probabilit y P cc for signals ha ving CNRs of 3dB, 6dB and 9dB. T able 2 shows a comparison of classiﬁcation p erformance b et ween compressive cyclic cum ulant based classiﬁer ( γ = 0 . 1) and cyclic cum ulant based classiﬁer ( γ = 1) for v arious n umber of pro- cessed symbols. W e note that when γ = 0 . 1, pro cessing 3250 sym b ols corresp onds to using 3277 lo w-rate nonuniform samples to estimate compressiv e cyclic cum ulants. On the other hand, when γ = 1, pro cessing 400 symbols corresponds to using 4096 uniform samples to estimate cyclic cumu- lan ts. Eviden tly , the classiﬁcation p erformance of compressive cyclic cumulan t based classiﬁers are similar to cyclic cumulan t based classiﬁers when making pairwise comparisons, e.g., 3250 symbols v ersus 400 sym b ols, 6500 sym b ols v ersus 800 sym b ols, 13000 sym b ols versus 1625 symbols. Let us remind the reader that when γ = 0 . 1, the av erage sampling rate commensurates with the infor- mation rate of the signal. Due to the compressibilit y of the lag pro duct of the signal, compressiv e cyclic cumulan t estimates can give reliable approximations of their cyclic cumulan t counterparts at a fraction of the original sampling rate. The results of T able 2 also provide the follo wing insights regarding the use of low-rate nonuni- 20 (a) 13000 Symbols: n = 4 , q = 0 (b) 13000 Sym b ols: n = 4 , q = 2 (c) 13000 Sym b ols: n = 6 , q = 3 (d) 1625 Symbols: n = 4 , q = 0 (e) 1625 Symbols: n = 4 , q = 2 (f ) 1625 Symbols: n = 6 , q = 3 (g) 800 Symbols: n = 4 , q = 0 (h) 800 Symbols: n = 4 , q = 2 (i) 800 Sym b ols: n = 6 , q = 3 (j) 400 Symbols: n = 4 , q = 0 (k) 400 Sym b ols: n = 4 , q = 2 (l) 400 Symbols: n = 6 , q = 3 Figure 4: NMSE of compressive cyclic cumulan ts against cyclic cumulan ts. Across each row, the plots show NMSE of compressiv e cyclic cumulan ts v ersus cyclic cumulan ts for select v alues of n and q as a function of γ for select pro cessed data length under noiseless conditions. 21 (a) 13000 Symbols: n = 4 , q = 0 (b) 13000 Sym b ols: n = 4 , q = 2 (c) 13000 Sym b ols: n = 6 , q = 3 (d) 1625 Symbols: n = 4 , q = 0 (e) 1625 Symbols: n = 4 , q = 2 (f ) 1625 Symbols: n = 6 , q = 3 (g) 800 Symbols: n = 4 , q = 0 (h) 800 Symbols: n = 4 , q = 2 (i) 800 Sym b ols: n = 6 , q = 3 (j) 400 Symbols: n = 4 , q = 0 (k) 400 Sym b ols: n = 4 , q = 2 (l) 400 Symbols: n = 6 , q = 3 Figure 5: NMSE of compressive cyclic cumulan ts against cyclic cumulan ts. Across each row, the plots show NMSE of compressiv e cyclic cumulan ts v ersus cyclic cumulan ts for select v alues of n and q as a function of γ for select pro cessed data length when CNR=9dB. 22 (a) 13000 Symbols: n = 4 , q = 0 (b) 13000 Sym b ols: n = 4 , q = 2 (c) 13000 Sym b ols: n = 6 , q = 3 (d) 1625 Symbols: n = 4 , q = 0 (e) 1625 Symbols: n = 4 , q = 2 (f ) 1625 Symbols: n = 6 , q = 3 (g) 800 Symbols: n = 4 , q = 0 (h) 800 Symbols: n = 4 , q = 2 (i) 800 Sym b ols: n = 6 , q = 3 (j) 400 Symbols: n = 4 , q = 0 (k) 400 Sym b ols: n = 4 , q = 2 (l) 400 Symbols: n = 6 , q = 3 Figure 6: NMSE of compressive cyclic cumulan ts against cyclic cumulan ts. Across each row, the plots show NMSE of compressiv e cyclic cumulan ts v ersus cyclic cumulan ts for select v alues of n and q as a function of γ for select pro cessed data length when CNR=6dB. 23 (a) 13000 Symbols: n = 4 , q = 0 (b) 13000 Sym b ols: n = 4 , q = 2 (c) 13000 Sym b ols: n = 6 , q = 3 (d) 1625 Symbols: n = 4 , q = 0 (e) 1625 Symbols: n = 4 , q = 2 (f ) 1625 Symbols: n = 6 , q = 3 (g) 800 Symbols: n = 4 , q = 0 (h) 800 Symbols: n = 4 , q = 2 (i) 800 Sym b ols: n = 6 , q = 3 (j) 400 Symbols: n = 4 , q = 0 (k) 400 Sym b ols: n = 4 , q = 2 (l) 400 Symbols: n = 6 , q = 3 Figure 7: NMSE of compressive cyclic cumulan ts against cyclic cumulan ts. Across each row, the plots show NMSE of compressiv e cyclic cumulan ts v ersus cyclic cumulan ts for select v alues of n and q as a function of γ for select pro cessed data length when CNR=3dB. 24 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (a) 13000 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (b) 13000 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (c) 13000 Symbols: γ = 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (d) 1625 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (e) 1625 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (f ) 1625 Sym b ols: γ = 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 1 0 49 (g) 800 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (h) 800 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 1 0 49 (i) 800 Symbols: γ = 1 × + # 2 × 50 0 0 0 + 0 49 0 1 # 0 0 50 0 2 0 2 1 47 (j) 400 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 3 0 47 (k) 400 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (l) 400 Symbols: γ = 1 Figure 8: Confusion matrix using minim um distance classiﬁer based on b F s for signals under noiseless setting where × , +, # and 2 represent 2PSK, 4PSK, 8PSK and 16QAM respectively . form sampling (in the CTHOCS case) versus the use of high rate uniform sampling (in the THOCS case). F or a ﬁxed observ ation time, CTHOCS-based classiﬁcation p erformance degrades with de- creasing γ . On the other hand, for an equiv alent num b er of samples, non uniform sampling and uniform sampling hav e comparable classiﬁcation performance. While uniform sampling w ould al- lo w a shorter observ ation time, a higher rate sampler is required whic h could b e exp ensive or hav e limited resolution. Therefore, non uniform sampling can viewed as a means to safely sample slo wer than Nyquist rate without risking aliasing (whic h is not supported by uniform sampling). If the degradation in CTHOCS-based classiﬁcation p erformance is not acceptable (due to decreasing γ ), an increase in observ ation in terv al w ould ultimately impro ve CTHOCS-based classiﬁcation p erfor- mance such that it is comparable to THOCS-based classiﬁcation p erformance for the same num b er of measurements. W e note that when acquiring signals of high symbol rates where a lo w-rate non uniform sampling is truly useful, the eﬀect of increasing the observ ation in terv al is mitigated since man y symbols can still b e observed in a relatively short time span for suc h signals. 25 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (a) 13000 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (b) 13000 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (c) 13000 Symbols: γ = 1 × + # 2 × 50 0 0 0 + 0 49 0 1 # 0 0 50 0 2 0 1 0 49 (d) 1625 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (e) 1625 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (f ) 1625 Sym b ols: γ = 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 6 0 44 (g) 800 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (h) 800 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 3 0 47 (i) 800 Symbols: γ = 1 × + # 2 × 49 0 0 1 + 0 42 0 8 # 0 0 50 0 2 0 12 1 38 (j) 400 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 6 0 44 (k) 400 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 1 0 49 (l) 400 Symbols: γ = 1 Figure 9: Confusion matrix using minimum distance classiﬁer based on b F s for signals having CNR of 9dB. 3250 400 6500 800 13000 1625 Sym b ols 0.1 1 0.1 1 0.1 1 γ 3277 4096 6554 8192 13107 16384 Samples CNR 3dB 96% 93.5% 96.5 % 97% 99.5% 100% 6dB 100% 99% 99.5% 98.5% 100% 100% 9dB 100% 99.5% 100% 98.5% 100% 100% T able 2: Comparison of av erage P cc for select v alues of pro cessed symbol length b etw een γ = 0 . 1 and γ = 1 for signals ha ving CNRs of 3dB, 6dB and 9dB. 26 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (a) 13000 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (b) 13000 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (c) 13000 Symbols: γ = 1 × + # 2 × 50 0 0 0 + 0 47 0 3 # 0 0 50 0 2 0 3 0 47 (d) 1625 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (e) 1625 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (f ) 1625 Sym b ols: γ = 1 × + # 2 × 50 0 0 0 + 0 42 0 8 # 0 0 50 0 2 0 9 4 37 (g) 800 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (h) 800 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 3 0 47 (i) 800 Symbols: γ = 1 × + # 2 × 49 1 0 0 + 0 39 0 11 # 0 0 50 0 2 0 14 1 35 (j) 400 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 47 0 3 # 0 0 50 0 2 0 5 0 45 (k) 400 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 2 0 48 (l) 400 Symbols: γ = 1 Figure 10: Confusion matrix using minimum distance classiﬁer based on b F s for signals having CNR of 6dB. 27 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (a) 13000 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (b) 13000 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (c) 13000 Symbols: γ = 1 × + # 2 × 50 0 0 0 + 0 42 0 8 # 0 0 50 0 2 0 6 1 43 (d) 1625 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 47 0 3 # 0 0 50 0 2 0 3 0 47 (e) 1625 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 0 0 50 (f ) 1625 Sym b ols: γ = 1 × + # 2 × 50 0 0 0 + 1 31 15 3 # 0 0 50 0 2 1 21 9 19 (g) 800 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 45 0 5 # 0 0 50 0 2 0 5 0 45 (h) 800 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 50 0 0 # 0 0 50 0 2 0 6 0 44 (i) 800 Symbols: γ = 1 × + # 2 × 44 6 0 0 + 6 33 3 8 # 0 0 50 0 2 2 15 9 24 (j) 400 Symbols: γ = 0 . 1 × + # 2 × 50 0 0 0 + 0 43 0 7 # 0 0 50 0 2 0 12 0 38 (k) 400 Symbols: γ = 0 . 5 × + # 2 × 50 0 0 0 + 0 44 0 6 # 0 0 50 0 2 0 7 0 43 (l) 400 Symbols: γ = 1 Figure 11: Confusion matrix using minimum distance classiﬁer based on b F s for signals having CNR of 3dB. Average P cc (a) CNR:3dB Average P cc (b) CNR:6dB Average P cc (c) CNR:9dB Figure 12: Comparison of av erage P cc for select v alues of pro cessed data length when CNR=3dB, 6dB and 9dB. 28 5 Conclusion In extending the theory of THOCS to CTHOCS, w e ha ve prop osed the use of a low-rate non uniform sampling proto col and a lo w complexit y technique for estimating compressiv e cyclic cumulan ts from these nonuniform samples. W e also considered sp eciﬁcally the signal selectivity qualit y of compressiv e cyclic cumulan ts and sa w promising classiﬁcation results for a four class scenario even at sampling rates which commensurate with the information rate of the signal. Our w ork could p oten tially enable mo dulation classiﬁcation for muc h higher bandwidth signals than is currently p ossible due to sampling hardw are limitations. There are op en questions not yet addressed b y our pap er. In our signal mo del, we ha v e not accoun ted for p oten tial oscillator drifts common in most comm unication radios and sensing systems. Hence, the impact of carrier frequency drifts and timing oﬀset drifts on CTHOCS remains to b e seen. In our simulations, w e ha ve also not accoun ted for drifts in the sampling time instants of the non uniform sampler, although in other con texts CS tests on actual non uniform sampling hardw are ha ve b een promising [24]. F urther theoretical analysis of CTHOCS could also allow a b etter characterization of ho w small one can actually choose γ in practice. In this pap er, we also did not thoroughly discuss the detection asp ect of sp ectral p eaks (com- pressiv e cyclic moments) whic h is fundamen tal for the accurate estimation of compressive cyclic cum ulants or cyclic cumulan ts. Though the p oten tial lo cations of the sp ectral p eaks were assumed to b e kno wn in our sim ulations ( α = ( n − 2 q )∆ f c + k T ), a p eak detection scheme was emplo yed to detect the p eaks since the n umber of p eaks ( k ) v aries for eac h signal conﬁguration. That is, in our sim ulations we knew wher e to lo ok for p eaks but we did not assume the presence or absence of any p eak was kno wn a priori. In practice, it would b e necessary to compute the p o w er of these sp ectral p eaks (given v arious signal parameters) and accordingly , design robust p eak detection algorithms when noise is present in the incoming signal for a giv en lev el of detection probabilit y , false alarm rate and compression rate ( γ ). W e ha ve used the standard Cell-av eraging Constan t F alse Alarm Rate (CA-CF AR) algorithm to detect the sp ectral p eaks in our sim ulations. A Pro of of Lemma 3.1 Our argumen t uses the following result. Lemma A.1 ( [5]) . L et u , v ∈ R L b e given, and supp ose that a matrix A satisﬁes the RIP of or der max ( k u + v k 0 , k u − v k 0 ) with isometry c onstant δ . Then |h Au, Av i − h u, v i| ≤ δ k u k 2 k v k 2 . (31) Let e k ∈ R L denote the k th cardinal basis v ector, equal to 1 in its k th en try and 0 elsewhere. Letting v = e k in Lemma A.1, w e hav e that |h Au, Ae k i − h u, e k i| ≤ δ k u k 2 , (32) as long as A satisﬁes the RIP of order k u k 0 +1 with isometry constant δ . W e ﬁrst derive an upp er b ound on the estimation error b η lag [ k ] − η lag [ k ]: b η lag [ k ] − η lag [ k ] = e H k b η lag − e H k η lag = e H k A H w lag − e H k η lag 29 = h w lag , Ae k i − h η lag , e k i = h Aη lag , Ae k i − h η lag , e k i = h Aη lag ,s , Ae k i − h η lag ,s , e k i + h Aη lag ,t , Ae k i − h η lag ,t , e k i ≤ δ k η lag ,s k 2 + h Aη lag ,t , Ae k i − h η lag ,t , e k i = δ k η lag ,s k 2 + m − 1 X j =1 h Aη lag ,t,j , Ae k i − h η lag ,t,j , e k i ≤ δ k η lag ,s k 2 + δ k η lag ,t, 1 k 2 + δ k η lag ,t, 2 k 2 + · · · + δ k η lag ,t,m − 1 k 2 ≤ δ k η lag ,s k 2 + δ √ s k η lag ,t k ∞ + δ √ s k η lag ,t k ∞ + · · · + δ √ s k η lag ,t k ∞ = δ k η lag ,s k 2 + δ ( m − 1) √ s k η lag ,t k ∞ . Here, the ﬁfth line follo ws from the linearity of the inner pro duct, the sixth line follows from direct application of (32), the sev enth line follows from breaking η lag ,t in to m − 1 comp onen ts η lag ,t, 1 , η lag ,t, 2 , . . . , η lag ,t,m − 1 , each of sparsity s , the eighth line follows from direct application of (32), and the nin th line follows from a standard b ound relating the ` 2 norm of an s -sparse vector to its ` ∞ norm. Using similar arguments, one can deriv e a lo wer b ound for the error estimate as b η lag [ k ] − η lag [ k ] ≥ − δ k η lag ,s k 2 − δ ( m − 1) √ s k η lag ,t k ∞ . This completes the pro of of Lemma 3.1. B Pro of of Lemma 3.2 The result follo ws from the following equalities: 8 R x ncs ( t, τ ) n,q = ˆ E [ x lag ncs ( t, τ ) n,q ] = ˆ E " n Y i =1 x ( ∗ ) i i ncs ( t, τ i ) # = ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) ! n Y i =1 ∞ X m = −∞ a m p s/h ( t − mT s ) !# = ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) # ˆ E " n Y i =1 ∞ X m = −∞ a m p s/h ( t − mT s ) # = ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) # ˆ E " ∞ X m = −∞ n Y i =1 a m p s/h ( t − mT s ) # = ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) # · ˆ E " ∞ X m = −∞ n Y i =1 a m ! n Y i =1 p s/h ( t − mT s ) !# = ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) # · 8 Using tec hniques found in the FOT literature suc h as [11, 20], Lemma 3.2 also holds using the usual Dirac impulse train as the sampling op erator and one would arrive at the same conclusion. 30 " ∞ X m = −∞ E n Y i =1 a m ! n Y i =1 p s/h ( t − mT s ) !# = ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) # " ∞ X m = −∞ γ n Y i =1 p s/h ( t − mT s ) !# = γ ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) # " n Y i =1 ∞ X m = −∞ p s/h ( t − mT s ) # = γ ˆ E " n Y i =1 x ( ∗ ) i ( t − τ i ) ∞ X m = −∞ p s/h ( t − mT s ) # = γ ˆ E " n Y i =1 x ( ∗ ) i i s/h ( t, τ i ) # = γ ˆ E  x lag s/h ( t, τ )  = γ R x s/h ( t, τ ) n,q , where the fourth equality follo ws due to the indep endence of the time series, the ﬁfth equality follo ws since the s/h pulses do not ov erlap, the seven th equality follo ws since the probabilistic temp oral moment function is equal to the FOT temp oral moment function, the eigh th equalit y follo ws since E ( a n m ) = γ , and the ninth equalit y follows since the s/h pulses do not ov erlap. C Pro of of Lemma 3.3 The result follo ws from the following equalities:     b R α,c x ncs ( t, τ ) n,q    2  t = * 1 γ 2 b T 2 Z t + b T 2 − τ ∗ 0 2 t − b T 2 + τ ∗ 0 2 Z t + b T 2 − τ ∗ 0 2 t − b T 2 + τ ∗ 0 2 x lag ncs ( r , τ ) n,q  x lag ncs ( s, τ ) n,q  ∗ e − j 2 π α ( r − s ) dr ds + t = * 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 n Y i =1 x ( ∗ ) i i ncs ( r + t, τ i )  x ( ∗ ) i i ncs ( s + t, τ i )  ∗ e − j 2 π α ( r − s ) dr ds + t = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i i ncs ( r + t, τ i )  x ( ∗ ) i i ncs ( s + t, τ i )  ∗ e − j 2 π α ( r − s ) + t dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i i ncs ( r + t, τ i )  x ( ∗ ) i i ncs ( s + t, τ i )  ∗ + t e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i ) ∞ X m = −∞ a m p s/h ( r + t − mT s ) × x ( ∗ ) i ( s + t − τ i ) ∞ X m = −∞ a m p s/h ( s + t − mT s ) ! ∗ + t e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t 31 × * n Y i =1 ∞ X m = −∞ a m p s/h ( r + t − mT s ) ∞ X m = −∞ a m p s/h ( s + t − mT s ) !+ t e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t × *" ∞ X m = −∞ a n m p n s/h ( r + t − mT s ) # " ∞ X m = −∞ a n m p n s/h ( s + t − mT s ) #+ t e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t × * ∞ X m = −∞ a 2 n m p n s/h ( r + t − mT s ) p n s/h ( s + t − mT s ) + X u 6 = v a n u a n v p n s/h ( r + t − uT s ) p n s/h ( s + t − v T s ) + t e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t ×   * ∞ X m = −∞ a 2 n m p n s/h ( r + t − mT s ) p n s/h ( s + t − mT s ) + t + * X u 6 = v a n u a n v p n s/h ( r + t − uT s ) p n s/h ( s + t − v T s ) + t   × e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t ×   * ∞ X m = −∞ E  a 2 n m  p n s/h ( r + t − mT s ) p n s/h ( s + t − mT s ) + t + * X u 6 = v E ( a n u a n v ) p n s/h ( r + t − uT s ) p n s/h ( s + t − v T s ) + t   × e − j 2 π α ( r − s ) dr ds = 1 γ 2 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t ×   * ∞ X m = −∞ γ p n s/h ( r + t − mT s ) p n s/h ( s + t − mT s ) + t + * X u 6 = v γ 2 p n s/h ( r + t − uT s ) p n s/h ( s + t − v T s ) + t   × e − j 2 π α ( r − s ) dr ds = 1 b T 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 Z b T 2 − τ ∗ 0 2 − b T 2 + τ ∗ 0 2 * n Y i =1 x ( ∗ ) i ( r + t − τ i )  x ( ∗ ) i ( s + t − τ i )  ∗ + t ×   1 γ * ∞ X m = −∞ p n s/h ( r + t − mT s ) p n s/h ( s + t − mT s ) + t + * X u 6 = v p n s/h ( r + t − uT s ) p n s/h ( s + t − v T s ) + t   × e − j 2 π α ( r − s ) dr ds, (33) where the second equality follo ws from changes of v ariables, the third equality follows from using the Dominated Con vergence Theorem assuming the time series considered here are uniformly b ounded from ab o ve, the fourth equality follows since the complex exp onen tial term is indep endent of t , the 32 sixth equality follows from the indep endence of the time series and that both a m and p s/h ( t ) are real, the sev enth equality follo ws since the s/h pulses do not o verlap, the nin th equality follows from linearit y of the time-av eraging op erator, the tenth equality follows since the probabilistic temp oral momen t function is equal to the FOT temp oral moment function for b oth terms within the square brac kets, and the eleven th equality follows since E  a 2 n m  = γ and E ( a n u a n v ) = γ 2 ∀ u 6 = v . References [1] E. Cand ` es, J. Rom b erg, and T. T ao. Robust uncertain ty principles: Exact signal reconstruc- tion from highly incomplete frequency information. IEEE T r ans. Inf. The ory , 52(2):489–509, F ebruary 2006. [2] E. J. Cand` es. The Restricted Isometry Prop erty and its implications for compressed sensing. Compte R endus de l’A c ademie des Scienc es, Paris , 346:589–592, 2008. [3] A.V. Danda wate and G.B. Giannakis. Asymptotic theory of mixed time av erages and k th- order cyclic-moment and cumulan t statistics. IEEE T r ans. Inf. The ory , 41(1):216–232, January 1995. [4] M. A. Dav enp ort, P . T. Boufounos, M. B. W akin, and R. G. Baraniuk. Signal pro cessing with compressiv e measurements. IEEE J. Sel. T opics Signal Pr o c ess. , 4(2):445–460, April 2010. [5] M.A. Dav enp ort and M.B. W akin. Analysis of orthogonal matc hing pursuit using the restricted isometry prop ert y . IEEE T r ans. Inf. The ory , 56(9):4395–4401, Septem b er 2010. [6] O. A. Dobre, A. Ab di, Y. Bar-Ness, and W. Su. Higher-order cyclic cumulan ts for high order mo dulation classiﬁcation. In Pr o c. MILCOM , v olume 1, pages 112–117, Octob er 2003. [7] D. Donoho. Compressed sensing. IEEE T r ans. Inf. The ory , 52(4):1289–1306, April 2006. [8] W. A. Gardner. Intr o duction to cyclostationary signals . IEEE Press, 1994. [9] W. A. Gardner and W. A. Bro wn. F raction-of-time probabilit y for time-series that exhibit cyclostationary . Signal Pr o c ess. , 23:273–292, June 1991. [10] W. A. Gardner and C. M. Sp o oner. The cum ulant theory of cyclostationary time-series, Part I: Foundation. IEEE T r ans. Signal Pr o c ess. , 42(12):3387–3407, Decem b er 1994. [11] L. Izzo and A. Nap olitano. Sampling of generalized almost-cyclostationary signals. IEEE T r ans. Signal Pr o c ess. , 51(6):1546–1556, June 2003. [12] J. Le ´ sko w and A. Nap olitano. Conﬁdence b ound calculation and prediction in the fraction-of- time probabilit y framework. T echnical rep ort, 2002. [13] C. W. Lim and M. B. W akin. CHOCS: A framework for estimating compressiv e higher order cyclostationary statistics. In SPIE Defense Se curity Sens. , April 2012. [14] P . Marchand, J. L Lacoume, and C. Le Martret. Multiple h yp othesis mo dulation classiﬁcation based on cyclic cum ulants of diﬀeren t orders. In IEEE Int. Conf. on A c oust., Sp e e ch and Signal Pr o c ess. , volume 4, pages 2157–2160, 1998. 33 [15] P . Marc hand, C. Le Martret, and J. L Lacoume. Classiﬁcation of linear mo dulations by a com bination of diﬀerent orders cyclic cumulan ts. In IEEE Signal Pr o c essing Workshop on Higher-Or der Statistics , pages 47–51, 1997. [16] D. Needell and J. A. T ropp. CoSaMP: Iterativ e signal recov ery from incomplete and inaccurate samples. Appl. Comput. Harmonic Anal. , 26(3):301–321, Ma y 2009. [17] J. Reichert. Automatic classiﬁcation of communication signals using higher order statistics. In IEEE Int. Conf. A c oust., Sp e e ch, Signal Pr o c ess. , v olume 5, pages 221–224, Marc h 1992. [18] M. Rudelson and R. V ershynin. On sparse reconstruction from Fourier and Gaussian mea- suremen ts. Comm. Pur e Appl. Math. , 61(8):1025–1045, August 2008. [19] C. M. Sp o oner. The ory and applic ation of higher or der cyclostationarity . PhD thesis, Dept. Elec., Comput. Eng., Univ. of CA, Davis, 1992. [20] C. M. Sp ooner and W. A. Gardner. The cumulan t theory of cyclostationary time-series, Part I I: Developmen t and applications. IEEE T r ans. Signal Pr o c ess. , 42(12):3409–3429, December 1994. [21] C.M. Spo oner. Classiﬁcation of co-c hannel comm unication signals using cyclic cumulan ts. In Pr o c. Asilomar Conf. , volume 1, pages 531–536, 1995. [22] C.M. Spo oner. On the utility of sixth-order cyclic cumulan ts for rf signal classiﬁcation. In Pr o c. Asilomar Conf. , volume 1, pages 890–897, 2001. [23] C.M. Sp ooner, W.A. Bro wn, and G.K. Y eung. Automatic radio-frequency en vironment anal- ysis. In Pr o c. Asilomar Conf. , v olume 2, pages 1181–1186, 2000. [24] M. W akin, S. Beck er, E. Nak amura, M. Grant, E. Sov ero, D. Ching, J. Y o o, J. Romberg, A. Emami-Neyestanak, and E. Cand´ es. A non-uniform sampler for wideband sp ectrally-sparse en vironments. IEEE J. Emer ging Sel. T opics Cir cuits Syst. , 2(3):516–529, September 2012. 34

Technical Report: Compressive Temporal Higher Order Cyclostationary Statistics

Original Paper

Comments & Academic Discussion

Leave a Comment

Original Paper

Related Papers

Comments & Academic Discussion

Leave a Comment