Pseudospectrum using eigenvector method
[S,w] = peig(x,p)
[S,w] = peig(x,p,w)
[S,w] = peig(...,nfft)
[S,f] = peig(x,p,nfft,fs)
[S,f] = peig(x,p,f,fs)
[S,f] = peig(...,'corr')
[S,f] = peig(x,p,nfft,fs,nwin,noverlap)
[...] = peig(...,freqrange)
[...,v,e] = peig(...)
peig(...)
[S,w] = peig(x,p)
implements
the eigenvector spectral estimation method and returns S
,
the pseudospectrum estimate of the input signal x
,
and w
, a vector of normalized frequencies (in rad/sample)
at which the pseudospectrum is evaluated. The pseudospectrum is calculated
using estimates of the eigenvectors of a correlation matrix associated
with the input data x
, where x
is
specified as either:
A row or column vector representing one observation of the signal
A rectangular array for which each row of x
represents
a separate observation of the signal (for example, each row is one
output of an array of sensors, as in array processing), such that x'*x
is
an estimate of the correlation matrix
You can use the output of corrmtx
to
generate such an array x
.
You can specify the second input argument p
as
either:
A scalar integer. In this case, the signal subspace
dimension is p
.
A twoelement vector. In this case, p(2)
,
the second element of p
, represents a threshold
that is multiplied by λ_{min},
the smallest estimated eigenvalue of the signal's correlation matrix.
Eigenvalues below the threshold λ_{min} × p(2)
are assigned
to the noise subspace. In this case, p(1)
specifies
the maximum dimension of the signal subspace.
If the inputs to peig
are real sinusoids,
set the value of p
to double the number of input
signals. If the inputs are complex sinusoids, set p
equal
to the number of inputs.
The extra threshold parameter in the second entry in p
provides
you more flexibility and control in assigning the noise and signal
subspaces.
S
and w
have the same
length. In general, the length of the FFT and the values of the input x
determine
the length of the computed S
and the range of the
corresponding normalized frequencies. The following table indicates
the length of S
(and w
) and
the range of the corresponding normalized frequencies for this syntax.
S Characteristics for an FFT Length of 256 (Default)
Real/Complex Input Data  Length of S and w  Range of the Corresponding Normalized Frequencies 

Realvalued  129  [0, π] 
Complexvalued  256  [0, 2π) 
[S,w] = peig(x,p,w)
returns
the pseudospectrum in the vector S
computed at
the normalized frequencies specified in vector w
,
which has two or more elements
[S,w] = peig(...,nfft)
specifies
the integer length of the FFT nfft
used to estimate
the pseudospectrum. The default value for nfft
(entered
as an empty vector []
) is 256.
The following table indicates the length of S
and w
,
and the frequency range for w
for
this syntax.
S and Frequency Vector Characteristics
Real/Complex Input Data  nfft Even/Odd  Length of S and w  Range of w 

