1. 從概率預測到決策函式之間的鴻溝
在統計學領域，統計學家常常不考慮能獲得多少，而是考慮損失了多少，我們認為損失的負數就是獲得。如何衡量損失是一件值得深入思考的事情。
0x1：例子：只考慮預測結果不考慮不同預測結果背後的實際損失帶來的問題
考慮如下例子：
假設一個氣象學家正在預測颶風襲擊他所在城市的概率是多少。他估計颶風不襲擊該城市的可能性為 99%-100%，且這一概率的可信度為 95%。氣象學家非常滿意他的計算精度，並建議大家不用疏散。不幸的是颶風襲擊了城市，城市被淹沒了。
這個典型的例子說明了單純使用精度去度量結果時存在的缺陷，很顯然，這是我們很多開發朋友包括筆者自己平時在專案中無時不刻遇到的情況。
我們很幸苦收集了大量打標樣本，設計了一套特徵工程以及演算法模型後，又經過了很久的訓練得到了一個預測模型，這時候，預測模型可以得到幾千到幾萬不等的預測結果，其中置信度概率從0.8 ~ 0.99不等，這個時候問題來了，我們要拿哪些結果進行輸出，例如輸出給下游業務方或者輸出給客戶？
是 >= 0.99？那麼是否意味著漏報了很多資料？
是 >= 0.9？那麼是否意味著產生很多誤報？
退一步說， >= 0.99 就一定是100%沒誤報的嗎？我們設定的這個所謂的 D（>= 0.99）的決策函式，100%不會出問題嗎？
造成這個問題的根本原因是什麼呢？
筆者認為根本原因是：我們對損失函式的定義還不夠完善，如本篇blog的tilte所說，損失函式的本質是”連線統計推斷和問題域的橋樑“，所以損失函式就不能單純只考慮統計推斷的概率結果，還要考慮問題域的實際情況，即不同的決策結果帶來的實際現實世界的損失，損失函式必須是承上啟下的。
0x2：大致的正確比精確的錯誤更好
從貝葉斯推斷的世界觀來看，我們不能太過於強調精確度的度量，儘管它往往是一個有吸引力和客觀的度量，但會忽視了執行這項統計推斷的初衷，那就是：推斷的結果。
此外，我們希望有一種強調“決策收益的重要性”的方法，而不僅僅是估計的精確度。
0x3：以尋找最優結果為目標的損失評價函式
在機器學習專案中，我們更在意的是得到最佳的預測結果，而不僅僅是在所有可能的引數下達到最佳精度。
尋找貝葉斯行動等價於尋找一些引數，這些引數優化的不是引數精度，而是任意某種表現。不過我們需要先定義表現（損失函式、AUC、ROC、準確率 / 召回率）
Relevant Link:
https://github.com/CamDavidsonPilon/Probabilistic-Programming-and-Bayesian-Methods-for-Hackers
2. 損失函式
在開始討論這章之前，筆者希望先丟擲一句話：損失函式是連線統計推斷和問題所在領域的橋樑。我們整篇文章會圍繞這個觀點展開討論。
損失函式不僅僅侷限於我們在機器學習和深度學習書籍裡看到的最小二乘損失、交叉熵損失、0-1損失等等。實際上，損失函式是可以根據實際情況，在具體專案中進行自定義的定義，它的核心作用是將現實世界裡的問題通過數學方式數字化表達，讓損失這個抽象的邏輯概念可以數值化並通過優化演算法進行運算和優化。
0x1：損失函式的定義
我們先來學習統計學和決策理論中的損失函式。損失函式是一個關於真實引數以及對該引數的估計的函式：
損失函式的重要性在於，他們能夠衡量我們的估計的好壞：損失越大，那麼損失函式得到的估計就越差。
一個簡單而普遍的例子是平方誤差損失函式，這是一種典型的，與誤差的平方成正比的損失函式。線上性迴歸、無偏差統計量計算、機器學習的許多領域都有廣泛的應
0x2：常見的損失函式
1. 0-1損失函式
0-1損失非常直觀和簡單，非對即錯，是一種階躍函式，在實際專案很少使用到。
2. 平方損失函式 - 非線性損失評價函式
平方損失誤差損失的缺點在於它過於強調大的異常值。這是因為隨著估計值的偏離，損失是平方增加的，而不是線性增加的。
對3個單位偏離的懲罰要小於對5個單位偏差的懲罰；但是對3個單位偏離的懲罰卻不比對1個單位偏離的懲罰大很多，雖然在一些場景下，這兩種情況下偏離的差值是相同的。
均方損失的這種不線性，意味著較大的誤差會導致糟糕的結果。在實際專案中，是否可以使用需要資料科學家自己根據經驗去權衡。
如果我們的問題是迴歸分析的場景，我們是希望擬合模型能儘量地概括資料的分佈，離得越遠就意味著迴歸擬合能力越差，因為懲罰也越大，因此可能均方損失的這種非線性特性反而是合理的。
因此要選擇損失函式的時候，需要我們對問題場景要非常的瞭解，知道選擇不同的損失函式對最終獲得的結果意味著什麼。
3. 對數損失函式（交叉熵損失函式）
交叉熵損失函式基於資訊理論中的熵值概念，將損失看成是當前系統不確定性的度量，預測值和目標值差距越大，意味著不確定性越大，熵值也就越大，損失也就越大。
4. 絕對值損失函式 - 強調偏離度的線性損失評價函式
一種更穩健的損失函式是誤差的絕對值函式，公式如下：
可以看到，絕對值損失函式不考慮預測值和目標值是”左偏“還是”右偏“，它只是線性地衡量”預測值和目標值的距離“。
絕對值損失函式是機器學習和穩定統計中經常使用的損失函式。
筆者插入：可以很容易地看到，絕對損失函式沒有對偏離的正負進行區分。but！這個重要嗎？截止目前為止，這個問題筆者很難作出回答，我自己也在摸索中，從筆者目前為止為數不多的專案經驗來看，兩種情況都會遇到，你需要自己靈活判斷。
假如你需要的預測的一臺機器的惡意網路發包行為的可疑程度，我們對網路established行為進行建模，進行迴歸分析，迴歸模型擬合的值其實是一個閾值，低於或者高於這個迴歸值的行為模式都被認為是有一定可疑程度的。但是很顯然，從安全業務場景裡去定義這個問題時，低於這個閾值被認為是可疑程度相對較低的，高於這個閾值的被認為是可疑程度相對較高的，這個時候，絕對值損失函式就不能100%客觀地表徵我對這個業務場景的理解了，我可能需要將絕對值損失函式進行改造，改為分段函式，並在負分段和正分段分別賦予不同的懲罰因子。
我們在下個小節會討論的分段損失函式，可以緩解這個問題。
0x3：面向具體業務場景的自定義損失函式
從筆者的理解上來看，目前資料科學家和資料工程師們選擇損失函式主要基於如下兩方面：
1. 數學計算上的便利性，由於計算機在數學上的便利性，我們可以自由地設計自己的損失函式； 2. 應用上的穩健性，即損失函式可以客觀的對損失進行度量。損失函式確實是客觀的，它們都是估計值和真值之差的函式，不管該誤差是正還是負，和估計的收益物件無關； 3. 損失函式是”可優化的“，這涉及到凸優化理論，理論上說，只有凸函式可以通過計算進行優化，但是值得注意的是，現在對非凸函式的優化也有大量的研究，有興趣的朋友可以google相關資料；
我們文章開頭說過，損失函式的定義需要承上啟下，既要能夠對概率問題進行數值化，也要考慮實際問題域的情況。這就要求我們對常規損失函式進行定製。
再次考慮以上颶風的例子，統計學家也可以這樣預測：
颶風襲擊的概率在 0%到 1%之間，但是如果發生颶風的損失是巨大的。99%的概率沒有颶風，沒有颶風帶來的收益有限。基於這種損失評價體系，統計學家最後給出的建議可能截然不同。
很顯然，這種改動意味著我們把關注的重心從更精確的概率估計轉到概率估計帶來的現實結果上來，即損失函式要更加關注決策結果和業務收益。
1. 分段區間非對稱平方誤差損失函式
設想一種場景，預測值低於目標值造成的損失，小於預測值高於目標值。我們可以用一種非對稱的均方誤差損失函式來描述這種關係：
上式表明，這種損失函式對預測值小於目標值的情況，懲罰更嚴重。在這種損失函式的驅動下，模型會盡量避免讓預測值小於目標值。
這種損失函式形式更適合於估計下個月的網路流量，因為為了避免伺服器資源分配不足情況，在對伺服器資源分配進行估計時要多分配一些。
從這個例子可以看到，損失函式並不一定都是簡單的說去衡量模型的預測值和目標值之間的數值差異，而是更加貼近我們的業務目標。損失函式衡量的是我們的模型預測結果對我們的業務目標的貢獻程度。
2. 強調引數估計結果更接近 0 或者 1 的損失函式
該損失函式更傾向於估計靠近0或者1位置的引數 ?。或者理解為噹噹目標值為靠近0或者1時，對預測值的一點點偏離都會產生比較大的懲罰。
因為如果真值?靠近 0 或者 1 時，損失函式將變的非常大，除非估計值?̂也 接近 0 或者 1。
基於這種損失函式進行引數優化，模型得到的引數是傾向於靠近 0 或者 1 的。
這種損失函式更適合政治評論家使用，因為他們更傾向於給出 “是/否”的答案（即目標值為0或者1）。
3. 傾向於誤差小的估計的損失函式
該損失函式的取值範圍為 0 到 1，這種損失函式表明其使用者不太關心誤差大的估計結果，即預測值和目標值偏離程度並不會顯著導致損失增加，這和 0-1 損失函式比較類似，但是在真值附件的懲罰程度又沒有 0-1 損失函式那麼強烈。0-1損失的懲罰太階躍了。
4. 帶有對某個結果傾向性的損失函式
天氣預報人員使用的損失函式常常是帶有明顯的傾向性的。氣象預報員有很強的動力盡量準確地預報下雨的概率，也同樣有動力去錯誤地暗示可能有雨。為什麼會出現這種傾向性？
人們更喜歡有所準備，即使可能不會下雨，也不喜歡下雨造成的措手不及。出於這個原因，預報員傾向於人為地增加降雨概率和誇大這一估計，因為這能帶來更大“收益”。
這個和我們在開頭討論的洪水預測有類似的意思，正是中國那句老話：寧可信其有，不可信其無。
3. 貝葉斯推斷理論體系下的損失函式
0x1：樣本集資料驅動下的經驗損失函式
需要明白的是，我們之前在討論損失函式時，都假設我們知道了真值（目標值），但是實際情況下，損失函式理論下的引數的真實值基本上是不可能已知的，否則，如果我們已經知道真實值，也就沒有必要再去估計了。討論和真值之間的損失只是為我們提供一個理論的上界，作為一個benchmark參考系。
因此不管是在傳統機器學習/深度學習中，還是貝葉斯推斷中，損失函式都不是針對真是引數值的，而是針對訓練樣本集的，即經驗損失函式。
在這一點上，貝葉斯統計推斷和機器學習是一樣的，本質上都是訓練樣本資料驅動的。
0x2：貝葉斯後驗分佈估計與貝葉斯點估計 - 從概率分佈和最大後驗估計的角度看後驗分佈結果
在貝葉斯推斷中，把未知引數看成是與先驗和後驗分佈的隨機變數。貝葉斯推斷不會給出一個具體的估計值，而是該出一個估計值的後驗分佈。
就後驗分佈而言，從後驗分佈得到的樣本都可能是真實的引數值（真實的引數值可能就是分佈中的某個點）。
而貝葉斯點估計是什麼呢？我們可以理解為在貝葉斯後驗估計的概率分佈空間中取一個點作為輸出結果。這麼講有點抽象，舉一個例子說明。
一維高斯函式：
上圖中的曲線可以理解為對均值和方差【u，σ】的貝葉斯後驗估計，它是一個後驗概率分佈。
而我們根據最大後驗似然估計的原則，選擇【u，σ】作為這個分佈的最終輸出，得到的【u，σ】就是一個貝葉斯估計點。
0x3：貝葉斯估計的期望損失 - 貝葉斯估計損失的理論值
當我們知道了未知引數的後驗分佈，就可以基於樣本集，計算某種引數估計相關的損失函式了。
但是問題是後驗分佈的損失函式還是一個分佈，我們無法做出決策，我們需要一個確定的結果或者實值，例如在很多專案中，當我們通過ML得到一個預測概率值的時候，這時我們面臨一種不確定性時，我們往往會通過決策函式（例如一個經驗閾值）來將不確定性提煉成一個單獨的確定性決策動作（0/1）。
我們需要提煉我們的後驗分佈為單一的值（或向量（在多變數的情況下））。如果後驗分佈的值取得合理，可以避免頻率學派無法提供不確定性的缺陷，同時這種結果的資訊更豐富。
回到貝葉斯估計得到的後驗估計話題，比起貝葉斯後驗估計，我們更關心的是在某種估計下的期望損失，即一個確定的實值。
這麼做的原因有很多，首先期望損失有利於解釋貝葉斯點估計。現實世界中的系統和機器無法把先驗分佈作為輸入引數。實際上最終關心的是估計結果，後驗分佈只是計算估計結果的中間步驟， 如果直接給出後驗分佈是沒有太大作用的。
與之類似，我們需要使後驗分佈變得陡峭，理想情況是變成一個點。如果後驗分佈的值取得合理，可以避免頻率學派無法提供不確定性的缺陷，這種結果的資訊更豐富。
從貝葉斯後驗中選擇點， 就是貝葉斯點估計，也即求期望損失。
假設 ?(?|?) 是在觀測資料 X 的條件下 ? 的後驗概率，? 的期望損失 
為
從後驗分佈得到 ? 個樣本 ?? , ? = 1, ⋯ , ?， 損失函式為?，可以通過大數定理，利用 ?̂ 的估計值近似計算出期望損失
筆者插入：這裡要注意，期望損失計算的物件貝葉斯估計的理論值，在實際工程中，我們是沒有上帝視角，我們只有訓練樣本，因此基於訓練樣本得到的貝葉斯後驗分佈以及貝葉斯估計結果，本身就包含了一定的偏置（bias），我們將基於樣本得到的貝葉斯估計結果成為經驗損失。
0x2：貝葉斯點估計相比於頻率派估計的優勢
注意，與 MAP 相比，計算損失函式的期望值利用了後驗分佈中的更多資訊，因為在 MAP 中僅僅是尋找後驗分佈的最大值，而忽略了後驗分佈的形狀。
MAP 中忽略掉的資訊可能使你暴露在長尾部分的風險當中，像颶風這種可能性很小但是存在的風險卻很巨大，MAP 將導致估計結果無視引數的未知性。
類似地，頻率學派的目標只是最小化化誤差，不考慮與誤差結果相關的損失。頻率派方法幾乎可以保證永遠不會絕對地準確（我們都是預測模型幾乎很少能得到100%的預測結果）。
貝葉斯點估計方法是通過提前計劃解決這個問題：如果你的估計將是錯誤的，那還不如以正確的姿態犯錯 --- 模糊地正確性勝過精確的錯誤。
4. 貝葉斯統計例項：優化“價格競猜”遊戲的展品出價
這個小節我們使用貝葉斯推斷中的核心思想：使用基於面向最終結果收益的損失函式，進行策略的優化。
0x1：業務問題場景
我們的目標是在一個價格競猜的遊戲中取得最大收益，“價格競猜”的遊戲規則如下：
1. 比賽雙方爭奪競猜站臺上的商品的價格。 2. 每位參賽者看到的商品都不一樣，價格都是獨一無二的，所以並不存在價格各項干擾的情況。 3. 觀看後，每位參賽者被要求給出對於自己那套獎品的投標價格。 4. 如果投標價格超過實際價格，投標者被取消獲獎資格。 5. 如果投標價格低於真正的價格，且差距在250$以內，投標者獲得兩套獎品。
遊戲的難度在於如何平衡價格的不確定性，保持你的出價足夠低，以便不超過實際價格，同時又不能太低，需要儘量接近真實價格。
0x2：對問題進行數學建模
很明顯，這個問題是一個迴歸預測問題，解決這個問題的方法有很多種，我們可以使用機器學習迴歸分析，也可以構建一個深度學習網路進行有監督學習預測，同時，我們也可以使用本章介紹的主題：貝葉斯概率，將我們的先驗知識、訓練樣本、後驗估計結果全部概率化，在概率分佈的框架內對問題進行建模分析。
1. 確定隨機變數
不管是機器學習/深度學習/概率統計方法，在建模的時候都暗含了一個要求，即確認隨機變數。
所以接下來的問題是，在這次程式設計建模中，隨機變數是什麼？
1. 基於歷史價格對對待預測價格的先驗估計：是輸入x； 2. 後驗估計分佈：是預測值（y）； 3. 對其的信心分佈引數：是函式模型的權重引數； 4. 在訓練中觀察到兩個獎項更新後的真實價格：是Y目標值；
2. 確定先驗分佈
先驗分佈的選擇，是一個非常複雜的話題，在PAC可學習理論中，先驗假設是對假設搜尋空間的一種約束，先驗越強，這種約束就越強。先驗分佈的好處是使搜尋結果更快且更靠近先驗約束，缺點就是如果先驗引入的不好，可能會導致最後的搜尋結果欠擬合。
這裡我們假設我們記錄了之前的“價格競猜”比賽，我們將歷史的競猜結果作為真實價格的先驗分佈。
我們假設它遵循一個正態分佈：真實價格 ~ Normal（?p, ?p）
通過之前的記錄，得到的先驗假設?p = 35000，?p = 7500。
筆者思考_1：注意，這裡我們的先驗不是一個確定的實值，而是一個概率分佈，在貝葉斯概率程式設計體系中，讀者朋友一定要深刻去理解這種”一切皆有可能“的思想，不論是先驗知識，還是後驗估計，所有的東西都是一個概率分佈，世界不是確定的，而是概率的。
筆者思考_2：這個問題是一個憑空預測價格的場景，並不是提供商品的某些特徵然後基於這些特徵來預測價格，所以這裡其實暗含了一個前提，所有待預測價格的商品，價格都是接近的，收斂在某個價格區間中，所以不論是先驗分佈還是信念分佈，訓練過程本質上都是在儘量地擬合（靠近）那個價格區間，如果沒有這個前提，這個問題是無解的，或者至少不能用簡單的一個信念函式來解，請讀者朋友仔細體察這點。
3. 決策函式的選擇 - 函式模型F() - 信念函式
我們把先驗分佈理解為機器學習中的輸入x，決策函式相當於F(x)的形式，決策函式負責將先驗分佈轉化為後驗概率估計。同時它也是決定我們應該怎麼玩遊戲的策略判斷依據。
我們現在有一個關於價格的大概想法（即先驗分佈估計），但這種猜測可能與真實價格顯著不同（在舞臺上壓力會成倍增加，你可以看到為什麼有些出價這麼瘋狂）。所以我們需要定義一個信念函式，這個信念函式決定我們有多大可能性會相信自己的先驗判斷。
我們假設我們關於獎品價格的信念也是正態分佈：
Prize i~ Normal（??, ??），i = 1，2。這裡 i = 1，2 表示每個參賽者自己的那套獎品只有兩個獎品（但這可以擴充套件到任何數量）。
該套獎品的價格，可以由 Prize1 + Prize2 + w 決定，其中 w 是一個誤差項。
訓練的結果是得到一個最匹配訓練樣本（歷史價格分佈）的模型引數，即信念函式，這個信念函式是一個連續正態函式，根據不同的價格給出不同的信念，基本上來說，越偏離信念均值中心的價格，信念就越低。
4. 確定後驗估計函式 - 對真實價格估計的模型
首先明確一點，貝葉斯估計是面向概率分佈的一種演算法，演算法得到的後驗估計不是一個確定的實值，而是一個分佈。
我們選擇高斯分佈作為後驗估計的模型，理由是高斯分佈非常適合進行決策。高斯分佈的均值中心天然就包含了最佳後驗概率估計的特性。
同時因為獎品有2個，所有最終的結果是兩個分佈累加的結果，即後驗估計模型函式F(x)是一個複合函式，F(x) = F(x1) + F(x2)。
0x3：基於PyMC預測獎品價格例項
1. 先驗價格分佈
讓我們取一些具體的值，假設套裝有兩個獎品：一趟奇妙的加拿大多倫多之旅，一個可愛的新吹雪機。
我們對這些獎品的真實價格有一些猜測（先驗假設），但是我們也非常不確定。按照上面所討論的，我們可以用正態分佈來表達這種不確定性。
即吹雪機 ~ Normal（3000，500），多倫多之旅 ~ Normal（12000，3000）
這個先驗假設包含了一個事實，我們相信前往多倫多旅行的真實價格為 12000 美元，但有68.2% 的概率價格會下降1個標準差，也就是說，有 68.2% 的概率行程價格在 [9000,15000]區間中。
# -*- coding: utf-8 -*- import scipy.stats as stats import numpy as np import matplotlib.pyplot as plt norm_pdf = stats.norm.pdf plt.subplot(311) x = np.linspace(0, 60000, 200) # 歷史價格先驗分佈 sp1 = plt.fill_between(x, 0, norm_pdf(x, 35000, 7500), color="#348ABD", lw=3, alpha=0.6, label="historical total prices") p1 = plt.Rectangle((0, 0), 1, 1, fc=sp1.get_facecolor()[0]) plt.legend([p1], [sp1.get_label()]) plt.subplot(312) x = np.linspace(0, 10000, 200) # 吹雪機價格先驗分佈 sp2 = plt.fill_between(x, 0, norm_pdf(x, 3000, 500), color="#A60628", lw=3, alpha=0.6, label="snowblower price guess") p2 = plt.Rectangle((0, 0), 1, 1, fc=sp2.get_facecolor()[0]) plt.legend([p2], [sp2.get_label()]) plt.subplot(313) x = np.linspace(0, 25000, 200) # 多倫多之旅價格先驗分佈 sp3 = plt.fill_between(x, 0, norm_pdf(x, 12000, 3000), color="#7A68A6", lw=3, alpha=0.6, label="Trip price guess") plt.autoscale(tight=True) p3 = plt.Rectangle((0, 0), 1, 1, fc=sp3.get_facecolor()[0]) plt.legend([p3], [sp3.get_label()]) plt.show()
2. PYMC訓練及預測後驗概率估計
# -*- coding: utf-8 -*- import scipy.stats as stats import numpy as np import matplotlib.pyplot as plt import pymc as pm norm_pdf = stats.norm.pdf plt.subplot(311) x = np.linspace(0, 60000, 200) # 歷史價格先驗分佈 sp1 = plt.fill_between(x, 0, norm_pdf(x, 35000, 7500), color="#348ABD", lw=3, alpha=0.6, label="historical total prices") p1 = plt.Rectangle((0, 0), 1, 1, fc=sp1.get_facecolor()[0]) plt.legend([p1], [sp1.get_label()]) plt.subplot(312) x = np.linspace(0, 10000, 200) # 吹雪機價格先驗分佈 sp2 = plt.fill_between(x, 0, norm_pdf(x, 3000, 500), color="#A60628", lw=3, alpha=0.6, label="snowblower price guess") p2 = plt.Rectangle((0, 0), 1, 1, fc=sp2.get_facecolor()[0]) plt.legend([p2], [sp2.get_label()]) plt.subplot(313) x = np.linspace(0, 25000, 200) # 多倫多之旅價格先驗分佈 sp3 = plt.fill_between(x, 0, norm_pdf(x, 12000, 3000), color="#7A68A6", lw=3, alpha=0.6, label="Trip price guess") plt.autoscale(tight=True) p3 = plt.Rectangle((0, 0), 1, 1, fc=sp3.get_facecolor()[0]) plt.legend([p3], [sp3.get_label()]) # plt.show() # 吹雪機和多倫多之旅的價格猜測先驗分佈 data_mu = [3e3, 12e3] data_std = [5e2, 3e3] # 對價格的信念分佈就是歷史價格的先驗分佈 mu_prior = 35e3 std_prior = 75e2 true_price = pm.Normal("true_price", mu_prior, 1.0 / std_prior ** 2) prize_1 = pm.Normal("first_prize", data_mu[0], 1.0 / data_std[0] ** 2) prize_2 = pm.Normal("second_prize", data_mu[1], 1.0 / data_std[1] ** 2) price_estimate = prize_1 + prize_2 @pm.potential def error(true_price=true_price, price_estimate=price_estimate): return pm.normal_like(true_price, price_estimate, 1 / (3e3) ** 2) mcmc = pm.MCMC([true_price, prize_1, prize_2, price_estimate, error]) mcmc.sample(50000, 10000) # train price_trace = mcmc.trace("true_price")[:] # predict x = np.linspace(5000, 40000) plt.plot(x, stats.norm.pdf(x, 35000, 7500), c="k", lw=2, label="prior dist. of suite price") _hist = plt.hist(price_trace, bins=35, normed=True, histtype="stepfilled") plt.title("Posterior of the true price estimate") plt.vlines(mu_prior, 0, 1.1 * np.max(_hist[0]), label="prior's mean", linestyles="--") plt.vlines(price_trace.mean(), 0, 1.1 * np.max(_hist[0]), label="posterior's mean", linestyles="-.") plt.legend(loc="upper left") plt.show()
從上圖可以看到，基於吹雪機和多倫多旅行的價格先驗，總體價格從歷史先驗的35000下降到了20000（圖中藍色色塊）。
3. 結合實際業務場景的損失函式
從最大後驗的角度來看，我們此時已經解決問題了，直接給出最終出價即可。但是在貝葉斯統計領域，我們有更多的資訊，我們應該將損失納入我們的最終出價中。
我們需要使用結合業務場景的損失函式來找到最佳出價，即對於我們的損失函式來說是最優解。
def showcase_loss(guess, true_price, risk=80000): if true_price < guess: return risk elif abs(true_price - guess) <= 250: return -2*np.abs(true_price) else: return np.abs(true_price - guess - 250)
其中 riks 是一個引數，表明如果你的猜測高於真正的價格的損失程度。
這裡選擇了一個值: 80000，風險較低意味著你更能忍受出價高於真實價格的情況。
如果我們出價低於真實價格，並且差異小於250，我們將獲得兩套獎品；否則，當我們出價比真實價格低，我們要儘可能接近，因此其損失是猜測和真實價格之間距離的遞增函式。
接下來的問題是，這個 risk 引數如何選擇呢？
對於每個可能的出價，我們計算與該出價相關聯的期望損失，通過改變 risk 引數來看它如何影響我們的最終損失。
# -*- coding: utf-8 -*- import scipy.stats as stats import numpy as np import matplotlib.pyplot as plt import pymc as pm norm_pdf = stats.norm.pdf # plt.show() # 吹雪機和多倫多之旅的價格猜測先驗分佈 data_mu = [3e3, 12e3] data_std = [5e2, 3e3] # 對價格的信念分佈就是歷史價格的先驗分佈 mu_prior = 35e3 std_prior = 75e2 true_price = pm.Normal("true_price", mu_prior, 1.0 / std_prior ** 2) prize_1 = pm.Normal("first_prize", data_mu[0], 1.0 / data_std[0] ** 2) prize_2 = pm.Normal("second_prize", data_mu[1], 1.0 / data_std[1] ** 2) price_estimate = prize_1 + prize_2 @pm.potential def error(true_price=true_price, price_estimate=price_estimate): return pm.normal_like(true_price, price_estimate, 1 / (3e3) ** 2) mcmc = pm.MCMC([true_price, prize_1, prize_2, price_estimate, error]) mcmc.sample(50000, 10000) # train price_trace = mcmc.trace("true_price")[:] def showdown_loss(guess, true_price, risk=80000): loss = np.zeros_like(true_price) ix = true_price < guess loss[~ix] = np.abs(guess - true_price[~ix]) close_mask = [abs(true_price - guess) <= 250] loss[close_mask] = -2 * true_price[close_mask] loss[ix] = risk return loss guesses = np.linspace(5000, 50000, 70) risks = np.linspace(30000, 150000, 6) expected_loss = lambda guess, risk: \ showdown_loss(guess, price_trace, risk).mean() for _p in risks: results = [expected_loss(_g, _p) for _g in guesses] plt.plot(guesses, results, label="%d" % _p) plt.title("Expected loss of different guesses, \nvarious risk-levels of \ overestimating") plt.legend(loc="upper left", title="Risk parameter") plt.xlabel("price bid") plt.ylabel("expected loss") plt.xlim(5000, 30000) plt.show()
最大限度地減少我們的損失是我們的目標，這對應於上面圖中的每條曲線的最小值點。更正式地說，我們希望通過尋找以下公式的解來減少我們的期望損失
期望損失的最小值被稱為貝葉斯行動。我們可以基於SciPy的fmin api進行智慧搜尋，尋找任一單變數或多變數函式的極值（不一定是全域性極值）
# -*- coding: utf-8 -*- import scipy.stats as stats import numpy as np import matplotlib.pyplot as plt import pymc as pm import scipy.optimize as sop norm_pdf = stats.norm.pdf # plt.show() # 吹雪機和多倫多之旅的價格猜測先驗分佈 data_mu = [3e3, 12e3] data_std = [5e2, 3e3] # 對價格的信念分佈就是歷史價格的先驗分佈 mu_prior = 35e3 std_prior = 75e2 true_price = pm.Normal("true_price", mu_prior, 1.0 / std_prior ** 2) prize_1 = pm.Normal("first_prize", data_mu[0], 1.0 / data_std[0] ** 2) prize_2 = pm.Normal("second_prize", data_mu[1], 1.0 / data_std[1] ** 2) price_estimate = prize_1 + prize_2 @pm.potential def error(true_price=true_price, price_estimate=price_estimate): return pm.normal_like(true_price, price_estimate, 1 / (3e3) ** 2) mcmc = pm.MCMC([true_price, prize_1, prize_2, price_estimate, error]) mcmc.sample(50000, 10000) # train price_trace = mcmc.trace("true_price")[:] ax = plt.subplot(111) def showdown_loss(guess, true_price, risk=80000): loss = np.zeros_like(true_price) ix = true_price < guess loss[~ix] = np.abs(guess - true_price[~ix]) close_mask = [abs(true_price - guess) <= 250] loss[close_mask] = -2 * true_price[close_mask] loss[ix] = risk return loss guesses = np.linspace(5000, 50000, 70) risks = np.linspace(30000, 150000, 6) expected_loss = lambda guess, risk: \ showdown_loss(guess, price_trace, risk).mean() for _p in risks: _color = next(ax._get_lines.prop_cycler) _min_results = sop.fmin(expected_loss, 15000, args=(_p,),disp = False) _results = [expected_loss(_g, _p) for _g in guesses] plt.plot(guesses, _results , color = _color['color']) plt.scatter(_min_results, 0, s = 60, \ color= _color['color'], label = "%d"%_p) plt.vlines(_min_results, 0, 120000, color = _color['color'], linestyles="--") print("minimum at risk %d: %.2f" % (_p, _min_results)) plt.title("Expected loss & Bayes actions of different guesses, \n \ various risk-levels of overestimating") plt.legend(loc="upper left", scatterpoints=1, title="Bayes action at risk:") plt.xlabel("price guess") plt.ylabel("expected loss") plt.xlim(7000, 30000) plt.ylim(-1000, 80000) plt.show()
當我們降低風險閾值（不那麼擔憂出價過高），我們可以提高我們的出價，想要更接近真實價格。
這個例子中，因為維度較低，我們可以肉眼找到極值點，但是在高維函式中，肉眼是無法找到的極值的。
這裡列舉幾個損失函式的貝葉斯行動可以用顯式公式表達的例子：
1. 使用均方損失: 貝葉斯行動是後驗分佈的均值； 2. 當後驗分佈的中位數將絕對期望損失函式最小化時，用樣本的中位數來接近是非常準確的； 3. MAP估計是某個損失函式收縮到0-1損失的解；
Relevant Link:
https://www.zhihu.com/question/21134457 http://www.xuyankun.cn/2017/05/13/bayes/ http://mindhacks.cn/2008/09/21/the-magical-bayesian-method/
5. 貝葉斯行動預測的另一個例子 - 金融預測
0x1：問題場景
我們需要設計一個模型，對某隻股票未來的價格進行預測，我們的利潤和損失將直接依賴於基於預測價格的行動。
很顯然，如果預測未來價格是上行的(價格上漲)，而實際也上漲，那不論是精確預測還是預測結果高於實際價格，至少是盈利的；但是如果股票價格實際下行了，預測卻是上行，那麼無論預測結果是靠近還是偏離，都是虧損的，區別只是虧多虧少而已，當然，虧少肯定比虧多要好。
同時還有一個約束條件，作為金融機構來說，穩健的收益回報比大比例的盈利和虧損更被看重，所以即使是預測上行對了，也儘量會避免過大的偏離實際值。
0x2：模型設計
1. 確定損失評價函式
如何衡量模型的預測結果和未來實際價格之間差距的損失？假如我們使用平方誤差損失函式，那麼對正負號是不加以區分的，預測值 -0.01 和 0.03 和實際價格 0.02 之間的懲罰相同:
(0.01 - (-0.01))2 = (0.01 - 0.03)2 = 0.0004
如果我們基於採用了這種損失策略的模型的預測下了堵住，那麼 0.03 的預測值會使你賺錢，而 -0.01 的預測值會使你賠錢，但損失函式無法提示這一點。我們需要一個更好的損失函式，即考慮了預測價格的正負號和真正的價值的損失函式。
# -*- coding: utf-8 -*- import scipy.stats as stats import numpy as np import matplotlib.pyplot as plt import pymc as pm import scipy.optimize as sop def stock_loss(true_return, yhat, alpha=100.): if true_return * yhat < 0: # opposite signs, not good return alpha * yhat ** 2 - np.sign(true_return) * yhat \ + abs(true_return) else: return abs(true_return - yhat) true_value = .05 pred = np.linspace(-.04, .12, 75) plt.plot(pred, [stock_loss(true_value, _p) for _p in pred], label="Loss associated with\n prediction if true value = 0.05", lw=3) plt.vlines(0, 0, .25, linestyles="--") plt.xlabel("prediction") plt.ylabel("loss") plt.xlim(-0.04, .12) plt.ylim(0, 0.25) true_value = -.02 plt.plot(pred, [stock_loss(true_value, _p) for _p in pred], alpha=0.6, label="Loss associated with\n prediction if true value = -0.02", lw=3) plt.legend() plt.title("Stock returns loss if true value = 0.05, -0.02") plt.show()
上圖中我們可以看到幾個關鍵的特徵：
1. 損失函式反應了使用者並不想猜測符號，猜錯符號的損失遠大於猜對符號； 2. 不管是否猜對符號，使用者都不想大幅度猜錯，即大幅度偏離真實值。這也體現了金融機構應對下行風險(預測方向是錯的，量級很大)和上線風險(預測方向是正確的，量級很大)的態度是相似的。兩者都被視為危險的行為而不被孤立。因此，當我們進一步遠離真實價格，我們的損失會增加，但在正確方向上的極端損失相對錯誤方向上的較小。
2. 確定函式模型F()
我們將會在被認為能夠預測未來回報的交易訊號上做一個迴歸。資料是人工虛構的，因為絕大多數時候金融資料都不是線性的。
下圖中，我們沿著最小方差線畫出了資料分佈情況：
# -*- coding: utf-8 -*- import scipy.stats as stats import numpy as np import matplotlib.pyplot as plt import pymc as pm import scipy.optimize as sop # Code to create artificial data N = 100 X = 0.025 * np.random.randn(N) Y = 0.5 * X + 0.01 * np.random.randn(N) ls_coef_ = np.cov(X, Y)[0, 1] / np.var(X) ls_intercept = Y.mean() - ls_coef_ * X.mean() plt.scatter(X, Y, c="k") plt.xlabel("trading signal") plt.ylabel("returns") plt.title("Empirical returns vs trading signal") plt.plot(X, ls_coef_ * X + ls_intercept, label="Least-squares line") plt.xlim(X.min(), X.max()) plt.ylim(Y.min(), Y.max()) plt.legend(loc="upper left") plt.show()
線性迴歸的函式形式如下：
對一個特定的交易訊號 x，回報的分佈有如下形式：
對於給定的損失函式，我們希望找到：
下圖中我們對比了在平方損失函式和我們自定義損失函式條件下，不同的交易訊號的貝葉斯行動。
從上圖中可以看出：
1. 當交易訊號接近為0時，正負回報都可能出現，我們最好的行為是預測接近為0，也就是說處於中立； 2. 只有當我們都非常自信時，我們才進場下注，當對不確定性感到不安時，選擇不作為； 3. 當訊號變得越來越極端時，我們對正負回報的預測越來越自信，自定義損失函式的模型將會收斂到最小二乘；
具有以上特點的模型被稱為係數預測模型。稀疏模型不是想要千方百計地去擬合資料，稀疏模型試圖找打針對我們定義的股票損失的最好預測，而最小二乘只試圖找到相對於平方損失下的資料的最佳擬合。
Relevant Link: