Re: [R] best way to fit a model

From: Achim Zeileis <Achim.Zeileis_at_wu-wien.ac.at>
Date: Sat 10 Sep 2005 - 04:25:48 EST

On Fri, 9 Sep 2005 14:55:04 -0300 Ronaldo Reis-Jr. wrote:

> Hi,
>
> I have some data that have this behaviour:
>
> |
> |*******
> | *
> | *
> | *
> | *
> |----------------
>
> What is the best and simpler way to fit this in R?

If the changepoint is known, then this is straightforward using lm:

## generate example data
set.seed(20050909)
x <- seq(0, 10, by = 0.25)
y.mean <- ifelse(x >= 5, x, 5)
y <- y.mean + rnorm(41)
plot(y ~ x)
lines(y.mean ~ x)

## fit linear model with break
fm <- lm(y ~ I((x-5) * (x >= 5)))
fm
y.fit1 <- fitted(fm)
lines(y.fit1 ~ x, col = 2)

If it is unknown, it can be estimated using Vito Muggeo's segmented package:

## estimate change point in x
library("segmented")
sfm <- segmented(lm(y ~ x), x, 6)
sfm
y.fit2 <- fitted(sfm)
lines(y.fit2 ~ x, col = 3)

This fits a continuous mean function. Alternatively, breakpoints() in strucchange can be used to estimate a break point:

## estimate break point in x
library("strucchange")
bp <- breakpoints(y ~ x)
summary(bp)
y.fit3 <- fitted(bp)
lines(y.fit3 ~ x, col = 4)

This does not enforce that the line is continuous, hence the jump in the fitted mean. Of course, the estimated breakpoint could be used to fit a continuous line model, but this is not what is optimized in breakpoints().
Z

> Thanks
> Ronaldo
> --
> Ela pilotava um Continenal 2001 com ignição automática Magiclic...
> --
> |> // | \\ [***********************************]
> | ( õ õ ) [Ronaldo Reis Júnior ]
> |> V [UFV/DBA-Entomologia ]
> | / \ [36570-000 Viçosa - MG ]
> |> /(.''`.)\ [Fone: 31-3899-4007 ]
> | /(: :' :)\ [chrysopa@insecta.ufv.br ]
> |>/ (`. `'` ) \[ICQ#: 5692561 | LinuxUser#: 205366 ]
> | ( `- ) [***********************************]
> |>> _/ \_Powered by GNU/Debian Woody/Sarge
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide!
> http://www.R-project.org/posting-guide.html
>



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sat Sep 10 04:28:39 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:40:09 EST