Semireplicated + REML in R [🇷 for BE/BA]

posted by Helmut Homepage – Vienna, Austria, 2020-07-10 21:03 (1555 d 08:15 ago) – Posting: # 21670
Views: 20,478

Hi ElMaestro,

❝ I am writing a little optimizer based on REML in R for use with semi-replicated designs. I am starting from the bottom and working upwards, since the various packages do not really help when R is replicated and T is not.


Wow! Nobody mastered that in R so far. I don’t want to dampen your enthusiasm but I think that is not possible at all.

❝ Does anyone have a dataset with well-established REML reference values for the variance components for the TRR/RRT/RTR design? A little more than the CVref of EMA dataset II would be appreciated :-)


Here you are (for DS II, of course). Results in Phoenix for the FDA’s covariance structure „Banded No-Diagonal Factor Analytic (f=2)” which is in SAS-lingo FA0(2). After 7 iterations (-2REML -48.2, AIC -26.2):
lambda(1,1) 0.19030501
lambda(1,2) 0.23981134
lambda(2,2) 0.074350441
s2wR        0.013246493
s2wT        0.0074570998
G Matrix*
0.03621600 0.04563730
0.04563730 0.06303747
As usual:
WARNING: Newton's algorithm converged with modified Hessian. Output is suspect.
Model may be over-specified. A simpler model could be tried.
s2wT is nonsense, since T is not replicated and the model is over-specified.
The lambdas are called in SAS-lingo FA(x,y). After 7 iterations (-2REML -30.7, AIC -20.7):
Cov Parm     Subject    Group    Estimate
FA(1,1)      Subj                  0.1903
FA(2,1)      Subj                  0.2398
FA(2,2)      Subj                  0.1072
Residual     Subj       trt R      0.01325
Residual     Subj       trt T      0.001497

    Estimated G Matrix
Row  trt    Col1      Col2
 1    R   0.03622   0.04564
 2    T   0.04564   0.06900
Also:
NOTE: Convergence criteria met but final hessian is not positive definite.
NOTE: Asymptotic variance matrix of covariance parameter estimates has been found to be
      singular and a generalized inverse was used. Covariance parameters with zero variance do
      not contribute to degrees of freedom computed by DDFM=SATTERTH.

No problems with convergence with f=1 (or FA0(1) in SAS). FA0(1) helps in general but not always.
lambda(1,1) 0.19030501
lambda(1,2) 0.23981134
s2wR        0.013246493
s2wT        0.012985088
Different s2wT but crap as well.

[image]You can get anything you like. The same data set throwing warnings in one software but not the other; same with convergence. In Certara’s forum a user posted a data set which failed to converge both in SAS and Phoenix, even with FA0(1). Then one is in deep shit.

Even if a data set converges without warnings, you might get different estimates (except the first one and s2wR) in SAS and Phoenix. Not surprising cause the optimizer walks around in parameter-space trying to find the maximum of a surface which is essentially flat. Increasing the number of iterations and/or decreasing the tolerance rarely helps.

I still wait for someone explaining why he/she is using such a wacky design. Beyond me.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна! [image]
Helmut Schütz
[image]

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

Complete thread:

UA Flag
Activity
 Admin contact
23,249 posts in 4,885 threads, 1,666 registered users;
76 visitors (0 registered, 76 guests [including 8 identified bots]).
Forum time: 05:19 CEST (Europe/Vienna)

I believe there is no philosophical high-road in science,
with epistemological signposts. No, we are in a jungle
and find our way by trial and error,
building our road behind us as we proceed.    Max Born

The Bioequivalence and Bioavailability Forum is hosted by
BEBAC Ing. Helmut Schütz
HTML5