[R] Proper syntax for using varConstPower in nlme

Michael A. Gilchrist mikeg at utk.edu
Thu Oct 15 21:05:50 CEST 2009


Hello,

Excuse me for posting two questions in one day, but I figured it would be 
better to ask my questions in separate emails.  I will again give the caveat 
that I'm not a statistician by training, but have a fairly decent 
understanding of probability and likelihood.

As before, I'm trying to fit a nonlinear model to a dataset which has two main 
factors using nlme. Within the dataset there are two Type categories and four 
Tissue categories, thus giving me 8 datasets in total.  The dataset is in a 
groupedData object with
 	       formula= Count ~ Time|Type/Tissue
and there are three basic model parameters that I am trying to fit: T0, aN, 
and aL.

Calling nlsList gives
------------------------------------------
>nlsList(Count ~ quad.PBMC.model(aL, aN, T0), data = tissueData, start = 
list(T0 = 1000, aN = exp(-2), aL = exp(-2)))
Call:
   Model: Count ~ quad.PBMC.model(aL, aN, T0) | Type/Tissue
    Data: tissueData

Coefficients:
                     T0         aN         aL
Naive/CLN    1530.0088 0.26508876 0.04730525
.
.
.
Memory/IngLN 2114.1868 0.05298126 0.12795589
Memory/MLN   1050.0811 0.02277018 0.13560749

Degrees of freedom: 96 total; 72 residual
Residual standard error: 271.4194
>
---------------------------------------------

I find that T0 varies greatly with each treatment, so I'm going to leave that 
alone for now. However, aN and aL values don't seem to vary w/in a Type or 
between Types.  As a result, I would like a mixed effects model using nlme.

Further, looking at the residuals, I find that they are heteroscedastic. As a 
result, I would like to try and model the variance in the data using 
varConstPowerFun within nlme.  I've been trying to understand how to use this 
option by reading Pinheiro and Bates's book on mixed effects models.  Based on 
this, I've tried using the syntax,

---------------------------------------------
>  nlme(Count ~ quad.PBMC.model(aL, aN, T0),
+   data = tissueData,
+   weights = varConstPower(form =~ Count),
+   start = list( fixed = c(rep(1000, 8), -2, -2) ),
+   fixed = list(T0 ~ TypeTissue-1, aL ~ 1, aN ~ 1),
+   random = aL + aN ~ 1|Tissue
+   )
Error in MEestimate(nlmeSt, grpShrunk) :
   Singularity in backsolve at level 0, block 1
>
>
---------------------------------------------
The above command clearly this doesn't work, but if I comment out the "weights 
= ..." line, it executes with no problem.

I'd greatly appreciate any guidance on how to properly set up varConstPower.


Sincerely,

Mike


-----------------------------------------------------
Department of Ecology & Evolutionary Biology
569 Dabney Hall
University of Tennessee
Knoxville, TN 37996-1610

phone:(865) 974-6453
fax:  (865) 974-6042

web: http://eeb.bio.utk.edu/gilchrist.asp




More information about the R-help mailing list