Estimation of signal parameters via rotational invariance techniques

Last updated June 27, 2024

Estimation of signal parameters via rotational invariant techniques (ESPRIT), is a technique to determine the parameters of a mixture of sinusoids in background noise. This technique was first proposed for frequency estimation.^[1] However, with the introduction of phased-array systems in everyday technology, it is also used for angle of arrival estimations.^[2]

One-dimensional ESPRIT

At instance ${\textstyle t}$ , the ${\textstyle M}$ (complex -valued) output signals (measurements) ${\textstyle y_{m}[t]}$ , ${\textstyle m=1,\ldots ,M}$ , of the system are related to the ${\textstyle K}$ (complex -valued) input signals ${\textstyle x_{k}[t]}$ , ${\textstyle k=1,\ldots ,K}$ , as $y_{m}[t]=\sum _{k=1}^{K}a_{m,k}x_{k}[t]+n_{m}[t],$ where ${\textstyle n_{m}[t]}$ denotes the noise added by the system. The one-dimensional form of ESPRIT can be applied if the weights have the form ${\textstyle a_{m,k}=e^{-j(m-1)\omega _{k}}}$ , whose phases are integer multiples of some radial frequency $\omega _{k}$ . This frequency only depends on the index of the system's input, i.e., ${\textstyle k}$ . The goal of ESPRIT is to estimate $\omega _{k}$ 's, given the outputs ${\textstyle y_{m}[t]}$ and the number of input signals, ${\textstyle K}$ . Since the radial frequencies are the actual objectives, ${\textstyle a_{m,k}}$ is denoted as ${\textstyle a_{m}(\omega _{k})}$ .

Collating the weights ${\textstyle a_{m}(\omega _{k})}$ as $\mathbf {a} (\omega _{k})=[\,1\,\ e^{-j\omega _{k}}\,\ e^{-j2\omega _{k}}\,\ ...\,\ e^{-j(M-1)\omega _{k}}\,]^{\top }$ and the ${\textstyle M}$ output signals at instance ${\textstyle t}$ as $\mathbf {y} [t]=[\,y_{1}[t]\,\ y_{2}[t]\,\ ...\,\ y_{M}[t]\,]^{\top }$ , $\mathbf {y} [t]=\sum _{k=1}^{K}\mathbf {a} (\omega _{k})x_{k}[t]+\mathbf {n} [t],$ where $\mathbf {n} [t]=[\,n_{1}[t]\,\ n_{2}[t]\,\ ...\,\ n_{M}[t]\,]^{\top }$ . Further, when the weight vectors $\mathbf {a} (\omega _{1}),\ \mathbf {a} (\omega _{2}),\ ...,\ \mathbf {a} (\omega _{K})$ are put into a Vandermonde matrix $\mathbf {A} =[\,\mathbf {a} (\omega _{1})\,\ \mathbf {a} (\omega _{2})\,\ ...\,\ \mathbf {a} (\omega _{K})\,]$ , and the $K$ inputs at instance ${\textstyle t}$ into a vector $\mathbf {x} [t]=[\,x_{1}[t]\,\ ...\,\ x_{k}[t]\,]^{\top }$ , we can write $\mathbf {y} [t]=\mathbf {A} \,\mathbf {x} [t]+\mathbf {n} [t].$ With several measurements at instances ${\textstyle t=1,2,\dots ,T}$ and the notations ${\textstyle \mathbf {Y} =[\,\mathbf {y} [1]\,\ \mathbf {y} [2]\,\ \dots \,\ \mathbf {y} [T]\,]}$ , ${\textstyle \mathbf {X} =[\,\mathbf {x} [1]\,\ \mathbf {x} [2]\,\ \dots \,\ \mathbf {x} [T]\,]}$ and ${\textstyle \mathbf {N} =[\,\mathbf {n} [1]\,\ \mathbf {n} [2]\,\ \dots \,\ \mathbf {n} [T]\,]}$ , the model equation becomes $\mathbf {Y} =\mathbf {A} \mathbf {X} +\mathbf {N} .$

Dividing into virtual sub-arrays

