7

I am current using the package (lme4) for a mixed effect model with random effects.

My model takes the form:

mod <- lmer(response ~ b1(predict1, bp) + b2(predict2, bp) + (1 | site), data = data)

The functions b1 and b2 come from the helpful advice in:

Estimating the break point in a broken stick / piecewise linear model with random effects in R [code and output included]

b1 <- function(x, bp) ifelse(x < bp, bp - x, 0)
b2 <- function(x, bp) ifelse(x < bp, 0, x - bp)

I am wondering, how would I set up my model so that the second part (line) of the broken stick model has a constant gradient?

So for example, I would like something like:

If X < breakpoint, Y = m1X + C1

If X > breakpoint, Y = C2

Hope I am making sense.

Matt.

mjburns
  • 1,077
  • 3
  • 12
  • 16

2 Answers2

2

I may have found a solution to my problem, but could someone please confirm what I've done?

I slightly modified the above functions b1 and b2. I have called them b3 and b4:

b3 <- function(x, bp) ifelse(x < bp, bp - x, 0)   #Y = mx + c
b4 <- function(x, bp) ifelse(x < bp, 0, 1)        #Y = C

The mixed model is similar, but b3 and b4 are used:

mod2 <- lmer(SIGNAL ~ b3(dci, bp) + b4(dci,bp) + (1 | gauge), data = cand_bug_data)
mjburns
  • 1,077
  • 3
  • 12
  • 16
1

I see your post is pretty old, but I'm working on the same issue -- and I found a slightly different solution than yours.. figured I'd post for others out there.

b1 <- function(x, bp) ifelse(x < bp, x, bp)

#Wrapper for Mixed effects model with variable break point
foo <- function(bp)
 {
   mod <- lmer(y ~ b1(x, bp) + (1|Site), data = dat) 
   REMLcrit(mod)
 }

mod <- lmer(y ~ b1(x, bp) + (1 | Site),  data = dat)

And everything else as appears on the original post

Sarah
  • 81
  • 6