如何用Python正确滤波白噪声生成一阶高斯-马尔可夫过程？

我是Python新手，希望从高斯白噪声合成一阶Gauss-Markov过程。根据信号处理理论，可通过设计合适的噪声整形滤波器实现。一阶Gauss-Markov过程有两个参数：sigma（过程的标准差）和时间常数beta。

整形滤波器的传递函数如下：

我的代码如下：

import scipy.signal as dsp
import numpy as np

Nsamples = 2000
fs = 100

time = np.arange(Nsamples) / fs

rng = np.random.default_rng()
gaussianNoise = rng.standard_normal(size=time.shape)

wgn = (gaussianNoise - np.mean(gaussianNoise)) / np.std(gaussianNoise)
print('\n\n\nWGN MEAN: ', np.mean(wgn))
print('WGN STD: ', np.std(wgn))


beta = 0.01
sigma = 0.1
b = np.array([np.sqrt(2 * beta * sigma**2)])
a = np.array([1, beta])

gaussMarkovNoise = dsp.lfilter(b, a, whiteGaussianNoise)

但结果有问题：Gauss-Markov噪声的自相关应呈指数衰减，但按此方式滤波后，其自相关仍像白噪声那样仅在原点有尖峰。请问我忽略了什么？

1. 变量名语法错误 #

代码最后一行使用了未定义的whiteGaussianNoise，应替换为前面定义的wgn，否则代码无法正常执行。

2. 连续域到离散域的滤波器转换错误 #

你使用的传递函数是连续时间域的，但scipy.signal.lfilter需要离散时间域的滤波器系数，直接套用连续域系数会导致滤波逻辑完全错误。

连续时间传递函数为：
$$ H(s) = \frac{\sqrt{2\beta}\sigma}{s + \beta} $$

针对一阶系统，采用冲激响应不变法离散化更直观，结合采样间隔 $T = 1/fs$，离散时间的递归关系为：
$$ x[n] = e^{-\beta T} x[n-1] + \sigma \sqrt{1 - e^{-2\beta T}} w[n] $$

对应的滤波器系数：

• 分子 $b = [\sigma \sqrt{1 - e^{-2\beta T}}]$
• 分母 $a = [1, -e^{-\beta T}]$

这里的关键错误是分母系数符号：连续域的$s+\beta$对应离散化后递归项的负系数，你之前写的+beta完全颠倒了逻辑。

3. 冗余的白噪声归一化 #

rng.standard_normal()生成的已经是均值为0、标准差为1的高斯白噪声，无需再做(gaussianNoise - np.mean(gaussianNoise)) / np.std(gaussianNoise)的归一化操作。

修正后的代码 #

import scipy.signal as dsp
import numpy as np
import matplotlib.pyplot as plt

Nsamples = 2000
fs = 100
T = 1 / fs  # 采样间隔

time = np.arange(Nsamples) / fs

rng = np.random.default_rng()
wgn = rng.standard_normal(size=time.shape)
print('WGN MEAN: ', np.mean(wgn))
print('WGN STD: ', np.std(wgn))


beta = 0.01  # 连续时间的时间常数
sigma = 0.1  # 过程的标准差

# 计算离散化后的滤波器系数
beta_discrete = np.exp(-beta * T)
b = np.array([sigma * np.sqrt(1 - beta_discrete**2)])
a = np.array([1, -beta_discrete])

gaussMarkovNoise = dsp.lfilter(b, a, wgn)

# 计算并绘制自相关函数验证结果
def compute_autocorrelation(x):
    result = np.correlate(x, x, mode='full')
    return result[len(result)//2:] / result[len(result)//2]

wgn_ac = compute_autocorrelation(wgn)
gm_ac = compute_autocorrelation(gaussMarkovNoise)

plt.figure(figsize=(10, 6))
plt.plot(time, wgn_ac[:Nsamples], label='高斯白噪声自相关')
plt.plot(time, gm_ac[:Nsamples], label='一阶Gauss-Markov过程自相关')
plt.xlabel('时间 (s)')
plt.ylabel('归一化自相关')
plt.title('自相关函数对比')
plt.legend()
plt.grid(True)
plt.show()

验证结果 #

运行修正后的代码，会看到Gauss-Markov过程的自相关函数呈现明显的指数衰减，完全符合理论预期——这是因为离散化后的滤波器正确实现了一阶马尔可夫过程的递归逻辑：当前时刻的值由前一时刻的值衰减后，叠加缩放后的白噪声项构成。

内容的提问来源于stack exchange，提问作者eljamba