pxx = pyulear(x,order)
returns the power spectral density estimate, pxx of
a discrete-time signal vector, x, using the Yule-Walker
method. pxx is the distribution of power per
unit frequency. The frequency is expressed in units of radians/sample.
order is the order of the autoregressive (AR)
model used to produce the PSD estimate. pyulear uses
a default DFT length of 256.

pxx = pyulear(x,order,nfft) uses nfft points
in the discrete Fourier transform (DFT). For real x, pxx has
length (nfft/2+1) if nfft is
even, and (nfft+1)/2 if nfft is
odd. For complex–valued x, pxx always
has length nfft. If empty, the default nfft is
256.

[pxx,w]
= pyulear(___) returns the vector of normalized
angular frequencies, w, at which the PSD is
estimated. w has units of radians/sample. For
real–valued signals, w spans the interval
[0,π] when nfft is even and [0,π)
when nfft is odd. For complex–valued
signals, w always spans the interval [0,2π).

[pxx,f]
= pyulear(___,fs) returns
a frequency vector, f, in cycles per unit time.
The sampling frequency, fs, is the number of
samples per unit time. If the unit of time is seconds, then f is
in cycles/sec (Hz). For real–valued signals, f spans
the interval [0,fs/2] when nfft is
even and [0,fs/2) when nfft is
odd. For complex-valued signals, f spans the
interval [0,fs).

[pxx,w]
= pyulear(x,order,w) returns
the two-sided AR PSD estimates at the normalized frequencies specified
in the vector, w. The vector, w,
must contain at least 2 elements.

[pxx,f]
= pyulear(x,order,f,fs) returns
the two-sided AR PSD estimates at the frequencies specified in the
vector, f. The vector, f,
must contain at least 2 elements. The frequencies in f are
in cycles per unit time. The sampling frequency, fs,
is the number of samples per unit time. If the unit of time is seconds,
then f is in cycles/sec (Hz).

[___] = pyulear(x,order,___,freqrange) returns
the AR PSD estimate over the frequency range specified by freqrange.
Valid options for freqrange are: 'onesided', 'twosided',
or 'centered'.

[pxx,f,pxxc]
= pyulear(___,'ConfidenceLevel',probability) returns
the probabilityx100% confidence intervals for
the PSD estimate in pxxc.

pyulear(___) with no output
arguments plots the AR PSD estimate in dB per unit frequency in the
current figure window.

Create a realization of an AR(4) wide-sense
stationary random process. Estimate the PSD using the Yule-Walker
method. Compare the PSD estimate based on a single realization to
the true PSD of the random process.

Create an AR(4) system function. Obtain the frequency
response and plot the PSD of the system.

A = [1 -2.7607 3.8106 -2.6535 0.9238];
[H,F] = freqz(1,A,[],1);
plot(F,20*log10(abs(H)),'b','linewidth',2);
xlabel('Hz'); ylabel('dB/Hz');
title('True Power Spectral Density of AR(4) System Function')

Create a realization of the AR(4) random process. Set
the random number generator to the default settings for reproducible
results. The realization is 1000 samples in length. Assume a sampling
frequency of 1. Use pyulear to estimate the PSD
for an 4-th order process. Compare the PSD estimate with the true
PSD.

rng default;
x = randn(1000,1);
y = filter(1,A,x);
[Pxx,F] = pyulear(y,4,1024,1);
hold on;
plot(F,10*log10(Pxx),'r'); hold on;
legend('True Power Spectral Density','PSD Estimate')

Number of DFT points, specified as a positive integer. For a
real-valued input signal, x, the PSD estimate, pxx has
length (nfft/2+1) if nfft is
even, and (nfft+1)/2 if nfft is
odd. For a complex-valued input signal,x, the
PSD estimate always has length nfft. If nfft is
specified as empty, the default nfft is used.

Sampling frequency specified as a positive scalar. The sampling
frequency is the number of samples per unit time. If the unit of time
is seconds, the sampling frequency has units of hertz.

Normalized frequencies for Goertzel algorithm, specified as
a row or column vector with at least 2 elements. Normalized frequencies
are in radians/sample.

Cyclical frequencies for Goertzel algorithm, specified as a
row or column vector with at least 2 elements. The frequencies are
in cycles per unit time. The unit time is specified by the sampling
frequency, fs. If fs has
units of samples/second, then f has units of
Hz.

Frequency range for the PSD estimate, specified as a one of 'onesided', 'twosided',
or 'centered'. The default is 'onesided' for
real-valued signals and 'twosided' for complex-valued
signals. The frequency ranges corresponding to each option are

'onesided' — returns the
one-sided PSD estimate of a real-valued input signal, x.
If nfft is even, pxx will
have length nfft/2+1 and is computed over the
interval [0,π] radians/sample. If nfft is
odd, the length of pxx is (nfft+1)/2
and the interval is [0,π) radians/sample. When fs is
optionally specified, the corresponding intervals are [0,fs/2]
cycles/unit time and [0,fs/2) cycles/unit time
for even and odd length nfft respectively.

'twosided' — returns the
two-sided PSD estimate for either the real-valued or complex-valued
input, x. In this case, pxx has
length nfft and is computed over the interval
[0,2π) radians/sample. When fs is optionally
specified, the interval is [0,fs) cycles/unit
time.

'centered' — returns the
centered two-sided PSD estimate for either the real-valued or complex-valued
input, x. In this case, pxx has
length nfft and is computed over the interval
(-π,π] radians/sample for even length nfft and
(-π,π) radians/sample for odd length nfft.
When fs is optionally specified, the corresponding
intervals are (-fs/2, fs/2]
cycles/unit time and (-fs/2, fs/2)
cycles/unit time for even and odd length nfft respectively.

Coverage probability for the true PSD, specified as a scalar
in the range (0,1). The output, pxxc, contains
the lower and upper bounds of the probabilityx100%
interval estimate for the true PSD.

PSD estimate, specified as a real-valued, nonnegative column
vector. The units of the PSD estimate are in squared magnitude units
of the time series data per unit frequency. For example, if the input
data is in volts, the PSD estimate is in units of squared volts per
unit frequency. For a time series in volts, if you assume a resistance
of 1 ohm and specify the sampling frequency in Hz, the PSD estimate
is in watts/Hz.

Normalized frequencies, specified as a real-valued column vector.
If pxx is a one-sided PSD estimate, w spans
the interval [0,π] if nfft is even and
[0,π) if nfft is odd. If pxx is
a two-sided PSD estimate, w spans the interval
[0,2π). For a DC-centered PSD estimate, f spans
the interval (-π,π] radians/sample for even length nfft and
(-π,π) radians/sample for odd length nfft.

Cyclical frequencies, specified as a real-valued column vector.
For a one-sided PSD estimate, f spans the interval
[0,fs/2] when nfft is even
and [0,fs/2) when nfft is
odd. For a two-sided PSD estimate, f spans the
interval [0,fs). For a DC-centered PSD estimate, f spans
the interval (-fs/2, fs/2]
cycles/unit time for even length nfft and (-fs/2, fs/2)
cycles/unit time for odd length nfft .

Confidence bounds, specified as an N-by-2 matrix with real-valued
elements. The row dimension of the matrix is equal to the length of
the PSD estimate, pxx. The first column contains
the lower confidence bound and the second column contains the upper
confidence bound for the corresponding PSD estimates in the rows of pxx.
The coverage probability of the confidence intervals is determined
by the value of the probability input.