3
$\begingroup$

I have a dataset that I'm currently smoothing using a 6th order least-squares polynomial which models the initial sharp curve well but as the data tapers out some noise at the end causes a downwards curve at the end. It's probably easiest to show a plot of my current data followed by the results of fitting extrapolated by a few extra points:

Polynomial fitting

I'm expecting the data to taper out as it moves along so my initial thoughts were to maybe to break the data into thirds and apply lower order polynomials as it goes along. But I don't have a good background in statistics and I'm guessing this problem has some elegant solution I'm not aware of. Because the dataset is fairly short I thought I'd include a copy in case anyone was inspired to demonstrate some techniques:

0   0.632696209
1   0.954847196
2   1.211956632
3   1.282479631
4   1.348572899
5   1.36480867
6   1.411463353
7   1.413356854
8   1.443252208
9   1.5086791
10  1.519056572
11  1.515427082
12  1.545547955
13  1.567970692
14  1.554558155
15  1.539432376
16  1.662380739
17  1.630363821
18  1.627295705
19  1.574494008
20  1.651035552
21  1.59834609
22  1.612328401
23  1.634990079
24  1.668922299
25  1.622188079
26  1.759410212
27  1.849993968
28  1.78653232
29  1.7872419
30  1.639521297
31  1.70151144
32  1.868248296
33  1.707963085
34  1.703480801
35  1.811524105
36  1.790136202
37  1.79970863
38  1.944272907
39  1.818103305
40  1.750744767
43  1.77700684
44  1.80647612
$\endgroup$
6
  • $\begingroup$ Is this time series data with equal intervals between measurements ? $\endgroup$
    – IrishStat
    Jul 3 '19 at 12:50
  • $\begingroup$ @IrishStat It's not a time series as such, the X is how many times an event has occurred over a period that is the same for all samples $\endgroup$
    – PeterJ
    Jul 3 '19 at 13:02
  • $\begingroup$ can you please detail precisely what .632696209 is at X=0 $\endgroup$
    – IrishStat
    Jul 3 '19 at 13:05
  • $\begingroup$ @IrishStat it's a probability, I don't want to be vague but someone asked me not to disclose full details but it means based on the X value the event has a 0.63 chance of occurring. $\endgroup$
    – PeterJ
    Jul 3 '19 at 13:10
  • $\begingroup$ Do you have an idea of when the break in slope occurs, e.g., at what horizontal coordinate? $\endgroup$ Jul 3 '19 at 13:12
5
$\begingroup$

My equation search on your data turned up a simple two-parameter power equation, "y = a * pow(x, b)", with parameters a = 1.0769014059561925E+00 and b = 1.4153886539395866E-01 yielding R-squared = 0.789 and RMSE = 0.112

plot

$\endgroup$
1
  • $\begingroup$ Note that parameter "a" is almost equal to 1.0, and I tried the one-parameter equation "y = pow(x, b)" and this two-parameter version is slightly better. If you have a large number of data sets to be fitted, it is possible that the one-parameter version is all that you need. $\endgroup$ Jul 3 '19 at 17:52
3
$\begingroup$

You have a specific kind of structural change in your data generating model. Specifically, your straight regression line seems to change in slope at a particular point. In order that the two line segments meet up, the intercept needs to change in a specific way at that same point.

You should look at "piecewise linear regression". Our tags and may also be helpful. For instance, here is a tutorial in R.

If you already know when (at what predictor value) the change occurs, then you can set up a specific design matrix with constraints to model this. However, if you also need to detect the change point, you need specific tools. For instance, in R, the segmented package will fit a piecewise linear regression where you specify the number of breaks (and optionally guesses for where they occur). Here is how segmented deals with your data, where we specify either one or two change points:

segmented

If you need to decide whether you have one, two or more break points, I would recommend something like cross validation.

xx <- c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 
18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 
34, 35, 36, 37, 38, 39, 40, 43, 44)
yy <- c(0.632696209, 0.954847196, 1.211956632, 1.282479631, 1.348572899, 
1.36480867, 1.411463353, 1.413356854, 1.443252208, 1.5086791, 
1.519056572, 1.515427082, 1.545547955, 1.567970692, 1.554558155, 
1.539432376, 1.662380739, 1.630363821, 1.627295705, 1.574494008, 
1.651035552, 1.59834609, 1.612328401, 1.634990079, 1.668922299, 
1.622188079, 1.759410212, 1.849993968, 1.78653232, 1.7872419, 
1.639521297, 1.70151144, 1.868248296, 1.707963085, 1.703480801, 
1.811524105, 1.790136202, 1.79970863, 1.944272907, 1.818103305, 
1.750744767, 1.77700684, 1.80647612)

library(segmented)
model_1 <- segmented(lm(yy~xx),npsi=1)
model_2 <- segmented(lm(yy~xx),npsi=2)

xx_fit <- seq(min(xx),max(xx),by=.01)

plot(xx,yy)
lines(xx_fit,predict(model_1,newdata=data.frame(xx=xx_fit)))
lines(xx_fit,predict(model_2,newdata=data.frame(xx=xx_fit)),col="red")
$\endgroup$
3
  • $\begingroup$ Thanks you've brought up a few useful techniques / keywords I wasn't aware of which I'll check out further but it seems like a good way to go. $\endgroup$
    – PeterJ
    Jul 3 '19 at 13:39
  • 4
    $\begingroup$ This is a good application for a regression spline, e.g., piecewise cubic polynomial. The join points for the various polynomials are called knot. Put a couple of knots near where you think rapid changes may be, from prior knowledge. $\endgroup$ Jul 3 '19 at 13:48
  • $\begingroup$ A simple power equation fits this data well, see my answer to this question. $\endgroup$ Jul 3 '19 at 15:36

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.