How to find a split point in a curve?

classic Classic list List threaded Threaded
5 messages Options
Reply | Threaded
Open this post in threaded view
|

How to find a split point in a curve?

Luigi
Dear all,
I am trying to find a turning point in some data. In the initial phase, the
data increases more or less exponentially (thus it is linear in a nat log
transform), then reaches a plateau. I would like to find the point that
marks the end of the exponential phase.
I understand that the function spline can build a curve; is it possible
with it to find the turning points? I have no idea of how to use spline
though.
Here is a working example.
Thank you

```
Y = c(259, 716, 1404, 2173, 3944, 5403, 7140, 9121,
      11220, 13809, 16634, 19869, 23753, 27447,
      30590, 33975, 36627, 39600, 42067, 44082,
      58190, 63280, 65921, 67929, 69977, 71865,
      73614, 74005, 74894, 75717, 76365, 76579,
      77087, 77493, 77926, 78253, 78680, 79253,
      79455, 79580, 79699, 79838, 79981, 80080,
      80124, 80164, 80183, 80207, 80222, 80230,
      80241, 80261, 80261, 80277, 80290, 80303,
      80337, 80376, 80422, 80461, 80539, 80586,
      80653, 80708, 80762, 80807, 80807, 80886,
      80922, 80957, 80988, 81007, 81037, 81076,
      81108, 81108, 81171, 81213, 81259, 81358,
      81466, 81555, 81601, 81647, 81673, 81998,
      82025, 82041, 82053, 82064, 82094, 82104,
      82110, 82122, 82133, 82136, 82142, 82164,
      82168, 82180, 82181, 82184, 82187, 82188,
      82190, 82192, 82193, 82194)
Y = log(Y)
X = 1:length(Y)
plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="zoomed in"))
abline(lm(Y[1:3] ~ X[1:3]))
abline(lm(Y[1:5] ~ X[1:5]), lty=2)
text(7, 6, "After third or fifth point, there is deviance", pos=3)
text(2.5, 10, "Solid line: linear model points 1:3", pos =3)
text(2.5, 9, "Dashed line: linear model points 1:5", pos =3)
plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="overall"))
abline(lm(Y[1:3] ~ X[1:3]))
```

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: How to find a split point in a curve?

Rui Barradas
Hello,

Are you looking for a segmented regression?

fit <- lm(Y ~ X)
seg <- segmented::segmented(fit, seg.Z = ~X)
seg$psi[, 'Est.']
#[1] 29.21595

plot(X, Y)
plot(seg, add = TRUE)


Hope this helps,

Rui Barradas


Às 16:12 de 14/05/20, Luigi Marongiu escreveu:

> Dear all,
> I am trying to find a turning point in some data. In the initial phase, the
> data increases more or less exponentially (thus it is linear in a nat log
> transform), then reaches a plateau. I would like to find the point that
> marks the end of the exponential phase.
> I understand that the function spline can build a curve; is it possible
> with it to find the turning points? I have no idea of how to use spline
> though.
> Here is a working example.
> Thank you
>
> ```
> Y = c(259, 716, 1404, 2173, 3944, 5403, 7140, 9121,
>        11220, 13809, 16634, 19869, 23753, 27447,
>        30590, 33975, 36627, 39600, 42067, 44082,
>        58190, 63280, 65921, 67929, 69977, 71865,
>        73614, 74005, 74894, 75717, 76365, 76579,
>        77087, 77493, 77926, 78253, 78680, 79253,
>        79455, 79580, 79699, 79838, 79981, 80080,
>        80124, 80164, 80183, 80207, 80222, 80230,
>        80241, 80261, 80261, 80277, 80290, 80303,
>        80337, 80376, 80422, 80461, 80539, 80586,
>        80653, 80708, 80762, 80807, 80807, 80886,
>        80922, 80957, 80988, 81007, 81037, 81076,
>        81108, 81108, 81171, 81213, 81259, 81358,
>        81466, 81555, 81601, 81647, 81673, 81998,
>        82025, 82041, 82053, 82064, 82094, 82104,
>        82110, 82122, 82133, 82136, 82142, 82164,
>        82168, 82180, 82181, 82184, 82187, 82188,
>        82190, 82192, 82193, 82194)
> Y = log(Y)
> X = 1:length(Y)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="zoomed in"))
> abline(lm(Y[1:3] ~ X[1:3]))
> abline(lm(Y[1:5] ~ X[1:5]), lty=2)
> text(7, 6, "After third or fifth point, there is deviance", pos=3)
> text(2.5, 10, "Solid line: linear model points 1:3", pos =3)
> text(2.5, 9, "Dashed line: linear model points 1:5", pos =3)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="overall"))
> abline(lm(Y[1:3] ~ X[1:3]))
> ```
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: How to find a split point in a curve?

