*******************************************************************;
   * The first line specifies the input data set and a convergence    *;
   * criterion that is stricter than default.                        *;
   *******************************************************************;
PROC NLMIXED data=cloc.data1 gconv=1e-12;
   *******************************************************************;
   * The next statement names the parameters and gives them starting *;
   * values.                                                         *;
   *******************************************************************;
PARMS gammaB=11, gammaL=16, gammaF=9, S=.8, s2e=50, 
      s2L=116, s2F=70, sLF sBL=30, sBF=30, s2B=60;
   *******************************************************************;
   * The next three lines specify the level-2 model.                 *;
   *******************************************************************;
 B=gammaB+zetaB;
 L=gammaL+zetaL;
 F=gammaF+zetaF;
   *******************************************************************;
   * The next line specifies the level-1 model.                      *;
   *******************************************************************;
 Y=d*B+(1-d)*(F+(L-F)*exp(-S*t));
   *******************************************************************;
   * The next line states that the residuals of the level-1 model    *;
   * are assumed to be normally distributed with variance s2e.       *;
   *******************************************************************;
MODEL cesd ~ normal(Y,s2e);
   *******************************************************************;
   * The next line states that the level-2 random effects are        *;
   * multivariate normal and specifies the lower triangle of the     *;
   * variance-covariance matrix.                                     *;
   *******************************************************************;
   * The out=rfx option outputs estimates of zetaB, zetaL, and zetaF *;
   * for each person to a new data set called rfx. These estimates    *;
   * can be used to plot individual curves, as in Figure 4.          *;
   *******************************************************************;
RANDOM zetaB zetaL zetaF ~ normal([0,0,0],[s2B,sBL,s2L,sBF,sLF,s2F])
       subject=id out=rfx;
   *******************************************************************;
   * The next three lines request that estimates of the random       *;
   * effect correlations (as well as standard errors and tests) be   *;
   * displayed. ESTIMATE statements can be used flexibly to estimate *;
   * any expression that is a combination of model parameters.       *;
   *******************************************************************;
ESTIMATE 'corr(B,L)' sBL/sqrt(s2B)/sqrt(s2L);
ESTIMATE 'corr(B,F)' sBF/sqrt(s2B)/sqrt(s2F);
ESTIMATE 'corr(L,F)' sLF/sqrt(s2L)/sqrt(s2F);
RUN;