如何在 R 中模拟粉红噪声

发布于 2024-12-24 01:16:46 字数 57 浏览 2 评论 0原文

我知道可以通过将 rnorm() 的输出视为时间序列来实现白噪声。关于如何模拟粉红噪声有什么建议吗?

I know that white noise can be achieved by treating the output of rnorm() as a timeseries. Any suggestions on how to simulate pink noise?

如果你对这篇内容有疑问,欢迎到本站社区发帖提问 参与讨论,获取更多帮助,或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论, 你可以免费 注册 一个本站的账号。

评论(2

你的背包 2024-12-31 01:16:46

软件包 tuneR 具有 noise 函数,可以生成白噪声或粉红噪声的波对象:

require(tuneR)
w <- noise(kind = c("white"))
p <- noise(kind = c("pink"))
par(mfrow=c(2,1))
plot(w,main="white noise")
plot(p,main="pink noise")

编辑:我意识到上面的方法不会生成矢量(doh) 。将其转换为矢量的残酷方法是添加以下代码:

writeWave(p,"p.wav")#writes pink noise on your hard drive
require(audio)#loads `audio` package to use `load.wave` function
p.vec <- load.wave("path/to/p.wav")#this will load pink noise as a vector

在此处输入图像描述

Package tuneR has noise function which can generate a wave object that is either white or pink noise:

require(tuneR)
w <- noise(kind = c("white"))
p <- noise(kind = c("pink"))
par(mfrow=c(2,1))
plot(w,main="white noise")
plot(p,main="pink noise")

EDIT: I realized that the method above doesn't generate the vector (doh). Brutal way to convert it into the vector is to add the code below:

writeWave(p,"p.wav")#writes pink noise on your hard drive
require(audio)#loads `audio` package to use `load.wave` function
p.vec <- load.wave("path/to/p.wav")#this will load pink noise as a vector

enter image description here

迷迭香的记忆 2024-12-31 01:16:46

正如@mbq所说,您可以使用 p@left 来获取向量,而不是保存和读取 wav 文件。另一方面,您可以直接使用在tuneR中生成时间序列的函数:

TK95 <- function(N, alpha = 1){ 
    f <- seq(from=0, to=pi, length.out=(N/2+1))[-c(1,(N/2+1))] # Fourier frequencies
    f_ <- 1 / f^alpha # Power law
    RW <- sqrt(0.5*f_) * rnorm(N/2-1) # for the real part
    IW <- sqrt(0.5*f_) * rnorm(N/2-1) # for the imaginary part
    fR <- complex(real = c(rnorm(1), RW, rnorm(1), RW[(N/2-1):1]), 
                  imaginary = c(0, IW, 0, -IW[(N/2-1):1]), length.out=N)
     # Those complex numbers that are to be back transformed for Fourier Frequencies 0, 2pi/N, 2*2pi/N, ..., pi, ..., 2pi-1/N 
     # Choose in a way that frequencies are complex-conjugated and symmetric around pi 
     # 0 and pi do not need an imaginary part
    reihe <- fft(fR, inverse=TRUE) # go back into time domain
    return(Re(reihe)) # imaginary part is 0
}

这​​非常有效:

par(mfrow=c(3,1))
replicate(3,plot(TK95(1000,1),type="l",ylab="",xlab="time"))

< img src="https://i.sstatic.net/Xd7EI.png" alt="在此处输入图像描述">

As said by @mbq, you can just use p@left to get the vector, instead of saving and reading the wav file. On the other hand, you could directly use the function generating the time serie in tuneR:

TK95 <- function(N, alpha = 1){ 
    f <- seq(from=0, to=pi, length.out=(N/2+1))[-c(1,(N/2+1))] # Fourier frequencies
    f_ <- 1 / f^alpha # Power law
    RW <- sqrt(0.5*f_) * rnorm(N/2-1) # for the real part
    IW <- sqrt(0.5*f_) * rnorm(N/2-1) # for the imaginary part
    fR <- complex(real = c(rnorm(1), RW, rnorm(1), RW[(N/2-1):1]), 
                  imaginary = c(0, IW, 0, -IW[(N/2-1):1]), length.out=N)
     # Those complex numbers that are to be back transformed for Fourier Frequencies 0, 2pi/N, 2*2pi/N, ..., pi, ..., 2pi-1/N 
     # Choose in a way that frequencies are complex-conjugated and symmetric around pi 
     # 0 and pi do not need an imaginary part
    reihe <- fft(fR, inverse=TRUE) # go back into time domain
    return(Re(reihe)) # imaginary part is 0
}

and this works perfectly :

par(mfrow=c(3,1))
replicate(3,plot(TK95(1000,1),type="l",ylab="",xlab="time"))

enter image description here

~没有更多了~
我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站,即表示您同意使用 Cookies 和您的相关数据。
原文