# Kernel Density Estimate for cdf

**URL:** <https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790>\
**Category:** General Usage\
**Tags:** question, statistics\
**Created:** [September 25, 2022, 10:56pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790 "2022-09-25T22:56:14Z")\
**Posts on this page:** 12\
**Page:** 1

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [September 25, 2022, 10:56pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/1 "2022-09-25T22:56:15Z")

</div>

Hi all,

I would like to use kernel density estimation to estimate the cdf from data. One way to get the cdf is to numerically integrate the estimated pdf with `quadgk`, but the drawback is it is somewhat inefficient. Is there a more efficient method?

---

<div class="post-metadata">

**Author:** ![nsajko](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nsajko/32/221187_2.png) [@nsajko](https://discourse.julialang.org/u/nsajko)\
**Post date:** [September 26, 2022, 1:02am UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/2 "2022-09-26T01:02:34Z")

</div>

For estimating the cdf from a sample you probably want to use the ecdf (empirical distribution function):

> **[Empirical distribution function](https://en.wikipedia.org/wiki/Empirical_distribution_function)**
>
> In statistics, an empirical distribution function (commonly also called an empirical cumulative distribution function, eCDF) is the distribution function associated with the empirical measure of a sample. This cumulative distribution function is a step function that jumps up by 1/n at each of the n data points. Its value at any specified value of the measured variable is the fraction of observations of the measured variable that are less than or equal to the specified value.
> The empirical distri...

It’s trivial to implement, and much more foolproof to use than KDE, because there’s no bandwidth parameter necessary.

You can find an implementation of mine here (it’s a package but not yet registered):

> <https://gitlab.com/nsajko/UnivariateProbabilityDistributionVisualizationCalculation.jl/-/blob/main/src/UnivariateProbabilityDistributionVisualizationCalculation.jl#L92-L120>

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [September 26, 2022, 12:23pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/3 "2022-09-26T12:23:44Z")

</div>

Thank you for your reply. I was considering `ecdf` in StatsBase.jl. The main downside is that it can be noticeably jagged, particularly in the tails. I don’t know if there is some smoothing that would be reasonable to apply. The main constraint would be that the smoothed curve is monotonically increasing.

---

<div class="post-metadata">

**Author:** ![nsajko](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nsajko/32/221187_2.png) [@nsajko](https://discourse.julialang.org/u/nsajko)\
**Post date:** [September 26, 2022, 12:46pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/4 "2022-09-26T12:46:40Z")

</div>

Would linear interpolation be enough?

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [September 26, 2022, 12:46pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/5 "2022-09-26T12:46:53Z")

</div>

> [@Christopher\_Fisher](#):
>
> I would like to use kernel density estimation to estimate the cdf from data.

Usually one estimates a distribution (which parametric, semi-parametric, or non-parametric methods), which then has a CDF.

Why do you need the CDF specifically? If you want to sample from the estimated distribution, there is usually a more direct way. Eg for KDE, there is a simple two-step method using the kernel.

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [September 26, 2022, 1:19pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/6 "2022-09-26T13:19:51Z")

</div>

What I am ultimately trying to do is estimate the cumulative hazard function of a model without a closed-form cdf. The cumulative hazard function is `-log(1 - F(x))`, where `F` is the cdf. What I have found with the empirical cdf is that it can be somewhat unstable in the tails. I could use something like linear interpolation like nsajko suggested or perhaps run more monte carlo simulations of the model to improve stability.

But I have been interested in exploring options with KDE. Can you tell me more about the two step method?

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [September 26, 2022, 1:33pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/7 "2022-09-26T13:33:36Z")

</div>

See eg [page 5 of these lecture notes](https://www.stat.cmu.edu/~cshalizi/350/lectures/28/lecture-28.pdf).

That said, if your tails have few observations (as tails usually do 😉) but are relevant for your results, I would recommend a parametric approach to impose some structure on it. You can make it a mixture to fit the data better: eg start with a simple function, estimate, simulate to see where the discrepancies are, then extend until it fits.

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [September 26, 2022, 1:44pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/8 "2022-09-26T13:44:20Z")

</div>

Much appreciated. Thank you both for your recommendations. I will compare the lecture notes to the [implementation](https://github.com/statsmodels/statsmodels/blob/825581cf17f1e79118592f15f49be7ad890a7104/statsmodels/nonparametric/kde.py#L202) in the Python library called StatisticsModels, but it looks like that just uses numerical integration.

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [September 26, 2022, 2:01pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/9 "2022-09-26T14:01:09Z")

</div>

Here is one more option in case someone is interested. In [this](https://rdrr.io/cran/spatstat.core/src/R/quantiledensity.R) R library, a three step process is used: estimate kernel density, use kernel density to find cdf to approximate cdf with a selected step size, apply linear interpolation.

---

<div class="post-metadata">

**Author:** ![dlakelan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlakelan/32/8491_2.png) [@dlakelan](https://discourse.julialang.org/u/dlakelan)\
**Post date:** [September 26, 2022, 2:10pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/10 "2022-09-26T14:10:07Z")

</div>

The KDE is basically a convolution of the kernel with the discrete distribution of the data. You can randomly select data points and then add a random value from the kernel to bulk up your dataset, then use the ECDF on the denser dataset to get a smoother ECDF.

```julia
bulkedup = let data = rand(37), # small dataset
     bulkdata = Float64[]
     for i in 1:10000
         push!(bulkdata,rand(data) + rand(Normal(0,1)))
     end
     bulkdata
end

```

something like that, I’m using a Normal(0,1) kernel, but use whatever kernel you like

---

<div class="post-metadata">

**Author:** ![Christopher\_Fisher](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/christopher_fisher/32/26132_2.png) [@Christopher\_Fisher](https://discourse.julialang.org/u/Christopher_Fisher)\
**Post date:** [September 26, 2022, 5:23pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/11 "2022-09-26T17:23:25Z")

</div>

After some digging, I came across the [Nelson-Aalen estimator](https://juliastats.org/Survival.jl/latest/na/#Nelson-Aalen-Estimator) for cumulative hazard functions. My original motivation was to use the cdf as a simple way to compute the cumulative hazard function indirectly. This might be useful to someone in the future.

---

<div class="post-metadata">

**Author:** ![bertschi](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bertschi/32/33462_2.png) [@bertschi](https://discourse.julialang.org/u/bertschi)\
**Post date:** [September 26, 2022, 7:34pm UTC](https://discourse.julialang.org/t/kernel-density-estimate-for-cdf/87790/12 "2022-09-26T19:34:39Z")

</div>

In case of normal kernels, the main parameter that is fitted by KDE is actually the bandwidth b, i.e., standard deviation of the kernels. The actual density is then a mixture distribution of the form

p\_{KDE}(x) = \sum\_i \frac{1}{N} \mathcal{N}(x | x\_i, b)

Thus, the following might work:

```julia
using Distributions
using KernelDensity
using Plots

x = rand(Normal(0, 1), 100)
b = KernelDensity.default_bandwidth(x)

fit1 = kde(x)
fit2 = MixtureModel(Normal.(x, b))

scatter(pdf(fit1, x), pdf.(fit2, x)) # these should be basically identical

x = sort(x)
plot(x, cdf.(fit2, x)) # Mixture model has cdf method ... note that broadcasting

```