2

EDIT

enter image description here

I have this data which I have plotted against each other. Y is measuring cumulative planted area (so it goes from 0% to 100%) while X is measuring % cumulative rainfall with respect to that location's climatology. So I am interested in know how much rainfall needs to be accumulated in order to plant 50% of the area 

   y <- c(0.5,3.0,22.2,46.0,77.3,97.0,98.9,100.0)
   x <- c(0.96,10.68,17.55,32.46,41.04,47.51,60.98,80.99)

   plot(x,y, pch = 19, xlab = "X",ylab = "Y")

I am looking for a method to say know what it the value of Y for a given value of X. I am specially interested in understanding what is the value of Y for X = 50. Could you please help me how to fit a curve and in R?

Nick Cox
  • 48,377
  • 8
  • 110
  • 156
89_Simple
  • 751
  • 1
  • 9
  • 23
  • 3
    If you were to explain what Y is measuring and how it was measured, that would create an opportunity to assess the uncertainty in the prediction you ask for. Otherwise you will be left either with an arbitrary model (as described in one answer so far) or an interpolator (as described in another answer). The former will have no evident meaning while the latter will not be able to give you any information about the prediction error. – whuber Apr 19 '18 at 19:32
  • Y is measuring cumulative planted area (so it goes from 0% to 100%) while X is measuring cumulative rainfall. So I am interested in know how much rainfall needs to be accumulated in order to plant 50% of the area – 89_Simple Apr 20 '18 at 08:36
  • 2
    Do you have other datasets? I am not likely to revise my answer (quite a lot of work and I am reluctant to throw more time at it). But knowing now what I just guessed before -- the curve should be monotonic and is bounded by 0 and 100% -- casts a more negative light on polynomial regression and on any interpolation method that doesn't automatically respect these properties. – Nick Cox Apr 20 '18 at 09:10
  • I have more datasets. They are like this only for different locations. For each location, I want to determine the percentage accumulated rainfall till the 50% of planting. Thank you for your answer. – 89_Simple Apr 20 '18 at 09:19
  • 1
    There is a pchip implementation in R. I used my own Stata implementation which was a fairly straight translation of the original MATLAB code. – Nick Cox Apr 20 '18 at 10:52
  • You can use a monotone spline fit, there is an example in: https://stats.stackexchange.com/questions/206073/looking-for-function-to-fit-sigmoid-like-curve/316446#316446 – kjetil b halvorsen Jul 02 '18 at 14:01

3 Answers3

5

How to do it in R is off-topic here, but if R is half as good as they say you should have no difficulty finding implementations.

I'd suggest not thinking of curve fitting, but of interpolation. Here are three methods of interpolation: linear (familiar since childhood), cubic spline and piecewise cubic Hermite. The aspect ratios for the graphs exaggerate the differences, which is needed for comparison. Cubic spline in particular waggles locally but pchip respects local maxima and minima and thus is monotonic for these data. (If you have subject-specific knowledge on whether the relationship should be monotonic, do tell.)

To do the interpolation I added a grid at $x =$ 1(1)80.

enter image description here

For $x =$ 50 I get 

  +-----------------------------------+
  |    linear      spline       pchip |
  |-----------------------------------|
  | 97.351225   100.56174   97.639561 |
  +-----------------------------------+

For comparison I give here polynomial fits up to 6

enter image description here

and a logistic equation

enter image description here

The latter predicts 93.12383 at $x =$ 50. The polynomials show the usual sad progression from underfitting to overfitting. What extra spin the logistic may deserve from relevance to underlying process I can't say. There might be a case for taking 100 as known upper limit and leaving only two parameters to estimate.

