inst/doc/lm.br.R

In lm.br: Linear Model with Breakpoint

### R code from vignette source 'lm.br.Rnw'

###################################################
### code chunk number 1: lm.br.Rnw:161-189
###################################################
library(lm.br)
log_odds <- c( -1.194, -2.023, -2.285, -1.815, -1.673, -1.444, -1.237, -1.228 )
year <- c(1998.92, 2001.25, 2002.29, 2003.37, 2004.37, 2005.71, 2006.71, 2007.71)
VarCov <- matrix(  c(   0.0361, 0, 0, 0, 0, 0, 0, 0,
          0, 0.0218, 0.0129, 0, 0, 0, 0, 0,
          0, 0.0129, 0.0319, 0, 0, 0, 0, 0,
          0, 0, 0, 0.0451, 0.0389, 0, 0, 0,
          0, 0, 0, 0.0389, 0.0445, 0, 0, 0,
          0, 0, 0, 0, 0, 0.0672, 0.0607, 0.0607,
          0, 0, 0, 0, 0, 0.0607, 0.0664, 0.0607,
          0, 0, 0, 0, 0, 0.0607, 0.0607, 0.0662 ) ,  
            nrow = 8, ncol = 8 )
dd <- lm.br( log_odds ~ year, w = VarCov, inv = TRUE, var.known = TRUE )
bounds <- dd$cr( CL=0.90, out='v')
n <- length(dd$x1)
nbd <- nrow(bounds)
title <- "Exact 90% confidence region for changepoint"
x <- y <- matrix( NA, max(n,nbd), 4 )
x[1:n,1] <- dd$x1
y[1:n,1] <- dd$y
x[1:nbd,2:3] <- bounds[,1]
y[1:nbd,2:3] <- bounds[,2:3]
x[1:n,4] <- dd$x1
y[1:n,4] <- c( -1.4571, -1.60646, -1.65565, -1.67545, -1.53292,
                   -1.76094, -1.55349, -1.54487 )
matplot( x, y, 
   type=c('p','l','l','p'), pch=c(15,0), lty='solid', col='black',
   lwd=2, cex=1.5, main=title, xlab=dd$x1nm, ylab=dd$ynm )


###################################################
### code chunk number 2: lm.br.Rnw:201-216
###################################################
log_odds <-  c( -1.194, -2.023, -2.285, -1.815, -1.673, 
  -1.444, -1.237, -1.228 )
year <-  c( 1998.92, 2001.25, 2002.29, 2003.37, 2004.37, 
  2005.71, 2006.71, 2007.71 )
VarCov <-  matrix(  c(0.0361,  0,  0,  0,  0,  0,  0,  0,
  0,  0.0218,  0.0129,  0,  0,  0,  0,  0,
  0,  0.0129,  0.0319,  0,  0,  0,  0,  0,
  0,  0,  0,  0.0451,  0.0389,  0,  0,  0,
  0,  0,  0,  0.0389,  0.0445,  0,  0,  0,
  0,  0,  0,  0,  0,  0.0672,  0.0607,  0.0607,
  0,  0,  0,  0,  0,  0.0607,  0.0664,  0.0607,
  0,  0,  0,  0,  0,  0.0607,  0.0607,  0.0662), nrow=8, ncol=8)
dd <- lm.br( log_odds ~ year, w= VarCov, inv= T, var.known= T )
dd$ci( )
dd$ci( method = "AF" )


###################################################
### code chunk number 3: lm.br.Rnw:242-243
###################################################
lm.br( mpg ~ wt + hp, data = mtcars )