[R] HELP with NLME

[Spencer Graves]

	  Have you tried fitting the model with "verbose=TRUE" and possibly 
also "control=nlmeControl(msVerbose=TRUE)"?

	  Also, have you consulted Pinheiro and Bates (2000) Mixed-Effects 
Models in S and S-PLUS (Springer)?  That book contains several examples 
of use of lme and nlme.  R script files named "ch01.R", "ch02.R", ..., 
"ch06.R", "ch08.R" can be found in "~R\library\nlme\scripts", where 
"~R\" is your R installation directory.  If you work through these files 
one line at a time until you think you understand it all, you might get 
an answer to your question.

	  If this fails to produce answer your question, I would try to find a 
much simpler and self contained example that still returns your error 
message.  I might do this by removing terms from the model, removing 
parameters from the model, etc., trying to write a very few lines of 
code that would produce a data set that will still generate the error 
message you see.

	  Hope this helps.
	  Spencer Graves
p.s.  Is it feasible for you to upgrade to R2.3.1?  R2.1 is by now 
rather old.

Loki Natarajan wrote:
> Hi,
> 
> I was very much hoping someone could help me with the following.
> I am trying to convert some SAS NLMIXED code to NLME in R (v.2.1),
> but I get an error message. Does anyone have any suggestions?
> I think my error is with the random effect "u" which seems to be
> parametrized differently in the SAS code. In case it's helpful, 
> what I am essentially trying to do is estimate parameters using ML in a 
> measurement error setting with some validation data (indicated by 
> vs.flag). Any help would be greatly appreciated.I apologize for the
> clumsiness of the R code.
> Many thanks in advance.
> 
> Sincerely,
> Loki
> #################################################################
> SAS Code:
> 
> proc nlmixed data=repdat
> parms b0 -3 b1 -.135 a0 3 a1 4 sigsq 0.25;
> 
> if vs.flag = 1 then do;
>     eta1 = b0 + b1*lnbldT;
>     llbin = anybc.cens.ind*eta1 - log(1+exp(eta1));
>     eta2 = a0 + a1*lnndsTs;
>     llnorm = -1/(2*sigsq)*(lnbldT - eta2)**2 - .5*log(sigsq);
>     ll = llbin + llnorm;
>     end;
> 
> else do;
>     eta2 = a0 + a1*lnndsTs;
>     eta1 = b0 + b1*eta2 + u;
>     llbin = anybc.cens.ind*eta1 - log(1+exp(eta1));
>     ll = llbin;
>     end;
> 
> sigma2 = sigmasq*b1**2 /*variance of random effect;
> model anybc.cens.indic ~ general(ll);
> random u ~ normal(0,sigma2) subject = CaseID; */;
> run;
> ####################################################################
> R Code with error message:
> 
> me.km.nlme <- nlme(model = anybc.cens.indic ~ 
> vs.flag*((anybc.cens.indic*(b0+b1*lnbldT) - log(1 + exp(b0+b1*lnbldT))) + 
> (((-1/(2*sigsq))*(lnbldT -a0 -a1*lnndsTs)^2) - 0.5*log(sigsq))) + 
> (1-vs.flag)*((anybc.cens.indic*(b0+b1*(a0 + a1*lnndsTs + u))) - log(1 + 
> exp(b0+b1*(a0 + a1*lnndsTs + u)))), fixed = 
> list(a0~1,a1~1,b0~1,b1~1,sigsq~1),na.action=na.omit, data=rc.df, 
> method="ML", random=u~1|CaseID,
> start = c(a0=0, a1=1, b0=-3, b1 = -0.135, sigsq = 0.25))
> 
> + Error: Singularity in backsolve at level 0, block 1
> In addition: There were 16 warnings (use warnings() to see them)
>> warnings()
> Warning messages:
> 1: Singular precision matrix in level -1, block 5
> 2: Singular precision matrix in level -1, block 5
> 3: Singular precision matrix in level -1, block 5
> 4: NA/Inf replaced by maximum positive value
> 5: Singular precision matrix in level -1, block 5
> 6: Singular precision matrix in level -1, block 5
> 7: Singular precision matrix in level -1, block 5
> 8: NA/Inf replaced by maximum positive value
> 9: NaNs produced in: log(x)
> 10: NaNs produced in: log(x)
> 11: NaNs produced in: log(x)
> 12: NaNs produced in: log(x)
> 13: NaNs produced in: log(x)
> 14: NaNs produced in: log(x)
> 15: NaNs produced in: log(x)
> 16: NaNs produced in: log(x)
> #####################################################################
> 
> 
> Loki Natarajan
> Associate Professor of Biostatistics
> Moores UCSD Cancer Center
> 3855 Health Sciences Drive #0901
> La Jolla, CA 92093-0901
> 
> 
> phone:     858 822 4763
> Fax:       858 822 6897
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.