Using interaction term for describing non-linear relationship in regression model

Question

I have created a Rscript for illustrating my question:

library(epicalc)

temp.df <- data.frame(x=1:100,y=log(1:100))
set.seed(1)
temp.df[,"y"] <- temp.df[,"y"]*sample((10000:10100)/10000,size=100,replace=TRUE)
temp.df[temp.df[,"x"]<=15,"x_gp"] <- "<=15"
temp.df[temp.df[,"x"]>15,"x_gp"] <- ">15"
temp.df[,"x_gp"] <- factor(temp.df[,"x_gp"],levels=c("<=15",">15"))

glm.null <- glm(y~NULL,data=temp.df)
glm.01 <- glm(y~x,data=temp.df)
glm.02 <- glm(y~log(x),data=temp.df)
glm.03 <- glm(y~x*x_gp,data=temp.df)

lrtest(glm.null,glm.01)
lrtest(glm.null,glm.02)
lrtest(glm.null,glm.03)
lrtest(glm.01,glm.02)
lrtest(glm.02,glm.03)
lrtest(glm.01,glm.03)

(logLik(glm.02)*-2)-(logLik(glm.03)*-2)

plot(temp.df[,c("x","y")],col="black",pch=15)
points(temp.df[,"x"],predict(glm.01),col="blue",pch=20)
points(temp.df[,"x"],predict(glm.02),col="red")
points(temp.df[,"x"],predict(glm.03),col="green",pch=16)

The dataframe is describing a log-linear relationship. But in reality, maybe it is not log-linear, even if it is log-linear, it may not be preferred to use log, as it may be more difficult to understand.

I thought of "bending" the linear regression into two "segments" (glm.03) using the interaction term. I am aware that it is not as good as log in my example but already a lot better than using a straight line (glm.02) alone.

I wonder if my method will be too unorthodox in the field of public health, so can anyone suggest other ways for describing a log-linear-like relation, without using log? Thanks.

Answer 1

Yeah, how about a nonlinear regression model. There are a few models that are used all of the time in biostats. I used the Michaelis-Menten model for somebody one time: http://en.wikipedia.org/wiki/Michaelis%E2%80%93Menten_kinetics.

Looking at the plot in your script, this model would be a great fit. Would this be accepted in public health?

The function is smooth - it looks piecewise because I only plotted it at the fitted values. I used the nls function in R to estimate this model:

  model1 = nls(DATA$R ~ ((a * DATA$CONC) / (b + DATA$CONC)), start = list(a = 100, b = -1), algorithm = "port")
  plot(DATA$CONC, DATA$R)
  lines(DATA$CONC, fitted.values(model1), lty = 1)