表示连续概率分布

2023-12-29

我有一个涉及连续概率分布函数集合的问题,其中大部分是根据经验确定的(例如出发时间、过境时间)。我需要的是某种方法来获取其中两个 PDF 并对它们进行算术运算。例如。如果我有两个值 x 取自 PDF X,y 取自 PDF Y,我需要获取 (x+y) 的 PDF,或任何其他操作 f(x,y)。

分析解决方案是不可能的,因此我正在寻找允许此类操作的 PDF 的某种表示形式。一个明显(但计算成本昂贵)的解决方案是蒙特卡罗:生成大量 x 和 y 值,然后仅测量 f(x, y)。但这需要太多的CPU时间。

我确实考虑过将 PDF 表示为范围列表,其中每个范围具有大致相等的概率,从而有效地将 PDF 表示为均匀分布列表的并集。但我不知道如何将它们结合起来。

有人对这个问题有什么好的解决办法吗?

Edit:目标是创建一种用于操作 PDF 的迷你语言(也称为领域特定语言)。但首先我需要理清底层的表示和算法。

Edit 2:dmckee 建议使用直方图实现。这就是我通过统一分布列表得到的结果。但我不知道如何将它们组合起来创建新的发行版。最终,在 P(x

Edit 3:我有一堆直方图。它们分布不均匀,因为我是根据出现数据生成它们的,所以基本上,如果我有 100 个样本并且我想要直方图中有 10 个点,那么我会为每个条形分配 10 个样本,并使条形宽度可变但面积恒定。

我发现要添加 PDF,你需要对它们进行卷积,并且我已经为此加强了数学知识。当您对两个均匀分布进行卷积时,您会得到一个具有三个部分的新分布:较宽的均匀分布仍然存在,但每侧都有一个三角形,其宽度为较窄的分布。因此,如果我对 X 和 Y 的每个元素进行卷积,我将得到一堆这些元素,全部重叠。现在我试图弄清楚如何将它们全部相加,然后得到一个最接近它的直方图。

我开始怀疑蒙特卡洛到底是不是一个坏主意。

Edit 4: 这张纸 http://www.heldermann-verlag.de/eqc/eqc01_16/eqc16002.pdf详细讨论了均匀分布的卷积。一般来说,您会得到“梯形”分布。由于直方图中的每个“列”都是均匀分布的,因此我希望可以通过对这些列进行卷积并对结果求和来解决该问题。

然而,结果比输入复杂得多,并且还包括三角形。Edit 5:[错误内容已删除]。但是,如果这些梯形近似于具有相同面积的矩形,那么您就会得到正确的答案,并且减少结果中的矩形数量看起来也非常简单。这可能是我一直在寻找的解决方案。

Edit 6:解决了!这是这个问题的最终 Haskell 代码:

-- | Continuous distributions of scalars are represented as a
-- | histogram where each bar has approximately constant area but
-- | variable width and height.  A histogram with N bars is stored as
-- | a list of N+1 values.
data Continuous = C {
      cN :: Int,
      -- ^ Number of bars in the histogram.
      cAreas :: [Double],
      -- ^ Areas of the bars.  @length cAreas == cN@
      cBars :: [Double]
      -- ^ Boundaries of the bars.  @length cBars == cN + 1@
    } deriving (Show, Read)


{- | Add distributions.  If two random variables @vX@ and @vY@ are
taken from distributions @x@ and @y@ respectively then the
distribution of @(vX + vY)@ will be @(x .+. y).

This is implemented as the convolution of distributions x and y.
Each is a histogram, which is to say the sum of a collection of
uniform distributions (the "bars").  Therefore the convolution can be
computed as the sum of the convolutions of the cross product of the
components of x and y.

When you convolve two uniform distributions of unequal size you get a
trapezoidal distribution. Let p = p2-p1, q - q2-q1.  Then we get:


>   |                              |
>   |     ______                   |
>   |     |    |           with    |  _____________
>   |     |    |                   |  |           |
>   +-----+----+-------            +--+-----------+-
>         p1   p2                     q1          q2
> 
>  gives    h|....... _______________
>            |       /:             :\
>            |      / :             : \                1
>            |     /  :             :  \     where h = -
>            |    /   :             :   \              q
>            |   /    :             :    \
>            +--+-----+-------------+-----+-----
>             p1+q1  p2+q1       p1+q2   p2+q2

However we cannot keep the trapezoid in the final result because our
representation is restricted to uniform distributions.  So instead we
store a uniform approximation to the trapezoid with the same area:

>           h|......___________________
>            |     | /               \ |
>            |     |/                 \|
>            |     |                   |
>            |    /|                   |\
>            |   / |                   | \
>            +-----+-------------------+--------
>               p1+q1+p/2          p2+q2-p/2

-}
(.+.) :: Continuous -> Continuous -> Continuous
c .+. d = C {cN     = length bars - 1,
             cBars  = map fst bars, 
             cAreas = zipWith barArea bars (tail bars)}
    where
      -- The convolve function returns a list of two (x, deltaY) pairs.
      -- These can be sorted by x and then sequentially summed to get
      -- the new histogram.  The "b" parameter is the product of the
      -- height of the input bars, which was omitted from the diagrams
      -- above.
      convolve b c1 c2 d1 d2 =
          if (c2-c1) < (d2-d1) then convolve1 b c1 c2 d1 d2 else convolve1 b d1 
d2 c1 c2
      convolve1 b p1 p2 q1 q2 = 
          [(p1+q1+halfP, h), (p2+q2-halfP, (-h))]
               where 
                 halfP = (p2-p1)/2
                 h = b / (q2-q1)
      outline = map sumGroup $ groupBy ((==) `on` fst) $ sortBy (comparing fst) 
$ concat
                [convolve (areaC*areaD) c1 c2 d1 d2 |
                 (c1, c2, areaC) <- zip3 (cBars c) (tail $ cBars c) (cAreas c),
                 (d1, d2, areaD) <- zip3 (cBars d) (tail $ cBars d) (cAreas d)
                ]
      sumGroup pairs = (fst $ head pairs, sum $ map snd pairs)

      bars = tail $ scanl (\(_,y) (x2,dy) -> (x2, y+dy)) (0, 0) outline
      barArea (x1, h) (x2, _) = (x2 - x1) * h

其他运算符留给读者作为练习。


不需要直方图或符号计算:如果采取正确的观点,一切都可以在语言级别以封闭形式完成。

[我将交替使用术语“度量”和“分布”。另外,我的 Haskell 已经生锈了,请您原谅我在这方面的不精确。]

概率分布确实是codata.

令 mu 为概率测度。您可以对度量做的唯一一件事是将其与测试函数集成(这是“度量”的一种可能的数学定义)。请注意,这是您最终要做的:例如,针对身份进行集成将取平均值:

mean :: Measure -> Double
mean mu = mu id

另一个例子:

variance :: Measure -> Double
variance mu = (mu $ \x -> x ^ 2) - (mean mu) ^ 2

另一个例子,计算 P(mu

cdf :: Measure -> Double -> Double
cdf mu x = mu $ \z -> if z < x then 1 else 0

这表明了一种二元性的方法。

方式Measure因此应表示类型(Double -> Double) -> Double。这允许您对 MC 仿真、针对 PDF 的数值/符号求积等结果进行建模。例如,函数

empirical :: [Double] -> Measure
empirical h:t f = (f h) + empirical t f

返回 f 对 an 的积分经验测量通过例如获得。 MC 采样。还

from_pdf :: (Double -> Double) -> Measure
from_pdf rho f = my_favorite_quadrature_method rho f

根据(常规)密度构建测量值。

现在,好消息。如果 mu 和 nu 是两个度量,则卷积 mu ** nu是(谁)给的:

(mu ** nu) f = nu $ \y -> (mu $ \x -> f $ x + y)

因此,给定两个度量,您可以针对它们的卷积进行积分。

另外,给定一个随机变量 X 的规律mu,a * X 的定律由下式给出:

rescale :: Double -> Measure -> Measure
rescale a mu f = mu $ \x -> f(a * x)

此外,phi(X) 的分布由下式给出图像测量 http://en.wikipedia.org/wiki/Pushforward_measurephi_* X,在我们的框架中:

apply :: (Double -> Double) -> Measure -> Measure
apply phi mu f = mu $ f . phi

因此,现在您可以轻松地制定出用于测量的嵌入式语言。这里还有很多事情要做,特别是关于实线以外的样本空间、随机变量之间的依赖性、条件,但我希望你明白这一点。

特别是,前推是函数式的:

newtype Measure a = (a -> Double) -> Double
instance Functor Measure a where
    fmap f mu = apply f mu

它也是一个 monad(练习 -- 提示:这看起来很像延续 monad。什么是return?什么是类似物call/cc ?).

此外,与微分几何框架相结合,这可能会变成自动计算贝叶斯后验分布的东西。

在一天结束时,你可以写类似的东西

m = mean $ apply cos ((from_pdf gauss) ** (empirical data))

计算 cos(X + Y) 的平均值,其中 X 具有 pdfgaussY 已通过 MC 方法采样,其结果为data.

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

表示连续概率分布 的相关文章

  • 加权图的 BFS 算法 - 寻找最短距离

    我看过很多帖子 即 post1 https stackoverflow com questions 30409493 using bfs for weighted graphs post2 https cs stackexchange co
  • 算法 - 树中所有节点的最大距离

    所以 找到树中两个节点之间的最长路径相当容易 但我想要的是找到从节点出发的最长路径x到树中的另一个节点 对于所有x 这个问题也可以用以下方式表达 计算从给定的树中可以生成的所有有根树的高度 One way of course is to j
  • 为什么我不能声明推断类型?

    我有以下内容 runcount Eq a Num b gt a gt b runcount runcountacc 0 runcountacc Eq a Num b gt b gt a gt b runcountacc n runcount
  • Python 比编译的 Haskell 更快?

    我有一个用 Python 和 Haskell 编写的简单脚本 它读取包含 1 000 000 个换行符分隔的整数的文件 将该文件解析为整数列表 对其进行快速排序 然后将其写入已排序的不同文件中 该文件与未排序的文件具有相同的格式 简单的 这
  • 使用 Parsec 解析正则表达式

    我正在尝试通过实现一个小型正则表达式解析器来学习秒差距 在 BNF 中 我的语法类似于 EXP EXP LIT EXP LIT 我尝试在 Haskell 中实现这一点 expr try star lt gt try litE lt gt l
  • 并行 Haskell - GHC GC 火花

    我有一个正在尝试并行化的程序 带有可运行代码的完整粘贴here http lpaste net 101528 我进行了分析 发现大部分时间都花在findNearest这本质上是一个简单的foldr超过一个大Data Map findNear
  • n的渐近增长选择下限(n/2)

    如何找到 n select Floor n 2 的渐近增长 我试过 使用扩展并得到它等于 n n 1 floor n 2 1 n floor n 2 知道我该如何从那里去吗 感谢任何帮助 更喜欢提示而不是答案 我同意上面的答案 但想提供更多
  • 这是 unsafeCoerce 的安全使用吗?

    我遇到的情况是 我目前正在使用极其可怕的函数 unsafeCoerce 幸运的是 这并不是为了任何重要的事情 但我想知道这是否是该函数的安全使用 或者是否有其他方法可以解决其他人知道的这个特定问题 我的代码类似于以下内容 data Toke
  • 检索 Haskell 项目中所有导入的列表

    因此 我的最终目标是通过确保项目导入的所有实体都存在于其声称可以使用的版本中 来评估 cabal 文件中依赖项的准确性 一个好的开始是找到单个源文件使用的所有导入实体的列表 可选地包含有关它们来自何处的信息 我愿意暂时忽略类实例的情况 因为
  • CSR 矩阵 - 矩阵乘法

    我有两个方阵A and B 我必须转换B to CSR Format并确定产品C A B csr C 我在网上找到了很多关于CSR 矩阵 向量乘法 http www mathcs emory edu cheung Courses 561 S
  • 让电脑实现360度=0度,旋转炮塔

    我正在制作一个游戏 其中有一个计算机控制的炮塔 炮塔可360度旋转 它使用 trig 找出枪瞄准所需的角度 obj deg 并将枪的当前角度存储在 gun deg 下面的代码以设定的速度旋转枪 if objdeg gt gundeg gun
  • 优化 Haskell 内循环

    仍在 Haskell 中进行 SHA1 实现 我现在已经有了一个有效的实现 这是内部循环 iterateBlock Int gt Word32 gt Word32 gt Word32 gt Word32 gt Word32 gt Word3
  • Haskell数据类型转换问题

    我目前正在学习 Haskell 并且一直在编写一些非常简单的程序来练习 我的程序之一是 import System IO main do putStrLn Give me year y lt getLine let res show cal
  • 如何从列中创建对称矩阵?

    例如 我想转动以下列 90 175 600 650 655 660 代入矩阵 90 175 600 650 655 660 175 600 650 655 660 655 600 650 655 660 655 650 650 655 66
  • 如何与更高级别的类型合作

    玩弄教堂的数字 我遇到了无法指导 GHC 类型检查器处理高阶类型的情况 首先我写了一个版本 没有任何类型签名 module ChurchStripped where zero z z inc n z s s n z s natInteger
  • R 中 if-else 中的逻辑运算符

    我有一个名为 mat 的下表 5 列和 3 行 AC CA RES 1 0 2 2 1 3 0 0 0 1 正在执行的操作是mat 1 mat 1 mat 2 我正在测试以下内容 1 如果一行的两列都为零 则结果应为 NA 2 如果一行中只
  • 我该如何解决? KnapSack - 值完全相同,但每个对象都有三个权重

    我在解决我的练习时遇到问题 我读到了动态规划和算法 我认为我的练习是 特定背包问题 我用暴力法解决了它 但我无法用动态规划解决它 我有一艘重300吨的船 背包 有些晶体本身含有 3 种物质 X Y Z 每种物质都有重量 并且所有晶体都具有相
  • 单词预测算法

    我确信有一篇关于此问题的帖子 但我找不到提出这个确切问题的帖子 考虑以下 我们有字典可供使用 我们收到了许多单词段落 我希望能够根据此输入预测句子中的下一个单词 假设我们有几个句子 例如 你好 我的名字是汤姆 他的名字是杰瑞 他去了没有水的
  • 具有最小刻度的图表的漂亮标签算法

    我需要手动计算图表的刻度标签和刻度范围 我知道漂亮刻度的 标准 算法 参见 我也知道这个Java实现 http erison blogspot nl 2011 07 algorithm for optimal scaling on char
  • 不同类型的列表?

    data Plane Plane point Point normal Vector Double data Sphere Sphere center Point radius Double class Shape s where inte

随机推荐

  • R中的预分配列表

    在 R 中 在循环中扩展数据结构效率很低 我如何预分配list具有一定的尺寸 matrix通过以下方式可以轻松做到这一点ncol and nrow论据 如何在列表中做到这一点 例如 x lt list for i in 1 10 x i l
  • 如何找到实体框架的水晶报表?

    如何将 Crystal Reports 绑定到实体框架实体 我确实还没有找到解决方法 而且我还没有足够的积分来对现有问题进行投票 至于我 我阅读了下面的链接 http aspalliance com 2049 Use LINQ to Ret
  • 使用IoC时单元测试的策略应该是什么?

    在读完有关依赖注入和 IoC 的所有内容后 我决定尝试在我们的应用程序中使用 Windsor Container 它是一个 50K LOC 多层 Web 应用程序 所以我希望它不是一个矫枉过正的东西 我使用了一个简单的静态类来包装容器 并在
  • 修复我的网络活动指示器

    我的网络活动指示器有一个问题 有时它会在不应该显示的情况下继续显示 我为它编写了自己的管理器 并将其替换为使用NSAssert像这样的声明 void setNetworkActivityIndicatorVisible BOOL setVi
  • 如何在 Eclipse 中的 Android 库项目中引用外部 jar

    哦 安卓 我多么喜欢你的言辞 我有一个工作区 里面有一些项目 App1和App2是Android应用程序 Common是一个Android库项目 App1 和 App2 依赖于 Common 通过 Android 选项卡链接 Common
  • 更新命令行应用程序状态

    我有一个命令行应用程序 当前打印增加的百分比 1 2 3 4 输出是连续的 但我见过命令行工具显示更改 就好像它是内联更新一样 1 2 与第一个位置相同 3 与第一个位置相同 4 与第一个位置相同 我怎样才能做到这一点 我正在使用 Java
  • 将 3 列文件转换为矩阵格式

    我有一个如下例所示的文件格式 显示了 5 个人 包括他们自己 之间的关系 1 1 1 0 2 1 0 5 3 1 0 1 4 1 0 3 5 1 0 1 2 2 1 0 3 2 0 5 4 2 0 2 5 2 0 3 3 3 1 0 4 3
  • 静态属性和单例有什么区别?

    使用 C 实现的单例可能类似于 public class Singleton private static Singleton instance private Singleton public static Singleton Insta
  • 如何使用 NSUserDefaults 注册用户默认值而不覆盖现有值?

    我有一个 AppDelegate 类 void initialize我用来注册一些默认值的方法 这是我使用的代码 void initialize NSDictionary defaults NSDictionary dictionaryWi
  • 可移植的有符号/无符号字节转换,C++

    我正在使用有符号到无符号字节 int8 t 转换来打包字节 uint32 t uint8 t byte lt lt n 这可以在 Intel Linux 上使用 GCC 运行 对于其他平台 编译器 例如 PowerPC 是否可移植 有更好的
  • JavaFX:使用自定义节点作为 TreeView 的折叠/展开分支开关

    是否可以替换a的展开和折叠箭头TreeView与定制Node 形状 而不是image https stackoverflow com q 15521806 1725096 The fx shape箭头的 css 属性提供了基本的 SVG 形
  • 我如何 - 长轮询和调度程序?

    我正在尝试安排一个长轮询机制 我想知道是否可以利用调度程序来实现这一点 到目前为止 这就是我一直在想的 通过计时器进行调度 但仅在上一次迭代已经完成的情况下才将下一次迭代排入队列 当上一次迭代完成时 将下一次迭代排入队列 我一直在查看现有的
  • 检测图像中的图案并检索其位置

    我有一张 1920x1080 的图像 我需要获取图像中每个矩形的位置 最好为 2 个点 左上 右下 我对 Python 还很陌生 我想也许使用 opencv 模块 但如果你能给我一些指示 那将会非常有帮助 Thanks 我建议使用OpenC
  • 如何在不安装Python Egg文件的情况下直接运行它们?

    是否可以像使用 Java 运行 jar 文件一样直接运行 Python Egg 文件 例如 使用 Java 您可能会执行以下操作 java jar jar file A 蟒蛇蛋 http peak telecommunity com Dev
  • C# 按原始顺序获取 FieldInfos/PropertyInfos?

    如何按照它们在类中的排列顺序获取类型 FieldInfos PropertyInfos 作为 MemberInfo 数组 class Test public bool First get set public int Second publ
  • 如何使新的 numpy 数组与给定数组大小相同并用标量值填充它

    我有一个带有一些随机数的 numpy 数组 如何创建一个具有相同大小的新数组并用单个值填充它 我有以下代码 A np array 2 2 2 2 B np copy A B B fill 1 我想要一个与 A 大小相同但填充 1 的新数组
  • 在 Oracle 中通过 SELECT 查询声明变量并设置其值

    在 SQL Server 中我们可以使用这个 DECLARE variable INT SELECT variable mycolumn from myTable 我怎样才能在 Oracle 中做同样的事情 我目前正在尝试以下操作 DECL
  • 使用 bookworm 在 Docker 中为 PHP8.1 添加 Memcached 支持时出现问题

    我有一个Dockerfile依靠PHP 8 1 apache 运行了几个月 Once PHP 8 1 apache开始使用 Debian bookworm memcached 客户端在构建镜像时开始出错 The Dockerfile涉及的行
  • 如何制作八角形盒子?

    我想使用 CSS 创建下面给出的框 我进行了很多搜索并阅读了非常好的文章 例如https css tricks com examples ShapesOfCSS https css tricks com examples ShapesOfC
  • 表示连续概率分布

    我有一个涉及连续概率分布函数集合的问题 其中大部分是根据经验确定的 例如出发时间 过境时间 我需要的是某种方法来获取其中两个 PDF 并对它们进行算术运算 例如 如果我有两个值 x 取自 PDF X y 取自 PDF Y 我需要获取 x y