Probability

允許計算黃金比例的統計實驗有哪些例子?

  • March 15, 2021

有一些非常簡單的經驗可以由一個孩子在家裡完成,其結果可以讓一個人在統計上接近著名的數字,例如 π 或者 e .

一個例子 π 出現可能是同類中最著名的一個。在Buffon 的針問題中,我們在地板上畫條並放下一根針。針穿過兩條帶之間的線的概率包括 π . 多次重複這個過程可以讓我們接近 π 準確地說,我們是否願意將經驗重複足夠多次。

一個例子 e 出現包括繪製一個隨機大小的樣本 n 從一定規模的人口中替換 n . 未選擇總體成員的概率為 p=(11/n)n . 如果 n 然後 p1/e .

我的問題是,有哪些實驗示例可以讓人們在統計上接近黃金比例的值 Φ=(1+5)/2=1.618033 ? 或者換句話說,如何接近 Φ 通過蒙特卡羅模擬。

(有一個條件是不能對經驗進行微調以獲得結果。例如,如果我們在地板上畫一個輪廓,並以某種黃金比例將其分成兩部分,然後隨機向其扔石頭,我們顯然可以通過計算每個部分落下的石頭的數量來恢復黃金比例。我要求以更意想不到的方式出現結果的例子。)

麻煩先生只有一個。事實上,除非他有一個兒子可以繼承他的姓氏,否則他將是最後一個 三胞胎先生。可悲的是,他那個時代和地方的社會習俗甚至不允許這樣一個精緻的姓氏通過女係而生存。

Tripletoddletrouble 先生有一種罕見且在數學上很方便的遺傳條件,Tripletoddletroubles 的任何後代都將繼承這一遺傳條件:如果他有任何兒子,而且他有 50% 的機會,這將是一組三胞胎。因此,在家族樹的每一步,我們同樣可能看到三個兒子或一個兒子。

對名稱學家來說,壞消息 是,一個美妙的名字有 50% 的機會在下一代中滅絕。好消息是,每個三胞胎先生平均有 1.5 個兒子——因為這安全地高於一個,三胞胎麻煩的預期人口呈指數增長,而且他們的姓氏有可能永遠存在。

Tripletoddletrouble 姓氏最終滅絕的概率是多少?

這是一個快速R模擬。

set.seed(123)

nsims <- 1e5
ngens <- 20

simulate_extinction <- function(ngens) {
 nsurvivors <- c(1, rep(NA, ngens - 1))
 for (gen in seq_len(ngens - 1)) {
   nsurvivors[gen + 1] <- 3 * rbinom(1, nsurvivors[gen], 0.5)
 }
 extinct <- (!is.na(nsurvivors) & nsurvivors == 0) # rbinom gives NA if population huge
 return(extinct)
}

pextinct <- rowMeans(replicate(n = nsims, simulate_extinction(ngens)))

plot(pextinct, xlab = "Generation number", ylab = "Probability of extinction")
abline(h = (sqrt(5) - 1)/2, col = "red")
sprintf("Estimated probability of extinction = %f", pextinct[ngens])

圖中的紅線位於

φ1=φ1=5120.618034

姓氏滅絕概率圖接近 1/phi

[1] "Estimated probability of extinction = 0.618150"

這是一個關於分支流程的問題。事實上,對他們的隨機行為的最早調查之一起源於維多利亞時代對不尋常姓氏滅絕的擔憂。生成的Galton-Watson 工藝文件可在線獲取:

**“3個或0個後代,同樣可能”可以說是最簡單的分支過程,具有非平凡的滅絕概率。**如果機會發揮作用,我們至少需要兩個結果,包括零後代才能使滅絕成為可能。“非一子一零”顯然是注定的:沒有枝條萌芽,第一次沒有兒子的家族就滅絕了。“2 或 0 個後代,同樣可能”給出了恰好一個後代來替換每個個體的平均值。當命運在這個刀刃上取得平衡時,從長遠來看,即使家譜成功發芽幾次,滅絕也是必然的。我們可以調整後代分佈以產生一系列所需的滅絕概率,但只能通過引入不相等的概率或兩個以上的結果。這種設置並沒有人為地“調整”到黃金比例中的鞋拔。

