# Power spectral density of $\left(x(t)\right)^2$?

The relation between $x(t)$ and output $y(t)$ of a non-linear device is expressed as $$y(t) = (x(t))^2$$ Let $x(t)$ be zero-mean stationary Gaussian random process with auto-correlation $$R_x(\tau) = e^{-a\lvert \tau \rvert} , \quad a>0$$ Find the output PSD?

I tried using the formula of output $\textrm{PSD}_y = \lvert H(f)\rvert^2 \cdot \textrm{PSD}_x$ , but stuck in process of finding $H(f)$ of a non linear system squarer. Not able to think of any other idea.

• Hint: $R_y(\tau) = E[Y(t)Y(t+\tau)] = E[X^2(t)X^2(t+\tau)]$ where $X(t)$ and $X(t+\tau)$ are jointly Gaussian zero-mean unit-variance random variables with correlation coefficient $\rho = \exp(-a|\tau|)$. So start by ignoring fancy wordings and notation and figure out whether you can compute the expected value you need: $E[A^2B^2]$ where $A$ and $B$ are jointly Gaussian zero-mean unit-variance random variables with correlation coefficient $\rho$ for some $\rho \in (0,1)$. Sep 21 '16 at 19:51
• @DilipSarwate: I just saw your comment after finishing my answer. I was wondering how a beginner is supposed to figure out $E[A^2B^2]$ other than by looking up the formula (1) in my answer. Maybe I'm missing some simple shortcut ... Sep 21 '16 at 20:00
• @MattL. My hint was intended to get the OP off the track of trying to find out what "$H(f)$" is for a nonlinear system (cf. first sentence of your answer) , and look at the problem without the thicket of notation of $t$ and $\tau$ and random process and what not. Once the unnecessary stuff is removed, what is really needed becomes quite clear. Sep 21 '16 at 20:44
• hey @DilipSarwate, haven't seen you around recently. i would stick with small-case for functions of time (like $x(t)$) and leave the large-case names for functions of frequency (like $X(f)$). Sep 21 '16 at 21:22

Since the question has been raised as to whether the hint that I had given to the OP in a comment on the original question was appropriate for a newcomer to signal processing, here goes.

Stripped of extraneous baggage and notation, the question is whether it is possible to determine the value of $E[X^2Y^2]$ straightforwardly where $X$ and $Y$ are zero-mean unit variance jointly Gaussian random variables with correlation coefficient $\rho$, or is it necessary to resort to advanced methods or results such as the one cited in the answer by @MattL or generalizations thereof such as Price's Theorem (cf. the book cited in MattL's answer). One straightforward method involves conditional expectations and the result that is known to be a LIE:

$$E[X^2Y^2] = E\big[E[X^2Y^2\mid Y]\big]$$

Sorry, I couldn't resist: LIE is an acronym for the Law of Iterated Expectation which asserts that the expectation of $A$ is the same as the expectation of the random variable $E[A\mid B]$. Note that $E[A\mid B]$ is a function of $B$, not $A$ as one might naively believe. We already squeezed out all of the randomness in $A$ when we took its expectation. The LIE is that the expected value of $E[A\mid B]$, a function of $B$) works out, by a miracle of modern mathematics, to be the same as the expected value of $A$.

Now, for any choice of real number $y$, given that $Y=y$,

$$E[X^2Y^2\mid Y = y] = E[X^2y^2\mid Y = y] = y^2 E[X^2\mid Y = y].$$ But, given that $Y = y$, we know that the conditional distribution of $X$ given $Y=y$ is a Gaussian distribution with mean $\rho y$ and variance $1-\rho^2$. Consequently, $$E[X^2\mid Y = y] = \operatorname{var}(X\mid Y=y) + \left(E[X\mid Y = y]\right)^2 = (1-\rho^2) + (\rho y)^2$$ leading to $$E[X^2Y^2\mid Y = y] = y^2((1-\rho^2) + (\rho y)^2)) = y^2(1-\rho^2) + \rho^2y^4$$ giving us that the random variable $E[X^2Y^2\mid Y]$ is just $Y^2(1-\rho^2) + \rho^2Y^4$. Since $E[Y^2] = 1$ and $E[Y^4] = 3$, we get that $$E[X^2Y^2] = E\big[E[X^2Y^2\mid Y]\big] = E[Y^2(1-\rho^2) + \rho^2Y^4] = (1-\rho^2) + 3\rho^2 = 1 +2\rho^2.$$

Finally, noting that $\rho = e^{-a|t|}$, we can write $$R_y(t) = 1 + 2e^{-2a|t|}$$ just as has been explained in MattL's answer.

• +1 for creative solution. Your hint was of course good. My point was just that finding an appropriate expression for $E[X^2Y^2]$ is probably beyond the skills of most undergraduate students, and possibly many graduate students. Sep 23 '16 at 8:50

For this problem you can't use the formula involving $|H(f)|^2$ because it only applies to linear time-invariant (LTI) systems, and a squarer is obviously a non-linear system.

The only way to solve this problem that I can think of is to use the formula

$$E\{x^2y^2\}=E\{x^2\}E\{y^2\}+2E^2\{xy\}\tag{1}$$

which is valid for jointly Gaussian and zero mean random variables $x$ and $y$. I don't know stuff like this by heart but I look it up (e.g., Probability, Random Variables, and Stochastic Processes by Papoulis, 3rd ed., Eq. (7-36)).

We know that

$$R_y(\tau)=E\{y(t)y(t+\tau)\}=E\{x^2(t)x^2(t+\tau)\}\tag{2}$$

Since $x(t)$ and $x(t+\tau)$ are jointly Gaussian with zero mean, we can use $(1)$ to obtain

\begin{align}R_y(\tau)&=E\{x^2(t)\}E\{x^2(t+\tau)\}+2E^2\{x(t)x(t+\tau)\}\\&=R^2_x(0)+2R_x^2(\tau)\\&=1+2e^{-2a|\tau|}\tag{3}\end{align}

The PSD of $y(t)$ is now easily obtained by computing the Fourier transform of $(3)$.

• It is Equation 6-199 in the 4th edition of Papoulis and Pillai Sep 21 '16 at 21:05
• i didn't want to answer the question, but i would have used the fact that the spectrum of $x_1(t) \cdot x_2(t)$ is the convolution of the two spectra $X_1(f) \circledast X_2(f)$. Sep 21 '16 at 21:20
• @robertbristow-johnson Your comment is correct but not relevant to the problem being discussed. Sep 21 '16 at 22:09
• really @DilipSarwate? even if $x_1(t)=x_2(t)$? Sep 21 '16 at 22:33
• @robertbristow-johnson Yes, even if. We are not just multiplying two signals in the time domain, we are multiplying two random variables and then computing the expected value of the product, Thus, the convolution of the Fourier transforms is irrelevant. Sep 21 '16 at 22:44