У меня проблема с набором функций непрерывного распределения вероятностей, большинство из которых определяются эмпирически (например, время отправления, время прохождения). Мне нужен какой-то способ взять два из этих PDF-файлов и сделать арифметику с ними. Например. если у меня есть два значения x, взятых из PDF X, и y, взятых из PDF Y, мне нужно получить PDF для (x + y) или любую другую операцию f (x, y).
Аналитическое решение невозможно, поэтому я ищу какое-то представление PDF-файлов, которое позволяет такие вещи. Очевидным (но вычислительно дорогим) решением является метод Монте-Карло: генерировать множество значений x и y, а затем просто измерять f (x, y). Но это занимает слишком много процессорного времени.
Я действительно думал о представлении PDF в виде списка диапазонов, где каждый диапазон имеет примерно равную вероятность, эффективно представляя PDF как объединение списка равномерных распределений. Но я не вижу, как их объединить.
Есть ли у кого-нибудь хорошие решения этой проблемы?
Редактировать: Цель состоит в том, чтобы создать мини-язык (также известный как предметно-ориентированный язык) для работы с PDF-файлами. Но сначала мне нужно разобраться в базовом представлении и алгоритмах.
Редактировать 2: dmckee предлагает реализацию гистограммы. Это то, к чему я шел с моим списком унифицированных дистрибутивов. Но я не вижу, как их объединить для создания новых дистрибутивов. В конечном итоге мне нужно найти такие вещи, как P (x
Редактировать 3: У меня есть куча гистограмм. Они распределяются неравномерно, потому что я генерирую их из данных о вхождении, поэтому, в основном, если у меня есть 100 выборок и мне нужно десять точек на гистограмме, я выделяю 10 выборок для каждого столбца и делаю столбцы переменной ширины, но постоянной площади.
Я понял, что для добавления PDF-файлов вы их сверните, и для этого я взялся за математику. Когда вы сворачиваете два равномерных распределения, вы получаете новое распределение с тремя разделами: более широкое равномерное распределение все еще там, но с треугольником, прикрепленным с каждой стороны шириной более узкого. Так что, если я сверну каждый элемент X и Y, я получу кучу таких, все перекрываются. Теперь я пытаюсь выяснить, как сложить их все, а затем получить гистограмму, которая является наилучшим приближением к ней.
Я начинаю задумываться, не была ли Монте-Карло такой плохой идеей.
Редактировать 4: В этой статье обсуждаются свертки равномерных распределений в некоторых деталях. В общем, вы получаете «трапециевидное» распределение. Поскольку каждый «столбец» в гистограммах представляет собой равномерное распределение, я надеялся, что проблему можно решить путем свертки этих столбцов и суммирования результатов.
Однако результат значительно сложнее, чем входы, а также включает в себя треугольники. Редактировать 5: [Неверный материал удален]. Но если эти трапеции приближаются к прямоугольникам с одинаковой площадью, тогда вы получите правильный ответ, и уменьшение количества прямоугольников в результате выглядит довольно просто. Это может быть решением, которое я пытался найти.
Редактировать 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
Другие операторы оставлены в качестве упражнения для читателя.