The weight vector ${\textstyle \mathbf {a} (\omega _{k})}$ has the property that adjacent entries are related. $[\mathbf {a} (\omega _{k})]_{m+1}=e^{-j\omega _{k}}[\mathbf {a} (\omega _{k})]_{m}$ For the whole vector $\mathbf {a} (\omega _{k})$ , the equation introduces two selection matrices $\mathbf {J} _{1}$ and $\mathbf {J} _{2}$ : $\mathbf {J} _{1}=[\mathbf {I} _{M-1}$ and $\mathbf {J} _{2}=[\mathbf {0} \quad \mathbf {I} _{M-1}]$ . Here, $\mathbf {I} _{M-1}$ is an identity matrix of size ${\textstyle (M-1)}$ and $\mathbf {0}$ is a vector of zero.

The vectors $\mathbf {J} _{1}\mathbf {a} (\omega _{k})$ $[\mathbf {J} _{2}\mathbf {a} (\omega _{k})]$ contains all elements of $\mathbf {a} (\omega _{k})$ except the last [first] one. Thus, $\mathbf {J} _{2}\mathbf {a} (\omega _{k})=e^{-j\omega _{k}}\mathbf {J} _{1}\mathbf {a} (\omega _{k})$ and $\mathbf {J} _{2}\mathbf {A} =\mathbf {J} _{1}\mathbf {A} \mathbf {H} ,\quad {\text{where}}\quad {\mathbf {H} }:={\begin{bmatrix}e^{-j\omega _{1}}&\\&e^{-j\omega _{2}}\\&&\ddots \\&&&e^{-j\omega _{K}}\end{bmatrix}}.$ The above relation is the first major observation required for ESPRIT. The second major observation concerns the signal subspace that can be computed from the output signals.

Signal subspace

The singular value decomposition (SVD) of ${\textstyle \mathbf {Y} }$ is given as $\mathbf {Y} =\mathbf {U} \mathbf {\Sigma } \mathbf {V} ^{\dagger }$ where ${\textstyle \mathbf {U} \in \mathbb {C} ^{M\times M}}$ and ${\textstyle \mathbf {V} \in \mathbb {C} ^{T\times T}}$ are unitary matrices and ${\textstyle \mathbf {\Sigma } }$ is a diagonal matrix of size ${\textstyle M\times T}$ , that holds the singular values from the largest (top left) in descending order. The operator ${\textstyle \dagger }$ denotes the complex-conjugate transpose (Hermitian transpose).

Let us assume that ${\textstyle T\geq M}$ . Notice that we have ${\textstyle K}$ input signals. If there was no noise, there would only be ${\textstyle K}$ non-zero singular values. We assume that the ${\textstyle K}$ largest singular values stem from these input signals and other singular values are presumed to stem from noise. The matrices in SVD of ${\textstyle \mathbf {Y} }$ can be partitioned into submatrices, where some submatrices correspond to the signal subspace and some correspond to the noise subspace. ${\begin{aligned}\mathbf {U} ={\begin{bmatrix}\mathbf {U} _{\mathrm {S} }&\mathbf {U} _{\mathrm {N} }\end{bmatrix}},&&\mathbf {\Sigma } ={\begin{bmatrix}\mathbf {\Sigma } _{\mathrm {S} }&\mathbf {0} &\mathbf {0} \\\mathbf {0} &\mathbf {\Sigma } _{\mathrm {N} }&\mathbf {0} \end{bmatrix}},&&\mathbf {V} ={\begin{bmatrix}\mathbf {V} _{\mathrm {S} }&\mathbf {V} _{\mathrm {N} }&\mathbf {V} _{\mathrm {0} }\end{bmatrix}}\end{aligned}},$ where ${\textstyle \mathbf {U} _{\mathrm {S} }\in \mathbb {C} ^{M\times K}}$ and ${\textstyle \mathbf {V} _{\mathrm {S} }\in \mathbb {C} ^{N\times K}}$ contain the first ${\textstyle K}$ columns of ${\textstyle \mathbf {U} }$ and ${\textstyle \mathbf {V} }$ , respectively and ${\textstyle \mathbf {\Sigma } _{\mathrm {S} }\in \mathbb {C} ^{K\times K}}$ is a diagonal matrix comprising the ${\textstyle K}$ largest singular values.

Thus, The SVD can be written as $\mathbf {Y} =\mathbf {U} _{\mathrm {S} }\mathbf {\Sigma } _{\mathrm {S} }\mathbf {V} _{\mathrm {S} }^{\dagger }+\mathbf {U} _{\mathrm {N} }\mathbf {\Sigma } _{\mathrm {N} }\mathbf {V} _{\mathrm {N} }^{\dagger },$ where ${\textstyle \mathbf {U} _{\mathrm {S} }}$ , ⁣ ${\textstyle \mathbf {V} _{\mathrm {S} }}$ , and ${\textstyle \mathbf {\Sigma } _{\mathrm {S} }}$ represent the contribution of the input signal ${\textstyle x_{k}[t]}$ to ${\textstyle \mathbf {Y} }$ . We term ${\textstyle \mathbf {U} _{\mathrm {S} }}$ the signal subspace. In contrast, ${\textstyle \mathbf {U} _{\mathrm {N} }}$ , ${\textstyle \mathbf {V} _{\mathrm {N} }}$ , and ${\textstyle \mathbf {\Sigma } _{\mathrm {N} }}$ represent the contribution of noise ${\textstyle n_{m}[t]}$ to ${\textstyle \mathbf {Y} }$ .

Hence, from the system model, we can write ${\textstyle \mathbf {A} \mathbf {X} =\mathbf {U} _{\mathrm {S} }\mathbf {\Sigma } _{\mathrm {S} }\mathbf {V} _{\mathrm {S} }^{\dagger }}$ and ${\textstyle \mathbf {N} =\mathbf {U} _{\mathrm {N} }\mathbf {\Sigma } _{\mathrm {N} }\mathbf {V} _{\mathrm {N} }^{\dagger }}$ . Also, from the former, we can write $\mathbf {U} _{\mathrm {S} }=\mathbf {A} \mathbf {F} ,$ where ${\mathbf {F} }=\mathbf {X} \mathbf {V} _{\mathrm {S} }\mathbf {\Sigma } _{\mathrm {S} }^{-1}$ . In the sequel, it is only important that there exist such an invertible matrix ${\textstyle \mathbf {F} }$ and its actual content will not be important.

Note: The signal subspace can also be extracted from the spectral decomposition of the auto-correlation matrix of the measurements, which is estimated as $\mathbf {R} _{\mathrm {YY} }={\frac {1}{T}}\sum _{t=1}^{T}\mathbf {y} [t]\mathbf {y} [t]^{\dagger }={\frac {1}{T}}\mathbf {Y} \mathbf {Y} ^{\dagger }={\frac {1}{T}}\mathbf {U} {\mathbf {\Sigma } \mathbf {\Sigma } ^{\dagger }}\mathbf {U} ^{\dagger }={\frac {1}{T}}\mathbf {U} _{\mathrm {S} }\mathbf {\Sigma } _{\mathrm {S} }^{2}\mathbf {U} _{\mathrm {S} }^{\dagger }+{\frac {1}{T}}\mathbf {U} _{\mathrm {N} }\mathbf {\Sigma } _{\mathrm {N} }^{2}\mathbf {U} _{\mathrm {N} }^{\dagger }.$

Estimation of radial frequencies

We have established two expressions so far: ${\textstyle \mathbf {J} _{2}\mathbf {A} =\mathbf {J} _{1}\mathbf {A} \mathbf {H} }$ and ${\textstyle \mathbf {U} _{\mathrm {S} }=\mathbf {A} \mathbf {F} }$ . Now, ${\begin{aligned}\mathbf {J} _{2}\mathbf {A} =\mathbf {J} _{1}\mathbf {A} \mathbf {H} \implies \mathbf {J} _{2}(\mathbf {U} _{\mathrm {S} }\mathbf {F} ^{-1})=\mathbf {J} _{1}(\mathbf {U} _{\mathrm {S} }\mathbf {F} ^{-1})\mathbf {H} \implies \mathbf {S} _{2}=\mathbf {S} _{1}\mathbf {P} ,\end{aligned}}$ where $\mathbf {S} _{1}=\mathbf {J} _{1}\mathbf {U} _{\mathrm {S} }$ and $\mathbf {S} _{2}=\mathbf {J} _{2}\mathbf {U} _{\mathrm {S} }$ denote the truncated signal sub spaces, and $\mathbf {P} =\mathbf {F} ^{-1}\mathbf {H} \mathbf {F} .$ The above equation has the form of an eigenvalue decomposition, and the phases of eigenvalues in the diagonal matrix $\mathbf {H}$ are used to estimate the radial frequencies.

Thus, after solving for $\mathbf {P}$ in the relation $\mathbf {S} _{2}=\mathbf {S} _{1}\mathbf {P}$ , we would find the eigenvalues ${\textstyle \lambda _{1},\ \lambda _{2},\ \ldots ,\ \lambda _{K}}$ of $\mathbf {P}$ , where $\lambda _{k}=\alpha _{k}\mathrm {e} ^{j\omega _{k}}$ , and the radial frequencies $\omega _{1},\ \omega _{2},\ \ldots ,\ \omega _{K}$ are estimated as the phases (argument) of the eigenvalues.

Remark: In general, $\mathbf {S} _{1}$ is not invertible. One can use the least squares estimate ${\mathbf {P} }=(\mathbf {S} _{1}^{\dagger }{\mathbf {S} _{1}})^{-1}\mathbf {S} _{1}^{\dagger }{\mathbf {S} _{2}}$ . An alternative would be the total least squares estimate.

Algorithm summary

Input: Measurements ${\textstyle \mathbf {Y}$ , the number of input signals ${\textstyle K}$ (estimate if not already known).

Compute the singular value decomposition (SVD) of ${\textstyle \mathbf {Y} =\mathbf {U} \mathbf {\Sigma } \mathbf {V} ^{\dagger }}$ and extract the signal subspace ${\textstyle \mathbf {U} _{\mathrm {S} }\in \mathbb {C} ^{M\times K}}$ as the first ${\textstyle K}$ columns of ${\textstyle \mathbf {U} }$ .
Compute ${\textstyle \mathbf {S} _{1}=\mathbf {J} _{1}\mathbf {U} _{\mathrm {S} }}$ and ${\textstyle \mathbf {S} _{2}=\mathbf {J} _{2}\mathbf {U} _{\mathrm {S} }}$ , where $\mathbf {J} _{1}=[\mathbf {I} _{M-1}\quad \mathbf {0} ]$ and $\mathbf {J} _{2}=[\mathbf {0} \quad \mathbf {I} _{M-1}]$ .
Solve for ${\textstyle \mathbf {P} }$ in ${\textstyle \mathbf {S} _{2}=\mathbf {S} _{1}\mathbf {P} }$ (see the remark above).
Compute the eigenvalues ${\textstyle \lambda _{1},\lambda _{2},\ldots ,\lambda _{K}}$ of ${\textstyle \mathbf {P} }$ .
The phases of the eigenvalues ${\textstyle \lambda _{k}=\alpha _{k}\mathrm {e} ^{j\omega _{k}}}$ provide the radial frequencies ${\textstyle \omega _{k}}$ , i.e., ${\textstyle \omega _{k}=\arg \lambda _{k}}$

Notes

Choice of selection matrices

In the derivation above, the selection matrices $\mathbf {J} _{1}=[\mathbf {I} _{M-1}\quad \mathbf {0} ]$ and $\mathbf {J} _{2}=[\mathbf {0} \quad \mathbf {I} _{M-1}]$ were used. However, any appropriate matrices ${\textstyle \mathbf {J} _{1}}$ and $\mathbf {J} _{2}$ may be used as long as the rotational invariance ${\textstyle (}$ i.e., $\mathbf {J} _{2}\mathbf {a} (\omega _{k})=e^{-j\omega _{k}}\mathbf {J} _{1}\mathbf {a} (\omega _{k})$ ${\textstyle )}$ , or some generalization of it (see below) holds; accordingly, the matrices ${\textstyle \mathbf {J} _{1}\mathbf {A} }$ and ${\textstyle \mathbf {J} _{2}\mathbf {A} }$ may contain any rows of ${\textstyle \mathbf {A} }$ .

Generalized rotational invariance

The rotational invariance used in the derivation may be generalized. So far, the matrix $\mathbf {H}$ has been defined to be a diagonal matrix that stores the sought-after complex exponentials on its main diagonal. However, $\mathbf {H}$ may also exhibit some other structure.^[3] For instance, it may be an upper triangular matrix. In this case, ${\textstyle \mathbf {P} =\mathbf {F} ^{-1}\mathbf {H} \mathbf {F} }$ constitutes a triangularization of $\mathbf {P}$ .

Related Research Articles

In physics, the cross section is a measure of the probability that a specific process will take place in a collision of two particles. For example, the Rutherford cross-section is a measure of probability that an alpha particle will be deflected by a given angle during an interaction with an atomic nucleus. Cross section is typically denoted $σ$ (sigma) and is expressed in units of area, more specifically in barns. In a way, it can be thought of as the size of the object that the excitation must hit in order for the process to occur, but more exactly, it is a parameter of a stochastic process.

In particle physics, the Dirac equation is a relativistic wave equation derived by British physicist Paul Dirac in 1928. In its free form, or including electromagnetic interactions, it describes all spin-1/2 massive particles, called "Dirac particles", such as electrons and quarks for which parity is a symmetry. It is consistent with both the principles of quantum mechanics and the theory of special relativity, and was the first theory to account fully for special relativity in the context of quantum mechanics. It was validated by accounting for the fine structure of the hydrogen spectrum in a completely rigorous way. It has become vital in the building of the Standard Model.

<span class="mw-page-title-main">Navier–Stokes equations</span> Equations describing the motion of viscous fluid substances

The Navier–Stokes equations are partial differential equations which describe the motion of viscous fluid substances. They were named after French engineer and physicist Claude-Louis Navier and the Irish physicist and mathematician George Gabriel Stokes. They were developed over several decades of progressively building the theories, from 1822 (Navier) to 1842–1850 (Stokes).

In mathematics, the exterior algebra or Grassmann algebra of a vector space $is an associative algebra that contains which has a product, called exterior product or wedge product and denoted with, such that for every vector in The exterior algebra is named after Hermann Grassmann, and the names of the product come from the "wedge" symbol and the fact that the product of two elements of are "outside"$

<span class="mw-page-title-main">Hooke's law</span> Physical law: force needed to deform a spring scales linearly with distance

In physics, Hooke's law is an empirical law which states that the force needed to extend or compress a spring by some distance scales linearly with respect to that distance—that is, $F s = kx$ , where $k$ is a constant factor characteristic of the spring, and $x$ is small compared to the total possible deformation of the spring. The law is named after 17th-century British physicist Robert Hooke. He first stated the law in 1676 as a Latin anagram. He published the solution of his anagram in 1678 as: ut tensio, sic vis. Hooke states in the 1678 work that he was aware of the law since 1660.

In physics and mathematics, the Lorentz group is the group of all Lorentz transformations of Minkowski spacetime, the classical and quantum setting for all (non-gravitational) physical phenomena. The Lorentz group is named for the Dutch physicist Hendrik Lorentz.

In quantum field theory, the Dirac spinor is the spinor that describes all known fundamental particles that are fermions, with the possible exception of neutrinos. It appears in the plane-wave solution to the Dirac equation, and is a certain combination of two Weyl spinors, specifically, a bispinor that transforms "spinorially" under the action of the Lorentz group.

Stellar dynamics is the branch of astrophysics which describes in a statistical way the collective motions of stars subject to their mutual gravity. The essential difference from celestial mechanics is that the number of body

In physics, the Hamilton–Jacobi equation, named after William Rowan Hamilton and Carl Gustav Jacob Jacobi, is an alternative formulation of classical mechanics, equivalent to other formulations such as Newton's laws of motion, Lagrangian mechanics and Hamiltonian mechanics.

The Rayleigh–Ritz method is a direct numerical method of approximating eigenvalues, originated in the context of solving physical boundary value problems and named after Lord Rayleigh and Walther Ritz.

In physics, the Majorana equation is a relativistic wave equation. It is named after the Italian physicist Ettore Majorana, who proposed it in 1937 as a means of describing fermions that are their own antiparticle. Particles corresponding to this equation are termed Majorana particles, although that term now has a more expansive meaning, referring to any fermionic particle that is its own anti-particle.

In numerical linear algebra, the method of successive over-relaxation (SOR) is a variant of the Gauss–Seidel method for solving a linear system of equations, resulting in faster convergence. A similar method can be used for any slowly converging iterative process.

<span class="mw-page-title-main">Prony's method</span> Method to estimate the components of a signal

Prony analysis was developed by Gaspard Riche de Prony in 1795. However, practical use of the method awaited the digital computer. Similar to the Fourier transform, Prony's method extracts valuable information from a uniformly sampled signal and builds a series of damped complex exponentials or damped sinusoids. This allows the estimation of frequency, amplitude, phase and damping components of a signal.

The Jaynes–Cummings model is a theoretical model in quantum optics. It describes the system of a two-level atom interacting with a quantized mode of an optical cavity, with or without the presence of light. It was originally developed to study the interaction of atoms with the quantized electromagnetic field in order to investigate the phenomena of spontaneous emission and absorption of photons in a cavity.

Stokes' theorem, also known as the Kelvin–Stokes theorem after Lord Kelvin and George Stokes, the fundamental theorem for curls or simply the curl theorem, is a theorem in vector calculus on $. Given a vector field, the theorem relates the integral of the curl of the vector field over some surface, to the line integral of the vector field around the boundary of the surface. The classical theorem of Stokes can be stated in one sentence: The line integral of a vector field over a loop is equal to the surface integral of its curl over the enclosed surface. It is illustrated in the figure, where the direction of positive circulation of the bounding contour \partialΣ, and the direction n of positive flux through the surface Σ, are related by a right-hand-rule. For the right hand the fingers circulate along \partialΣ and the thumb is directed along n .$

In the geometry of numbers, the Klein polyhedron, named after Felix Klein, is used to generalize the concept of continued fractions to higher dimensions.

In mathematics, the Lindström–Gessel–Viennot lemma provides a way to count the number of tuples of non-intersecting lattice paths, or, more generally, paths on a directed graph. It was proved by Gessel–Viennot in 1985, based on previous work of Lindström published in 1973.

In physics, particularly in quantum field theory, the Weyl equation is a relativistic wave equation for describing massless spin-1/2 particles called Weyl fermions. The equation is named after Hermann Weyl. The Weyl fermions are one of the three possible types of elementary fermions, the other two being the Dirac and the Majorana fermions.

Symmetries in quantum mechanics describe features of spacetime and particles which are unchanged under some transformation, in the context of quantum mechanics, relativistic quantum mechanics and quantum field theory, and with applications in the mathematical formulation of the standard model and condensed matter physics. In general, symmetry in physics, invariance, and conservation laws, are fundamentally important constraints for formulating physical theories and models. In practice, they are powerful methods for solving problems and predicting what can happen. While conservation laws do not always give the answer to the problem directly, they form the correct constraints and the first steps to solving a multitude of problems. In application, understanding symmetries can also provide insights on the eigenstates that can be expected. For example, the existence of degenerate states can be inferred by the presence of non commuting symmetry operators or that the non degenerate states are also eigenvectors of symmetry operators.

Generalized pencil-of-function method (GPOF), also known as matrix pencil method, is a signal processing technique for estimating a signal or extracting information with complex exponentials. Being similar to Prony and original pencil-of-function methods, it is generally preferred to those for its robustness and computational efficiency.

References

↑ Paulraj, A.; Roy, R.; Kailath, T. (1985), "Estimation Of Signal Parameters Via Rotational Invariance Techniques - Esprit", Nineteenth Asilomar Conference on Circuits, Systems and Computers, pp. 83–89, doi:10.1109/ACSSC.1985.671426, ISBN 978-0-8186-0729-5, S2CID 2293566
↑ Volodymyr Vasylyshyn. The direction of arrival estimation using ESPRIT with sparse arrays.// Proc. 2009 European Radar Conference (EuRAD). – 30 Sept.-2 Oct. 2009. - Pp. 246 - 249. -
↑ Hu, Anzhong; Lv, Tiejun; Gao, Hui; Zhang, Zhang; Yang, Shaoshi (2014). "An ESPRIT-Based Approach for 2-D Localization of Incoherently Distributed Sources in Massive MIMO Systems". IEEE Journal of Selected Topics in Signal Processing. 8 (5): 996–1011. arXiv: 1403.5352 . Bibcode:2014ISTSP...8..996H. doi:10.1109/JSTSP.2014.2313409. ISSN 1932-4553. S2CID 11664051.