1、Required Knowledge

1.1、概率密度函数

用于描述连续型随机变量在不同取值上的概率密度，记作$f(x)$。
如随机变量$X$的分布为正态分布，则其概率密度函数为：
$$f(x)=\frac{1}{\sigma \sqrt{2\pi }}{e}^{-\frac{(x-\mu {)}^2}{2{\sigma }^2}}$$
记作：
$$X\sim \mathcal{N}(\mu ,{\sigma }^2)$$

1.2、联合概率密度

用于描述两个或多个连续型随机变量联合分布的概率密度。它表明这些随机变量在不同取值上的概率密度，是多个变量的联合概率分布的具体表现形式。
例如，记往坐标轴上打靶的位置为随机事件$A$，该事件受到随机变量$X$和$Y$的影响，那么打靶的位置，即二维随机变量$(X,Y)$的概率密度叫作$X$和$Y$的联合概率密度，记作$f(x,y)$。

1.3、边缘分布概率密度

二维随机变量$(X,Y)$有概率密度函数，而$X$和$Y$都是随机变量，各自也有各自的概率密度，分别记作$f_X(x)$ 和 $f_Y(y)$，有：
$$f_X(x)={\int }_{-\infty }^{\infty }{f}_{X,Y}(x,y)dy$$
$$f_Y(y)={\int }_{-\infty }^{\infty }{f}_{X,Y}(x,y)dx$$

1.4、独立性

A、B是两事件，如果满足
$$P(A\cap B)=P(A)P(B)$$
等同于：
$$P(AB)=P(A)P(B)$$
则称事件 $A$，$B$ 相互独立。事件 $A$ 和事件 $B$ 同时发生的概率等于事件A发生的概率乘以事件B发生的概率，即：事件 $A$，$B$ 是否发生不受另一事件的影响。

2、独立二维随机变量的联合概率密度

如二维随机变量$(X,Y)$相互独立，则其联合概率密度$f(x,y)=f(x)f(y)$。现有变量$X,Y$均服从标准正态分布，其联合分布的概率密度为如图 $A$ 所示的曲面。

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.stats import norm# 假设 X 和 Y 都是独立的标准正态分布 N(0, 1)
mu_X, sigma_X = 0, 1
mu_Y, sigma_Y = 0, 1# 定义网格范围和步长
x = np.linspace(-4, 4, 200)
normal_x = norm.pdf(x, mu_X, sigma_X)
y = np.linspace(-4, 4, 200)
normal_y = norm.pdf(y, mu_Y, sigma_Y)
x, y = np.meshgrid(x, y)# 计算联合概率密度 f(x, y)
f_X = norm.pdf(x, mu_X, sigma_X)
f_Y = norm.pdf(y, mu_Y, sigma_Y)
f_XY = f_X * f_Y  # 由于独立性，联合密度等于单独密度的乘积# 计算 Z = X + Y 的概率密度
z = x + y  # Z 的值
f_Z = norm.pdf(z, mu_X + mu_Y, np.sqrt(sigma_X**2 + sigma_Y**2))# 创建一个三维图形对象
fig = plt.figure(figsize=(18, 6))# 绘制联合分布的三维图
ax1 = fig.add_subplot(131, projection='3d')
surf = ax1.plot_surface(x, y, f_XY, cmap='hot', edgecolor='k', alpha=0.5)
ax1.set_title('Joint Probability Density of X and Y')
ax1.set_xlabel('X')
ax1.set_ylabel('Y')
ax1.set_zlabel('Density')
ax1.view_init(elev=30, azim=45)  # 调整视角# 分别在x=4、y=4的平面上绘制原始的X、Y标准正态分布图
ax1.plot(x[0, :], np.full_like(x[0, :], 4), normal_x, color='green', lw=2, label='Normal X')
ax1.plot(np.full_like(y[:, 0], 4), y[:, 0], normal_y, color='black', lw=2, label='Normal Y')# 添加图例
ax1.legend()plt.tight_layout()
plt.show()

如图所示，绿色和黑色分别为 $X$ 和 $Y$ 的标准正态概率密度函数图，而曲面是联合概率密度函数图像。可以看出 $(0,0)$ 这个点的概率密度最大，符合标准正态分布的特性。

3、边缘分布的几何学解释

