piece-wise linear regression nls function
Instead of reinventing the wheel, why not use the "segmented" package?
Having stored your data in a data frame "X" I did:
require(segmented)
m1 <- lm(FM ~ BMIJS,data=X)
m2 <- segmented(m1,seg.Z=~BMIJS,psi=list(BMIJS=35))
which worked instantaneously. The break point is estimated as 41.580, with
a standard error of 2.094 I then did:
with(X, plot(BMIJS,FM))
plot(m2,add=TRUE)
which seems to look as good as one can expect.
I must say however that the plot of your data does not look to me
as though a broken-stick model is appropriate. Why not just a straight
line?
cheers,
Rolf Turner
On 01/10/2013 02:33 PM, John Sorkin wrote:
windows 7, R 2.12 I am trying to run a piecewise linear regression with a single knot, i.e. a regression composed of two straight lines where the two lines intersect at an x value given by the variable knot. I wish to estimate the slope of both lines, the value of knot, the x value where the two lines intersect, and an intercept. I am using the nls code below, and get the following error message: Error in nls(FM ~ blow * BMIJS + bhi * sapply(BMIJS - knot, max, 0), start = list(knot = 25, : singular gradient nls code: test <- nls(FM~blow*BMIJS+bhi*sapply(BMIJS-knot,max,0),start=list(knot=25,blow=1,bhi=1),data=FeMaleData) summary(test) greatly shortened version of my data (the full data set has 450 records) FM BMIJS 2 55.878 40.57273 4 34.270 27.76939 5 20.123 21.73818 6 19.320 19.71203 9 49.701 43.55356 10 51.188 37.84742 11 46.753 37.71003 13 65.079 37.23438 14 37.097 36.81806 15 30.625 29.92783 17 50.617 42.42754 18 63.954 48.78709 20 29.790 26.97648 21 36.558 34.79373 22 41.275 33.03063 24 27.682 27.24508 26 37.968 35.41399 28 24.878 27.20250 30 47.513 35.77961 31 51.315 37.46032 33 41.944 36.40212 34 38.150 32.83818 35 60.719 42.48594 36 42.643 34.29355 38 40.728 32.42817 42 34.814 30.57573 43 32.896 29.32912 44 30.430 25.44183 46 48.986 37.90910 49 47.485 36.34642 52 46.312 38.64647 54 45.228 33.08783 55 45.391 35.86965 59 37.256 32.66507 60 27.367 28.49880 63 38.663 34.34131 64 34.527 29.57858 67 58.368 38.97266 68 13.473 17.35397 69 22.456 20.80958 71 28.829 25.50056 73 15.487 20.22202 76 18.313 21.38991 77 41.535 36.85707 78 56.124 40.51978 80 52.587 40.77256 81 24.991 25.48543 83 56.327 39.97214 84 70.836 36.52915 85 62.294 42.45244 86 39.689 35.18527 87 35.006 35.15136 88 47.378 37.54779 89 18.149 23.99236 90 33.041 28.10476 91 28.884 26.74443 92 37.670 32.25230 94 55.410 43.72364 99 34.461 35.05930 101 59.727 42.83035 102 41.913 35.64677 104 66.644 41.01642 105 55.250 43.86426 107 45.196 31.78370 108 36.476 33.45537 109 34.386 29.08402 110 39.277 36.98500 111 53.789 45.54654 112 33.077 29.09559 116 57.246 39.98031 120 52.546 40.12191 122 34.409 29.70977 123 31.188 28.75295 126 54.567 38.15226 129 19.193 22.71878 133 39.322 33.45712 134 41.415 31.28980 136 57.616 36.94016 140 28.162 24.40219 142 37.524 29.92673 143 29.611 29.15452 144 26.780 26.53462 146 47.219 35.14919 147 35.341 28.68955 148 44.827 37.68317 149 54.180 41.12226 150 41.636 30.00930 151 33.626 28.00164 156 34.334 29.64970 160 36.317 30.12031 161 46.823 35.64603 163 39.506 34.27740 164 61.619 39.20019 169 48.984 35.77558 171 66.467 41.59008 172 70.144 42.79996 173 37.324 31.56521 174 66.882 46.04938 182 54.239 38.21065 184 48.800 32.01630