1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
|
test_that("pooled parameters", {
skip_if_not_installed("mice")
data("nhanes2", package = "mice")
set.seed(123)
imp <- mice::mice(nhanes2, printFlag = FALSE)
models <- lapply(1:5, function(i) {
lm(bmi ~ age + hyp + chl, data = mice::complete(imp, action = i))
})
pp <- pool_parameters(models)
expect_equal(pp$df_error, c(9.2225, 8.1903, 3.6727, 10.264, 6.4385), tolerance = 1e-3)
expect_snapshot(print(pp))
})
test_that("pooled parameters", {
skip_if_not_installed("mice")
skip_if_not_installed("datawizard")
data("nhanes2", package = "mice")
nhanes2$hyp <- datawizard::slide(as.numeric(nhanes2$hyp))
set.seed(123)
imp <- mice::mice(nhanes2, printFlag = FALSE)
models <- lapply(1:5, function(i) {
glm(hyp ~ age + chl, family = binomial, data = mice::complete(imp, action = i))
})
pp1 <- pool_parameters(models)
expect_equal(pp1$df_error, c(Inf, Inf, Inf, Inf), tolerance = 1e-3)
pp2 <- pool_parameters(models, ci_method = "residual")
m_mice <- with(data = imp, exp = glm(hyp ~ age + chl, family = binomial))
pp3 <- summary(mice::pool(m_mice))
expect_equal(pp2$df_error, pp3$df, tolerance = 1e-3)
})
skip_on_cran()
test_that("pooled parameters, glmmTMB, components", {
skip_if_not_installed("mice")
skip_if_not_installed("glmmTMB")
sim1 <- function(nfac = 4, nt = 10, facsd = 0.1, tsd = 0.15, mu = 0, residsd = 1) {
dat <- expand.grid(fac = factor(letters[1:nfac]), t = 1:nt)
n <- nrow(dat)
dat$REfac <- rnorm(nfac, sd = facsd)[dat$fac]
dat$REt <- rnorm(nt, sd = tsd)[dat$t]
dat$x <- rnorm(n, mean = mu, sd = residsd) + dat$REfac + dat$REt
dat
}
set.seed(101)
d1 <- sim1(mu = 100, residsd = 10)
d2 <- sim1(mu = 200, residsd = 5)
d1$sd <- "ten"
d2$sd <- "five"
dat <- rbind(d1, d2)
set.seed(101)
dat$REfac[sample.int(nrow(dat), 10)] <- NA
dat$x[sample.int(nrow(dat), 10)] <- NA
dat$sd[sample.int(nrow(dat), 10)] <- NA
impdat <- suppressWarnings(mice::mice(dat, printFlag = FALSE))
models <- lapply(1:5, function(i) {
glmmTMB::glmmTMB(
x ~ sd + (1 | t),
dispformula = ~sd,
data = mice::complete(impdat, action = i)
)
})
out <- pool_parameters(models, component = "conditional")
expect_named(
out,
c(
"Parameter", "Coefficient", "SE", "CI_low", "CI_high", "Statistic",
"df_error", "p"
)
)
expect_equal(out$Coefficient, c(187.280225, -87.838969), tolerance = 1e-3)
out <- suppressMessages(pool_parameters(models, component = "all", effects = "all"))
expect_named(
out,
c(
"Parameter", "Coefficient", "Effects", "SE", "CI_low", "CI_high",
"Statistic", "df_error", "p", "Component"
)
)
expect_equal(
out$Coefficient,
c(187.280225, -87.838969, 3.51576, -1.032665, 0.610992, NaN),
tolerance = 1e-3
)
out <- pool_parameters(models, component = "all", effects = "fixed")
expect_named(
out,
c(
"Parameter", "Coefficient", "SE", "CI_low", "CI_high",
"Statistic", "df_error", "p", "Component"
)
)
expect_equal(
out$Coefficient,
c(187.280225, -87.838969, 3.51576, -1.032665),
tolerance = 1e-3
)
})
test_that("pooled parameters, glmmTMB, zero-inflated", {
skip_if_not_installed("mice")
skip_if_not_installed("glmmTMB")
skip_if_not_installed("broom.mixed")
data(Salamanders, package = "glmmTMB")
set.seed(123)
Salamanders$cover[sample.int(nrow(Salamanders), 50)] <- NA
Salamanders$mined[sample.int(nrow(Salamanders), 10)] <- NA
impdat <- suppressWarnings(mice::mice(Salamanders, printFlag = FALSE))
models <- lapply(1:5, function(i) {
glmmTMB::glmmTMB(
count ~ mined + cover + (1 | site),
ziformula = ~mined,
family = poisson(),
data = mice::complete(impdat, action = i)
)
})
out <- pool_parameters(models, ci_method = "residual")
expect_named(
out,
c(
"Parameter", "Coefficient", "SE", "CI_low", "CI_high", "Statistic",
"df_error", "p", "Component"
)
)
expect_equal(
out$Coefficient,
c(0.13409, 1.198551, -0.181912, 1.253029, -1.844026),
tolerance = 1e-3
)
# validate against mice ---------------
m_mice <- suppressWarnings(with(data = impdat, exp = glmmTMB::glmmTMB(
count ~ mined + cover + (1 | site),
ziformula = ~mined,
family = poisson()
)))
mice_summ <- summary(mice::pool(m_mice, dfcom = Inf))
expect_equal(out$Coefficient, mice_summ$estimate, tolerance = 1e-3)
expect_equal(out$SE, mice_summ$std.error, tolerance = 1e-3)
expect_equal(out$p, mice_summ$p.value, tolerance = 1e-3)
out <- pool_parameters(models, component = "all", effects = "all")
expect_named(
out,
c(
"Parameter", "Coefficient", "Effects", "SE", "CI_low", "CI_high",
"Statistic", "df_error", "p", "Component"
)
)
expect_equal(
out$Coefficient,
c(0.13409, 1.198551, -0.181912, 1.253029, -1.844026, 0.158795),
tolerance = 1e-3
)
out <- pool_parameters(models, component = "conditional", effects = "fixed")
expect_named(
out,
c(
"Parameter", "Coefficient", "SE", "CI_low", "CI_high",
"Statistic", "df_error", "p"
)
)
expect_equal(
out$Coefficient,
c(0.13409, 1.198551, -0.181912),
tolerance = 1e-3
)
})
|
