The NLMOD Procedure

Example 14.1 Segmented Model

Suppose you are interested in fitting a model that consists of two segments that connect in a smooth fashion. For example, the following model states that the mean of Y is a quadratic function in x for values of x less than and that the mean of Y is constant for values of x greater than :

In this model equation, , , and are the coefficients of the quadratic segment, and c is the plateau of the mean function. The NLMOD procedure can fit such a segmented model even when the join point, , is unknown.

Suppose you also want to impose conditions on the two segments of the model. First, the curve should be continuous—that is, the quadratic and the plateau section need to meet at . Second, the curve should be smooth—that is, the first derivative of the two segments with respect to x needs to coincide at .

The continuity condition requires that

The smoothness condition requires that

If you solve for and substitute your result in the expression for c, the two conditions jointly imply that

Although there are five unknowns, the model contains only three independent parameters. Together the continuity and smoothness restrictions completely determine two parameters, given the other three.

The following DATA step creates the data table for this example in your CAS session:

data mycas.A;
   input y x @@;
   datalines;
.46 1  .47  2 .57  3 .61  4 .62  5 .68  6 .69  7
.78 8  .70  9 .74 10 .77 11 .78 12 .74 13 .80 13
.80 15 .78 16
 ;

The following PROC NLMOD statements fit this segmented model:

proc nlmod data=mycas.A out=mycas.B;
   parms alpha=.45 beta=.05 gamma=-.0025;

   x0 = -.5*beta / gamma;

   if (x < x0) then
        yp = alpha + beta*x  + gamma*x*x;
   else
        yp = alpha + beta*x0 + gamma*x0*x0;

   model y ~ residual(yp);

   estimate 'join point' -beta/2/gamma;
   estimate 'plateau value c' alpha - beta**2/(4*gamma);
   predict 'predicted' yp pred=yp;
   predict 'response' y pred=y;
   predict 'x' x pred=x;
run;

The parameters of the model are , , and , which are represented in the PROC NLMOD statements by the variables alpha, beta, and gamma, respectively. To model the two segments, a conditional statement assigns the appropriate expression to the mean function, depending on the value of . The ESTIMATE statements compute the values of and c. The PREDICT statement computes predicted values for plotting and saves them to data table b.

The results from fitting this model are shown in Output 14.1.1 through Output 14.1.3. The iterative optimization converges after six iterations (Output 14.1.1). Output 14.1.2 shows the estimated parameters. Output 14.1.3 indicates that the join point is 12.7477 and the plateau value is 0.7775.

Output 14.1.1: Nonlinear Least Squares Iterative Phase

Quadratic Model with Plateau

The NLMOD Procedure

Iteration History
IterationEvaluationsObjective
Function
ChangeMaximum
Gradient
050.0035144531 7.184063
120.00073527160.002779182.145337
220.00062927510.000106000.032551
320.00062912610.000000150.002952
420.00062912440.000000000.000238
520.00062912440.000000000.000023
620.00062912440.000000002.313E-6

Convergence criterion (GCONV=1E-8) satisfied.


Output 14.1.2: Least Squares Analysis of the Quadratic Model

Analysis of Variance
SourceDFSum of
Squares
Mean
Square
F ValueApprox
Pr > F
Model20.17690.0884114.22<.0001
Error130.01010.000774  
Corrected Total150.1869   

Parameter Estimates
ParameterEstimateStandard
Error
DFt ValueApprox
Pr > |t|
Approximate 95% Confidence
Limits
LowerUpper
alpha0.39210.02671314.70<.00010.33450.4497
beta0.06050.00842137.18<.00010.04230.0787
gamma-0.002370.00055113-4.300.0009-0.00356-0.00118


Output 14.1.3: Additional Estimates for the Quadratic Model

Additional Estimates
LabelEstimateStandard
Error
DFt ValueApprox
Pr > |t|
AlphaApproximate Confidence
Limits
LowerUpper
join point12.74771.2781139.97<.00010.059.986415.5089
plateau value c0.77750.01231363.11<.00010.050.75090.8041


The following statements produce a graph, shown in Output 14.1.4, of the observed and predicted values along with reference lines for the join point and plateau estimates:

data B;
   set mycas.B;
run;

proc sort data = B;
   by x;
run;
proc sgplot data=B noautolegend;
   yaxis label='Observed or Predicted';
   refline 0.7775  / axis=y label="Plateau"    labelpos=min;
   refline 12.7477 / axis=x label="Join point" labelpos=min;
   scatter y=y  x=x;
   series  y=yp x=x;
run;

Output 14.1.4: Observed and Predicted Values for the Quadratic Model

 Observed and Predicted Values for the Quadratic Model


Last updated: December 21, 2018