[R] Piecewise Regression with a known slope

[vmuggeo at dssm.unipa.it]

Dear Jarrett,

If I understand correctly your post, your constraints may be achieved
straightforwardly in segmented. See the code below.

> The error is also most likely gamma distributed..[SNIP]
The 'error' component can be specified in the 'initial' model by means of
the family argument in the glm() function

> I have attempted to use the segmented package, specifying something
> close to the visually estimated breakpoint as the value of psi, but
> it continues to return fits with a breakpoint that occurs somewhere
> in the middle of the linear portion of the line, ..[SNIP]

Mhmm..you can specify different starting values and to assess the
differences in the 'final' estimates. However if you says "..it continues
to return fits with a breakpoint that occurs somewhere in the middle.."
probably the ML estimate of the breakpoint is really in the middle.

Hope this helps you,

best,
vito



####simulate data
n<-50
x<-1:n/n
y<- 0-pmin(x-.5,0)+rnorm(50)*.03
plot(x,y) #This should be your scatterplot..
abline(0,0,lty=2)

o<-lm(y~x) #or glm(y~x,family=..)
o1<-segmented(o,seg.Z=~x,psi=list(x=.3))
slope(o1) #get the slope
points(x,fitted(o1),col=2)

#a parsimonious modelling: constrain right slope=0
o<-lm(y~1)
xx<- -x
o2<-segmented(o,seg.Z=~xx,psi=list(xx=-.3))
slope(o2)
points(x,fitted(o2),col=2)

#now constrain \hat{\mu}(x)=0 for x>psi

o<-lm(y~0)
o3<-segmented(o,seg.Z=~xx,psi=list(xx=-.3))
slope(o3)
points(x,fitted(o3),col=3)

> Hey, all.  I'm working on a data set with a broken stick linear
> regression where I know one of the two slopes.  It is a negative
> linear function until the line intersects with the x-axis, at which
> point it becomes 0.  It is not a nonlinear asymptotic function, and,
> indeed, using negative exponential or logistic types of fits as an
> approximation has tended to lead to an under or overestimation of
> values.  I am also very interested to know just what the breakpoint
> in the data is.
>
> Essentially
>
> if x<psi y = a + bx + error, where b is negative
> else y=0+error
>
> and I want to know the value of psi, as well as a and b.  The error
> is also most likely gamma distributed, as values<0 are not possible
> (nutrient concentrations).
>
> I have attempted to use the segmented package, specifying something
> close to the visually estimated breakpoint as the value of psi, but
> it continues to return fits with a breakpoint that occurs somewhere
> in the middle of the linear portion of the line, and the slope and
> intercept of the second half of the regression is not 0.
>
> Is there either a package that exists that will allow me to estimate
> such a model, or a function that I can use for optim or nlm (I admit,
> I am a rank novice at coding such functions).
>
> Thanks so much!
>
> -Jarrett
>
> ______________________________________________
> R-help at 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
> and provide commented, minimal, self-contained, reproducible code.
>
>