LLHelper人種鑒定器演算法解析

01-28

1 引言

用過LLhelper的都知道，LLhelper裡面有一個人種鑒定器的功能，輸入你抽卡的成果，可以算出你是歐洲人，亞洲人還是非洲人。判斷人種的機制我在以前的一個答案中介紹過：

Love Live里歐洲人亞洲人非洲人到底該怎麼判定呢？ - 知乎用戶的回答

簡而言之，人種鑒定的過程實際上是計算抽出UR數量的概率分布，而這個分布由所輸入的11連、單抽、歐洲券等分布加和而成。翻譯成數學的語言，就是：

已知離散隨機變數 $X_i$ 的概率分布函數為 $P(X_i) (0leq ileq N_i)$ ，求 $X_1+X_2+...+X_i$ 的累積分布函數在 $x_{UR}-1$ 的值 $F(x_{UR}-1)$ 和概率分布函數的在 $x_{UR}$ 值 $p(x_{UR})$ 。

特別的，當 $x_{UR}=0$ ，累積分布函數不再有意義。

在以前，大家都認為保底機制是當11連出現了11R，會把一張R替換成SR，而這不影響出UR的概率。所以按照老演算法，所有的11連和單抽都滿足二項分布 $B(n, 0.01)$ ，而歐洲券滿足二項分布 $B(n,0.2)$ 。求兩個二項分布的和分布相對比較方便，逐項求和就行了。這就是原先的人種鑒定器演算法。

然而在今年夏天LLhelper的開發者們發現，保底機制不是11R換成10R 1SR，而是11R中的一張有10%的概率變成UR，90%的概率變成SR。這就導致保底期間11連和非保底期間11連的出貨率有著明顯的區別，直接導致原來人種鑒定器作廢。

而更致命的是，保底期間11連出UR的分布不是二項分布，在求總概率分布時，還要涉及保底11連、非保底11連加單抽、歐洲券三種分布的加和。所以人種鑒定器重新上線後，採用了全新的演算法。

2 特徵函數與概率分布

所以為了解決這個問題，我們需要一種能夠很方便地計算多個隨機變數之和分布的方法。這個方法的關鍵就是：特徵函數。

在概率論中，任何隨機變數的特徵函數(縮寫：ch.f,複數形式：ch.fs)完全定義了它的概率分布。在實直線上，它由以下公式給出，其中X是任何具有該分布的隨機變數：
$phi_X(t)=E(e^{itX})$

(引自特徵函數 (概率論) - 維基百科)

打開來寫，就是這個樣子(X的取值為0,1,2...N-1)

$phi_X(t)=sum_{x=0}^{N-1}{p(x)e^{itx}}$

(注意：大寫的 $P(X=a)$ 表示概率，而小寫的 $p(x)$ 表示概率分布函數。)

可以看出，特徵函數實際上就是概率分布函數的傅里葉變換。傅里葉變換具有這樣的性質：兩個函數的卷積的變換，等於它們變換的乘積。而求隨機變數之和的分布，正是卷積運算。所以特徵函數也具有這樣的性質：

多個隨機變數之和的分布的特徵函數，等於它們的特徵函數的乘積。

又因為傅里葉變換是可逆的，只要把特徵函數之積做逆變換，就能得到我們所求的分布。逆變換的做法是：

$p(X)=frac{1}{N}sum_{k=0}^{N-1}{phi_X(frac{2pi k}{N})e^{-frac{2pi ki}{N}}}$

所以，LLhelper人種鑒定器的核心思想就是，通過計算每次抽卡分布的特徵函數並求乘積，得到抽出UR總數的累積分布函數在所實際輸入UR數的值。

3 數學問題

我們現在可以把人種鑒定翻譯成純數學問題了。已知：

一次保底期間11連抽出UR個數 $X_1$ 的概率分布為 $P(X_1)(0 leq X_1 leq 11)$ ，進行保底期間11連的次數為 $n_1$ ；
一次單抽出UR個數 $X_2$ 的概率分布為 $P(X_2=0)=0.99, P(X_2=1)=0.01$ ，單抽和非保底期間11連的總計抽卡數為 $n_2$ ；
一次歐洲券出UR個數 $X_3$ 的概率分布為 $P(X_3=0)=0.8, P(X_3=1)=0.2$ ，抽歐洲券次數為 $n_3$ ；
通過以上途徑一共抽出的UR數為 $x$ 。

求：

$P(n_1X_1+n_2X_2+n_3X_3=x)$ 和 $P(n_1X_1+n_2X_2+n_3X_3<x)$ 。