由联合分布的概率密度函数可知：
$$f_X(x)={\int }_{-\infty }^{\infty }{f}_{X,Y}(x,y)dy$$
$$f_Y(y)={\int }_{-\infty }^{\infty }{f}_{X,Y}(x,y)dx$$
如何理解这个公式呢？

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.stats import norm# 假设 X 和 Y 都是独立的标准正态分布 N(0, 1)
mu_X, sigma_X = 0, 1
mu_Y, sigma_Y = 0, 1# 定义网格范围和步长
x = np.linspace(-4, 4, 200)
normal_x = norm.pdf(x, mu_X, sigma_X)
y = np.linspace(-4, 4, 200)
normal_y = norm.pdf(y, mu_Y, sigma_Y)
x, y = np.meshgrid(x, y)# 计算联合概率密度 f(x, y)
f_X = norm.pdf(x, mu_X, sigma_X)
f_Y = norm.pdf(y, mu_Y, sigma_Y)
f_XY = f_X * f_Y  # 由于独立性，联合密度等于单独密度的乘积# 计算 Z = X + Y 的概率密度
z = x + y  # Z 的值
f_Z = norm.pdf(z, mu_X + mu_Y, np.sqrt(sigma_X**2 + sigma_Y**2))## 创建一个三维图形对象
fig = plt.figure(figsize=(18, 6))# 绘制联合分布的三维图
ax1 = fig.add_subplot(131, projection='3d')
surf = ax1.plot_surface(x, y, f_XY, cmap='hot', edgecolor='k', alpha=0.5)
ax1.set_title('Joint Probability Density of X and Y')
ax1.set_xlabel('X')
ax1.set_ylabel('Y')
ax1.set_zlabel('Density')
ax1.view_init(elev=30, azim=45)  # 调整视角# 分别在x=4、y=4的平面上绘制原始的 X、Y 标准正态分布图
ax1.plot(x[0, :], np.full_like(x[0, :], 4), normal_x, color='green', lw=2, label='Normal X')
ax1.plot(np.full_like(y[:, 0], 4), y[:, 0], normal_y, color='black', lw=2, label='Normal Y')# 计算边缘分布
# 计算边缘分布
marginal_x = np.sum(f_XY, axis=1) * (y[1, 0] - y[0, 0])  # 对 y 积分得到 x 的边缘分布
marginal_y = np.sum(f_XY, axis=0) * (x[0, 1] - x[0, 0])  # 对 x 积分得到 y 的边缘分布# 分别在x=4.2、y=4.2的平面上绘制计算 X、Y 边缘分布图
ax1.plot(x[0, :], np.full_like(x[0, :], 4.2), marginal_x, color='red', lw=2, label='Marginal X')
ax1.plot(np.full_like(y[:, 0], 4.2), y[:, 0], marginal_y, color='blue', lw=2, label='Marginal Y')# 添加图例
ax1.legend()plt.tight_layout()
plt.show()

可以看到，原始的 X、Y 分布与基于联合分布函数计算得到的边缘分布是一样的。

以 X 为例，根据代码 marginal_x = np.sum(f_XY, axis=1) * (y[1, 0] - y[0, 0]) 可知，这就是对公式的代码化表达。

要想绘制 $X$ 的边缘分布，就是对每个 $X$ 的取值计算其概率密度。不妨令 $X=0$ 来理解这个公式。当 $X=0$ 时，联合概率密度如下图：

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm# 定义网格范围和步长
y = np.linspace(-4, 4, 200)  # y 的取值范围
x = np.linspace(-4, 4, 200)  # x 的取值范围
y, x = np.meshgrid(y, x)# 假设 X 和 Y 都是独立的标准正态分布 N(0, 1)
mu_X, sigma_X = 0, 1
mu_Y, sigma_Y = 0, 1# 计算联合概率密度 f(x, y)
f_X = norm.pdf(x, mu_X, sigma_X)
f_Y = norm.pdf(y, mu_Y, sigma_Y)
f_XY = f_X * f_Y  # 由于独立性，联合密度等于单独密度的乘积# 仅保留 x = 0 的值
tolerance = 0.05  # 容差，用于确定接近0的范围
mask = np.abs(x) < tolerance
f_XY_filtered = np.where(mask, f_XY, 0)  # 将不满足条件的设为 NaN# 创建一个三维图形对象
fig = plt.figure(figsize=(10, 8))# 绘制 f(x,y) （仅当 x = 0 时）
ax = fig.add_subplot(111, projection='3d')
surf = ax.plot_surface(x, y, f_XY_filtered, cmap='viridis', edgecolor='k', alpha=0.7)
ax.set_title('f(x,y) (with $x = 0$)')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Density')plt.tight_layout()
plt.show()

