[R] plotmath vector problem; full program enclosed

Charles C. Berry cberry at tajo.ucsd.edu
Tue Jul 6 21:54:36 CEST 2010


On Tue, 6 Jul 2010, David Winsemius wrote:

>
> On Jul 6, 2010, at 1:41 PM, Duncan Murdoch wrote:
>
>> On 06/07/2010 10:54 AM, Paul Johnson wrote:
>> > Here's another example of my plotmath whipping boy, the Normal 
>> > distribution.
>> > 
>> > A colleague asks for a Normal plotted above a series of axes that
>> > represent various other distributions (T, etc).
>> > 
>> > I want to use vectors of equations in plotmath to do this, but have
>> > run into trouble.  Now I've isolated the problem down to a relatively
>> > small piece of working example code (below).  If you would please run
>> > this, I think you will see the problem.  When plotmath meets one
>> > vector of expressions, it converts all but one to math, so in the
>> > figure output i get, in LaTeX speak
>> > 
>> > b1         $\mu-1.0 \sigma$    $\mu$
>> > 
>> > All values except the first come out correctly.
>> > 
>> > This happens only when I try to use bquote or substitute to get
>> > variables to fill in where the 1.96, 1.0, and so forth should be.  In
>> > the figure output, you should see a second axis where all of the
>> > symbols are resolved correctly.
>> > 
>> > As usual, thanks in advance for your help, sorry if I've made an
>> > obvious mistake or overlooked a manual.
>> > 
>> > ### Filename: plotMathProblem.R
>> > ### Paul Johnson July 5, 2010
>> > ### email me <pauljohn at ku.edu>
>> > 
>> >  sigma <- 10.0
>> >  mu <- 4.0
>> > 
>> >  myx <- seq( mu - 3.5*sigma,  mu+ 3.5*sigma, length.out=500)
>> > 
>> >  myDensity <- dnorm(myx,mean=mu,sd=sigma)
>> > 
>> >  ### xpd needed to allow writing outside strict box of graph
>> >  ### Need big bottom margin to add several x axes
>> >  par(xpd=TRUE, ps=10, mar=c(18,2,2,2))
>> > 
>> > plot(myx, myDensity, type="l", xlab="", ylab="Probability Density ",
>> > main=myTitle1, axes=FALSE)
>> >  axis(2, pos= mu - 3.6*sigma)
>> >  axis(1, pos=0)
>> > 
>> >  lines(c(myx[1],myx[length(myx)]),c(0,0)) ### closes off axes
>> > 
>> > 
>> >  addInteriorLine <- function(x, m, sd){
>> >    for (i in 1:(length(x))){
>> >      lines( c(x[i],x[i]), c(0, dnorm(x[i],m=m,sd=sd)), lty= 14, lwd=.2)
>> >    }
>>> }
>> > 
>> > 
>> >  dividers <- c(qnorm(0.025), -1, 0, 1, qnorm(0.975))
>> >  addInteriorLine(mu+sigma*dividers, mu,sigma)
>> > 
>> >  # bquote creates an expression that text plotters can use
>> >  t1 <-  bquote( mu== .(mu))
>> >  mtext(bquote( mu == .(mu)), 1, at=mu, line=-1)
>> > 
>> > 
>> >  addInteriorLabel <- function(pos1, pos2, m, s){
>> >    area <- abs(100*( pnorm(m+pos1*s,m,s)-pnorm(m+pos2*s, m,s)))
>> >    mid <- m+0.5*(pos1+pos2)*s
>> >    text(mid, 0.5*dnorm(mid,m,s),label=paste(round(area,2),"%"))
>>> }
>> > 
>> > 
>> >  addInteriorLabel(dividers[1],dividers[2],  mu, sigma)
>> >  addInteriorLabel(dividers[2],dividers[3],  mu, sigma)
>> >  addInteriorLabel(dividers[3],dividers[4],  mu, sigma)
>> >  addInteriorLabel(dividers[4],dividers[5],  mu, sigma)
>> > 
>> > 
>> > 
>> > 
>> > ### Following is problem point: axis will
>> > ### end up with correct labels, except for first point,
>> > ### where we end up with "b1" instead of "mu - 1.96*sigma".
>> >   b1 <- substitute( mu - d*sigma, list(d=round(dividers[1],2)) )
>> >   b2 <- substitute( mu - sigma )
>> >   b3 <- substitute( mu )
>> >   b4 <- substitute( mu + sigma )
>> >   b5 <- substitute( mu + d*sigma, list(d=round(dividers[5],2)) )
>> > ##   plot(-20:50,-20:50,type="n",axes=F)
>> >   axis(1, line=4,at=mu+dividers*sigma,
>> > labels=c(expression(b1),b2,b3,b4,b5), padj=-1)
>> > 
>> 
>> You want "as.expression(b1)", not "expression(b1)".  The latter means "the 
>> expression consisting of the symbol b1".  The former means "take the object 
>> stored in b1, and convert it to an expression.".
>> 
>> It's not perfect, because you'll end up with "mu - -1.96sigma" (i.e. two 
>> minus signs), but it's closer than what you had.
>
> Easily addressed in this case with "~" instead of "-". The value of "d" 
> provides the minus:
>
> b1 <- substitute( mu ~ d*sigma, list(d=round(dividers[1],2)) )
>

But if d >= 0, there will be no sign so maybe use

b1 <- substitute( mu ~ d*sigma, list(d=sprintf("%+.2f", dividers[1]))) 
b5 <- substitute( mu ~ d*sigma, list(d=sprintf("%+.2f", dividers[5])))

etc.

Chuck

>> 
>> Duncan Murdoch
>> 
>> > 
>> > 
>> > 
>> > ### This gets "right result" but have to hard code the dividers
>> >  b1 <- expression( mu - 1.96*sigma )
>> >  b2 <- expression( mu - sigma )
>> >  b3 <- expression( mu )
>> >  b4 <- expression( mu + sigma )
>> >  b5 <- expression( mu + 1.96*sigma )
>> >  axis(1, line=8,at=mu+dividers*sigma, labels=c(b1,b2,b3,b4,b5), padj=-1)
>> > 
>> > 
>> > 
>> > 
>> > 
>> 
>> ______________________________________________
>> R-help at 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.
>
> David Winsemius, MD
> West Hartford, CT
>
> ______________________________________________
> R-help at 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.
>

Charles C. Berry                            (858) 534-2098
                                             Dept of Family/Preventive Medicine
E mailto:cberry at tajo.ucsd.edu	            UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/  La Jolla, San Diego 92093-0901



More information about the R-help mailing list