新冠肺炎的超声波——深度学习方法
关于一个利用超声波自动检测新冠肺炎的开放倡议
您将从本文中学到什么:
- CT 和 X 射线是新冠肺炎诊断过程中常用的评估方法，但医学界也提倡超声成像
- 超声波是优选的，因为它便宜、便携、易于消毒、无辐射且随处可得
- 我们的深度学习方法首次证明了在超声上自动检测新冠肺炎可能是可行的
- 有了可解释的机器学习方法，我们向医生的决策支持工具迈进了一步
- 在 GitHub 上开始了一项开放存取数据的倡议——贡献吧！
https://github.com/jannisborn/covid19_ultrasound
迄今为止，已有超过 1 . 4 亿人感染了新冠肺炎病毒，300 多万人已经死亡。这种疾病仍在全球蔓延，几乎所有国家都受到影响。努力控制病毒的一个关键因素是快速可靠的检测。检测病毒 RNA 的新冠肺炎 RT-聚合酶链式反应测试(逆转录聚合酶链反应)被认为是最可靠的，但它的灵敏度在各国之间差异很大，有时假阴性高达 30%[1]，在发展中国家几乎没有。此外，PCR 测试需要几个小时，这是有问题的，例如在分诊的情况下，医生必须立即决定是否隔离病人。
这就是医学成像的用武之地。迄今为止最突出的是 CT，它已经在疫情早期用于新冠肺炎的快速诊断。研究表明，在 CT 上可以看到明显的不同模式，例如，空气间隙实变和所谓的(多灶性)毛玻璃阴影[2]。从报告来看，CT 似乎是一个非常有前途的工具，有时甚至在聚合酶链式反应测试失败时也很敏感[3]。但是它有很大的缺陷:CT 辐射大，价格昂贵，难以消毒。显然，这些缺点阻碍了 CT 在诊断中的广泛应用。作为一种替代方法，研究人员考虑了 X 射线，但预测能力较差。但是还有另一种医学成像工具，尽管它很受欢迎，但在新冠肺炎的环境中却没有得到太多的关注。我们在谈论超声波。
护理点超声波(POCUS)是
- 便宜:虽然一次 x 光检查估计要花费 370 美元左右，CT 从 500 美元到 3000 美元不等，但超声波(US)很便宜，只需大约 140 美元。此外，该设备本身很便宜，因此易于销售，便携式设备的起价为 2000 美元。
- 易于使用:几乎所有的医生都知道如何进行超声波检查。没有辐射那样的安全措施，而且设备很方便。
- 快速:使用一台设备，每小时可以进行 4 到 5 次肺部筛查
- 便携式:“护理点”说明了一切。病人不必移动，这节省了大量的时间和精力。
- 安全:我们不使用任何辐照元素。句号。任何 X 射线或 CT 检查都会略微增加患癌症的终生风险，尤其是对年轻患者而言。
尽管有这些优点，但超声仅在过去几年中被整合到肺部疾病的诊断过程中。可以分析某些特征病理模式，例如所谓的 B 线、A 线和条形码标记，以诊断所谓的胸腔积液、肺泡实变、间质综合征和气胸。由于新冠肺炎以类似的方式改变肺部超声模式，超声对新冠肺炎的适用性现在在几个出版物中进行了研究，其敏感性与 CT 进行了比较[4]。结果很清楚:根据[5，6，7]，[8](发表在权威杂志*《柳叶刀呼吸医学】* 上)的作者主张 us 在新冠肺炎诊断中发挥更突出的作用，并提供证据表明 US 检测新冠肺炎的灵敏度与 ct 非常相似。
那么这是否意味着超声将很快取代 CT 作为诊断工具呢？不幸的是，事情没那么简单。一个值得注意的缺点是，在美国，医生必须经过培训才能识别 COVID，并且观察 COVID 特定的模式不是一项容易的任务，而是需要一些经验。在目前的情况下，进行广泛培训的时间有限。
我们相信，一个用于医生的辅助系统给出了一个解决方案，该系统利用计算机视觉技术对超声记录进行自动分类。这种系统可以支持医生的决策，并对病人被感染的可能性进行初步评估。在我们的研究中，我们朝着这个系统迈出了第一步，从各种在线来源收集了一组美国记录，对它们进行预处理，并训练一个神经网络对图像进行分类。我们在 arXiv“poco vid-Net:从新的肺部超声成像数据集(POCUS)中自动检测新冠肺炎”的预印本实际上是 COVID 在超声方面的第一项工作。在我们发表在应用科学杂志“使用可解释的超声图像分析加速肺部病理检测”的文章中，我们随后提供了一项对用于对肺部 US 数据进行分类的(可解释的)深度学习方法的综合研究。
让我们开始更详细地讨论我们的工作:
一个新的 POCUS 数据集
一个主要的贡献是收集了目前********(92 个新冠肺炎，73 个细菌性肺炎和 90 个健康对照)的 250 多个数据集。该数据集包含来自各种来源的数据，包括我们的合作者在诺森比亚(英国)的一家医院收集的未发表的临床数据，以及在 Neuruppin(德国)扫描的健康志愿者的记录，以及来自其他出版物和教育网站的数据。请注意，预处理数据需要做更多的工作——简而言之，在从视频中选择帧来创建图像数据集之前，我们手动裁剪了视频并移除了伪影。我们甚至在团队中一名医生的帮助下检查了数据的质量，他可以对每个视频中的可见模式提供宝贵的意见。最后，这产生了一个超过 3000 张图像的干净数据集。当然，这绝不是详尽无遗的，还有许多数据需要发现和处理。然而，我们认为在我们的 GitHub 页面上提供数据和预处理管道是我们工作的主要部分，在这里您可以为我们的开放访问计划贡献数据。****
一名新冠肺炎患者的肺部超声波检查(取自 http://www.thepocusatlas.com/covid19(第 4 天)，根据许可归属-非商业性 4.0 国际版提供)
深度学习方法
但是这些数据有用吗，神经网络可以学习从 US 图像中检测新冠肺炎吗？为了回答这个问题，我们训练了各种神经网络来区分新冠肺炎、肺炎和健康患者。我们的最佳性能模型利用预训练的 VGG-16 模型作为主干，后面是一个完全连接的层。诸如旋转、移动和翻转之类的增强有助于防止数据集上的过度拟合，对于具有两百万可训练参数的模型来说，数据集仍然相当小。
结果非常有希望——我们的模型实现了 88%的准确率来分类为新冠肺炎、肺炎和健康患者。下面你可以看到混淆矩阵，分别沿每个轴标准化。这些结果是在分层的 5 折叠交叉验证中获得的，确保来自同一视频的帧不在不同的折叠中。如下图所示，不同类别之间的性能是平衡的，这意味着对新冠肺炎和细菌性肺炎的敏感性和特异性都很高。这也表明，深度学习方法不仅对新冠肺炎来说，而且对一般的肺部病理学来说，都可能是一项值得努力的工作。
结果:新冠肺炎的特异性为 0.9，敏感性为 0.88
可解释人工智能——一个决策支持系统
当我们提出这些结果时，我们经常被问及我们是否试图用这些方法取代医生。明确的答案是否定的！相反，我们希望提供第二种意见，一种可以在诊断过程中帮助医生的系统。例如，如果一个软件可以在你可以在超声波上看到有趣的东西的特定时间点停止视频，然后在图像上突出显示异常，这不是很好吗？
********
这样做的一种方法称为类激活映射(CAM)[9]，它突出显示图像中对神经网络预测最有决定性的部分。我们为我们的分类网络实现了 CAM，并观察到有趣的结果:通常，热图集中在超声上可见的所谓实变，这表明肺部有液体。此外，如示例所示，通常健康肺部的超声水平“A 线”也经常被摄像头拾取。
现在，为了找出这些激活图是否真的对临床实践有用，我们请两位医学博士对覆盖在一段视频上的 50 个激活图例子进行评价。这些凸轮总的来说是有用的。因此，我们认为这可能是一个非常有趣的研究方向，以改善特定病理模式的检测，并利用医生可以使用的实际软件的实施。
下一步是什么？
那么如何将这项初步但有前景的工作转化为实际应用呢？还有很长的路要走，但是为此，我们需要你的帮助。作为第一步，我们需要更多的数据来验证我们的模型——如果你有任何对我们有用的数据:通过我们的 GitHub 页面贡献给我们的数据集。我们也渴望支持任何旨在用我们的数据集改进机器学习模型的工作。
有了这些初步结果，重要的一步是在真实的临床数据上评估该方法。为此，我们目前正与一家医院合作进行一项临床研究，以验证我们的方法。进一步的研究可以在临床研究中系统地比较 CT、X 射线和超声波以及它们各自的自动检测方法。这将是一个巨大的飞跃。
但是，如果你不知道任何数据来源，也不是机器学习专家，仍然有许多可能性可以做出贡献。一件非常简单的事情就是提高对使用 POCUS 作为诊断工具的可能性的认识。分享我们的出版物，观看我们的宣传视频或查看我们的开发帖子，并访问我们的 GitHub 登录页面，其中提供了所有代码。如果您想投稿或有问题，请联系我们，欢迎任何反馈。感谢您阅读我们的工作！
参考
[1]杨，杨，等，“2019-nCoV 感染的病毒脱落的实验室诊断和监测。”(2020) MedRxiv 。
[2]杰弗里·P·凯恩。中国武汉 2019 年新型冠状病毒(2019-ncov)感染的胸部 ct 发现(2020 年)放射科医师的要点。
[3]艾，陶，等.“冠状病毒病 2019(新冠肺炎)胸部 CT 与 RT-PCR 检测相关性:1014 例报告”(2020) 放射学 : 200642。
[4]菲亚拉，M. J .“超声在新冠肺炎:与 CT 相关的超声发现时间表。”(2020) 临床放射学。
[5] Smith，M. J .等人，“COVID‐19 患者的床旁肺部超声检查-叙述性综述。”(2020) 麻醉。
[6]Sofia，Soccorsa 等，“胸部超声和 SARS-新冠肺炎:一篇图片论文。”(2020) 超声杂志:1–5。
[7] Soldati，Gino 等人，“肺部超声在 COVID‐19 疫情期间有作用吗？."(2020) 医学超声杂志。
[8]布昂森索、达尼洛、达维德·帕塔和安东尼奥·基亚雷蒂。“新冠肺炎爆发:少用听诊器，多用超声波.”(2020) 柳叶刀呼吸内科。
[9]周，，等.“学习深层特征用于鉴别性定位”IEEE 计算机视觉和模式识别会议论文集。2016.
UMAP 和 K-意味着对角色进行分类——以及它为什么有用(英雄联盟)
利用降维和聚类算法将英雄联盟“冠军”进行分类。
截至今天(2020 年 10 月)的，在线游戏《英雄联盟》中有 151 个“冠军”(可玩角色)，每个人都提供独特的个人游戏风格，这使其成为世界上最受欢迎的电子竞技。尽管这种多样性提供了一个引人入胜的竞争环境，但同时也给新玩家和分析师带来了复杂性。如果每个冠军都是独一无二的，一个新手怎么可能理解每场比赛的复杂性，或者一个分析师怎么可能总结一个球员的表现？
2016 年 4 月，Riot(游戏背后的公司)试图通过引入“冠军课程”来帮助解决这个问题。这由 12 个子类组成，分属 6 个类。这些是 Riot 的开发团队手工制作的，为新玩家习惯游戏提供了一个很好的起点。你开始玩几个游戏，看到一个叫“Thresh”的角色，快速检查会发现他的职业是“Catcher”，另一个角色也是:“Blitzcrank”，你以前和他交过手。虽然你错过了一些细节，但你明白总体思路是避免被它们抓住。
然而，这些课程很快就被搁置了，现在很少在任何内容中提到它们。这促使我回答了一个问题:AI 会如何班每个冠军？这与暴乱的解释一致吗？如果是这样，我们还能了解到什么？为此，我们将该方法分为四个阶段。
- 数据收集、清理和特征创建。
- 通过提取关键信号来降低数据的维数。
- 使用聚类算法将冠军分成几类。
- 分析课程以确定趋势和主题。
关于游戏背景的预先说明(如果您熟悉游戏，请忽略):
如果你不熟悉《英雄联盟》,你会漏掉一些内容。如果你想变得更熟悉，那么你可以阅读这个游戏的 介绍 。否则，最起码要知道这是一场 5v5 游戏，每个队都有一名球员扮演以下角色之一:顶级球员、中层球员、Jungler、ADC、支援。
如果你希望在文章中加入商业元素，那么将“冠军”转化为“顾客”，将“死亡”等数据转化为“购买”和“网站反弹”可能会有所帮助。你可以跳到文章的底部，获得一些在“真实世界”中使用的例子。
数据收集、清理和特征创建
通常，第一种方法是启动一个 Riot API 会话，并开始查询他们庞大的数据湖。 注意: *我打算创建一个完整的指南来说明如何做到这一点，但同时请查看 GitHub 以了解详细信息。*我们首先决定将结果限制在钻石(接近顶级)玩家，因为它提供了大量被认为质量相对较高的游戏。总的来说，我们拉了 68，000 多场比赛，其中每一行都包含关于冠军的信息和一系列关键统计数据，从杀戮，死亡和助攻到最大的致命一击，抑制剂和龙的安全。
除了“标准”数据点，我们还可以考虑从可用数据中提取的替代特征。例如，可以查看一个名为“时间线事件”的数据集，并查看每次杀戮的详细信息。这包括地点、时间、凶手和被杀者。从这里，我们可以创建两个新功能:单杀:在一场游戏中有多少次玩家在没有帮助的情况下获得了一次杀戮& 早期漫游:在一场游戏的前 10 分钟里有多少次玩家在不是自己的车道上获得了一次杀戮或帮助。
在没有过多阐述单调的情况下，数据被清理、转换，然后聚合到每一行都包含每个角色中冠军的平均统计数据(如果有的话)。这意味着 Ekko 在中间车道和丛林中有单独的一行，两者都包含每个统计数据的相关平均值。最终特征如下:
分析的每场比赛的所有统计列表。
通过提取关键信号降低数据的维度
我们现在的问题是，有太多的数据，很容易对每个冠军进行分类。重要的是，数据中存在大量的多重共线性*，其中许多列高度相关。例如，一个冠军平均杀死的数量和玩家连续杀死的次数，或者他们减少的伤害和他们受到的伤害之间有明显的关系。*
幸运的是，降维提供了一个可靠的解决方案。无需在这一部分花费太多精力，这种技术采用所有已知的变量，并试图将数据“缩小”成更少的特征。为了帮助描绘一幅心理图画，想象你只有两个变量:杀死和获得的黄金。如果你把这些标绘成二维空间，你会得到类似如下的东西:
用回归线绘制的获得的黄金和杀死的总数的二维相关图。
现在，注意贯穿数据的最佳拟合线。如果你把它分离出来，并把它翻转过来，你会得到一条单一的一维线！实质上，1D 线是两个变量综合结果。左边是“穷不杀人”，右边是“富不杀人”。尽管您丢失了粒度数据，但您保留了它们在单一维度中的本质。虽然这是对幕后发生的事情的一种过于简单化的描述，但这可能有助于你的直觉。
至于降维工具包，有三个标准选项:
为了这个项目，我决定和 UMAP 一起去。我有两个原因:当试验各种技术和参数的组合时，我发现 UMAP 提供的类与我的直觉一致(是的，这使它有偏见！);在我以前的项目中，使用有正确或错误答案的标记数据(想想:癌症检测或数字识别)，UMAP 提供了比 t-SNE 和 PCA 更好的准确性结果。
应用这项技术非常简单，只用三行代码就可以完成。唯一值得一试的属性是数据将被简化成多少个部分。
UMAP 包的 Python 实现。
但是，请注意，在应用 UMAP 变换之前，数据应该标准化，为此，请在数据帧上使用 SKLearn 函数 StandardScaler。GitHub 可以看到更多关于这方面的信息。
*作为总结，本节将 39 个原始特性转化为 5 个核心组件。每个组件都“包含”来自原始数据集的信息，但是已经根据它们的底层关系进行了分组。如果这有助于你的理解，想象一下一个组件被称为“杀死”，虽然包含了所有特征的一部分，但它与和杀死相关的变量关系最大(想想:黄金，杀死条纹，高伤害，等等…).**你可能会问，5 从何而来？*实验和测试！它似乎提供了最可行的结果(再次，偏见警告)。
使用聚类算法将冠军分成几类
我们现在有了每个角色中每个冠军的 5 维数据集，是时候尝试将他们分为四类了。4 也是相对任意的，实际上可以使用任意数量的类。然而，从实验中我发现，当分成四个类时，简单性和多样性得到了很好的平衡。尽管这再次在过程中产生了偏差，但鉴于我们正在处理无监督学习，没有“正确”的答案或基于结果的方法来优化。在实践中，这些决定最终留给了那些具有领域知识的人(或者更常见的是，基于需求)。
为了创建这四个类，我们可以使用一种叫做 k 均值聚类 的技术。该技术着眼于将数据放入围绕中心点的聚类中，其中每个数据点将根据它们与这些中心中的哪一个最近而被分类。为了做到这一点，该算法首先随机选择四个质心，并根据最接近的质心对每个冠军进行分类。一旦它们都被分类，它就找到这些新分配的簇的中心，然后使用它们作为更新的质心。这样一直持续下去，直到这些中心点稳定下来。
这里有一个可视化的算法来帮助你:
k-means 算法的可视化示例，从随机的三个点开始，然后慢慢向每个类的最终质心移动。
分析课程以确定趋势和主题
一旦 k-means 在减少的数据上运行，我们现在将所有冠军分成四个单独的集群。由于我们是在一个 5 维空间中操作的，所以不可能可视化输出来帮助我们理解这些类，所以我们可以基于这四个类聚集所有的原始统计数据。因此，我们将有一个表，其中包含四个类以及我们从 Riot 数据湖中获取的每个要素的平均值。以下是摘录:
显示 ADC 角色中四个类的平均聚合统计信息的表。
现在是研究这些的时候了，以便更好地感受它们之间的不同。为此，我运行了一个脚本，将每个统计数据的类平均值与车道平均值进行比较，然后根据百分比差异对它们进行排序。
例如，在《丛林》中，我发现一个职业使得在他们的游戏中有 18.23%的第一次杀人，而的角色平均只有 12.77%。同样的职业也比平均水平少清理了 6%的丛林营地。被归类到这个集群的冠军？沙可，李辛，伊莉斯，赵信。这帮助我描绘了一个团队，我后来将它命名为“3 营地 Gankers”，参考这些冠军采取的共同丛林路径，他们在那里做红 Buff，蓝 Buff，Gromp，然后立即攻击一条车道或入侵敌人的丛林营地(感谢 jung.gg 提供的数据)。
您可以在这里找到每条赛道所有级别的列表、简短描述、哪些 meta Champions 适合它们，以及超过平均水平的 3 个关键统计数据和表现不佳的 3 个关键统计数据:
关键点:我在这里强调的一点很关键，提供这些优点和缺点是为了给每个职业一个快速的直觉。对于个人冠军而不是来说，拥有和我选择的人一样的实力是很常见的。他们仍然属于这个班级的原因是，平均而言，在所有的统计数据中，他们比任何其他人都更接近于那个，但他们不一定是这个班级的完美孩子！**
例如，中间车道的冠军 Kassadin 被归类为漫游者他们的优势在于他们经常在游戏早期攻击敌人的车道。卡萨丁本人并没有分享这种特殊的品质，但是他和这个职业分享了比其他职业更多的其他属性。**
那又怎样？
正如我经常说的，为人工智能而人工智能是浪费时间。将冠军分段成班有什么好处*？这是一个非常有效的问题，也是在研究任何技术、方法或项目时都应该提出的问题。这里有几个例子(虽然有些夸张和虚构)可能有助于回答这个问题(*警告:除非你熟悉这个游戏，否则你可能很难理解，如果不熟悉——第三个例子适用于商业问题):
你是一家领先的欧洲职业电子竞技团队的分析师。在一次前所未有的事件转折中，对方球队的首发中场球员在走上舞台时戏剧性地摔倒受伤，并被一名第一次上场的新秀取代。离冠军评选还有 10 秒钟，你的主管蔻驰跑过来对你说“我需要了解他什么？!"。根本就没有时间开始列出他的共同冠军和独奏队列性能统计数据。你需要一句俏皮话来给教练提供足够的信息，让他有优势。幸运的是，你已经有了早期人工智能生成的冠军职业，所以你对玩家最近 50 场单人队列游戏进行了快速分析，并说:“他最常见和表现最好的冠军职业是安全定标员！”。你已经向教练简单介绍过了，他知道你的意思。他的目标是 bans Azir 和 Corki，然后首先选择 Orianna。这个家伙窒息了，因为他要么被迫玩一个较低层次的安全定标器，要么被迫离开他的舒适区，在那里他不断被抓到并喂食。
这是一个过于感性化的例子，但却传达了一个重要的信息:现在的总结有时比以后的细节更好。
一个更快且(稍微)不太引人注目的例子是，如果你在第 10 季中以毛凯的身份统治第一车道，却在 10.22 补丁中读到他将被遗忘，他的所有基础数据都减少到-40。你喜欢他的踢球风格，想找个人来代替他。你会发现毛凯被认为是一个推动者和参与者，就像马尔菲特和奥恩一样。Ornn 目前有 98%的胜率，所以你决定他将是你最后几周攀登的合适人选。
寄语:知道如何宽泛地描述你喜欢的东西有助于你找到你可能喜欢的其他东西。事实上，我写了一篇文章使用类似的技术来推荐冠军。
我还将提供一个业务示例，因为我认为欣赏这些技术的“真实世界”使用是很重要的:
假设您是一家小配件公司的营销主管，最近您获得了一笔预算，要创建五个营销活动发送给当前的客户。如果你执行上述相同的步骤，但使用人口统计、购买和网站互动等数据，你可以将你的客户分为五个“类别”。通过研究这些，您将能够针对特定类型的客户量身定制您的活动，从而提高其有效性。例如，您可能会发现一组客户是“重度研究人员”,他们不定期购物，在做出决定之前会在网站上花很长时间，通常会将物品放在篮子里几个星期。另一方面，有些“自发”顾客点击了 Instagram 链接，并在 30 秒内结账。第一类中的所有客户都会收到一封电子邮件，其中包含关于您最新产品的详细报告、优点、缺点和更详细的信息。二班收到一段 Instagram 视频，视频中有一只可爱的狗狗在玩你的产品。
信息？个体层次的分析是强大的，但是对于用例来说通常是不可行的。在这种情况下，为了完成工作，可以“扔掉”一些更精细的细节。
你已经看到文章的结尾了！我叫 Jack J，是一名将人工智能应用于竞技游戏和电子竞技的专业数据科学家。我是iTero的创始人。GGT5【和jung . GG。你可以在 Twitter 上关注我，加入 iTero Discord 或者给我发邮件到 jack@itero.gg 。下一场见。
结合 Knime 和 Python 的嵌入和聚类
实践教程
KNIME 中聚类 MNIST 数据库的 UMAP 降维和 DBSCAN
聚类。橄榄和树叶。形状和颜色。(图片由作者提供)
Knime 是一个免费的开源数据分析、报告和整合平台。KNIME 通过其模块化的数据流水线概念，集成了用于机器学习和数据挖掘的各种组件。对于像我这样没有很强编码背景的人来说，Knime 是端到端数据科学体验的大门。在过去两年与 Knime 的合作中，我已经能够理解数据科学的全球图景，我将我的恐惧抛在身后，现在我在充满挑战和令人惊叹的数据科学世界中越走越深。
几天前，我在《走向数据科学》上读到了一篇文章，作者是罗莎丽亚·西里波和米莎·利索夫伊，展示了一个非常好的 Knime 工作流，其中涉及应用于 MNIST 数据集的 t-SNE 。我曾经和 UMAP 一起做过类似的工作流程，效果不错。
受到 Linkedin 上的 Rosaria Silipo 的鼓励，我决定发表我的第一个故事。
这个故事的提议是什么？
本文的原始资料是闫乐存的从 0 到 9 的手写数字的 OpenML 的 MNIST 数据集。该数据集包含 70000 个 28x28 像素的手写数字，灰度等级为 0-255。也就是说，我们有 70000 个向量，维数是 784。
该数据在 Knime 内的 Python 脚本中读取，并在该 Python 脚本上应用 UMAP 降维算法，因此我们可以将数据集的每个图像映射到低维空间的一个点上。出于表现结果的目的，逻辑上，我使用 2D 和 3D 作为目标尺寸。
之后，为了更有趣，由于生成的 UMAP 嵌入数据忽略了标签，我在 Knime 中应用了一个 DBSCAN 节点来对结果进行聚类。我们稍后将检查结果。
UMAP:降维的一致流形逼近和投影
均匀流形近似和投影(UMAP) 是一种降维技术，可用于类似于 t-SNE 的可视化，也可用于一般的非线性降维。该算法基于对数据的三个假设
-数据在黎曼流形上均匀分布；
-黎曼度量是局部常数(或者可以近似为常数)；
-歧管是局部连接的。
什么是黎曼流形？等等！！什么是流形？流形是 n 维的拓扑对象，在一个点的邻域上，它的局部行为类似于 n 维欧几里得空间。一个好的外行例子是地球。Rieamannian 流形是光滑的可微流形，其中可以定义函数的长度、曲线、角度、体积、梯度等概念。更多关于流形或黎曼流形的信息，试试这里的和这里的。
维基百科中的一个黎曼流形。(斯特拉特，公共领域，通过[维基共享](http://StuRat, Public domain, via Wikimedia Commons))
UMAP 有一些参数来控制其性能。
n _ neighbors*:n _ neighbors 越大，算法寻找的邻域越大，并且趋向于全局整体结构。n_neighbors 越低，UMAP 将关注的细节结构越多。*
最小距离 :定义算法将考虑在最终缩减空间中将点打包在一起的最小距离。这意味着较小的最小距离值将导致较小的聚类嵌入。另一方面，最小距离的高值将导致 UMAP 算法指向全局结构。
度量 :我们可以定义算法在原始空间中计算距离的度量。
难以消化！我知道。对工程师来说近乎神奇。从黑暗中走出来，明确地(简单地！！)，UMAP 通过 n 维点(n 低至 2)对每个高维对象进行建模，使得相似的对象通过附近的点进行建模，而不相似的对象通过远处的点进行建模的概率很高。
基于密度的噪声应用空间聚类
DBSCAN， 基于密度的具有噪声的应用的空间聚类 ，是在无监督的 lerning 方法中考虑的一种聚类算法。基于三个输入参数，该算法能够通过将紧密地打包在一起的点(具有一定数量的附近邻居的点)分组在一起，将单独且远离地位于低密度区域(其最近的邻居太远)中的那些点标记为噪声，来将一组数据群集化。
可以想象，支配 DBSCAN 算法的性能的参数与决定一个点是否彼此靠近的距离 epsilon 有关，并且还与认为一个点是“邻域”或核心点所必需的附近点的最小数目最小点。集群是所有核心点及其邻居的总和。最后，在 Knime 扩展 KNIME 距离矩阵中包含的 DBSCAN 的 KNIME 节点中，我们可以看到需要一个度量输入，因此算法具有必要的距离模型来计算数据集上点之间的 de 距离。
正在配置 DBSCAN 节点。(作者图片)
总结一下:
epsilon :指定点之间的距离，因此它们可以被认为是集群的一部分。如果两个点之间的距离小于 epsilon，则这两个点将被视为属于同一集群的邻居。
最小点数 :从点 p 可达的最小点数，使得点 p 可被视为核心点。
距离模型 如前所述，需要用距离模型通知 DBSCAN 节点对数据集中的点间距离进行评估。为此，我们使用 Mumeric Distance 节点。
DBSCAN 算法的显著优点是:
- DBSCAN 不需要用户指定群集的数量。
- DBSCAN 管理噪声概念，并丢弃和标记噪声点
Knime 工作流
这两个 Knime 工作流均在 Knime Hub 中发布和提供:
UMAP DBSCAN MNIST 2D 巨蟒。UMAP 2D 维嵌入
UMAP DBSCAN MNIST 3D 巨蟒。UMAP 三维嵌入
我将只关注 2D 解释，因为两个工作流实际上是相同的。不过，我会在最后展示两者的比较结果。
遵循你可以看到整个工作流程，稍后我会解释它逐块。
2D UMAP 嵌入的工作流程。(图片由作者提供)
Python 脚本
接下来，您可以找到第一个数据块的详细图片，使用该数据块从 Knime Python 扩展的 Python 源节点中提取数据并计算 UMAP。
Python 抓取和 UMAP。(图片由作者提供)
中间的橙色块是 Python 源节点，包括用于从 OpenML 获取数据和执行 UMAP 技术的脚本。两个变量节点“Metric”和“Number of rows”被合并并像 UMAP 算法的输入一样被传递。对于最终的结果，我使用了欧几里德度量，关于“行数”整数输入节点，我只是在测试时用它来改变低数值。只是时间问题。在右边，我使用列重命名节点在工作流中用可理解的名称命名列。
在这里，您可以找到 Python 源节点中的代码:
from pandas import DataFrame
from sklearn.datasets import fetch_openml
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns# Dimension reduction and clustering libraries
import umapin_batch = flow_variables[‘input_batch’]
umap_metric = flow_variables[‘metric’]#Varias metricas posibles
# mnist_784
# Fashion-MNIST
mnist = fetch_openml(‘mnist_784’, version=1, as_frame=False)
mnist.target = mnist.target.astype(int)standard_embedding = umap.UMAP(
n_neighbors=30,
min_dist=0.0,
n_components=2,
metric=umap_metric,
random_state=42
).fit_transform(mnist.data[:in_batch, 1:])#Añadimos los target
array2D_1=standard_embedding
array2D_2=mnist.target[:in_batch]
array2D_2=np.reshape(array2D_2,(in_batch,1))
array_to_save=np.concatenate((array2D_1,array2D_2),axis=1)# Output table for the node
output_table = DataFrame(array_to_save)
第一次绘图
在执行了 UMAP 算法之后，我们有了 MNIST 数字的 2D 映射。为了比较结果，我使用了两种不同的绘图方法。
策划 UMAP。(图片由作者提供)
UMAP 输出散点图(Plotly)结合色彩管理器
来自的散点图 Plotly Knime 扩展是一个非常好的节点。结果是交互式的，输出颜色可通过颜色管理器节点进行配置，因此缩减的 2D 空间中的数据集点可根据其在嵌入前在原始空间中的位数进行标注。
用散点图绘制的 UMAP 算法定义的 2D 簇。(图片由作者提供)
或者，在 3D 中，我们可以使用散点图节点，而不是散点图节点，结果如下。
用散点图 3D 绘制的 UMAP 算法定义的 3D 聚类。(图片由作者提供)
请注意，UMAP 算法已经在 2D 和 3D 中检测和绘制了相同的 10 个清楚区分的簇。
用 Python 视图输出 UMAP
使用 Python 源节点和一段代码，也可以找到类似的结果。
用 Python 绘制的 UMAP 算法定义的 2D 聚类。(图片由作者提供)
用 Python 绘制的 UMAP 算法定义的三维聚类。(图片由作者提供)
接下来，您可以在 Knime 的 Python 视图节点中找到执行这两个绘图的代码片段。
为了 2D
from io import BytesIO
import matplotlib as mplt
import matplotlib.pyplot as plt
import numpy as nplabel_array=input_table['Col2']
colors = ['red', 'blue', 'green', 'yellow', 'orange', 'peru', 'lime', 'pink', 'steelblue', 'olive']plt.scatter(input_table['Col0'], input_table['Col1'],
c=label_array, alpha=0.5, marker='o',
cmap=mplt.colors.ListedColormap(colors),
s=0.1,
label="Label")
plt.xlabel("X")
plt.ylabel("Y")
plt.legend(loc='best')buffer=BytesIO()
plt.savefig(buffer,format="svg")
output_image = buffer.getvalue()
对于 3D:
from io import BytesIO
import matplotlib as mplt
import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d import Axes3Dlabel=input_table['Col3']colors = ['red', 'blue', 'green', 'yellow', 'orange', 'peru', 'lime', 'pink', 'steelblue', 'olive']fig = plt.figure(figsize=(12,12))
ax = fig.add_subplot(111, projection='3d')ax.scatter(input_table['Col0'],input_table['Col1'] ,input_table['Col2'],
c=input_table['Col3'],
cmap=mplt.colors.ListedColormap(colors),
marker='o',
s=1)ax.set_xlabel('X Label')
ax.set_ylabel('Y Label')
ax.set_zlabel('Z Label')buffer=BytesIO()
plt.savefig(buffer,format="svg")
output_image = buffer.getvalue()
DBSCAN 节点和 sorroundings
因此，到目前为止，我们有一组 2D(或 3D)点，代表一组 700000 个书面图像。我们知道哪个原始数字对应于缩减空间的每个点，因为原始标签已经添加到 UMAP 输出中，正如您在此处看到的:
第 0 列:X，第 1 列:Y，第 2 列:对应于原始数字的标签。(图片由作者提供)
但是想象一下我们没有原始标签。我们会有大量没有进一步信息的点，如下所示:
未标记的点，在灰色物质下面有什么结构吗？。(图片由作者提供)
现在是时候使用 DBSCAN 并尝试找到在这个分散的无意义的图下面是什么结构了。
DBSCAN 节点。(图片由作者提供)
我们将 UMAP 的输出驱动到 DBSCAN 节点。此外，我们为其余的 DBSCAN 输入提供变量，尽管它们可以在节点内部配置。
DBSCAN 的配置屏幕。参见注释“Eps 和 Minpts”由变量控制。(图片由作者提供)
对于 2D DBSCAN，我使用了ε= 0.15和最小点数=50 。数值距离节点中的公制被配置为欧几里得。
经过一些字符串操作后，输出是一个包含聚集数据的表。
请参见包含集群数据的最后一列。请参阅高亮显示的行，其中一个点被认为是噪点。(图片由作者提供)
如果我们在 2D 或 3D 的缩减空间内绘制点，并用 DBSCAN 产生的聚类进行标记，我们会发现什么？
让我们看看。
总体结果
在下图中，您可以看到 DBSCAN 识别的群集与 UMAP 根据原始映像交付的群集相同。这是非常令人满意和美好的。
2D 数据库扫描聚类。(图片由作者提供)
3D DBSCAN 聚类。(图片由作者提供)
DBSCAN 的结果是否足够好与原始标签有关？
与原始标签相关的 DBSCAN 聚类的准确率为 95.9%。
所有这些工作都是在 Knime 中完成的，结合了原生 Knime 节点，并在其工作流中嵌入了 Python 代码，从而展示了 Knime environement 的功能。
我正在思考的下一个故事，将是关于如何执行一个类似的分析，但用一个神经自动编码器达到嵌入状态。会随时通知你。
数据集成 UMAP
生命科学的数理统计和机器学习
基于图形的单细胞组学与 UMAP 的集成
作者图片
这是来自 生命科学的数理统计和机器学习 专栏的第二十篇文章，我试图以简单的方式解释生物信息学和计算生物学中使用的一些神秘的分析技术。数据整合是下一个重要的步骤，通过结合多种信息来源，利用协同效应来提高分析准确性。在计算生物学和生物医学中，数据集成在单细胞研究领域取得了特别的进展。去年， 自然 认定单细胞多模态组学集成为 2019 年度 方法。在本文中，我将展示如何通过使用图交集方法和图和 UMAP 来执行单细胞组学( scOmics )数据的基于图的跨模态集成。
数据集成背后的理念
当组合来自不同统计分布的数据时，简单地将它们连接在一起而不考虑它们各自的分布可能不是一个好主意，因为它们中的一些可能是二元的、分类的或连续的。使用 scOmics 数据，人们经常处理连续的 scRNAseq (转录组)/ scProteomics (蛋白质组)和二元的 scBSseq (甲基化)/ scATAseq (开放染色质区域)数据。假设个体数据类型能够以 78%、83%和 75%的准确率(数字是编造的)区分细胞与患病个体和健康个体。然后，将多个 scOmics 组合在一起的一个期望将是预测精度的提升(例如高达 96%)，这是由于特定技术噪声的降低和 scOmics 信号一致性的增强。
作者图片
组合单个 scOmics 可以有不同的方法。一种可能的方法是将 scOmics 转换到一个公共的非参数空间，在那里它们失去了关于其原始技术的记忆。这就是人工神经网络(ANN)如何能够压缩多个信息源，我在之前的一篇文章中详细描述了这一点。【SNF】相似度网络融合也属于这一类方法。另一种方法是明确地对单个数据分布建模，并使用贝叶斯规则将它们组合成一个 联合概率。最后，人们可以尝试提取共同的变化(即理清数据中特定技术和共享的变化)并将其分解(分解成正交分量)以便更好地解释，这是 PLS 、 CCA 和因子分析使用的方法。
基于图形的数据集成
基于图表的数据整合属于上一节中的第一种方法，即当转换为图表时，单个数据集表示数据点之间的成对连接，而不需要任何关于生成单个数据集的统计过程的“记忆”。在图形空间中，通过保持各个数据集的图形上的数据点之间的边缘一致，可以直接找到各个数据集的各个图形之间的交点。这种方法的优势在于，当将原始数据转换为图形时，可以应用适当的 距离度量，即，在处理二进制数据时，可能希望使用汉明距离来计算数据点之间的成对连接，而在处理连续数据时，应用欧几里德距离可能就足够了。
为了展示基于图的数据整合的思想，我们将把它应用于 scNMT 单细胞组学数据集，该数据集包含:1)基因表达(scRNAseq)，2) DNA 甲基化(scBSseq)，以及 3)来自相同生物细胞的开放染色质区域(scATACseq)信息。我们将从读取 scNMTseq 数据开始，即 scRNAseq、scATACseq 和 scBSseq，过滤数据，例如通过去除 scRNAseq 中低表达的基因，应用对数变换以使数据看起来更正态分布，并检查数据分布。
作者图片
单细胞基因表达(scRNAseq)数据在某种程度上是零膨胀的，并且明显是连续的，因此使用欧几里得距离计算细胞之间的成对相似性可能是合适的。然而，当我们绘制 scBSseq 和 scATACseq 数据分布时，很明显它们遵循双峰分布，因此使用欧几里得距离来测量使用这种类型的二进制数据的细胞之间的相似性是不合适的。因此，我们将对 scBSseq 和 scATACseq 使用汉明距离，而骰子和 Jaccard 相似度也是一个可能的选择。
作者图片
如果我们查看 scBSseq 数据矩阵，元素似乎在 0%和 100%之间变化，这意味着某个 CpG 位点可以是甲基化的(100%)或未甲基化的(0%)。实际上，并非矩阵的所有元素都是严格的 0%或 100%，这可能是由于技术采样和排序的原因。但是，我们可以将 scBSseq 数据(同样适用于 scATACseq 数据)视为二进制数据。为了使数据显式为二进制，我们将把所有低于 50%的计数编码为 0，高于 50%的计数编码为 1。然后对 scATACseq 数据应用类似的二值化过程。您不一定要这样做，但为了解释方法，使用图相交方法演示连续数据 scRNAseq 如何与完全二进制数据 scBSseq 和 scATACseq 集成是很有趣的。
现在，我们将使用来自 scRNaseq 数据集的原始基因表达数据构建一个 KNN 图。我们将在单元之间应用欧几里得相似性度量，并为每个单元保留 30 个最近的邻居。事实上，当从全连通图开始时，scikit-learn 中的 KNN 方法将检查每一对数据点的边的权重，并为每个数据点选择 30 个权重最大的邻居。结果我们会得到一个 KNN 邻接矩阵 ，其中 1 表示单元格连通，0 表示不连通。要了解图形构造的其他技巧，请查看来自瑞本菲塔斯的本精彩笔记本。
作者图片
这里，我们选择了一对随机的点，例如 ESC_F04 和 ESC_F07 ，它们似乎在 scRNAseq 数据集中是相连的。我们将在 scBSseq 和 scATACseq 数据集中跟踪这种联系。下面我们将使用 MDS 布局可视化 KNN 图，并通过红边突出显示 ESC_F04 和 ESC_F07 单元格之间的连接。我们现在已经可以观察到 ESC 细胞和 EB 细胞彼此分开的某种聚集。现在，我们将对 scBSseq 数据进行相同的图形构建，详情请查看完整的 jupyter 笔记本。在这种情况下，由于 scBSseq 是一个二进制数据，因此像我们对 scRNAseq 数据所做的那样，用欧几里德距离来度量细胞之间的相似性不是最佳的。这里我们选择汉明距离，尽管 Dice 和 Jaccard 距离可能也适用于分类/二进制数据。
作者图片
作者图片
太好了！现在，我们已经从 scRNAseq、scBSseq 和 scATACseq 构建了 3 个 KNN 图，并准备相交这些图。正如在Wolfram Math Worldpublic resource 中完美展示的那样，图交集的思想是保持多个 scOmics 数据集之间的边一致。我们将使用在 igraph 库中实现的交集函数。在我们执行了图形交集之后，一些节点将被隔离。这里，为了简单起见，我们将忽略/删除这些节点，并且在稍后运行基于 Leiden 图的集群时不将它们视为单独的集群。
来自Wolfram Math World Public Domain的图形交集插图
让我们将通过与相交获得的共有图可视化，这 3 个图来自单独的 scOmics，并证明 ESC_F04 和 ESC_F07 细胞之间的连接/边缘仍然存在，我们将用红色突出显示该边缘。现在我们可以运行 Leiden 基于图的聚类算法(对 Louvain 算法的修改)，该算法在共识图上找到 4 个聚类。
作者图片
我们可以看到，EB 细胞似乎形成了一个单独的(绿色)簇，而 ESC 细胞现在分裂成至少两个(蓝色和红色)簇。有一小组黄色细胞被 Leiden 发现为第四个群，但是这可能应该用 Leiden 算法的分辨率参数来调整，因为黄色细胞在视觉上似乎不是一个单独的群。另一方面，这也可能是布局的问题。尝试另一种图形布局可能会证明黄色细胞确实形成了一个独特的集群。
单细胞组学整合的 UMAP
事实证明，手动构建 KNN 图，然后按照我们在上一节中所做的方式对各个 scOmics 的图进行交集，并不是真正必要的，因为这个选项似乎基本上包含在 UMAP 降维算法中。的确，正如原 UMAP 文章中所说:
利用数据 UMAP 的模糊单纯集表示可以潜在地扩展到支持具有不同数据类型的数据集的降维。每种数据类型都可以被视为底层结构的另一种视图，每种数据类型都有不同的关联度量，例如，分类数据可能使用 Jaccard 或 Dice 距离，而顺序数据可能使用哈曼距离。每个视图和度量都可以用来独立地生成模糊单纯集，然后可以将这些模糊单纯集相交在一起，以创建单个模糊单纯集用于嵌入。
作者图片
因此，UMAP 允许快速构建单个 scOmics 数据集的模糊单纯集(类似于构建 KNN 图)及其交集，只需几行代码，参见关于使用 UMAP 混合不同数据类型的有趣讨论此处。
为了演示 UMAP 如何用于整合 scOmics 数据，我们将交叉来自单细胞 RNAseq 和来自 8617 个血液单核细胞(CBMCs)的蛋白质水平数据集的图，这些数据是通过CITE-seq 技术测序的，我们将制作一个一致性 scOmics UMAP 图。同样，想法是将两个 scOmics 数据集转换到非参数图形空间(模糊单纯集),其中它们忘记它们的技术差异，并在该非参数空间中相交图形。产生的共识/整合图将用于构建 UMAP 嵌入。
然而，我们将从显示两个单独的 scOmics 数据集的 UMAP 嵌入开始。作为最近邻的最佳数量，我们将选择 sqrt(N) = sqrt(8617) = 93，这一选择背后的动机在这里中已经涉及，并且与格子上的随机漫步的最佳大小有关，超过这个大小，随机漫步的“记忆”就会丢失。
作者图片
作者图片
UMAP 维数缩减图上的颜色对应于血液单核细胞(CBMCs)的 12 个细胞群。现在我们将构建两个图/ 模糊单纯集的交集，并可视化由此产生的嵌入。
作者图片
在生成的嵌入图中，从两个 scOmics 的相交中，可以得到一些有趣的观察结果。例如，先前在 scRNAseq UMAP 图上同质的灰色聚类分裂成三个子聚类。另一方面，以前在 scRNAseq UMAP 图上不可区分的紫色细胞群现在形成了与蓝色细胞分开的明显的簇。这两个结果似乎都是通过将蛋白质水平的信息纳入到所得的共有 UMAP 嵌入中来驱动的。
摘要
在本文中，我们了解到基于图形的数据集成是一种有趣的方式，它可以将多个信息源结合起来，以提高分析的准确性。将原始数据转换到具有适当的个体相似性度量的非参数图形空间中，允许仅保留成对连接信息，同时忽略数据中的特定于技术的噪声。以下图形的交集(或 UMAP 术语中的模糊单纯集)提供了数据集成的直观可视化和评估。UMAP 似乎是数据集成的理想之选，因为它只需几行代码就可以指定各个距离度量和图形的交集，即混合不同种类的数据。
在下面的评论中，让我知道哪些来自生命科学的分析技术对你来说似乎是 T2 特别神秘的，我会在以后的文章中介绍它们。在我的 Github 上查看帖子中的代码。在 Medium 关注我，在 Twitter @NikolayOskolkov 关注我，在 Linkedin 关注我。在下一篇文章中，我们将讨论如何在 UMAP 空间聚集，敬请关注。
PyTorch 中多任务学习的不平衡数据加载
在多个不平衡数据集上训练多任务模型的 PyTorch 实用指南
解决多任务学习 (MTL)问题需要独特的培训设置，主要是在数据处理、模型架构和绩效评估指标方面。
在这篇文章中，我将回顾数据处理部分。具体来说，如何在多个数据集上训练多任务学习模型，以及如何处理高度不平衡数据集的任务。
我将分三步描述我的建议:
- 将两个(或更多)数据集合并成一个 PyTorch *数据集。*该数据集将成为 PyTorch 数据加载器的输入。
- 修改批制备过程，以在每批中产生一个任务，或者在每批中混合来自两个任务的样品。
- 通过使用批量取样器作为数据加载器的一部分来处理高度不平衡的数据集。
我只审查了数据集和数据加载器的相关代码，忽略了其他重要模块，如模型、优化器和指标定义。
为了简单起见，我使用一个通用的两个数据集的例子。然而，数据集的数量和数据的类型不应该影响主设置。我们甚至可以使用同一个数据集的几个实例，以防同一组样本有多组标签。例如，具有对象类别和空间位置的图像数据集，或者具有每个图像的面部情感和年龄标签的面部情感数据集。
一个 PyTorch 数据集 类需要实现__getitem__()函数。该函数处理给定索引的样本获取和准备。当使用两个数据集时，就可能有两种不同的方法来创建样本。因此，我们甚至可以使用单个数据集，获得具有不同标签的样本，并改变样本处理方案(输出样本应该具有相同的形状，因为我们将它们作为批量张量进行堆叠)。
首先，让我们定义两个数据集:
我们定义两个(二进制)数据集，一个有 10 个 1 的样本(平均分布)，第二个有 55 个样本，50 个数字 5 的样本和 5 个数字 5 的样本。这些数据集仅用于说明。在真实的数据集中，你应该既有样本又有标签，你可能会从数据库中读取数据或者从数据文件夹中解析数据，但是这些简单的数据集足以理解主要概念。
接下来，我们需要定义一个 数据加载器 。我们用我们的 concat_dataset 提供它，并设置加载器参数，比如批量大小，以及是否重排样本。
这一部分的输出如下所示:
tensor([ 5., 5., 5., 5., -5., 5., -5., 5.])
tensor([5., 5., 5., 5., 5., 5., 5., 5.])
tensor([-1., -5., 5., 1., 5., -1., 5., -1.])
tensor([5., 5., 5., 5., 5., 5., 5., 5.])
tensor([ 5., 5., 5., 5., -5., 1., 5., 5.])
tensor([ 5., 5., 5., 1., 5., 5., 5., -1.])
tensor([ 5., 5., 5., 5., -1., 5., 1., 5.])
tensor([ 5., -5., 1., 5., 5., 5., 5., 5.])
tensor([5.])
每一批都是来自我们的 concat_dataset 的 8 个样本的张量。顺序是随机设置的，样本是从样本池中选择的。
直到现在，一切都相对简单。这些数据集被组合成一个数据集，并从两个原始数据集中随机选取样本来构建小批量。现在让我们试着控制和操作每批样品。我们希望在每个小批量中只从一个数据集获取样本，每隔一批在它们之间切换。
这是一个 BatchSchedulerSampler 类的定义，它创建了一个新的样本迭代器。首先，通过为每个内部数据集创建一个 RandomSampler 。第二种方法是从每个内部数据集迭代器中抽取样本(实际上是样本索引)。因此，构建一个新的样本索引列表。使用 8 的批量意味着我们需要从每个数据集获取 8 个样本。
现在让我们使用一个新的数据加载器运行并打印样本，它将我们的 BatchSchedulerSampler 作为输入采样器(使用采样器时，shuffle 不能设置为 True )。
输出现在看起来像这样:
tensor([-1., -1., 1., 1., -1., 1., 1., -1.])
tensor([5., 5., 5., 5., 5., 5., 5., 5.])
tensor([ 1., -1., -1., -1., 1., 1., -1., 1.])
tensor([5., 5., 5., 5., 5., 5., 5., 5.])
tensor([-1., -1., 1., 1., 1., -1., 1., -1.])
tensor([ 5., 5., -5., 5., 5., -5., 5., 5.])
tensor([ 1., 1., -1., -1., 1., -1., 1., 1.])
tensor([5., 5., 5., 5., 5., 5., 5., 5.])
tensor([-1., -1., -1., -1., 1., 1., 1., -1.])
tensor([ 5., -5., 5., 5., 5., 5., -5., 5.])
tensor([-1., 1., -1., 1., -1., 1., 1., -1.])
tensor([ 5., 5., 5., 5., 5., -5., 5., 5.])
tensor([ 1., -1., -1., 1., 1., 1., 1., -1.])
tensor([5., 5., 5., 5., 5., 5., 5.])
万岁！！！
对于每个小批量，我们现在只获得一个数据集样本。为了对更重要的任务进行降采样或升采样，我们可以尝试这种调度方式。
我们批次中剩下的问题来自第二个高度不平衡的数据集。这是 MTL 经常出现的情况，有一个主任务和一些其他的卫星子任务。一起训练主任务和子任务可能会提高性能，并有助于整体模型的一般化。问题是子任务的样本通常非常稀疏，只有几个正(或负)样本。让我们使用我们以前的逻辑，但是也强制在每个任务中关于样本分布的平衡批处理。
为了处理不平衡的问题，我们需要用一个ImbalancedDatasetSampler替换 BatchSchedulerSampler 类中的随机采样器(我使用了这个资源库中的一个很好的实现)。这个类处理数据集的平衡。我们也可以混合使用 RandomSampler 用于一些任务，使用不平衡数据采样器用于其他任务。
我们首先创建ExampleImbalancedDatasetSampler，它继承了 ImbalancedDatasetSampler ，只修改了 _get_label 函数来适应我们的用例。
接下来，我们使用BalancedBatchSchedulerSampler，它类似于前面的 BatchSchedulerSampler 类，只是将 RandomSampler 用于不平衡任务的用法替换为examplembanceddatasetsampler。
让我们运行新的数据加载器:
输出如下所示:
tensor([-1., 1., 1., -1., -1., -1., 1., -1.])
tensor([ 5., 5., 5., 5., -5., -5., -5., -5.])
tensor([ 1., 1., 1., -1., 1., -1., 1., 1.])
tensor([ 5., -5., 5., -5., -5., -5., 5., 5.])
tensor([-1., -1., 1., -1., -1., -1., -1., 1.])
tensor([-5., 5., 5., 5., 5., -5., 5., -5.])
tensor([-1., -1., 1., 1., 1., 1., -1., -1.])
tensor([-5., 5., 5., 5., 5., -5., 5., 5.])
tensor([ 1., -1., 1., 1., 1., -1., 1., -1.])
tensor([ 5., 5., 5., -5., 5., -5., 5., 5.])
tensor([-1., -1., -1., -1., 1., 1., 1., 1.])
tensor([-5., 5., 5., 5., 5., 5., -5., 5.])
tensor([-1., 1., -1., 1., 1., 1., 1., 1.])
tensor([-5., -5., 5., 5., -5., -5., 5.])
不平衡任务的小批量现在更加平衡了。
这个设置还有很大的发挥空间。我们可以以平衡的方式组合任务，通过将 samples_to_grab 设置为 4，这是批量大小的一半，我们可以获得一个混合的小批量，每个任务中有 4 个样本。为了对更重要的任务产生 1:2 的比率，我们可以为第一个任务设置 samples_to_grab=2 ，为第二个任务设置 samples_to_grab=6 。
就是这样。完整的代码可以从我的库下载。
疏通 Power BI——向 Power BI 部署渠道(BIOps)的飞跃(Ep。6)
使用 Power BI 自动化报告生命周期(BIOps)
Power BI 最近在 preview 中推出了一项名为部署管道的新功能，以自动化 Power BI 报告和内容创建(Ft。Narisorn Limpaswadpaisarn
上一集
- 解除 Power BI 封锁—Power BI 报告服务器(内部部署)可以导出 CSV 和 Excel 格式的数据吗？(Ep。1)
- 一天内解除电力 BI-仪表板阻塞(Ep。2)
- 疏通 Power BI —如何获得 Power BI 桌面本地与云？(Ep。3)
- 疏通 Power BI —让我们明确 Power BI 产品组合！(Ep。4)
- 解除 Power BI-Power BI Premium 作为分析服务的超集(Ep。5)
Power BI Service 有一个工作区的概念，您可以在其中创建仪表板、报告和分页报告的集合，并与工作区内的任何人进行协作。点击了解有关 Power BI 工作区[的更多详情。在过去，我听说过客户需要为 BI 项目生命周期管理分离出开发、测试和生产工作空间的需求。我们将这个概念称为“分阶段部署”，这是一个“隔离适合流程每个阶段的数据、用户和内容”的框架。根据](http://Create the new workspaces in Power BI)【掌握权力毕】，分阶段部署的步骤包括:
- 创建开发应用工作区，并将 Power BI 内容构建到工作区中。
- 一个应用程序被发布或更新，并提供给少数用户审阅。
- BI 经理或项目负责人审查正在开发的内容的状态，并向开发人员提供反馈。
- Power BI REST API 用于将完整的内容从开发工作区迁移到测试工作区。
- 发布或更新测试应用程序，并提供给少量用户进行审查。
- 用户验收测试(UAT)用户或团队审查与需求相关的内容，并提供反馈。
- Power BI REST API 用于将批准的内容从测试工作区迁移到生产工作区。
- 发布或更新生产应用程序，并提供给用户组使用。
- 业务用户组通过生产应用程序从任何设备访问和使用仪表盘和报告。
下图和九步流程描述了分阶段部署生命周期的主要元素。来源:掌握权力 BI(布雷特鲍威尔)
然而，我对 Powell 方法的经验是，高级 BI 管理员需要创建个人工作区来定义开发、测试和生产，这在技术上是可能的，但还不够直观。我很高兴看到“ Power BI 部署管道的公开预览版，这是 BI 团队在 Power BI Premium 中管理内容生命周期的一种新的改进方法。”简而言之，我将这一特性称为“ BIOps ”，因为 Power BI 可以自动化商业智能报告生命周期，并以前所未有的速度提供更新。
先决条件
您必须:
- 使用 Power BI Premium 许可证。
- 成为新工作空间体验的管理员。简单来说就是不支持经典工作空间。
要开始使用，请前往https://app.powerbi.com/并使用您的 Power BI 帐户登录。创建一个工作空间并获取一些 Power BI 样本。许多样品都很容易获得:
选择一些样本到选定的工作区。
创建部署管道(在预览中)
- 创建管道
- 分配您的工作空间
- 开发和测试你的内容
- 与您的用户分享
这是“创建部署管道(预览)”的登录页面
对于部署管道，不需要为开发、测试和生产分别创建 3 个工作区。相反，Power BI 将为您自动完成这一过程。只需单击“部署到测试”和“部署到生产”。单元测试相应地运行。
我只需要为 1 个管道附加 1 个工作空间。当您确认部署到测试或部署到生产时，将自动创建后面的阶段。您的新工作区将保持相同的命名约定，除了末尾的标签指定它是开发工作区、测试工作区还是生产工作区。
观察
- 工作空间 必须 是新的工作空间体验。经典的工作经验是行不通的，没有例外。
- 仅为 1 个部署管道分配 1 个工作空间。
- 只有在您愿意的情况下，您才能部署数据集。如果存在依赖关系，报告可能无法部署。我建议部署一组数据集和报告。
- 在撰写和发布本文时， Power BI 数据流尚不支持部署管道。
- 此功能仍在公开预览中；因此，还不建议将其用于生产。我会密切关注这个空间，因为它绝对很棒。
参考
免责声明:以下内容并非本人雇主官方背书。本文表达的观点和意见仅代表作者的观点和意见，不一定反映当前或以前的雇主、组织、委员会、其他团体或个人的官方政策或立场。本文中进行的分析是基于有限日期的开源信息。分析中所做的假设并不反映任何前任或现任雇主的立场。
用近似梯度法解套索正则化
迭代软阈值算法
我首先告诉你，我不是数学家，因为这个话题需要一些数学知识(振作起来！)我想说明的是，它不会完全严谨。
本文旨在总结应用 ISTA(迭代软阈值算法)解决梯度下降优化套索正则化问题的整个过程。
当然，我假设如果你正在阅读，至少你对我们将要谈论的内容有所暗示，我希望我的帖子会有所帮助。
让我们首先构建一个数学难题，它将引导我们理解如何使用梯度下降计算 lasso 正则化，即使成本函数不是完全可微的，如 Lasso 的情况。
子梯度
首先要掌握的概念是凸函数的定义。
一个完全可微的函数 f 被说成是凸的如果:
换句话说，其局部 线性逼近 (泰勒展开)总是低估 f 。这种凸函数的一个例子是我们钟爱的均方误差。
让我们更进一步，从凸函数直接得出次梯度，当我们处理一个在每一点都不可微的函数时，这个工具是必要的，而且由于套索使用了我们变量(绝对值)的 L1 范数 我们知道这个函数在零点 不可微 。
次梯度 是一种概括凸函数在不可微点处梯度的方法。凸函数 f 在 x 点的次梯度是任意的 g 这样:
该表达式也是对 f 的线性近似，其中我们使用次梯度来代替梯度，次梯度在满足上述条件的值的 集合 内。
凸函数 f(x)的次梯度。点 x1 是可微的，我们有一个唯一的梯度，而在点 x2(不可微的)我们有一组符合上述条件的可能的子梯度值
g 1 ，g 2 ，g3是* f 在x1*，x* 2 的子渐变。***
在函数可微的点(x 1 )上，次梯度是 唯一的 确切地说是:g=∇f(x)、 相反，在 x 2 上，我们可能有无限数量的次梯度(在上图中，我们只能看到其中的几个)
在 f(x)=|x|的情况下，x = 0 处的子梯度 g 如下图所示，它们的范围在[-1，1]中包含的所有可能值之间。任何一条斜率在此范围* 内穿过 x = 0 的直线，都会对函数 进行下界。***
x=0 时 f(x)=|x|的一些次梯度
凸 f 的所有次梯度的集合称为次微分:
次微分定义可以推广到各种凸函数，无论可微与否。主要的一点是，当一个函数在它的定义域中的每一点都是完全可微的，那么次微分是由整个梯度集合(分别是唯一的值)组成的，相反，当一个函数不是处处可微时，次梯度就来了。
再次考虑单变量 x 的 L1 范数:
绝对值函数(左)及其次微分∂f(x)作为 x 的函数(右)
f(x) = |x|的次微分；这种情况下 k=1，2，3
给定次微分，则任意 f (可微或不可微)的最优性条件为:
*x 是极小元当且仅当 0 是fat***x **的次梯度。这叫做 次梯度最优性条件 。
为什么？简单，g= 0 是次梯度意味着对于所有 y:
注意一个凸且可微的函数 f 的含义:
近端标测
现在我们已经熟悉了次梯度，我们需要了解的其他工具有:近似算子和软阈值算子。
近端操作符 定义:
这里我们搜索的是点*【x 】，它最小化了一个通用凸函数 f*****但是 同时保持接近一个参考点 我们通过α选择接近多少作为这些术语之间的相对权重或权衡参数*****
软阈值算子 定义:
或者更简洁的形式:
例如给定λ= 1，如果 x = 3，sign(3) = 1，max((| 3 |–1)，0) = max(2，0)= 2->1 * 2 = 2；如果 x = -0.5，sign(-0.5) = -1，max((|-0.5 |–1)，0) = max(-0.5，0)= 0->-1 * 0 = 0**；如果 x = -3，sign(-3) = -1，max((|-3 |–1)，0) = max(2，0) = 2 - > -1 * 2 = -2**
**好吧，我知道，你在想: ***为什么是他们？？？*我要啤酒不要操作员！
我理解你的兄弟姐妹。但是在你喝啤酒之前，你需要经历最艰难的部分。继续战斗！
那么，让我试着解释一下所有这些，并开始应用近似运算符。
定义 f 为局部线性逼近在xk对于一个通用函数，我们假设* 凸且可微，我们把它代入最接近的算子:***
已知，点 x 是局部极小点如果:*
我敢打赌这提醒了你一些事情……等等，我用一种更熟悉的符号来写:
现在好些了吗？是的，这是一个梯度下降步骤，仅用线性近似的近似映射来计算。事实上，梯度下降和近似算子是非常亲密的朋友。
下一轮我们将使用不同的函数:
同样，我们可以将它插入到 prox 操作符中，我们得到:
从 L1 和 L2 范数的定义我们得到:
这些是可分离的函数，因此我们可以将它们作为组件式标量问题来解决:
由于 h 是凸的，但在 x = 0 处不可微，我们用它的次微分来工作并回忆:
于是我们得到了最后的结果:
或者用上面介绍的软阈值运算符重新排列:
最后这个结果可能有点棘手，让我解释一下是怎么回事。根据绝对值的次微分和上面获得的最优性条件，我们有 3 种可能的情况:
所以最后我们可以概括为:
****
图片来源皮埃尔·阿布林
开始吧，最后一轮…阿德里亚娜！
现在考虑:
作为一个一般的凸可微函数 g(即最小二乘)加上一个凸不可微函数 h 的和，在我们的情况下是 L1 范数，我们想最小化它。为此，我们首先将一个 局部二次近似 应用于 g，并保留 h 不变。
****
在 x=1 时，具有不同 alphas 的一般 g(x)的局部二次近似
注意，我们在没有 Hessian 的情况下获得的局部二次近似只是我们之前计算出的线性近似的近似算子。
为了我们的目的，我们还需要重新排列 g:
在去掉不依赖于 x 的常数项并加上 h 后，我们可以最小化整个 f:
经过大量的数学计算，最终我们可以将 ISTA (迭代收缩阈值算法)定义为:
对于我们的 n 变量可分函数的每个分量。
是的，看起来很糟糕，但并不像看起来那么糟糕。您只需像往常一样计算梯度下降，唯一的区别是在将值传递给下一次迭代之前，使用软阈值操作符应用进一步的变换。
然后，对于梯度下降，你可以定义一个规则来停止迭代并非常接近函数的最小值。
值得一提的是，这个结果不依赖于 g，它可以像我们想要的那样复杂，我们只需要能够计算它的梯度。
搞定！从现在开始，不再有数学，只有啤酒。我们应得的。
打开解析的盒子！
图片来源:Pixabay
句法分析是分析句子中所有单词并确定它们的句法和/或语义角色的过程。就像潜意识里我们用母语中的语法来组成和理解一个句子。类似地，我们在 NLP 中使用语法来破译一个句子是如何构造的，并理解它的句法和语义角色。主要有两种解析: 选区解析 和 依存解析 。
依存解析 顾名思义就是侧重于将一个句子分解成单词或记号或短语，然后分析它们之间的二元依存关系，了解它们的语法作用。它有助于解释句子中的语义关系。
依赖关系分析能够识别的一些二进制关系如下所示:
依存分析通常使用有向弧来产生这种关系，其中弧的一端是中心词，弧指向的一端是词本身。短语(一组词)中的中心词仅仅是一个词，如果去掉这个词，短语就失去了意义。
一个短语(一组词)中的中心词就是一个词，一旦这个词被去掉，这个短语就失去了意义。
例如:
短语: 《好邻居》
在这里，如果你去掉“邻居”这个词，你就失去了这个短语的任何可能的逻辑意义，因为“好”本身没有任何意义。因此，这个短语的标题是“邻居”。
现在，有向弧和节点不会让你想起有向图吗？如果是的话，那么你是绝对正确的！我们在执行依存解析之后获得的依存树只不过是一个有向图，其中所有的单词都是顶点，并且除了根顶点之外，每个顶点或单词都只有一条引入弧。这是因为，只有依存关系树的根顶点(根单词)是没有引入弧的节点。此外，从根到图中的任何其他节点总是有唯一的路径。
现在，让我给你看一个依赖解析的例子。
罗斯取消了去迪拜的航班
现在，看看这个句子的依存解析。
这里的小句关系如 nsubj 和 dobj 把主语标识为“Rose”，把谓语“cancelled”的直接宾语标识为“flight”。请注意，“cancelled”是这两个短语中的中心词，它是图中唯一没有引入弧的节点，因此它是根。像 prep 这样的关系表明介词关系“到”的中心是“飞行”。关系 poss 是一个所有格修饰语，表明“her”是中心词“flight”的所有格修饰语，而 pobj 表明介词“Dubai”的宾语与其中心词“to”相连。
简而言之，我们现在可以看到给定句子中的每个单词是如何在语义上构造的，并且我们对给定句子中所有可能的语义上可推导的二元关系的信息是有限的。因此，依赖解析的一个巨大优势是它提供了谓词和它们的参数之间的语义关系的近似[1]。依存句法分析的另一个优势在于分析语序更灵活的语言的句子。换句话说，在语言中，人们可以更自由地放置某些词类，依存句法分析帮助我们获得关于文本的重要语义信息。
请注意，这只是一个包含少量依赖项的简短示例。如果你想了解更多关于这种关系的知识，请看看斯坦福大学的类型依赖手册。
现在，我们理解了依赖解析的作用，下一个重要的问题是理解它可以用在哪里。这种句法分析技术对于理解句子中的语义关系至关重要的应用，如共指消解、问题回答、信息抽取等，具有重要意义。
下面给出了一个使用 spacy 查找句子依存解析的例子。
依存性分析说明
选区解析 另一方面涉及到考虑关于句子的句法信息。它将一个句子分解成短语，并根据预先确定的语法和词典(或词汇)组合短语以形成一棵树，称为解析树。这些短语可以是名词短语(NP)、动词短语(VP)、状语短语、形容词短语和/或介词短语(PP)。因此，选区分析具有作为语法分析器的角色。
在我们继续之前，让我快速地告诉你什么是 NLP 中的语法。自然语言处理中的 语法 是任何语言中所有规则的正式规范。在本文中，我们将使用英语作为我们的主要语言，以帮助简化解释！
例如:
名词短语(NP): “牛逼的车”
但是这个短语不也是由其他词类构成的吗？
“the”: 冠词/限定词；“牛逼”:形容词；【汽车】:名词
我们如何用语法规则来表达这一点呢？
我们说， NP → Det Adj 名词
在这里，名词、形容词和名词被称为词尾，因为它们直接表示词典或词汇中的单词。在词典中可能是这样的:
Det → the
Adj →牛逼
名词→汽车
因此，我们英语语法中的所有规则都可以用类似的方式来表示。那么，语法类型是只有一种还是有多种呢？答案是肯定的！乔姆斯基的分类主要给了我们 4 种语法。
乔姆斯基的分类
-类型 0:无限制语法
格式:x → y，其中 x 和 y 可以是非终结符和/或终结符。完全没有限制。
-类型 1:上下文敏感语法
格式:xAy → xzy，其中 x，y，z 可以是任意值，A 在上下文 xAy 中是非终结符。
-类型 2:上下文无关语法
格式:A→ x，其中 A 是非终结符，x 可以是终结符或非终结符的序列。
-类型 3:常规语法
格式 1 : A→ bt，其中 A 和 b 为非终结点，t 为终结点。
或者
格式 2 : A→ t，其中 A 为非终结符，t 为终结符。
因此，我们可以根据需要创建属于任何类型的语法，并使用它来解析句子。为了方便起见，让我们假设在本文中使用常规语法，除非另有说明。
现在，让我们回到理解选区解析上来！
选区分析器为每个句子创建一个分析树。让我们看一个例子。
句子 : “玫瑰是一种花”
语法:
- S → NP VP
- NP → NNP
- 副总裁→ VBZ NP
- NP → DT NN
词汇:
- NNP →玫瑰，百合，兰花
- VBZ →是，是
- DT → a，the，an
- NN →花
这个句子的解析树如下所示。
您可以在上面的图片中看到，解析树包含所有单词作为叶或终端，S 作为根。现在，主要有两种解析方式: 自下而上和 。在 自下而上解析 中，解析器首先在内部执行词性标注，开始确定每个单词的句法角色。然后，它使用预先确定的语法来建立两个单词之间的语法关系，并递归地组合短语，直到到达词根:S。然而，在 自上而下解析 中，解析器从 S 开始，并使用语法规则来解开每个短语，直到到达句子中的所有单词。
请记住，我们在这里讨论的是一棵 有效解析树 和一棵 有效解析树 总是包含 S 作为根，并且只包含给定句子的单词作为叶。例如，如果你的解析树包含的单词比你的句子中出现的单词多，那么它就是一个无效的解析树。
一个 有效解析树 总是包含 S 作为根，并且只包含给定句子的单词作为叶。
另一个需要记住的重要事情是，一个句子可以有多个解析树，但是其中只有一个在语法和语义上都有效。让我们看一个例子:
**句子:**我拍了一头穿着睡衣的大象。
你的预定语法和词汇或词汇如下:
L1 微型英语语法和词汇
现在，您的句子可以有两个可能的解析树，如下所示。
这两棵树看起来在语法上是有效的，但在语义上并不都是有效的。评估语义有效性时要记住的重要事情是看短语是如何组合的。在左边的树上，由于名词短语“穿着睡衣的大象”与“shot”相连，听起来作者想说他或她射杀了一头穿着睡衣的大象。然而，如果你仔细观察右边的解析树，动词短语(VP)“射杀一头大象”与介词短语(PP)“穿着我的睡衣”组合在一起，这传达了作者穿着他或她的睡衣射杀一头大象的意思。因此，句子的正确含义仅由右边的解析树传达，并且只有该解析树才是正确的解析树。
在 NLP 文献中，这个问题被称为结构歧义问题，其中一个句子可以有多个解析树，但其中只有一个是准确的。
现在，让我们看看如何使用 python 来实现这一点。我使用 nltk 的 CoreNLP 解析器，它的后端有斯坦福解析器。请注意，为了让这段代码工作，你必须在你的本地计算机上下载并启动斯坦福 NLP 服务器。
选区分析说明
您一定已经注意到，我只讨论了单独利用语法和词典来形成解析树的解析技术。使用其他可能给我们带来更好结果的方法怎么样？使用自底向上或自顶向下解析的动态编程的选区解析的两个例子是 CKY 解析和厄尔利解析。
CKY 语法分析器使用上下文无关语法和自下而上的动态编程方法为任何给定的句子构建一个语法分析树。
而 Earley Parser 使用任意上下文无关语法和动态编程自顶向下的方法为任何给定的句子构建解析树。
现在，你一定想知道为什么不把概率的概念引入我们的方法中。我们都知道，一些语法规则比其他规则出现的频率更高，因此，与它们的对应规则相比，我们可以给它们更高的概率，从而形成一个概率语法，然后可以用来执行解析。现在，我们可以简单地通过乘以用于构建所述树的所有规则的概率来计算最终解析树的概率，并且可以选择具有最高概率的解析树作为给定句子的精确解析结构。
概率解析 使用出现概率较高的语法规则进行解析，找到最可能的解析树。
这个想法仍然有它的缺陷。这种解析对它展开的词的种类不敏感。
例如:
第一句:我吃了健怡可乐炸鸡。
**第二句:**和朋友吃了炸鸡。
在第一句中，“用健怡可乐”是 NP(名词短语)，由另一个 NP:“健怡可乐”和 PP(介词短语):“用”。换句话说:
规则 1 : NP → NP PP
规则 1 的概率: 0.5。(让我们假设)
在第二句中，“with friends”是 VP(动词短语)，由另一个 NP:“friends”和 PP(介词短语):“with”。换句话说:
规则二: VP → VP PP (0.5)
规则 2 的概率: 0.5(假设)
现在，想象你是一台计算机，你正试图找到哪条规则来解析你的句子。因为两个规则的概率是相等的，所以您可以选择其中任何一个并得到一个语法上有效的解析树。但是，第一句和第二句的意思不一样。第一句暗示你吃了炸鸡和健怡可乐，而第二句暗示你和你的朋友吃了炸鸡。因此，概率分析器仍然缺乏区分两个规则的概率是否相等的能力，从而使它对单词的词汇信息不敏感。
因此，下一步是使用一个考虑了概率解析器忽略的词汇信息的解析器。答案是使用词汇化的概率分析器。
词汇化概率解析器使用关于单词的词汇信息(例如，使用关于每个短语的头部的信息)以及规则出现的概率来解析句子。
还有其他种类的解析器，它们利用自然语言的不同方面来解析一个句子，但是在本文中，我只关注最流行的几种。本文中使用的 python 代码可以在这里找到。
我希望这篇文章能帮助您更好地理解解析！
感谢您的阅读！
附言:如果你想让我写一些特定的话题，请在下面评论！谢谢大家！😃
参考文献:
不确定性感知强化学习
窥探知道自己在做什么的建筑代理
基于模型的强化学习(RL)最受样本效率的青睐。它对期望的输入量要求不高，并且对我们期望模型实现的目标有一个上限。
这个模型不太可能完美地再现环境。当通过训练有素的代理与现实世界互动时，我们可能会遇到与训练中看到的不同的状态和奖励。为了让基于模型的 RL 发挥作用，我们需要克服这个问题。这至关重要。这有助于我们的特工知道自己在做什么。
第一，无模型 RL 怎么样？无模型 RL 总是在训练和测试代理时使用地面实况环境转换。除非有我们引入的偏移，比如模拟到真实的转移，在这种情况下，我们不能责怪算法。因此，不确定性不是一个大问题。对于像 Q 函数 *Q(s，a)*这样优化动作的东西，我们可以尝试在动作选择上集成确定性意识。但是既然它工作得很好，就目前而言，闭上我们的眼睛，假装我们没有看到，这并没有什么坏处。
内容
1.基于模型的 RL
2 的不确定性来源。不确定性意识的好处
3。构建不确定性感知模型
-什么看起来可行
-什么可行
4。结论
不确定性的来源
模型不确定性是由模型在测试期间看到的数据和用于训练模型的数据之间的分布不匹配造成的。我们在不同于培训期间看到的分布上测试代理。
确切地说，不确定性意识会产生什么样的影响？
在训练开始时，模型p*(sₜ₊₁|sₜ、* a 、ₜ已经暴露了自己相当小的真实世界数据。我们希望函数不要过度适合这个小数量，因为我们需要它足够有表现力来捕捉后面时间步骤中的转换。然后，将积累真实数据来学习精确的模型。
这在基于模型的 RL 中很难实现。为什么？RL 的简单目标是最大化未来的累积回报。计划者以此为目标，试图遵循模型预测高回报的计划。因此，如果模型高估了特定行动序列将获得的回报，规划者将很高兴遵循这一闪闪发光但错误的估计。在现实世界中选择这样的行为会导致滑稽的行为。简而言之，计划者有动机利用模型的积极错误。
(我们可以把计划器看作是我们在给定世界状态的情况下选择最佳行动的方法)。
这还能更糟吗？在高维空间中——比如输入是一幅图像——由于潜在变量，模型会犯更多的错误。在基于模型的 RL 中，通常通过使用策略上的数据收集来缓解分布不匹配，在现实世界中观察到的转换被添加到训练数据中，并用于重新规划和纠正模型中的偏差。然而，在这种情况下，错误太多，政策上的修正无法赶上失去的模式。大量的错误可能会导致我们每次重新计划时都改变策略，结果，模型可能永远不会收敛。
我们可以选择为模型可能犯的每一个错误收集数据，但如果我们能够检测到模型可能出错的地方，那么基于模型的规划者就可以避免可能导致严重后果的行动，这不是更好吗？
估计不确定性
首先，让我们用一个简单的故事来表达我们所知道的。
一对恩爱的夫妇得到了一个婴儿和一个机器人的祝福——不一定同时得到。作为一名保姆，机器人的目标是让小朱丽叶开心。虽然实现这一目标受到奖励的激励，但机器人避免任何损坏或可能伤害婴儿的事情也是可取的。
现在，婴儿变得喜欢在指着虫子的时候哭泣——因为好孩子会这样做——机器人的最佳奖励计划变成压扁虫子，让婴儿朱丽叶看着它把害虫喂给猫。
不过，为了改变，假设机器人在电视上指着一些可怕的东西时遇到了婴儿哭泣——这是一种不熟悉的状态，似乎接近婴儿朱丽叶的哭泣指行为。不确定动态，机器人的最佳计划可能是压扁电视，喂猫。我们不确定这是否会让宝宝开心，但这肯定会造成伤害。
然而，如果模特因为不自信而回避了这个动作，那么最好不要碰电视，以一个悲伤的朱丽叶为代价来避免伤害。
一个不确定性感知模型可以让代理人知道哪里出现不良结果的可能性高，哪里需要更加小心。但是如果模型对采取行动后的结果没有信心，那么利用它来达到目标可能是好的。
如果我们的机器人确信泡菜能让朱丽叶平静下来，而不会带来任何风险，那么它可能会考虑跑到厨房，让她咀嚼一个，因为这样一来，它就能达到让她开心的目的。
一个能够获得其不确定性的精确估计的模型，给予基于模型的计划者避免有不小可能导致不希望的结果的行动的能力。渐渐地，模型将学会做出更好的估计。不确定性意识也将告知模型它需要探索更多的状态。
有什么解决方案吗？
利用熵
我们知道熵是随机性的一种度量，或者随机变量的概率分布的扩散程度。
离散随机变量的熵 X
当每一个结果(【xᵢ】)以相等的概率发生时，熵将达到峰值，即最大不确定性，并且当存在具有高概率 p (几乎为 1)的单个结果时，熵将最小，而其余结果共享接近于零的概率( 1 — p )，即最大确定性。
这里，A 比 b 具有更低的不确定性，对于 A 来说 单峰 高斯分布 ，熵和方差(围绕均值扩散)倾向于在相同的相对方向上移动。
然而，问题就在这里——数据的不确定性会影响我们模型的不确定性吗？为了更好地表达这个问题，这里有两个图表:
你将会看到手工制作的数据的可视化图表。把注意力放在我们要用它们传达的概念上，而不是它们代表什么，好吗？谢了。
左 : 模型预测为线性函数——模型拟合数据，但是“太好了”。——————*——右 :拟合多项式回归函数
这两个情节有助于我们区分两个相关的概念。在它们之间，你认为哪一个是不确定的？抓住一个轻松的时刻再次观察它们，并尝试想出一个答案。
你选好了吗？让我给你一些思考的空间。
在第一个图中，模型对几千个样本进行了线性回归(LR)。数据似乎是确定的。过度拟合显示确定性，对吗？当然，模型对数据很有信心。红线预测的方差将非常接近于零。但是，这个估计是否给出了关于观察结果的最佳解释呢？
为此，模型需要通过在数据和可能的协变量之间建立良好的线性关系来解释观察结果。但是不确定包括哪些和省略哪些，它包括了(几乎)所有的东西！这影响了模型的不确定性，因为我们对模型没有信心。
在多项式回归(PR)中，数据中存在噪声。它在观测值中产生了不确定性。样本显示了与平均值的显著偏差，但模型相对来说拟合得很好——我们仍然可以在没有提供牛奶的情况下从心理学课上获得一些笑声，这似乎足够真实。与通过过度拟合线性函数接收的误差相比，这将给出更大的 MSE。
所以回到我们最初的问题，数据中的不确定性不会揭示我们的模型-不确定性。这使得熵，作为模型不确定性的度量，并不总是有效的，因为当数据的方差接近零时，即使模型仍然不确定，熵也是低的，如 LR 图所示。
不确定性感知模型——有效的解决方案
a)学习预测不良行为的函数
考虑无人机在雨林中学习飞行。我们希望它能学会在环境中导航，同时避免与树木相撞。在 RL 中，代理通过尝试一个动作来学习该动作的结果。所以要学会躲避树，它必须经历几次撞击。但是高速撞击肯定会造成破坏。
我们可以通过让无人机经历温和、低速的撞击来训练它，这样它就能学习森林环境。当它遇到训练分布中不存在的森林部分时，它需要了解其策略的不确定性，以便在收集新的训练数据时能够与该部分进行安全交互。一旦对该部分有信心，它就可以在未来高速飞行。这是一个安全探险的例子。
为了实现这一点，我们在 RL 成本函数 c(st，at) 中集成了击中树的成本，以使 c(st，at) + C_bad。 C_bad 是分配给导致不良行为(冲突)的行为的新成本。它影响着无人机何时可以快速飞行，何时应该小心翼翼。
为了估计 C_bad，我们使用一个不良行为预测神经网络 P ，带权重θ**。它将无人机的当前状态 st 作为输入，它是观测值oₜt19】加上一系列动作aₜ*，aₜ₊₁…ah*无人机计划执行 e 并估计碰撞发生的概率。**
动作序列由模型预测控制 (MPC)在从当前时间步长 t 直到 t + H 的后退时间范围内选择和优化。不良行为模型P**θ输出一个伯努利分布(二进制 0 或 1)，指示在该时间范围内是否发生碰撞。**
为每个水平 h 记录碰撞标签。这意味着对于标签 1，在时间步长 t 和t+h之间的子序列中发生了不良行为。利用以上述输入为条件的这个概率标签，不良行为模型可以简单地表示为:
使用 P_theta 模型估计不良行为 C_bad 的成本
类似地，一个简单的实现应该是这样的:
然而，你可能已经注意到Pθ输出的是不良行为的概率分布，而不是该行为的实际费用。因此，实际的不良行为成本将乘以这个概率 p 得到 pC_bad 。最后，我们用一个标量 λ 来调整它，这个标量决定了与实现目标相比，代理避免风险结果的重要性。**
RL 成本由目标成本和代理的不良行为加权成本构成****
值得注意的是，虽然我们想要一个预测不安全行为的函数，但一个接受输入并给出安全估计的判别模型可能并不总是让我们高兴——它的预测在不熟悉的状态下可能毫无意义。最好是，在其预测中加入模型不确定性是有益的。
b)贝叶斯网络
**神经网络可以被称为条件模型 P(y|x，w)，，其给定输入 *x，*使用权重 w，将概率分配给每个可能的输出 y 。对于贝叶斯神经网络，不是每个神经元具有单个输出，而是将权重表示为可能值上的概率分布。
右图:贝叶斯神经网络。代替每个权重的单个固定标量，概率分布被分配给所有可能的值( 源 )
这是如何工作的？
使用一组训练样本 D ，我们在模型权重上找到一个后验，以这些样本 P ( w |D)为条件。为了预测特定标签或的分布，由后验分布缩放的权重的每个可行组合对相同的输入 x 进行预测。
如果一个单位对观测值不确定，这将在输出中表示出来，因为不确定性越高的权重在预测中引入的可变性越大。这在模型发现极少数据或没有数据的地区很常见，这将鼓励勘探。随着观察越来越多，模型做出的决策也越来越确定。
对权重P(w| D)的后验分布进行近似。这是通过试图找到权重 q(w|θ)上不同分布的参数θ来完成的，通过使其尽可能接近真实的后验分布*P(w|| D)。*这是变分推论；有点超出我们目前的范围:)。
c)辍学
在 RL 辍学是一个坏主意。不过，这不是我们要冒的风险。还记得我们说过，一个有区别的模型不会总是让我们快乐，除非它能在不良行为预测中包含不确定性？辍学是一个简单的方法。
丢弃是一种正则化技术，它以概率 p、随机丢弃神经网络中的一个单元，或者以概率 1 — p 保留它。它经常在训练中使用，以防止神经元过度依赖彼此。这在每次训练迭代期间创建了新的但相关的神经网络。
在实践中，已知仅在训练期间应用漏失，并在测试时移除，以实现高测试精度。然而，通过保留测试时的漏失，我们可以通过寻找不同前向传递的样本均值和方差来估计不确定性。这是估算不确定性的简单方法。
它的警告是，作为一种变分推理方法，由于变分下限而严重低估了不确定性。
那是什么意思？
为了理解这一点，我们需要引入 KL 散度——一种对同一随机变量的两个概率分布之间差异的度量。
有时，在大的实值分布上寻找真实概率是昂贵的。因此，取而代之的是使用该分布的近似值，并且两者之间的 KL 散度(差)被最小化。
p(x)与 q(x)之间 KL 发散的一个例证——DKL(p | | q)(出处 )
在上图中， q(x) 是精确分布 p(x) 的近似值。这种近似的目的是在 p(x) 具有高概率的情况下放置高发生概率。在插图上，注意 q(x) 是单高斯，而 p(x) 是两个高斯的混合？为了在 p(x) 的概率高的地方放置高概率， q(x) 使 p 中的两个高斯分布相等，以在两者上相等地放置高概率质量。
**类似地，在输入 x 和标签 y 的条件下，漏失在模型的权重 w 上具有真实的后验概率 p(w| x，y)。q(w)被用作该后验概率的近似分布。然后我们降低 q(w) 和实际后验概率 p(w| x，y) 之间的 KL 散度，使它们尽可能接近。然而，这样做将对 q(w) 在 p(w) 没有概率质量的情况下放置概率质量进行惩罚，但是忽略 q(w) 在 p(w) 实际上具有高概率的情况下没有放置高概率质量。这就是低估模型不确定性的原因。
c)引导程序集合
训练多个独立的模型，并对它们的预测进行平均。如果这些模型近似一个几乎相似的输出，这将表明他们同意，表明他们的预测的确定性。
平均 N 个集合模型
为了使模型相互独立，每个模型的权重θᵢ用数据的子集训练，数据用来自训练集的替换采样。然而，已知训练期间权重的随机初始化θᵢ和随机梯度下降使它们足够独立。
作为不确定性的一种度量，漏失可以被假定为集合方法的廉价近似，其中每个采样漏失充当不同的模型。贝叶斯神经网络也有集成概念——通过对权重*(w| D)取后验分布下的期望，贝叶斯网络变得相当于无限个集成——越多意味着越好。*
d)好奇 iLQR
将好奇心视为解决代理环境中不确定性的灵感。让我们看看如何将好奇行为添加到代理的控制循环中。
一些 LQR 的背景
在 RL 中，线性二次型调节器(LQR)输出一个用于开发模型的线性控制器。当处理非线性动力学时，我们使用线性回归在每个时间步拟合模型(sₜ₊₁|sₜ，aₜ)。这种迭代控制过程被称为迭代 LQR (iLQR)，是差分动态规划(DDP)的一种形式。**
系统动力学由以下等式表示:
系统动态学
f 代表学习到的动力学模型而 xₜ ₊₁ 是下一时间步的状态，表示为当前状态加上模型对当前状态的预测变化 x ₜ 当采取动作 u ₜ 时。例如，如果状态是机器人的速度， x ₜ 将是当前速度，而 f(x ₜ ，uₜ)δt将是选择 u ₜ 时的预测变化，从而产生新的速度
让人好奇
对于上述系统动力学中不确定性的积分，它被写成一个高斯分布，由一个均值向量μ和一个协方差矩阵σ表示。
高斯策略具有将状态和动作对映射到状态平均变化的神经网络。这种状态的变化就是模型μ(f)表示的均值向量。**
我们可以通过从正态分布绘制模型 f 来将系统动力学实现为高斯过程(GP)，其中我们试图学习使成本函数最小化的最佳均值向量 μ 。GP 随后使用以下等式进行预测:
作为高斯过程的系统动力学
其中 f(x ₜ ，u ₜ ) 为可训练动力学函数表示的均值向量，σₜ₊₁为当前状态和动作下 GP 预测的协方差矩阵。****
这个 GP 等同于一个普通的非好奇的 LQR 随机动力学方程。有什么不同？在非好奇的 iLQR 中，由于高斯的对称性，我们将忽略方差参数σₜ₊₁。然而，好奇 iLQR 需要预测分布的协方差来确定模型不确定性。高模型不确定性等于高方差。σₜ₊₁代表当前状态和动作下模型对预测 xₜ ₊₁ 的不确定性( x ₜ ，u ₜ).
来自 GP 模型的这种不确定性随后被用于恳求代理采取行动来解决模型在这种状态下的未来不确定性。简而言之，鼓励代理人选择减少模型方差的行动。这是通过奖励代理人包含一定程度的不确定性的行为，同时仍然最大化特定目标的奖励来实现的。
山地车实验:目标是到达山顶。由 开始寻求负回报 (向后移动到对面的山上)汽车获得足够的动力到达顶部。看到只有好奇的特工才知道 探究 这个。( 来源 )
理解 基于模型的 RL 中的 LQR 优化 可能有点杂耍，但对于掌握好奇号算法是如何导出的是必不可少的。不过，我们现在不要用这些方程来吓唬自己。
奖励好奇的行为使代理比使用标准的 iLQR 更快地达到目标。防止模型陷入局部最优，在短时间内找到更好的解。
结论
不确定性意识可用于影响代理的学习策略，这取决于它如何添加到培训成本中。它可以鼓励探索或悲观行为，在这种情况下，代理人避免可能有风险的结果。后者是 AI 安全感兴趣的。为了更加稳健，增加不确定性意识的解决方案可以一起使用。
参考文献
[1] A. Nagabandi，G. Kahn，S. Fearing，S. Levine，无模型微调的基于模型的深度强化学习的神经网络动力学 (2018)，ICRA 2018。
[2] S. Daftry，S. Zeng，J. A. Bagnell，M. Hebert，内省知觉:学习预测视觉系统中的故障 (2016)，IROS，2016。
[3] C. Blundell，J. Cornebise，K. Kavukcuoglu，D. Wierstra，神经网络中的权重不确定性 (2015)，ICML 2015。
Y. Gal 和 Z. Ghahramani。作为贝叶斯近似的辍学:表示深度学习中的模型不确定性。ICML，2016。
[5]李彦宏和李彦宏。具有α发散度的贝叶斯神经网络中的脱落推理，2017。
[6]基于模型的强化学习，深度 RL 决策与控制 (2019)， Berkley 。
[7]卡恩、维拉弗洛、庞和莱文。碰撞避免的不确定性感知强化学习，2017。
[8] S .贝克特尔，a .拉伊，y .林，l .里盖蒂，f .迈耶。好奇 iLQR :解决基于模型的 RL 中的不确定性。ICML，2019。
插补的不确定性
你在预测中考虑到它了吗？
几周前，我用 DataCamp 开发的处理 R 中缺失数据的课程上线了。从目前为止我从学生那里得到的反馈来看，他们从课程中受益最大的部分是关于如何将插补的不确定性整合到建模中。在本文中，我将在使用机器学习模型进行预测的背景下简要讨论它。
为了演示这个主题，我们将分析来自fivethirtyeight R 包的活检数据集。我们将使用原始数据的一个小的经过处理的版本。让我们来看看吧。
传记片数据包含一些传记电影的信息。这些特征包括电影的生产国、发行年份、票房收入记录、电影中出现的主题数量、主题类型或认可原因以及它们的性别。这里的任务是根据剩余的特征预测对数收益。我们正在考虑拍摄我们自己的电影。这部电影将于明年在美国上映，讲述一个单一角色的故事:一名女运动员。我们想知道我们可能挣多少钱。
不幸的是，数据集并不完整:在每个特性列中都有一些缺失值。
一种方法是以某种方式估算(即填充)缺失的数据，一旦数据完整，就将其输入模型。这其实是大部分人都会做的，那我们也来做吧。我们将使用VIM包中的热卡插补对数据进行插补，并使用ranger包训练一个随机森林。如果你对 hot-deck 的工作原理、其他可用的插补方法以及如何调整这些方法使其更好地工作感兴趣，请随时查看本文底部链接的课程。
0.01210784
我们得到了 0.01210784 的单一预测。一旦我们反转对数变换，我们将会看到我们的电影有望获得 exp(0.01210784)≈100 万美元。我们是百万富翁！但是我们能有多确定呢？
上述预测完全忽略了随机森林模型适合估算数据的事实。
值得注意的是，估算值并不是一成不变的。它们只是估计值，而估计值带有不确定性。
用统计学术语来表达就是:估算值是根据观测值估算出来的，后者是从某个人群中抽取出来的；因此，存在双重不确定性:与数据采样相关的不确定性和来自插补方法的不确定性。根据这样一个不确定的预测采取行动是不明智的！
虽然没有任何方法可以神奇地将统计学中的不确定性转化为确定性，但有一些方法可以捕捉这种不确定性，并在根据模型的预测采取行动时将其考虑在内。
在我们的例子中，一种方法是自举。我们需要做的是从原始的、不完整的数据中抽取一个 bootstrap 样本。这意味着我们正在对替换的数据行进行采样。然后，我们对 bootstrap 样本进行估算，对估算的数据建立模型，最后进行预测。我们重复这个过程很多次，比如说 1000 次，随机抽取不同的样本。因此，我们获得了我们计划的电影的预测收入分布，这种分布捕捉了插补的不确定性。
由于有了boot包，这在 R 中很容易实现。我们只需要定义一个接收原始数据和要采样的行的索引作为输入的函数。然后，该函数应该选择这些行，估算 bootstrap 样本，训练模型，并返回预测值。剩下要做的就是将这个自定义函数传递给boot函数。
该输出告诉我们，与自举分布的平均值相比，我们最初的单个预测值 0.01210784 有些偏差，并且该分布具有相当大的标准误差，这意味着相当大的不确定性。
现在让我们来计算我们预测的 95%置信区间。
因此，我们有 95%的把握我们的电影将获得-1.02 到 1.18 之间的收入。把它转换回美元，我们就有了$362𝑘和$3.2𝑚.的范围相当大的差距。让我们将预测的分布、单个预测和置信区间都绘制在一个图表中。
0.01 处的橙色线表示我们根据仅估算一次的数据训练单个模型得到的预测。这一枪连最有可能的结果都没给我们！这是由于从可视化中变得清晰的偏差:最有可能的收益位于分布的峰值，比我们最初的预测稍微偏右。然后:你看分布多广！
拍摄这部电影可能仍然是值得的，但取决于你的风险厌恶，你可能会想采取预防措施，以防它没有回报，这不是不可能的。
感谢阅读！我希望你已经学到了对你的项目有益的东西🚀
你可以找到这篇文章的数据集和所有代码，包括绘图，在这里。如果您有兴趣了解更多关于缺失数据、各种插补方法以及其他将插补不确定性纳入建模的方法，请查看我的课程。
如果你喜欢这篇文章，试试我的其他文章。不能选择？从这些中选择一个:
揭秘著名的竞赛获奖算法。
towardsdatascience.com](/boost-your-grasp-on-boosting-acf239694b1) [## 线性分类器:综述
本文讨论了四个流行的线性函数的数学性质和 Python 的实际应用
towardsdatascience.com](/linear-classifiers-an-overview-e121135bd3bb) [## 用有限混合模型和 EM 算法估计经济状态
towardsdatascience.com](/estimating-the-state-of-the-economy-with-finite-mixture-models-and-the-em-algorithm-d975d280dbc6)
用贝叶斯推理评估预测的不确定性
计算贝叶斯统计导论。
Maksym Kaharlytskyi 在 Unsplash 上的照片
介绍
频数统计 有几个局限性，比如在预测中缺少**(我们通常只对期望建模)，没有内置的 正则化 ，或者没有包含先验知识。一个可以从过去的信息中获得，例如以前的实验，但是也可以从一个训练有素的 主题专家 的纯粹主观的评估中产生。**
贝叶斯统计 是统计学领域的一种替代理论，基于概率的贝叶斯解释，其中概率表示关于事件的信念或信息(知识)的程度。这与频率主义者的解释不同，后者认为概率是经过多次试验后某一事件相对频率的极限。许多贝叶斯方法需要大量的计算来完成，这就是为什么在上个世纪广泛使用的大多数方法都是基于频率主义者的解释。然而，随着功能强大的计算机和新算法的出现，贝叶斯方法在 21 世纪的统计中得到了越来越多的应用。
贝叶斯定理
贝叶斯统计围绕着 贝叶斯定理 的用法，以【托马斯·贝叶斯】命名，以及 和 概率的乘积规则 p(y) 在获得更多证据 p(x|y )后，产生一个 后验概率分布 p(y|x )。
方程式 1。
等式 1 中分布 p(x，y) = p(x|y)p(y) 的变量 Y 的所有可能配置的总和导致 X 的 边际分布 ，即p(X)；这种计算叫做边缘化。在连续变量的情况下，等式 1 的分母中的和变成了通常难以处理的积分。此外， p(x) 是归一化常数，其缩放后验分布，使得其在整个空间上的积分再次变为 1——这是 概率分布 的属性之一。(但是，如果省略归一化，等式 1 就变成了 比例 。)**
**原则上，任何参数模型(线性回归、逻辑回归，甚至是**神经网络)都可以通过陈述参数 β 的先验来公式化为相应的贝叶斯版本。使用训练数据( X ， y )和贝叶斯定理，参数的后验概率分布计算如下
方程式 2。
一般来说，通过分析得出后验分布是不可能/不切实际的。但是，通过使用专用算法可以获得数值结果，这将在下面的章节中解释。
蒙特卡罗方法
首先需要介绍几个关键概念，一个是 伪随机数生成 。随机数通过采样 算法产生，例如拒绝采样。这些算法的目的是生成作为一个集合遵循某种特定分布(例如高斯分布)的数。然而，拒绝采样依赖于一个更简单的构建模块，该模块在 0 和 1(U(0，1))之间的区间内生成均匀分布的随机数，例如 线性同余生成器 。拒绝抽样适用于维数很少的非常简单的分布；在高维中，概率质量集中在很小的体积内，导致大量样本被拒绝(效率低下)。对于这样的场景，已经开发了更先进的算法；然而，它们仍然利用更简单的采样算法，例如上面提到的那些。
另一个重要概念是 随机数值积分 利用 蒙特卡罗方法 (或称蒙特卡罗实验)。蒙特卡罗方法是一种广泛的计算算法，它依靠随机数的生成来获得数值结果。例如，通过从期望的分布中生成随机样本，等式 2 的分母中的积分变成了正则期望(平均值的最大似然估计)，其精度与抽取的样本数成比例。
一类被称为 马尔可夫链蒙特卡罗 的采样算法已经被开发出来，用于从复杂的概率分布中抽取样本。一个 马尔可夫链 是一个可能事件(结果)的序列，其中每个事件的概率只取决于前一个事件达到的状态。它是满足 马尔可夫性质 (马尔可夫过程)的随机过程的结果，即过程未来状态的条件概率分布仅取决于最当前的状态，而不取决于之前的事件序列。在马尔可夫链蒙特卡罗算法类中，有几种变体，例如，【Gibbs 采样】 或Metropolis-Hastings；早期的变体之一，Metropolis 算法，将在下一节简要描述，因为它非常容易理解。
贝叶斯计算的目标是从后验分布中获得一组独立的绘图，以合理的精度估计感兴趣的量。更准确地说，这种随机抽取用于总结模型参数的后验分布，例如，通过报告抽样分布的 2.5%、25%、50%、75%和 97.5%的点。然而，绘图也可以用散点图或直方图来可视化。
最后，一旦模拟了模型参数的后验分布，就可以通过忽略条件后验分布中的参数来进行新样本的预测。等式 3 本质上是另一个积分，可以使用数值方法通过随机采样来计算。
方程式 3。
大都会算法
如上所述，从多维(可能不是归一化的)概率分布 p(z) 中抽取样本具有挑战性。采样算法的思想是从更简单的建议分布中反复抽取样本，例如，一维高斯分布的乘积。一种这样的蒙特卡罗方法被称为 Metropolis 算法 ，其生成马尔可夫链并最终收敛到期望的分布 p(z) 。更精确地说，随机数是使用建议分布生成的，新样本是被接受还是被丢弃，取决于它们是比前一个样本更有可能还是更不可能——并且仅仅是前一个样本(因此具有马尔可夫特性)。在实践中，由于计算能力增加，并行模拟几个链以检查 收敛 和 混合 。
图 1: Metropolis 算法。(来源:作者)
具有低方差的建议分布导致小的步长和高度相关的样本；相反，高方差导致较大的步长，但接受率较低，因此许多样本将被拒绝，从而减慢算法。由于这个原因，已经开发了解决这个问题的更高级的算法。关于这个话题的更多理论，可以推荐 Andrew Gelman 的书《贝叶斯数据分析》,这本书有免费的 PDF 文件。
这是本书的 pdf 格式，可以下载用于非商业目的。
www.stat.columbia.edu](http://www.stat.columbia.edu/~gelman/book/)
贝叶斯逻辑回归
让我们看一个使用葡萄酒数据集的具体例子。逻辑回归是一种对变量对二元类标签( y =0， y =1)的影响进行建模的算法，即解决监督分类问题。后验概率分布给出了模型的每个权重或输出的区间估计。在数据分析中，不仅要提供一个好的模型，还要提供结论的不确定性估计，这一点至关重要。
我们收集了一批酒瓶；我们希望为一种特定类型的葡萄酒开发一个分类模型，这种葡萄酒由于价格高而经常被模仿。此外，只有当我们(不)足够确信这是(不是)模仿时，我们才会采取行动。为了表征一款酒，我们用它的化学成分(由 高效液相色谱 测定)相对于 黄酮类脯氨酸 。
图 2: 不同葡萄酒的小提琴图和散点图(深色部分)。(来源:作者)
使用 Python 库 PyMC 通过采样开发和推断后验概率分布 p(y=1|x，β) 的线性模型。从上图可以看出，线性决策边界可能就足够了。此外，由于还没有任何先验知识，我们必须设置 无信息先验 ，这相当于逻辑回归的一个 集合 。请注意，变量在分析前没有标准化，因此参数的大小没有任何意义。
使用葡萄酒数据集的 PyMC 示例。
下面的图中给出了产生的链，显示了收敛和混合，但应评估相关的统计数据。请注意，每条链的样本总量(通常总共四个)是 3000，但是每条链的前三分之一被丢弃。
图 3: 马尔可夫链的推理程序；trace 显示了良好的混合和收敛(至少在视觉上)。(来源:作者)
为了评估采样算法的成功，开发了几个统计数据。下表是推理过程的总结。
最重要的统计数据是 R-hat ，它是一种诊断方法，通过比较多个链之间的方差和每个链内的方差来测试缺乏收敛。如果收敛已经实现，链间和链内的差异应该是相同的。为了最有效地检测不收敛的证据，每个链应该已经被初始化为相对于目标分布分散的起始值。在这种情况下， R-hat 为 1，因此实现了收敛和混合。
另一个统计是有效样本大小，它基本上是由序列内的的量校正的样本总量(样本乘以链数)，即总量除以校正因子。如果这个指标(太)低，可能需要更多的样本。
表 1: 模型和推理总结。如指标所示，链已经收敛。
让我们来看看我们的估计参数。94%贝叶斯可信区间不包括零，因此我们似乎有信心两个变量都与模型相关。
图 4: 估计参数的联合分布。(来源:作者)
让我们看看变量的部分相关性。两者都是随着值的增加而增加概率。然而，由于线的分散性较大，黄酮类化合物的含量有更多的不确定性。
****
图 5: 通过不透明度和带宽可视化预测的不确定性。与脯氨酸相比，黄酮类化合物的分类包含更多的不确定性。(来源:作者)
现在，我们有两瓶新酒。我们有多确定这是(或者不是)我们要找的酒？让我们做一个预测，评估后验分布。
****
图 6: 两个新样本后验概率的带宽预测不确定性。(来源:作者)
第一瓶 94%可信区间[0.90，1.00]，第二瓶[0.50，1.00]。对于第二个，我们不太有信心，事实上，模型在推断过程中从未见过这样的瓶子。
变分贝叶斯
值得一提的是，还有另一种解决整个问题的方法，称为 变分贝叶斯 ，其中后验概率用一些更简单的参数分布来近似，即所谓的变分分布，由此可以直接计算积分。该分布通常被公式化为单个单变量概率分布(平均场近似)的乘积，并且其参数通过最小化真实后验分布和变分分布之间的 Kullback-Leibler 散度 来找到。变分法精度较低，但速度更快，可以作为基于模拟的方法的起点。
摘要
贝叶斯统计是一个强大的框架，它解决了频繁主义方法不能解决的问题，如预测中的不确定性。这种不确定性是由于训练数据不足或测试数据中的异常样本造成的。此外，它还是质量控制的一种形式，通过 PyMC3 和 TensorFlow Probability 等库，这些方法已经变得人人都可以使用。
概率模型不确定性预测的逐步指南。
towardsdatascience.com](/bayesian-neural-networks-with-tensorflow-probability-fbce27d6ef6)******
不相关与独立随机变量——定义、证明和例子
它们不是一回事…
莎伦·麦卡琴在 Unsplash 上的照片
关于技术知识，我通常支持对正在使用的方法有基础的理解。我通常不喜欢记忆任何东西，并尽可能避免这样做。相反，我专注于发展强有力的概念基础，从这些基础中我可以用数学方法推导出我可能需要的任何东西。
在概率论和数理统计领域，利用方法/定理通常依赖于通用的数学假设和约束条件。两个这样的数学概念是随机变量(RVs)是“不相关的”，以及 RVs 是“独立的”。我已经看到了很多关于这些概念的混乱(包括在媒体平台上)。这些是数学上明确规定的术语，它们并不意味着同一件事。对于任何对统计学、数据科学或机器学习感兴趣的人来说，这些概念绝对值得理解。
本文旨在:
- 从数学上规定 RVs 不相关和 RVs 独立的定义
- 证明独立的 RVs 根据定义也是不相关的
- 证明 RVs 可以不相关但不独立(举例)
1.数学定义:
出于演示的目的，让我们假设我们有连续的 RVs X 和 Y ，它们都是在一些实值支持上定义的。
不相关的 RVs:
两个 RVs X 和 Y 是不相关的，如果它们的联合分布的期望值等于它们各自的边际分布的期望值的乘积。从数学上来说，这是:
这相当于说明 X 和 Y 之间的协方差为零:
独立的房车 :
两个 RVs X 和 Y 是独立的，如果它们的联合分布值等于它们各自的边际分布值的乘积，对于沿着它们各自的支撑的 X 和 Y 的任何可能范围。从数学上来说，这是:
这相当于陈述了 X 和 Y (互为条件)的条件分布等价于它们各自的边际分布。换句话说，拥有关于 X 的信息不会提供关于 Y 的额外信息，反之亦然:
2.根据定义，独立的 RVs 是不相关的:
如上所述，根据定义，独立的 RVs 也是不相关的。让我们使用我们的示例 RVs X 和 Y 以及上面的数学定义来证明这一点:
3.RVs 可以是不相关的，但不是独立的:
如上所述，RVs 可以是不相关的，但不是独立的。下面是一个经典的玩具例子来说明:
附录:标准正态分布 RV 的奇矩为零的证明
在第 3 节展示 RVs 不相关的例子时，我利用了标准正态分布的一个特殊属性(奇数阶矩都为零)。对于那些感兴趣的人来说，这个属性的证明如下所示:
最后的想法
希望以上推导和举例有见地。我认为没有足够的人花时间去做这些类型的练习。上述概念的知识对于理解高斯-马尔可夫定理、恢复有效方差估计量的条件以及该领域的更多概念至关重要。
希望以上有见地。正如我在以前的一些文章中提到的，我认为没有足够的人花时间去做这些类型的练习。对我来说，这种基于理论的洞察力让我在实践中更容易使用方法。我个人的目标是鼓励该领域的其他人采取类似的方法。我打算在未来写一些基础作品，所以请随时在 LinkedIn 上与我联系，并在 上关注我在 上的更新！
用统计模型揭示政府偏见
照片由 Unsplash 上的 Aditya Joshi 拍摄
对澳大利亚“体育丑闻”的数据驱动分析
你可以跑，但你躲不了(从统计数据上)。澳大利亚自由党很快就会发现这一点。
最近几周，自由党被指控在 2019 年大选前使用 1 亿美元的体育拨款来赢得选票。作为一名数据科学家，我决定通过分析数据来验证这一指控。因此，我着手研究我们能从公开的数据中了解到什么，这些数据与现在众所周知的“体育流氓”丑闻有关。
我的想法是，如果自由党根据选民的投票行为来分配补助金，那么我们应该能够使用统计模型来揭示他们的策略。
在我继续讲下去之前，这里有一些背景资料，以防你还没有看新闻…
背景
澳大利亚自由党自 2013 年联邦选举以来一直执政。在去年选举之前，政府宣布了一项社区体育基础设施计划，为当地体育组织提供高达 50 万美元的拨款。作为促进全国体育参与的一项举措，各组织必须提交一份申请，由澳大利亚体育局进行独立评估和评分。体育部长应该用这些分数来决定批准哪些拨款。
今年，澳大利亚国家审计署发布了一份审计报告，该报告暗示，但没有断言，自由党优先考虑“边缘”和“目标”选民的申请，而不是遵循澳大利亚体育局的建议。本质上，自由党被指控为目标地区提供更多拨款，以赢得 2019 年选举的选票。
澳大利亚选区地图内 684 个体育资助的位置。
回到行动
我分析的最终目的是看是否有证据表明拨款是根据上次选举中选民的投票行为分配的。本质上，有没有一个潜在的策略？如果是的话，自由党针对的是什么样的投票模式？
为了进行这一分析，我需要结合三个关键的数据源，所有这些数据源都是公开的:
- 关于批准的每项资助的信息(地点和金额)
- 上一次选举(2016 年)中每个选民的选举结果，以及
- 控制其他可能起作用的社会人口因素的人口普查数据。
有了这些信息，我将能够生成一个数据集，其中包含分配给每个选民的拨款总额，以及上一次选举的相应投票结果和人口普查的其他控制。
(注意所有的数据收集、争论和分析都是使用 *R* )
获取数据
澳大利亚体育局在他们的网站上发布了一份 684 项批准拨款的列表，所以我用rvest抓取了这份列表。接下来，这个列表被提供给一个 Google Maps API，以便获得每笔赠款的经度坐标(使用来自 dkahle 的 GitHub 的 T2)。这适用于几乎所有的资助。一些快速的谷歌搜索有助于找到这些组织的其余部分的位置。
为了将每个体育组织分配给其相应的选民，我需要 2016 年联邦选举的选民地图。方便的是，这可以在eechidna包中获得，还有选举投票数据。使用来自sp和rgeos的一些漂亮的函数来遍历选举人，每笔拨款都被成功分配。
现在，补助金已经汇总，因此得到的数据集包含 150 个联邦选区中每个选区的补助金总额。然后，这将与 2016 年选举中的两党优先投票和摇摆投票以及 2016 年人口普查中的社会人口变量集合结合起来，这些变量也是从eechidna中获得的。两党首选投票是倾向于自由党而不是工党的投票百分比，摇摆投票是 2013 年至 2016 年两党首选投票的百分点变化。
这是结果数据集的快照。
LiberalVote 为 2016 年两党优先投票百分比， TotalGrantAmount 为给予选民的拨款总额， Swing 为 2013 年至 2016 年的百分点摇摆投票。其他变量是 2016 年人口普查的社会人口特征。
现在是分析的时候了
我真正想揭示的是选民在 2016 年大选中对自由党的两党偏好投票与获得的拨款金额之间的关系。首先，让我们看看拟合了线性回归线的散点图是什么样子(使用ggplot2)。
赠款总额(TotalGrantAmount)和两党优先投票(LiberalVote)散点图，有回归线和 95%置信区间。观察到正相关。
这里有一些证据表明，更支持自由党的选民更有可能获得体育拨款。但实际上，如果自由党真的实施了有针对性的策略，这两个变量之间的关系可能不是线性的。向最安全的地区提供最多资金的策略不会是一个非常好的策略——因为在 2019 年的选举中，安全的自由派选民不太可能容易易手。
因此，我没有拟合线性回归线，而是绘制了一条黄土曲线——一条平滑的非参数回归曲线，为局部数据点赋予更高的权重。这有助于更好地了解 2016 年投票和拨款金额之间的关系。
拟合黄土曲线而不是线性回归线。在两党的优先投票中，拨款明显飙升，略高于 50%。
哇！黄土曲线揭示了一种重要的关系，这种关系在最初的散点图中并不明显。
在 2016 年，自由党赢得的选民获得的拨款金额明显增加，但没有以巨大优势获胜。这些人被称为“边缘”自由党选民。这一阴谋本身就是资助计划中自由主义偏见的令人信服的证据。
为了更进一步，我决定拟合一个广义加性模型(GAM ),以拨款额作为响应，两党偏好投票(从 2016 年开始)、摇摆投票(从 2013 年到 2016 年两党投票的百分点变化)和人口普查特征作为协变量。每个协变量都作为平滑项包含在模型中。在这种情况下，我使用了惩罚回归样条——一种类似黄土的灵活的非参数回归。
GAM 与mcgv拟合，造型效果的剧情由visreg产生。
两党优先投票如何影响每个选民的总拨款的拟合效果。高于 50%的峰值显然是显著的。
可以清楚地看到，即使在考虑了摇摆投票和人口普查变量后，边缘自由党选民的拨款金额仍有巨大的飙升！
在其他条件相同的情况下，一个获得 54%自由党支持的选民和一个获得 54%工党支持的选民所获得的资金差异是惊人的 50 万美元。
GAM 回归输出证实了这是我着手发现的决定性证据。在 p 值< 0.001 的情况下，数据证明两党优先投票对分配给每个选民的总拨款具有统计上的显著影响。
为了检验我的发现对异常值的稳健性，我在去除了两个拥有巨额赠款(180 万美元和 200 万美元——给布斯比和道森的选民)的选民后，重新估计了模型，发现观察到了同样的自由主义偏差。此外，模型假设似乎是有效的，因为残差是近似正态的(使用分位数-分位数图观察)。对于那些想知道的人来说，完全指定的 GAM 解释了 46.3%的偏差——在这种情况下这是一个很好的解释。
最后
我的统计模型提供了经验证据，证明自由党利用体育拨款试图巩固对 2019 年大选的支持。数据清楚地表明，边缘自由选民是政府的目标，因此获得了更多的拨款。
虽然这一分析从政治角度来看非常有趣，但它也展示了如何利用公开数据产生强有力的见解。
对于任何对探索澳大利亚选举和人口普查数据感兴趣的人，我强烈推荐查看eechidna R包(我是该包的作者和维护者)。该软件包可以轻松访问 2001 年至 2019 年的澳大利亚选举和人口普查数据(包括地图)，可在 CRAN 上获得。
所有用于分析的代码都可以在这个 GitHub repo 中找到。如果你有任何想法，想要联系，这是我的 LinkedIn 。
滚动跨期分析揭示动量效应
作者图片
扩展滚动窗口预测以发现经验策略的洞察力
这篇文章展示了跨期方法，该方法扩展和概括了滚动时间序列技术的范围，用于导出过渡过程和经验策略的模型。该方法是在解释动量溢价的背景下说明的，动量溢价是一个长期的持续挑战。
Jegadeesh 和 Titman [1]在 1993 年对动量效应进行了记录，表明动量效应从 1965 年到 1989 年为美国普通股带来了异常的正回报。自那以后，传统的横截面动量策略被概括为根据过去一年的回报滞后一个月对资产进行排序，并做多过去赢家资产的子集，同时做空过去输家资产的子集。多篇论文指出，在广泛的资产类别、市场和时间范围内存在动量，这种观点认为动量因素几乎无所不在。然而，动量溢价的来源仍然是一个悬而未决的问题，没有一个单一的模型来主导叙事。为了解释这一现象，人们提出了许多基于风险的行为模型，其中 Clifford Asness 和 Eugene Fama 对这一主题进行了一次值得注意的讨论[2]。
最近，动量效应被揭示动量论文[3]中提出的滚动跨期分析所解决。在这种方法中，选定的十分位数在时间范围内向前滚动，同时在组合的 11 个月排名、1 个月滞后和相对于组合形成日期的持有间隔期间收集每个月的组合回报。类似于滚动窗口预测，排名周期可以与样本内区间相关联，而持有周期构成样本外区间。下面是揭示动量文件的图表，展示了对 2006 年 8 月至 2017 年 8 月牛市状态的顶部动量十分位数运行滚动跨期方法的结果。
图:样本内和样本外动量最高的 10%月度回报的箱线图。
箱线图突出显示了样本内和样本外间隔之间的凸起。先前的研究主要致力于分析保持期。这些图把我们的注意力转移到了转变过程上。这个凸起的左侧可以用随机抽样模型来解释。根据这一模型，动量效应(样本外)应该为零，因此定义了评估基本动量理论和模型的标准。本文扩展了跨期方法，将动量效应解释为表现强劲的股票组合。该实现基于 Quantopian 平台，包括三个步骤:运行 Quantopian pipeline、选择特性和实施跨期方法。
量子管道
Quantopian 提供了一个数据科学平台，用于在交互式 Python 研究环境中开发和测试量化金融投资策略。这个框架的主要组件之一是一个充当不同异构数据集的公共接口的因子，因此可以在复合管道中统一处理这些数据集。特色应用需要技术和基本特征。
如下面的代码片段所示，技术因素，如动量和已实现的波动性，可以通过 CustomFactor 类使用股票价格来实现。
import numpy as np
from quantopian.pipeline.factors import CustomFactor, DailyReturns
from quantopian.pipeline.data.builtin import USEquityPricingclass Momentum(CustomFactor):
inputs = [USEquityPricing.close]
window_length = 252
def compute(self, today, assets, out, prices):
out[:] = (prices[-21] - prices[-252])/prices[-252]
class RealizedVolatility(CustomFactor):
inputs = [DailyReturns()]
window_length = 126
def compute(self, today, assets, out, returns):
ret_r2 = returns**2
out[:] = np.sqrt(np.nansum(ret_r2, axis=0))
基本面因素是通过晨星和 FactSet 数据集的 Quantopian 内置包装器定义的。
from quantopian.pipeline.factors import MarketCap, Latest
import quantopian.pipeline.data.factset as fset
import quantopian.pipeline.data.morningstar as mstarcap = MarketCap()
bp = 1/Latest([mstar.Fundamentals.pb_ratio]) # book-to-price ratio
ebit = fset.Fundamentals.ebit_oper_af.latest
shldrs_eq = fset.Fundamentals.shldrs_eq.latest
op = ebit/shldrs_eq # operating profitability
ag = fset.Fundamentals.assets_gr_af.latest # asset growth
此外，Quantopian 平台提供了一种机制，用于根据一组约束条件(如市值、美元交易量等)定义不同的股票领域。具体来说，QTradableStocksUS 宇宙包含大约 1600-2100 只股票，其市值大于 5 亿美元，在过去 200 天内的日均美元交易量中值大于或等于 250 万美元，收盘价高于每股 5 美元，以及过去 200 天的可用数据。然后，宇宙用于选择股票，以计算管道内的系数。
from quantopian.pipeline.filters import QTradableStocksUS
universe = QTradableStocksUS()
一旦定义了因素和领域，所有因素都可以合并到一个管道中，并在 Quantopian 引擎中运行。对于这种应用，这些因素被额外扩展用于生成四分位数等级。
from quantopian.pipeline import Pipeline
from quantopian.research import run_pipelinedef make_pipeline():
pipe = Pipeline()
pipe.add(Momentum(mask=universe), "mom")
pipe.add(RealizedVolatility(mask=universe).quartiles(), "rv")
pipe.add(cap.quartiles(mask=universe), "cap")
pipe.add(bp.quartiles(mask=universe), "bp")
pipe.add(op.quartiles(mask=universe), "op")
pipe.add(ag.quartiles(mask=universe), "ag")
pipe.set_screen(universe)
return pipestart, end = '2010-01-01', '2019-06-30'
pipe_data = run_pipeline(make_pipeline(), start, end)pipe_data.replace(-1, np.nan, inplace=True)
pipe_data.dropna(inplace=True)
最后，使用 Alphalens 模块将远期收益和因子合并到一个多指数数据框架中，该框架具有两个指数级别:日期和资产。
import pandas as pd
from alphalens.utils \
import get_clean_factor_and_forward_returns as add_returnsassets = pipe_data.index.levels[1]
prices = get_pricing(assets, start, end, fields='open_price')ad = add_returns(pipe_data['mom'], prices, quantiles=10)
ad.rename(columns={'factor_quantile': 'mom_decile'}, inplace=True)factor_names = ['rv', 'cap', 'bp', 'op', 'ag']
asset_data = pd.merge(ad, pipe_data[factor_names], \
left_index=True, right_index=True)
筛选功能
本文假设动量效应是表现强劲的公司的副产品，而不是投资者的行为。因此，作为下一步，我们需要筛选不同的公司特征(由晨星和 FactSet 数据集提供)，以捕捉动量十分位数内的分布异常。由于可用特征的数量很大，本文将范围缩小到 Fama-French 五因素模型[4]中使用的特征:规模、市净率、营业利润率和总资产增长率。根据之前的动量研究，这个基本面列表进一步增加了一个技术特征，即已实现波动性。
asset_flat = asset_data.reset_index()df_ranks = pd.DataFrame(index=range(0,4))
df_ranks.index.name = 'rank'm10 = asset_flat.query('mom_decile == 10')
m10_count = len(m10)
for n in ['rv', 'cap', 'bp', 'op', 'ag']:
q_counts = \
[len(m10.query('{}=={}'.format(n, q))) for q in range(0,4)]
df_ranks[n] = [q_count/m10_count for q_count in q_counts]
df_ranks
表:动量上十分位数中每个因素四分位数的分数。
上表显示了动量前十分位数中每个因素的概率分布。基于输出，筛选发现了两个主要异常:已实现波动率(rv)和市净率(bp)。对于全面的研究，特征重要性分析可以通过机器学习方法(例如，分类树的平均减少杂质)进行扩展，López de Prado 在“资产管理器的机器学习”[5]中对此进行了概述。交互影响可以通过“HHP 多维模型的交互分析”[6]中介绍的技术进一步增加。
滚动跨期分析方法
根据上述特征筛选分析，动量上十分位数的 83%和 76%分别由已实现波动率的上半部分和市净率四分位数的下半部分组成。然后，这些结果促使我们对从动量前十分之一和整个股票宇宙中选择的相应样本所表现出的跨期行为进行比较。
下面的代码片段实现了每 21 天运行一次的函数，并在投资组合形成日之前和之后收集样本内和样本外窗口(wdw1 和 wdw2)的每日等权重投资组合回报。
def rolling_inout_rtns(df_input, # factors and returns
wdw1, wdw2, # in- and out-of-sample windows
query): # query of selected factor deciles
days = []; dfs = []
ds = df_input.index.levels[0].unique() for t in range(wdw1*21, len(ds) - wdw2*21, 21):
df_t = []
# select assets in the formation date t
xs_t = df_input.xs(ds[t])
assets = xs_t.query(query).index
df_t.append(len(assets))
# calculate means of selected assets for [-wdw1, wdw2]
for dt in range(-wdw1*21, wdw2*21+1):
xs_dt = df_input.xs(ds[t+dt]).loc[assets.values]
df_t.append(xs_dt['1D'].mean())
days.append(ds[t])
dfs.append(df_t)
columns = []
columns.append('count')
for dt in range(-wdw1*21, wdw2*21+1):
columns.append(str(dt))
df_output = pd.DataFrame(dfs, index=days, columns=columns)
df_output.index.rename('date', inplace=True)
return df_output
在定义了样本内和样本外窗口之后，该过程可以随后应用于具有高波动性账面价格比股票的样本。
wdw1 = 3 # in-sample months
wdw2 = 1 # out-of-sample monthsquery = '2 <= rv <= 3 and 0 <= bp <= 1'
asset_rtns = rolling_inout_rtns(asset_data, wdw1, wdw2, query)query += ' and mom_decile == 10'
mom10_rtns = rolling_inout_rtns(asset_data, wdw1, wdw2, query)
最后，绘制结果。
import datetime
from pandas.tseries.offsets import BDay
import matplotlib.pyplot as pltmkt_start = datetime.datetime.strptime(start, '%Y-%m-%d') - BDay(1)
mkt_prices = get_pricing('SPY', mkt_start, end, fields='open_price')
mkt_rtns = mkt_prices.pct_change(periods=1)[1:]t0 = -wdw1*21
t1 = wdw2*21 + 1mkt_means = [mkt_rtns.mean() for dt in range(t0, t1)]mom10_means = [mom10_rtns[str(dt)].mean() for dt in range(t0, t1)]
asset_means = [asset_rtns[str(dt)].mean() for dt in range(t0, t1)]mom10_sma = pd.Series(mom10_means).rolling(window=21).mean()
asset_sma = pd.Series(asset_means).rolling(window=21).mean()x = range(t0, t1)plt.figure(figsize=(9,4))plt.plot(x, mkt_means, label='mkt', color='y')asset_label = 'rv[2,3] & bp[0,1]'
plt.plot(x, asset_means, label=asset_label, color='grey')
plt.plot(x, asset_sma, label=asset_label+' SMA', color='green')mom10_label = 'mom 10 & ' + asset_label
plt.plot(x, mom10_means, label=mom10_label, color='grey')
plt.plot(x, mom10_sma, label=mom10_label+' SMA', color='red')plt.xlim(-40,20)
plt.ylim(-0.003,0.006)
plt.xlabel('Days')
plt.ylabel('Average Daily Returns')
plt.legend(loc='upper right')
图:动量赢家和高波动性成长股的样本内和样本外结果。
如图所示，样本内和样本外区间由三个时期组成:排名、一个月跳票、持有。动量赢家(红色)和高波动性增长四分位数(绿色)之间的排名周期差异是由动量过程采样的高特质回报决定的。在跳跃区间，当排名周期在图上的第 21 天结束时，动量赢家接近高波动成长股的平均水平。在形成日期之后，它们都经历了类似的短期反转。
有趣的是，虽然滞后一个月之前被记录为经历短期反转，但图中的结果显示，跳过最近一个月并因此将其从常规排名区间中排除不再是必要的。然而，在这项研究中，将月份排除在排名之外，有助于通过捕捉动量赢家和高波动性成长股之间的一致性，深入了解动量效应的主要驱动因素。
76%的动量赢家是成长股，因此价值股的贡献并不大。尽管如此，同样的过程也可以应用于相应的高波动性值子集，并解释这种影响。
摘要
本文提供了基于滚动跨期方法的应用程序的分步演示，以解释动量策略的优势。该方法通过连接用于导出过渡过程模型的样本内和样本外区间，扩展了滚动窗口预测技术。具体来说，对于 2010 年至 2019 年，文章揭示了动量策略作为对高波动性成长股进行采样的程序，并通过排名、跳跃和持有区间来跟踪它们的转变。这个应用补充了两篇基于量子理论的关于动量的论文，包括波动时机和揭示动量。如果感兴趣，可以在 Quantopian 论坛的以下主题中找到更多信息和笔记本:
参考
[1] N. Jegadeesh 和 S. Titman，回报买入赢家和卖出输家:对股票市场效率的影响 (1993)，金融杂志
[2] C .阿斯内斯，法玛论动量 (2016)，AQR 见解
[3] Y. Malitskaia，揭开势头 (2019)，SSRN 工作文件
[4] E. Fama 和 K. French，一个五因素资产定价模型 (2015)，金融经济学杂志
[5] M.M. López de Prado，资产管理者的机器学习 (2020)，剑桥量化金融要素
[6] Y. Malitskaia、W.H. Greene 和 Y. Malitsky，HHP 多维模型的相互作用分析 (2012)，ERCIM 2012
网站
Quantopian 是一个免费的在线平台和社区，用于教育和投资算法的创建。量子乌托邦…
www.quantopian.com](https://www.quantopian.com) [## 晨星|让投资者成功
我们的独立研究、评级和工具正在帮助整个投资生态系统的人们撰写他们自己的…
www.morningstar.com](https://www.morningstar.com) [## FactSet |面向投资专业人士的集成金融数据和软件解决方案
了解 FactSet 与众不同的最佳方式是自己去看。我们与投资部门合作过…
www.factset.com](https://www.factset.com)
揭开宇宙:寻找系外行星候选者和其他异常的机器学习方法
演示如何使用机器学习快速准确地搜索大量天文数据
图片来源:spaceengine.org
灯塔的目的是在黑暗中给别人指路。当我们观察浩瀚的夜之海洋而没有人点燃它时，会有一些令人不安的暗示。—伊萨克·阿瑟
简介
T 天文学领域正日益成为一个大数据问题，随着更大的调查和更精确的仪器，数据进入了万亿字节的范围并呈上升趋势。在更广阔的物理领域更是如此。大型强子对撞机每小时获得的数据比脸书一整年获得的数据还多。因此，反馈可能会很慢，在发现之前，发现可能会在数据集中存在多年。
我概述了一种基于机器学习的方法，它允许我们:
- 天文时间序列数据的快速目标搜索
- 自动检测最异常的测量值
- 这样做有合理的计算要求
然后，我将这种方法应用于从最近的卫星观测中获得的真实天文数据。
光曲线
我们目前的分辨率还不足以从视觉上分辨邻近太阳系的细节。然而，我们对物理学有很好的理解，可以通过间接测量来理解正在发生的事情。其中之一是恒星的亮度。
如果一种现象(自然的或其他的)正在恒星周围发生，我们可以预期它的间接影响会随着时间的推移出现在恒星的亮度中。如果一个物体挡住了恒星的去路，那么恒星就会变暗，变暗的程度与该物体相对于恒星的表观尺寸成正比。根据亮度曲线的轮廓，我们可以推断出物体的属性，例如它的大小和它离恒星有多近。这是用来寻找系外行星的主要方法，被称为凌日法。
中转方法。图片来源:NASA
在分析这类数据时，有几件事要记住。尤其是寻找像行星这样的自然现象。
- 你需要观察恒星足够长的时间才能看到变暗的重复。例如，如果一个外星文明正在使用这种方法观察太阳，他们将不得不观察它一年，以观察地球重复造成的变暗。理想情况下，你会希望它发生三次。
- 你需要恒星轨道平面向我们倾斜。人们倾向于认为所有的行星都以与银道面相同的角度运行，但事实并非如此。例如，我们自己的太阳系与银道面呈 60 度倾斜，这就是为什么我们在南半球能更好地看到银河系中心的原因。不幸的是，如果一颗恒星的两极直接指向我们，我们将不会因为它的行星的运动而观察到亮度的任何变化。
- 1)和 2)的组合以及多年来光度计灵敏度的变化导致了测量偏差。如果你只是在短时间内观察恒星，并让现象在这段时间内重复出现，那么你多半是在识别恒星周围的短年天体。现在，我们可以敏感地注意到地球大小的物体造成的亮度变化，但很长一段时间以来，我们所知道的唯一行星是木星大小的物体，其轨道比水星更近(因为内轨道的年数更短)。
最终目标是分析来自 transit 方法测量的时间序列信息，以找到感兴趣的对象。
挑战
尽管受限于与太阳系处于同一平面的物体，银河系中恒星的数量是如此之多，以至于仅从一次观测中仍会产生数十万或数百万条光变曲线。每个测量还会有数千个数据点和附加信息。仪表板是为公民科学家开发的，以帮助天文学家处理数据。现在著名的 KIC 8462852(虎斑星)就是这样被发现的。这种方法确实意味着数据可以存储在后端，将计算需求缩小到一个地方。但是仍然存在以下问题:
- 这仍然是一种蛮力方法
- 寻找特定的现象仍然是耗时的，并且可能仍然需要不能通过前端完成的编程方法
- 如果影响微弱，现象可能会被忽略
- 很难快速获得对数据的整体理解
给我们快速反馈的方法非常有用，尤其是当我们想用其他仪器进行后续测量时。如果调查中的数据量呈上升趋势，上述问题的解决方案需要廉价、快速且可扩展。
这里可以引入几个有用的机器学习概念。
相关机器学习概念
大多数非该领域的人都理解机器学习模型的概念，即进行分类预测，无论是二元(是或否)还是多分类。不幸的是，在许多实际情况下，这些标签可能不可用或不可靠。即使您有时间和资源，数据的正确标记可能也没有那么有用，并且不是每个对数据感兴趣的从业者都会对相同的标记感兴趣。我们想要的是利用深度学习架构的潜力来理解复杂的结构，并以一种无人监管的方式返回见解。
时间序列编码
请考虑以下情况。如果我问一个人两幅画是否相似，通常不会(只是)有一些不同的特征使它们如此。我们的语言有时不足以描述为什么它们是相似的，但你知道它们是一种手动波浪形的方式，这与风格等有关。这是因为我们对我们所看到的东西有一种复杂的更高层次的理解(或感觉),而不是我们必须写下来或交流的。
深度学习模型通过理解数字阵列或向量表示形式的输入来工作。即使我们做一个分类问题，我们所做的只是将网络的最后一层从这个向量表示转换成一个显式输出。由于网络中的所有非线性，这些表示通常不是人类可以解释的。
对于我们的光曲线，我们不想做预测，而是让模型将光曲线变成这些矢量表示中的一个，捕捉它的所有结构。当我们需要洞察力或进行搜索时，我们最终将只是与这种表征一起工作。
既然我们有了一种以相同方式表示所有光线曲线复杂性的方法，我们可以执行聚类并评估相似度。我们还可以决定我们希望这个表示有多大，从而使它有利于压缩。
类似
如果我们将每个向量想象成空间中的一个点，其中每个数字对应一个坐标*，如纬度或经度*，我们可以通过其位置以及相对于其他点的距离和方向来提取信息。这本质上是相似的基础。
余弦相似度
通常，两个数据点相对于彼此的方向比它们的欧几里德(L2)距离更能说明问题。如果我们有数据的矢量表示，两个矢量的方向就是它们的归一化(按大小)点积。
余弦相似度。图片来源:维基百科
由于操作简单，即使在高级语言中计算也非常快，因为它们经常在 C 中运行算术运算。
对于纯软件思维的人来说，有一种诱惑，即逻辑上认为是一组 For 循环的总和。但处理器实际上在非常低的级别上有特殊的算术运算，称为融合乘加，它将运算结合在一起，比单独做乘法和加法更快。因此，对于相似性搜索，如果我们想最大化速度，我们希望利用这种算法。这样做比用高级语言手工实现要快几个数量级。这是我们想要将数据转换成潜在向量表示的另一个原因。理解这一点是用合理的硬件进行数十亿规模的相似性搜索的关键。
然而，这一指标有一些警告需要讨论，因为时间序列中相似性的概念不像其他数据类型那样简单明了。有相似性度量可以证明更适合时间序列相似性，包括:
- 复杂性不变相似性
- 基于相关性的相似性
- 基于自相关的相似度
在以后的工作中记住这一点，我稍后证明余弦相似性对于我们试图解决的问题仍然是足够的。
建筑与培训
为了获得一个给我们矢量表示的模型，我们需要训练一个神经网络。即使我们是在无人监督的情况下做这件事，我们仍然需要给模型一个任务来解决。
一个自动编码器是一种由两部分组成的架构——编码器和解码器。编码器将数据还原为我们选择的矢量表示，解码器尝试从该表示中重建信号。因为它们被一起训练，所以编码器将试图找到最有意义的表示，使得解码器可以创建最准确的重建。
由于我们处理的是时间序列，因此有两种类型的输入重构可用于训练:
插值 —给定一个值在区间【t₀，tₙ】内的时间序列，我们以概率为 10–50%的点子集为条件，然后尝试在相同范围内重建完整信号。
外推 —给定一个值在区间[t₀，t₂ₙ]we]的时间序列，编码信号的前半部分[t₀，tₙ]并尝试重建信号的后半部分[tₙ，t₂ₙ].]像插值一样，我们可以以第一个区间上的点的子集为条件。
为了获得光曲线的编码，我选择了插值方法，因为我们的目标是比较信号，而不一定是预测。
ODE-RNNs
递归神经网络(RNNs)是时间序列的主要结构类型，因为它们具有捕捉长期相关性的能力。但如果时间序列是不规则采样，一般会采用各种变通方法。直到最近，处理该问题的一种常用方法是将时间步长信息作为输入要素与信号一起传递，在缺失点处估算值，或者以大于最大间隙的分辨率进行固定时间离散化并聚合值。
神经常微分方程领域的最新发展使我们能够通过将 RNNs 转化为连续时间模型来更好地拟合缺失时间点。这在论文中探讨了不规则采样时间序列的潜在常微分方程，并被称为常微分方程-RNN。
奥德-RNN 插值缺失数据。图片来源:潜颂
处理不规则采样时间序列的能力对于我们的编码器来说是非常有用的能力，因为我们的向量表示将更能指示光曲线的真实动态。如果有几个间隔时间的后续测量，这是特别有用的。我们可以用它作为我们的网络编码器。
凌日系外行星调查卫星
TESS 于 2018 年发射，是首个在整个天空中搜寻系外行星的太空任务。它通过将天空划分为每个半球的 13 个观测扇区来实现这一点。每个扇区观察 27 天，在重新定向以对下一个扇区重复该过程之前，每 2 分钟记录一次测量。重点是附近的恒星，每个扇区观测大约 200，000 颗恒星。在编写本报告时，主要任务的 26 个预定部门中的 17 个已经完成，并向公众提供了数据。每个扇区的时间序列数据大约为 40GB，因此我将分析重点放在最近的两个扇区(16 和 17)。
TESS 相机的视野和半球分区。图片来源:TESS 天文台指南。
方法
TESS 数据包含对象 id 和各种带有相关时间戳的测量值。我将时间戳和流量(亮度)按对象 id 分组。然后实施以下预处理步骤:
- 离散化——尽管每两分钟测量一次，但飞船上的时钟分辨率很高。为了批量求解常微分方程，我们必须将所有时间步长的并集输入到模型中。如果时钟分辨率很高，这会产生许多独特的时间点，即使实际上它们非常接近。因为我们希望减少所有时间步长的并集大小，所以我将时间步长四舍五入到两位小数。从检查来看，这对于光变曲线没有明显的差异，并且聚集没有显著改变每个样本的观察数量。
- 归一化方案 —乍一看，通过每次测量最大值进行归一化似乎是合适的，因为这是人们通常描述亮度下降的方式，例如 0.99 相当于亮度下降 1%。然而，该模型很难同时适用于亮和暗的物体，因为亮和暗物体的量值变化可能要低几个数量级。为了使样本比例不变，I 使用每个样本的平均值和标准偏差对每个样本进行 z 归一化(非全局)。本文后面的所有结果都显示了标准化的光线曲线。
我训练了一个多层 ODE-RNN，但是对潜在 ODE 库做了以下修改，因为我发现它们为我的任务提供了更好的结果:
- 损失函数 —我把标准的 L2 损失换成了 L1 损失。我发现它总体上更好地捕捉了局部动态，而不是试图适应大峰和大谷。
- 微分方程解算器 —默认使用欧拉法。尽管这种方法训练速度更快，但随着步数的增加，误差会累积。我把它改成了自适应龙格-库塔解算器，这是 torchdiffeq 包中的标准配置。
- 禁用偏置 —对于所有层，我禁用了偏置参数，因为我发现它给出了更好的结果，并且训练更快。
使用 10 的批量大小对来自总数据集的 5000 个随机子集进行训练。Adamax 优化器以 1e-3 的固定学习速率使用。选择的潜在维度为 40。在具有 12GB GDDR5 RAM 的 NVIDIA GTX 1070 上进行训练。
推论
训练之后，目标是将我们所有的光线曲线转换成我们的潜在表现形式。
虽然一起训练编码器和解码器，但我们在推断时只使用编码器计算每条光曲线的潜在表示。这只需要做一次，并且在评估结果时比训练时快得多，因为我们不需要更新任何梯度，也不需要将所有梯度保存在内存中。尽管如此，评估的速度将是编码器有多少参数的函数。
在将结果保存到磁盘之前，我们只需一次将一批数据加载到内存中。在我的情况下，我能够在笔记本电脑 CPU 上每秒评估大约 2500 个潜在表示。完成后，我们就可以进行主要的分析了。
可视化结果
尽管将每条光变曲线简化为 40 个数字比成千上万个数据点需要处理的信息要少得多，但乍看之下，给定值的含义肯定不明显。能够一次看到所有的数据给了我们一个大概的概念，我们可以期望看到不同类型的测量。
t-分布式随机邻居嵌入(t-SNE)是一种有用的算法，可以将数据的维数降低到我们选择的大小。这个算法是而不仅仅是高维空间到低维空间的投影。简而言之，通过估计一个点选择另一个点作为其邻居的条件概率来赋值。通过这样做，你基本上保存了尽可能多的来自更高维度的信息。
选择二维输出意味着我们现在可以将完整的光变曲线测量视为 x-y 平面上的单点。彼此靠近的点是彼此相似的测量值。出现的是可以研究的集群。
苔丝 16 和 17 区的每颗恒星编码的 t-SNE 可视化。每个点都是一条光曲线。灯光曲线是 z 规格化的。
从我们的 t-SNE 图中我们可以看到，有清晰的星团，你可以自己在它们周围画一个圈。有些似乎包含常见的测量值，而另一些则比较罕见。在小岛上磨砺显然有利于快速识别更独特的观测类型。然而，从计算的角度来看，我们可能希望完全跳过可视化，直接从潜在表示进入集群标签。
使聚集
尽管 t-SNE 在将测量值组合在一起方面做得很好(以及提供可视化的明显好处)，但我们有时确实必须小心它的解释。减少到二维可能仍然隐藏测量的真实变化。通常，直接在 40 维潜在空间向量上聚类是有意义的。
有几种不同类型的聚类算法，每种算法处理问题的方式都非常不同。不幸的是，大多数类型的聚类算法对相似性度量、聚类大小或数据集的几何形状做出假设。令人烦恼的是，只有少数算法不需要指定集群的总数。
如果我们不得不做出任何妥协，那么这几乎就违背了我们一开始就试图做这件事的目的。这是因为我们不一定知道我们期望看到什么。由于数据集很大，我们希望尽快获得有用的反馈。
考虑到这些挑战，我们对天文观测的有用聚类算法的标准应该是能够:
- 自动检测簇的总数
- 创建具有不同大小和分布的集群
- 识别无关的数据点，并且不强制它们进入聚类
- 提供一致的结果，对超参数调节不敏感
尽管有几十种聚类算法，但基于层次密度的带噪声应用空间聚类 (HDBSCAN)算法几乎没有做出妥协，并且达到了上述所有标准。
在 TESS 数据集上应用该算法返回了 300 个聚类。有些似乎与微妙的差异密切相关，但总的来说是非常明显的。每次的计数从几个到几百个不等。下面是来自三个不同集群的三条光变曲线的示例。
不同聚类中的光变曲线样本。每行是一个簇。
如果我们想对一个候选名单进行后续测量，从聚类水平开始可以使这一决定更容易，因为我们现在可以量化一种类型的观察因其聚类的大小而有多常见。
使用快速人工智能相似性搜索(Faiss)进行快速搜索
如果我们利用前面讨论的相似性度量，可以实现快速搜索相似亮度曲线的能力。Faiss 是一个有用的库。
Faiss 是一个为快速高效的相似性搜索而优化的库。它是用 C++和 python 包装器编写的，因此使它易于访问和使用。它假设相似性可以通过 L2 距离或点积来评估，并且具有 GPU 支持。
Faiss 提供多种索引选择，以满足内存消耗、速度和精度要求。在这种情况下，我可以将所有向量放入内存中，并希望得到精确的结果，因此我选择了一个称为“扁平内积”的索引。平面索引不会简化计算，也不会对数据库的结构做任何假设。后者不是必需的，因为搜索是彻底的。
内存消耗是每个向量 4 x d 字节，其中 d 是我们潜在表示的维度(在我的例子中是 40)。这种方法的一个附带好处是，通过选择正确的 d 并进行适当的权衡，可以设计一个在任何硬件上工作的完整流水线。Faiss 文档中有一个关于索引选择的有用指南这里是。但是作为补充说明，如果我们要在此分析中包括整个 TESS 数据集，它可能不适合内存，并且可能需要(但不一定)使用 GPU 来提高速度，并选择一个索引来优化内存。
Faiss 只会做两个向量的点积。为了获得余弦相似性，我们需要在将样本添加到数据库时对每个样本进行 L2 归一化。然后，我们也对输入向量执行此操作，这将导致两个 L2 归一化向量的点积，从而给出它们的余弦相似性。
# Obtain a vector representation of our curve and normalize
input_vec = input_vec/np.linalg.norm(input_vec) # L2 normalization
input_vec = input_vec.reshape(1,-1)# Search our faiss db with the vectors from all TESS observations
probs, ixs = db.search(input_vec, k=5) # return top 5 results
因为目标是快速获得结果，所以注意这一过程的时间复杂性是很重要的。在整个数据库中找到前 5 个最相似的光变曲线用了 不到 1 毫秒 。此外，无论我们想要多少结果，这都可以很好地扩展。切换到前 1000 个结果不会对查询时间产生很大影响。
%timeit db.search(input_vec, k=5)
**762 µs** ± 32.1 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)%timeit db.search(input_vec, k=1000)
**1.41 ms** ± 264 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
基本用法和结果
获得查询结果后，可以通过查找对象 id 或使用模型的解码器组件从潜在表示中重建原始光曲线来绘制结果。我手动选择了一些表面上看起来有点独特的曲线，并进行了搜索。我在下面显示了输入曲线和前 3 个结果。尽管我上面提到了余弦相似性的警告，但我们可以看到，潜在的表示捕捉了许多我们用来将两个时间序列视为相似的属性。
输入光曲线(红色)和查询结果(黑色)。结果是 z 归一化的。
定向搜索
值得注意的是，我们的输入向量不一定来自数据库内部。例如，我们可能已经有了一个特定现象的理论模型，但是无法在数据库中找到一个现有的向量来找到类似的对象。如果我们知道如何用数学方法表达预期的光变曲线或者从模拟中得到一个结果，我们只需将时间序列输入到编码器中就可以得到表示。这成为我们查询数据库的输入向量。
顺便提一下，提到这如何能帮助寻找外星智慧是很有用的。通常搜索寻找不寻常的无线电或光学激光脉冲发射。然而，考虑到太阳帆、遮阳篷、轨道栖息地或任何故意设计来宣布它们存在的结构的潜在大小，这开辟了凌日方法作为这种搜索的有效方法。关于这些不自然的光变曲线看起来会是什么样子，已经有了研究和建议。
吕克的论文。F. A. Arnold 题为人造物体的凌日光变信号研究了我们期望从不期望自然存在的物体中看到的东西，包括 1)等边三角形 2)两个屏幕和 3)类似百叶窗的六个屏幕。它提醒我们，模拟存在单个对象、多个对象以及它们是否旋转的情况也很重要。值得注意的是，如果测量值低于一定的分辨率，仍然有可能将一些物体与环形行星等事物混淆。因此，给出的结果是与同等大小的球形物体相比的幅度差，而不是原始光曲线。
人造物体的例子(左)和亮度大小差异(右)。图像来源:人造物体的过境光变特征(2005)
外星文明可能有意制造一个明显是人造的重复信标。这种信号的候选者有很多，但一个例子是二进制中圆周率的前几个数字。然而，这可能不太可能在可见光谱中实现，因为这不是能量有效的。但有争议的是，一个试图吸引注意力的文明可能会通过在其他人的常规天文观测中试图被他们看到来做到这一点。如果是这种情况，那么在凌日法观测过程中，一系列产生明显不自然光变曲线的物体就值得寻找。
用二进制信号表示的圆周率的前几个数字。
要点是，给定一个理论化的光变曲线(自然的或人工的)，我们可以使用这种方法在可接受的时间内快速找到候选者。
最不相似的物体
除了将一个向量与整个数据库进行比较，还可以有一批输入向量。可以预见，这比一次做一个要快，并且伸缩性好。
发现奇怪光曲线的一种方法是在整个数据库中进行完全的成对相似性搜索，这是可能的，因为我们可以进行批量输入。由于 Faiss 返回顶部结果的相似性分数，我们可以推断最奇怪的光曲线是那些顶部相似性分数仍然很低的(找到最小的最大值)。这可以通过 thresh-holding 来完成，或者简单地根据数据库批处理查询的得分从最低到最高排序，返回 k =1 个最相似的结果。
下面是 TESS 扇区 16 和 17 的前 20 个最不相似的光变曲线。该操作耗时 1.1 秒:
Tess 部门 16 和 17 中最不相似的 20 个结果
关于这个结果的一些观察结果——有 4 或 5 个观察结果可能是在测量开始时具有大尖峰的小故障。除此之外，我们还有各种短周期和高变化曲线的组合。
更广阔的电磁频谱和未来的工作
对于那些不熟悉物理的人来说，我们用眼睛看到的只是更广泛的电磁(EM)光谱的一小部分。本质上，我们看到的是不同波长(或能级)的光子，每个波长有不同的颜色。在长波长，我们有无线电频谱，而可见光，红外线，X 射线和伽马射线组成了较短的波长。用专用的硬件来观察不可见的波长是可能的，这在天文学中是很常见的。这是因为电磁辐射揭示了许多物理过程，有些物体只有在这些波长下才可见。下面是对猎户座的可见光和红外观测的比较。
可见光(左)和红外(右)中的同一个星座。图片来源:NASA
TESS 卫星在可见光和近红外范围(600 纳米-1000 纳米)进行观测，因为这被认为是寻找类地系外行星的良好范围。上述方法的真正潜力在于更广泛的电磁频谱。
机器学习方法的一个优点是，我们可以训练编码器接受许多输入，并为每个对象提供一个组合的矢量表示。使用这种编码器的矢量数据库可以按照相同的过程建立。如上所述，我们可以根据我们对整个电磁光谱的预期来查询预期的光变曲线。例如，一颗红外辐射很强但没有可见光的恒星可能被尘埃包围，或者被一个结构包围。
我们可以将多个输入组合成相同的表示，并遵循相同的流程
另一个好处是，如果某个源的观测测量值不可用，还可以处理某些输入中的缺失数据。如果输入来自不同的观测站，而这些观测站的访问受到限制，或者观测是在不同时间协调的，就可能出现这种情况。尽管输入的不同部分存在差异，但网络可以概括为具有足够的表达能力，ODE-RNNs 有可能很好地处理这一点。
虽然可能需要更长的训练时间，但这不应该降低推理时间的性能或增加内存消耗，因为我们在可视化、聚类和搜索过程中仍然使用相同大小的向量。
摘要
我已经证明，如果我们使用自动编码器将凌日光变曲线测量结果编码成潜在空间表示，就有可能以合理的硬件要求对大型数据集进行快速有效的分析。神经常微分方程的最新发展使我们能够以精确的方式处理不规则采样的时间序列和可变长度的测量。
即使我们的硬件有所改进，我们开始精确地积累数十亿(或数万亿)的对象，我们肯定有能力快速查询数据。在编码器和解码器架构以及聚类算法方面可以继续进行大量的改进。
此外，在时间序列分析是关键的更广泛的物理学和宇宙学领域，有很大的机会应用这样的方法。
如果你想联系我，请随时在 glenn.kroegel@gmail.com 联系我，或者在 LinkedIn 上给我发消息。
揭开魔法:解读机器学习黑盒模型
你是否曾经开发出一个机器学习模型，它具有很高的准确性和惊人的 AUC，但却不知道如何解释这些预测？
在预测能力和可解释性之间进行权衡是使用黑盒模型时面临的一个常见问题，尤其是在必须向非技术受众解释结果的业务环境中。可解释性对于质疑、理解和信任人工智能和人工智能系统至关重要。它还为数据科学家和工程师提供了更好的方法来调试模型，并确保它们按预期工作。
动机
本教程旨在介绍在黑盒模型中处理模型解释的不同技术。
**免责声明:**本文旨在向普通数据科学家介绍可解释机器学习领域的一些有用技术，并推动其采用。他们中的大多数已经从 Christoph Molnar 的这本非常值得推荐的书中总结出来:可解释的机器学习。
本文中使用的完整代码可以在我的 GitHub 中找到
内容
- 可解释性方法的分类
- 数据集和模型培训
- 全球重要性
- 本地重要性
1.可解释性方法的分类
- **内在还是事后？**该标准区分可解释性是通过限制机器学习模型的复杂性(内在)还是通过在训练后应用分析模型的方法(事后)来实现。
- **特定于模型还是与模型无关？**线性模型有特定于模型的解释，因为回归权重的解释特定于该类模型。同样，决策树拆分也有其特定的解释。另一方面，模型不可知的工具可以用在任何机器学习模型上，并且在模型已经被训练(后期)之后被应用。
- **本地还是全球？**局部可解释性是指解释单个预测，而全局可解释性与解释预测任务中的模型一般行为相关。这两种类型的解释都很重要，并且有不同的工具来处理它们。
2.数据集和模型培训
本文使用的数据集是来自 UCI 机器学习知识库的 成人普查收入 。预测任务是确定一个人的年收入是否超过 5 万美元。
由于本文的重点不在 ML 管道的建模阶段，所以执行了最小特征工程，以便用 XGBoost 对数据进行建模。
为该模型获得的性能指标如下:
图 1:训练和测试装置的接收操作特性(ROC)曲线。
图 2: XGBoost 性能指标
这个模型的性能似乎还可以接受。
3.全球重要性
用于评估模型整体行为的技术包括:
3.1 -特征重要性(由 XGBoost 模型和 SHAP 评估)
3.2 -汇总图(SHAP)
3.3 -排列重要性(ELI5)
3.4 -部分相关图(PDPBox 和 SHAP)
3.5 -全局代理模型(决策树和逻辑回归)
3.1 -功能重要性
- XGBoost(特定型号)
feat_importances = pd.Series(clf_xgb_df.feature_importances_, index=X_train.columns).sort_values(ascending=True)
feat_importances.tail(20).plot(kind='barh')
图 3: XGBoost 特性的重要性
使用 XGBoost 时，在解释特性的重要性时必须小心，因为结果可能会产生误导。这是因为该模型计算了几个具有不同解释的重要性指标。它创建了一个重要性矩阵，这是一个表，第一列包括提升树中实际使用的所有特征的名称，另一列是用不同的度量标准计算的结果“重要性”值*(增益、覆盖、频率)*。在这里可以找到更全面的解释。
增益是解释每个特征的相对重要性(即准确性的提高)的最相关属性。
- SHAP
总的来说，SHAP 库被认为是一个解决可解释性的模型不可知的工具(我们将在本地重要性章节中讨论 SHAP 的直觉)。但是，对于基于树的机器学习模型，如决策树、随机森林和梯度提升树，该库有一种特定于模型的方法。
explainer = shap.TreeExplainer(clf_xgb_df)
shap_values = explainer.shap_values(X_test)shap.summary_plot(shap_values, X_test, plot_type = 'bar')
图 4: SHAP 特征重要性
XGBoost 特征重要性用于评估训练数据集模型输出中预测值的相关性，SHAP 特征重要性用于评估测试数据集模型输出中预测值的相关性，以评估最重要的特征在两种方法和集合中是否相似。
据观察，模型中最重要的变量得到了保留，尽管其重要性顺序不同(年龄似乎在 SHAP 方法的测试集中具有更大的相关性)。
3.2 剧情概要(SHAP)
SHAP 汇总图是评估模型特性的一个非常有趣的图，因为它提供了比传统的特性重要性更多的信息:
- 特征重要性:变量按重要性降序排列。
- 对预测的影响:横轴上的位置表示每个特征的数据集实例的值对模型输出的影响是多还是少。
- 原始值:颜色表示每个特性的高值或低值(在每个特性的范围内)。
- 相关性:特征与模型输出的相关性可以通过评估其颜色(其取值范围)和对水平轴的影响来分析。例如，观察到年龄与目标正相关，因为对输出的影响随着特征值的增加而增加。
shap.summary_plot(shap_values, X_test)
图 5: SHAP 汇总图
3.3 -排列重要性(ELI5)
评估预测值全局重要性的另一种方法是随机排列数据集中每个要素的实例顺序，并使用训练好的模型进行预测。如果通过在顺序中进行这种干扰，评估度量没有实质上改变，则该特征不那么相关。相反，如果评估度量受到影响，那么该特征在模型中被认为是重要的。这个过程是为每个特征单独完成的。
为了评估训练的 XGBoost 模型，ROC 曲线的曲线下面积(AUC)将被用作性能度量。将在训练和测试中分析排列重要性:
# Train
perm = PermutationImportance(clf_xgb_df, scoring = 'roc_auc', random_state=1984).fit(X_train, y_train)
eli5.show_weights(perm, feature_names = X_train.columns.tolist())# Test
perm = PermutationImportance(clf_xgb_df, scoring = 'roc_auc', random_state=1984).fit(X_test, y_test)
eli5.show_weights(perm, feature_names = X_test.columns.tolist())
图 6:训练集和测试集排列重要性。
尽管最重要的功能的顺序发生了变化，但看起来最相关的功能保持不变。有趣的是，与 XGBoost 特征重要性不同，训练集中的年龄变量具有相当强的影响(如测试集中的 SHAP 特征重要性所示)。此外，根据排列重要性的 6 个最重要的变量被保持在训练和测试中(顺序的差异可能是由于每个样本的分布)。
估算全球重要性的不同方法之间的一致性使模型输出的解释更加可信。
3.4 -部分相关图(PDPBox 和 SHAP)
部分相关图(PDP) 表示某一特性对预测输出的边际影响。为此，在其他条件不变的情况下，修改特征，并观察平均预测的变化。进行的过程如下:
1)选择特征
2)定义网格值
3)为每个网格值:
3.1)用网格值替换特征
3.2)平均预测值
4)绘制曲线
PDP 可以指示特性和输出之间的关系是线性的、单调的还是更复杂的。值得注意的是观察到的关系是与预测的关系，而不是与目标变量的关系。然而，根据模型的性能，可以产生目标对所评估特征的依赖性的直觉。
PDP 的优势是非常 容易实现 并且 非常直观 :如果所有数据点都被强制采用每个特定值，那么特定特征中的函数代表平均预测。
我们将使用 PDPBox 和 SHAP 分析部分相关图。
- PDPBox
例如，将分析两个最相关的观察特征的 PDP:
# Create the data that we will plot
pdp_education = pdp.pdp_isolate(model=clf_xgb_df, dataset=X_test, model_features=X_test.columns, feature='education.num')
pdp_age = pdp.pdp_isolate(model=clf_xgb_df, dataset=X_test, model_features=X_test.columns, feature='age')# Plot it
pdp.pdp_plot(pdp_education, 'education.num',figsize=(12, 6))
pdp.pdp_plot(pdp_age, 'age', figsize=(12, 6))
plt.show()
图 7:教育的部分相关图
看起来受教育年限(从 7 年开始)和收入超过 5 万美元的概率之间存在线性关系。这一特征在模型输出中的影响被证明高达 0.6(满分 1)。
图 8:年龄的部分相关性图
人们似乎更有可能在 50 多岁的时候赚到超过 5 万美元。
- SHAP 依赖情节
将使用 SHAP 方法生成相同的 PDP。该库除了指示特征对模型输出的边际影响之外，还通过颜色指示与最常交互的特征的关系。
shap.dependence_plot('education.num', shap_values, X_test)
shap.dependence_plot('age', shap_values, X_test)
图 9:教育的部分相关图
尽管 y 轴刻度不同于 PDPBox 图*(我们将在本地可解释性部分中了解原因)*，但“education.num”的趋势似乎与之前的方法相同。此外，SHAP 已经识别出特征【已婚 _1】是与它交互最多的特征(这意味着，对于该模型来说，具有高教育年限的已婚人士更有可能收入超过 5 万美元)。
图 10:年龄的部分相关性图
该方法中的年龄趋势与 PDPBox 方法一致。和它互动最多的功能是“education.num”。
在陈述了 PDP 对可解释性领域的优势后，也值得(也公平地)陈述一下劣势:
- 它没有考虑特征 的分布:它可能会产生误导，因为它可能会错误地解释数据非常少的区域(通过给所有数据点分配这些值，特征值被过度表示，这可能会导致错误的结论)。
- 特性独立性假设 :这是 PDP 最大的弊端之一。假设为其计算部分相关性的特征与其余预测值不相关。
- 异质效应可能被隐藏 :这是因为只计算平均边际效应。在极限情况下，PDP 可能是一条水平线，值均匀分布在上下，表明该特征对预测没有影响。
为了克服 PDP 的一些缺点，可以使用个体条件期望(ICE) 和累积局部效应(ALE) 图。尽管这些方法的实现没有包含在本文中，但是我们将简要地解释它们，以展示它们如何改进 PDP 方法。
个人条件期望(ICE)
ICE 图是单个数据点的 PDP 等效图。该图为数据集的每个实例显示一条线，指示该实例的预测如何随着要素值的变化而变化。PDP 是 ICE 图中所有线条的平均值。ICE 图允许可视化边际效应的变化，能够检测异质效应。
累积局部效应(ALE)图
当一个特征与其他预测因子高度相关时，PDP 会出现严重的问题，因为实际上不太可能发生的合成实例预测会被平均化(例如，年龄为 16 而教育程度同时为 10 是不太可能的)。这可能会在估计该特征的效果时产生显著的偏差。ALE 图除了计算速度更快之外，还是计算特征对模型预测的影响的无偏解决方案，因为它们评估的是其条件分布。也就是说，对于网格的值 x1，他们仅使用值与 x1 相似的实例的预测进行估计，从而避免使用现实中不可能的实例。
此外，为了估计特征对预测的影响，他们不使用平均值(将特征的影响与所有相关预测因子的影响相混合)，而是计算预测之间的差异。
冰和啤酒的区别
ICE 图解决了 PDP 可能呈现的异质效应的问题，但是没有解决由于相关特征引起的偏差。相反，ALE 图解决了偏倚问题，考虑了特征的条件分布及其与其余预测因子的相关性。
3.5 -全球代理模型
全局代理模型是一种可解释的模型，它被训练成近似黑盒模型的预测。我们可以通过解释代理模型得出关于黑箱模型的结论。用 Christoph Molnar 的话说:“通过使用更多的机器学习来解决机器学习的可解释性！”
我们将尝试使用逻辑回归和决策树作为全局代理模型来近似 XGBoost。
- 逻辑回归
# Train
log_clf = LogisticRegression().fit(X_train, y_train)# Predictions
y_pred_train_log = log_clf.predict(X_train)
y_proba_train_log = log_clf.predict_proba(X_train)[:, 1]
y_pred_test_log = log_clf.predict(X_test)
y_proba_test_log = log_clf.predict_proba(X_test)[:, 1]# R-squared
print('R-squared Train RL-XGB: ', r2_score(y_proba_train_log, y_proba_train))
print('R-squared Test RL-XGB: ', r2_score(y_proba_test_log, y_proba_test))
图 11:逻辑回归和 XGBoost 预测之间的 R 平方。
训练集和测试集的 R 平方都是负的。当拟合比简单地使用平均值更差时，就会发生这种情况。所以得出结论 Logistic 回归不是一个好的替代模型。
- 决策树
# Train
tree_clf = tree.DecisionTreeClassifier(random_state=0, max_depth=4).fit(X_train, y_train)# Predictions
y_pred_train_tr = tree_clf.predict(X_train)
y_proba_train_tr = tree_clf.predict_proba(X_train)[:, 1]
y_pred_test_tr = tree_clf.predict(X_test)
y_proba_test_tr = tree_clf.predict_proba(X_test)[:, 1]# R-squared
print('R-squared Train DT-XGB: ', r2_score(y_proba_train_tr, y_proba_train))
print('R-squared Test DT-XGB: ', r2_score(y_proba_test_tr, y_proba_test))# Metrics
clf_metrics(y_pred_train_tr, y_proba_train_tr, y_train, y_pred_test_tr, y_proba_test_tr, y_test)
图 12:决策树和 XGBoost 预测之间的 R 平方，以及前者的性能度量。
XGBoost 模型预测中的方差与决策树非常接近，因此它可以作为解释主模型的代理模型。事实上，性能指标也非常接近原始模型。
值得注意的是，虽然 XGBoost 预测的方差可以由决策树 很好地解释，但不能保证 后者会以与前者相同的方式使用这些特性。可能发生的情况是，树在输入空间的某些区域中正确地逼近 XGBoost，但是在其他区域中表现得完全不同。
将对结果树进行分析，以便评估所使用的特征是否对应于迄今为止已经检测到的最重要的特征:
# Plot tree
fig, ax = plt.subplots(figsize=(30, 10))
tree.plot_tree(tree_clf, feature_names= X_train.columns.to_list(), ax=ax, filled=True)
plt.show()
图 13:经过训练的决策树
该树用于估算收入的 5 个特征，按重要性排序为:
1.已婚 _1
2。教育编号
3。资本收益。资本损失
5。年龄
这些特征对应于其他方法检测到的最重要的特征。
4.本地重要性
局部代理模型是用于解释黑盒机器学习模型的个体预测的可解释模型。
4.1 -局部可解释的模型不可知解释(LIME)
LIME 分析当输入数据发生变化时，模型预测中会发生什么。它生成一个新的数据集，其中包含置换样本及其来自原始模型的相应预测。在这个合成集上，训练可解释的模型(逻辑回归、决策树、套索回归等)。)，然后通过采样实例与感兴趣实例的接近度对其进行加权。
对实例 X 的解释将是最小化损失函数(代理模型的预测和原始模型的预测之间的性能测量，例如 MSE)的代理模型的解释，保持模型的复杂度低。
# Generate explainer
explainer = lime.lime_tabular.LimeTabularExplainer(X_train.values, mode='classification',feature_names=X_train.columns.tolist(), discretize_continuous=False, random_state=1984)# Generate explanation for instance i
i = 546
exp = explainer.explain_instance(X_test.values[i], clf_xgb_array.predict_proba)# Plot
fig = exp.as_pyplot_figure();
图 14:每个特征在个体预测中的相对重要性。
# Generate explanation for sample
sp_obj = submodular_pick.SubmodularPick(explainer, X_test.values, clf_xgb_array.predict_proba, sample_size=3, num_exps_desired=3)[exp.show_in_notebook() for exp in sp_obj.sp_explanations]
图 15:对抽样实例的单独解释。
据观察，在所有区分阶级的个人解释中，最有影响力的特征是资本收益。接下来，根据具体情况，最相关的预测值是已婚、教育编号、年龄和性别。这些特征与在全局重要性算法中识别的特征相同。
4.2 -沙普利添加剂解释(SHAP)
SHAP 是一种解释基于沙普利值计算的个人预测的方法，一种来自联合博弈论的方法。它试图回答问题**“相对于平均预测，每个特征值对预测的贡献有多大？”**。为此，Shapley 值根据“玩家”对“总支付”的贡献将“支付”分配给“玩家”。玩家在联盟中合作，并为此类合作获得一定的奖励。
在机器学习上下文中，“游戏”是数据集实例的预测任务。 “总支付额”是对该实例的预测，减去对整个数据集的平均预测。 “玩家”是实例的特性值，它们联合起来接受“报酬”(预测)。Shapley 值是所有可能组合的特征值的平均边际贡献。它表示“总支出”(预测)如何在所有“玩家”(特征值)之间分配。
SHAP 带来的一项创新是，Shapley 值解释被表示为一种附加特征归因方法，即线性模型。就这样， SHAP 把石灰的好处与沙普利值 联系起来。
# Create explainer
explainer = shap.TreeExplainer(clf_xgb_df, model_output='probability', feature_dependence='independent', data=X_test)# Generate explanation for instance i
i= 150
data_for_prediction = X_test.iloc[i, :]shap_values = explainer.shap_values(data_for_prediction)
shap.initjs()
有几种方法可以把 SHAP 的解释形象化。我们将在本文中讨论其中的两个:力图和决策图。
力图
shap.force_plot(explainer.expected_value, shap_values, data_for_prediction)
图 16:单个实例的 SHAP 力图解释
力图显示了每个特征对预测的影响。有两个相关值需要注意:输出值**(实例的模型预测)和基本值**(整个数据集的平均预测)。更大的条意味着更高的影响，颜色指示特征值是否将预测从基础值向 1 (红色)或 0 (蓝色)移动。
决定情节
**shap.decision_plot(explainer.expected_value, shap_values, data_for_prediction)**
图 17:单个实例的 SHAP 决策图解释
决策图显示了与力图基本相同的信息。灰色垂直线是基值，红线表示每个特征是否将输出值移动到比平均预测值更高或更低的值。
这个图可以比前一个图更清晰和直观一点，尤其是当有许多特征要分析的时候。在力图中，当预测因子的数量很高时，信息可能看起来很浓缩。
结论
这篇文章旨在帮助数据科学家更好地理解他们的机器学习模型是如何工作的，并能够以一种清晰的方式解释这些结果。它对于调试模型和确保它们按预期工作也很有用。
我们已经介绍了可解释性方法的不同分类(内在的/事后的，模型特定的/模型不可知的，局部的/全局的)并且我们使用了几个库和技术来评估全局和局部的重要性。
总之，使用的库和技术有:
- XGBoost :特性重要性
- ELI5 :排列重要性
- PDPBox :部分依赖图
- 全局代理模型:逻辑回归，决策树
- 石灰:当地的重要性
- SHAP :特征重要性、概要图、部分依赖图、局部重要性
那么，哪一个库是解决 ML 模型可解释性的最好的库呢？ 在我看来，使用几个库和技术有助于建立模型输出的可信度(前提是结果一致)。然而，如果我必须选择一个，我肯定会去 SHAP。**
SHAP 在可解释机器学习领域做出了巨大贡献。这是因为在这里全局解释与个体解释是一致的，因为 Shapley 值是全局解释的“原子单位”(在博弈论中有坚实的理论基础)。例如，如果 LIME 用于局部解释，而 PDP 或排列重要性用于全局解释，则这两种方法之间没有共同的理论基础。
我希望这篇文章能达到它的目的，作为破解黑盒模型的一般指南。完整代码可以在我的 GitHub 中找到
在我的下一篇文章中，我将讨论模型公平性、，它在过去的几年中已经获得了越来越多的关注。该领域旨在评估该模型在处理数据中预先存在的偏见时的公平性: 职位匹配系统在 CEO 面试中偏向男性候选人，因为这与历史数据相匹配，这公平吗？
敬请期待！
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
