NoiseLang：N = 5 是一个狄拉克δ函数

# NoiseLang：当 N=5 是狄拉克δ函数时——Manu Martínez-Almeida 来源：https://manualmeida.dev/articles/noiselang/ 在我的电信工程学位（https://www.tel.uva.es/en/studies/degrees/itec.htm）期间，我修了一门信号与噪声的课程，花了很多个晚上手写概率：期望、方差、两个随机变量落在某个区域的概率。这总是很糟糕，当我试图在计算机上运行时，有太多样板代码。 这个愿望变成了 NoiseLang（https://noiselang.com/）。我大约九年前开始了它，然而，我从未完成它。直到最近，借助 AI 工具以及一些远比第一次独自构建时更雄心勃勃的东西，我才重新拾起它。 ## 一切都是分布 整个语言建立在一个想法之上，即**每个值都是一个概率分布**。一个普通数字是一个狄拉克脉冲，一个将所有权重集中在一个值上的分布。由于常数和随机变量是同一类对象，语言中的每个运算符都将分布映射到分布。 一个名称总是指向一个固定的节点，就像在整页数学中 `X` 是同一个 `X` 一样。所以 `X + X` 是 `2X`，而 `X - X` 正好是 `0`。如果你想要变量独立性，你可以分别写抽取，即多次使用 `~`，或者使用 `~[N]` 来抽取 `N` 个独立变量到一个向量中。 ``` X ~ unif_int(1, 6) Y ~ unif_int(1, 6) X + Y # 两个独立的骰子，一个真正的 2d6 分布 ``` 在你询问某些结果之前，什么都不会运行，例如 `P(X + Y < 10)`，此时它会强制运行时运行数百万次模拟（如果可能的话跨所有核心），并返回一个带有标准误差的估计值。 ``` Bday = unif_int(1, 365) days ~[23] Bday # 房间里有 23 个人 P(has_duplicates(days)) # 生日悖论，约 0.507 ``` 这比描述起来更容易观察，所以我构建了一些很酷的演示！ #### 每个值都是一个分布 ``` A = 5 D1 ~ unif_int(1, 6) D2 ~ unif_int(1, 6) S = A + D1 + D2 E(S) ``` 绘制 **0** **A** ≈ **—** 图 1. 通过四个步骤的同一张直方图：一个常数是一个单峰，波浪号将其扩散成一个骰子，添加两个骰子使其向钟形曲线弯曲，最后查询读出均值。 #### 一个数字是一个分布 `A = 5` 看起来像一个普通常数，确实如此，但 Noise 将其视为一个每次抽取都落在 5 上的概率分布。统计学家称之为狄拉克δ函数，一个无限细的单峰。 #### 波浪号使其扩散 `D1 ~ unif_int(1, 6)` 绑定了一个公平骰子，一个离散均匀分布，所有六个面等可能。单峰扩散成六个平坦条，所以 `D1` 现在是不确定的，但你仍然像任何其他变量一样写它。 #### 连接它们，钟形曲线出现 用 `S = A + D1 + D2` 将常数和两个**独立**的骰子相加。平坦的形状**卷积**成一个在 12 处达到峰值的三角形，已经向钟形曲线弯曲，如果你再添加几个，它们会锐化成一个真正的高斯分布。这就是中心极限定理的免费体现！ #### 一次查询折叠数百万次抽取 在你询问之前，什么都不会运行。`E(S)`，期望值，触发一次蒙特卡洛传递，将**数百万**次投掷平均成一个数字，根据大数定律，运行中的均值稳定在 `12`。你写了几行数学代码，而一个专家级的内核在底层运行。 ## 为什么它搁置了九年 设计从来都不是难点，因为为这种语言编写一个解析器和树遍历解释器只是一个周末的工作量。问题在于其他所有事情：编写一个高效的蒙特卡洛运行时，而不是一个朴素的解释器、条件贝叶斯推理等等。 当前版本是一个编译器，一个 JIT（使用出色的 Cranelift (https://cranelift.dev/))、一个 WASM (https://webassembly.org/) 后端，以及一堆精心编写的数值代码，所以对于一个可爱的周末项目来说，它始终遥不可及。 ## 使用智能体构建雄心勃勃的版本 在我的日常工作和副业项目中，我正在探索当今 AI 智能体能够做什么的边界。例如，我正在将我 15 年前为 iOS 构建的一款游戏，从过时的 Objective-C 移植到现代游戏引擎（相对成功）。 对于 NoiseLand，我意识到，AI 非常擅长构建 JIT 部分、运行时部分、数值部分，但它不擅长提出好的语言设计想法，很多时候会为了不同的目的重写现有语言特性，或者为非正交特性提出不同的语法。 ## 一个 IR，三个后端 在底层，`~` 和分布构造函数构建了一个只追加的 DAG，称为 RvGraph。这个图是唯一的真实来源，之后会被转换成三种不同的代码路径： - 一个列式批处理解释器，可以在任何地方工作，并充当正确性仲裁者； - 一个 Cranelift JIT，将整个表达式融合成一个原生内核； - 一个 WASM 发射器，为浏览器做同样的事情。 回退 回退 X ~ unif(1, 6) RvGraph 批处理解释器 Cranelift JIT WASM 发射器 一个共享模块定义了图形的含义，因此两个代码生成器保持纤薄，并且不会偏离。任何后端无法成功编译的内容都会回退到解释器，并且结果在后端和核心数量之间保持一致。所有测试都在所有三个代码路径中运行，并进行比较以确保位一致。 ## 使蒙特卡洛循环廉价 所有性能工作都围绕一个循环：抽取几百万个样本，对每个样本评估表达式，并汇总结果。少数技术承担了大部分工作，同时保持结果确定性（这是困难的部分）。 内核融合将每个中间值保存在寄存器中，因此算术密集型表达式可以留在寄存器中。PRNG（xoshiro256++ (https://prng.di.unimi.it/)）编译到内核中，而 `ln`、`sin` 和 `cos` 变成内联多项式近似，将内核速度提高了 2 倍。 我最喜欢的技巧在 RNG 中。生成随机数是一个串行依赖链，所以与其与之对抗，不如让内核同时运行四个独立流，并让乱序核心重叠它们。这个技巧最终击败了手写 SIMD 内核！ 在我的 14 核 M4 Pro 上，一行 `P(...)` 可以达到每秒约 58 亿个样本，并且从单核扩展到全核时速度提升约 9.6 倍。每个核心，生成的内核运行速度大约是 LLVM 编译的手写 Rust 的 1.15 倍。相同的融合循环，作为 WASM 发出，在 V8 内部运行速度大约是原生速度的一半到四分之三。 ``` X ~ rand::unif(-1, 1) Y ~ rand::unif(-1, 1) C = X^2 + Y^2 < 1 4 * P(C) ``` 飞镖 **0** 内部 **0** π ≈ **—** 图 4. 蒙特卡洛 π。圆覆盖了正方形的 π/4 面积，因此随着飞镖的积累，青色比例的四倍收敛于 π。 #### 两次随机抽取 这和之前一样，使用了两次 `~`。`X ~ unif(-1, 1)` 和 `Y ~ unif(-1, 1)` 各自在 [-1, 1] 范围内抽取一个数字，它们一起在正方形中选择一个随机点。 #### 问一个是/否问题 `C = X^2 + Y^2 < 1` 恰好当点落在单位圆内时为真。在 Noise 中，这个比较本身就是一个随机变量，一个每次抽取结果为真或假的伯努利分布。 #### 现在投掷飞镖 每个飞镖是 `(X, Y)` 的一次抽取，落在圆内时为青色，否则为灰色。它们均匀分布，因为均匀抽取不会偏爱任何点。 #### 四倍的比例就是 π 圆覆盖了正方形面积的 π/4，因此青色比例在 π/4 附近徘徊，这使得 `4 * P(C)` 成为 π 的估计值。你投掷的飞镖越多，估计值就越精确。 ## Noise 的定位 NoiseLang 是一个玩具语言，**你很可能不应该将其用于任何严肃的事情**，然而，我希望这个语言在我的大学时代就存在。 对于一个语言迷来说，它是一个小的、静态的随机变量代数，带有前向蒙特卡洛、基于表达式和基于拒绝的条件推理语言。 你可能会问，它与 NumPy 或 Stan 相比如何？NumPy 让你自己编写模拟，Stan 让你声明一个模型并等待采样器。Noise 让你像写数学一样写概率，同时在底层运行蒙特卡洛来获得答案。 表 1. Noise 的定位。Stan 和 PyMC 在将后验拟合到大量数据方面优于它，NumPy 在原始数组计算方面优于它，但 Noise 可以比它们任何一个都更快地将你从概率问题带到可见答案。Stan 和 PyMC 在它们为 Built for 的事情上优于 Noise：使用它们的 HMC/NUTS 采样器将后验拟合到大量连续数据，NumPy 在原始数组计算方面优于它。Noise 中的条件是基于拒绝的，因此它适用于少量离散观测，但对于一万个连续测量来说变得毫无用处，并且还没有状态模拟（还没有马尔可夫链）。Noise 的胜出之处在于：当你有一个概率问题并且想毫不费力地知道答案时。 因此，在问题的白板阶段使用 Noise，当你想要运行你刚刚写的数学时，并在你需要真正的后验时转向 Stan 或 PyMC，或者在你需要投入生产时转向 NumPy 和 JAX。 ## 回到信号与噪声 回到我的大学时代，有一门课叫做“Señales Aleatorias Y Ruido”，这是“随机信号与噪声”的西班牙语翻译。 教科书《Señales Aleatorias y Ruido》的封面 这本教科书，以其所有的荣耀。 这门课是许多学生的噩梦，包括我自己，事实上，我挂科了。说实话，第一年我没有投入足够的努力，但第二年我不得不再次修这门课时，情况发生了变化。教授很棒，让这门课变得有趣。让我印象深刻的是，他能够解释为什么 FM 能在有噪声的信道中生存而 AM 却不能。 所以，这是我向这门课致敬的方式，一个单屏的 Noise 程序，模拟了为什么 FM 能在有噪声的信道中生存而 AM 却不能。 ``` msg = 0.3 * signal::sine(3); am_modulate(m) = 1 + m; fm_modulate(m, swing) = exp(i * swing * m); static ~ signal::noise_white_complex(0.4); rx_am = am_modulate(msg) + static; rx_fm = fm_modulate(msg, 3) + static; am_demodulate(c) = abs(c) - 1; fm_demodulate(c, swing) = arg(c) / swing; rec_am = am_demodulate(rx_am); rec_fm = fm_demodulate(rx_fm, 3); Print("AM error", E(mse(rec_am, msg))); Print("FM error", E(mse(rec_fm, msg))); ``` 滚动浏览步骤 → 图 3. AM vs FM 端到端。载波是一个相量，AM 将消息写入其长度，FM 写入其角度，因此相同的静态会破坏 AM 读取的内容，而几乎不影响 FM 读取的内容。 #### 消息 在 Noise 中，一个值可以携带整个信号，这里它是一个温和的音调，一个慢速的三周期正弦波。这是我们想要通过噪声信道并完好无损地恢复的内容。 #### AM 将消息放入长度 载波是一个复数，一个旋转的箭头。**AM** 的调制器将消息放入箭头的*长度*中，因此尖端向内和向外滑动，随着消息的上升和下降穿过单位圆。 #### FM 将其放入角度 **FM** 将相同的消息放入*角度*中，因此尖端沿着单位圆移动，而其长度保持不变。信息是相同的，只是编码在旋转中。 #### 添加相同的静态 现在相同的白噪声冲击两个尖端。它恰好沿着 AM 读取的*半径*涂抹 AM，但只是轻轻地推动 FM *围绕*圆，几乎不改变角度。你可以从两个云的形状中看到这一点。 #### 恢复 解调：对 AM 读取*长度*，对 FM 读取*角度*，并叠加在原始消息上。你已经可以看到，AM 轨迹恢复得参差不齐，而 FM 轨迹紧贴原始轨迹。 #### 相同的风暴，非常不同的损伤 `mse` 在波形上平均误差，而 `E` 在静态的多次抽取上平均该误差。相同的信号能量，相同的噪声能量，然而 FM 恢复得干净约 5 倍，并且在低噪声极限下，优势接近 swing² = 9。这就是它花费在角度上的带宽带来的回报。 ## 在浏览器中运行 NoiseLang 以上所有内容都在 `@noiselang/core` (https://www.npmjs.com/package/@noiselang/core) 上运行，这得益于编译为 WebAssembly 的 Rust 引擎。 ``` npm install @noiselang/core ``` ``` import { run } from "@noiselang/core"; const result = await run(` X ~ rand::unif(-1, 1); Y ~ rand::unif(-1, 1); 4 * P(X^2 + Y^2 < 1) `); console.log(result.value); // "3.1415..." — 最后一个语句的值 console.log(result.output); // 所有 Print(...) 发出的内容 ``` `run` 从不抛出异常，失败会通过 `result.error` 返回，并带有源代码跨度。还有 `runWithIntrospection`，这是 noiselang.com (https://noiselang.com/) 上变量检查器背后的 API。 NoiseLang 可以在浏览器中在 noiselang.com (https://noiselang.com/) 上玩耍。打开它，输入 `X ~ unif(-1, 1); Y ~ unif(-1, 1); 4 * P(X^2 + Y^2 < 1)`，然后观看几百万次抽取从你的浏览器选项卡中估计 π。