[R] Help resolving issues with generating a chi-squared density plot from scratch

[Rolf Turner]

I am far too lazy to go through all your rather complicated code, but my 
basic impression is that you are re-inventing quite a few wheels.

Just to get a plot of the density function you can simply do, e.g.:

	plot(function(x){dchisq(x,1)},0,qchisq(0.999,1))

I can't see why you are fooling about with the ylim values; just let the 
plot() function choose ylim.

As to why you can't get things to work when df=1 ... well don't try to 
set the upper y limit equal to infinity! You have a finite plotting 
region, after all.

I have no idea what you mean by "The y-axis is numbered only 
relatively"; this makes no sense at all.  What *do* you want?  The 
y-axis labelling that you get will be/is the appropriate labelling.

The arrows will go where you tell them to go, so if they are "going 
diagonal" you are telling them to go diagonal.

cheers,

Rolf Turner

On 10/03/14 13:30, Levi Robinson wrote:
> I wrote the code to graph a chi-squared density function, shade the
> percentile, and point to the CV, but it has a few issues I can't seem to
> resolve
>
> 1. It won't work at all for DF = 1 due to ylim going to infinity, but I
> haven't been able to resolve this still after hours of trying.
> 2) The y-axis is numbered only relatively; I'd prefer they were the actual
> prob densities, but again, I fiddled with a few things, but it just
> wouldn't get me what I wanted.
> 3) For low degrees of freedom and higher percentiles, the arrow pointing to
> CV seems to end up going diagonal instead of straight down
>
> Here's the code below and here's the URL for R fiddle for the code which
> might make it easier to fix:
> http://www.r-fiddle.org/#/fiddle?id=ChFi0dyJ&version=4
>
>
> chi.dens = function(x = NULL, df = NULL, cv = NULL) {
>
>          # x = percentile/quantile
>          # df = degrees of freedom
>          # cv = critical value
>
>          if(x>1 ||x<0) stop("Percentile must be between 0 and 1") #
> Error-handling
>
>
>          qend = qchisq(x, df)
>          perc = x
>
>          qt0=qchisq(0.5, df)   # Defining for arrows later
>          dy0=dchisq(0.45, df)  # Defining for arrows later
>
>          xrange = qchisq (0.999, df)
>          x = seq(0, xrange)
>          y = dchisq(x, df)
>          yheight = max(dchisq(x, df))
>
>          # Creating plot
>          plot(x,y,type = "l", ylim=c(0,yheight),axes=FALSE)
>
>          title( "Chi-squared Density Curve with")
>          mtext(bquote(paste("DF = ", .(df), " and Percentile = ", .(perc))),
> side = 3, line = 0) # Input information
>
>          # Shading left tail-region
>
>          qt = signif(qend,5)
>          x0=seq(0, qt)
>          y0=dchisq(x0, df)
>          polygon(c(0, x0, qt), c(0, y0, 0), col = "lightblue")
>          xtks=signif(seq(0,xrange,length=10),3)
>          axis(1, pos=0, at=xtks, labels=xtks)
>          y.unit=max(dchisq(x, df))/5
>          y.pos=seq(0,5*y.unit, length=5)
>          y.lab=c(0, 1, 2, 3, 4) # Y axis numbers only reflect relative
> densities to each other.
>          axis(2,pos=0, at=y.pos, labels=y.lab)        # set up y-axis
>
>          # "Normal" CV less than the 99.9 Percentile:
>
>    if(qt <= xrange){
>       if(length(cv)==0) text(0.75*xrange,3*y.unit, bquote(paste("CV = ",
> .(qend))), cex=1.2, col = "red",adj=0.5)
>       #
>       if(length(perc)==1)  text(0.35*xrange,3*y.unit, paste("Area = ",
> perc), cex=1.2, col="darkblue", adj=0.5)
>       if(perc>0.5){
>          arrows(0.35*xrange, 2.5*y.unit, qt0, 0.3*dy0,  length=0.1,
> angle=20)     # pointing to the shaded region
>          }
>       if(perc<=0.5){
>          arrows(0.35*xrange, 2.5*y.unit, qt/2, 0.3*dy0,  length=0.1,
> angle=20)     # pointing to the shaded region
>        }
>      arrows(0.75*xrange, 2.5*y.unit, qt, 0, length=0.1, angle=20)
>      points(qt,0,pch=17, col="red")
>     }
>
>          # When CV is greater than the 99.9 Percentile:
>
>   if(qt > xrange){
>
>       print("CV may be too far to the right to be shown on graph due to the
> high percentile")
>
>       if(length(cv)==0) text(0.75*xrange,3*y.unit, bquote(paste("CV = ",
> .(qend))), cex=1.2, col = "red",adj=0.5)
>
>       if(length(perc)==1)  text(0.35*xrange,3*y.unit, paste("Area = ",
> perc), cex=1.2, col="darkblue", adj=0.5)
>       if(perc>0.5){
>          arrows(0.35*xrange, 2.5*y.unit, qt0, 0.3*dy0,  length=0.1,
> angle=20)     # pointing to the shaded region
>          }
>       if(perc<=0.5){
>          arrows(0.35*xrange, 2.5*y.unit, qt/2, 0.3*dy0,  length=0.1,
> angle=20)     # pointing to the shaded region
>        }
>      arrows(0.75*xrange, 2.5*y.unit, qt, 0, length=0.1, angle=20)
>      points(qt,0,pch=17, col="red")
>    }
>
> }