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
______________________________________________
R-help@r-project.org 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.
