[R] how to change one of the axis to normal probability scale

[Charles Annis, P.E.]

Dear C. Chang:

It's very helpful to draw your own grid in situations like this.  I have
included the code for mine below.

First some disclaimers:

1) The probability axis goes from p=1.e-7 to p=1-1.e-7, ridiculous
extrapolations by most accounts.  However, the engineers with whom I work do
this extrapolation, often unthinkingly.  At least plotting on this grid
forces them to consider what they are doing and how much (or little) data
they have to substantiate it.  You can fix this easily in your routine, by
only considering, say, p=0.001, to p=0.999.

2) The grid looks great as drawn but scales poorly because the location of
labels etc. are absolute rather than relative or some cobbled combination of
these.  Since it does what I need, I use it.  In your grid you may want to
make it easier to scale.

3) Notice that the plot is set up like this:
plot(NA, NA, xlim=c(3, 7), ylim=c(z.min, z.max), ...
so that you will be plotting on a Cartesian grid in x-units of log10() and
y-units of qnorm()

Anyway here's the code.  Watch out for superfluous carriage-returns inserted
by the mail handler.

#####  R code to draw a lognormal CDF grid.   #####
lognormal.CDF.fn <- function(x.axis.title="Nf, Cycles"){
#   lognormal CDF grid
# First draw the grid with no points.  
# The x-axis is log cycles.  The y-axis, although labled probability, is
ploted as normal z units.
z.min <- -5.5
z.max <- -z.min
z.norm <- seq(z.min, z.max, length=101)
plot(NA, NA, xlim=c(3, 7), ylim=c(z.min, z.max), type="n", xaxt="n",
yaxt="n", frame = FALSE, xlab = "", ylab = "Cumulative Probability (Fraction
less than N)")
#
# Draw the x-axis
axis(side = 1, labels = FALSE, at = c(3, 4, 5, 6, 7), line = 0., tick=TRUE,
outer = FALSE)
text(x=c(3, 4, 5, 6, 7),      y=rep(z.min-1.0, 5), rep("10", 5), xpd = NA,
cex=1.)
text(x=c(3, 4, 5, 6, 7)+ 0.1, y=rep(z.min-0.8, 5)+0.08, c("3", "4", "5",
"6", "7"), xpd = NA, cex=0.8)
text(x=5.05, y=z.min-1.6, x.axis.title, xpd=NA, cex=1.)
# Draw the interior log tick marks
for (i in 3:6){
axis(side=1, at=log10(c(2, 3, 4, 5, 6, 7, 8, 9))+ i,   cex=1, labels=NA,,
line = 0, tck = -0.01)
text(x=log10(c(2, 3, 4, 5, 6, 7, 8))+ i, y=rep(z.min-0.7, 7), c("2", "3",
"4", "5", "6", "7", "8"), xpd = NA, cex=0.7)
}
#
# Draw the y-probability axis
probs <- c(1e-7, 1e-6, 1e-5, 1e-4, 0.001, 0.002, 0.005, 0.01, 0.02, 0.05,
0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.95, 0.98, 0.99, 0.995, 0.998,
0.999, 0.9999, 0.99999, 0.999999, 0.9999999)
z.vals <- qnorm(probs)
axis(side=2, at=z.vals , labels=NA, line = 0, tck = -0.01)
text(x=rep(3.-0.4, 4),  y= -0.1+qnorm(c(1e-7, 1e-6, 1e-5, 1e-4)), cex=0.8,
xpd = NA, labels=rep("10", 6))
text(x=rep(3.-0.3  , 4),  y=  0.1+qnorm(c(1e-7, 1e-6, 1e-5, 1e-4)), cex=0.6,
xpd = NA,labels=c("-7", "-6", "-5", "-4"))
text(x=rep(3.-0.25, 6),  y= qnorm(c(0.001, 0.01, 0.1, 0.5, 0.9, 0.99,
0.999)), cex=0.7,  xpd = NA,labels=c("0.001", "0.01", "0.1", "0.5", "0.9",
"0.99", "0.999"), adj=1)
text(x=rep(3.-0.4, 4),  y= -0.1-qnorm(c(1e-7, 1e-6, 1e-5, 1e-4)), cex=0.8,
xpd = NA, labels=rep("1-10", 6))
text(x=rep(3.-0.27  , 4),  y=  0.1-qnorm(c(1e-7, 1e-6, 1e-5, 1e-4)),
cex=0.6,  xpd = NA,labels=c("-7", "-6", "-5", "-4"))
box()
}
### Test it
lognormal.CDF.fn()
###################################################################


Charles Annis, P.E.
 
Charles.Annis at StatisticalEngineering.com
phone: 561-352-9699
eFax:  503-217-4849
http://www.StatisticalEngineering.com

-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of C Chang
Sent: Tuesday, February 03, 2004 3:51 PM
To: r-help at stat.math.ethz.ch
Subject: [R] how to change one of the axis to normal probability scale

Hi,
 
I would like to plot a graph with one axis in log scale and the other in
normal probability scale using R. I cannot change the axis scale to
probability scale. What can be a solution? Thanks.
 
CT

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide!
http://www.R-project.org/posting-guide.html