在給定的最小值和最大值下生成服從指數分佈的隨機樣本

服從指數分佈的隨機樣本可以通過使用指數分佈的分位數函數通過逆採樣技術生成：

x=F−1(u)=−1λln(u)

在哪裡 u 是從單位區間上的均勻分佈中抽取的樣本 (0,1) .

在OpenFOAM軟件中，可以使用稱為exponential( here ) 的分佈模型來生成指數分佈的隨機樣本，據推測，它的用戶可以在隨機數生成之前為指數分佈的樣本選擇最小值和最大值。

本軟件實現的控製表達式如下：

x=−1λln[exp(−λtmin)+uexp(−λtmax)−exp(−λtmin)]

在哪裡 tmin 和 tmax 分別是用戶定義的最小值和最大值。

出現了兩個（相當廣泛的）問題：

1. 您是否曾經在文獻中遇到過上述表達式（或類似的表達式），用於生成具有給定最小值-最大值的指數分佈隨機樣本？或者你認為這個表達式看起來像（或者是）一個啟發式解決方案？
2. 如果是啟發式的，您會建議一種對此表達式進行驗證測試的方法，以測試該表達式是否產生服從 [min,max] 內的指數分佈的樣本？（繪製歸一化直方圖（即在 bin 中的計數除以觀察次數乘以 bin 寬度）並將其與分析指數分佈進行比較似乎是有問題的，因為有最小/最大限制）。

您將截斷描述為間隔。我會詳細說明。

認為 X 是任何隨機變量（例如指數變量）並且讓 FX 是它的分佈函數，

FX(x)=Pr(X≤x).

對於一個區間 [a,b], 截斷限制_ X 到那個區間。這減少了一些概率 X, 即機會 X 要么小於 a 或大於 b. 剩下的機會是

Pr(X∈[a,b])=Pr(X≤b)−Pr(X≤a)+Pr(X=a)=FX(b)−FX(a)+Pr(X=a).

因此，要使總概率為 1, 截斷的分佈函數 X 必須為零時 x<a, 1 什麼時候 x≥b, 否則是

FX(x;a,b)=Pr(X∈[a,x])Pr(X∈[a,b])=FX(x)−FX(a)+Pr(X=a)FX(b)−FX(a)+Pr(X=a).

當你可以計算分佈函數的倒數時——這幾乎總是意味著 X 是一個連續變量——生成樣本很簡單：繪製一個統一的隨機概率 U （從區間 [0,1], 當然）並找到一個數字 x 為此 FX(x)=U. 這個值被寫入

x=F−1X(U).

F−1X 稱為“百分點函數”或“逆分佈函數”。

例如，當 X 有一個帶速率的指數分佈 λ>0,

U=FX(x)=1−exp(−λx),

我們可以解決得到

F−1X(U)=−1λlog(U).

這稱為“反轉分佈”或“應用百分比函數”。

事實證明——這就是這篇文章的重點——當你可以反轉 FX, 您還可以反轉截斷分佈。 給定 U, 這相當於解決

U=FX(x;a,b)=FX(x)−FX(a)FX(b)−FX(a),

因為（因為我們現在假設 X 是連續的）條款 Pr(X=a)=0 退出。解決方案是

x=F−1X(U;a,b)=F−1X(FX(a)+[FX(b)−FX(a)]U).

也就是，唯一的變化就是畫完之後 U, 您必須重新縮放並移動它以使其值介於 FX(a) 和 FX(b), 然後你反轉它*。*

這產生了問題中的第二個公式。

一個等效的過程是繪製一個統一的值 V 從區間 [FX(a),FX(b)] 併計算 F−1X(V). 這是有效的，因為縮放和移位的版本 U 在這個區間內分佈均勻。我在下面的代碼中使用了這種方法。

該圖說明了該算法的結果 λ=1/2 並截斷到區間 [2,7]. 我認為僅此一項就可以很好地驗證程序。

R代碼是通用的：替換（ff實現 FX ) 和f.inv(實現 F−1X ) 與任何連續隨機變量的相應函數。

#
# Provide a CDF and its percentage point function.
#
lambda <- 1/2
ff <- function(x) pexp(x, lambda)
f.inv <- function(q) qexp(q, lambda)
#
# Specify the interval of truncation.
#
a <- 2
b <- 7
#
# Simulate data and truncated data.
#
n <- 1e6
x <- f.inv(runif(n))
x.trunc <- f.inv(runif(n, ff(a), ff(b)))
#
# Draw histograms.
#
dx <- (b - a) / 25
bins <- seq(a - ceiling((a - min(x))/dx)*dx, max(x)+dx, by=dx)

h <- hist(x.trunc, breaks=bins, plot=FALSE)
hist(x, breaks=bins, freq=FALSE, ylim=c(0, max(h$density)), col="#e0e0e0",
    xlab="Value", 
    main="Histogram of X and its truncated version")
plot(h, add=TRUE, freq=FALSE, col="#2020ff40")

abline(v = c(a,b), lty=3, lwd=2)
mtext(c(expression(a), expression(b)), at = c(a, b), side=1, line=0.25)

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