的分佈X4(X1−X3)+X5(X2−X1)X4(X1−X3)+X5(X2−X1)x_4(x_1-x_3)+x_5(x_2-x_1)獨立同居X一世∼N(0,1)X一世∼ñ(0,1)x_i sim N(0,1)

問題

給定的是 5 個獨立的標準正態變量 x1,x2,x3,x4,x5 .

什麼是pdf x4(x1−x3)+x5(x2−x1) ?

我知道的

x1−x3∼N(0,√2)

x2−x1∼N(0,√2)

x4(x1−x3)∼1π√2K0(|z|√2)

x5(x2−x1)∼1π√2K0(|z|√2)

在哪裡 K0 是 Bessel 函數，eq.(3,4) 是正態乘積分佈。仍然要添加 2 個加法器 x4(x1−x3) 和 x5(x2−x1) 但它們是依賴的，不能卷積。

文獻中的相關問題

如果兩個變量被轉換為 x1(ax1+bx2) 和 a∈R,b∈R+ 那麼pdf是已知的 [1] 但這在這裡不適用。

模擬

仿真表明，分佈 x4(x1−x3)+x5(x2−x1) 可以用帶參數的拉普拉斯分佈來近似 (0,1.34) . 但正確的答案不是拉普拉斯分佈。

---

[1] R. Gaunt：關於零均值相關正態隨機變量乘積分佈的註釋，Statistica Neerlandica，2018

線性代數展示

2(x4(x1−x3)+x5(x2−x1))=(x2−x3+x4+x5)2/4−(−x2+x3+x4+x5)2/4+√3(−√1/3x1+√1/12(x2+x3)+(1/2)(−x4+x5))2−√3(−√1/3x1+√1/12(x2+x3)+(1/2)(x4−x5))2.

每個平方項是獨立標準正態變量的線性組合，按比例縮放以具有方差 1, 這些方塊中的每一個都有一個 χ2(1) 分配。這四個線性組合也是正交的（正如快速檢查所證實的那樣），因此不相關；因為它們是不相關的聯合隨機變量，所以它們是獨立的。

因此，分佈是 (a) 兩個iid差異的一半 χ2(1) 變量加(b) √3 獨立iid差異的一半 χ2(1) 變量。

（iid的差異 χ2(1) 變量具有拉普拉斯分佈，因此這等效地是具有不同方差的兩個獨立拉普拉斯分佈的總和。）

因為一個特徵函數 χ2(1) 變量是

ψ(t)=1√1−2it,

這個分佈的特徵函數是

ψ(t/2)ψ(−t/2)ψ(t√3/2)ψ(−t√3/2)=[(1+t2)(1+3t2)]−1/2.

這不是任何拉普拉斯變量的特徵函數——也不是任何標準統計分佈的 cf。我一直無法找到它的傅里葉逆變換的封閉形式，它與 pdf 成正比。

這是疊加在估計值上的​​公式圖（紅色） ψ 基於 10,000 個值的樣本（實部為黑色，虛部為灰點）：

協議非常好。

---

編輯

仍然存在關於 PDF 是什麼的問題 f 好像。 它可以通過計算對傅里葉變換進行數值反轉來計算

f(x)=12π∫Re−ixtψ(t),dt=12π∫Re−ixt√(1+t2)(1+3t2),dt.

順便說一句，這個表達式完全回答了原來的問題。本節其餘部分的目的是表明它是一個實用的答案。

數值積分一旦成為問題 |x| 超過 10 或者 15, 但稍有耐心就可以準確計算出來。

根據https://stats.stackexchange.com/a/72486/919上對 Gamma 變量差異的分析，很容易通過混合兩個拉普拉斯分佈來近似結果。分佈中間附近的最佳近似值約為 0.4 次拉普拉斯 (1) 更多的 0.6 次拉普拉斯 (√3). 然而，這個近似的尾部有點太重了。

該圖中的左側圖是 100,000 個實現的直方圖 x4(x1−x3)+x5(x2−x1). 在它上面疊加（黑色）的數值計算 f 然後，用紅色表示它的混合近似值。近似值非常好，它與 f. 然而，它並不完美，如右側的相關情節所示。這情節 f 及其在對數尺度上的近似值。尾部近似的精度下降是顯而易見的。

這是一個R用於計算由其特徵函數指定的 PDF 值的函數。它適用於任何數值表現良好的CF（尤其是快速衰減的CF）。

cf <- Vectorize(function(x, psi, lower=-Inf, upper=Inf, ...) {
 g <- function(y) Re(psi(y) * exp(-1i * x * y)) / (2 * pi)
 integrate(g, lower, upper, ...)$value
}, "x")

作為其使用示例，以下是圖中黑色圖形的計算方式。

f <- function(t) ((1 + t^2) * (1 + 3*t^2)) ^ (-1/2)
x <- seq(0, 15), length.out=101)
y <- cf(x, f, rel.tol=1e-12, abs.tol=1e-14, stop.on.error=FALSE, subdivisions=2e3)

該圖是通過連接所有這些來構建的 (x,y) 價值觀。

這個計算為 101 的值 |x| 之間 0 和 15 大約需要一秒鐘。它是大規模可並行化的。

為了獲得更高的準確性，請增加subdivisions參數——但預計計算時間會成比例地增加。（使用的數字subdivisions=1e4。）

引用自：https://stats.stackexchange.com/questions/555736