為了解這個問題，顯然我們需要利用特徵函數的性質。

設 $X_i$ 的特徵函數為 $phi_i(t)$ ，則所有途徑抽出的總UR數的特徵函數為

$phi(t)=phi_1^{n_1}(t)phi_2^{n_2}(t)phi_3^{n_3}(t)$

根據逆變換的的公式，我們可以從特徵函數反過來求概率。

首先考慮x在總分布中的概率：

$p(x)=frac{1}{N}sum_{k=0}^{N-1}{phi(frac{2pi k}{N})e^{-frac{2pi kix}{N}}}$ ，其中 $N=11n_1+n_2+n_3+1$

這個概率可以直接通過特徵函數通過一次求和算出來。

為了求得人品，我們還需要知道累積分布函數在(x-1)的值，即概率分布函數的加和

$F(x-1)=p(0)+p(1)+...+p(x-1)$

為了簡化計算，我們可以把每一個概率的表達式都打開，然後化簡。把兩個求和號交換順序，因為特徵函數和t無關，可以提到對t求和的求和號外面：

$F(x-1)=frac{1}{N}sum_{t=0}^{x-1}{sum_{k=0}^{N-1}{phi(frac{2pi k}{N})e^{-frac{2pi kit}{N}}}} =sum_{k=0}^{N-1}{[phi(frac{2pi k}{N})sum_{t=0}^{x-1}{e^{-frac{2pi kit}{N}}}]}$

不難看出，對t求和的符號裡面是一個等比數列。令 $w=2pi k/N$ ，我們可以直接套求和公式。在k=0時，有

$phi(0) sum_{t=0}^{x-1}{e^{0}}=1cdot x=x$ ，而k>0時，有

$sum_{t=0}^{x-1}{e^{-frac{2pi kit}{N}}}=sum_{t=0}^{x-1}{e^{-iwt}}=frac{1-e^{-iwx}}{1-e^{-iw}}$

根據複數的性質，這個分式還可以進一步化簡。

考慮到 $e^{itheta}-e^{-itheta}=2isintheta$ ，有

$sum_{t=0}^{x-1}{e^{-iwt}}=frac{1-e^{-iwx}}{1-e^{-iw}}=frac{e^{-frac{1}{2}iwx}(e^{frac{1}{2}iwx}-e^{-frac{1}{2}iwx})}{e^{-frac{1}{2}iw}(e^{frac{1}{2}iw}-e^{-frac{1}{2}iw})}=e^{-frac{1}{2}iw(x-1)}cdotfrac{sin frac{1}{2}wx}{sin frac{1}{2}w}$

這樣，我們也能一次求和求出累積分布了。

$F(x-1)=frac{1}{N}sum_{k=0}^{N-1}{[phi(frac{2pi k}{N})e^{-frac{1}{2}i(2pi k/N)(x-1)}cdotfrac{sin (pi kx/N)}{sin (pi k/N)} ]}$

4 演算法實現

現在我們已經得到了計算概率和累積分布的表達式，就可以把它寫成程序實現了。需要計算的兩個值都含有一個求和號，所以顯而易見，程序的主幹是一個for循環，而時間複雜度也自然地是 $O(n)$ 。另外，需要把k=0的第一次循環結果先算出來以避免分母為0，還要用二元數組表示實部和虛部。於是我們有：

這個程序的輸出是：

你的人品已經擊敗了f的玩家，有p的玩家和你人品一樣。

而人種按照如下標準：

rank < 0.01 部落酋長
0.01 <= rank < 0.15 非洲人
0.15 <= rank < 0.35 偏黑亞洲人
0.35 <= rank < 0.65 亞洲人
0.65 <= rank < 0.85 偏白亞洲人
0.85 <= rank < 0.99 歐洲人
rank > 0.99 官托

現在我們需要完成主程序里的幾個函數。主程序里用了這四個函數：

特徵函數獲取：輸入概率分布數組prob和自變數t，返回特徵函數在t的值。

根據特徵函數的表達式 $phi_X(t)=sum_{x=0}^{N-1}{p(x)e^{itx}}$ ，而N為概率分布自變數可取值的總數，即分布數組的長度。我們可以寫出這個函數。

function chf(prob, t) {n if (t==0) {n return [1, 0];n }n var real = 0, imag = 0;n for (var i = 0; i < prob.length; i++) {n real += prob[i] * Math.cos(i*t);n imag += prob[i] * Math.sin(i*t);n }n return [real, imag];n}n

複數相乘：輸入兩個複數，返回它們的乘積。

