Transforms the density function of a continuous random variable into a discrete probability distribution with minimal distortion using the Lloyd-Max algorithm.
Arguments
- density_fn
probability density function.
- n_levels
cardinality of the set of all possible outcomes.
- eps
convergence threshold for the algorithm.
Value
A list containing:
- prob
discrete probability distribution.
- endp
endpoints of intervals that partition the continuous domain.
- repr
representation points of the intervals.
- dist
distortion measured as the mean-squared error (MSE).
Details
The function addresses the problem of transforming a continuous random variable \(X\) into a discrete random variable \(Y\) with minimal distortion. Distortion is measured as mean-squared error (MSE): $$ \text{E}\left[ (X - Y)^2 \right] = \sum_{k=1}^{K} \int_{x_{k-1}}^{x_{k}} f_{X}(x) \left( x - r_{k} \right)^2 \, dx $$ where:
- \(f_{X}\)
is the probability density function of \(X\),
- \(K\)
is the number of possible outcomes of \(Y\),
- \(x_{k}\)
are endpoints of intervals that partition the domain of \(X\),
- \(r_{k}\)
are representation points of the intervals.
This problem is solved using the following iterative procedure:
- \(1.\)
Start with an arbitrary initial set of representation points: \(r_{1} < r_{2} < \dots < r_{K}\).
- \(2.\)
Repeat the following steps until the improvement in MSE falls below given \(\varepsilon\).
- \(3.\)
Calculate endpoints as \(x_{k} = (r_{k+1} + r_{k})/2\) for each \(k = 1, \dots, K-1\) and set \(x_{0}\) and \(x_{K}\) to \(-\infty\) and \(\infty\), respectively.
- \(4.\)
Update representation points by setting \(r_{k}\) equal to the conditional mean of \(X\) given \(X \in (x_{k-1}, x_{k})\) for each \(k = 1, \dots, K\).
With each execution of step \((3)\) and step \((4)\), the MSE decreases or remains the same. As MSE is nonnegative, it approaches a limit. The algorithm terminates when the improvement in MSE is less than a given \(\varepsilon > 0\), ensuring convergence after a finite number of iterations.
This procedure is known as Lloyd-Max's algorithm, initially used for scalar quantization and closely related to the k-means algorithm. Local convergence has been proven for log-concave density functions by Kieffer. Many common probability distributions are log-concave including the normal and skew normal distribution, as shown by Azzalini.
References
Azzalini, A. (1985). A class of distributions which includes the normal ones. Scandinavian Journal of Statistics 12(2), 171–178.
Kieffer, J. (1983). Uniqueness of locally optimal quantizer for log-concave density and convex error function. IEEE Transactions on Information Theory 29, 42–47.
Lloyd, S. (1982). Least squares quantization in PCM. IEEE Transactions on Information Theory 28 (2), 129–137.
Examples
discretize_density(density_fn = stats::dnorm, n_levels = 5)
#> $prob
#> [1] 0.1063142 0.2444957 0.2983803 0.2444957 0.1063142
#>
#> $endp
#> [1] -Inf -1.2463707 -0.3831349 0.3831349 1.2463707 Inf
#>
#> $repr
#> [1] -1.7264716 -0.7662698 0.0000000 0.7662698 1.7264716
#>
#> $dist
#> [1] 0.07994185
#>
discretize_density(density_fn = function(x) {
2 * stats::dnorm(x) * stats::pnorm(0.5 * x)
}, n_levels = 4)
#> $prob
#> [1] 0.1647781 0.3380772 0.3359810 0.1611637
#>
#> $endp
#> [1] -Inf -0.5537368 0.3597728 1.2808536 Inf
#>
#> $repr
#> [1] -1.04423728 -0.06323628 0.78278182 1.77892538
#>
#> $dist
#> [1] 0.1026681
#>