R语言如何通过绘制各种优化器直观地评估混合效应模型的收敛性

cbwuti44 于 2023-06-27 发布在其他

关注(0)|答案(2)|浏览(155)

为了评估收敛警告是否会导致结果无效，或者相反，尽管有警告，结果仍然可以被视为有效，Bates et al. (2023)建议使用各种优化器重新调整受收敛警告影响的模型。作者认为，如果不同的优化器产生实际等效的结果，结果是有效的。'lme 4'包中的allFit函数允许使用多个优化器来改装模型。要使用上面列出的七个优化器，必须安装两个额外的软件包：'dfoptim'和'optimx'（参见lme 4手册）。allFit（）的输出包含每个优化器拟合的固定和随机效应的若干统计信息。

library(lme4)
#> Loading required package: Matrix
library(dfoptim)
library(optimx)

# Create data using code by Ben Bolker from 
# https://stackoverflow.com/a/38296264/7050882

set.seed(101)
spin = runif(600, 1, 24)
reg = runif(600, 1, 15)
ID = rep(c("1","2","3","4","5", "6", "7", "8", "9", "10"))
day = rep(1:30, each = 10)
testdata <- data.frame(spin, reg, ID, day)
testdata$fatigue <- testdata$spin * testdata$reg/10 * rnorm(30, mean=3, sd=2)

# Model
fit = lmer(fatigue ~ spin * reg + (1|ID),
           data = testdata, REML = TRUE)

# Refit model using all available algorithms
multi_fit = allFit(fit)
#> bobyqa : [OK]
#> Nelder_Mead : [OK]
#> nlminbwrap : [OK]
#> nmkbw : [OK]
#> optimx.L-BFGS-B : [OK]
#> nloptwrap.NLOPT_LN_NELDERMEAD : [OK]
#> nloptwrap.NLOPT_LN_BOBYQA : [OK]

# Show results 
summary(multi_fit)$fixef
#>                               (Intercept)      spin       reg  spin:reg
#> bobyqa                          -2.975678 0.5926561 0.1437204 0.1834016
#> Nelder_Mead                     -2.975675 0.5926559 0.1437202 0.1834016
#> nlminbwrap                      -2.975677 0.5926560 0.1437203 0.1834016
#> nmkbw                           -2.975678 0.5926561 0.1437204 0.1834016
#> optimx.L-BFGS-B                 -2.975680 0.5926562 0.1437205 0.1834016
#> nloptwrap.NLOPT_LN_NELDERMEAD   -2.975666 0.5926552 0.1437196 0.1834017
#> nloptwrap.NLOPT_LN_BOBYQA       -2.975678 0.5926561 0.1437204 0.1834016

创建于2023-06-24带有reprex v2.0.2
然而，数字对普通人来说很难。我们在视觉上更好。那么，是否有一种方法可以可视化allFit（）的输出，并查看各种优化器（例如，bobyqa，Nelder-Mead等）的一组预测器的参数？

来源：https://stackoverflow.com/questions/76546377/how-to-visually-assess-the-convergence-of-a-mixed-effects-model-by-plotting-vari

2条答案

按热度按时间

mxg2im7a1#

只需在multi_fit对象上执行sapplyfixef，然后使用Map执行plot。

pframe <- as.data.frame(t(sapply(multi_fit, fixef)))

layout(matrix(c(1:5, 5), 3, 2, byrow=TRUE), heights=c(2, 2, 1.5))
par(mar=c(4, 4, 3, 3))

Map(plot, as.list(pframe), ylab=colnames(pframe), col=list(seq_len(nrow(pframe))), pch=3)
plot.new(); legend('top', ncol=2, pch=3, legend=rownames(pframe), col=seq_len(nrow(pframe)), cex=.9)

赞(0）回复(0）举报 2023-06-27

olhwl3o22#

是的，有一种方法，使用一个名为plot.fixef.allFit的自定义函数，在https://pablobernabeu.github.io/2023/a-new-function-to-plot-convergence-diagnostics-from-lme4-allfit中进行了描述。

library(lme4)
#> Loading required package: Matrix
library(dfoptim)
library(optimx)

# Create data using code by Ben Bolker from 
# https://stackoverflow.com/a/38296264/7050882

set.seed(101)
spin = runif(600, 1, 24)
reg = runif(600, 1, 15)
ID = rep(c("1","2","3","4","5", "6", "7", "8", "9", "10"))
day = rep(1:30, each = 10)
testdata <- data.frame(spin, reg, ID, day)
testdata$fatigue <- testdata$spin * testdata$reg/10 * rnorm(30, mean=3, sd=2)

# Model
fit = lmer(fatigue ~ spin * reg + (1|ID),
           data = testdata, REML = TRUE)

# Refit model using all available algorithms
multi_fit = allFit(fit)
#> bobyqa : [OK]
#> Nelder_Mead : [OK]
#> nlminbwrap : [OK]
#> nmkbw : [OK]
#> optimx.L-BFGS-B : [OK]
#> nloptwrap.NLOPT_LN_NELDERMEAD : [OK]
#> nloptwrap.NLOPT_LN_BOBYQA : [OK]

# Show results 
summary(multi_fit)$fixef
#>                               (Intercept)      spin       reg  spin:reg
#> bobyqa                          -2.975678 0.5926561 0.1437204 0.1834016
#> Nelder_Mead                     -2.975675 0.5926559 0.1437202 0.1834016
#> nlminbwrap                      -2.975677 0.5926560 0.1437203 0.1834016
#> nmkbw                           -2.975678 0.5926561 0.1437204 0.1834016
#> optimx.L-BFGS-B                 -2.975680 0.5926562 0.1437205 0.1834016
#> nloptwrap.NLOPT_LN_NELDERMEAD   -2.975666 0.5926552 0.1437196 0.1834017
#> nloptwrap.NLOPT_LN_BOBYQA       -2.975678 0.5926561 0.1437204 0.1834016

# Read in function from GitHub
source('https://raw.githubusercontent.com/pablobernabeu/plot.fixef.allFit/main/plot.fixef.allFit.R')

plot.fixef.allFit(multi_fit, 
                  
                  select_predictors = c('spin', 'reg', 'spin:reg'), 
                  
                  # Increase padding at top and bottom of Y axis
                  multiply_y_axis_limits = 1.3,
                  
                  y_title = 'Fixed effect (*b*)')
#> Loading required package: dplyr
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
#> Loading required package: reshape2
#> Loading required package: stringr
#> Loading required package: scales
#> Loading required package: ggplot2
#> Loading required package: ggtext
#> Loading required package: patchwork

# Alternative using plot-specific Y axes and other modified settings

plot.fixef.allFit(multi_fit, 
                  
                  select_predictors = c('spin', 'spin:reg'), 
                  
                  # Use plot-specific Y axis limits
                  shared_y_axis_limits = FALSE,
                  
                  decimal_points = 7, 
                  
                  # Move up Y axis title
                  y_title_hjust = 4.5,
                  
                  y_title = 'Fixed effect (*b*)')

创建于2023-06-26带有reprex v2.0.2

赞(0）回复(0）举报 2023-06-27

我来回答

R语言如何通过绘制各种优化器直观地评估混合效应模型的收敛性

2条答案

相关问题

热门标签

最新问答

R语言 如何通过绘制各种优化器直观地评估混合效应模型的收敛性

2条答案

相关问题

热门标签

最新问答

R语言如何通过绘制各种优化器直观地评估混合效应模型的收敛性