function times(c1,c2) {n real = c1[0] * c2[0] - c1[1] * c2[1];n imag = c1[0] * c2[1] + c1[1] * c2[0];n return [real, imag];n}n

複數乘方：輸入一個複數和指數n，得到這個複數的n次方。用複數的指數形式計算會比較方便。

function power(c, n) {n if (n == 0) {n return [1,0];n }n logmode = Math.log(c[0] * c[0] + c[1] * c[1]);n mode = Math.exp(logmode * 0.5 * n);n angle = Math.atan2(c[1], c[0]) * n;n return [mode * Math.cos(angle), mode * Math.sin(angle)];n}n

複數指數，輸入幅角theta，得到 $e^{itheta}$ 的值。

function iexp(theta) {n return [Math.cos(theta), Math.sin(theta)];n}n

以上就是LLhelper人種鑒定器的全部演算法。

本文中的代碼和人種鑒定器網頁的javascript代碼不完全一樣，但基本思路和程序結構是幾乎相同的。LLhelper的實際代碼里，各變數的名稱可能跟本文有些區別，並且會有一些數學處理上的差異。比如沒有用Math.log函數，而是math.log1p函數等。另外LLhelper代碼里還有一些冗餘的預算，比如在函數中使用了分布的期望值，而期望值在運算過程中會被抵消。筆者也不明白進行這些運算的原因，並且基於原程序，去掉這些運算之後，結果和運行時間不受影響。

實際上，本文提供的代碼可以看作是人種鑒定器的重構。把本文的代碼替換入人種鑒定器的網頁(輸入輸出變數名也要替換)，仍然可以執行出正確的結果。為了讓演算法的思路更清晰，我沒有進一步再優化語句。所以執行的時間會比LLhelper的網站里略長。

5 演算法優化

顯然，這個演算法對於總抽卡數量N，需要進行N次循環。那麼我們有沒有辦法讓循環的次數減少，從而使計算快一些呢？

有。我們需要用到特徵函數的對稱性。回到特徵函數的表達式：

$phi_X(t)=sum_{x=0}^{N-1}{p(x)e^{itx}}$

我們知道， $e^{itx}=cos tx + isin tx$ 。在主程序的循環里，我們特徵函數的變數t按照等間隔取遍了一個2π的周期。在一個周期內，顯然三角函數具有對稱性。我們不妨求一下 $phi(2pi-t)$ ：

$phi(2pi-t)=sum_{x=0}^{N-1}{p(x)e^{it(2pi-x)}}=sum_{x=0}^{N-1}{[p(x)e^{2pi it}e^{-itx}]}$

在t為整數的時候，有 $e^{2pi it}=1, e^{-itx}=(e^{itx})^*$ (星號表示共軛複數)

而再看逆變換中的一次循環，同樣可把循環節中w用(2π-w)代入，有

$phi(2pi-w)e^{-(2pi - w)x}=phi(w)^* (e^{-wx})^* =(phi(w)e^{-wx})^*$

也就是說，對關於w=π對稱的兩個w，得到每一個求和項的實部是一樣的！對於累積分布函數，因為每一個概率都取共軛，其結果也自然是共軛。那麼我們就可以把求和的項數減少一半了。主程序中的for循環只要取前一半，然後再單獨把w=π的情況補上就行了。這樣優化後的主程序代碼如下：

for (var k = 1; k <= N / 2; k++) {n var w = (2*Math.PI*k)/N;n var c = [1,0];n if (n1 > 0) {n var c1 = power(chf(prob1, w), n1); /*n1次保底11連的特徵函數*/n c = times(c, c1);n }n if (n2 > 0) {n var c2 = power(chf(prob2, w), n2);/*n2次單抽的特徵函數*/n c = times(c, c2);n }n if (n3 > 0) {n var c3 = power(chf(prob3, w), n3);/*n2次單抽的特徵函數*/n c = times(c, c3);n }n if (k == N/2) {n f += times(c, iexp(-0.5*w*(x-1) ))[0]*Math.sin(0.5*w*x)/Math.sin(0.5*w);n p += times(c, iexp(-1*w*x))[0]; n }ttn else {n f += 2*times(c, iexp(-0.5*w*(x-1) ))[0]*Math.sin(0.5*w*x)/Math.sin(0.5*w);n p += 2*times(c, iexp(-1*w*x))[0]; n } n}n

使用上述優化後的演算法，可以顯著減少人種鑒定器的計算用時。

6 總結

現在我們知道了新版的人種鑒定器裡面，你的人種是如何得出的原理。然而即使我費了很大功夫讀懂了人種鑒定器，也改變不了我非洲人的命運。