Bioequivalence and Bioavailability Forum 13:12 CET

Main page Policy/Terms of Use Abbreviations Latest Posts

 Log in |  Register |  Search

Diagnostics: R [Study As­sess­ment]

posted by zizou - Plzeň, Czech Republic, 2016-05-22 19:07  - Posting: # 16348
Views: 12,507

Dear Helmut.
» Estimates and their SEs are exactly the same. CIs are not (due to different DFs?).
Exactly, different DFs.

# equivalent of your code:
library(lme4)
Subj   <- c(1, 2, 4, 5, 6, 4, 5, 6, 7, 8, 9, 7, 8, 9)
Dose   <- c(25, 25, 50, 50, 50, 250, 250, 250, 75, 75, 75, 250, 250, 250)
AUC    <- c(326.40, 437.82, 557.47, 764.85, 943.59, 2040.84, 2989.29,
            4107.58, 1562.42, 982.02, 1359.68, 3848.86, 4333.10, 3685.55)
Cmax   <- c(64.82, 67.35, 104.15, 143.12, 243.63, 451.44, 393.45,
            796.57, 145.13, 166.77, 296.90, 313.00, 387.00, 843.00)
resp   <- data.frame(Subj, Dose, Cmax, AUC)
muddle <- lmer(log(Cmax) ~ log(Dose) + (1|Subj), data=resp)
summary(muddle)


DF, p-values, CIs not reported in lmer for the reasons stated there.

For more decimal places:
print(muddle, digits=7, ranef.comp=c("Var","Std.Dev."))
# REML criterion at convergence: 6.2435 (same as in SPSS 6.243548)
# note: -2 (Restricted) Log Likelihood from R lme: (-2)*(-3.121774) = 6.243548


Some reference for lmer can be found in this PDF.

For DFs, p-values, ...: library lmerTest can be used, see this PDF (page 2 - description, details, references to SAS).
library(lmerTest) # Library includes modification of function lmer (if I am not mistaken).
muddle <- lmer(log(Cmax) ~ log(Dose) + (1|Subj), data=resp)
summary(muddle)


Results:
Linear mixed model fit by REML t-tests use Satterthwaite approximations to
  degrees of freedom [lmerMod] Formula: log(Cmax) ~ log(Dose) + (1 | Subj)
   Data: resp

REML criterion at convergence: 6.2

Scaled residuals:
     Min       1Q   Median       3Q      Max
-1.07548 -0.35579 -0.03301  0.45089  0.91854

Random effects:
 Groups   Name        Variance Std.Dev.
 Subj     (Intercept) 0.11205  0.3347 
 Residual             0.01456  0.1207 
Number of obs: 14, groups:  Subj, 8

Fixed effects:
            Estimate Std. Error      df t value Pr(>|t|)   
(Intercept)  1.94139    0.24314 9.19600   7.985 1.98e-05 ***
log(Dose)    0.76174    0.04728 5.89600  16.111 4.24e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
          (Intr)
log(Dose) -0.863

I think there is a little bug in df visualization - rounded on 3 decimal places, visible on 5 (actually I have not the latest version of R).
summary(muddle)$coefficients["(Intercept)","df"]  # [1] 9.195607

For calculation of 90% CIs: (below is only example for Intercept lower and upper limit)
alpha=0.05
summary(muddle)$coefficients["(Intercept)","Estimate"] - qt(1-alpha,summary(muddle)$coefficients["(Intercept)","df"]) * summary(muddle)$coefficients["(Intercept)","Std. Error"]
summary(muddle)$coefficients["(Intercept)","Estimate"] + qt(1-alpha,summary(muddle)$coefficients["(Intercept)","df"]) * summary(muddle)$coefficients["(Intercept)","Std. Error"]


Best regards,
zizou

Complete thread:

Activity
 Mix view
Bioequivalence and Bioavailability Forum |  Admin contact
19,032 posts in 4,059 threads, 1,299 registered users;
online 13 (0 registered, 13 guests [including 10 identified bots]).

When a distinguished but elderly scientist states that
something is possible, he is almost certainly right.
When he states that something is impossible,
he is very probably wrong.    Arthur C. Clarke

The BIOEQUIVALENCE / BIOAVAILABILITY FORUM is hosted by
BEBAC Ing. Helmut Schütz
HTML5 RSS Feed