2
$\begingroup$

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.

$\endgroup$

1 Answer 1

4
$\begingroup$

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?

enter image description here

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)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.