# Shading interval in density plot

**URL:** <https://discourse.julialang.org/t/shading-interval-in-density-plot/100260>\
**Category:** New to Julia\
**Tags:** question, plotting\
**Created:** [June 13, 2023, 2:26am UTC](https://discourse.julialang.org/t/shading-interval-in-density-plot/100260 "2023-06-13T02:26:10Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![eljfe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eljfe/32/50688_2.png) [@eljfe](https://discourse.julialang.org/u/eljfe)\
**Post date:** [June 13, 2023, 2:26am UTC](https://discourse.julialang.org/t/shading-interval-in-density-plot/100260/1 "2023-06-13T02:26:10Z")

</div>

Hi,  
Per the title of the post, trying to figure out how to shade a density plot interval `using StatsPlots`

I don’t really have any usable example code. A couple of days playing with this and none the wiser.

Trying to reproduce Statistical Rethinking Chapter 3 plots.

Tried the formula outlined [here](https://discourse.julialang.org/t/shade-area-under-the-curve-between-two-points-in-plots-jl/76190) which uses calls to `plot(...)` but couldn’t get it to work for `density(...)`. No clue how to approach this problem using `StatsPlots`.

The goal,

 ![Screen Shot 2023-06-12 at 10.22.56 PM](https://global.discourse-cdn.com/julialang/original/3X/7/e/7e25bbfea4ea7f00b1febe7121b42361ef91e49a.png)

thx,

---

<div class="post-metadata">

**Author:** ![nilshg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nilshg/32/2283_2.png) [@nilshg](https://discourse.julialang.org/u/nilshg)\
**Post date:** [June 13, 2023, 12:49pm UTC](https://discourse.julialang.org/t/shading-interval-in-density-plot/100260/2 "2023-06-13T12:49:26Z")

</div>

It will work the same, i.e. you can do (roughly):

```julia
density(something)
plot!(x, y, fillrange = zero(x))

```

but you will need the `y` values for the range of `x` you want to fill, so something like `pdf.(distribution, x)`

---

<div class="post-metadata">

**Author:** ![eljfe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eljfe/32/50688_2.png) [@eljfe](https://discourse.julialang.org/u/eljfe)\
**Post date:** [June 17, 2023, 1:58am UTC](https://discourse.julialang.org/t/shading-interval-in-density-plot/100260/3 "2023-06-17T01:58:45Z")

</div>

Very big thanks to you Nils. That set me on the right path.

FWIW I was able to replicate the _Statistical Rethinking_ plots after considerable homework.

The big missing piece for me was that `denstiy()` plots were based on kernel density calculations. That was novel enough. Enter `KernelDensity.jl`… Then to make it all work `using Plots.jl`.

It’s no wonder the graphic appears in an introductory chapter of _SR_ without explanation.

Here is the code, and the plot itself.

```julia
using Distributions
using KernelDensity
using StatsPlots

p_grid = collect(range(0,1,10^3))
prior = ones(10^3)
likelihood = pdf.(Binomial.(9, p_grid), 6)
posterior = likelihood .* prior
posterior = posterior ./ sum(posterior)
samples = sample( p_grid, Weights(posterior), 10^4)

ksamp = kde(samples)
x = collect(ksamp.x)
y = ksamp.density

xi = findall( i -> i<0.5, x)
d11 = density(samples, xlabel="proportion water (p)", ylabel="density", label="kde", xticks=0.0:0.25:1, legend=:topleft)
plot!(x[xi], y[xi], fillrange=zeros(1), fc=:blues, label="p<0.5")

xi = findall(i -> (0.5<i) & (i<0.75), x)
d12 = density(samples, xlabel="proportion water (p)", ylabel="density", label="kde", xticks=0.0:0.25:1, legend=:topleft)
plot!(x[xi], y[xi], fillrange=zeros(1), fc=:blues, label="0.5<p<0.75")

qa, qb = quantile(samp, (0.0, 0.8))
xi = findall(i -> (qa<=i) & (i<=qb), x)
d21 = density(samples, xlabel="proportion water (p)", ylabel="density", label="kde", xticks=0.0:0.25:1, legend=:topleft)
plot!(x[xi], y[xi], fillrange=zeros(1), fc=:blues, label="lower 80%")

qa, qb = quantile(samp, (0.1, 0.9))
xi = findall(i -> (qa<=i) & (i<=qb), x)
d22 = density(samples, xlabel="proportion water (p)", ylabel="density", label="kde", xticks=0.0:0.25:1, legend=:topleft)
plot!(x[xi], y[xi], fillrange=zeros(1), fc=:blues, label="middle 80%")

plot(d11, d12, d21, d22, layout=(2, 2), size=(800,700))

```

 ![Screen Shot 2023-06-16 at 9.55.06 PM](https://global.discourse-cdn.com/julialang/original/3X/9/7/97fef96734103476851ac51525ba6442efbfc907.png)
