如何选择我们的损失函数 - 连接统计推断和问题所在领域的桥梁
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:
