# A probability density function for pink noise

I would like to understand pink (or $$1/f$$) noise better. However, clearly written resources are difficult to come by, and are usually concerned with its Fourier spectrum, or qualitative descriptions of it as a time-dependent process.

I hope that someone can describe how a pink noise process $$P(t)$$ depends on its values at earlier times $$t' < t$$, in the form of a probability density function on $$P(t) - P(t')$$, in the same way that we can do for Brown noise (i.e. Brownian motion).

Preliminaries

As a 'process', white noise can be modelled as a collection of i.i.d. random variables $$W(t)$$, where for each $$t$$ the distribution of $$W(t)$$ is uniform over some interval $$[-a,a]$$. Brown or red noise is the time-integral of white noise, in the sense that $$B(t) = \int_0^t W(s) \,\mathrm d s \;.$$ Note that $$B(t)$$ is not independent of $$B(t')$$ for $$t' < t$$, in that $$B_t - B_{t'}$$ is normally distributed with mean $$0$$ and variance proportional to $$(t - t')$$. Thus we can talk about the expected properties of the trajectory at discrete time intervals in a straight-forward way.

If we take the Fourier transform of a white noise process over a time interval, we obtain a function whose expected modulus-square is constant over all frequencies, albeit with increasing variance for higher frequencies. If we do the same for Brown noise, we obtain a function whose expected modulus-square varies as $$1/f^2$$ with respect to frequency, again with increasing variance for higher frequencies.

'Defining' pink noise

Intermediate to these is pink noise, which is a process whose Fourier transform has a modulus-square proportional to $$1/f$$ (so that it is often called $$1/f$$ noise). For a discrete-time process, where we define $$W(t)$$ on the integers $$t \in [0,T]$$ for some $$T$$, we may mock up a pink noise process by setting \begin{aligned} \hat P(0) &= \hat W(0) \\ \hat P(k) &= \hat W(k) \big/ \sqrt{\lvert k \rvert}, \quad k \in [-T{+}1,T] \setminus \{ 0 \} \end{aligned} and then defining $$P$$ to be the inverse Fourier transform of $$\hat P$$: this is what the clearest sources that I can find suggest.

A small amount of numerical experimentation on my part suggests that the resulting function $$P(t)$$ has a much smaller variance than the original white-noise process $$W(t)$$. The time-integral $$C(t)$$ (i.e., a cumulative sum) of $$P(t)$$ describes a trajectory which is less volatile than the Brown noise process $$B(t)$$ obtained from the same white-noise function, but otherwise the two time-integrated noise processes follow a similar trajectory. Both of these observations are not very surprising, given that $$P$$ arises by decreasing the importance of the high-frequency components of $$W$$.

Question.

There are some qualitative descriptions of pink noise, in terms of how $$P(t)$$ depends on $$P(t')$$ for $$t' < t$$, such as that the autocorrelation of $$P(t)$$ with $$P(t')$$ decreases exponentially towards some non-zero value and stays there. However, I would like a quantitative result similar to what we have for Brown noise — specifically: is there a precise description for the probability density function of $$P(t) - P(t')$$?

(Failing this, is there a better description for how to simulate a pink noise process than the one that I've outlined?)

• What you call pink noise would be called log-correlated noise in mathematics. See for example arxiv.org/abs/1711.00427 for a non-stationary construction. The construction cannot be made stationary on any space of functions / distributions, but it can be made stationary on the space of Schwartz distributions quotiented by constants. – Martin Hairer Aug 16 '19 at 15:41
• For my own future reference, I'll note that in the above preprint, Ref. [12] –– B. Duplantier, R. Rhodes, S. Sheffield, and V. Vargas. Log-correlated Gaussian fields: an overview. [arXiv:1407.5605], 2014 –– seems a promising introductory resource. – Niel de Beaudrap Aug 16 '19 at 15:58

The spectral density of $$1/f$$ noise is $$S(\omega)=c/\omega$$ for $$\omega_1<\omega<\omega_2$$. The corresponding autocorrelation function in the time domain $$P(t)=\lim_{T\rightarrow \infty}\frac{1}{T}\int_0^T x(t')x(t'+t)\,dt'$$ is given by the Wiener-Khinchine theorem, $$P(t)=\frac{c}{2\pi}\int_{\omega_1}^{\omega_2}\frac{\cos\omega t}{\omega}\,dt=\frac{c}{2\pi}\left[C(\omega_1 t)-C(\omega_2 t)\right],$$ with $$C(x)$$ the cosine integral. The variance of the noise process is $$\mathbb{E}[x(t)^2]=P(0)=\frac{c}{2\pi}\log(\omega_2/\omega_1).$$
• I'm having some difficulties connecting your post to what I have read. Your remark on spectral density (or the power distribution) is of course standard. But it isn't clear to me why the autocorrelation function should only involve a cosine contribution, or how I should relate the autocorrelation which you denote by $P(t)$ to the expected value (resp. the variance) of the noise process, which in the notation of my post would be $\langle P(t) \rangle$ (resp. $\langle P(t)^2 - \langle P(t)^2 \rangle \rangle$), as $P(t)$ differs from $0$ for finite $t$ (and converges to $0$ as $t \to \infty$). – Niel de Beaudrap Aug 16 '19 at 16:30
• the Fourier transform of $S(\omega)$ has only a $\cos\omega t$ contribution because $P(t)=P(-t)$ for a stationary process. – Carlo Beenakker Aug 16 '19 at 18:04
• This is much clearer, thanks. I must admit that in my numerical experiments (where I am considering a discrete-time variant in which I suppose that I have $\omega_2/\omega_1 = 10^5$), I seem to have a variance a few orders of magnitude smaller than what I would expect from the spectral density, given your analysis. However, it may be that I should more carefully consider my analysis. – Niel de Beaudrap Aug 19 '19 at 10:45