5
$\begingroup$

I am trying to plot model fit curve. I am not able to get a smoother fit curve. Any advice?

XY <- data.frame(cbind(Values = c(91.8, 95.3,   99.8,   123.3,  202.9,  619.8,  1214.2, 1519.1, 1509.2, 1523.3, 1595.2, 1625.1),
            Concn = c(1000, 300,    100,    30, 10, 3,  1,  0.3,    0.1,    0.03,   0.01,   0)))
nls.fit <- nls(Values ~ (ymax* Concn / (ec50 + Concn)) + Ns*XY$Concn + ymin, data=XY,
				   start=list(ymax=max(XY$Values), ymin = min(XY$Values), ec50 = 3, Ns = 0.2045514))
plot(XY$Values ~ XY$Concn , data = XY, col = 4,
     main = "XY Std curve", log = "x")
lines(XY$Concn, predict(nls.fit))

enter image description here

$\endgroup$
  • $\begingroup$ Try estimating a cubic relationship, i.e., using OLS and adding squared and cubed versions of your independent variable. $\endgroup$ – rolando2 Mar 12 '17 at 17:07
  • $\begingroup$ You might do better on a programming site if you do not get any useful responses here. $\endgroup$ – mdewey Mar 12 '17 at 17:10
  • $\begingroup$ I'm doubtful about my suggestion. Any curve I get is still jagged. $\endgroup$ – rolando2 Mar 12 '17 at 17:25
  • 3
    $\begingroup$ @mdewey This is not a problem that will be solved by writing a better program. It's also not going to be solved by throwing polynomial terms into the mix. The method I describe at stats.stackexchange.com/a/83613/919 will work very well to fit the data, provided a more suitable functional form is adopted. In other words, I recommend two fundamental changes: replace the rational function with one that is supported by theory (most likely predicting an exponential decay) and replace nls by a loss function that forces a close fit near the smallest and largest responses. $\endgroup$ – whuber Mar 12 '17 at 17:28
5
$\begingroup$

This seems to be a case of dose-response modelling. There is an excellent paper by Ritz et al. (2015) that describes how these analyses can be performed using R. An introduction is provided here.

Using your data and the R package drc (which is the package described in the paper by Ritz et al.), I fitted a 4-parameter log-logistic function to the data.

# Load package
library(drc)

# The data
XY <- structure(list(Values = c(91.8, 95.3, 99.8, 123.3, 202.9, 619.8, 
1214.2, 1519.1, 1509.2, 1523.3, 1595.2, 1625.1), Concn = c(1000, 
300, 100, 30, 10, 3, 1, 0.3, 0.1, 0.03, 0.01, 0)), .Names = c("Values", 
"Concn"), class = "data.frame", row.names = c(NA, -12L))    

# Fit a four-parameter log-logistic function to the data
fit <- drm(Values~Concn, data = XY, fct = LL.4())

# Plot the fit
plot(fit, type = "confidence", broken = TRUE, col = "grey50", lwd = 2)
plot(fit, type = "obs", broken = TRUE, pch = 1, lwd = 2, col = "blue", add = TRUE)

Fitted Dose-Response Curve

The fit looks pretty good to me.

$\endgroup$
  • 2
    $\begingroup$ Before fitting a curve, one has to specify the goal. If the goal is to fit a theoretically valid (or reasonable) model to get back parameters that make sense, then this answer is the only way to go. If the goal is just to get a curve that sort of kind of looks smooth, with no regard to theory or parameters, then there are lots of ways to get acceptable curves. $\endgroup$ – Harvey Motulsky Mar 12 '17 at 19:43
  • $\begingroup$ @HarveyMotulsky Thanks Harvey, that's a good point. As the goal was not specifically stated by the OP, I just conjectured that the goal was of the first type that you describe (i.e. fit a model to estimate parameters). I leave the answer to the second goal to others (I think that BenOgorek's answer is along those lines). $\endgroup$ – COOLSerdash Mar 12 '17 at 19:50
  • 1
    $\begingroup$ The goal was to use this Standard Curve to get back parameters and use it to plot drug response curves. @COOLSerdash has helped me create the curve that is an aesthetically appealing plot as well. Thanks a lot. $\endgroup$ – RanonKahn Mar 13 '17 at 1:14
3
$\begingroup$

You've plotted your curve as a linear interpolation between eleven data points. Even if the true underlying curve smooth, drawing it by taking eleven sample points and interpolating linearly is going to look pointy.

You need more sample points when drawing the curve. Create a sequence of x-values to use as sample points:

x <- seq(from = min(XY$Values), to = max(XY$Values), length.out = 250)

My go to is to create 250 sample points. Then feed these into your predict function and you will get a (within human perception) smooth rendering of the curve.

$\endgroup$
  • 1
    $\begingroup$ The problem isn't the "pointiness"; it's the curving of the fit near the most extreme values of the response, especially for low responses. So long as the functional form is retained (which forces an eventual linear increase in responses), this problem will persist. $\endgroup$ – whuber Mar 12 '17 at 17:29
3
$\begingroup$

I also think the original plot is smooth to some extent already. Another approach is to have ggplot do the smoothing for you. Not sure if this is what you want.

ggplot(XY, aes(x=log(Concn), y = Values)) + geom_smooth(method="loess")

enter image description here

$\endgroup$
  • 3
    $\begingroup$ While elegant I think your suggestion hides so much under the carpet that is almost damaging. Think of someone asking the OP what is shown or why is this smoothed curve shown and not another? In addition if anything the choice of a LOESS smoother here introduces some very strong edge artefacts.. The smooth line suggests that actually the trend is reversing back towards central values as we go outside of the current sample support; this is something implausible given the the raw data shown. One should be more careful with "canned solutions"... $\endgroup$ – usεr11852 Mar 12 '17 at 21:14
1
$\begingroup$

Your function looks pretty smooth to me, and what I think you want is more flexibility.

So here's a spline with 8 degrees of freedom:

library(splines)
my_glm <- glm(Values ~ ns(Concn, df = 8), data = XY)
plot(XY$Values ~ XY$Concn , data = XY, col = 4,
     main = "XY Std curve", log = "x")
lines(XY$Concn, predict(my_glm))

[Can't get the picture to upload at the moment]

Now that's flexible and goes right through all the points. But then I see you've got this nice theoretical model that must have come from your subject matter:

Values ~ (ymax* Concn / (ec50 + Concn)) + Ns*XY$Concn + ymin

and I'm kind of jealous. Is there anything else that might explain bias in this theoretical model (if the discrepancy we're seeing is indeed bias)? Even if you can't explain away the bias, it sure is simple and does a pretty nice job explaining the variation. In fact, the R-squared is above 99.9%:

cor(XY$Values, predict(my_glm))^2
$\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.