1

I'm trying to plot 3 non-linear models in ggplot2. It's working in automatic scale but not in log10 scale where I get the "singular gradient error". What could be the problem(s)?

The data I'm trying to fit (df1):

x   y
4.17    0.55
10.08   0.48
40.25   0.38
101.17  0.32
400.33  0.24

The code I tried:

plot <- ggplot(df1, aes(x=x, y=y))+
  stat_smooth(method="nls",
              formula=y~I(a*x^(-n)),
              data=df1,
              start=list(a=1,n=1),
              se=FALSE,
              colour="red")+
  stat_smooth(method="nls",
              formula=y~m*(x+m^(1/n))^(-n),
              data=df1,
              start=list(m=0.7, n=0.17),
              se=FALSE,
              colour="blue")+
  stat_smooth(method="nls",
              formula=y~m*(x+m^(1/n))^(-n)+b,
              data=df1,
              start=list(m=0.7, n=0.17, b=1),
              se=FALSE,
              colour="green")+
  geom_point()+
  scale_x_log10()+
  scale_y_log10()+
  theme_bw()
plot
user2165907
  • 1,401
  • 3
  • 13
  • 28

1 Answers1

3

The problem seems to be that when you specify scale_x_log10 or scale_y_log10, the values of your data are transformed before being passed along to the different stats or geoms. This means while your nls may work on the untransformed data, it does not work on the log-transformed data.

#OK
nls(y~m*(x+m^(1/n))^(-n), df1, start=list(m=0.7, n=0.17))
#NOT OK
nls(y~m*(x+m^(1/n))^(-n), log10(df1), start=list(m=0.7, n=0.17))

There doesn't seem to be much you can do in ggplot2 to fix this. Instead, you could fit the NLS models ahead of time on the untransformed scale and just plot the results with ggplot2. For example

mods<-list(
    list(y~I(a*x^(-n)), list(a=1,n=1)),
    list(y~m*(x+m^(1/n))^(-n), list(m=0.7, n=0.17)),
    list(y~m*(x+m^(1/n))^(-n)+b, list(m=0.7, n=0.17, b=1))
)

fits<-lapply(mods, function(x, xr) {
    mod<-nls(x[[1]], data=df1, start=x[[2]])
    xx<-seq(xr[1], xr[2], length.out=100)
    yy<-predict(mod, newdata=data.frame(x=xx))
    data.frame(x=xx, y=yy)
}, xr=range(df1$x))

library(ggplot2)
ggplot(df1, aes(x=x, y=y))+
  geom_line(data=fits[[1]], color="red") +
  geom_line(data=fits[[2]], color="blue") +
  geom_line(data=fits[[3]], color="green") +
  geom_point()+
  scale_x_log10()+
  scale_y_log10()+
  theme_bw()

will produce

enter image description here

MrFlick
  • 195,160
  • 17
  • 277
  • 295
  • Changing 'length.out' in the definition of 'xx' changes the plotted model... How do you choose the "best" value for this parameter? – user2165907 Jul 28 '14 at 07:28
  • Bigger is better, but usually you just end up plotting points you'll never see. Trial and error I guess. – MrFlick Jul 28 '14 at 12:22