讓我們找到最終滅絕的概率, θ ,代數。直觀地說,這個概率分為兩部分:或者最初的三胞胎先生沒有兒子,他的血統立即滅絕,或者他成功地擁有了三個兒子,但他們的三個血統中的每一個最終都滅絕了。由於一個兒子和原來的 Mr Tripletoddletrouble 處於同一位置,他們的譜係也各有滅絕概率 θ . 由於我們只關心直系男性後裔,因此每一系的命運都是獨立於其他的。鑑於有三個兒子,因此姓氏滅絕的概率為 θ3 .

姓氏滅絕分支過程的樹形圖。 要么有 0 個兒子(立即滅絕),概率為 0.5,要么有 3 個兒子,概率為 0.5。 如果是後者,所有 3 行滅絕的概率是 theta^3,或者至少有一個倖存的概率是 1-theta^3。

從樹形圖中,我們看到了滅絕的概率 θ 必須服從方程

θ=12+12θ3

我們可以解決(很快就會解決)。首先,讓我們將其與一些更廣泛的分支過程理論聯繫起來。任何個體的後代數量是具有概率分佈的隨機變量 p0=p3=0.5pn=0 否則,其概率生成函數為:

Π(s)=npnsn=12+12s3

看起來很熟悉?沒有巧合。以後再說。通過評估 pgf 在 s=1 ,我們得到後代的平均數。這個號碼 R0=Π(1) 在人口生態學和人類人口學中很重要,它被稱為淨繁殖率(它通常被定義為每個女性生育的女兒的平均數量,而不是每個男性生育的兒子的平均數量 - 母性比父性更容易追踪,並且在許多物種中雌性可以通過孤雌生殖進行繁殖),而在流行病學中,它是基本繁殖數(在完全易感人群中,一個受感染個體直接產生的平均感染數)。如果 Π(1)1 那麼最終的滅絕是肯定的。如果 Π(1)>1 滅絕的概率低於一。我們有

Π(s)=3s22Π(1)=1.5>1

所以姓氏的生存概率為正。每一代能存活*多少個三胞胎麻煩?*以一個人為“零代”,讓 Zn 是後裔的數量 n 幾代人。 Zn 是一個隨機變量,其概率分佈可以從其 pgf 的係數中讀出 Πn(s) , 我們通過迭代應用找到 Π , 後代 pgf, n 次:

Πn(s)=Π(Π(Π(s)))

為什麼? Zn 是後代的總和 Zn1 上一代的倖存者。每個倖存者的後代數量是獨立的,同分佈(iid)隨機變量,帶有 pgf Π , 我們加起來的數量有 pgf Πn1 ,因此根據隨機數 iid 變量之和的 pgf 規則(此答案中的證明), Zn 有 pgf Πn(s)=Πn1(Π(s)) . 例如,經過兩代

Π2(s)=Π(Π(s))=12+12(12+s32)3=916+3s316+3s616+s916

所以有一個 116 九個後代的機會,但 Π2(0)=916 滅絕已經發生的可能性。 E(Z2) ,兩代後的預期後代數,由 Π2(1)=2.25 . 這不是巧合 1.52 .

單個個體的後代數量的均值和方差分別為 μ=Π(1)σ2=Π(1)+μμ2 . 你可以通過歸納證明 E(Zn)=μn . 現在很明顯為什麼最終滅絕是肯定的 μ<1 . 和 μ=1.5 儘管我們很可能提前滅絕,但我們平均看到了指數級增長。從本質上講,姓氏傳播鏈要么失敗,要么崩潰,而且 μ=1.5 保證有足夠的機會爆發滅絕並非不可避免。Tripletoddletroubles 的好消息;如果我們將上下文從姓氏切換到傳染病 R0>1 . 感染鏈可以隨機“變大或回家”的方式,而不是遵循諸如“每個病例恰好感染兩個易感者”之類的確定性規則,這與由於聚集或超級傳播事件導致的過度分散的流行病學思想有關。後後代數量的方差 n 代數可能相當大,因為 Zn 可能是巨大的或零。再次通過歸納,我們發現:

