Chung-hong Chan <chainsawtiney <at> gmail.com> writes: > I have plot the PDF of t distribution with df = 74. > curve(dt(x,df=74),from=-4, to=4) > > how can I shade the area under curve (for example, col="red") from t=+- 1.996?
If you could used lattice plots instead, below is a modified version of the standard lattice densityplot. You can use it with a call like: densityplot(form, data=DHL, shadelimit=shadelimit, panel=panel.shadeddensityplot, shadecol='#FFEECC') I am sure Deepayan will have a more elegant solution, but this worked for me. Dieter ------ # Densityplot with limits panel.shadeddensityplot = function (x, darg = list(n = 100), plot.points = "jitter", ref = FALSE, groups = NULL, jitter.amount = 0.01 * diff(current.panel.limits()$ylim), type = "p",shadelimit=0, shadecol="red", ...) { x=na.omit(x) if (ref) { reference.line <- trellis.par.get("reference.line") panel.abline(h = 0, col = reference.line$col, lty = reference.line$lty, lwd = reference.line$lwd) } plot.line <- trellis.par.get("plot.line") superpose.line <- trellis.par.get("superpose.line") if (!is.null(groups)) { panel.superpose(x, darg = darg, plot.points = plot.points, ref = FALSE, groups = groups, panel.groups = panel.densityplot, jitter.amount = jitter.amount, type = type, ...) } else { if (sum(!is.na(x)) > 1) { h <- do.call("density", c(list(x = x), darg)) lim <- current.panel.limits()$xlim id <- h$x > min(lim) & h$x < max(lim) panel.lines(x = h$x[id], y = h$y[id], ...) idp = h$x> min(lim) & h$x < shadelimit xp = h$x[idp] xp = c(xp[1],xp,xp[length(xp)]) yp = h$y[idp] yp = c(0,yp,0) panel.polygon(x = xp, y = yp,col=shadecol, ...) grid.text(x=0.01,y=0.99,default.units="npc", gp=gpar(fontsize=10),label= paste("n = ",length(x[x<=shadelimit]),"/",length(x),sep=""), just=c("left","top")) } switch(as.character(plot.points), "TRUE" = panel.xyplot(x = x, y = rep(0, length(x)), type = type, ...), rug = panel.rug(x = x, start = 0, end = 0, x.units = c("npc", "native"), type = type, ...), jitter = panel.xyplot(x = x, y = jitter(rep(0, length(x)), amount = jitter.amount), type = type, ...)) } } ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.