成功的数据科学之路在于你的学习能力。但是,应该学习什么呢?
成功交付数据科学项目的机会最大时,你需要不断学习,但要专注于哪些方面却不总是明确的
·发布于 Towards Data Science ·7 分钟阅读·2023 年 6 月 29 日
–
图片由 Tamara Menzi 提供,来源于 Unsplash
在过去的十年中,数据科学取得了许多伟大的进展,但尽管如此,许多项目却从未实现。作为数据科学家,我们不仅要展示强大的技术能力,还要理解业务背景,有效地与利益相关者沟通,并将他们的问题转化为可以推动业务价值的可操作性建议。这是否合理,还是说业务在寻找新的独角兽?在这篇博客中,我将描述业务如何随着岁月的推移而变化,这将为你提供一个更好的视角,以了解你可能需要学习哪些内容才能成功交付数据科学项目。
简短介绍
十多年前,公司认识到挖掘数据集可以产生增加收入、优化流程和降低(生产)成本的信息。这催生了一个新领域和新角色;数据科学领域及数据科学家。但企业的需求多年来不断变化。因此,了解企业的需求对于了解作为数据科学家需要学习什么非常重要。 在接下来的部分中,我们将首先回顾描述数据科学领域在过去十年中的演变。这可以帮助你了解;1. 什么 曾经 是重要的学习内容,2. 现在 是什么 重要的学习内容,以及 3. 什么 可能是 未来工作中重要的学习内容。接下来的部分让我们一起回顾过去。
首先是科学程序员。
数据科学有许多学科,其基础建立在源自数十年(学术)研究和开发的统计学和数学上。许多最初的核心算法在文本挖掘、图像识别、传感和时间序列等学科中形成了基础。早期,这些方法的发布并没有附带代码。为了应用这些方法,公司会聘请科学程序员来完成具有挑战性和耗时的任务。 但在编写一行代码之前,通常会有一个思考过程,以了解为何要付出这些努力,以及可以期待什么样的结果。在过去十年中,这一情况发生了剧变,因为谷歌、Meta 等公司开始开源它们的库。此外,社区也开始开发开源包,如 sklearn、scipy 等。现在安装只需一行代码。
数据科学领域发展迅速,但企业到底需要什么?
现在,科学程序员已经变成了数据科学家。不过,情况发生了一些变化。企业也需要能够与利益相关者有效沟通、识别商业机会,并将技术洞察转化为推动业务价值的可操作建议的数据科学家。这催生了一种新型的数据科学家;应用型数据科学家。
应用型数据科学家与基础型数据科学家
“数据科学家”这个词经常被用作数据科学领域不同角色的统称。从数据管理员、机器学习工程师、数据工程师、统计学家以及更多角色。但当我们谈论真正的数据科学家时,基本上有两种类型;基础型数据科学家和应用型数据科学家。
-
基础数据科学家具有强大的统计和机器学习技术知识,以分析复杂的数据集并得出见解。此人可以告诉你所有关于底层数据分布的内容,并能轻松创建(修改)解决问题的算法/方法。这些人通常在研究和开发环境以及学术机构中表现出色。
-
应用数据科学家专注于应用现有技术和方法来解决特定业务问题或开发数据驱动的产品和/或服务。通常,这些角色在某一领域内加速发展,例如文本挖掘、图像识别、传感器或时间序列。创新通常通过在领域相关的数据上应用新方法来实现,而不是创建新算法或方法。
两个角色各有其优缺点,我有三个提示可以帮助成功交付数据科学项目。
提示 1:学习编程基础。
感谢 Coursera、Udemy、YouTube 和 Medium 等优秀平台,有大量材料可以学习编程基础。
-
以已知风格编写代码,如 PEP8。
-
编写内联注释;说明你做的事情及其原因。
-
编写文档字符串。
-
使用合理的变量名。
-
降低代码复杂性(尽量降低)。
-
编写单元测试。
-
编写文档。
-
保持代码整洁。
编程是数据科学领域的主要挑战之一。它被严重低估,但却是决定数据科学项目能否进入生产的关键因素之一。想想看,你是否愿意维护一个没有文档、没有单元测试且写成意大利面式风格的模型或编程代码?我不这么认为。
每个数据科学项目都需要可重复的代码,而进入生产的步骤则需要可维护的代码。最终,每个项目都不过是一堆代码行,需要有人将其投入生产。保持整洁。保持有序。
提示 2:项目的成功不仅仅依赖于机器学习解决方案。
数据科学项目通常开始时充满热情,但这种热情可能会迅速消退,因为项目需要的远不止一个机器学习解决方案。在最近发表的一篇文章[1]中,非常清楚地描述了数据科学项目中最重要的技术步骤。然而,要将一个想法转化为生产,你需要的不仅仅是技术技能。以下是可以帮助提高项目成功率的步骤总结:
-
从目标开始。在项目开始时,了解如何将项目落地到组织或公司中。数据治理、伦理和隐私是开始时重要的。
-
检查哪个平台或基础设施可供协作。 例如,这可以是带有 CI/CD 管道的 git 和 cookie-cutter 模板。
-
理解领域。 在进行任何分析之前,需要对你所工作的领域有基本的理解。你需要知道如何根据你工作的领域和背景来处理数据。没有一种通用的数据科学解决方案。
-
***正确进行数据分析。***这看似简单,但知道如何用
pip install
安装一个包并不等于你就是专家。进行自己的研究和阅读文章。避免那些你无法解释的(复杂的)机器学习解决方案。使用训练-测试-验证集。将你的结果与基准比较。与有经验的科学家和具有领域知识的人讨论你的想法和结果。 -
***报告你的结果。***保持透明。讲述基于事实的故事。不要将故事泛化超出数据范围。描述过程比模型输出的那个单一数字更为重要。
-
编写可重复和可维护的代码。 证明结果是可重复的,代码是可维护的。
-
交付结果。 如果所有步骤都完成了,结果或产品需要以客户能够使用的方式交给他们。仅仅将带有工作代码的个人笔记本电脑交给客户并不是解决方案。
如果你仔细观察这些步骤,只有一个步骤(#4)涉及数据分析和模型创建。让这一点深入人心。
提示 3:聪明地学习并重复。
数据科学是一个高度复杂且快速发展的领域,不同的专业领域在这里汇聚。每个数据科学家都有不同的背景,持续学习是这一过程的一部分。这意味着个性化的学习和成长路径可以带来很大益处,这取决于你的学位/起点、经验、领域知识、数学、统计、编程、工程、沟通和演示技能的背景。与同事讨论你可以改进的地方,并制定个人的学习和成长路线图。请注意,参加随机的数据科学课程可能很有趣,但可能与公司的使命甚至你的个人成长路径不一致。
学习的能力是每个人都应该不断练习的技能,终身学习可能是你能给自己最好的礼物。
总有更多的东西可以学习。
成功的道路不是某个特定的在线课程,而可能需要多年甚至几十年的奉献、努力和奋斗。投资自己,学习基础知识,超越肤浅的了解,专业化,并认识到成功是许多小步骤的积累,其中建模部分只是整个过程中的一步。
让我尝试将这分解成几个部分。首先,沟通非常重要。也许你能创造出最天才的方法,但你需要有效地将复杂的技术概念阐述给技术和非技术的利益相关者。问题解决:你应该能够用结构化和系统化的思维方法处理复杂问题。批判性思考,从多个角度分析问题,并提出有效的解决方案。你可以通过在Stack Overflow等网站上帮助社区来轻松练习。当你在职业生涯和职位上不断成长时,你应该能够指导和辅导开发人员。提供指导,分享最佳实践,并帮助提升他们的技术技能。适应能力。你不应固守于你所知道的单一技术,而应拥抱新的技术、方法和工具。你应该能够快速学习并适应不断变化的项目要求或行业趋势。时间管理。有效管理你的时间。优先处理任务,遵守截止日期,平衡各种需求。专注于交付高质量的工作。
保重。保持警惕。
Cheers E.
如果你觉得这篇文章对你有帮助,欢迎 关注我 ,因为我写了更多关于贝叶斯因果学习的内容。如果你考虑购买 Medium 会员,你可以通过我的 推荐链接来支持我的工作。这与购买一杯咖啡的价格相同,但可以让你每月阅读无限的文章。
让我们联系吧!
参考文献
-
Michael A. Lones, 如何避免机器学习陷阱:学术研究者指南,arXiv: 2108.02497
-
Tessa Xie, 数据科学职业错误避免, 2021,
-
数据科学家是否变成了过时的职业? 数据科学中央
人工智能监管之路
原文:
towardsdatascience.com/the-path-towards-ai-regulation-7199eb803040
意见
《人工智能法案》的贡献及其全球影响
·发表于Towards Data Science ·7 分钟阅读·2023 年 3 月 10 日
–
2021 年中期,欧盟委员会提出了针对欧盟的人工智能(AI)新监管框架。因此,它代表了对人工智能实施某种程度监管的首次尝试。
你可以阅读AI 法案提案的完整官方文件(24 种语言版本)
目前,它仍然是最先进且被广泛接受的相关监管提案。虽然它仍处于初步提案阶段,但欧盟正在努力推动短期内的强烈议程以实现采纳。
目前对人工智能系统的定义非常宽泛和肤浅。
这一倡议受到全球的密切关注,因为它源于对保护个人基本权利的需求,因此可能对整个社会产生影响。
这一法律框架从使用和风险的角度考虑人工智能系统,同时接受了对人工智能及其系统概念的任何定义的局限性。基本上,目前对人工智能系统的定义非常宽泛和肤浅,因为它止步于将人工智能系统定义为一种软件技术,包括机器学习、逻辑和知识基础系统以及统计方法。
人工智能法案根据风险等级对人工智能系统进行分类,并根据风险等级规定了限制和义务。
图片由janilson furtado提供,来源于Unsplash
来自 Towards Data Science 编辑的说明: 虽然我们允许独立作者根据我们的 规则和指南**发布文章,但我们不对每位作者的贡献表示认可。你不应在未寻求专业建议的情况下依赖作者的作品。详细信息请参见我们的 读者条款*。
AI 风险分类
规章提案并不限制技术或方法。相反,它采用了一种基于风险的方法,提出了一个根据风险等级对 AI 系统进行分类的方案,并根据风险等级指定了不同的要求、限制和义务。特别是:
-
呈现‘不可接受’风险的 AI 系统将被禁止。 这些包括涵盖操控性‘潜意识技术’、或利用特定易受害群体、社会评分目的以及实时远程生物特征识别系统的 AI 系统。
-
‘高风险’ AI 系统将被授权,但需遵守非常严格的要求和义务,才能进入欧盟市场。这些包括作为关键产品安全组件的 AI 系统、生物特征识别、关键基础设施的运行、教育工具;就业和工人管理;获取必要服务的工具;执法;移民、庇护和边境控制;以及司法和民主过程的管理。
-
仅呈现‘有限风险’的 AI 系统将仅需遵守基本的透明度义务。这些包括与人互动的系统(即聊天机器人)、情感识别、生物特征分类,以及生成或操控文本、图像、音频或视频内容(包括深度伪造)的系统。
注意这里的(遗漏)大象:军事和防御应用不在规章范围之内。实际上,《AI 法案》明确指出
专门为军事用途开发或使用的 AI 系统应排除在本规章范围之外。
同样也在
第 2.3 条:本规章不适用于专门为军事目的开发或使用的 AI 系统。
难道它们不应该受到监管吗?
图片由 Tingey Injury Law Firm 提供,发布于 Unsplash
AI 合规执行
从某种意义上说,与这种基于风险的视角一起,这个框架还涵盖了治理和执行机构的定义:它建立了一个欧洲人工智能委员会(监督过程并确保各国的一致实施)和一个每个成员国的国家监督机构(监督法规的应用和实施),同时国家市场监督机构将负责通过广泛的保密访问(包括源代码)来监控参与者的合规性。
执法机构将采取适当措施限制、禁止、召回或撤回不符合规则的人工智能系统。这包括最高达 3000 万欧元或全球年总营业额 6% 的行政罚款,具体取决于违规的严重性。
人工智能法案及其相关执行机制采取了产品责任法的形式。
照片由 Tingey Injury Law Firm 提供,发布在 Unsplash。
从务实的角度来看,基于风险的方法是完全有意义的。然而,这对未来的可持续性和普遍性提出了若干挑战。特别是:
-
人工智能的定义是否足够广泛和精确?
-
我们如何确保我们将应用程序分类到正确的风险类别?
-
对于我们现在无法预见的新解决方案和应用,怎么办?
此外,请注意,人工智能法案及其相关执行机制基本上采取了产品责任法的形式。的确,这种框架显然是在公司和产品层面制定规则和权利,而不是个人层面。这对生产和服务基于人工智能解决方案的公司有深远的影响。公司需要声明并展示其解决方案和应用的风险水平。它们还需要保证遵守各自风险水平的规则。
考虑到整个平台始于保护个人基本权利的需求,法案所采取的方法似乎与基本原则不太一致。此外,这种方法还隐含着实际风险:如果出现任何问题,究竟由谁负责?是开发公司、供应商,还是任何部署、配置或定制解决方案的中介?
人工智能监管框架还涵盖了其他维度:
-
标准化:与标准化机构如 ISO 和 IEEE 的讨论应确保技术和定义的采用尽可能广泛且达成一致。北约本身现在也在推进相关的标准化倡议。
-
卓越:框架不应(仅仅)限制人工智能的潜力。一方面,它将限制危险的应用,另一方面,它应在基础和应用层面上促进和推动创新与卓越。
如果不同时对数据进行监管和保护,那么对人工智能的监管可能没有任何意义。
至少在欧洲层面上,对这一倡议的重要性和必要性存在广泛的共识。其他国家也开始了类似的进程:例如,美国在 2022 年制定了 人工智能权利法案 的蓝图,而巴西制定了初步的人工智能法律。许多国家(包括美国和许多欧盟成员国)也在 2020 年代初期制定了国家人工智能倡议。欧洲联盟可能会像它在数据保护领域所做的那样,为未来全球人工智能的监管铺平道路。
确实,对不对数据进行监管和保护的情况下对人工智能进行监管似乎没有任何意义。确实,人工智能法规的应用在已经拥有良好建立的数据保护和监管法律的国家中会更为可行。
在这方面,欧洲联盟再次在全球范围内引领潮流,推出了GDPR,该法规目前在更广泛的层面上产生了广泛的影响。
欧洲联盟可能会为未来全球人工智能系统的监管铺平道路,就像它在数据保护方面所做的那样。
AI 法案的详细描述可以在这个欧盟议会简报中找到,而 AI 法案提案的完整官方文件发布在这里(提供 24 种语言)。关于此事的讨论和观点也来源于IdeasLab会议和由布鲁塞尔欧洲政治研究中心(CEPS)于 2023 年 2 月 28 日组织的智库。
[AI 法案提案的官方标题]
总结:AI 法案可以代表国际层面上对人工智能系统的首次尝试性监管。一旦获得批准,它可能会影响并在某种程度上迫使其他几个国家,包括美国,采取类似的举措,可能会与这一提案保持一致或受到其启发。
国际间的协调可能会带来好处,这与近年来的数据监管情况相对立。在过去,欧盟通过 GDPR 法规引领了数据监管,但美国没有匹配 GDPR 要求的保障水平,从而在国际层面上产生了混乱、法律障碍以及数据传输和存储的限制。
另一方面,这可能代表了软件系统和模型将首次被明确和直接限制使用和实施。人们可能会好奇,为什么这对 AI 系统适用,而过去对于其他软件或数学模型没有出现类似需求。这可能会成为未来技术监管的重要先例。
为什么我们需要对 AI 系统进行监管,而过去对其他软件模型和算法没有进行监管?
你怎么看?应该对 AI 进行监管吗?如何进行?
你有自己的看法要分享吗?
基于风险的方法在长远中可能暴露哪些陷阱?
这样的监管会影响(你的)业务、创新过程或数据科学实践吗?
照片由Christopher Burns拍摄,发布在Unsplash上
让我们来看看“混乱的数据科学家的 PATH 变量:如何管理它”
了解 PATH 是什么,以及如何在 Windows 和类 Unix 系统中添加路径
·发表于 Towards Data Science ·6 分钟阅读·2023 年 6 月 26 日
–
图片由我使用 Midjourney 制作
有时我觉得 StackOverflow 上的人好像头顶上装了指南针。他们总是能找到 PATH。
我认为它不在你的 PATH 中。
你可能搞乱了你的 PATH。
你是否将它添加到你的 PATH 中?
检查可执行文件是否在你的 PATH 中。
“我的意思是,PATH 到底是什么?”每当我阅读这样的句子时,我总是会这么说,脸红着试图修复错误。现在,我在数据科学的旅程中已经三年多了,我完全知道它是什么了。差不多。
在这篇文章中,我打算教你如何在 Windows 和类 Unix 系统中管理这个令人困惑的环境变量。
开始吧!
命令也有路径
你使用最多的终端命令是什么?毫无疑问,我的是 git
,因为我在写文章时经常进行提交。
我问这个问题的原因是,因为大多数终端命令在操作系统中也有自己的路径。要找到那个路径,你只需要运行 which valid_command_name
。例如,这是我 git
实例的路径:
$ which git # In windows, use `where`
/usr/bin/git
你可能会认识到上面的是一个绝对路径。但它到底指向什么呢?
一个可执行的二进制文件!(如果你使用的是 Windows,你会得到一个 .exe
文件路径)
它包含了当 git
命令被执行时应该做什么的说明。所有终端命令都是以某种形式的 shell 脚本或可执行文件,执行特定的任务。
例如,你可以对内置终端命令如 clear
、cat
或 touch
使用 which
命令,你会得到绝对路径(甚至是 which
命令本身的路径)。
一旦你得到可执行文件的路径,你可以尝试用其完整路径运行命令,如下所示:
$ /usr/bin/git status -s # the same as `git status -s`
M 2023/6_june/7_path/dump.ipynb
并且它确实有效!但有趣的是——如何使得词语 git
相当于绝对路径 /usr/bin/git
?我意思是,git
只是一个词。
你的计算机如何知道如此多命令的完整可执行路径,并且仅用一个关键词就能运行它们?
(这是 Eureka 时刻!)
使用 PATH 变量!
看啊——PATH!
PATH 是一个变量,保存了执行几乎所有终端命令的绝对路径,无论是内置命令还是用户安装的命令。
当你输入任何命令,如 conda env create
时,操作系统会在你的 PATH 变量中寻找一个名为 conda
的可执行文件。如果找不到可执行文件,它会返回一个 command not found
错误:
$ unknown_command
unknown_command: command not found
你可以用一个命令打印出 PATH 中的所有路径。
在类 Unix 系统中:
$ echo $PATH
/home/bexgboost/.local/bin:/home/bexgboost/anaconda3/bin:...
在 Windows 中:
$ ECHO %PATH%
C:\Python39\Scripts\;C:\Python39\;C:\WINDOWS\system32;...
在 Unix 中,PATH 中的每个路径由冒号分隔,而在 Windows 中由分号分隔。
正是这些路径使得你不需要记住系统中安装的可执行文件的数十个绝对路径。
在 Windows 中处理 PATH
编辑 PATH 变量在类 Unix 系统和固执的 Windows 中有所不同。本节内容讲述在 Windows 中如何管理你的 PATH 变量。
你最常进行的 PATH 操作是添加新路径。例如,要将 C:\Users\user\Desktop\articles
目录添加到 PATH,你应以管理员权限打开命令提示符或 Powershell。然后,运行以下命令:
$ setx PATH "%PATH%;C:\Users\user\Desktop\articles"
但,别急——如果不小心,setx
命令可能会对你的 PATH 造成很大伤害。因此,在 Windows 中,安全的方法是通过系统属性窗口添加路径。
这里是一步一步的指示和一个 GIF,用于将新路径添加到 PATH:
-
点击开始按钮,搜索“编辑系统环境变量”,然后点击第一个结果。
-
从弹出窗口中点击“环境变量”。
-
决定新路径是添加到当前用户还是整个系统。
-
在“用户变量 for your_username”或“系统变量”下,向下滚动,直到找到“Path”变量。
-
选择它并点击 编辑。
-
在新的弹出窗口中,点击“新建”并粘贴新路径,例如
C:\Users\user\Desktop\articles
-
点击确定三次。
GIF 由我制作。
在最后一个弹出窗口中,你也可以删除或更改路径的顺序(见最后一节关于顺序的重要性)。
在 Unix 中处理 PATH
在类 Unix 系统中,添加新的路径到 PATH 要容易得多。
例如,要添加 /home/bexgboost/articles
路径,你需要运行以下命令:
$ export PATH="/home/bexgboost/articles:$PATH"
新路径将被添加到 PATH 的最前面。要将其添加到末尾,你需要改变 $PATH
的位置:
$ export PATH="$PATH:/home/bexgboost/articles"
但新路径只会在当前终端会话期间有效。要永久添加到你的 PATH 环境变量中,你需要将 export
命令追加到你的 shell 配置文件中。
具体的配置文件取决于你使用的 shell:
-
对于 Bash:将
export
命令添加到~/.bashrc
或~/.bash_profile
中(你也可以考虑.profile
,有关区别请参见这里)。 -
对于 Zsh:将
export
命令添加到~/.zshrc
中。 -
对于 Fish:将
export
命令添加到~/.config/fish/config.fish
中。
一旦你从列表中选择了自己的文件(我会选择.bash_profile
),使用以下echo
命令:
$ echo 'export PATH="$PATH:/home/bexgboost/articles"' >> ~/.bash_profile
然后,通过source ~/.bash_profile
重新加载配置文件,操作完成。
路径的顺序
PATH 变量中的路径顺序决定了执行命令或程序时目录被搜索的优先级。
当你在终端中输入一个命令时,操作系统会从左到右依次检查 PATH 中的每一个路径来查找其可执行文件。如果多个目录中包含相同名称的可执行文件,则会使用第一个。
例如,如果你的系统中有多个 Python 版本,则当你在终端中运行python
时,会使用 PATH 中最左边的 Python 可执行文件。
结论
在数据科学和机器学习的世界里,终端常常被忽视。相反,在线课程和训练营将初学者的注意力吸引到像 VSCode 或 JupyterLab 这样更闪亮的 IDE 上。
但总有一天,不可避免地,他们必须在自己的机器上打开生锈的 shell。这时,他们会在 PATH 中遇到一个巨大的障碍。
我希望这篇文章能帮助你击碎那些障碍。
感谢阅读!
喜欢这篇文章,并且,直白点说,它那古怪的写作风格?想象一下访问更多类似的文章,全部由一位才华横溢、迷人、机智的作者(顺便说一下,就是我 😃。
仅需 4.99 美元的会员费,你不仅可以访问我的故事,还可以获取来自 Medium 上最优秀的思想者的知识宝库。如果你使用我的推荐链接,你将获得我超级 nova 般的感激和对支持我工作的虚拟击掌。
[## 使用我的推荐链接加入 Medium - Bex T.
作为 Medium 会员,你的会员费的一部分会分给你阅读的作者,你可以完全访问每一个故事……
ibexorigin.medium.com](https://ibexorigin.medium.com/membership?source=post_page-----b469bfb45785--------------------------------)
Jupyter 的完美文本编辑器:一个完整的 Python IDE
原文:
towardsdatascience.com/the-perfect-text-editor-for-jupyter-a-complete-python-ide-608dd466560b
从语法高亮到代码补全,一个完整的 Python IDE 在 Jupyter 中
·发表于 Towards Data Science ·阅读时间 5 分钟·2023 年 3 月 9 日
–
本文是系列文章的第四部分。查看完整系列:第一部分,第二部分,第三部分。
在过去几天里,我们在 Jupyter 中构建了一个完整的 Python IDE。在本文中,我们将进行最后的调整,并将一切打包成 Docker 镜像,以创建一个便于数据科学家和机器学习工程师使用的便携工作环境。
Jupyter 并不是一个完整的 IDE。它甚至不是许多人认为的 IPython 用户界面。我认为 Jupyter 是一个开发平台,你可以按照自己的喜好创建工作区。由于它包含一个终端模拟器,你可以做你想做的任何事。因此,我们将利用它来创建一个功能丰富的 Python IDE。
因此,在之前的文章中,我们安装了 Neovim,并将其配置为与现在最流行的文本编辑器:Visual Studio Code 相似。在本文中,我们将安装一个全面的工具:代码补全、代码格式化、git 集成、拼写检查和代码对齐。让我们开始吧!
学习速率 是一个面向对 ML 和 MLOps 世界感兴趣的人的新闻通讯。如果你想了解更多类似的主题,请在这里订阅。每个月的最后一个星期天,你会收到我的更新和对最新 MLOps 新闻和文章的见解!
Conquer of Completion
许多工具承诺在 Vim/Neovim 的使用上让开发者的生活更轻松,但我个人偏好的是Conquer of Completion。
Conquer of Completion (CoC) 是一个流行的 Vim 和 Neovim 插件,提供强大的自动补全引擎。以下是有关 CoC 的一些关键点:
-
CoC 支持多种语言:它可以为包括 Python、JavaScript、C、Go 等在内的各种编程语言提供自动补全建议。它使用语言服务器根据你的代码上下文提供智能推荐。
-
CoC 是高度可配置的:它允许用户以类似 VS Code 的方式自定义自动补全引擎。用户可以设置触发补全的自定义映射,调整补全源的优先级等。
-
CoC 拥有庞大的用户社区:这是一个流行的插件,拥有大量用户社区,这意味着有许多资源可以用来排除故障和进行自定义。该插件也在积极维护和更新。此外,社区提供了几个我们将使用的优秀扩展。
总的来说,Conquer of Completion 是一个强大的插件,可以显著提升你在 Vim 或 Neovim 中的自动补全体验,尤其是当你处理多种编程语言时。
安装
CoC 需要一些设置。你需要安装和配置一些额外的组件,包括 NodeJS 和 npm
包管理器。根据情况,你可能还需要为你喜欢的编程语言安装一个语言服务器。
在我们的情况下,首先安装 NodeJS。
NodeJS 安装
NodeJS 和 npm
是安装 Neovim 的 CoC 插件并启用多种编程语言自动补全的硬性要求。
要安装这两个软件包,请执行以下命令:
- 安装 PPA 以访问其软件包:
curl -sL https://deb.nodesource.com/setup_18.x -o nodesource_setup.sh
sudo bash nodesource_setup.sh
2. 安装 Node.js:
sudo apt install nodejs
3. 验证:
node -v
如果你看到类似 v18.15.0
(在撰写本文时的最新版本)的版本号,你就可以开始使用了。即使你安装了其他版本,也应该可以继续使用。
安装 CoC
安装 CoC 插件与其他插件一样简单。只需将以下行放入你的 init.vim
配置文件中:
Plug 'neoclide/coc.nvim', {'branch': 'release'} " Auto complete Python
重启编辑器并运行 :PlugInstall
Vim 命令。vim-plug
,我们正在使用的插件管理器,会为你安装 CoC。为完成安装过程,请再次重启编辑器。
要详细查看 Neovim 插件的安装过程,请参考本系列的第二篇博客。
我们准备安装所需的扩展,将编辑器转变为一个功能齐全的 IDE。运行 :CocInstall coc-pyright
Vim 命令以安装 Python 支持。安装完成后,退出编辑器的所有窗口,然后打开一个 Python 文件,如下所示:
nvim test.py
自动补全应该启用;然而,它并不是很实用。因此,让我们配置 CoC。在你的 init.vim
配置文件中添加以下几行:
" CoC customization
" Applying code actions to the selected code block
vmap <leader>a <Plug>(coc-codeaction-selected)<cr>
xmap <leader>a <Plug>(coc-codeaction-selected)<cr>
" Format action on <leader>
vmap <leader>f <Plug>(coc-format-selected)
xmap <leader>f <Plug>(coc-format-selected)
" GoTo code navigation
nmap <silent> gd <Plug>(coc-definition)
nmap <silent> gv :vsp<CR><Plug>(coc-definition)<C-W>L
nmap <silent> gy <Plug>(coc-type-definition)
nmap <silent> gi <Plug>(coc-implementation)
nmap <silent> gr <Plug>(coc-references)
" Symbol renaming
nmap <leader>rn <Plug>(coc-rename)
" Refactor
nmap <silent> <leader>re <Plug>(coc-codeaction-refactor)
xmap <silent> <leader>r <Plug>(coc-codeaction-refactor-selected)
nmap <silent> <leader>r <Plug>(coc-codeaction-refactor-selected)
" Confirm selection by pressing Tab
inoremap <expr> <Tab> pumvisible() ? coc#_select_confirm() : "<Tab>"
" Open autocomplete menu by pressing Ctrl + Space
inoremap <silent><expr> <c-space> coc#refresh()
" Open documentation by pressing K
nnoremap <silent> K :call ShowDocumentation()<cr>
function! ShowDocumentation()
if CocAction('hasProvider', 'hover')
call CocActionAsync('doHover')
else
call feedkeys('K', 'in')
endif
endfunction
我在每一行上留下了评论,以便你理解它的作用。现在 CoC 已配置好,你可以使用这些快捷键来充分利用它。
然而,事情并没有就此结束。正如我们所说,CoC 有多个扩展(如 coc-pyright
)可以安装。它们包括:
-
coc-git
-
coc-docker
-
coc-yaml
-
coc-json
-
coc-prettier
-
coc-pairs
-
coc-spell-checker
要安装它们中的任何一个,只需运行 :CocInstall
命令并提供扩展的名称。有关扩展的完整列表,请访问 coc-extensions
页面:github.com/topics/coc-extensions
最后,你可以像在 VS Code 中一样配置 CoC。运行 :CocConfig
命令,coc-settings.json
文件将弹出。我在其中设置了以下配置:
{
"python.venvPath": "/home/dimpo/.pyenv/versions/",
"python.linting.flake8Enabled": true,
"python.formatting.provider": "autopep8",
"pyright.inlayHints.variableTypes": false,
"pyright.inlayHints.functionReturnTypes": false,
"git.addGBlameToVirtualText": true,
}
每一行都是自解释的,但如果你想详细了解可以配置的内容,请参考每个 CoC 扩展的文档。
像往常一样,使用 docker commit
命令保存你的工作(如果你不知道这意味着什么,请参考本系列以前的文章)。你可以在 DockerHub 上找到我的镜像。
总结
在这一系列文章中,我们在 JupyterLab 环境中构建了一个完整的 Python IDE。我们从一个简单的 Neovim 安装开始,并将其配置为像今天最受欢迎的文本编辑器 VS Code 一样。
这篇文章结束了这个故事,但如果你想看到更多内容,如如何使用调试器或进一步配置特定插件,请留言,我会尽力而为!在此之前,祝编程愉快。
关于作者
我的名字是 Dimitris Poulopoulos,我是一名在 Arrikto 工作的机器学习工程师。我为包括欧盟委员会、欧洲统计局、国际货币基金组织、欧洲中央银行、经济合作与发展组织和宜家在内的主要客户设计和实施了 AI 和软件解决方案。
如果你对机器学习、深度学习、数据科学和数据操作等更多帖子感兴趣,请在 Medium、LinkedIn 或 Twitter 上关注我。
所表达的观点仅代表我个人,并不代表我的雇主的观点或意见。
平滑你的噪声数据的完美方法
使用 Whittaker-Eilers 方法进行极快且可靠的平滑和插值。
·
关注 发表在 Towards Data Science · 13 min read · 2023 年 10 月 25 日
–
现实世界的数据从来不会干净。无论你是在进行调查、测量降雨量还是接收来自太空的 GPS 信号,噪声数据总是存在。处理这些数据是数据科学家的主要工作内容。它并非全是光鲜的机器学习模型和人工智能——而是清理数据,以便尽可能提取有意义的信息。如果你当前看到的图表上有太多的曲折线而无法使用,那么,我有你需要的解决方案。
Whittaker-Eilers 平滑法
Whittaker-Eilers 平滑器 [1] 是我在从事地球观测工作时接触到的。对数千张高分辨率卫星图像进行像素级分析需要极其快速的算法来清理数据。某一天可能会有云层覆盖,一些像素可能被烟雾遮挡,或者传感器可能会出现伪影。列举的情况可以一直延续下去。你得到的是数 TB 的相当嘈杂、存在空缺的时间序列数据,需要进行平滑和插值——这正是 Whittaker 平滑器的强项。
平滑
让我们从一个图表开始,展示 1880 到 2022 年的 全球温度异常 [2]。橙色表示测量数据,绿色表示使用 Whittaker-Eilers 方法平滑的数据。
图 1) Whittaker-Eilers 平滑的全球温度异常
如你所见,原始的时间序列数据相当嘈杂。我们想提取的结论并不是关于逐年的波动,而是过去一个世纪的数据总体趋势。平滑数据提供了一种直接的方法来突出趋势,更好的是,它只需要 4 行代码就能运行。
from whittaker_eilers import WhittakerSmoother
temp_anom = [-0.17, -0.09, -0.11, -0.18] # and so on...
whittaker_smoother = WhittakerSmoother(
lmbda=20, order=2, data_length=len(temp_anom)
)
smoothed_temp_anom = whittaker_smoother.smooth(temp_anom)
如果你不喜欢 Python,这也只需要 4 行 Rust 代码。
use whittaker_eilers::WhittakerSmoother;
temp_anom = vec![-0.17, -0.09, -0.11, -0.18] // and so on...
let whittaker_smoother =
WhittakerSmoother::new(20, 2, temp_anom.len(), None, None)
.unwrap();
let smoothed_data = whittaker_smoother.smooth(&temp_anom).unwrap();
免责声明,这个包是我编写的 — Python: [*pip install whittaker-eilers*](https://pypi.org/project/whittaker-eilers/)
或 Rust: [*cargo add whittaker-eilers*](https://crates.io/crates/whittaker-eilers)
插值
现在让我们在数据中引入一些空缺。与其每年均匀收集数据,不如减少为每隔一年一次,并增加几个较长的空缺。
图 2) Whittaker-Eilers 平滑并插值的全球温度异常
Whittaker 方法处理起来非常轻松。你只需要给每个测量值分配一个权重。在插值时,存在的测量值权重为 1,而你希望获取值的测量用一个虚拟值,比如 -999,权重为 0。Whittaker 方法会处理剩下的部分。它也可以用来根据测量的不确定性进行相应的加权,而且只需几行 Python 代码即可完成。
from whittaker_eilers import Whittaker
temp_anom = [-0.17, -0.09, -0.11, -0.18, -0.3] # and so on...
weights = [1.0, 0.0, 1.0, 0.0, 1.0]
whittaker_smoother = WhittakerSmoother(
lmbda=150, order=2, data_length=len(temp_anom), weights=weights
)
smoothed_temp_anom = whittaker_smoother.smooth(temp_anom)
除了数据中的空缺,Whittaker 平滑器还可以轻松处理不均匀间隔的测量。如果你不在意插值,只需提供一个包含测量时间或地点的 x_input
,例如 [1880, 1885, 1886, 1900]
,就可以获得平滑后的值。
配置
通过一个单一的参数 λ(lambda),你可以持续控制数据的平滑度。lambda 越大,数据越平滑。
图 3) Whittaker-Eilers 平滑数据,λ 值不同
平滑器的阶数 d 也是可以控制的。阶数越高,Whittaker 在计算时间序列的平滑程度时会考虑更多的邻近元素。核心要点是,插值数据将是阶数为 2d 的多项式,而外推数据将是阶数为 d 的多项式;不过,我们稍后会详细探讨这些数学原理。
图 4)Whittaker-Eilers 平滑数据的多阶次
与其他方法的比较
为了展示 Whittaker-Eilers 方法的优越性,让我们将其与几种不同的技术进行比较:
-
高斯核平滑(也称为 RBF 核)
-
Savitzky-Golay 滤波器
-
局部回归 (LOWESS)
第一个是核平滑器,本质上是一种对邻近点进行加权平均的复杂方法。第二种是 Savitzky-Golay 滤波器,它与 Eilers 在 2003 年的原始论文相对立,并且使用非常广泛。它通过对相邻数据子集进行多项式的最小二乘拟合来平滑数据。最后,局部回归是 NASA 选择的平滑温度异常数据的方法。这种方法在连续的数据子集上执行迭代加权线性拟合,并根据残差重新加权点。
我选择不包括任何实时平滑的方法,因为这些方法在仅评估过去数据时会必然引入信号延迟。如果你对这种方法感兴趣,务必查看移动平均、指数平滑和卡尔曼滤波器。
平滑
那么让我们再看看温度异常时间序列,这一次使用额外的方法进行平滑。
图 5)不同平滑技术的比较
如你所见,变化不大——至少在你调节好参数后如此。让我们尝试一个例子,既可以揭示一些差异,也可以进行性能基准测试。我们将生成 4 个基本的正弦波,范围在 0 到 2π 之间,并添加一些高斯噪声。每个信号的样本大小将增加,从 40 开始,直到 10,000 结束。对于 Savitzky-Golay、高斯核和 LOWESS 方法来说,已经提出了挑战。它们期望一个窗口长度——每个连续子集中应包含的数据长度,以便进行拟合/平均。如果你在相同的时间段内增加测量数量,你需要确保你的窗口长度随着整体数据长度的变化而变化,以获得每个数据集的最佳平滑效果。下面,每种方法大致考虑了整体数据的十分之一作为每个窗口。
图 6) 不同测量采样频率下平滑的正弦波
然而,Whittaker 不需要如此动态的参数调整。你只需提供测量位置 x,设置一次λ,它就会每次都给出平滑的数据。这种测量位置的知识使其能够处理不等间距的数据。LOWESS 也可以在这种数据上运行,但另外两种方法要求数据间隔相等。总体而言,我发现 Whittaker 最易于使用。
Whittaker 平滑器的另一个好处是其对边界条件的适应性。在数据边缘,Savitzky-Golay 滤波器的行为会发生变化,因为它不能在没有填充或镜像数据的情况下拟合超出数据集末尾的多项式。另一种选择是让它拟合最终窗口中的多项式并使用这些结果。高斯核平滑器的情况更糟,因为没有未来的测量值进入平均值,它开始表现出来自以前值的大偏差,这在下面的图表中可以清楚地看到。
图 7) 不同测量采样频率下平滑正弦波的边界
插值
Whittaker 的一个关键特性是其内建的数据插值能力。但与其他方法相比,它的效果如何呢?这个问题并不那么明确。Savitzky-Golay 平滑可以进行插值,但仅限于数据间隙小于其窗口大小的情况,LOWESS 平滑也是如此。高斯核平滑根本没有插值的能力。解决此问题的传统方法是先对数据应用线性插值,然后再进行平滑。因此,我们将这种方法应用于其他三种技术,并将结果与 Whittaker 进行比较。
每种方法将与其自身的平滑基线进行比较,该基线取自本节开始时的图(图 5)。我删除了每个其他点并引入了两个大的间隙,创建了一个与文章开始处插值示例(图 2)中的数据集完全相同的数据集。基线和插值运行的参数保持不变。
图 8) 线性插值联合使用技术的比较
通过线性插值填补空白后,各方法的表现总体较好。通过计算均方根误差(RSME)来比较平滑后的数据(没有间隙)和有间隙的平滑数据,我们得到了以下结果。
-
线性插值 + Savitzky-Golay: 0.0245 °C
-
Whittaker : 0.0271 °C
-
线性插值 + 高斯核: 0.0275 °C
-
线性插值 + LOWESS: 0.0299 °C
线性插值结合 Savitzky-Golay 方法最接近原始平滑数据,其次是 Whittaker,两者之间差别不大!
我想简单提一下,我以这种方式执行了插值基准测试,对比了各自平滑的基线,以避免调整参数。我本可以使用添加噪声的正弦波,删除一些数据并尝试将其平滑回原始信号,但这会让我头疼于为每种方法找到最佳参数。
基准测试
让我们重新审视正弦波数据,以生成这些方法的速度基准。我选择了每种方法中最流行的 Python 实现。Savitzky-Golay 和高斯核滤波器使用SciPy
实现,LOWESS 使用statsmodels
实现,Whittaker 使用我基于 Rust 的 Python 包实现。下图显示了每种方法平滑正弦波所需的时间以及数据长度的变化。报告的时间是平滑每个数据集 50 次所需的总时间。
图 9) 平滑不同长度的时间序列所需时间 50 次
迄今为止最快的方法是 Whittaker。它可以在不到一秒的时间内平滑 50 个每个长度为 100,000 数据点的时间序列,比高斯滤波器快 10 倍,比 Savitzky-Golay 滤波器快 100 倍。最慢的是 LOWESS,尽管它被配置为不迭代地重新加权每个线性回归(这是一个昂贵的操作)。值得注意的是,通过调整窗口长度可以加快这些方法的速度,但这会牺牲数据的平滑度。这是 Whittaker 的一个非常好的特点——它的计算时间随着数据长度线性增加(O(n)),而且你无需担心窗口大小。此外,如果你的数据中有间隙,你将以无速度损失的方式进行插值,而其他方法则需要某种形式的预处理!
数学原理
现在我们已经覆盖了主要内容,让我们深入探讨 Whittaker-Eilers 平滑器背后的数学,并了解它为何对噪声数据是如此优雅的解决方案 [2] [3]。
想象一下你的噪声数据y。存在一些序列z,你认为它具有对你的y的最佳平滑度。z变得越平滑,它与原始数据y之间的残差越大。Whittaker-Eilers 方法在这些残差和数据的平滑度之间找到最佳平衡。残差计算为标准的平方差之和,
数据平滑度的度量可以使用相邻测量值之间的平方差之和来计算,
S 和 R 是我们需要平衡的两个属性。但我们还希望让用户控制合适的平衡点,我们通过引入λ来调整平滑度。
现在我们的目标变成了找到能最小化 Q 的序列 z,因为这是平滑度度量和残差都达到最小的地方。让我们展开方程 3 并尝试解出 z。
此时,理想的做法是用向量替换我们的求和,
然后我们可以使用一个巧妙的技巧,将 Δz 表示为矩阵和向量,
其中 m 是数据的长度。如果你将矩阵 D 与一个向量进行比较,你会发现它给出了相邻元素之间的差异——这正是我们想要的。现在我们面临一个最小二乘问题。为了找到 Q 的最小值,我们将其梯度设为 0,
其中 I 是单位矩阵(来自分解 z,一个向量)。我们知道 I、D、λ 和 y,所以我们剩下一个简单的线性方程,
这可以通过你喜欢的任何矩阵分解方法来解决,从而获得平滑的数据系列 z。
插值
上述解决方案仅适用于数据均匀分布且所有测量值都可用的情况。如果需要插值怎么办?你需要为每个测量值应用权重。
就像重新审视方程 1 并为每个残差应用权重,并将其表示为对角矩阵一样简单,
然后进行与之前相同的计算,
再次,这可以通过简单的矩阵分解来解决,返回平滑和插值的数据。所有需要做的事是在需要插值的地方将 y 填充为虚拟值,例如 -999,并将这些测量的权重设置为 0,然后观察奇迹发生。数据如何插值取决于滤波器的阶数。
滤波器阶数
Whittaker-Eilers 平滑器的阶数是在配置部分提到的内容。现在我们有了描述平滑器的数学框架,它可能更有意义。在创建 R(我们的平滑度度量)时,我们首先选择了“一级”差异。我们可以相当容易地计算二阶差异,这样不仅根据相邻数据点对平滑器进行惩罚,还可以根据一阶差异的变化来惩罚它,就像计算导数一样。
这可以扩展到三阶、四阶、五阶差分,等等。通常表示为d,实现起来并不太复杂,因为唯一改变的是矩阵D,如下所示:
使其在与z相乘时展开为方程 17。可以实现一个简单的函数来生成这个矩阵,给定一个通用的d。
稀疏矩阵
这可以并且已经按照 Eilers 的建议使用稀疏矩阵来实现。矩阵I和D非常稀疏,如果以稀疏矩阵的形式存储,将极大地节省内存和计算资源。上述所有数学计算都可以通过稀疏矩阵包轻松处理,包括 Cholesky 分解(及其他)。如果不使用稀疏矩阵实现,这个算法在处理较长时间序列时可能会非常缓慢,比我比较的其他方法慢得多。
总结与进一步阅读
这是一个了不起的算法,我简直不敢相信它没有被更广泛地应用。加权平滑和插值被封装成快速、高效的矩阵运算。有什么不好的吗?
我在whittaker-eilers包的 repo 中包含了我用于进行基准测试和插值测试的 Python 脚本。还有许多示例展示了如何在 Python 或 Rust 中入门,以及与 Eilers 原始 MATLAB 算法的测试[1]。但如果你不在意那么详细的信息,
Python: [*pip install whittaker-eilers*](https://pypi.org/project/whittaker-eilers/)
或 Rust: [*cargo add whittaker-eilers*](https://crates.io/crates/whittaker-eilers)
尽管这是一篇长文章,我仍然没有涵盖所有内容。Eilers 2003 年的论文还涉及了平滑不均匀间隔数据的数学原理以及如何使用交叉验证来找到最佳λ。如果你想深入了解算法背后的数学,我推荐查看这篇论文。我还建议阅读 Sophocles J. Orfanidis 的《应用最优信号处理》,它提供了关于所有信号处理内容的深入数学指南。感谢阅读!请务必查看我的个人网站上的这篇文章和其他文章。
我之后写了一篇跟进文章,解释了如何最佳地调整 Whittaker!
利用 Whittaker-Eilers 平滑和留一交叉验证最大化数据价值
[towardsdatascience.com
本文中的所有图像均由作者制作。
参考文献
[1] Paul H. C. Eilers, 完美平滑器, Analytical Chemistry 2003 75 (14), 3631-3636, DOI: 10.1021/ac034173t
[2] NASA/GISS,《全球温度》,NASA 戈达德太空研究所(GISS)。网址:climate.nasa.gov/vital-signs/global-temperature/
[3] 索福克勒斯·J·奥尔法尼迪,《应用最优信号处理》,罗格斯大学,网址:www.ece.rutgers.edu/~orfanidi/aosp
使我获得数据科学家职位的作品集
原文:
towardsdatascience.com/the-portfolio-that-got-me-a-data-scientist-job-513cc821bfe4
剧透:制作作品集意外地简单(而且免费)
·发表于Towards Data Science ·10 分钟阅读·2023 年 3 月 24 日
–
照片由John Schnobrich拍摄,发布在Unsplash
获得数据科学家职位很难。
现在已经不是 2015 年了:仅仅知道几个 pandas 函数并在简历上写上“Big Data”已经不够了。顶级职位的竞争非常激烈。在最近浏览 LinkedIn 职位板时,我发现很难找到申请人数少于 100 人的伦敦数据科学家职位。
好消息是,这种竞争并不是因为工作岗位不足。即使在 2023 年,数据科学仍然是一个快速增长的领域,美国劳工统计局估计数据科学家职位的数量将在 2021 年至 2031 年间增长 36%[1]。
因此,挑战不在于数据科学工作岗位的缺乏——实际上有很多!而在于大量的人试图进入这个行业,使得在众多候选人中脱颖而出并获得那个丰厚的首个职位变得越来越困难。
解决方案:建立数据科学作品集
在当前的就业市场中,我相信,区别于竞争对手的最佳方法之一就是通过建立个人数据科学项目的作品集来展示你的技能和经验。如果你没有大量的数据科学工作经验,这一点尤其重要。借用数据科学家 Will Stanton 的话[2]:
如果你没有数据科学家的经验,那么你绝对必须做独立项目。
在这篇文章中,我将带你了解我创建个人项目作品集的过程,并尝试为你提供一些灵感,帮助你自己创建一个作品集。
如果创建作品集的想法让你感到有些畏惧,不用担心:建立一个作品集其实出乎意料的简单(而且免费的),而且你不需要任何网页开发知识就能开始。
我的背景(也就是我创建我的作品集的原因)
我的数据科学之旅相对较晚才开始。我刚从学校/大学毕业时的第一份工作是在项目管理、销售和市场营销等领域。当我决定在 2020 年转向数据科学时,我几乎没有实际操作数据的经验,也没有计算机科学或统计学等领域的技术本科学位。
哇,这不是一个很强的起点。
在我之前的文章中,我分享了我如何通过参加一系列在线课程、攻读数据科学硕士学位和做几次数据科学实习来解决这个问题。在此过程中,我参与了许多不同的项目,并建立了广泛的数据相关技能。
然而,对于我申请的前几份工作,我发现我很难获得面试机会:感觉招聘人员甚至没有仔细查看我的申请,尽管我符合所有最低要求。那时,这非常令人沮丧。我相信我有很多可以提供的,但似乎对话还没开始就被中断了。
更糟糕的是,我了解到招聘人员通常会在过度申请的职位发布上花费仅 6 到 8 秒来查看简历[3]。这让我相信我需要一种更快的方式来传达我的技能并脱颖而出,特别是考虑到我作为数据科学家的经验并不丰富。
为了解决这个问题,我决定创建一个在线作品集,以一种视觉上引人注目且易于消化的格式展示我的经验。
我的数据科学作品集:mattschapman.github.io/
。图片来源:作者自制。
然后,一旦我建立了我的作品集,我在简历的顶部包含了一个链接,以便招聘人员可以轻松找到:
我的简历/履历的标题。我喜欢有一个可以点击的链接或二维码,直接带雇主到我的作品集,这样可以帮助我从人群中脱颖而出。图片来源:作者自制。
我是如何创建我的作品集的
在创建我的作品集时,我遵循了一些基本原则:
-
包括一些真正独特且有趣的项目
-
保持简短和简单
-
让它变得美观
-
不要在“网页开发”方面浪费时间
在本文的其余部分,我将简要讨论每一项原则,并解释它们如何帮助我创建一个成功的作品集。
包括一些真正独特且有趣的项目
这可能听起来显而易见,但这是一个如此重要的点,以至于值得重复:
不要包含别人也在包含的内容。
如果你作品集的全部目的就是区分自己,那为什么还要包括那些已经被别人做得死去活来的通用项目呢?这是 Jeremie Harris 在 数据科学家不被聘用的 4 种最快方式 [4] 中提到的一个观点,他说:
很难想到比在突出个人项目中展示你在琐碎的概念验证数据集上的工作更快地将你的简历扔进‘绝对不考虑’堆的方法。
如果有疑问,以下这些项目会让你受害多于受益:
使用 Titanic 数据集 进行生存分类。
使用 MNIST 数据集 进行手写数字分类。
使用 iris 数据集 进行花卉物种分类。
这意味着,如果你对计算机视觉非常感兴趣,就包括几个相关的项目。或者,如果你对时间序列问题情有独钟,就加上一些相关主题的项目。说实话,无论你对贝叶斯优化还是编织篮子感兴趣,关键是确保它是你真正感兴趣并认为能给别人带来价值的东西。就我个人而言,我非常关注自然语言处理和社会科学相关的问题,因此我尝试包括了这些领域的项目。但说实话,只要你能令人信服地展示为什么它们能带来价值,内容的具体细节其实并不是那么重要。
保持简洁明了
没有人愿意在无聊或冗长的网页上花时间,这在在线作品集中尤为真实。你提供的是一个快照,而不是写一篇博士论文。我诚心认为,如果你花太多时间编写冗长的项目描述和结果部分,你会浪费时间。
要理解这一点,我们必须回到作品集的初衷,其核心是一个简单的想法:
作品集的目的是为了让你有机会进入面试,而不是直接获得工作。
换句话说,你绝不会仅凭作品集就收到工作邀请。现实是,即使招聘人员喜欢你作品集的外观,你仍然需要通过他们的面试和评估。作品集的目的是提供你技能的快速快照,展示你能做什么。
在我的作品集中,例如,我提到了我之前参与的七个项目。在这七个项目中,六个项目的描述都不超过两句简短的话。
是的,你没看错:两句话。
这是我作品集中一个项目描述的例子:
我作品集中的一个示例项目。图片来源:作者自有。
如你所见,我没有浪费时间设置背景和解释结果的细微差别。当然,这些内容非常重要。但是如果潜在的雇主感兴趣,他们会在面试中询问你。在这个阶段,你只需要一个简短的描述,能够引起他们的兴趣并促使进一步讨论。如果你愿意,你甚至可以包括一个指向更多描述的链接;这也是我通过在项目描述下方包含代码链接所做的。但在你作品集的主页上,你只需要提供高层次的细节。其他的都是浪费时间。
让它变得漂亮
我的下一个建议是使你的作品集具有视觉吸引力。
设身处地为一个潜在的招聘者考虑:你已经阅读了数百(甚至数千)份平淡无奇的简历,甚至都快忘记绿色是什么样子,因为已经很久没见过不是黑白的东西了。然后,突然出现一份简历/作品集——简直是开门见山——其实挺好看的。
你不需要是火箭科学家(甚至不需要是数据科学家)就能明白招聘者会喜欢看你的作品集。
要让你的作品集更漂亮,我的首要建议是包含大量有趣的图表/图形。例如,这里是我在作品集中包含的另一个项目的描述,其中包含了一些关于 Covid-19 限制和封锁期间人口流动的图表:
我作品集中的一个项目描述。图片来源:作者自有。
图表向潜在雇主展示了我知道如何以视觉方式传达想法并制作漂亮的图像。这可能看起来微不足道,但实际上并非如此:作为数据科学家,你会经常需要向非技术性利益相关者传达你的发现。通过在作品集中包含一些图表,你向读者展示了你认识到沟通的重要性,并具备将其付诸实践的技能。换句话说,这是展示你技能并在竞争激烈的领域中脱颖而出的绝佳方式。
不要在网页开发上浪费不必要的时间
我的最后一个建议是不要花太多时间实际“构建”包含你作品集的网站。记住:你申请的是数据科学家的职位,而不是网页开发人员。你的评价标准不会是 HTML 和 CSS 技能,而是你的数据科学技能。
因此,除非你在网页开发方面有丰富的经验,否则我建议选择低代码/无代码选项来构建你的网站,例如使用像WordPress或Webflow这样的工具,它们允许你“拖放”作品集的各个部分,并从一系列预制模板中选择。
数据科学作品集的选项范围。图像来源:作者自制。
就我而言,我选择了一些不同的东西,并遵循了 这篇 由 Ivanna Kacewica 撰写的优秀指南,她使用了 Minimal Jekyll 主题 来为 GitHub Pages 构建网站的基本框架 [5]。我选择这个选项的主要原因是 (a) 它是免费的,(b) 它不要求你在 URL 中包含广告或“wordpress”字样。
Ivanna 方法的一个额外好处是它允许你在 GitHub Pages 上免费托管你的作品集,地址为username.github.io,其中username是你在 GitHub 上的用户名(或组织名称)。由于你的 GitHub 用户名包含在 URL 中,别人很容易找到你的 GitHub 账户,这对于你在 GitHub 上存储的其他项目非常有用。
我不会详细介绍如何使用 Ivanna 的主题(因为她解释得非常好),但这种方法的基本理念是通过 Jekyll 和 markdown 添加文本和图像。如果你不熟悉 Jekyll,它是一个简单的构建网站的框架,不需要编程或 HTML/网页开发知识。
使用 Jekyll,你建立网站的基本结构,然后用 markdown 填充文本和图像。如果你以前没用过 markdown,不用担心——markdown 只是一个非常简单的文本格式化方式。当我说“简单”时,我是指非常简单:例如,要写标题,你只需在文本前面加上一个井号。要制作项目符号,你只需在文本前面加上一个星号。就是这么简单。
如何使用 markdown 添加内容。图像来源:作者自制。
如果你想查看我如何自定义主题的详细信息,可以查看我 GitHub 上的 index.md
文件 repository。
结论
构建在线作品集是在人满为患的数据科学就业市场中脱颖而出的好方法。在这篇文章中,我提供了一些关于如何构建我的个人作品集的见解以及构建你自己作品集的最佳建议。
告诉我你的进展情况,如果你想获得关于你的作品集的反馈或建议,可以在评论中留下链接,我一定会查看。
如果你想获取更多关于如何构思作品集的建议,你可能会想看看我的其他文章:
忘掉泰坦尼克号和 MNIST:选择一个独特的项目来提升你的技能,并让你从人群中脱颖而出
towardsdatascience.com
哦,还有一件事 —
我开设了一个免费的通讯,AI in Five,每周分享 5 个要点,涵盖最新的 AI 新闻、编码技巧和数据科学家/分析师的职业故事。没有炒作,没有“数据是新石油”的废话,也没有来自埃隆的推文——只有实际的技巧和见解,帮助你在职业生涯中发展。如果这对你来说合适,点击这里订阅!
[## AI in Five | Matt Chapman | Substack
最新的新闻、职业故事和编码技巧,来自数据科学和 AI 领域,浓缩为 5 个要点……
来源
[1] 美国劳工统计局。职业前景手册:数据科学家。(2022 年 9 月 8 日)。www.bls.gov/ooh/math/data-scientists.htm
[2] Will Stanton. 创建一个出色的数据科学博客。(2015 年 7 月 15 日)。will-stanton.com/creating-a-great-data-science-resume/
[3] Standout CV. 招聘人员查看你的简历的时间有多长?(2022 年 11 月)。standout-cv.com/how-long-recruiters-spend-looking-at-cv#:~:text=Research%20shows%20that%20recruiters%20spend,15%20minutes%20reviewing%20a%20CV
。
[4] Jeremie Harris. 作为数据科学家被拒绝的 4 种最快方式。(2018 年 6 月 12 日)。towardsdatascience.com/the-4-fastest-ways-not-to-get-hired-as-a-data-scientist-565b42bd011e
[5] Ivanna Kacewica. 使用 Github Pages 在不到 10 分钟内建立你的作品集网站。(2019 年 3 月 1 日)。medium.com/@evanca/set-up-your-portfolio-website-in-less-than-10-minutes-with-github-pages-d0efa8ff56fd
使用蒙特卡罗模拟传播误差的力量和简单性
在数据分析和模型拟合中掌握不确定性,包括动手代码和示例
LucianoSphere (Luciano Abriata, PhD)
·发表在Towards Data Science ·15 分钟阅读·2023 年 7 月 21 日
–
关于蒙特卡罗方法的内容在Towards Data Science上有相当多的介绍,但关于它在误差传播中非常重要且有用的应用介绍却不多,除了Shuai Guo的精彩介绍以及其他几篇文章:
演示蒙特卡罗模拟
towardsdatascience.com
在这里,我想提出一些具体的数值应用和代码,供你实际尝试,并亲身感受蒙特卡罗方法在几乎任何类型的计算中如何极其有用且易于实现。
我将从一个非常简单的应用开始,传播减法操作中的误差,然后举例说明如何用基本相同的思路来传播几乎任何类型的数值例程中的误差,从简单的线性回归到复杂的拟合过程,这些过程在解析上很难处理。
蒙特卡罗模拟中的误差传播
误差传播是数据分析和科学计算中的一个基本概念。当你有带有不确定性的测量值时,对这些值进行数学运算会导致最终计算结果中的误差传播。对于简单的算术操作,误差传播可以使用公式进行解析。如果你对解析误差传播感兴趣,可以查看这个资源:
版权 © 2000 年 7 月 1. 系统误差与随机误差 2. 确定随机误差 (a) 仪器误差极限…
然而,对于涉及多个变量和非线性函数的复杂操作,或者对于数据拟合或神经网络执行等大型计算程序,解析误差传播可能很快变得不切实际甚至不可能。
蒙特卡洛模拟提供了一种替代的误差传播方法,特别是在解析解决方案具有挑战性或未知的情况下。蒙特卡洛方法涉及使用随机抽样来模拟大量场景,每次计算所需的量,然后分析结果的分布。这种统计方法提供了对最终结果不确定性的估计,使我们能够在复杂计算中传播误差。
为了演示蒙特卡洛模拟在误差传播中的应用,我们考虑一个非常简单的减法操作a - b,其中a和b都有一些相关的不确定性da和db。
(注:本文及所有后续示例均使用 Matlab 编写,但可以轻松移植到任何其他脚本语言。)
clear; close all; clc
a = 8; da = 2; % Example a = 8 with uncertainty 2
b = 4; db = 1; % Example b = 4 with uncertainty 1
c = a — b; dc = sqrt(da*da + db*db); % Analytical error propagation
disp(‘Expected’)
disp([num2str(c), ‘ +/- ‘, num2str(dc)]) % Print
fprintf(‘\n’)
% =============== Start here Monte Carlo simulation
c_mc = []; % Will hold all results
for i=1:100000
a_mc = normrnd(a, da); % Get a random number with a Gaussian distribution around a
b_mc = normrnd(b, db); % Same around b
c_mc = [c_mc; a_mc — b_mc]; % Run calculation and store
end
disp(‘From MonteCarlo with 10 iterations’)
disp([num2str(mean(c_mc(1:10))), ‘ +/- ‘, num2str(std(c_mc(1:10)))]) % Print result for 10 iterations
fprintf(‘\n’)
disp(‘From MonteCarlo with 100 iterations’)
disp([num2str(mean(c_mc(1:100))), ‘ +/- ‘, num2str(std(c_mc(1:100)))]) % Result for 100 iterations
fprintf(‘\n’)
disp(‘From MonteCarlo with 1000 iterations’)
disp([num2str(mean(c_mc(1:1000))), ‘ +/- ‘, num2str(std(c_mc(1:1000)))]) % Result for 1000 iterations
fprintf(‘\n’)
disp(‘From MonteCarlo with 10000 iterations’)
disp([num2str(mean(c_mc(1:10000))), ‘ +/- ‘, num2str(std(c_mc(1:10000)))]) % Result for 10000 iterations
fprintf(‘\n’)
disp('From MonteCarlo with 100000 iterations')
disp([num2str(mean(c_mc(1:100000))), ' +/- ', num2str(std(c_mc(1:100000)))]) % Result for 100000 iterations
fprintf('\n')
在这个示例中,不确定性既通过解析方法传播,也通过蒙特卡洛程序传播,然后对其结果进行比较。首先,将c上的不确定性解析传播为a和b的平方不确定性之和的平方根。
接下来是蒙特卡洛模拟。我们从以各自值为中心的高斯分布中随机抽取a和b,其标准差等于它们的不确定性。样本值存入临时变量a_mc和b_mc中,这些变量的差值存储在数组c_mc中。我们重复这个过程很多次,这里是 100,000 次,然后计算这个数组的平均值和标准差,以估计c的中心值及其相关的不确定性。这个示例程序报告了在 10 次、100 次、1000 次、10000 次和 100000 次迭代后的平均值和标准差,从中我们可以评估估计的收敛情况。
以下是一个示例运行的结果:
预期
4 +/- 2.2361
从蒙特卡洛模拟(10 次迭代)中得出的结果
3.2695 +/- 3.1021
从蒙特卡洛模拟(100 次迭代)中得出的结果
3.8985 +/- 2.5303
从蒙特卡洛模拟(1000 次迭代)中得出的结果
3.9573 +/- 2.2683
从蒙特卡洛模拟(10000 次迭代)中得出的结果
4.0182 +/- 2.2139
从蒙特卡洛模拟(100000 次迭代)中得出的结果
3.9984 +/- 2.2439
如你所见,计算结果和从解析传播中预期的误差(2.2361)经过蒙特卡洛过程的迭代后均合理收敛。特别注意误差如何下降并非常接近解析传播值。同时,即使是 1000 次迭代也已经产生了一个非常合理的误差估计。
请注意,如果你自己运行代码,当然会得到不同的数字,但趋势应该始终存在:随着迭代次数的增加,你的误差应该收敛到类似于示例中报告的值(以及解析推导值)。
再提供两个例子,现在估计通过数据拟合获得的参数的不确定性
在上一节中,我们探讨了蒙特卡洛模拟如何用于传播简单算术操作中的误差。现在,让我们深入研究更实际的例子,应用蒙特卡洛模拟来估计通过数据拟合获得的参数的不确定性。
在这个例子中,我们将考虑线性拟合过程。给定一组数据点 (x, y),我们希望找到形式为 y = ax + b 的最佳拟合直线。然而,现实世界的数据通常包含测量误差或其他不确定性,这可能会影响拟合参数 (a 和 b) 的准确性。为了考虑这些不确定性并获得更可靠的估计,我们可以使用蒙特卡洛模拟。
在这个例子中,我们生成一些样本数据,并添加以 0 为中心、宽度为 +/- 1 的高斯噪声以模拟测量误差。这成为了接下来的实验向量 y 值。
clear; close all; clc
a = 5; b = 4; % Slope and intercept
x = [1:10]; % x from 1 to 10
y = a * x + b + normrnd(0, 1, 1, 10); % Simulate data and add Gaussian noise centered at 0 and with a width of +/- 1
首先,我们使用 Matlab 的 fitlm 函数进行线性回归,该函数计算最佳拟合直线并提供解析推导的斜率和截距及其不确定性:
analyticfit = fitlm(x,y) % fitlm performs a linear fit and displays the results, which include
% the analytically-derived slope and intercept and their uncertainties
该过程的输出结果如下:
analyticfit =
Linear regression model:
Estimated Coefficients:
Estimate SE tStat pValue
________ _______ ______ __________
(Intercept) 3.6469 0.94449 3.8612 0.0048006
x1 5.0534 0.15222 33.198 7.3956e-10
Number of observations: 10, Error degrees of freedom: 8
Root Mean Squared Error: 1.38
R-squared: 0.993, Adjusted R-Squared: 0.992
F-statistic vs. constant model: 1.1e+03, p-value = 7.4e-10
因此,你可以看到在这一轮运行和生成器提供的随机数中,解析拟合得出的 斜率 = 5.0534 +/- 0.15222 和 截距 = 3.6469 +/- 0.94449
现在,让我们看看 Monte Carlo 程序的发现。只需按照我们传播减法操作不确定性的相同思路,我们必须取 x 和 y 向量,随机添加噪声(这里仅对y进行,因为我们假设x没有误差,但我们也可以完美地对其添加噪声!),然后运行拟合程序 10,000 次。这次我们使用polyfit,一个用于拟合多项式(这里是 1 次)的 Matlab 函数,它不返回拟合参数的不确定性——这些我们将通过 Monte Carlo 获得。
所以,就像在前面的示例中一样,我们添加噪声,运行线性拟合,并将参数存储在专门用于斜率(a_mc)和截距(b_mc)的向量中:
a_mc = []; % Vector with slopes from Monte Carlo run
b_mc = []; % Vector with intercepts from Monte Carlo run
for i=1:10000 % 10k iterations, and in each:
yy = normrnd(y,1); % Create temporal vector where noise is added to the experimental y values
P = polyfit(x,yy,1); % Fit that vector to x using linear regresion as a polynomial of degree 1
b_mc = [b_mc; P(1)]; % Store all a values
a_mc = [a_mc; P(2)]; %
end
fprintf(‘\n’)
disp(‘From Monte Carlo with 10 iterations’)
disp([‘Intercept = ‘, num2str(mean(a_mc(1:10))), ‘ +/- ‘, num2str(std(a_mc(1:10)))])
disp([‘Slope = ‘, num2str(mean(b_mc(1:10))), ‘ +/- ‘, num2str(std(b_mc(1:10)))])
fprintf(‘\n’)
disp(‘From Monte Carlo with 100 iterations’)
disp([‘Intercept = ‘, num2str(mean(a_mc(1:100))), ‘ +/- ‘, num2str(std(a_mc(1:100)))])
disp([‘Slope = ‘, num2str(mean(b_mc(1:100))), ‘ +/- ‘, num2str(std(b_mc(1:100)))])
fprintf(‘\n’)
disp(‘From Monte Carlo with 1000 iterations’)
disp([‘Intercept = ‘, num2str(mean(a_mc(1:1000))), ‘ +/- ‘, num2str(std(a_mc(1:1000)))])
disp([‘Slope = ‘, num2str(mean(b_mc(1:1000))), ‘ +/- ‘, num2str(std(b_mc(1:1000)))])
fprintf(‘\n’)
disp(‘From Monte Carlo with 10000 iterations’)
disp([‘Intercept = ‘, num2str(mean(a_mc(1:10000))), ‘ +/- ‘, num2str(std(a_mc(1:10000)))])
disp([‘Slope = ‘, num2str(mean(b_mc(1:10000))), ‘ +/- ‘, num2str(std(b_mc(1:10000)))])
fprintf(‘\n’)
这是我获得的一组示例结果:
来自 Monte Carlo,10 次迭代
截距 = 3.9435 +/- 0.97266
斜率 = 5.0174 +/- 0.16053
来自 Monte Carlo,100 次迭代
截距 = 3.581 +/- 0.69941
斜率 = 5.0716 +/- 0.11592
来自 Monte Carlo,1000 次迭代
截距 = 3.6461 +/- 0.69365
斜率 = 5.0551 +/- 0.11206
来自 Monte Carlo,10000 次迭代
截距 = 3.6482 +/- 0.68509
斜率 = 5.0533 +/- 0.10918
我们观察到,随着迭代次数的增加,Monte Carlo 对截距和斜率的估计逐渐接近通过fitlm获得的分析值。从 Monte Carlo 模拟中获得的不确定性(标准差)与fitlm提供的标准误差一致,确认了 Monte Carlo 方法在传播拟合参数误差方面的准确性。
然而,使用 Monte Carlo 方法的主要优点之一是我们可以选择输入的误差是什么,我们不仅可以考虑观测值y上的不确定性,还可以考虑自变量上的不确定性。方法的这些和其他特性允许更好、更受控且知识驱动的不确定性量化和传播。 (参见文末专门讨论这些和其他方法优点以及缺点的部分。)
一个更复杂的例子,来自实际应用
我们刚才用于线性程序的相同方法,实际上可以应用于几乎任何类型的计算,特别适用于复杂的拟合程序,它可以提供有关模型参数和预测可靠性的宝贵见解。
虽然上述示例主要是为了演示,但现在让我给你展示一个来自我自己研究的例子,在这个例子中,Monte Carlo 模拟是唯一实际可行的方法来将误差传播到拟合参数中。故事始于我需要拟合一个相当复杂的方程,其形式如下:
基础数据包括一个独立变量(我们的x),在这里是(1/tcp),假设没有噪声,以及一个被称为y的依赖变量,称为 R₂。我的实验数据集包括多个 x, y 值,报告了不同相关系统的结果,所有这些数据必须一起拟合,以生成 kₑₓ, pA, pB, R₂A⁰, R₂B⁰和 d𝓌的单一值。所有这些都在这篇文章中解释了,我使用了 Matlab 强大的nlinfit函数:
[## Matlab 的 nlinfit 进行高效非线性函数拟合
使用 Matlab 的 nlinfit 函数对复杂函数进行强大而多功能的拟合
pub.towardsai.net](https://pub.towardsai.net/efficient-nonlinear-function-fitting-with-matlabs-nlinfit-3a2948c0fde6?source=post_page-----9c8dcca9d90d--------------------------------)
在那篇文章中,我没有深入探讨误差估计,因为 Matlab 的nlinfit函数返回直接提供误差的矩阵……但是从数据中直接得到的,且没有办法告诉系统应该考虑输入数据中的具体误差。即,拟合参数的不确定性可以像这样一次性计算:
[b1,res,J] = nlinfit(x,y,’AuxiliaryFunctionForNlinfit’,guesses);
b1 = b1';
dof = length(x)-length(b1);
sfit1 = norm(res)/sqrt(dof);
sd1 = sfit1*sqrt(diag(inv(J’*J)));
现在,使用蒙特卡罗方法,我们可以进行更灵活和有意义的误差传播。例如,我们可以告诉程序,较大的y(R₂)具有更高的不确定性(实验固有的):
for mcIteration = 1:10000
mcIteration % Prints in what iteration the loop is
for j=1:numel(yexp) % Take experimental y values and add noise (more noise when higher)
yy(j) = normrnd(yexp(j), yexp(j) * 0.05); % In the example, the error is 5 % of the value
end
[b1,res,J]=nlinfit(x,yy,’AuxFunction’,beta0, opts);
predictedy = AuxFunction(b1,x);
mc_results = [mc_results; b1];
end
disp([‘-’])
disp([‘-’]) % Printe the results of the simulation
disp([‘================== MONTE CARLO SIMULATION’])
disp([‘R21 = ‘,num2str(mean(mc_results(:,1))), ‘ +/- ‘, num2str(std(mc_results(:,1)))])
disp([‘R22 = ‘,num2str(mean(mc_results(:,2))), ‘ +/- ‘, num2str(std(mc_results(:,2)))])
disp([‘kEX = ‘,num2str(mean(mc_results(:,3))), ‘ +/- ‘, num2str(std(mc_results(:,3)))])
disp([‘dw = ‘,num2str(mean(mc_results(:,4))), ‘ +/- ‘, num2str(std(mc_results(:,4)))])
disp([‘p1_0 = ‘,num2str(mean(mc_results(:,5))), ‘ +/- ‘, num2str(std(mc_results(:,5)))])
disp([‘p1_20 = ‘,num2str(mean(mc_results(:,6))), ‘ +/- ‘, num2str(std(mc_results(:,6)))])
disp([‘p1_40 = ‘,num2str(mean(mc_results(:,7))), ‘ +/- ‘, num2str(std(mc_results(:,7)))])
disp([‘p1_50 = ‘,num2str(mean(mc_results(:,8))), ‘ +/- ‘, num2str(std(mc_results(:,8)))])
你可以看到,这个过程非常简单,和之前完全一样:只需多次运行拟合,每次使用一些受噪声影响的输入数据,反映出其实际的不确定性。最后,通过计算均值和标准差来汇总结果,以估算最终值。
这是一个输出示例,比较了基于雅可比矩阵的误差估计(没有考虑每个测量的具体不确定性)和 10000 次蒙特卡罗模拟的输出:
================== 蒙特卡罗模拟(10,000 次运行)
R21 = 0.8034 +/- 0.11613
R22 = 32.5616 +/- 34.5635
kEX = 6983.9807 +/- 595.5142
dw = 519.5121 +/- 92.0618
p1_0 = 1.0092 +/- 0.0034707
p1_20 = 0.86244 +/- 0.026407
p1_40 = 0.69446 +/- 0.066818
p1_50 = 0.58735 +/- 0.096593
================== 使用雅可比矩阵误差的普通拟合
R21 = 0.31138 +/- 0.25225
R22 = 33.0533 +/- 4.7355
kEX = 6965.3075 +/- 273.3825
dw = 513.7739 +/- 21.398
p1–0 = 1.0013 +/- 0.0049281
p1–20 = 0.86324 +/- 0.015702
p1–40 = 0.71381 +/- 0.03781
p1–50 = 0.62331 +/- 0.053877
你可以看到,收敛值和误差相当相似,但使用基于雅可比矩阵的常规误差传播时,误差系统地较小。了解这个问题背后的物理学,蒙特卡罗传播的误差更有意义,特别是 R21 的相对误差较低,而 R22 的相对误差较高。
(有关问题的完整解释、物理学、实验和计算,请关注我的 Medium 账户,我将很快撰写一篇专门的文章)。
附加内容:使用蒙特卡罗方法估计数值过程输出的不确定性的一些实际优缺点
以下优点使得蒙特卡罗模拟在各种数据分析和建模场景中成为一种强大且灵活的方法:
-
处理复杂误差分布: 蒙特卡罗方法的一个重要优点是其处理复杂误差分布的能力。在许多实际情况中,观察数据(y 值)和输入变量(x 值)的不确定性可能不遵循简单的高斯分布,这在大多数分析误差传播推导中被假设。通过使用蒙特卡罗模拟,我们可以轻松定制误差分布,以更好地代表我们数据的特定特征。这种灵活性使我们能够更准确地建模不确定性,从而改进参数估计和提高预测的可靠性。
-
误差建模的灵活性:
与一些假设特定误差模型的直接分析方法不同,蒙特卡罗模拟允许我们纳入各种误差结构。例如,我们可以包括异方差误差,其中误差的变异性随着不同的 x 值变化。此外,我们还可以考虑相关误差、非高斯分布以及其他可能在分析上难以处理的复杂情况。这种适应性使得蒙特卡罗方法能够处理更广泛的实际数据场景。
-
处理非线性模型:
在许多情况下,数值模型可能涉及变量之间复杂的非线性关系。分析误差传播方法通常仅限于线性模型,使其在非线性情境下不够适用。相比之下,蒙特卡罗模拟天生适用于非线性模型,因为它们依赖于数值采样而不是显式的分析推导。这种能力使我们能够在更复杂和现实的模型中评估参数的不确定性。
-
对假设的鲁棒性:
分析误差传播依赖于高斯误差和线性假设,这在实际中可能并不总是成立。当这些假设被违反时,分析结果的准确性可能会受到影响。而蒙特卡罗模拟则不依赖于这些假设,可以处理更广泛的数据分布和模型复杂性。这种鲁棒性确保了即使在传统方法可能失效的情况下,也能充分捕捉不确定性。
-
纳入先验信息:
在贝叶斯蒙特卡罗方法中,我们可以在拟合数据之前整合关于参数的先验知识或信念。当有来自之前研究或专家意见的先验信息时,这一特性特别有用。通过将先验分布纳入模拟,我们可以有效更新我们的参数估计和不确定性评估,从而得到更具信息性和准确的建模结果。
-
直观量化不确定性:
蒙特卡罗模拟生成参数估计的分布,使我们能够直观地理解可能值的范围和相关的不确定性。这种方法通常比从分析方法中获得的单点估计更具信息性,因为单点估计可能无法传达不确定性的全貌。通过可视化不确定性分布,决策者可以基于不同结果的概率做出明智的选择。
然而,蒙特卡罗方法并非万无一失,因此了解这些缺点至关重要:
-
计算强度:
蒙特卡罗模拟通常涉及大量迭代以获得准确结果。随着迭代次数的增加,计算时间和所需资源也显著增加。对于复杂模型和大量数据点,蒙特卡罗模拟可能变得计算上昂贵且耗时,这可能会阻碍其在某些情况下的实用性。
-
收敛性问题:
对采样分布的敏感性:
-
在某些情况下,蒙特卡罗模拟可能会遇到收敛性问题,尤其是在处理高维参数空间或复杂模型时。确保收敛并从参数空间中获得代表性样本可能具有挑战性,结果的质量取决于所使用的采样策略。收敛不足可能导致偏差或不可靠的不确定性估计。
蒙特卡罗结果的准确性和可靠性在很大程度上依赖于用于生成模拟随机样本的采样分布的选择。如果所选分布不能充分代表数据中的真实不确定性,它可能导致不准确或误导的不确定性估计。正确选择采样分布需要仔细考虑和对数据特征的了解。
-
模型选择与验证:
蒙特卡罗模拟不能提供模型选择或验证的指导。虽然它们可以估计拟合参数的 uncertainties,但并不解决模型适用性或拟合优度相关的问题。验证所选择的模型并评估其对给定数据的适当性仍然是需要与蒙特卡罗不确定性分析一起认真考虑的独立任务。
-
先验分布的主观性:
在贝叶斯蒙特卡洛方法中,纳入先验信息涉及为模型参数指定先验分布。这些先验分布的选择可能引入主观性,因为不同的研究人员或分析师可能有不同的先验信念。这种主观性可能会影响最终的参数估计和不确定性量化,导致结果可能存在分歧。
-
解释复杂性:
解读和传达蒙特卡洛模拟结果可能比传统的分析方法更具挑战性。不确定性分布可能不如点估计易于解释,分布中的大量信息可能会使非专家感到不知所措,而简单的误差传播可能就足够了。
-
过拟合风险:
使用蒙特卡洛模拟进行不确定性估计时,有过拟合数据的风险。运行大量迭代并拟合多个模型可能会无意中找到在特定模拟数据集上表现良好的模型,但对新数据泛化效果不佳。需要仔细的验证和交叉验证策略来降低这一风险。当然,还需要在解决问题的背景下判断结果的含义。
进一步阅读
想了解更多关于蒙特卡洛方法的内容,请查看这三篇发表在TDS Editors上的文章,这些文章介绍了除帮助传播不确定性以外的引言和应用:
A Zero-Math Introduction to Markov Chain Monte Carlo Methods [## 零数学介绍:马尔科夫链蒙特卡洛方法
对许多人来说,贝叶斯统计学充其量是巫术,最糟糕的是完全主观的无稽之谈。在这些…
An Overview of Monte Carlo Methods [## 蒙特卡洛方法概述
蒙特卡洛(MC)方法是计算算法的一部分,使用重复随机抽样的过程来…
An Overview of Monte Carlo Methods Monte Carlo Simulation [## 赌场总是赢钱:蒙特卡洛模拟
赌场如何赚钱?诀窍很简单——你玩得时间越长,亏钱的概率就越高。让我们…
www.lucianoabriata.com 我撰写和拍摄关于我广泛兴趣范围内的所有事物:自然、科学、技术、编程等。
在这里打赏我 或 成为 Medium 会员 以访问所有故事(我会获得少量收入,而您无需支付任何费用)。 订阅以获取我的新故事 通过电子邮件。 关于小工作咨询 请查看我的 服务页面。您可以 在这里联系我。
贝叶斯因果推断的力量:揭示数据集中隐藏因果关系的库的比较分析。
使用最适合的贝叶斯因果推断库揭示数据集中的隐藏因果变量:与五个流行库的实操示例进行比较。
·发表于 Towards Data Science ·阅读时间 20 分钟·2023 年 5 月 22 日
–
图片由 Alexander Schimmeck 提供,来源于 Unsplash
理解系统或过程中的变量因果效应非常有价值。现在有许多 Python 库可以帮助确定因果关系。我将比较五个流行的因果推断库在功能、易用性和灵活性方面的表现。 每个库都配有实际操作示例。包含的库有 Bnlearn, Pgmpy, CausalNex, DoWhy,和 CausalImpact。通过阅读本文,你将更好地理解这五个因果推断库,并确定哪个最适合你的使用场景。
背景
因果推断是确定过程或系统中变量之间的因果关系。一般来说,我们可以将变量分为两类;驱动变量和乘客变量。驱动变量是那些直接影响结果或因变量的变量,而乘客变量则是那些没有直接影响但与结果变量相关联的变量。
识别驱动变量在预测维护或安全领域等项目中可能至关重要。驱动变量有助于解释预测变量和结果变量之间的因果关系。相反,乘客变量对结果变量没有直接影响。然而,它们仍然可以有用,因为它们可以提供额外的变化,从而理解数据收集的背景。例如,如果我们发现引擎故障与地点有强相关性,我们可能会怀疑有一个潜在的驱动变量导致特定地点的故障。 因果推断有助于通过识别要操作的变量和要监控的变量来做出更好的决策。
在因果推断中,我们不仅寻求确定一个事件是否会发生,还要理解它发生的机制。
然而,因果推断分析是一项具有挑战性的任务,因为它需要分离多个变量的影响,考虑混杂变量,并处理不确定性。幸运的是,Python 有几个库可以帮助数据科学家进行因果推断。在这篇文章中,我将比较 Python 中五个最受欢迎的因果推断库:Bnlearn、 Pgmpy、 DoWhy、 CausalNex*** 和 CausalImpact。
如果你觉得这篇关于贝叶斯因果学习的文章有帮助,欢迎 关注我 ,因为我会写更多关于贝叶斯统计的内容。如果你考虑订阅 Medium 会员,可以通过我的 推荐链接 来稍微支持我的工作。这和买一杯咖啡的价格一样,你可以每月无限阅读所有文章。
五个因果推断库的比较。
让我们逐一检查这五个软件包,并从解决问题的能力角度审视它们的功能。为了最一致地比较这些库,我将尽可能使用相同的数据集,即普查收入数据 集。这个多变量数据集包含 14 个节点,具有 48842 个实例(或样本)[1]。每个节点包含两个或更多可能的状态,可以用来估计拥有研究生学位是否对年收入超过 50k 很重要。我们可以使用Bnlearn 库加载数据集。
# Installation
pip install bnlearn
# Import libraries
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import bnlearn as bn
# Import data set and drop continous and sensitive features
df = bn.import_example(data='census_income')
drop_cols = ['age', 'fnlwgt', 'education-num', 'capital-gain', 'capital-loss', 'hours-per-week', 'race', 'sex']
df.drop(labels=drop_cols, axis=1, inplace=True)
df.head()
# workclass education ... native-country salary
#0 State-gov Bachelors ... United-States <=50K
#1 Self-emp-not-inc Bachelors ... United-States <=50K
#2 Private HS-grad ... United-States <=50K
#3 Private 11th ... United-States <=50K
#4 Private Bachelors ... Cuba <=50K
#
#[5 rows x 7 columns]
库 1:Python 的 Bnlearn。
Bnlearn 是一个适用于创建和分析贝叶斯网络的 Python 包,适用于离散、混合和连续数据集[2, 3]。它设计为易于使用,并包含了最受欢迎的贝叶斯管道,用于因果学习,包括结构学习、参数学习和推断。可以通过简单指定参数进行初始化来使用一系列统计测试。Bnlearn 还包含各种辅助函数,用于转换数据集、推导(整个)图的拓扑排序、比较两个图以及制作各种有洞察力的图表等。有关使用 bnlearn 进行结构学习的更多详细信息请参见这里:
## 使用贝叶斯结构学习在 Python 中检测因果关系的逐步指南。
入门指南,用于有效确定变量间的因果关系。
[towardsdatascience.com
Bnlearn 的一个重要功能是它可以仅基于数据集学习因果结构。为此任务实现了六种算法;hillclimbsearch, exhaustivesearch, constraintsearch, chow-liu, naivebayes 和 TAN
,并且可以与评分类型 BIC, K2, BDEU
结合使用。一些方法需要设置根节点,例如增强型朴素贝叶斯(TAN
),如果你知道结果(或目标)值的话,推荐使用。这也会显著降低计算负担,并且在特征很多的情况下推荐使用。此外,通过 independence test
,可以轻松地从模型中修剪掉虚假边缘。在下面的示例中,我将使用 hillclimbsearch
方法和评分类型 BIC
进行结构学习。在这个示例中,我们不会定义目标值,而是让Bnlearn 完全基于数据自身决定整个因果结构。
# Load library
import bnlearn as bn
# Structure learning
model = bn.structure_learning.fit(df, methodtype='hillclimbsearch', scoretype='bic')
# Test edges significance and remove.
model = bn.independence_test(model, df, test="chi_square", alpha=0.05, prune=True)
# Make plot
G = bn.plot(model, interactive=False)
# Make plot interactive
G = bn.plot(model, interactive=True)
# Show edges
print(model['model_edges'])
# [('education', 'salary'),
# ('marital-status', 'relationship'),
# ('occupation', 'workclass'),
# ('occupation', 'education'),
# ('relationship', 'salary'),
# ('relationship', 'occupation')]
要确定有向无环图(DAG),我们需要指定输入数据框,如上面的代码部分所示。拟合模型后,结果存储在 model
字典中,可以用于进一步分析。因果结构的交互式图表如图 1 所示。
图 1. 使用 Bnlearn 进行结构学习的人口普查收入数据 集的交互式图表。如果学习到 CPD,工具提示将描述估计的 CPD(图片由作者提供)。
通过学习到的有向无环图(图 1),我们可以估计条件概率分布(CPD,见下方代码部分),并使用do-calculus 进行推断。换句话说,我们可以开始向数据提出问题。
# Learn the CPDs using the estimated edges.
# Note that we can also customize the edges or manually provide a DAG.
# model = bn.make_DAG(model['model_edges'])
# Learn the CPD by providing the model and dataframe
model = bn.parameter_learning.fit(model, df)
# Print the CPD
CPD = bn.print_CPD(model)
问题 1: 在教育为博士的情况下,薪水>50k 的概率是多少:
P(salary | education=Doctorate)
从直观上看,我们可能会期望高概率,因为教育是“博士”。让我们从贝叶斯模型中找出后验概率。在下面的代码部分,我们得出了P=0.7093
的概率。这证实了当教育为博士时,与没有博士教育相比,薪水>50K 的概率更高。
# Start making inferences
query = bn.inference.fit(model, variables=['salary'], evidence={'education':'Doctorate'})
print(query)
+---------------+---------------+
| salary | phi(salary) |
+===============+===============+
| salary(<=50K) | 0.2907 |
+---------------+---------------+
| salary(>50K) | 0.7093 |
+---------------+---------------+
现在让我们询问一下较低教育是否也会导致薪水>50K 的概率较低。我们可以轻松地将教育更改为HS-grad
,然后再次提出这个问题。
问题 2: 在教育为 HS-grad 的情况下,薪水>50k 的概率是多少:
*P(salary | education=*HS-grad*)*
这导致了P=0.1615
的后验概率。根据这个数据集,学习对获得更高薪水非常有利。然而,请注意,我们没有使用任何其他约束作为可能影响结果的证据。
query = bn.inference.fit(model, variables=['salary'], evidence={'education':'HS-grad'})
print(query)
+---------------+---------------+
| salary | phi(salary) |
+===============+===============+
| salary(<=50K) | 0.8385 |
+---------------+---------------+
| salary(>50K) | 0.1615 |
+---------------+---------------+
直到这一部分,我们使用了单一变量,但 DAG 中的所有变量都可以用作证据。让我们做另一个更复杂的查询。
问题 3: 在教育为博士且婚姻状态为未婚的情况下,属于某个工作类别的概率是多少:
P(workclass| education=Doctorate, marital-status=never-married)
。
在下面的代码部分可以看到,这返回了每个工作类别的概率,其中工作类别为private
的概率最高:P=0.5639
。
# Start making inferences
query = bn.inference.fit(model, variables=['workclass'], evidence={'education':'Doctorate', 'marital-status':'Never-married'})
print(query)
+-----------------------------+------------------+
| workclass | phi(workclass) |
+=============================+==================+
| workclass(?) | 0.0420 |
+-----------------------------+------------------+
| workclass(Federal-gov) | 0.0420 |
+-----------------------------+------------------+
| workclass(Local-gov) | 0.1326 |
+-----------------------------+------------------+
| workclass(Never-worked) | 0.0034 |
+-----------------------------+------------------+
| workclass(Private) | 0.5639 |
+-----------------------------+------------------+
| workclass(Self-emp-inc) | 0.0448 |
+-----------------------------+------------------+
| workclass(Self-emp-not-inc) | 0.0868 |
+-----------------------------+------------------+
| workclass(State-gov) | 0.0810 |
+-----------------------------+------------------+
| workclass(Without-pay) | 0.0034 |
+-----------------------------+------------------+
总结
-
输入数据:输入数据可以是离散的、连续的或混合的数据集。
-
优势: 包含最受欢迎的贝叶斯管道,用于结构学习、参数学习和使用do-calculus进行推断。可以轻松创建图表并探索 CPD。适合初学者和不想自己构建管道的专家。
库 2:Pgmpy。
Pgmpy 是一个提供概率图模型工具的库。它包含了各种统计方法的实现,用于结构学习、参数估计、近似(基于采样)和精确推断。一个优点是核心函数是低级统计函数,这使得构建自己的因果块更加灵活。虽然这很棒,但它需要对贝叶斯建模有很好的了解,因此相比于像Bnlearn这样的库,对于刚开始学习贝叶斯建模的人来说可能更困难。
在功能方面,可以得到与 Bnlearn 相同的结果,因为 Bnlearn 也利用了 Pgmpy 的一些核心功能。然而,在 pgmpy 中,你确实需要自己构建整个管道,包括数据所需的转换步骤,以及结果的收集和绘图。在下面的代码部分,我将简要演示如何使用 hillclimbsearch
估计器和评分方法 BIC
学习因果结构。请注意,不同的估计器可能在管道中需要非常不同的步骤。
# Install pgmpy
pip install pgmpy
# Import functions from pgmpy
from pgmpy.estimators import HillClimbSearch, BicScore, BayesianEstimator
from pgmpy.models import BayesianNetwork, NaiveBayes
from pgmpy.inference import VariableElimination
# Import data set and drop continous and sensitive features
df = bn.import_example(data='census_income')
drop_cols = ['age', 'fnlwgt', 'education-num', 'capital-gain', 'capital-loss', 'hours-per-week', 'race', 'sex']
df.drop(labels=drop_cols, axis=1, inplace=True)
# Create estimator
est = HillClimbSearch(df)
# Create scoring method
scoring_method = BicScore(df)
# Create the model and print the edges
model = est.estimate(scoring_method=scoring_method)
# Show edges
print(model.edges())
# [('education', 'salary'),
# ('marital-status', 'relationship'),
# ('occupation', 'workclass'),
# ('occupation', 'education'),
# ('relationship', 'salary'),
# ('relationship', 'occupation')]
使用 do-calculus 对数据进行推断时,我们首先需要估计 CPDs,如下方代码部分所示。
vec = {
'source': ['education', 'marital-status', 'occupation', 'relationship', 'relationship', 'salary'],
'target': ['occupation', 'relationship', 'workclass', 'education', 'salary', 'education'],
'weight': [True, True, True, True, True, True]
}
vec = pd.DataFrame(vec)
# Create Bayesian model
bayesianmodel = BayesianNetwork(vec)
# Fit the model
bayesianmodel.fit(df, estimator=BayesianEstimator, prior_type='bdeu', equivalent_sample_size=1000)
# Create model for variable elimination
model_infer = VariableElimination(bayesianmodel)
# Query
query = model_infer.query(variables=['salary'], evidence={'education':'Doctorate'})
print(query)
总结
-
输入数据:输入数据集必须是离散的。
-
优点: 包含可以用来创建自己贝叶斯管道的基础构建模块。
-
缺点: 需要对贝叶斯建模有良好的了解。
库 3:CausalNex。
CausalNex 是一个用于从数据中学习因果模型、识别因果路径和估计因果效应的 Python 库 [5]。Causalnex 仅支持离散分布,而连续特征或具有大量类别的特征应在拟合贝叶斯网络之前进行离散化。描述中提到:“如果使用了许多变量,模型通常会拟合较差”。然而,提供了辅助功能以减少分类特征的基数。让我们使用 Census Income data set 来检查 Causalnex 的可能性。首先,我们需要将数据转换为数值,因为这正是底层 NOTEARS [7] 算法所期望的。我们可以通过标签编码非数值变量来实现(见代码部分)。
# Installation
pip install causalnex
# Load libraries
from causalnex.structure.notears import from_pandas
from causalnex.network import BayesianNetwork
import networkx as nx
from sklearn.preprocessing import LabelEncoder
le = LabelEncoder()
# Import data set and drop continous and sensitive features
df = bn.import_example(data='census_income')
drop_cols = ['age', 'fnlwgt', 'education-num', 'capital-gain', 'capital-loss', 'hours-per-week', 'race', 'sex']
df.drop(labels=drop_cols, axis=1, inplace=True)
# Next, we want to make our data numeric, since this is what the NOTEARS expect.
df_num = df.copy()
for col in df_num.columns:
df_num[col] = le.fit_transform(df_num[col])
我们现在可以应用 NOTEARS 算法来学习结构并使用绘图功能可视化因果模型。我们需要应用阈值处理以移除较弱的边缘,并防止图形完全连通。此外,为了避免伪边缘,可以对边缘施加约束。
# Structure learning
sm = from_pandas(df_num)
# Thresholding
sm.remove_edges_below_threshold(0.8)
# Use positions from Bnlearn
pos=G['pos']
# Make plot
plt.figure(figsize=(15,10));
edge_width = [ d['weight']*0.3 for (u,v,d) in sm.edges(data=True)]
nx.draw_networkx_labels(sm, pos, font_family="Yu Gothic", font_weight="bold")
nx.draw_networkx(sm, node_size=400, arrowsize=20, alpha=0.6, edge_color='b', width=edge_width)
# If required, remove spurious edges and relearn structure.
sm = from_pandas(df_num, tabu_edges=[("relationship", "native-country")], w_threshold=0.8)
图 3. 使用 Causalnex 对 Census Income data set 进行结构学习。没有节点的标签包含较弱的边缘,这些边缘通过阈值处理被移除。(图像作者提供。)
使用推导出的结构,我们可以学习条件概率分布(CPDs)并开始进行查询。然而,还需要进行一些额外的步骤,即减少分类特征的基数,并离散化数值特征。请注意,我们还需要以字典形式指定节点状态,以防止手动将数值转换为分类标签。尽管每个步骤都很重要,但为了简化起见,我将跳过这些步骤以避免代码行数过多。有关详细信息,建议阅读文档手册这里。它还演示了如何使用一些辅助方法来简化离散化过程。
# Step 1: Create a new instance of BayesianNetwork
bn = BayesianNetwork(structure_model)
# Step 2: Reduce the cardinality of categorical features
# Use domain knowledge or other manners to remove redundant features.
# Step 3: Create Labels for Numeric Features
# Create a dictionary that describes the relation between numeric value and label.
# Step 4: Specify all of the states that each node can take
bn = bn.fit_node_states(df)
# Step 5: Fit Conditional Probability Distributions
bn = bn.fit_cpds(df, method="BayesianEstimator", bayes_prior="K2")
# Return CPD for education
result = bn.cpds["education"]
# Extract any information and probabilities related to education.
print(result)
# marital-status Divorced ... Widowed
# salary <=50K ... >50K
# workclass ? Federal-gov ... State-gov Without-pay
# education ...
# 10th 0.077320 0.019231 ... 0.058824 0.0625
# 11th 0.061856 0.012821 ... 0.117647 0.0625
# 12th 0.020619 0.006410 ... 0.058824 0.0625
# 1st-4th 0.015464 0.006410 ... 0.058824 0.0625
# 5th-6th 0.010309 0.006410 ... 0.058824 0.0625
# 7th-8th 0.056701 0.006410 ... 0.058824 0.0625
# 9th 0.067010 0.006410 ... 0.058824 0.0625
# Assoc-acdm 0.025773 0.057692 ... 0.058824 0.0625
# Assoc-voc 0.046392 0.051282 ... 0.058824 0.0625
# Bachelors 0.097938 0.128205 ... 0.058824 0.0625
# Doctorate 0.005155 0.006410 ... 0.058824 0.0625
# HS-grad 0.278351 0.333333 ... 0.058824 0.0625
# Masters 0.015464 0.032051 ... 0.058824 0.0625
# Preschool 0.005155 0.006410 ... 0.058824 0.0625
# Prof-school 0.015464 0.006410 ... 0.058824 0.0625
# Some-college 0.201031 0.314103 ... 0.058824 0.0625
# [16 rows x 126 columns]
摘要
-
输入数据:输入数据集必须是离散的。连续数据不支持。
-
优点:可以学习因果网络并对数据进行推断。
-
缺点:需要大量预处理步骤,但提供了辅助功能以减少分类特征的基数。
库 4:DoWhy。
DoWhy是一个用于因果推断的 Python 库,支持因果假设的明确建模和测试[2, 3]。与 Bnlearn 和 Pgmpy 相比,在 DoWhy 库中,必须定义结果变量和处理变量。 处理变量是你想研究因果效应的变量,而结果变量是你想测量其效果的变量。此外,建议提供一个 DAG,即节点之间的因果关系。要为你的数据集创建因果图,你需要领域知识,或者可以将每个变量连接到目标和处理变量。注意,如果没有提供 DAG,后者会自动执行。对于这个数据集,我不包括已知的结构,而是让DoWhy在所有变量之间创建到结果和处理变量的边。第一个代码部分显示了所需的预处理步骤,第二个代码部分则是因果模型的创建。
# Installation
pip install dowhy
# Import libraries
from dowhy import CausalModel
import dowhy.datasets
from sklearn.preprocessing import LabelEncoder
# Import data set and drop continous and sensitive features
df = bn.import_example(data='census_income')
drop_cols = ['age', 'fnlwgt', 'education-num', 'capital-gain', 'capital-loss', 'hours-per-week', 'race', 'sex']
df.drop(labels=drop_cols, axis=1, inplace=True)
# Treatment variable must be binary
df['education'] = df['education']=='Doctorate'
# Next, we need to make our data numeric.
df_num = df.copy()
for col in df_num.columns:
df_num[col] = le.fit_transform(df_num[col])
# Specify the treatment, outcome, and potential confounding variables
treatment = "education"
outcome = "salary"
# Step 1\. Create a Causal Graph
model= CausalModel(
data=df_num,
treatment=treatment,
outcome=outcome,
common_causes=list(df.columns[~np.isin(df.columns, [treatment, outcome])]),
graph_builder='ges',
alpha=0.05,
)
# Display the model
model.view_model()
图 2。DoWhy 创建了这个 DAG 来建模人口普查收入数据 集,其中结果变量为“salary”,处理变量为“Education”(图像由作者提供)。
正如你在上面的代码部分中可能已经注意到的,处理变量必须是二元的(设置为 Doctorate),并且所有分类变量需要编码为数值。在下面的代码部分,我们将使用因果图的属性来识别因果效应。结果可能不会令人惊讶。它表明,拥有Doctorate
增加了获得良好salary
的机会。
# Step 2: Identify causal effect and return target estimands
identified_estimand = model.identify_effect(proceed_when_unidentifiable=True)
# Results
print(identified_estimand)
# Estimand type: EstimandType.NONPARAMETRIC_ATE
# ### Estimand : 1
# Estimand name: backdoor
# Estimand expression:
# d
# ────────────(E[salary|workclass,marital-status,native-country,relationship,occupation])
# d[education]
#
# Estimand assumption 1, Unconfoundedness: If U→{education} and U→salary then P(salary|education,workclass,marital-status,native-country,relationship,occupation,U) = P(salary|education,workclass,marital-status,native-country,relationship,occupation)
#
# ### Estimand : 2
# Estimand name: iv
# No such variable(s) found!
#
# ### Estimand : 3
# Estimand name: frontdoor
# No such variable(s) found!
# Step 3: Estimate the target estimand using a statistical method.
estimate = model.estimate_effect(identified_estimand, method_name="backdoor.propensity_score_stratification")
# Results
print(estimate)
#*** Causal Estimate ***
#
## Identified estimand
# Estimand type: EstimandType.NONPARAMETRIC_ATE
#
### Estimand : 1
# Estimand name: backdoor
# Estimand expression:
# d
# ────────────(E[salary|workclass,marital-status,native-country,relationship,occupation])
# d[education]
#
# Estimand assumption 1, Unconfoundedness: If U→{education} and U→salary then P(salary|education,workclass,marital-status,native-country,relationship,occupation,U) = P(salary|education,workclass,marital-status,native-country,relationship,occupation)
#
## Realized estimand
# b: salary~education+workclass+marital-status+native-country+relationship+occupation
# Target units: ate
#
## Estimate
# Mean value: 0.4697157228651772
# Step 4: Refute the obtained estimate using multiple robustness checks.
refute_results = model.refute_estimate(identified_estimand, estimate, method_name="random_common_cause")
摘要
-
输入数据:结果变量和处理变量。强烈建议提供一个因果 DAG。
-
要求: 处理变量必须是二元的,所有分类变量需要编码成数值。
-
(优)缺点: 输出包含许多细节,这些细节可能有利但也使理解变得复杂。无法从数据中学习因果 DAG。
库 5:CausalImpact。
CausalImpact 是一个用于使用贝叶斯结构 时间序列模型 进行因果推断的 Python 包 [4]。其主要目标是通过分析期望和观察到的时间序列数据之间的差异来推断给定干预的预期效果,如程序评估或处理效果分析。它假设响应变量可以通过线性回归精确建模。然而,响应变量不能受到已发生干预的影响。例如,如果一家公司想要推断某个营销活动对其“收入”的影响,则其每日“访问量”不能作为协变量使用,因为总访问量可能会受到该活动的影响。在以下示例中,我们将创建 100 个观测值,其中响应变量为 y
,预测变量为 x1
。请注意,在实际应用中,我们会有更多的预测变量。为了演示,干预效果 通过在时间点 71 之后将响应变量更改 10 个单位来创建。
为了估计因果效应,我们首先需要指定干预前期和干预后期。使用 run
函数,我们利用 MLE(默认)拟合时间序列模型,并估计因果效应。
# Installation
pip install causalimpact
# Import libraries
import numpy as np
import pandas as pd
from statsmodels.tsa.arima_process import arma_generate_sample
import matplotlib.pyplot as plt
from causalimpact import CausalImpact
# Generate samples
x1 = arma_generate_sample(ar=[0.999], ma=[0.9], nsample=100) + 100
y = 1.2 * x1 + np.random.randn(100)
y[71:] = y[71:] + 10
data = pd.DataFrame(np.array([y, x1]).T, columns=["y","x1"])
# Initialize
impact = CausalImpact(data, pre_period=[0,69], post_period=[71,99])
# Create inferences
impact.run()
# Plot
impact.plot()
# Results
impact.summary()
# Average Cumulative
# Actual 130 3773
# Predicted 120 3501
# 95% CI [118, 122] [3447, 3555]
# Absolute Effect 9 272
# 95% CI [11, 7] [326, 218]
# Relative Effect 7.8% 7.8%
# 95% CI [9.3%, 6.2%] [9.3%, 6.2%]
# P-value 0.0%
# Prob. of Causal Effect 100.0%
# Summary report
impact.summary(output="report")
平均列描述了处理后期的平均时间。累积列汇总了各个时间点的值,这在响应变量表示流量数量时(如查询、点击、访问、安装、销售或收入)是一个有用的视角。我们看到,在这个例子中,由于我们在数据中包含了这一点,因此因果效应的概率为 100%,P 值为 0,这是正确的。
图 4. CausalImpact 图(作者提供的图片)。
在图(图 4)中,我们可以看到三个面板,具体描述如下:在面板 1 中,我们可以看到数据以及专门针对处理后期的反事实预测。在第二个面板中,展示了观察数据与反事实预测之间的差异。这一差异表示了由模型确定的逐点因果效应。底部面板通过汇总前一个面板中的逐点贡献来描绘干预的累积效果。需要注意的是,这些推断的有效性依赖于协变量没有受到干预本身的影响。此外,模型假设在预期期间建立的协变量与处理时间序列之间的关系在整个后期保持一致。
总结
-
输入数据: 时间序列数据。
-
要求: 需要指定干预前和干预后的时间段。
-
优点: 可以确定时间序列数据上的因果效应。
图书馆之间的整体总结。
我们已经看到五个用于学习因果关系的库,每个库都有其自身的优缺点和解决特定问题的重点。总结如下:
BnLearn 库适用于具有离散、混合或连续数据集的情况,需要从数据中推导因果 DAG(结构学习)、进行参数学习(CPD 估计),或在你希望进行推理(do-calculus)的情况下。其优点在于许多预处理和后处理步骤都在管道中解决,从而减少了出错的可能性。它允许将类别(字符串)标签作为输入,并且可以使用统计测试(独立性测试)去除虚假的边。此外,还可以将边或节点列入白名单或黑名单。请注意,由于大多数处理步骤是在管道中执行的,因此无法完全自定义建模部分,这可能会降低用户的控制权。
Pgmpy 库包含低级统计函数,用于创建和组合各种结构学习、参数学习和推理方法。它在使用上比较灵活,但需要对贝叶斯建模有良好的知识。请注意,Bnlearn 库利用了来自 Pgmpy 的一些函数,因此功能有重叠。换句话说,如果你打算创建自定义管道,建议使用 Pgmpy。否则,我推荐使用 bnlearn。
CausalNex 库可以用来从数据中推导因果 DAG,并进行推理。请注意,它仅适用于离散数据集。在检测到因果 DAG 后,可以通过阈值来移除边。不幸的是,它无法处理类别标签,需要将其转换为数值。此外,还需要进行各种中间预处理步骤。例如,减少类别特征的基数,描述数值特征,并以字典形式指定每个节点可以采取的所有状态。尽管这给用户提供了更多的建模控制,但额外的建模增加了一层复杂性。
DoWhy 方法在你希望估计因果效应时是理想的,如果你有一个 结果变量 和一个 处理 变量。这非常适用于 A/B 测试的使用场景,例如确定酒店预订取消的原因。结果 可以设置为 “取消”,处理 可以设置为 “分配不同的房间”。另一个使用场景可能是估计会员奖励计划的效果。请注意,该库并不旨在从数据中学习因果 DAG。然而,它确实需要一个作为模型输入的因果 DAG,你需要手动提供。
CausalImpact 是一个用于因果推断的 Python 包,使用贝叶斯结构时间序列模型。虽然它也计算因果关系,但其功能与其他库不可比拟,因为它不打算创建因果 DAG 或对数据进行推断。目标是通过分析预期和观察到的时间序列数据之间的差异,推断给定干预的预期效果,例如程序评估或治疗效果分析。
最后的话。
对数据进行因果建模可能因多种原因而具有挑战性。首先需要选择合适的库,这些库与您的具体使用案例匹配。我们已经对五个流行的因果推断库进行了比较:Bnlearn、Pgmpy、CausalNex、DoWhy 和 CausalImpact。每个库都有其自身的优缺点,并且最适合某些特定的使用案例。我们使用相同的数据集比较了这些库的功能、能力和结果的可解释性。只有* CausalImpact*使用了不同的数据集,因为它只能建模连续值并且需要时间序列数据。
功能分组会将Bnlearn、Pgmpy 和 CausalNex 聚类在一起,因为它们包含结构学习、参数学习和推断的功能。第二组是DoWhy 和 CausalImpact,它们旨在测量处理变量对结果变量Y的因果效应,同时控制一组特征X。除了这里描述的五个因果库,还有更多库值得关注,如EconML、Pyro 和 CausalML,它们属于第二组(处理变量对结果变量的因果效应)。
保重。保持冷静。
干杯 E.
免责声明:我是 Python bnlearn 库的开发者。不管你对因果库的决定是什么,确保你完全理解建模部分并能够解释结果。
如果你觉得这篇关于贝叶斯因果学习的文章有帮助,欢迎 关注我 ,因为我会写更多关于贝叶斯统计的文章。如果你考虑加入 Medium 会员,你可以通过使用我的 推荐链接* 来稍微支持一下我的工作。价格和一杯咖啡相同,你可以每月无限阅读文章。*
软件
让我们联系吧!
与贝叶斯统计相关的工作
初学者指南,用于有效地确定变量之间的因果关系。
towardsdatascience.com ## 逐步指南:使用贝叶斯定理设计知识驱动模型
如果你没有数据但有专家知识。初学者指南将知识转化为计算机辅助…
towardsdatascience.com ## 寻找最佳增强模型的指南:使用贝叶斯超参数调优,但不包括…
增强型决策树算法可能会优于其他模型,但过拟合是一个现实的危险。使用…
towardsdatascience.com [## 分类模型的贝叶斯优化超参数调优实践指南
学习如何拆分数据、优化超参数、预防过拟合、选择最佳模型,并创建…
erdogant.medium.com](https://erdogant.medium.com/hands-on-guide-for-hyperparameter-tuning-with-bayesian-optimization-for-classification-models-2002224bfa3d?source=post_page-----d91e8306e25e--------------------------------)
参考文献
-
Census Income. (1996). UCI 机器学习库。
doi.org/10.24432/C5S595
。 -
Taskesen E, 使用 Python 中的贝叶斯结构学习检测因果关系的逐步指南,Medium,2021
-
Taskesen, E. (2020). 使用 bnlearn Python 包学习贝叶斯网络。 (版本 0.3.22) [计算机软件]。
-
Amit Sharma, Emre Kiciman. DoWhy: 一种用于因果推断的端到端库。2020 年。
arxiv.org/abs/2011.04216
-
Kay H. 等,使用贝叶斯结构时间序列模型推断因果影响,2015 年,《应用统计年鉴》(第 9 卷,第 247–274 页)
-
Xun Zheng 等,无泪的 DAG: 结构学习的连续优化,第 32 届神经信息处理系统会议(NeurIPS 2018),加拿大蒙特利尔。
-
AI4I 2020 预测维护数据集(2020)。UCI 机器学习库。
dbt 增量模型在大数据中的威力
在 BigQuery 上的实验
·
关注 发表在 Towards Data Science ·5 分钟阅读·2023 年 2 月 9 日
–
如果你在使用 dbt 模型处理几 MB 或 GB 的数据,这篇文章不适合你;你已经做得很好了!这篇文章是为那些需要在 BigQuery 中扫描 TB 级数据,以便计算某些计数、总和或滚动总计的可怜的灵魂准备的,这些计算可能是每天进行的,甚至更频繁。在这篇文章中,我将介绍一种针对“大数据”进行便宜的数据摄取和数据消费的技术。
图片由 Joshua Sortino 提供,来源于 Unsplash
让我们假设我们有以时间戳为粒度的原始数据,并且我们需要计算每个客户的总数。此外,我们希望进行按天和按月的分析,因为我们需要在不同粒度上汇总数据。我将以 Twitter 点赞为例,其中 user_id
是 Twitter 用户,post_id
是 Twitter 帖子,timestamp
是 user_id
点赞 post_id
的时间戳。原始数据以天为粒度存储为分区表。如果你的原始数据没有分区,这已经是你的第一个改进,因此考虑将原始数据迁移到分区表,参考 这个 示例。确保按需设置分区过滤器 cloud.google.com/bigquery/docs/managing-partitioned-tables#require-filter
,以避免扫描所有原始数据的意外情况。
Twitter 点赞的原始日志
这些数据实际上大约是 869.91 GB!所以如果我想知道每个用户所执行的点赞数量,我的[昂贵的]简单查询将是:
SELECT
user_id,
COUNT(*) as nr_likes
FROM twitter_likes
然而,Twitter 上的点赞是不可变的事件,无法被修改。这允许我逐步加载这些数据。我将遍历原始数据来加载当天的数据,检查每日汇总的数据以加载本月的数据,最后通过每日数据来更新发生变化的用户的点赞数。
对于这两种增量模型,我将使用静态分区技术,这被证明是性能最优的 discourse.getdbt.com/t/benchmarking-incremental-strategies-on-bigquery/981
。这意味着我将对分区进行操作而不是扫描数据,因此我会选择最近的两个日期分区(以防我有延迟到达的事件),对它们进行汇总,然后将其添加到我的 dbt 模型 twitter_likes_daily
中。
{% set partitions_to_replace = [
'current_date',
'date_sub(current_date, interval 1 day)'
] %}
{{
config(
materialized='incremental',
partition_by = { 'field': 'date', 'data_type': 'date' },
cluster_by = ["user_id"],
incremental_strategy = 'insert_overwrite'
partitions = partitions_to_replace
)
}}
SELECT
user_id,
DATE(timestamp) as date,
COUNT(*) as nr_likes
FROM {{ source ('twitter', 'twitter_likes')}}
-- I have a fake condition here because
-- my raw table does not let me query without a partition filter
WHERE DATE(timestamp) >= "1990-01-01"
{% if is_incremental() %}
-- this filter will only be applied on an incremental run
AND DATE(timestamp) in ({{ partitions_to_replace | join(',') }})
{% endif %}
GROUP BY 1, 2
twitter_likes_daily
仅为 55.52 GB,日加载扫描 536.6 MB!让我们创建每月汇总数据,为此,我们可以直接在 twitter_likes_daily
上工作,而不是原始的 twitter_likes
。
{% set partitions_to_replace = [
'date_trunc(current_date, month)',
'date_trunc(date_sub(current_date, interval 1 month), month)'
] %}
{{
config(
materialized='incremental',
partition_by = { 'field': 'month_year', 'data_type': 'date', "granularity": "month" },
cluster_by = ["user_id"],
incremental_strategy = 'insert_overwrite'
partitions = partitions_to_replace
)
}}
SELECT
user_id,
DATE_TRUNC(date, MONTH) as month_year,
SUM(nr_likes) as nr_likes
FROM {{ref('twitter_likes_daily')}}
{% if is_incremental() %}
-- this filter will only be applied on an incremental run
AND DATE_TRUNC(date, MONTH) in ({{ partitions_to_replace | join(',') }})
{% endif %}
GROUP BY 1, 2
twitter_likes_monthly
的大小为 14.32 GB,而一次加载仅扫描 2.5 GB!
现在让我们来看看总数。在这种情况下,我们可以决定运行一个 merge
操作,只更新过去 2 天内有新点赞的用户的总数。其余用户没有更新。因此,在这种情况下,我可以从 twitter_likes_daily
中挑选出今天有点赞的用户,计算他们在 twitter_likes_monthly
中的总数,最后将它们合并到包含总数的表 twitter_likes_total
中。请注意,该表使用 user_id
进行 clustered
,因为与 id 的集群合并比常规合并更高效 discourse.getdbt.com/t/benchmarking-incremental-strategies-on-bigquery/981
。
{{
config(
materialized='incremental',
unique_key = 'user_id',
cluster_by = ["user_id"],
incremental_strategy = 'merge'
)
}}
{% if is_incremental() %}
WITH users_changed as (
SELECT DISTINCT user_id
FROM {{ref('twitter_likes_daily')}}
WHERE date = current_date()
)
{% endif %}
SELECT
user_id,
sum(nr_likes) as total_likes
FROM {{ref('twitter_likes_monthly')}}
{% if is_incremental() %}
AND user_id in (select user_id from users_changed)
{% endif %}
group by 1
twitter_likes_total
是 1.6 GB,但扫描了 16 GB!所以从 869.91 GB 中,我们仅扫描了 16 GB,并生成了一些轻量且便宜的分析表。
然而,twitter_likes_total
的方法也可以是一个仅扫描 twitter_likes_monthly
的 dbt
表配置。如下所示:
{{
config(
materialized='table',
cluster_by = ["user_id"]
)
}}
SELECT
user_id,
sum(nr_likes) as total_likes
FROM {{ref('twitter_likes_monthly')}}
上述脚本将扫描 14.5 GB,因此实际上比 merge
过程更轻量!这到底发生了什么?重要的是要理解,merge
可能比完全重建表的效果差。鉴于我们扫描的表 twitter_likes_monthly
不大且为 14.5 GB,完全重建表比 merge
更高效。因此,当数据量不是很大时,merge
实际上可能是过度工程。
总结一下
原始不可变事件的扫描可以通过在 dbt
上使用增量模型高效处理,特别是使用静态插入+覆盖策略。改进可能是巨大的,从 TB 减少到每天几 GB 的 dbt 负载,但最重要的是,能够让数据分析师在原本成本高昂的原始表扫描上运行查询,转而使用一些轻量、预聚合的表。最后,虽然向表中追加数据在理论上可能比重建表便宜,但根据用例,有时使用 merge
可能比完全重建表的性能差,尤其是当源数据只有几 GB 时。
独立分量分析(ICA)在真实世界应用中的强大能力 — EEG 示例
独立分量分析(ICA)是一种强大的数据驱动工具,能够在数据中分离线性贡献
·
关注 发布在 Towards Data Science ·6 分钟阅读·2023 年 10 月 23 日
–
Photo by National Cancer Institute 在 Unsplash
- 介绍
独立成分分析(ICA)通常与降维任务相关。然而,这项技术最突出的应用是从数据中分离线性贡献,寻找统计独立的成分。例如,ICA 被广泛用于从音频中分离乐器轨迹。本文的目标是通过著名的 “鸡尾酒会” 示例介绍和激发 ICA 的兴趣,然后简要介绍 ICA 如何利用概率和信息理论的基础知识提取独立成分。接着,我们探索一个实际的例子——从脑电图(EEG)数据中识别和去除眼睛眨动的动作。
贡献与信号和 ICA 的分离
为了更好地理解从数据中分离线性贡献,标准的例子是 “鸡尾酒会”。这个场景设置在一个喧闹的鸡尾酒会中,两个人进行对话,分别由房间中的红色(Y1) 和蓝色(Y2) 扬声器表示。因此,两只麦克风被战略性地放置在附近,以捕捉这两个人的声音,由绿色(X1) 和黄色(X2) 麦克风表示。这些声音的可听度根据对话者与每个麦克风的距离而变化。
“鸡尾酒会” 的插图。信号由麦克风 X1 和 X2 捕获,我们希望将其分解为原始源 Y1 和 Y2。图片由作者提供。
挑战在于如何有效地将每个混合音频文件中的两个不同声音分离,以获得每个扬声器的孤立音频录音。这一难题可以通过熟练应用独立成分分析(ICA)来解决。在正确的数据假设下(我们将在下一节中看到),每个独立成分将是蓝色和红色扬声器的声音。
ICA 分解的插图。图片由作者提供。
ICA 的基本机制
独立成分分析(ICA)试图使用坐标系统 S 对数据 X 进行表示,在这个坐标系统中,成分 Si 是独立的。这种独立的坐标系统可以说是数据的一种自然表示。
这里最重要的一点是 ICA 假设成分是统计独立的,而不是无关的(如 PCA)[1]。
统计 独立 的定义。
主成分分析(PCA)假设数据和成分是正态分布的,当满足这一要求时,PCA 成分是统计独立的。然而,在数据分布非高斯的一般情况下,ICA 是寻找独立成分 Y (源) 代表数据 X (信号) 的强大方法。
ICA 变换由Y = AX给出,其中A的元素是定义组件的系数。因此,它仍然是一个线性方法,因为它假设源和信号之间的线性组合。
我们也可以进行反向操作,即根据 ICA 组件写入数据。
由于我们正在寻找统计上独立的组件,我们想要最小化组件之间共享的信息。
因此,信息论是帮助我们解决这个问题的自然选择。两随机变量X和Y之间共享的信息量被称为互信息。
因此,我们想要最小化组件之间共享的信息,
这是找到最小化I[Y] = I[ AˆT X ]的A。因此,结果的正交向量导致最独立的组件。
一旦我们理解了 ICA 的动机和基本机制,让我们在实际情况中应用它。
现实世界应用:去除 EEG 数据中的眼睛眨动
脑电图(EEG)数据受到各种运动伪影的污染,这些伪影由患者在图像采集过程中的自愿或非自愿运动引起。一个常见的伪影是眼睛眨动。因此,像 ICA 这样的工具是去除这些噪声的广泛使用技术[2]。如前所述,ICA 是一种分离信号独立贡献的方法。因此,我们直观地应该期待眼睛眨动对接近眼睛的通道有重要贡献。
EEG 测量跨皮层点的电位差随时间变化。因此,在 EEG 实验中,我们将从皮层的不同部位获得多条时间序列,这些序列包含关于神经回路的不同信息。
来自多个通道的 EEG 时间序列示例。图像由作者提供。
上述示例中的通道位置。图像由作者提供。
为了正确处理 EEG 数据并应用 ICA,我们将使用MNE(一个用于 EEG 数据处理和分析的库)。此外,我们将使用由 João Sato 教授在联邦 ABC 大学收集的数据,这些数据可在GitHub上找到。
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
import numpy as np
import mne as mne
sns.set_palette('pastel')
eeg = np.loadtxt('/Users/rodrigo/Post-Grad/Disciplines/Neural Signaling Processing/Aula_3/restEEG/OlhosAbertos.txt')
canais = pd.read_csv('/Users/rodrigo/Post-Grad/Disciplines/Neural Signaling Processing/Aula_3/restEEG/NOMEScanais.txt', dtype='str', header=None)
# Sampling frequency 250Hz:
Hz = 250
eeg = pd.DataFrame(eeg, columns=canais.iloc[:,0].values)
eeg['Time(s)'] = np.arange(0, (1/Hz)*len(eeg), 1/Hz)
创建 MNE EEG 对象。
# Creating MNE object
mne_info = mne.create_info(
list(eeg.columns[:-1]),
Hz,
ch_types=['eeg']*len(eeg.columns[:-1]))
mne_raw = mne.io.RawArray(eeg.iloc[:,:-1].values.T, mne_info)
# Defining the electrodes architecture
standard_1020 = mne.channels.make_standard_montage('standard_1020')
mne_raw.set_montage(standard_1020)
为了提取重要特征,我们需要过滤时间序列,以排除与伪影和噪声相关的高频(> 120 Hz)和低频(< 4 Hz),以及来自环境电源线的 60Hz 贡献。
# Set filter parameters
lowCut = 4 # Hz
highCut = 120 # Hz
freqNotch = 60 # Hz (eletric interference from the enviroment)
# Apply bandpass and notch filter
mne_raw.filter(lowCut, highCut, fir_design='firwin')
mne_raw.notch_filter(freqNotch, fir_design='firwin')
现在,数据准备好进行 ICA 处理。
ica_obj = mne.preprocessing.ICA(
n_components=0.95,
method='infomax',
max_iter="auto",
random_state=1,
fit_params=dict(extended=True)).fit(mne_raw)
ica = ica_obj.get_sources(mne_raw).get_data()
ica_obj.plot_components(picks=None,show=True, inst=mne_raw)
ICA 系数用作识别信号组件空间贡献的代理。
将值与每个成分的空间区域关联的过程示意图。图片由作者提供。
因此,可以为每个成分创建热图。
ICA 成分。图片由作者提供。
第一个成分是信号中的眼睑运动贡献。正如我们预期的那样,这一贡献在靠近眼睛的电极上更大。也可以更深入地查看该成分的时间/频率信息。
ICA 第一个成分的时间/频率信息。图片由作者提供。
现在,为了从数据中去除眼睑运动伪影,我们只需重建没有这个特定成分的信号。这意味着恢复没有与伪影相关的模式的时间序列。
结论
总结来说,独立成分分析(ICA)是一种强大的数据驱动工具,它假设数据为非高斯分布,并且源信号之间存在线性关系,以寻找统计上独立的成分。ICA 允许我们分离数据中的不同贡献,使其在需要隔离噪声或不需要的信号的场景中特别有用。在 EEG 数据的情况下,ICA 可以识别并去除如眼睑运动等运动伪影,这些伪影可能会污染神经信号。
备注
重要的是要指出,ICA 是一个强大的工具,但使用时需要谨慎。数据驱动的方法不假设对系统有任何先验知识。因此,混淆和虚假结论很常见,这使得检查要求和仔细评估结果至关重要。
致谢
这个项目的灵感来源于 João Ricardo Sato 教授。
本文的笔记本可用 在这里。
参考文献
[1] Hastie T, Tibshirani R, Friedman J, Hastie T, Tibshirani R, et al. (2009) 基础扩展与正则化。在:统计学习的元素,纽约,NY:Springer New York,Springer 统计系列。
[2] Mennes M, Wouters H, Vanrumste B, Lagae L, Stiers P. (2010). ICA 作为去除 EEG/ERP 中眼动伪影工具的验证。生理心理学。2010;47(6):1142–1150。doi:10.1111/j.1469–8986.2010.01015.x