[R] varFixed

Spencer Graves spencer.graves at pdf.com
Mon Dec 22 20:37:53 CET 2003


      If I had several days to work on this, I'd study Pinheiro and 
Bates (2000) Mixed-Effects Models in S and S-Plus (Springer) to see if I 
could use "lme or "gls" on this.  If I expected to encounter many 
similar problems in the future, I might modify "lme" to accept a prior 
distribution over the parameters to be estimated. 

      However, if I needed an answer today for a problem I might not see 
again for a while, I might just write a log likelihood and give it to 
"optim", something like the following:  The following looks to me like 
what you are describing, possibly simplified: 

      y[i,j] = mu + a[i] + e[i,j], where a[i] ~ N(0, s.a^2) and e[i,j] ~ 
N(0, 1), i = 1, ..., n, j = 1, ..., m[i].  

      If this is correct, then this is equivalent to the following: 

      y ~ N(mu*One, Sig), where y = vector of all observations, One = 
vector of all 1's, and Sig = s.a.^2*diag(sum(m)) + SigW, where SigW = 
block diagonal matrix with i-th block = m[i] x m[i] matrix of all 1's.  
Then the "deviance" = (-2)*log(likelihood) can be written as follows: 

      Deviance1 = 
sum(m)*log(2*pi)+log(det(Sig))+t(y-mu)%*%inverse(Sig)%*%(y-mu). 

      In another hour, I'd have a function written to compute 
"Deviance1" in terms of mu and log(s.a) -- not s.a directly.  I'd have 
the minimum deviance + hessian / information from "optim" and a contour 
plot of "Deviance1" in the regions.  From this I could get confidence 
regions using 2*log(likelihood ratio) is approximately chi-square, etc.  
In a couple of days, I could also do a Monte Carlo study to evaluate the 
accuracy of the normal approximations, etc. 

      hope this helps. 
      spencer graves

, esp. pp. 249-

Thomas Lumley wrote:

>On Sat, 20 Dec 2003, Harold Doran wrote:
>
>  
>
>>Dear List:
>> Earlier this week I posted a question and received no response, and I
>>continue to struggle with my model. My original question is pasted
>>below.
>> I am using lme and want to fix the variance of the within group
>>residual at 1 (e~n(0,1). I think the varFixed function should be used to
>>accomplish this, but I am struggling to figure out how to do this.
>>
>>    
>>
>
>I don't think this is possible.  I've tried to get this sort of effect
>(for meta-analyses) and have not been able to. I think there is always an
>overall scale factor for variances estimated.
>
>	-thomas
>
>______________________________________________
>R-help at stat.math.ethz.ch mailing list
>https://www.stat.math.ethz.ch/mailman/listinfo/r-help
>  
>




More information about the R-help mailing list