[S] Solving a segmented regression model in SPLUS

Bill Kreamer (kreamer@waisman.wisc.edu)
Thu, 19 Feb 1998 14:43:04 -0600


I am a novice Splus user and was wondering if someone had created a Splus
program to solve a segmented model. My data follows a quadratic decline of
a variable (k4neg) until reaching a constant value (plateau) over time
(hrbirth). I want to determine if the time points at which this model
reaches the constant (X0 "the hinge point") differs for different subjects
(subjid) and between treatment types.

I would eventually like to determine the "hinge-point" for each subject
and analyze this value for subjects grouped into different treatment types
(3). My data set containes 4 variables with the following structure:

1. Feedtype Subject's diet
2. Subjid Subject ID
3. hrbirth age of sample relative to subject's birth
4. k4neg negative of the variable value (needed for the SAS program)

I have used NLIN procedure in SAS to model the data and determine the
hinge points. I was hoping to do the same with SPLUS. (I included an
example of the SAS code for reference. The variable names used in the SAS
Program are enclosed with parenthesis above).

DATA TEMP1;
SET LISA.TOTAL2;

PROC SORT DATA=TEMP1;
BY SUBJID HRBIRTH;

*---------------------------------------------------------------------------
--*
| IF WE ASSUME THAT THE ELIMINATION OF ZCP FOLLOWS A SEGMENTED MODEL WERE
THE|
| INITIAL PHASE IS A PARABOLIC DECREASE FOLLOWED BY A CONSTANT LEVEL THEN
WE |
| CAN HYPOTHESIZE THAT:
|
|
|
| y = a + bx + cx**2 if x < x0
|
|
|
| y = p if x > x0
|
|
|
| THAT IS, FOR VALUES OF x LESS THAN x0 (THE HINGE-POINT), THE EQUATION
|
| RELATING y AND x IS QUADRATIC (A PARABOLA) AND, FOR VALUES OF x GREATER
|
| THAN x0, THE EQUATION IS CONSTANT (THE PLATEAU--A HORIZONTAL LINE).
PROC |
| NLIN (PAGES 575-606 IN SAS USERS GUIDE: STATISTICS VERSION 5 EDITION
AND |
| PAGES 1162-1165 IN SAS-STAT USERS GUIDE: VOLUME 2, GLM-VARCOMP VERSION 6
|
| 4TH EDITION) CAN FIT SUCH A SEGMENTED MODEL EVEN WHEN THE HINGE-POINT,
x0, |
| IS UNKNOWN. THE CURVE MUST BE CONTINOUS (THE TWO SECTIONS MUST MEET AT
|
| x0), AND THE CURVE MUST BE SMOOTH (THE FIRST DERIVATIVES WITH RESPECT TO
x |
| ARE THE SAME AT x0). THESE CONDITIONS IMPLY THAT:
|
|
|
| x0 = -b/2c
|
|
|
| p = a - (b**2/4c)
|
|
|
| THE SEGMENTED EQUATION INCLUDES ONLY THREE PARAMETERS, HOWEVER, THE
|
| EQUATION IS NON-LINEAR WITH RESPECT TO THESE PARAMETERS. THE FOLLOWING
|
| PROGRAM STATEMENTS USING PROC NLIN, CONDITIONALLY EXECUTES DIFFERENT
|
| SECTIONS OF CODE FOR THE TWO PARTS OF THE MODEL, DEPENDING ON WHETHER x
is |
| less than x0.
|
|
|
| A PUT STATEMENT IS USED TO PRINT THE CONSTRAINED PARAMETERS EVERY TIME
|
| THE PROGRAM IS EXECUTED.
|
*---------------------------------------------------------------------------
--*;

proc nlin DATA=TEMP1;
by subjid;
parms a=-8000 b=80 c=-0.1; ***ESTIMATES FOR INITIAL PARAMETERS***;
FILE PRINT;
x0=-.5*b/c; ***ESTIMATE OF HINGE-POINT***;
db=-.5/c; ***DERIVATIVE OF x0 WITH RESPECT TO b***;
dc=.5*b/c**2; ***DERIVATIVE OF x0 WITH RESPECT TO c***;
if hrbirth<x0 then do; ***QUADRATIC PART OF THE MODEL***;

model K4NEG = a + b*hrbirth + c*hrbirth*hrbirth;

der.a = 1;
der.b = hrbirth;
der.c = hrbirth*hrbirth;
end;
else do; ***PLATEAU PART OF THE MODEL***;

model K4NEG = a + b*x0 + c*x0*x0;

der.a = 1;
der.b = x0 + b*db + 2*c*x0*db;
der.c = b*dc + x0*x0 + 2*c*x0*dc;
end;
if _obs_=1 & _model_=1 then do; ***PRINT OUT PARAMETERS***;
plateau = a + b*x0 + c*x0*x0;
put x0= plateau=;
end;

output out=K4NEG predicted=yp;
id x0 plateau;
RUN;

*---------------------------------------------------------------------------
--*
| PLOT THE ESTIMATED VALUES AS PREDICTED BY THE CALCULATED MODEL (*), THE
|
| ACTUAL DATA (A-Z), AND THE PREDICTED HING-POINT (+) FOR EACH SUBJECT.
|
*---------------------------------------------------------------------------
--*;

proc plot DATA = K4NEG;

plot K4NEG*hrbirth yp*hrbirth='*' PLATEAU*x0='+' / overlay;

by subjid;
run;

proc sort DATA=K4NEG;
by feedtype;

*---------------------------------------------------------------------------
--*
| PLOT THE ESTIMATED VALUES AS PREDICTED BY THE CALCULATED MODEL (*), THE
|
| ACTUAL DATA (A-Z), AND THE PREDICTED HING-POINT (+) FOR EACH FEEDTYPE.
|
*---------------------------------------------------------------------------
--*;

proc plot DATA=K4NEG;

plot K4NEG*hrbirth yp*hrbirth='*' PLATEAU*x0='+' / overlay;

by feedtype;
run;

-----------------------------------------------------------------------
This message was distributed by s-news@wubios.wustl.edu. To unsubscribe
send e-mail to s-news-request@wubios.wustl.edu with the BODY of the
message: unsubscribe s-news