topic Re: A joint model using PROC NLMIXED in Statistical Procedures

A joint model using PROC NLMIXED

MunaSmith — Wed, 28 Jun 2017 16:50:06 GMT

Hello,

I have a problem in fitting the following joint model for binary and longitudinal outcome using NLMIXED, the concept that I want to predict the probability of developing a medical condition (which is HT in this model ) based on a set of weekly measuremnts of a biomarker (Y). the joint model consists of two parts the first one is a mixed model to summarize the biomarker and a logistic model to obtain the predicted probabilties .

proc nlmixed data=data qpoints=1;
parameters beta0= 2.5 beta1= -0.22
a11= 0.34 a12= -0.01 a22=0.03
alpha0=10 alpha1=10 alpha2=10 s2=0.03;

*the variance covariance matrix;
v11 = a11*a11;
v12 = a11*a12;
v22 = a12*a12 + a22*a22;

*The longitudinal model ;
linplong = (beta0 + u1) + (beta1 + u2)*time;
resid = (y-linplong);

lllong = -0.5*(1.837876 + resid**2 / s2 + log(s2));

*The binary model;
xb=alpha0+alpha1*u1+alpha2*u2;
prob = exp(xb)/(1+exp(xb));
liklhd = (prob**HT)*((1-prob)**(1-HT));
llbin = log(liklhd);

model HT ~ general(lllong + llbin);
random u1 u2 ~ normal([0, 0],[v11,v12,v22]) subject=id;

run;

The above model fails to converge and I get the following error msg

ERROR: No valid parameter points were found

any ideas of what might be wrong in the above code !

Thanks in advance

The first 5 subjects in the data

id y HT time
9 2,1 0 2
9 1,42 0 3
9 1,08 0 4
9 0,95 0 5
9 0,73 0 6
9 0,34 0 7
10 NA 0 2
10 NA 0 3
10 NA 0 4
10 NA 0 5
10 NA 0 6
10 2,15 0 7
11 NA 0 2
11 NA 0 3
11 NA 0 4
11 NA 0 5
11 NA 0 6
11 NA 0 7
12 NA 0 2
12 NA 0 3
12 1,67 0 4
12 1,6 0 5
12 1,51 0 6
12 1,04 0 7
13 NA 0 2
13 NA 0 3
13 NA 0 4
13 0,11 0 5
13 NA 0 6
13 NA 0 7
14 2,56 1 2
14 2,08 1 3
14 1,88 1 4
14 1,79 1 5
14 1,74 1 6
14 NA 1 7
15 NA 1 2
15 NA 1 3
15 NA 1 4
15 2,41 1 5
15 2,43 1 6
15 NA 1 7

Re: A joint model using PROC NLMIXED

Rick_SAS — Fri, 30 Jun 2017 13:44:00 GMT

First, check that the data set is created properly. One way is to plot the response for each subject over time:

proc sgplot data=data noautolegend;
series x=time y=y / group=id grouplc=HT;  /* GROUPLC= option requires recent version of SAS */
run;

Next, I think you need to use the BOUNDS statement to bound some of the variance parameters. Lastly, try different initial values for the regression parameters. For example:

proc nlmixed data=data;
parameters beta0= 2.5 beta1= -0.22
a11= 0.34 a12= -0.01 a22=0.03
alpha0=1 alpha1=1 alpha2=1 s2=0.03;
bounds 0 < a11, a12, a22, s2;

Re: A joint model using PROC NLMIXED

MunaSmith — Tue, 11 Jul 2017 13:27:21 GMT

Thank you for the reply

I tried your suggestion and tried a grid of starting values but still didn't work!

I was wondering if this likelihood is so complicated to evaluated by NLMIXED!

Re: A joint model using PROC NLMIXED

Rick_SAS — Tue, 11 Jul 2017 13:42:47 GMT

Most non-convergence issues are caused because the model is either specified incorrectly or the model does not fit the data.

If you post the NLMIXED code and the data AS A SAS DATA STEP, the issue will be easier to investigate.

Re: A joint model using PROC NLMIXED

MunaSmith — Wed, 12 Jul 2017 19:53:37 GMT

Thank you !

I uploaded the data set in attachemnt and here is the final code

Libname cls "C:\Users\user\Desktop";

proc import datafile="C:\Users\user\Desktop\book1.csv"
OUT= cls.data
DBMS=csv
REPLACE;
GETNAMES=YES;
RUN;*Please find the attached csv file;

proc nlmixed data=cls.data ;
parameters beta0= 2.5 beta1=-0.22
a11= 0.59 a12= -0.11 a22=0.18
alpha0=-1.38 alpha1=0.42 alpha2=21.32 s2=0.03;
bounds 0 < a11, a22, s2;

*The longitudinal model;
linplong = (beta0 + u1) + (beta1 + u2)*time;
resid = (y-linplong);
if (abs(resid) > 1.3E100) or (s2 < 1e-12) then do;
lllong = -1e20;
end; else do;
lllong = -0.5*(1.837876 + resid**2 / s2 + log(s2));
end;
*The logistic model;
xb=alpha0+alpha1*u1+alpha2*u2;
prob = exp(xb)/(1+exp(xb));
liklhd = (prob**HT)*((1-prob)**(1-HT));
llbin = log(liklhd);

model HT ~ general(lllong + llbin);
random u1 u2 ~ normal([0, 0],[a11,a12,a22]) subject=id;

run;
title1;