Matt Jacob
2017-Feb-16 01:14 UTC
[R] Validating Minitab's "Expanded Gage R&R Study" using R and lme4
I'm trying to validate the results of an "Expanded Gage R&R Study" in Minitab using R and lme4, but I can't get the numbers to match up in certain situations. I can't tell whether my model is wrong, my data is bad, or something else is going on. For instance, here's some data for which the results don't match: https://i.stack.imgur.com/5PCgm.png After running the gage study, these are the results according to Minitab: Study Var %Study Var %Tolerance Source StdDev (SD) (6 * SD) (%SV) (SV/Toler) Total Gage R&R 1.76277 10.5766 100.00 14.36 Repeatability 0.00000 0.0000 0.00 0.00 Reproducibility 1.76277 10.5766 100.00 14.36 B 0.00000 0.0000 0.00 0.00 A*B 1.76277 10.5766 100.00 14.36 Part-To-Part 0.00000 0.0000 0.00 0.00 A 0.00000 0.0000 0.00 0.00 Total Variation 1.76277 10.5766 100.00 14.36 But when I mimic Minitab's results by parsing the output from lmer() and doing the arithmetic in Excel, this is what I see: https://i.stack.imgur.com/EGg9F.png The raw output from lmer() was: Linear mixed model fit by REML ['lmerMod'] Formula: y ~ 1 + (1 | A) + (1 | B) + (1 | A:B) Data: d REML criterion at convergence: -100.1 Scaled residuals: Min 1Q Median 3Q Max -1.308e-07 -1.308e-07 -1.308e-07 -6.541e-08 1.308e-07 Random effects: Groups Name Variance Std.Dev. A:B (Intercept) 1.333e+00 1.154e+00 B (Intercept) 7.066e-04 2.658e-02 A (Intercept) 2.260e-03 4.754e-02 Residual 2.655e-14 1.629e-07 Number of obs: 8, groups: A:B, 4; B, 2; A, 2 Fixed effects: Estimate Std. Error t value (Intercept) 52.17 0.57 91.53 convergence code: 0 Model failed to converge with max|grad| = 0.422755 (tol = 0.002, component 1) Model is nearly unidentifiable: very large eigenvalue - Rescale variables? And the R code that produced that output is: library(lme4) A <- factor(c(1, 1, 2, 2, 2, 1, 2, 1)) B <- factor(c(1, 2, 1, 2, 1, 2, 2, 1)) y <- c(51.356124843620798, 51.356124843620798, 54.8816618912481, 51.356124843620798, 54.8816618912481, 51.356124843620798, 51.356124843620798, 51.356124843620798) d <- data.frame(y, A, B) fm <- lmer(y ~ 1 + (1|A) + (1|B) + (1|A:B), d) summary(fm) For a different measurement with a different response, it's a completely different situation! Given the following data: https://i.stack.imgur.com/cH0bO.png The resulting table from Minitab is: Study Var %Study Var %Tolerance Source StdDev (SD) (6 * SD) (%SV) (SV/Toler) Total Gage R&R 0.193649 1.16190 55.90 1.00 Repeatability 0.093541 0.56125 27.00 0.48 Reproducibility 0.169558 1.01735 48.95 0.88 B 0.132288 0.79373 38.19 0.68 A*B 0.106066 0.63640 30.62 0.55 Part-To-Part 0.287228 1.72337 82.92 1.49 A 0.287228 1.72337 82.92 1.49 Total Variation 0.346410 2.07846 100.00 1.79 And after plugging my R results into Excel, I get exactly the same thing: https://i.stack.imgur.com/jUEAP.png Which was produced by this R code: library(lme4) A <- factor(c(1, 1, 2, 2, 2, 1, 2, 1)) B <- factor(c(1, 2, 1, 2, 1, 2, 2, 1)) y <- c(-49.4, -49.8, -50.1, -50.1, -50.0, -49.9, -50.2, -49.6) d <- data.frame(y, A, B) fm <- lmer(y ~ 1 + (1|A) + (1|B) + (1|A:B), d) summary(fm) That generated the following lmer() summary: Linear mixed model fit by REML ['lmerMod'] Formula: y ~ 1 + (1 | A) + (1 | B) + (1 | A:B) Data: d REML criterion at convergence: -3.8 Scaled residuals: Min 1Q Median 3Q Max -0.7705 -0.6853 -0.1039 0.4379 1.4151 Random effects: Groups Name Variance Std.Dev. A:B (Intercept) 0.01125 0.10607 B (Intercept) 0.01750 0.13229 A (Intercept) 0.08250 0.28723 Residual 0.00875 0.09354 Number of obs: 8, groups: A:B, 4; B, 2; A, 2 Fixed effects: Estimate Std. Error t value (Intercept) -49.8875 0.2322 -214.9 Is the difference attributable to the warnings produced by lmer() about the model failing to converge and being nearly unidentifiable? What could Minitab be doing differently when the measurement data contains only two distinct values? Matt This question is cross-posted to http://stats.stackexchange.com/questions/262170/how-can-i-validate-minitabs-expanded-gage-rr-study-using-open-source-tools