带粉红噪声的时滞微分方程

Question

我需要解一组延迟微分方程，我想在这些方程中使用粉红噪声。在

我试着用Pydelay包for Python来做，但问题是要么生成噪声，然后将其传递给模拟，要么在模拟过程中生成噪声。第一个选项不起作用，因为在该解算器中指定的参数必须是常量，且噪波会随时间变化。第二个选项不起作用，因为Pydelay只支持生成每一个样本都独立于其他样本的噪声，而粉红噪声则不是这样。在

有人知道怎么做吗？在

这是我的代码（我需要I1和I2是粉红色的噪音，而不是像我的代码中那样是恒定的）：

# -*- coding: utf-8 -*-
"""
Created on Tue Aug  2 15:35:15 2016

@author: kasienka
"""
#!python
import pydelay
from time import time
import math
import numpy as np
from numpy import fft
from scipy import integrate
import matplotlib as mpl
mpl.use('Agg')
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.interpolate import griddata
import sys
import pylab as pl
from pydelay import dde23

# define the equations
eqns = {
    'y1' : 'I1 - y1(t-tau) + epsilon * pow(y2(t-tau), 1.1)',
    'y2' : 'I2 - y2(t-tau) + epsilon * pow(y1(t-tau), 1.1)'
    }
#define the parameters
params = {
    'I2' : 0.2,
    'I1' : 0.4,
    'tau': 1.31,
    'epsilon': 0.2
    }

# Initialise the solver
dde = dde23(eqns=eqns, params=params)

# set the simulation parameters
dde.set_sim_params(tfinal=50, dtmax=0.0001)

histdic = {
    'y1': lambda t: 0.2,
    'y2': lambda t: 0.4
    }
dde.hist_from_funcs(histdic, 1000)

# run the simulator
dde.run()
p1 = []
p2 = []
x1 = []
x2 = []
czasy = np.linspace(0, 50, 10000)
for el in times:
    x10 = dde.sol_spl(el)['y1']
    x20 = dde.sol_spl(el)['y2']
    prob1 = 1.0 / (1 + 3.14**(-(x10 - x20)) )
    prob2 = 1.0 / (1 + 3.14**(-(x20 - x10)) )
    p1.append(prob1)
    p2.append(prob2)
    x1.append(x10)
    x2.append(x20)
fig = pl.figure()
pl.plot(times, p1, label = '$p_1$')
pl.plot(times, p2, label = '$p_2$')
pl.xlabel('$time$')
pl.ylabel('$probabilities$')
pl.legend()
pl.savefig(sys.argv[1])

Answer 1

所以问题是为数据集中的所有元素生成粉色噪音？在

def voss(nrows, ncols=16):
    """Generates pink noise using the Voss-McCartney algorithm.

    nrows: number of values to generate
    rcols: number of random sources to add

    returns: NumPy array
    """
    array = np.empty((nrows, ncols))
    array.fill(np.nan)
    array[0, :] = np.random.random(ncols)
    array[:, 0] = np.random.random(nrows)

    # the total number of changes is nrows
    n = nrows
    cols = np.random.geometric(0.5, n)
    cols[cols >= ncols] = 0
    rows = np.random.randint(nrows, size=n)
    array[rows, cols] = np.random.random(n)

    df = pd.DataFrame(array)
    df.fillna(method='ffill', axis=0, inplace=True)
    total = df.sum(axis=1)

    return total.values

来源：ThinkDSP

Answer 2

我假设你想要动态噪声。在这种情况下，您首先应该认识到，有几种类型的微分方程具有根本不同的解算器：

• 延迟微分方程（DDE）通常采用以下任一方法求解：

  • 嵌入的、多步骤的Runge-Kutta方法，其中插入过去（这是Pydelay所做的），

  • 用积分步分割所有延迟的单步积分器。

• 随机微分方程（SDEs）的求解采用简单的单步方法（全部基于Euler方法）。多步方法仍然是一个热点问题，而嵌入式方法是最近才被提出的。在

  我所知的所有关于sde的论文（理论和方法，但它们并不多见），只考虑白噪声（Wiener过程）；认识到粉红噪声本身就是一个问题。我只需快速搜索就能找到a paper on simulating white noise with a white-noise SDE。

• 随机时滞微分方程（SDDEs）又需要特殊的解。我对它们不是很熟悉，但显然你至少继承了DDE和SDE的所有局限性。因此，最好的解算器应该是一个单步方法，用积分步骤分割所有延迟。快速浏览一下文献，这确实是我们所做的。请注意，这并不能解决粉红色噪音的问题。

如前所述，Pydelay使用带插值的多步方法。这不是为随机性设计的，也不能正确地处理随机性（如果是这样的话，它已经是SDE的一个出色的解决方案，也就是说，毫不拖延）。相反，它所做的是将多步方法视为一个单步方法，然后添加噪声。正如Pydelay的作者自己所说，这“相当粗糙”。此外，实际使用粉红噪声（即使使用粗略的方法）将需要您重新编程软件，并可能导致积分器误差估计固有的问题。我强烈建议不要这样。即使你成功了，使用高级DDE解算器（如自适应步骤）的所有优势也会丧失，因此从头开始编写一个新的、更简单的积分器会更容易。在

如果您真的需要这样做，我建议您先找出如何解决带有粉红噪声的SDEs，然后将该方法扩展到SDDEs（希望这是相当简单的）。在

---

^{¹如果您想要观察噪声，那很简单：只需将其添加到您的解决方案中。}