Realvalued  Even  (  [0, π] 
Realvalued  Odd  (  [0, π) 
Complexvalued  Even or odd 
 [0, 2π) 
[S,f] = peig(x,p,nfft,fs)
returns
the pseudospectrum in the vector S
evaluated at
the corresponding vector of frequencies f
(in Hz).
You supply the sampling frequency fs
in Hz. If
you specify fs
with the empty vector []
,
the sampling frequency defaults to 1 Hz.
The frequency range for f
depends on nfft
, fs
,
and the values of the input x
.
The length of S
(and f
) is the
same as in the S and Frequency Vector Characteristics above. The following
table indicates the frequency range for f
for this
syntax.
S and Frequency Vector Characteristics with fs Specified
Real/Complex Input Data  nfft Even/Odd  Range of f 

Realvalued  Even  [0, 
Realvalued  Odd  [0, 
Complexvalued  Even or odd  [0, 
[S,f] = peig(x,p,f,fs)
returns
the pseudospectrum in the vector S
computed at
the frequencies specified in vector f
, which has
two or more elements
[S,f] = peig(...,'corr')
forces
the input argument x
to be interpreted as a correlation
matrix rather than matrix of signal data. For this syntax x
must
be a square matrix, and all of its eigenvalues must be nonnegative.
[S,f] = peig(x,p,nfft,fs,nwin,noverlap)
allows
you to specify nwin
, a scalar integer indicating
a rectangular window length, or a realvalued vector specifying window
coefficients. Use the scalar integer noverlap
in
conjunction with nwin
to specify the number of
input sample points by which successive windows overlap. noverlap
is
not used if x
is a matrix. The default value for nwin
is
2 × p(1)
and noverlap
is nwin
– 1.
With this syntax, the input data x
is segmented
and windowed before the matrix used to estimate the correlation matrix
eigenvalues is formulated. The segmentation of the data depends on nwin
, noverlap
,
and the form of x
. Comments on the resulting windowed
segments are described in the following table.
Windowed Data Depending on x and nwin
Input data x  Form of nwin  Windowed Data 

Data vector  Scalar  Length is 
Data vector  Vector of coefficients  Length is 
Data matrix  Scalar  Data is not windowed. 
Data matrix  Vector of coefficients 

See the table, Eigenvector Length Depending on Input Data and Syntax, for related information on this syntax.
The arguments nwin
and noverlap
are
ignored when you include 'corr'
in the syntax.
[...] = peig(...,freqrange)
specifies
the range of frequency values to include in f
or w
. freqrange
can
be either:
'half'
— returns half the
spectrum for a real input signal, x
. If nfft
is
even, then S
has length nfft
/2 + 1 and is computed over the interval
[0,π]. If nfft
is odd,
the length of S
is (nfft
+ 1)/2 and the frequency interval is
[0,π). When your specify fs
,
the intervals are [0,fs
/2) and [0,fs
/2]
for even and odd nfft
, respectively.
'whole'
— returns the whole
spectrum for either real or complex input, x
. In
this case, S
has length nfft
and
is computed over the interval [0,2π). When
you specify fs
, the frequency interval is [0,fs
).
'centered'
— returns the
centered whole spectrum for either real or complex input, x
.
In this case, S
has length nfft
and
is computed over the interval (–π,π]
for even nfft
and (–π,π)
for odd nfft
. When you specify fs
,
the frequency intervals are (–fs
/2,fs
/2]
and (–fs
/2,fs
/2) for
even and odd nfft
, respectively.
You can put the arguments freqrange
or 'corr'
anywhere
in the input argument list after p
.
[...,v,e] = peig(...)
returns
the matrix v
of noise eigenvectors, along with
the associated eigenvalues in the vector e
. The
columns of v
span the noise subspace of dimension size(v,2)
.
The dimension of the signal subspace is size(v,1)size(v,2)
.
For this syntax, e
is a vector of estimated eigenvalues
of the correlation matrix.
peig(...)
with no output
arguments plots the pseudospectrum in the current figure window.
In the process of estimating the pseudospectrum, peig
computes
the noise and signal subspaces from the estimated eigenvectors v_{j} and
eigenvalues λ_{j} of
the signal's correlation matrix. The smallest of these eigenvalues
is used in conjunction with the threshold parameter p(2)
to
affect the dimension of the noise subspace in some cases.
The length n of the eigenvectors computed
by peig
is the sum of the dimensions of the signal
and noise subspaces. This eigenvector length depends on your input
(signal data or correlation matrix) and the syntax you use.
The following table summarizes the dependency of the eigenvector length on the input argument.
Eigenvector Length Depending on Input Data and Syntax
Form of Input Data x  Comments on the Syntax  Length n of Eigenvectors 

Row or column vector 


Row or column vector 


Row or column vector 
 2 × 
lbym matrix  If  m 
mbym nonnegative definite matrix 
 m 
You should specify nwin
> p(1)
or length(nwin)
> p(1)
if you want p(2)
> 1
to have any
effect.
The eigenvector method estimates the pseudospectrum from a signal
or a correlation matrix using a weighted version of the MUSIC algorithm
derived from Schmidt's eigenspace analysis method [1] [2].
The algorithm performs eigenspace analysis of the signal's correlation
matrix in order to estimate the signal's frequency content. The eigenvalues
and eigenvectors of the signal's correlation matrix are estimated
using svd
if you don't supply
the correlation matrix. This algorithm is particularly suitable for
signals that are the sum of sinusoids with additive white Gaussian
noise.
The eigenvector method produces a pseudospectrum estimate given by
$${P}_{ev}(f)=\frac{1}{{\displaystyle \sum _{k=p+1}^{N}{\text{v}}_{k}^{H}e(f){}^{2}/{\lambda}_{k}}}$$
where N is the dimension of the eigenvectors and v_{k }is the kth eigenvector of the correlation matrix of the input signal. The integer p is the dimension of the signal subspace, so the eigenvectors v_{k} used in the sum correspond to the smallest eigenvalues λ_{k} of the correlation matrix. The eigenvectors used span the noise subspace. The vector e(f) consists of complex exponentials, so the inner product
$${\text{v}}_{k}^{H}e(f)$$
amounts to a Fourier transform. This is used for computation of the pseudospectrum. The FFT is computed for each v_{k} and then the squared magnitudes are summed and scaled.
[1] Marple, S. Lawrence. Digital Spectral Analysis. Englewood Cliffs, NJ: PrenticeHall, 1987, pp. 373–378.
[2] Schmidt, R. O. “Multiple Emitter Location and Signal Parameter Estimation.” IEEE^{®} Transactions on Antennas and Propagation. Vol. AP34, March, 1986, pp. 276–280.
[3] Stoica, Petre, and Randolph L. Moses. Spectral Analysis of Signals. Upper Saddle River, NJ: Prentice Hall, 2005.