Nick Cox
  • 48,377
  • 8
  • 110
  • 156
  • 1
    +1. But in light of the additional information--that the curve should be monotonically increasing from (0,0) to (100,100)--one would be tempted to fit a GLM using something like a Beta CDF link (two parameters only). The logistic regression at least is the likeliest, out of the options you offer, to indicate the magnitude of the errors. The utility of this exercise is somewhat doubtful, unless the "prediction" is intended as a retrodiction: before the end of the season, how is one to know the total accumulated rainfall? – whuber Apr 20 '18 at 13:36
  • 1
    @whuber Indeed. I have a climatologist persona but left that on one side. My last sentence is key here. If the OP had revealed all in the question at the outset I would have treated 100 as a known constant. It's fairly reassuring that the model comes up with 103 or so as an estimate, but we can still do better. I learned from my namesake Sir David Cox years ago that good models always respect known limiting behaviour. – Nick Cox Apr 20 '18 at 13:52
  • I fit the data to a few hundred known, named equations with three or less parameters and logistic was near the top of the list. I note that the value should probably be clipped to a maximum of 100. Based on my analysis this equation is likely a good general solution to the remaining data sets. – James Phillips Apr 21 '18 at 21:23
1

You can use a polynomial regression to do this. It looks like a 3rd or 4th degree polynomial would fit your data. Try something like:

# This fits a 3rd degree polynomial to the data
fit <- lm(y ~ poly(x, 3, raw=TRUE))
predict(fit, data.frame(x=50))

Changing the 3 will change the degree of the polynomial. You can also put any vector into data.frame(x=c(...)) to get the predicted value of any x-values. You can plot the line by doing:

plot(x, y)
xx <- seq(0, 85, length.out=100)
lines(xx, predict(xx, data.frame(x=xx)))
JP Trawinski
  • 783
  • 7
  • 11
1

The fitting to a linear piecewise function might be sufficient to model your experiment. The figure below shows the result of fitting a three segments piecewise function to your data.

enter image description here

The equation of the function is :

enter image description here

The computed values of the parameters $a_1 , a_2 , p_1 , q_1 , p_2 , q_2 , p_3 , q_3 $ are shown on the figure.

The method of calculus is very simple (no iteration, no initial guessed value). See pages 30 and 31 in the paper : https://fr.scribd.com/document/380941024/Regression-par-morceaux-Piecewise-Regression-pdf

If you want the value of $x$ corresponding to a given value $y$ on the second segment : $$x=\frac{y-q_2}{p_2}$$ For example with $y=50$ we get $x=\frac{50-(-20.159)}{2.262}=31.02$

The logistic equation found by Nick Cox leads to $x=31.86$ which is close (difference = $0.84$).

On the other hand, at $x=50$ the piecewise function gives $y=2.262*50-20.159=92.94$

The N.Cox's logistic function gives $y=93.14$ very close again (difference = $0.20$)

On a statistical viewpoint, the logistic model might be more significant that the piecewise model. So, one could prefer the logistic model, even if the regression for fitting is more complicated than with the piecewise model.

JJacquelin
  • 551
  • 3
  • 8
  • 1
    I think on biological grounds one expects continuous growth and thus smooth change here. Thus this form is implausible in principle for all that in practice it gives similar predictions. With 8 parameters to be estimated it is not parsimonious! I also underline a point made elsewhere: the upper limit of 100% should be regarded as a known constant; that wasn't apparent in the question as first posted. The particular data example is not only just one example; it is too small to justify fine judgments about which method is best. – Nick Cox Jul 02 '18 at 17:36
  • @Nick It appears that the function is intended at least to be continuous, thereby making it a *linear spline* with "only" six rather than eight parameters. However, given that this is a cumulative function, fitting via least squares (as hinted at by the reference to "method of calculus") is inappropriate: the responses surely are strongly correlated and heteroscedastic. The latter considerations apply to any proposed solution, btw. – whuber Jul 02 '18 at 17:47
  • In retrospect I would delete "continuous growth and thus". – Nick Cox Jul 02 '18 at 17:56
  • @Nick Cox. Of course, that is what I pointed out in the last sentence of my answer. The aim of my answer is not to propose a physical model. It is a purely computational tool which can compute numerical results close to a physical model, but on a simpler way (without iterative process and without initial guess). That is a answer to Crop89's question of numerical calculus, without discussion about the physical model. The proposed piecewise function is not a physical model, it is a mathematical approximate of a physical model (which, by the way, is not clearly known). – JJacquelin Jul 02 '18 at 18:26
  • Glad we agree broadly, but you didn't give any reason _why_ a logistic might be preferred. Very minor point: "physical model" is indelicate wording when it's biology! A more serious point is that you can keep the same attitude and make the procedure simpler yet by resorting to linear interpolation. – Nick Cox Jul 02 '18 at 18:34