How to scale/transform graphics::plot() axes with any transformation, not just logarithmic (for Weibull plots)?
Asked Answered
T

1

14

I am building an R package to display Weibull plots (using graphics::plot) in R. The plot has a log-transformed x-axis and a Weibull-transformed y-axis (for lack of a better description). The two-parameter Weibull distribution can thus be represented as a straight line on this plot.

The logarithmic transformation of the x-axis is as simple as adding the log="x" parameter to plot() or curve(). How can I supply the y-axis transformation in an elegant way, so that all graphics-related plotting will work on my axis-transformed plot? To demonstrate what I need, run the following example code:

## initialisation ##
beta     <- 2;eta <- 1000
ticks    <- c(seq(0.01,0.09,0.01),(1:9)/10,seq(0.91,0.99,0.01))
F0inv    <- function (p) log(qweibull(p, 1, 1))
    # this is the transformation function
F0       <- function (q) exp(-exp(q))
    # this is the inverse of the transformation function
weibull  <- function(x)pweibull(x,beta,eta)
    # the curve of this function represents the weibull distribution 
    # as a straight line on weibull paper
weibull2 <- function(x)F0inv(weibull(x))

First an example of a Weibull distribution with beta=2 and eta=1000 on a regular, untransformed plot:

## untransformed axes ##
curve(weibull ,xlim=c(100,1e4),ylim=c(0.01,0.99))
abline(h=ticks,col="lightgray")

plot1

This plot is useless for Weibull analysis. Here is my currently implemented solution that transforms the data with function F0inv() and modifies the y-axis of the plot. Notice that I have to use F0inv() on all y-axis related data.

## transformed axis with F0inv() ##
curve(weibull2,xlim=c(100,1e4),ylim=F0inv(c(0.01,0.99)),log="x",axes=F)
axis(1);axis(2,at=F0inv(ticks),labels=ticks)
abline(h=F0inv(ticks),col="lightgray")

plot2

This works, but this is not very user-friendly: when the user wants to add annotations, one must always use F0inv():

text(300,F0inv(0.4),"at 40%")

I found that you can achieve a solution to my problem using ggplot2 and scales, but I don't want to change to a graphics package unless absolutely necessary since a lot of other code needs to be rewritten.

## with ggplot2 and scales ##
library(ggplot2)
library(scales)
weibull_trans <- function()trans_new("weibull", F0inv, F0)
qplot(c(100,1e4),xlim=c(100,1e4),ylim=c(0.01,0.99),
    stat="function",geom="line",fun=weibull) + 
    coord_trans(x="log10",y = "weibull") 

plot3

I think that if I could dynamically replace the code for applying the logarithmic transformation with my own, my problem would be solved.

I tried to find more information by Googling "R axis transformation", "R user coordinates", "R axis scaling" without useful results. Almost everything I have found dealt with logarithmic scales.

I tried to look into plot() at how the log="x" parameter works, but the relevant code for plot.window is written in C – not my strongest point at all.

Thapsus answered 8/4, 2013 at 12:29 Comment(2)
I will be interested if someone comes up with a better solution, but I think you've covered the relevant ground; I don't think you're going to do better than this -- outside of ggplot2 I don't know of any systems for generic axis transformation.Bauske
Thanks for taking the time to respond. In the mean time, I came to the same conclusion, it looks like I will be forced to use ggplot2.Thapsus
Y
1

While it doesn't appear to be possible in base graphics, you can make this function do what you want so that you can call it more simply:

F0inv    <- function (p) log(qweibull(p, 1, 1))
## this is the transformation function
F0       <- function (q) exp(-exp(q))

weibullplot <- function(eta, beta,
                        ticks=c(seq(0.01,0.09,0.01),(1:9)/10,seq(0.91,0.99,0.01)),
                        ...) {
  ## the curve of this function represents the weibull distribution 
  ## as a straight line on weibull paper
  weibull2 <- function(x)
    F0inv(pweibull(x, beta, eta))
  curve(weibull2, xlim=c(100, 1e4), ylim=F0inv(c(0.01, 0.99)), log="x", axes=FALSE)
  axis(1);
  axis(2, at=F0inv(ticks), labels=ticks)
  abline(h=F0inv(ticks),col="lightgray")
}

weibullplot(eta=1000, beta=2)
Yahiya answered 12/7, 2014 at 2:25 Comment(0)

© 2022 - 2024 — McMap. All rights reserved.