Skip to content

Plotting Confidence Intervals into a density plot

8 messages · Elysa Mitova, David Winsemius, Jim Lemon

Original
Elysa Mitova
# 
Hi,

I am desperately looking for a way to plot confidence intervals into a
density plot of only one variable (not a scatter plot etc.)

Have you any advice how to do this?

I've only found manual ways to do with "abline", but this is a rather
bothersome method and only works with ggplot (and not ggplot2).

Thank you!
David Winsemius
# 
This makes it appear that you expect this to be done in ggplot2 automagically. I suspect you must instead first find the right approach to construction of those upper and lower bounds before plotting. It's not clear what methods you expect to be needed. Your desperation is not a guide. Perhaps trying a bit of searching?

install.packages("sos")
library(sos)
findFn("confidence intervals density estimates")


Delivers quite a few results. Then searching on the text within that webpage you find 


208	2	27	54	nprobust	kdrobust	2016-11-14 16:41:50	27	Kernel Density Estimation with Robust Confidence Intervals
209	2	27	54	nprobust	lprobust	2016-11-14 16:41:50	27	Local-Polynomial Estimation with Robust Confidence Intervals

Is that what you seek?
I know you just subscribed, so now is the time to read the Posing Guide.

== 

David Winsemius
Alameda, CA, USA
Elysa Mitova
# 
Hi, thank you!

I've constructed the upper and lower bounds with

 a <- 2.505766
 s <- 0.7789832
 n <- 607
 error <- qnorm(0.975)*s/sqrt(n)
 left <- a-error
 right <- a+error
 left
right

Now, I have the numbers I need, but I have no idea how to plot them. I was
thinking of using a polygon, but somehow it doesn't work out, because my
y-axis shows only density and is in itself not a variable?

xx <- data

fit1 <- density(data,na.rm=TRUE)

fit2 <- replicate(10000, { x <- sample(xx, replace=TRUE);
	density(x, na.rm=TRUE, from=min(fit1$x), to=max(fit1$x))$y } )

fit3 <- apply(fit2, 1, quantile, c(0.025,0.975) )  - Probably herein
lies the problem?

plot(fit1, ylim=range(fit3))
polygon( c(fit1$x, rev(fit1$x)), c(fit3[1,], rev(fit3[2,])),
col='grey', border=F)
lines(fit1)

I tried working with this solution I found on the internet, but
somehow now the lines the shaded areas sporadically everywhere around
my density plot? I just want a polygon spreading from  2.44 to 2.57
along the x-axis.


Any tipps?




On Fri, Dec 2, 2016 at 1:24 AM, David Winsemius <dwinsemius at comcast.net>
wrote:

  
  


  


      

    

    

      
      

        


  

        

      Jim Lemon
    

    

      
      #
    

  


  

      
Hi Elysa,
I think you are going a bit off course in your example. Try this and
see if it is close to what you want:

data<-rnorm(100)+runif(100,0,15)
smu_data<-supsmu(1:100,data)
rollfun<-function(x,window=10,FUN=sd) {
 xlen<-length(x)
 xout<-NA
 forward<-window%/%2
 backward<-window-forward
 for(i in 1:xlen) {
  xstart<-i-backward
  if(xstart < 1) xstart<-1
  xend<-i+forward-1
  if(xend > xlen) xend<-xlen
  xout[i]<-do.call(FUN,list(x[xstart:xend],na.rm=TRUE))
 }
 return(xout)
}
mad_data<-rollfun(data,10,mad)
plot(data,ylim=c(0,17))
library(plotrix)
dispersion(smu_data$x,smu_data$y,mad_data,type="l",interval=TRUE,
 fill="lightgray")
lines(smu_data,lwd=2)
points(1:100,data)

Jim

        
On Fri, Dec 2, 2016 at 7:18 PM, Elysa Mitova <elysa.mitova at gmail.com> wrote:

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Elysa Mitova
    

    

      
      #
    

  


  

      
Thank you,

this seems to work, but it is not exactly what I need (it indeed looks
great, but a bit beyond my understanding)

I just need a shaded area between  2.44 to 2.57 along the x-axis - a
polygon inserted into my density plot (and not a confidence line along a
scatter plot like your suggested solution)

My x-axis is an index (a data frame), my y-axis is the automatically
constructed density

        
On Fri, Dec 2, 2016 at 10:01 AM, Jim Lemon <drjimlemon at gmail.com> wrote:

        

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

  
  


  


      

    

    

      
      

        


  

        

      Jim Lemon
    

    

      
      #
    

  


  

      
In order to display a polygon, you need x/y pairs for each point. If
you just want a rectangle, you only need four x/y pairs, e.g.:

plot(0,xlim=x(2.44,2.57),ylim=c(0,1),type="n")
polygon(c(2.44,2.57,2.57,2.44),c(0,0,1,1),col="lightgray")

Now if you have a series of x values and want to display a band of
constant width around it:

y_values<-runif(14)
plot(seq(2.44,2.57,by=0.01),y_values,ylim=c(-2,3))
dispersion(seq(2.44,2.57,by=0.01),y_values,ulim=rep(0.5,14),
 type="l",interval=TRUE,col="lightgray")
lines(seq(2.44,2.57,by=0.01),y_values)

Jim

        
On Fri, Dec 2, 2016 at 8:59 PM, Elysa Mitova <elysa.mitova at gmail.com> wrote:

            

      
      
      
    

        

      
      
    

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Jim Lemon
    

    

      
      #
    

  


  

      
Hang on, maybe you mean something like this:

erupt_dens<-density(faithful$eruptions)
plot(erupt_dens,ylim=c(0,0.65))
dispersion(erupt_dens$x,erupt_dens$y,ulim=erupt_dens$y/5,
 type="l",fill="lightgray",interval=TRUE)
lines(erupt_dens)

Jim

        
On Fri, Dec 2, 2016 at 9:36 PM, Jim Lemon <drjimlemon at gmail.com> wrote:

            

      
      
      
    

        

      
      
    

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Elysa Mitova
    

    

      
      #
    

  


  

      
Hi,

sadly it does not work either, because my index (x axis) is an atomic
vector.
Error Message: $ operator is invalid for atomic vectors
I think I have to stick to displaying the confidence intervals with
straight lines (ablines) , instead of a shaded area (polygon)

Thank you so much for your help!

        
On Fri, Dec 2, 2016 at 11:45 AM, Jim Lemon <drjimlemon at gmail.com> wrote: