# Sample from Kernel Density Estimator
**URL:** https://discourse.julialang.org/t/sample-from-kernel-density-estimator/50639
**Category:** Statistics
**Tags:** question
**Created:** [November 23, 2020, 4:29pm UTC](https://discourse.julialang.org/t/sample-from-kernel-density-estimator/50639 "2020-11-23T16:29:14Z")
**Posts on this page:** 4
**Page:** 1
### Author:  [@heliosdrm](https://discourse.julialang.org/u/heliosdrm)
#### Post date: [November 23, 2020, 4:29pm UTC](https://discourse.julialang.org/t/sample-from-kernel-density-estimator/50639/1 "2020-11-23T16:29:14Z")
From a sample of data I can create a KDE to estimate the probability density with [KernelDensity.jl](https://github.com/JuliaStats/KernelDensity.jl), e.g.:
```julia
density_estimator = kde(data)
```
Now, is there any package or custom code that I could use to retrieve a sample from a distribution with such a density, e.g. something that allows me to do `rand(xxx(density_estimator))`?
(Or maybe other implementations of KDE which have that feature?)
---
### Author:  [@konkam](https://discourse.julialang.org/u/konkam)
#### Post date: [November 23, 2020, 10:34pm UTC](https://discourse.julialang.org/t/sample-from-kernel-density-estimator/50639/2 "2020-11-23T22:34:10Z")
I do not know of an implementation, but you can write your own with few lines. Assuming a 1-D normal kernel density estimator, to get one sample you need to:
- pick a data point `x`, using for instance with something like `x = sample(data)`
- get the kernel bandwith `h`, be it the Silverman optimal bandwith or extracted from the object created by the function `kde_lscv(data)` from the package you mention
- generate a random sample with `rand(Normal(x,1/h))`
If you change the kernel beware of the parametrisation to make sure your bandwith is a proper scale parameter, and if you go multidimensional pay attention too.
But essentially you are sampling from a mixture model with equal weights, so first sample a component then sample from that component.
---
### Author:  [@heliosdrm](https://discourse.julialang.org/u/heliosdrm)
#### Post date: [November 24, 2020, 4:23pm UTC](https://discourse.julialang.org/t/sample-from-kernel-density-estimator/50639/3 "2020-11-24T16:23:27Z")
Thanks, that works, except that in the last step the standard deviation should be `h`, not `1/h`:
```julia
using Random, KernelDensity, Distributions
data = randn(1000)
h = KernelDensity.default_bandwidth(data)
newdata = [rand(Normal(rand(data) , 1/h)) for _=1:1000]
std(newdata) # => 4.317997200903775
```
but:
```julia
newdata = [rand(Normal(rand(data) , h)) for _=1:1000]
std(newdata) # => 1.0014533279921076 (as expected)
```
---
### Author:  [@konkam](https://discourse.julialang.org/u/konkam)
#### Post date: [November 25, 2020, 3:59pm UTC](https://discourse.julialang.org/t/sample-from-kernel-density-estimator/50639/4 "2020-11-25T15:59:04Z")
You are quite right, that “beware of the parametrisation” here is pretty ironic given I used a wrong parametrisation one line above 😅