Var(Zn)={μn1σ2(μn1)μ1,μ1\[2ex]nσ2,μ=1

一般來說,最終滅絕的概率是最小的正解θ , 等式 θ=Π(θ) . 這正是等式 (1) 我們在上面得出!但是我們怎麼知道要採取哪種解決方案呢?世代滅絕的概率 nΠn(0) ,因為這是常數或 s0 pgf 的期限 Zn , 因此表示 Pr(Zn=0) . 最終滅絕的概率一定是 limnΠn(0) 我們可以使用蛛網圖找到 y=Π(x)y=x 為了 0x1 . 自從 Π(0)=p0 ,一個人沒有後代的概率,我們可以假設 y 截距在 0<Π(0)1 (如果 p0=0 那麼滅絕顯然是不可能的)。所以 y=Π(x) 從上面開始 y=x ,並且第一次相交 y=x 必須從上面。自從 Π(x) 並且它的導數只有非負係數,它的圖是遞增和凸的 0x1 . 這意味著它可以相交 y=x 在此間隔內最多兩次:一次從上方,然後再次從下方。 Π(1)=pn=1 所以這些圖肯定在 (1,1) .

這個交叉點的性質取決於坡度 Π(1) , 代表後代的平均數 μ (生物學上, R0 )。如果 Π(1)>1 它必須比 y=x 所以 y=Π(x) 正在從下方擊中線,在這種情況下,一定有一個較早的交叉點 0<x<1 . 如果 Π(1)<1 它從上面打得更淺,而且沒有更早的根。如果 Π(1)=1 兩條曲線剛剛接觸 (1,1) , 但 y=Π(x) 之前一定更淺(它的平均斜率超過 0x11p0 所以低於一個),因此從上面接近這條線並且不可能有更早的根。這就是為什麼如果 μ=1p0>0 ,最終滅絕的概率為一。

尋找 limnΠn(0) 以圖形方式,從 y 截距處水平讀取 y=Π(0)y=x 線,現在在哪裡 x=Π(0) . 然後垂直讀取 y=Π(x) 圖,現在在哪裡 y=Π(Π(x))=Π2(x) . 水平讀取到該行 x=Π2(x) . 垂直於曲線讀取 y=Π(Π2(x))=Π3(x) . 請注意,所有水平讀數都向右,垂直讀數向上,因為 y=Π(x) 正在增加,因此每個垂直位置都高於前一個。這個過程必須收斂到第一個(即最小的正 x ) 交點 x=Π(x) , 在哪裡 y=Π(x) 命中 y=x 從上面。我們舉例說明三種情況 μ=1.5 ( p0=p3=12 ), μ=1 ( p0=23,p3=13 ) 和 μ=0.5 ( p0=56,p3=16 )。藍色虛線是切線 y=Π(x)(1,1) ,並顯示其斜率的作用 Π(1)=μ 在確定是否有較早的交叉路口。

後代 pgf 的蛛網圖顯示最終滅絕概率的收斂

我們需要最小的正解 θ(1) . 移動 θ 到右邊並加倍以清除分數,我們得到:

0=θ32θ+1=(θ1)(θ2+θ1)

解決方案是 φ<φ1<1 所以最小的正解是 θ=φ1 .

**是時候揭開“小提琴”的面紗了。**這個與黃金比例的聯繫不是我記得以前見過的結果,但我通過考慮最終方程所需的因式分解對其進行了逆向工程。自從 Π(1)=pn=1 ,我們總是有 θ=1 作為一個根 θ=Π(θ) , 所以 (θ1) 一旦我們將一側設置為零,它就必須作為一個因素出現。我也知道我想看到什麼二次方。之後,我回過頭來嘗試形成一個有效的 pgf。不允許使用負係數;正係數只需要歸一化,因此它們總和為單位。我希望由此產生的後代概率分佈會是一個“不錯”的分佈——我認為它是!

蜘蛛網圖的R代碼

ngens <- 100
par(mfrow=c(1, 3), pty = "s", xaxs = "i", yaxs = "i")

for(p0 in c(1/2, 2/3, 5/6)) {
 pgf <- function(x) {p0 + (1-p0)*x^3}
 mu <- p0*0 + (1-p0)*3
 plot(pgf, xlim = c(0,1), ylim = c(0,1), xlab = "", ylab = "", 
      main = paste0("Mean offspring = ", mu))
 segments(0, 0, 1, 1)
 abline(1 - mu, mu, col = "blue", lty = "dotted")
 pextinct <- c(0, rep(NA, ngens))
 for (n in seq_len(ngens)) {
   pextinct[n + 1] <- pgf(pextinct[n])
   segments(pextinct[n], pextinct[n], pextinct[n], pextinct[n + 1], col = "red")
   segments(pextinct[n], pextinct[n + 1], pextinct[n + 1], pextinct[n + 1], col = "red")
 }
 print(pextinct)
}

引用自:https://stats.stackexchange.com/questions/514045