Bert Gunter-2
In reply to this post by Luigi
You need to mathematically define 'turning point' first: "end of
exponential phase" is subjective and meaningless. Once you have a precise
mathematical formulation in hand, you can proceed.

Bert Gunter

"The trouble with having an open mind is that people keep coming along and
sticking things into it."
-- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )


On Thu, May 14, 2020 at 8:12 AM Luigi Marongiu <[hidden email]>
wrote:

> Dear all,
> I am trying to find a turning point in some data. In the initial phase, the
> data increases more or less exponentially (thus it is linear in a nat log
> transform), then reaches a plateau. I would like to find the point that
> marks the end of the exponential phase.
> I understand that the function spline can build a curve; is it possible
> with it to find the turning points? I have no idea of how to use spline
> though.
> Here is a working example.
> Thank you
>
> ```
> Y = c(259, 716, 1404, 2173, 3944, 5403, 7140, 9121,
>       11220, 13809, 16634, 19869, 23753, 27447,
>       30590, 33975, 36627, 39600, 42067, 44082,
>       58190, 63280, 65921, 67929, 69977, 71865,
>       73614, 74005, 74894, 75717, 76365, 76579,
>       77087, 77493, 77926, 78253, 78680, 79253,
>       79455, 79580, 79699, 79838, 79981, 80080,
>       80124, 80164, 80183, 80207, 80222, 80230,
>       80241, 80261, 80261, 80277, 80290, 80303,
>       80337, 80376, 80422, 80461, 80539, 80586,
>       80653, 80708, 80762, 80807, 80807, 80886,
>       80922, 80957, 80988, 81007, 81037, 81076,
>       81108, 81108, 81171, 81213, 81259, 81358,
>       81466, 81555, 81601, 81647, 81673, 81998,
>       82025, 82041, 82053, 82064, 82094, 82104,
>       82110, 82122, 82133, 82136, 82142, 82164,
>       82168, 82180, 82181, 82184, 82187, 82188,
>       82190, 82192, 82193, 82194)
> Y = log(Y)
> X = 1:length(Y)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="zoomed in"))
> abline(lm(Y[1:3] ~ X[1:3]))
> abline(lm(Y[1:5] ~ X[1:5]), lty=2)
> text(7, 6, "After third or fifth point, there is deviance", pos=3)
> text(2.5, 10, "Solid line: linear model points 1:3", pos =3)
> text(2.5, 9, "Dashed line: linear model points 1:5", pos =3)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="overall"))
> abline(lm(Y[1:3] ~ X[1:3]))
> ```
>
>         [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: How to find a split point in a curve?

Gabor Grothendieck
In reply to this post by Luigi
We can use nls2 to try each value in 10:100 as a possible split point
picking the one with lowest residual sum of squares:

library(nls2)
fm <- nls2(Y ~ cbind(1, pmin(X, X0)), start = data.frame(X0 = 10:100),
  algorithm = "plinear-brute")
plot(Y ~ X)
lines(fitted(fm) ~ X, col = "red")

> fm
Nonlinear regression model
  model: Y ~ cbind(1, pmin(X, X0))
   data: parent.frame()
     X0   .lin1   .lin2
18.0000  6.5570  0.2616
 residual sum-of-squares: 4.999

Number of iterations to convergence: 91
Achieved convergence tolerance: NA

On Thu, May 14, 2020 at 11:13 AM Luigi Marongiu
<[hidden email]> wrote:

>
> Dear all,
> I am trying to find a turning point in some data. In the initial phase, the
> data increases more or less exponentially (thus it is linear in a nat log
> transform), then reaches a plateau. I would like to find the point that
> marks the end of the exponential phase.
> I understand that the function spline can build a curve; is it possible
> with it to find the turning points? I have no idea of how to use spline
> though.
> Here is a working example.
> Thank you
>
> ```
> Y = c(259, 716, 1404, 2173, 3944, 5403, 7140, 9121,
>       11220, 13809, 16634, 19869, 23753, 27447,
>       30590, 33975, 36627, 39600, 42067, 44082,
>       58190, 63280, 65921, 67929, 69977, 71865,
>       73614, 74005, 74894, 75717, 76365, 76579,
>       77087, 77493, 77926, 78253, 78680, 79253,
>       79455, 79580, 79699, 79838, 79981, 80080,
>       80124, 80164, 80183, 80207, 80222, 80230,
>       80241, 80261, 80261, 80277, 80290, 80303,
>       80337, 80376, 80422, 80461, 80539, 80586,
>       80653, 80708, 80762, 80807, 80807, 80886,
>       80922, 80957, 80988, 81007, 81037, 81076,
>       81108, 81108, 81171, 81213, 81259, 81358,
>       81466, 81555, 81601, 81647, 81673, 81998,
>       82025, 82041, 82053, 82064, 82094, 82104,
>       82110, 82122, 82133, 82136, 82142, 82164,
>       82168, 82180, 82181, 82184, 82187, 82188,
>       82190, 82192, 82193, 82194)
> Y = log(Y)
> X = 1:length(Y)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="zoomed in"))
> abline(lm(Y[1:3] ~ X[1:3]))
> abline(lm(Y[1:5] ~ X[1:5]), lty=2)
> text(7, 6, "After third or fifth point, there is deviance", pos=3)
> text(2.5, 10, "Solid line: linear model points 1:3", pos =3)
> text(2.5, 9, "Dashed line: linear model points 1:5", pos =3)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="overall"))
> abline(lm(Y[1:3] ~ X[1:3]))
> ```
>
>         [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.



--
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.com

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: How to find a split point in a curve?

PIKAL Petr
In reply to this post by Luigi
Hi

Maybe I am wrong but it seems to me that it is cumulative data from recent
epidemy in some state.

If yes, instead of inventing wheel I would go to canned and proved solution
using tools from
https://www.repidemicsconsortium.org/

I found useful and enlightening this blog especially first part from
February 18th.
https://timchurches.github.io/blog/

Cheers
Petr

> -----Original Message-----
> From: R-help <[hidden email]> On Behalf Of Luigi Marongiu
> Sent: Thursday, May 14, 2020 5:12 PM
> To: r-help <[hidden email]>
> Subject: [R] How to find a split point in a curve?
>
> Dear all,
> I am trying to find a turning point in some data. In the initial phase,
the

> data increases more or less exponentially (thus it is linear in a nat log
> transform), then reaches a plateau. I would like to find the point that
> marks the end of the exponential phase.
> I understand that the function spline can build a curve; is it possible
> with it to find the turning points? I have no idea of how to use spline
> though.
> Here is a working example.
> Thank you
>
> ```
> Y = c(259, 716, 1404, 2173, 3944, 5403, 7140, 9121,
>       11220, 13809, 16634, 19869, 23753, 27447,
>       30590, 33975, 36627, 39600, 42067, 44082,
>       58190, 63280, 65921, 67929, 69977, 71865,
>       73614, 74005, 74894, 75717, 76365, 76579,
>       77087, 77493, 77926, 78253, 78680, 79253,
>       79455, 79580, 79699, 79838, 79981, 80080,
>       80124, 80164, 80183, 80207, 80222, 80230,
>       80241, 80261, 80261, 80277, 80290, 80303,
>       80337, 80376, 80422, 80461, 80539, 80586,
>       80653, 80708, 80762, 80807, 80807, 80886,
>       80922, 80957, 80988, 81007, 81037, 81076,
>       81108, 81108, 81171, 81213, 81259, 81358,
>       81466, 81555, 81601, 81647, 81673, 81998,
>       82025, 82041, 82053, 82064, 82094, 82104,
>       82110, 82122, 82133, 82136, 82142, 82164,
>       82168, 82180, 82181, 82184, 82187, 82188,
>       82190, 82192, 82193, 82194)
> Y = log(Y)
> X = 1:length(Y)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="zoomed in"))
> abline(lm(Y[1:3] ~ X[1:3]))
> abline(lm(Y[1:5] ~ X[1:5]), lty=2)
> text(7, 6, "After third or fifth point, there is deviance", pos=3)
> text(2.5, 10, "Solid line: linear model points 1:3", pos =3)
> text(2.5, 9, "Dashed line: linear model points 1:5", pos =3)
> plot(Y ~ X, ylab = "Ln(Y)", xlim=c(0,10, main="overall"))
> abline(lm(Y[1:3] ~ X[1:3]))
> ```
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.