当 $x=0$ 时，$f(x,y)$ 为图中平面。
$$f_X(0)={\int }_{-\infty }^{\infty }{f}_{X,Y}(0,y)dy$$
在 $f(x,y)$ 在 $y$ 轴上求积分，就是图中面的面积。就等于 $X$ 的边缘分布函数（标准正态分布）在 $x=0$ 时的值。

4、两个随机变量的函数的分布

之前讲了单个随机变量的概率密度函数，两个随机变量的概率密度函数，那么两个随机变量组成的新的随机变量，其概率密度函数是什么？

设 $(X,Y)$ 是二维连续型随机变量，它具有的概率密度 $f(x,y)$，则 $Z=X+Y$ 仍为连续型随机变量，其概率密度为：
$$f_Z(z)={\int }_{-\infty }^{\infty }f(z-y,y)dy$$
或
$$f_Z(z)={\int }_{-\infty }^{\infty }f(x,z-x)dx$$
若 X 和 Y 相互独立，X和Y的边缘概率密度分别为$f_X(x)$和$f_Y(y)$，上式变为：
$$f_Z(z)={\int }_{-\infty }^{\infty }{f}_X(z-y)f_Y(y)dy$$
$$f_Z(z)={\int }_{-\infty }^{\infty }{f}_X(x)f_Y(z-x)dx$$
这两个公式称为$f_X$和$f_Y$的卷积公式，记为$f_X\ast f_Y$

4.1、随机变量 Z 的概率密度

因为正态分布的随机变量的线性组合仍为正态分布（证明可见概率论与数理统计），所以随机变量 $Z$ 服从均值为 0，方差为 2 的正态分布。

随机变量 $Z$ 的概率密度的三维和二维图像如图 $B$ 和 $C$。
当 $Z=0$ 时，概率密度为 $0.25$ 左右，对应图 $B$ 中过原点，斜率为 0.5 的直线上对应的 $X$ 和 $Y$ 值。对应图 $C$ 中 $Z=0$ 时的概率密度，同样约为 $0.25$。

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.stats import norm# 定义网格范围和步长
x = np.linspace(-4, 4, 200)
y = np.linspace(-4, 4, 200)
x, y = np.meshgrid(x, y)# 假设 X 和 Y 都是独立的标准正态分布 N(0, 1)
mu_X, sigma_X = 0, 1
mu_Y, sigma_Y = 0, 1# 计算联合概率密度 f(x, y)
f_X = norm.pdf(x, mu_X, sigma_X)
f_Y = norm.pdf(y, mu_Y, sigma_Y)
f_XY = f_X * f_Y  # 由于独立性，联合密度等于单独密度的乘积# 计算 Z = X + Y 的概率密度
z = x + y  # Z 的值
f_Z = norm.pdf(z, mu_X + mu_Y, np.sqrt(sigma_X**2 + sigma_Y**2))# 计算边缘分布
marginal_x = np.sum(f_XY, axis=0) * (x[0, 1] - x[0, 0])
marginal_y = np.sum(f_XY, axis=1) * (y[1, 0] - y[0, 0])# 创建一个三维图形对象
fig = plt.figure(figsize=(18, 6))# 绘制联合分布的三维图
ax1 = fig.add_subplot(131, projection='3d')
ax1.plot_surface(x, y, f_XY, cmap='viridis', edgecolor='k', alpha=0.5)
ax1.set_title('Joint Probability Density of X and Y')
ax1.set_xlabel('X')
ax1.set_ylabel('Y')
ax1.set_zlabel('Density')
ax1.view_init(elev=30, azim=45)  # 调整视角# 绘制边缘分布图
ax1.plot(x[0, :], np.full_like(x[0, :], 4), marginal_x, color='red', lw=2, label='Marginal X')
ax1.plot(np.full_like(y[:, 0], 4), y[:, 0], marginal_y, color='blue', lw=2, label='Marginal Y')# 添加图例
ax1.legend()# 绘制 Z = X + Y 的三维图
ax2 = fig.add_subplot(132, projection='3d')
ax2.plot_surface(x, y, f_Z, cmap='plasma', edgecolor='k', alpha=0.5)
ax2.set_title('Probability Density of Z = X + Y')
ax2.set_xlabel('X')
ax2.set_ylabel('Y')
ax2.set_zlabel('Density')
ax2.view_init(elev=45, azim=0)  # 调整视角# 绘制 Z 的边缘概率密度图
z_values = np.linspace(-8, 8, 400)
f_Z_marginal = norm.pdf(z_values, mu_X + mu_Y, np.sqrt(sigma_X**2 + sigma_Y**2))ax3 = fig.add_subplot(133)
ax3.plot(z_values, f_Z_marginal, color='blue', lw=2)
ax3.fill_between(z_values, f_Z_marginal, color='blue', alpha=0.3)
ax3.set_title('Marginal Probability Density of Z = X + Y')
ax3.set_xlabel('Z')
ax3.set_ylabel('Density')plt.tight_layout()
plt.show()

