Description Details Author(s) References See Also Examples
Package simex
is an implementation of the SIMEX–algorithm by Cook
and Stephanski and the MCSIMEX–Algorithm by Küchenhoff, Mwalili and Lesaffre.
Package: | simex |
Type: | Package |
Version: | 1.8 |
Date: | 2019-07-28 |
License: | GPL 2 or above |
LazyLoad: | yes |
The package includes first of all the implementation for the SIMEX– and MCSIMEX–Algorithms. Jackknife and asymptotic variance estimation are implemented. Various methods and analytic tools are provided for a simple and fast access to the SIMEX– and MCSIMEX–Algorithm.
Functions simex()
and mcsimex()
can be used on models issued
from lm()
, glm()
with asymtotic estimation.
Models from nls()
, gam()
(package mgcv),
polr()
(package MASS),
lme()
, nlme()
(package nlme) and coxph()
(package survival) can also be corrected with
these algorithms, but without asymptotic estimations.
Wolfgang Lederer, Heidi Seibold, Helmut Küchenhoff
Maintainer: Wolfgang Lederer,wolfgang.lederer@gmail.com
Lederer, W. and Küchenhoff, H. (2006) A short introduction to the SIMEX and MCSIMEX. R News, 6/4, 26 – 31
and for functions generating the initial naive models:
lm
, glm
, nls
,
gam
, lme
, nlme
, polr
, coxph
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 | # See example(simex) and example(mcsimex)
## Seed
set.seed(49494)
## simulating the measurement error standard deviations
sd_me1 <- 0.3
sd_me2 <- 0.4
temp <- runif(100, min = 0, max = 0.6)
sd_me_het1 <- sort(temp)
temp2 <- rnorm(100, sd = 0.1)
sd_me_het2 <- abs(sd_me_het1 + temp2)
## simulating the independent variables x (real and with measurement error):
x_real1 <- rnorm(100)
x_real2 <- rpois(100, lambda = 2)
x_real3 <- -4*x_real1 + runif(100, min = -2, max = 2) # correlated to x_real
x_measured1 <- x_real1 + sd_me1 * rnorm(100)
x_measured2 <- x_real2 + sd_me2 * rnorm(100)
x_het1 <- x_real1 + sd_me_het1 * rnorm(100)
x_het2 <- x_real3 + sd_me_het2 * rnorm(100)
## calculating dependent variable y:
y1 <- x_real1 + rnorm(100, sd = 0.05)
y2 <- x_real1 + 2*x_real2 + rnorm(100, sd = 0.08)
y3 <- x_real1 + 2*x_real3 + rnorm(100, sd = 0.08)
### one variable with homoscedastic measurement error
(model_real <- lm(y1 ~ x_real1))
(model_naiv <- lm(y1 ~ x_measured1, x = TRUE))
(model_simex <- simex(model_naiv, SIMEXvariable = "x_measured1", measurement.error = sd_me1))
plot(model_simex)
### two variables with homoscedastic measurement errors
(model_real2 <- lm(y2 ~ x_real1 + x_real2))
(model_naiv2 <- lm(y2 ~ x_measured1 + x_measured2, x = TRUE))
(model_simex2 <- simex(model_naiv2, SIMEXvariable = c("x_measured1", "x_measured2"),
measurement.error = cbind(sd_me1, sd_me2)))
plot(model_simex2)
### one variable with increasing heteroscedastic measurement error
model_real
(mod_naiv1 <- lm(y1 ~ x_het1, x = TRUE))
(mod_simex1 <- simex(mod_naiv1, SIMEXvariable = "x_het1",
measurement.error = sd_me_het1, asymptotic = FALSE))
plot(mod_simex1)
## Not run:
### two correlated variables with heteroscedastic measurement errors
(model_real3 <- lm(y3 ~ x_real1 + x_real3))
(mod_naiv2 <- lm(y3 ~ x_het1 + x_het2, x = TRUE))
(mod_simex2 <- simex(mod_naiv2, SIMEXvariable = c("x_het1", "x_het2"),
measurement.error = cbind(sd_me_het1, sd_me_het2), asymptotic = FALSE))
plot(mod_simex2)
### two variables, one with homoscedastic, one with heteroscedastic measurement error
model_real2
(mod_naiv3 <- lm(y2 ~ x_measured1 + x_het2, x = TRUE))
(mod_simex3 <- simex(mod_naiv3, SIMEXvariable = c("x_measured1", "x_het2"),
measurement.error = cbind(sd_me1, sd_me_het2), asymptotic = FALSE))
### glm: two variables, one with homoscedastic, one with heteroscedastic measurement error
t <- x_real1 + 2*x_real2
g <- 1 / (1 + exp(-t))
u <- runif(100)
ybin <- as.numeric(u < g)
(logit_real <- glm(ybin ~ x_real1 + x_real2, family = binomial))
(logit_naiv <- glm(ybin ~ x_measured1 + x_het2, x = TRUE, family = binomial))
(logit_simex <- simex(logit_naiv, SIMEXvariable = c("x_measured1", "x_het2"),
measurement.error = cbind(sd_me1, sd_me_het2), asymptotic = FALSE))
summary(logit_simex)
print(logit_simex)
plot(logit_simex)
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.