R help - Oct 2009 - Proper syntax for using varConstPower in nlme

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

Michael A. Gilchrist wrote:
> 
> ---------------------------------------------
>>  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
>>
>>
> 
You could use varPower(form=~fitted()) (the default, also varPower()). In my
experience this runs into some limitation quickly with nlme, because some
boundary conditions make convergence fail.

Try varPower(fixed = 0.5) first and play with the number.

You should only use varConstPower when you have problems with values that
cover a large range, coming close to zero, which could make varPower go
havoc.

Always do a plot of the result; the default plot gives you residual, and
some indication how to proceed.

Dieter


-- 
View this message in context:
http://www.nabble.com/Proper-syntax-for-using-varConstPower-in-nlme-tp25914277p25927578.html
Sent from the R help mailing list archive at Nabble.com.

Hi Dieter,

Thanks for the reply.  I had played with the initial conditions, but 
apparently not enough.  I finally found some that avoided the singularity 
issue.  Many thanks.

More generally, I went over the documentation yet again in P&B and I still 
find it a bit confusing.  They talk about using form = ~fitted(.) when 
discussing varPower, but the rest of the documentation seems to suggest that 
"form = ~..." should be used to indicate which covariate you assume
the
variance changes with.

Could you or someone else provide some clarification?

Thanks.

Mike

On Fri, 16 Oct 2009, Dieter Menne wrote:
>
>
>
> Michael A. Gilchrist wrote:
>>
>> ---------------------------------------------
>>>  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
>>>
>>>
>>
>
> You could use varPower(form=~fitted()) (the default, also varPower()). In
my
> experience this runs into some limitation quickly with nlme, because some
> boundary conditions make convergence fail.
>
> Try varPower(fixed = 0.5) first and play with the number.
>
> You should only use varConstPower when you have problems with values that
> cover a large range, coming close to zero, which could make varPower go
> havoc.
>
> Always do a plot of the result; the default plot gives you residual, and
> some indication how to proceed.
>
> Dieter
>
>
> -- 
> View this message in context:
http://www.nabble.com/Proper-syntax-for-using-varConstPower-in-nlme-tp25914277p25927578.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at 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.
>