4.2、卷积公式的几何学解释

如上边缘分布的几何学解释，我们已知联合分布函数$(X,Y)$，有 $Z=X+Y$，那么Z的概率密度函数$f_Z(z)$ 如下，也就是卷积公式。

要想知道并理解 $Z$ 的概率密度公式。不妨令 $Z=0$ 来理解这个公式。当 $Z=0$ 时，联合概率密度如下图：

$$f_Z(z)={\int }_{-\infty }^{\infty }{f}_X(-y)f_Y(y)dy={\int }_{-\infty }^{\infty }{f}_X(x)f_Y(-x)dx$$

以 $f_Z(z)={\int }_{-\infty }^{\infty }{f}_X(-y)f_Y(y)dy$ 为例， $f_Y(y)$为标准正态分布， $f_X(-y)$也是标准正态分布，其中 $x=-y$。

即固定了 $y$，$x$也就为 $-y$，$f_X(-y)\ast f_Y(y)$ 的图像就是 $f(x,y)$ 的一个切面，如下图所示。平面方程为 $x+y=0$。

卷积公式就是对两个变量的其中一个做积分，积分结果就是 $Z=0$ 的概率密度。

对 $f_X(-y)f_Y(y)$ 函数在 $x$ 或 $y$ 方向上进行积分，可以理解为将该三维曲线投影到 $x$ 或 $y$ 平面上，然后计算投影曲线与坐标轴围成的面积。

由于 $45°$ 的原因，使得该曲线在 $x$ 或 $y$ 平面的投影面积相同。

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from scipy.stats import norm# 定义网格范围和步长
y = np.linspace(-4, 4, 200)
x = np.linspace(-4, 4, 200)
y, x = np.meshgrid(y, x)  # x 和 y 分别用于计算 f_X(-y) 和 f_Y(y)# 假设 X 和 Y 都是独立的标准正态分布 N(0, 1)
mu_X, sigma_X = 0, 1
mu_Y, sigma_Y = 0, 1# 计算 f_X(-y) 和 f_Y(y)
f_X = norm.pdf(x, mu_X, sigma_X)
f_Y = norm.pdf(y, mu_Y, sigma_Y)# 计算 f_X(-y) * f_Y(y)
f_XY = f_X * f_Y# 仅保留 x + y = 0 的值
tolerance = 0.05
mask = np.abs(x+y) < tolerance
f_XY_filtered = np.where(mask, f_XY, np.nan)
# f_XY_filtered = np.where(mask, f_XY, 0)# 创建一个三维图形对象
fig = plt.figure(figsize=(10, 8))# 绘制 f_X(-y) * f_Y(y) 的三维图（仅当 x + y = 0 时）
ax = fig.add_subplot(111, projection='3d')
surf = ax.plot_surface(x, y, f_XY_filtered, cmap='viridis', edgecolor='k', alpha=0.7)
ax.set_title(r'$f_X(-y) \times f_Y(y)$ (with $x + y = 0$)')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Density')plt.tight_layout()
plt.show()