资讯详情

正态随机过程仿真实验报告指南:从生成、验证到自相关与功率谱估计

📅 2026/10/11 16:46:31 | 华诺云谱 👁 阅读
正态随机过程仿真实验报告指南:从生成、验证到自相关与功率谱估计
简介这份资源是面向高校随机信号分析与处理课程学习者的实验报告聚焦「相关正态随机过程的仿真」这一典型实验适合正在完成随机过程仿真作业、需要参考完整实验流程与MATLAB实现思路的本科生。压缩包内仅含1个PDF文件约525KB完整收录实验目的、实验要求、程序代码、实验结果与实验体会五个部分结构清晰便于对照学习。报告围绕均匀分布随机数生成、白噪声正态序列构造、特定相关函数随机过程x(n)的递推生成、统计参数计算以及区间概率积分等任务展开并给出可直接参考的MATLAB代码片段与直方图分析结果。读者可借此理解正态分布与均匀分布随机过程之间的关系掌握均值、方差、相关函数的集合统计验证方法并对照理论值评估区间概率分布的一致性。目前已有183人学习下载适合作为实验报告撰写与代码调试的参考范例。1. 正态随机过程仿真到底在仿什么从一次实验报告翻车说起很多同学拿到“正态随机过程仿真”这个题目第一反应是打开软件生成一串randn画个图写句“服从正态分布”报告就交了。我带过几届做随机与信号实验的学生十份报告里有七份栽在同一个地方把“正态随机过程”当成了“正态随机变量”。前者是一族随时间变化的随机变量任意有限个时刻采样构成的向量服从多维正态分布后者只是单个时刻的分布。这个区别不搞清楚后面功率谱密度、自相关函数、各态历经性全是空中楼阁。这篇笔记面向正在做随机与信号实验、需要交出一份能自圆其说的仿真报告的人。我会把正态随机过程的仿真拆成可复现的步骤怎么生成、怎么验证分布、怎么估计自相关和功率谱、怎么判断各态历经以及报告里哪些图必须放、哪些参数必须交代。全程用常见工具实现不依赖任何特定平台你照着改参数就能跑。2. 生成正态随机过程从白噪声到有色谱的四步走2.1 先分清两种生成路线仿真正态随机过程有两条主流路线选错了后面全白做。第一条是直接法用独立同分布的高斯白噪声序列作为驱动经过一个线性时不变系统输出就是正态随机过程。因为高斯过程经过线性系统仍然是高斯过程这条路线数学上最干净适合验证自相关和功率谱的理论值。第二条是谱表示法在频域构造幅度谱乘上随机相位再做逆变换。这条路线适合生成具有指定功率谱密度形状的过程比如带限白噪声、一阶马尔可夫过程。缺点是相位随机化后单次实现的分布需要大样本才能逼近正态。我一般建议实验报告用直接法打底因为它的因果链条清晰白噪声→滤波器→有色噪声每一步都能单独验证。谱表示法作为补充用来展示功率谱控制的灵活性。2.2 用高斯白噪声驱动线性系统的最小实现下面这段 Python 代码生成一个一阶自回归过程也就是离散化的 RC 低通滤波输出它是典型的正态随机过程。import numpy as np import matplotlib.pyplot as plt # 参数设置 N 10000 # 样本点数 fs 1000 # 采样率 Hz rho 0.95 # 自回归系数决定相关时间 sigma_w 1.0 # 驱动白噪声标准差 # 生成高斯白噪声 np.random.seed(42) w np.random.normal(0, sigma_w, N) # 一阶自回归滤波x[n] rho * x[n-1] w[n] x np.zeros(N) for n in range(1, N): x[n] rho * x[n-1] w[n] # 理论方差sigma_x^2 sigma_w^2 / (1 - rho^2) sigma_x_theory sigma_w / np.sqrt(1 - rho**2) print(f样本方差: {np.var(x):.4f}) print(f理论方差: {sigma_x_theory**2:.4f})这段代码的核心逻辑是白噪声经过一阶递归滤波输出序列的任意有限维分布都是多维正态。rho越接近 1相关时间越长过程越“慢”。sigma_w控制驱动强度但输出方差由rho和sigma_w共同决定公式是 $\sigma_x^2 \sigma_w^2/(1-\rho^2)$。很多人只改sigma_w却发现方差对不上就是忘了这个关系。运行后你会看到样本方差和理论方差接近但不会完全相等因为这是单次实现。报告里要说明方差估计本身也是随机变量需要给出置信区间或多次实现的均值。2.3 参数怎么设相关时间、采样率与样本量的三角关系三个参数互相牵制设不好图就难看。相关时间$\tau_c$ 由rho决定。对于一阶自回归过程自相关函数是 $\rho^{|k|}$相关时间近似为 $-1/\ln(\rho)$ 个采样间隔。rho0.95时相关时间约 19.5 个采样点。如果你要观察自相关的衰减样本长度至少要是相关时间的 20 倍以上否则自相关估计的尾巴全是噪声。采样率$f_s$ 决定了你能看到的频率范围。归一化频率 $f/f_s$ 从 0 到 0.5。一阶自回归过程的功率谱是洛伦兹型半功率带宽约等于 $(1-\rho)/(1\rho)$ 乘以 $f_s/2$。rho0.95时带宽很窄谱线集中在低频这时候如果fs设得太小谱形状会被截断。样本量$N$ 影响估计的平滑度。功率谱估计的方差不会随 $N$ 增大而减小除非你做平均。单次周期图估计的标准差约等于均值本身也就是说谱估计的起伏和它的平均值一样大。报告里如果只放一张周期图审阅的人一眼就能看出你没做平滑。我通常这样配先定rho让相关时间在 10 到 50 个采样点之间再定fs让带宽占归一化频率的 5% 到 20%最后定N至少为相关时间的 100 倍。这样自相关和功率谱都能看得清楚。2.4 验证生成结果是否真的正态生成完了不能只看波形像不像要做分布检验。最直接的是画直方图叠加理论正态曲线再配合 Q-Q 图。from scipy import stats # 直方图与理论正态对比 fig, axes plt.subplots(1, 2, figsize(12, 4)) axes[0].hist(x, bins50, densityTrue, alpha0.7, label样本直方图) xx np.linspace(x.min(), x.max(), 200) axes[0].plot(xx, stats.norm.pdf(xx, 0, sigma_x_theory), r-, label理论正态) axes[0].legend() axes[0].set_title(分布对比) # Q-Q 图 stats.probplot(x, distnorm, plotaxes[1]) axes[1].set_title(Q-Q 图) plt.tight_layout() plt.show() # 偏度与峰度检验 skew stats.skew(x) kurt stats.kurtosis(x) print(f偏度: {skew:.4f} (理论 0)) print(f峰度: {kurt:.4f} (理论 0))偏度衡量对称性峰度衡量尾部厚度。正态分布的偏度为 0超额峰度为 0。样本量 10000 时偏度和峰度的标准差约分别为 $\sqrt{6/N}\approx0.024$ 和 $\sqrt{24/N}\approx0.049$。如果你的结果偏离超过三倍标准差要么样本量不够要么生成方法有问题。注意一阶自回归过程是相关序列直接做偏度和峰度检验时有效样本量小于名义样本量检验会偏乐观。更严格的做法是对不重叠的块取均值再检验或者用块自举法。实验报告里至少要用 Q-Q 图做定性判断不能只写“看起来像正态”。3. 自相关与功率谱估计把理论曲线和仿真曲线叠在一起3.1 自相关函数的无偏估计与边界效应自相关函数描述过程在不同时刻的相似程度。对于零均值过程理论自相关 $R[k]E[x[n]x[nk]]$。样本估计有两种有偏估计和无偏估计。def autocorr_biased(x, max_lag): N len(x) x x - np.mean(x) r np.correlate(x, x, modefull) r r[N-1:] / N return r[:max_lag] def autocorr_unbiased(x, max_lag): N len(x) x x - np.mean(x) r np.correlate(x, x, modefull) r r[N-1:] r r / (N - np.arange(len(r))) return r[:max_lag] max_lag 100 r_biased autocorr_biased(x, max_lag) r_unbiased autocorr_unbiased(x, max_lag) lags np.arange(max_lag) # 理论自相关 r_theory sigma_x_theory**2 * rho**lags plt.plot(lags, r_biased, label有偏估计) plt.plot(lags, r_unbiased, label无偏估计) plt.plot(lags, r_theory, k--, label理论值) plt.xlabel(滞后 k) plt.ylabel(自相关) plt.legend() plt.show()有偏估计除以固定的 $N$在滞后接近 $N$ 时估计值趋近于零偏差大但方差小。无偏估计除以 $N-k$均值是对的但滞后大时方差爆炸可能出现自相关大于理论最大值的情况。实验报告里我建议主图用有偏估计因为它的均方误差更小同时在图注里说明边界效应。理论曲线和仿真曲线叠在一起时你会发现小滞后处吻合很好大滞后处仿真曲线在零附近抖动。这个抖动范围大约是 $\pm 2/\sqrt{N}$ 乘以方差可以作为置信带画在图上。如果理论曲线超出这个带说明模型选错了。3.2 周期图法估计功率谱的方差问题功率谱密度是自相关函数的傅里叶变换。最直接的估计是周期图对样本做 FFT取模平方除以 $N$。from scipy.signal import periodogram, welch # 周期图估计 f, Pxx periodogram(x, fsfs, windowboxcar, nfft2048, return_onesidedTrue) # 理论功率谱一阶自回归过程 # S(f) sigma_w^2 / |1 - rho * exp(-j*2*pi*f/fs)|^2 f_theory np.linspace(0, fs/2, 500) S_theory sigma_w**2 / (1 rho**2 - 2*rho*np.cos(2*np.pi*f_theory/fs)) plt.semilogy(f, Pxx, label周期图) plt.semilogy(f_theory, S_theory, r--, label理论谱) plt.xlabel(频率 Hz) plt.ylabel(功率谱密度) plt.legend() plt.show()周期图的问题是方差大。单次周期图的每个频点近似服从指数分布标准差等于均值。也就是说谱估计的起伏和它的平均水平一样大。图上你会看到理论曲线穿过周期图的“毛刺”但毛刺幅度惊人。这不是代码错了是周期图的固有缺陷。解决办法是分段平均也就是 Welch 法。把长序列分成若干段每段做周期图然后平均。段数越多方差越小但频率分辨率越低。# Welch 法分段平均 f_w, Pxx_w welch(x, fsfs, windowhann, nperseg512, noverlap256, nfft1024) plt.semilogy(f_w, Pxx_w, labelWelch 估计) plt.semilogy(f_theory, S_theory, r--, label理论谱) plt.xlabel(频率 Hz) plt.ylabel(功率谱密度) plt.legend() plt.show()nperseg是每段长度noverlap是重叠点数。段数约等于 $N/(nperseg-noverlap)$。段数越多谱越平滑但频率分辨率约等于 $f_s/nperseg$。报告里要同时给出周期图和 Welch 图说明方差与分辨率的权衡。3.3 各态历经性怎么用仿真验证各态历经性是说时间平均等于集合平均。仿真只能验证不能证明。做法是生成多条独立实现计算每条的时间平均再看这些时间平均的分布是否集中在集合平均附近。M 50 # 独立实现次数 N_short 2000 # 每条实现的长度 means np.zeros(M) vars_ np.zeros(M) for m in range(M): w_m np.random.normal(0, sigma_w, N_short) x_m np.zeros(N_short) for n in range(1, N_short): x_m[n] rho * x_m[n-1] w_m[n] means[m] np.mean(x_m) vars_[m] np.var(x_m) print(f时间平均的均值: {np.mean(means):.4f}, 标准差: {np.std(means):.4f}) print(f时间方差的均值: {np.mean(vars_):.4f}, 标准差: {np.std(vars_):.4f}) print(f理论方差: {sigma_x_theory**2:.4f})如果过程是各态历经的时间平均的均值应该接近 0时间方差的均值应该接近理论方差。但注意一阶自回归过程是各态历经的前提是rho的绝对值小于 1。如果rho1过程变成随机游走就不是各态历经了。仿真时rho不要设成 1 或超过 1。时间平均的标准差随着 $N_short$ 增大而减小但减小的速度受相关时间影响。相关时间越长需要越长的 $N_short$ 才能让时间平均稳定。报告里可以画时间平均随 $N_short$ 变化的曲线展示收敛趋势。4. 避坑与排查正态随机过程仿真里最容易翻车的五个地方4.1 现象直方图不像正态尾部翘起原因驱动噪声不是高斯分布或者滤波器是非线性的。常见错误是用均匀分布随机数经过线性滤波输出虽然接近正态但尾部不对。另一个原因是rho太接近 1过程接近随机游走单次实现的分布严重偏离正态。解决确认驱动噪声用randn或等价的高斯生成器。检查滤波器是否线性。如果rho超过 0.99要么减小rho要么大幅增加样本量并做去趋势。4.2 现象自相关理论曲线和仿真曲线对不上原因最常见的是均值没去掉。理论自相关假设零均值如果样本均值不为零自相关会被一个常数项污染。另一个原因是理论公式用错了比如把连续时间的自相关套到离散序列上。解决估计自相关前先减去样本均值。离散一阶自回归的自相关是 $\rho^{|k|}$连续一阶马尔可夫的自相关是 $\exp(-|t|/\tau)$采样后是 $\exp(-|k|/(\tau f_s))$不要混用。4.3 现象功率谱估计在低频处远高于理论值原因样本均值没去掉相当于在过程上叠加了一个直流分量。周期图在零频附近会有一个巨大的尖峰泄漏到低频段。另一个原因是用了矩形窗旁瓣泄漏严重。解决去均值。换用 Hann 窗或 Hamming 窗。如果低频仍然偏高检查rho是否过大导致谱峰太窄周期图的分辨率不够。4.4 现象各态历经验证时时间方差波动很大原因每条实现的长度太短时间平均还没收敛。相关时间越长需要的实现长度越长。用rho0.95时相关时间约 20 个采样点实现长度 2000 只有 100 个相关时间时间方差的相对标准差约 $1/\sqrt{100}10%$。解决增加每条实现的长度或者增加实现次数看整体趋势。报告里不要只给一个数要给分布或误差棒。4.5 现象换一个随机种子结果差异巨大原因样本量不够或者过程的相关时间太长导致有效样本量很小。有效样本量约等于 $N$ 除以相关时间的两倍。N10000、相关时间 20 时有效样本量只有 250估计的波动自然大。解决增加样本量或者做多次实现取平均。报告里要固定随机种子以便复现同时说明单次结果的随机性范围。5. 让报告站得住脚从一张图到一套证据链实验报告的说服力不在于图多而在于每张图都在回答一个具体问题。我一般按这个顺序组织先给生成方法的框图说明再给一段典型波形接着是分布检验然后是自相关和功率谱的对比最后是各态历经的收敛曲线。每张图下面用一两句话说明“这张图验证了什么”。一个容易被忽略的技巧是把理论值算出来叠在仿真图上而不是只放仿真结果。理论自相关、理论功率谱、理论方差这些都有闭式解算出来不费事但能让审阅的人一眼看出你的仿真是否可信。如果理论曲线和仿真曲线在置信带内吻合你的工作就站住了。另一个技巧是给出参数表。rho、sigma_w、fs、N、nperseg、noverlap这些参数集中列在一张表里比散落在正文中更容易复现。我习惯在表里加一列“选择理由”比如“rho0.95 使相关时间约 20 个采样点便于观察自相关衰减”。最后各态历经的验证不要只写“时间平均接近集合平均”要给出时间平均的分布图或误差棒。如果时间平均的标准差是理论值的 5%就明确写出来。审阅的人想看到的是你对估计精度的量化意识而不是一句定性结论。我自己做这类实验时习惯在代码里留一个config字典把所有参数集中管理跑完自动生成报告所需的全部图表。这样换参数重跑只需要改一行不用满代码找数字。这个习惯帮我省下了大量重复劳动也避免了参数不一致导致的低级错误。希望帮到你。本文还有配套的精品资源点击获取
📝

华诺云谱内容团队

资深建站顾问 · 行业研究员

10年+企业数字化服务经验,专注智能建站、SEO优化与品牌营销,持续输出建站技巧、行业洞察与营销干货,已帮助5000+企业实现数字化增长。

你可能需要的服务

订阅华诺云谱资讯周报

每周一封,精选建站技巧、SEO与营销干货,直达邮箱。已有 8,000+ 企业主订阅,助你少走弯路。

↑