Type: | Package |
Title: | Tools for Model Specification in the Latent Variable Framework |
Version: | 2.0.3 |
Date: | 2024-02-22 |
URL: | https://github.com/bozenne/lavaSearch2 |
BugReports: | https://github.com/bozenne/lavaSearch2/issues |
Description: | Tools for model specification in the latent variable framework (add-on to the 'lava' package). The package contains three main functionalities: Wald tests/F-tests with improved control of the type 1 error in small samples, adjustment for multiple comparisons when searching for local dependencies, and adjustment for multiple comparisons when doing inference for multiple latent variable models. |
License: | GPL-3 |
Encoding: | UTF-8 |
Depends: | R (≥ 2.10), ggplot2, lava (≥ 1.6.4) |
Imports: | abind, doParallel, MASS, Matrix, methods, multcomp, mvtnorm, nlme, parallel, Rcpp, reshape2, sandwich, stats, utils |
Suggests: | clubSandwich, data.table, foreach, emmeans, lme4, lmerTest, numDeriv, pbapply, pbkrtest, R.rsp, riskRegression, survival, testthat |
LinkingTo: | Rcpp, RcppArmadillo |
VignetteBuilder: | R.rsp |
NeedsCompilation: | yes |
RoxygenNote: | 7.2.3 |
Packaged: | 2024-02-22 11:18:47 UTC; hpl802 |
Author: | Brice Ozenne |
Maintainer: | Brice Ozenne <brice.mh.ozenne@gmail.com> |
Repository: | CRAN |
Date/Publication: | 2024-02-23 09:10:02 UTC |
Tools for Model Specification in the Latent Variable Framework
Description
The package contains three main functionalities:
-
compare2
,summary2
: Wald tests/robust Wald tests/F-tests/robust F-tests with improved control of the type 1 error in small samples. -
glht2
: adjustment for multiple comparisons when doing inference for multiple latent variable models. -
modelsearch2
: searching for local dependencies with adjustment for multiple comparisons.
It contains other useful functions such as:
-
calibrateType1
: simulation study of the type 1 error of Wald tests. -
createContrast
: user-friendly function generating a contrast matrix. -
getVarCov2
: reconstruct the conditional variance covariance matrix. -
iidJack
: extract the jackknife iid decomposition.
Details
The latent variable models (LVM) considered in this package can be written
as a measurement model:
Y_i = \nu + \eta_i \Lambda + X_i K + \epsilon_i
and a structural model:
\eta_i = \alpha + \eta_i B + X_i \Gamma + \zeta_i
where \Sigma
is the variance covariance matrix of the residuals \epsilon
,
and \Psi
is the variance covariance matrix of the residuals \zeta
.
The corresponding conditional mean is:
\mu_i(\theta) = E[Y_i|X_i] = \nu + (\alpha + X_i \Gamma) (1-B)^{-1} \Lambda + X_i K
\Omega(\theta) = Var[Y_i|X_i] = \Lambda^{t} (1-B)^{-t} \Psi (1-B)^{-1} \Lambda + \Sigma
The package aims to provides tool for testing linear hypotheses on the model coefficients
\nu
, \Lambda
, K
, \Sigma
,
\alpha
, B
, \Gamma
, \Psi
.
Searching for local dependency enable to test whether the proposed model is too simplistic and if so to identify which additional coefficients should be added to the model.
Limitations
'lavaSearch2' has been design for Gaussian latent variable models. This means that it may not work / give valid results:
in presence of censored or binary outcomes.
with stratified models (i.e. object of class
multigroup
).
Compute the First Derivative of the Expected Information Matrix
Description
Compute the first derivative of the expected information matrix.
Usage
.dInformation2(
dmu,
dOmega,
d2mu,
d2Omega,
OmegaM1,
missing.pattern,
unique.pattern,
name.pattern,
grid.3varD1,
grid.2meanD1.1varD1,
grid.2meanD2.1meanD1,
grid.2varD2.1varD1,
name.param,
leverage,
n.cluster,
weights
)
Details
calc_dinformation
will perform the computation individually when the
argument index.Omega
is not null.
Add a New Link Between Two Variables in a LVM
Description
Generic interface to add links to lvm
objects.
Usage
addLink(object, ...)
## S3 method for class 'lvm'
addLink(
object,
var1,
var2,
covariance,
all.vars = lava::vars(object),
warnings = FALSE,
...
)
## S3 method for class 'lvm.reduced'
addLink(object, ...)
Arguments
object |
a |
... |
[internal] only used by the generic method and from |
var1 |
[character or formula] the exogenous variable of the new link or a formula describing the link to be added to the lvm. |
var2 |
[character] the endogenous variable of the new link. Disregarded if the argument |
covariance |
[logical] is the link is bidirectional? Ignored if one of the variables non-stochastic (e.g. exogenous variables). |
all.vars |
[internal] a character vector containing all the variables of the |
warnings |
[logical] Should a warning be displayed when no link is added? |
Details
The argument all.vars
is useful for lvm.reduce
object where the command vars(object)
does not return all variables. The command vars(object, xlp = TRUE)
must be used instead.
Arguments var1
and var2
are passed to initVarlink
.
Examples
library(lava)
set.seed(10)
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u~x1+x2
latent(m) <- ~u
m2 <- m
addLink(m, x1 ~ y1, covariance = FALSE)
addLink(m, y1 ~ x1, covariance = FALSE)
coef(addLink(m, y1 ~ y2, covariance = TRUE))
addLink(m2, "x1", "y1", covariance = FALSE)
addLink(m2, "y1", "x1", covariance = FALSE)
newM <- addLink(m, "y1", "y2", covariance = TRUE)
coef(newM)
2D-display of the Domain Used to Compute the Integral
Description
2D-display of the domain used to compute the integral.
Usage
## S3 method for class 'intDensTri'
autoplot(object, coord.plot = c("x", "y1"), plot = TRUE, ...)
Arguments
object |
output of the function |
coord.plot |
[character vector] the x and y coordinates. Can be |
plot |
[logical] should the plot be displayed? |
... |
[internal] Only used by the generic method. |
Value
A ggplot
object.
See Also
Graphical Display of the Bias or Type 1 Error
Description
Graphical display of the bias or type 1 error
for the output of calibrateType1
.
Usage
## S3 method for class 'calibrateType1'
autoplot(
object,
type = "bias",
plot = TRUE,
color.threshold = "red",
type.bias = "absolute",
alpha = 0.05,
nrow.legend = NULL,
name2label = NULL,
color = NULL,
keep.method = NULL,
...
)
Arguments
object |
output of |
type |
[character] if type equals |
plot |
[logical] should the plot be displayed? |
color.threshold |
[character] the color for the line representing the expected value(s). |
type.bias |
[character] if type.bias equals |
alpha |
[numeric, 0-1] the significance threshold to consider.
Only relevant when type equals |
nrow.legend |
[integer, >0] the number of rows for the legend.
Only relevant when type equals |
name2label |
[named character vector] the label for the legend.
The vector should contain the method names (see details).
Only relevant when type equals |
color |
[character vector] a vector of colours to be used to color the lines.
Only relevant when type equals |
keep.method |
[character vector] the methods names for which the type 1 error should be displayed.
Only relevant when type equals |
... |
[internal] Only used by the generic method. |
Details
Method names:
-
p.Ztest
-
p.Satt
-
p.KR
-
p.robustZtest
-
p.robustSatt
-
p.robustKR
Value
An list containing:
plot: a ggplot object.
data: the dataset used to generate the ggplot object.
Display the Value of a Coefficient across the Steps.
Description
Display the value of a coefficient across the steps.
Usage
## S3 method for class 'modelsearch2'
autoplot(
object,
param,
ci = TRUE,
step = 0:nStep(object),
conf.level = 0.95,
plot = TRUE,
add.0 = TRUE,
...
)
Arguments
object |
a |
param |
[character vector] the name of the coefficient(s) to be displayed. |
ci |
[logical] should the confidence intervals of the coefficient(s) be displayed. |
step |
[integer >0] the steps at which the coefficient value should be displayed. |
conf.level |
[numeric, 0-1] confidence level of the interval. |
plot |
[logical] should the graph be displayed? |
add.0 |
[logical] should an horizontal line representing no effect be displayed? |
... |
[internal] only used by the generic method. |
Value
A list containing
-
plot
: a ggplot object. -
data
: the data used to generate the ggplot object.
Examples
## Not run:
mSim <- lvm(Y~G+X1+X2+X3+X4+X5)
addvar(mSim) <- ~Z1+Z2
set.seed(10)
df.data <- lava::sim(mSim, 1e2)
mBase <- lvm(Y~G)
addvar(mBase) <- ~X1+X2+X3+X4+X5+Z1+Z2
e.lvm <- estimate(mBase, data = df.data)
res <- modelsearch2(e.lvm, method.p.adjust = "holm", alpha = 0.05)
autoplot(res, param = "Y~G")
autoplot(res, param = c("Y","Y~G"))
## End(Not run)
Adjust the p.values Using the Quantiles of the Max Statistic
Description
Adjust the p.values using the quantiles of the max statistic.
Usage
calcDistMaxIntegral(
statistic,
iid,
df,
iid.previous = NULL,
quantile.previous = NULL,
quantile.compute = lava.options()$search.calc.quantile.int,
alpha,
cpus = 1,
cl = NULL,
trace
)
calcDistMaxBootstrap(
statistic,
iid,
iid.previous = NULL,
quantile.previous = NULL,
method,
alpha,
cpus = 1,
cl = NULL,
n.sim,
trace,
n.repmax = 100
)
Arguments
statistic |
[numeric vector] the observed Wald statistic. Each statistic correspond to a null hypothesis (i.e. a coefficient) that one wish to test. |
iid |
[matrix] zero-mean iid decomposition of the coefficient used to compute the statistic. |
df |
[numeric] the degree of freedom defining the multivariate Student's t distribution.
If |
iid.previous |
[matrix, EXPERIMENTAL] zero-mean iid decomposition of previously tested coefficient. |
quantile.previous |
[numeric, EXPERIMENTAL] rejection quantiles of the previously tested hypotheses. If not |
quantile.compute |
[logical] should the rejection quantile be computed? |
alpha |
[numeric 0-1] the significance cutoff for the p-values. When the p-value is below, the corresponding link will be retained. |
cpus |
[integer >0] the number of processors to use. If greater than 1, the computation of the p-value relative to each test is performed in parallel. |
cl |
[cluster] a parallel socket cluster generated by |
trace |
[logical] should the execution of the function be traced? |
method |
[character] the method used to compute the p-values. |
n.sim |
[integer >0] the number of bootstrap simulations used to compute each p-values. Disregarded when the p-values are computed using numerical integration. |
n.repmax |
[integer >0] the maximum number of rejection for each bootstrap sample before switching to a new bootstrap sample. Only relevant when conditioning on a previous test. Disregarded when the p-values are computed using numerical integration. |
Value
A list containing
p.adjust: the adjusted p-values.
z: the rejection threshold.
Sigma: the correlation matrix between the test statistic.
correctedLevel: the alpha level corrected for conditioning on previous tests.
Examples
library(mvtnorm)
set.seed(10)
n <- 100
p <- 4
link <- letters[1:p]
n.sim <- 1e3 # number of bootstrap simulations
#### test - not conditional ####
X.iid <- rmvnorm(n, mean = rep(0,p), sigma = diag(1,p))
colnames(X.iid) <- link
statistic <- setNames(1:p,link)
r1 <- calcDistMaxIntegral(statistic = statistic, iid = X.iid,
trace = FALSE, alpha = 0.05, df = 1e6)
r3 <- calcDistMaxBootstrap(statistic = statistic, iid = X.iid,
method = "residual",
trace = FALSE, alpha = 0.05, n.sim = n.sim)
r4 <- calcDistMaxBootstrap(statistic = statistic, iid = X.iid,
method = "wild",
trace = FALSE, alpha = 0.05, n.sim = n.sim)
rbind(integration = c(r1$p.adjust, quantile = r1$z),
bootResidual = c(r3$p.adjust, quantile = r3$z),
bootWild = c(r4$p.adjust, quantile = r4$z))
#### test - conditional ####
## Not run:
Z.iid <- rmvnorm(n, mean = rep(0,p+1), sigma = diag(1,p+1))
seqQuantile <- qmvnorm(p = 0.95, delta = rep(0,p+1), sigma = diag(1,p+1),
tail = "both.tails")$quantile
r1c <- calcDistMaxIntegral(statistic = statistic, iid = X.iid,
iid.previous = Z.iid, quantile.previous = seqQuantile,
trace = FALSE, alpha = 0.05, df = NULL)
r3c <- calcDistMaxBootstrap(statistic = statistic, iid = X.iid,
iid.previous = Z.iid, quantile.previous = seqQuantile, method = "residual",
trace = FALSE, alpha = 0.05, n.sim = n.sim)
r4c <- calcDistMaxBootstrap(statistic = statistic, iid = X.iid,
iid.previous = Z.iid, quantile.previous = seqQuantile, method = "wild",
trace = FALSE, alpha = 0.05, n.sim = n.sim)
rbind(integration = c(r1c$p.adjust, quantile = r1c$z),
bootResidual = c(r3c$p.adjust, quantile = r3c$z),
bootWild = c(r4c$p.adjust, quantile = r4c$z))
## End(Not run)
Compute the Type 1 Error After Selection [EXPERIMENTAL]
Description
Compute the type 1 error after selection [EXPERIMENTAL].
Usage
calcType1postSelection(
level,
mu,
Sigma,
quantile.previous,
distribution,
df,
n = 10,
correct = TRUE,
...
)
Arguments
level |
[numeric 0-1] expected coverage. |
mu |
[numeric vector] the expectation of the joint distribution of the test statistics |
Sigma |
[matrix] the variance-covariance of the joint distribution of the test statistics. |
quantile.previous |
[numeric] significance quantile used at the previous step. |
distribution |
[character] distribution of the test statistics.
Can be |
df |
[integer > 0] the degree of freedom of the joint Student's t distribution.
Only used when |
n |
[integer > 0] number of points for the numerical integration |
correct |
[logical] if true, correct the level to account for previous testings. |
... |
arguments passed to |
Details
The number of tests at the current step (i.e. after selection) is assumed to be one less than the number of tests at the previous step (i.e. before selection).
Arguments mu
and Sigma
must contain the moments for the vector of test statistics
before and after selection (in that order).
Value
[numeric] the type 1 error.
Author(s)
Brice Ozenne
Examples
library(mvtnorm)
n <- 350
#### only 2 tests
Sigma <- rbind(c(1,0,0),c(0,1,1),c(0,1,1))
z2 <- qmvnorm(0.95, mean = rep(0,2), sigma = Sigma[1:2,1:2], tail = "both.tails")$quantile
## no selection since strong effect
mu <- c(10,0,0)
calcType1postSelection(0.95, quantile.previous = z2, distribution = "gaussian",
mu = mu, Sigma = Sigma, correct = TRUE)
## strong selection
## Not run:
mu <- c(0,0,0)
levelC <- calcType1postSelection(0.95, quantile.previous = z2, distribution = "gaussian",
mu = mu, Sigma = Sigma)
print(levelC) # more liberal than without selection
calcType1postSelection(levelC, quantile.previous = z2, distribution = "gaussian",
mu = mu, Sigma = Sigma, correct = FALSE)
## End(Not run)
#### 3 tests
Sigma <- diag(1,5,5)
Sigma[4,2] <- 1
Sigma[2,4] <- 1
Sigma[5,3] <- 1
Sigma[3,5] <- 1
z2 <- qmvnorm(0.95, mean = mu[1:3], sigma = Sigma[1:3,1:3], tails = "both.tails")$quantile
## no selection since strong effect
## Not run:
mu <- c(10,0,0,0,0)
calcType1postSelection(0.95, quantile.previous = z2, distribution = "gaussian",
mu = mu, Sigma = Sigma, correct = TRUE)
## strong selection
mu <- c(0,0,0,0,0)
levelC <- calcType1postSelection(0.95, quantile.previous = z2,
mu = mu, Sigma = Sigma, distribution = "gaussian")
calcType1postSelection(levelC, quantile.previous = z2, distribution = "gaussian",
mu = mu, Sigma = Sigma, correct = FALSE)
## End(Not run)
Simulation Study Assessing Bias and Type 1 Error
Description
Perform a simulation study over one or several sample size to assess the bias of the estimate and the type 1 error of the Wald test and robust Wald test
Usage
calibrateType1(object, param, n.rep, ...)
## S3 method for class 'lvm'
calibrateType1(
object,
param,
n.rep,
n,
correction = TRUE,
warmup = NULL,
null = NULL,
F.test = FALSE,
cluster = NULL,
generative.object = NULL,
generative.coef = NULL,
true.coef = NULL,
n.true = 1e+06,
round.true = 2,
bootstrap = FALSE,
n.bootstrap = 1000,
checkType1 = FALSE,
checkType2 = FALSE,
dir.save = NULL,
label.file = NULL,
seed = NULL,
cpus = 1,
trace = 2,
...
)
## S3 method for class 'lvmfit'
calibrateType1(
object,
param,
n.rep,
correction = TRUE,
F.test = FALSE,
bootstrap = FALSE,
n.bootstrap = 1000,
seed = NULL,
trace = 2,
cpus = 1,
...
)
Arguments
object |
a |
param |
[character vector] names of the coefficient whose value will be tested. |
n.rep |
[integer, >0] number of simulations per sample size. |
... |
[internal] Only used by the generic method. |
n |
[integer vector, >0] sample size(s) considered in the simulation study. |
correction |
[logical] should the type 1 error after correction be computed? |
warmup |
[list of lvm] a list of |
null |
[numeric vector] vector of null hypotheses, one for each model coefficient. By default a vector of 0. |
F.test |
[logical] should a multivariate Wald test be perform testing simultaneously all the null hypotheses? |
cluster |
[integer vector] the grouping variable relative to which the observations are iid.
Will be passed to |
generative.object |
[lvm] object defining the statistical model generating the data. |
generative.coef |
[name numeric vector] values for the parameters of the generative model.
Can also be |
true.coef |
[name numeric vector] expected values for the parameters of the fitted model. |
n.true |
[integer, >0] sample size at which the estimated coefficients will be a reliable approximation of the true coefficients. |
round.true |
[integer, >0] the number of decimal places to be used for the true value of the coefficients. No rounding is done if |
bootstrap |
[logical] should bootstrap resampling be performed? |
n.bootstrap |
[integer, >0] the number of bootstrap sample to be used for each bootstrap. |
checkType1 |
[logical] returns an error if the coefficients associated to the null hypotheses do not equal 0. |
checkType2 |
[logical] returns an error if the coefficients associated to the null hypotheses equal 0. |
dir.save |
[character] path to the directory were the results should be exported.
Can also be |
label.file |
[character] element to include in the file name. |
seed |
[integer, >0] value that will be set before adjustment for multiple comparisons to ensure reproducible results.
Can also be |
cpus |
[integer >0] the number of processors to use. If greater than 1, the simulations are performed in parallel. |
trace |
[integer] should the execution of the function be trace. Can be 0, 1 or 2. |
Value
An object of class calibrateType1
.
Author(s)
Brice Ozenne
See Also
link{autoplot.calibrateType1}
for a graphical display of the bias or of the type 1 error.
Examples
## Not run:
#### simulate data ####
m.Sim <- lvm(c(Y1[mu1:sigma]~1*eta,
Y2[mu2:sigma]~1*eta,
Y3[mu3:sigma]~1*eta,
eta~beta1*Group+beta2*Gender))
latent(m.Sim) <- ~eta
categorical(m.Sim, labels = c("M","F")) <- ~Gender
d <- lava::sim(m.Sim, 1e2)
#### calibrate type 1 error on the estimated model ####
m <- lvm(Y1~eta,
Y2~eta,
Y3~eta,
eta~Group+Gender)
e <- lava::estimate(m, data = d)
res <- calibrateType1(e, param = "eta~Group", n.rep = 100)
res <- calibrateType1(e, param = c("eta~Group","Y1~eta"), F.test = TRUE, n.rep = 100)
res <- calibrateType1(e, param = "eta~Group", n.rep = 100, cpus = 4)
summary(res)
## End(Not run)
Check that Validity of the Dataset
Description
Check whether the dataset can be used to fit the lvm
object.
Usage
checkData(object, data, trace)
## S3 method for class 'lvm'
checkData(object, data, trace = TRUE)
Arguments
object |
a |
data |
[data.frame] the dataset used to obtain the object. |
trace |
[logical] when |
Value
Invisible TRUE
or FALSE
.
Examples
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u~x
latent(m) <- ~u
d <- lava::sim(m,1e2)
try(checkData(m, data = d)) # return an error
checkData(m, data = d[,-4])
try(checkData(m, data = d[,-(3:4)])) # return an error
Simplify a lvm object
Description
Remove variables with no link.
Usage
clean(x, ...)
## S3 method for class 'lvm'
clean(x, rm.exo = TRUE, rm.endo = TRUE, rm.latent = TRUE, ...)
Arguments
x |
|
... |
additional arguments to lower level functions |
rm.exo |
should exogenous variables with no links be removed from the object? |
rm.endo |
should endogenous variables with no links be removed from the object? |
rm.latent |
should latent variables with no links be removed from the object? |
Examples
m <- lvm()
m <- regression(m, x=paste0("x",1:5),y="y1")
m <- regression(m, x=paste0("x",1:5),y="y2")
covariance(m) <- y1~y2
cancel(m) <- y1 ~ x1
cancel(m) <- y2 ~ x1
clean(m)
m <- lvm(y1 ~ eta + x1, y2 ~ eta, y3 ~ eta + x2)
latent(m) <- ~eta
clean(m)
m
cancel(m) <- y1 ~ eta
cancel(m) <- y2 ~ eta
cancel(m) <- y3 ~ eta
clean(m)
Model Coefficients With Small Sample Correction
Description
Extract the coefficients from a latent variable model.
Similar to lava::compare
but with small sample correction.
Usage
coef2(object, as.lava, ...)
## S3 method for class 'lvmfit'
coef2(object, as.lava = TRUE, ssc = lava.options()$ssc, ...)
Arguments
object |
a |
as.lava |
[logical] if |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the model coefficients.
Value
A numeric vector named with the names of the coefficients.
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
set.seed(10)
dW <- sampleRepeated(10, format = "wide")
set.seed(10)
dL <- sampleRepeated(10, format = "long")
dL$time2 <- paste0("visit",dL$time)
#### latent variable models ####
e.lvm <- estimate(lvm(c(Y1,Y2,Y3) ~ 1*eta + X1, eta ~ Z1), data = dW)
coef(e.lvm)
coef2(e.lvm)
coef2(e.lvm, as.lava = FALSE)
Extract the Coefficient by Type
Description
Extract specific types of coefficient from a lvm
object:
covariance coefficient(s) (coefCov
),
extra parameter(s) (coefExtra
),
position in the list of models for each coefficient (coefIndexModel
),
intercept coefficient(s) (coefIntercept
),
coefficient(s) that are used as reference (coefRef
),
regression coefficient(s) (coefReg
),
variance coefficient(s) (coefVar
).
Usage
coefCov(object, value, keep.var, ...)
## S3 method for class 'lvm'
coefCov(object, value = FALSE, keep.var = FALSE, ...)
## S3 method for class 'lvmfit'
coefCov(object, value = FALSE, keep.var = FALSE, ...)
## S3 method for class 'multigroup'
coefCov(object, value = FALSE, keep.var = FALSE, ...)
coefExtra(object, value, ...)
## S3 method for class 'lvm'
coefExtra(object, value = FALSE, ...)
## S3 method for class 'lvmfit'
coefExtra(object, value = FALSE, ...)
## S3 method for class 'multigroup'
coefExtra(object, value = FALSE, ...)
coefIndexModel(object, ...)
## S3 method for class 'lvm'
coefIndexModel(object, ...)
## S3 method for class 'lvmfit'
coefIndexModel(object, ...)
## S3 method for class 'multigroup'
coefIndexModel(object, ...)
## S3 method for class 'multigroupfit'
coefIndexModel(object, ...)
coefIntercept(object, value, ...)
## S3 method for class 'lvm'
coefIntercept(object, value = FALSE, ...)
## S3 method for class 'lvmfit'
coefIntercept(object, value = FALSE, ...)
## S3 method for class 'multigroup'
coefIntercept(object, value = FALSE, ...)
coefRef(object, value, ...)
## S3 method for class 'lvmfit'
coefRef(object, value = FALSE, ...)
coefReg(object, value, ...)
## S3 method for class 'lvm'
coefReg(object, value = FALSE, ...)
## S3 method for class 'lvmfit'
coefReg(object, value = FALSE, ...)
## S3 method for class 'multigroup'
coefReg(object, value = FALSE, ...)
coefVar(object, value, ...)
## S3 method for class 'lvm'
coefVar(object, value = FALSE, ...)
## S3 method for class 'lvmfit'
coefVar(object, value = FALSE, ...)
## S3 method for class 'multigroup'
coefVar(object, value = FALSE, ...)
Arguments
object |
a lvm model or a fitted lvm model |
value |
should the name of the coefficient be returned? Else return the coefficients |
keep.var |
should the variance coefficients be returned? |
... |
arguments to be passed to |
Value
A vector containing the names of the positions of the coefficients.
Examples
#### regression ####
m <- lvm(Y~X1+X2)
e <- estimate(m, lava::sim(m, 1e2))
coefCov(m)
coefCov(m, value = TRUE)
coefCov(m, keep.var = TRUE)
coefCov(m, value = TRUE, keep.var = TRUE)
coefIndexModel(m)
coefIndexModel(e)
coefIntercept(m)
coefIntercept(m, value = TRUE)
coefReg(m)
coefReg(m, value = TRUE)
#### LVM ####
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u~x1+x2
latent(m) <- ~u
covariance(m) <- y1~y2
m.Sim <- m
categorical(m.Sim, labels = c("a","b","c")) <- ~x2
e <- estimate(m, lava::sim(m.Sim, 1e2))
coefCov(m)
coefCov(m, value = TRUE)
coefCov(m, keep.var = TRUE)
coefCov(m, value = TRUE, keep.var = TRUE)
coefExtra(m)
coefIndexModel(m)
coefIndexModel(e)
## additional categorical variable
categorical(m, labels = as.character(1:3)) <- "X1"
coefExtra(m)
coefExtra(m, value = TRUE)
## additional categorical variable
categorical(m, labels = as.character(1:3)) <- "x1"
coefIntercept(m)
coefIntercept(m, value = TRUE)
coefIntercept(e)
coefReg(e, value = TRUE)
#### multigroup ####
m <- lvm(Y~X1+X2)
eG <- estimate(list(m,m), list(lava::sim(m, 1e2), lava::sim(m, 1e2)))
coefIndexModel(eG)
Extract the Type of Each Coefficient
Description
Extract the type of each coefficient of a lvm
object.
Usage
coefType(object, as.lava, ...)
## S3 method for class 'lvm'
coefType(object, as.lava = TRUE, data = NULL, ...)
## S3 method for class 'lvmfit'
coefType(object, as.lava = TRUE, ...)
## S3 method for class 'multigroup'
coefType(object, as.lava = TRUE, ...)
Arguments
object |
a |
as.lava |
[logical] export the type of coefficients mimicking |
... |
arguments to be passed to |
data |
[data.frame, optional] the dataset. Help to identify the categorical variables. |
Details
A lvm can be written as a measurement model:
Y_i = \nu + \Lambda \eta_i + K X_i + \epsilon_i
and a structural model:
\eta_i = \alpha + B \eta_i + \Gamma X_i + \zeta_i
where \Psi
is the variance covariance matrix of the residuals \zeta
and \Sigma
is the variance covariance matrix of the residuals \epsilon
.
coefType
either returns the Latin/Greek letter corresponding to the coefficients
or it groups them:
intercept:
\nu
and\alpha
.regression:
\Lambda
,K
,B
, and\Gamma
.covariance: extra-diagonal terms of
\Sigma
and\Psi
.variance: diagonal of
\Sigma
and\Psi
.
A link denotes a relationship between two variables. The coefficient are used to represent the strength of the association between two variable, i.e. the strength of a link. A coefficient may corresponds to the strength of one or several link.
Value
coefType
returns a data.frame
when as.lava=FALSE
:
name: name of the link
Y: outcome variable
X: regression variable in the design matrix (could be a transformation of the original variables, e.g. dichotomization).
data: original variable
type: type of link
value: if TRUE, the value of the link is set and not estimated.
marginal: if TRUE, the value of the link does not impact the estimation.
detail: a more detailed description of the type of link (see the details section)
lava: name of the coefficient in lava
When as.lava=TRUE
, coefType
returns a named vector containing the type of each coefficient.
Examples
#### regression ####
m <- lvm(Y~X1+X2)
e <- estimate(m, lava::sim(m, 1e2))
coefType(m)
coefType(e)
#### LVM ####
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u~x1+x2
latent(m) <- ~u
covariance(m) <- y1~y2
m.Sim <- m
categorical(m.Sim, labels = c("a","b","c")) <- ~x2
e <- estimate(m, lava::sim(m.Sim, 1e2))
coefType(m)
coefType(e)
## additional categorical variables
categorical(m, labels = as.character(1:3)) <- "X1"
coefType(m, as.lava = FALSE)
#### LVM with constrains ####
m <- lvm(c(Y1~0+1*eta1,Y2~0+1*eta1,Y3~0+1*eta1,
Z1~0+1*eta2,Z2~0+1*eta2,Z3~0+1*eta2))
latent(m) <- ~eta1 + eta2
e <- estimate(m, lava::sim(m,1e2))
coefType(m)
coefType(e)
#### multigroup ####
m <- lvm(Y~X1+X2)
eG <- estimate(list(m,m), list(lava::sim(m, 1e2), lava::sim(m, 1e2)))
coefType(eG)
Form all Unique Combinations Between two Vectors
Description
Form all unique combinations between two vectors (removing symmetric combinations).
Usage
.combination(..., levels = FALSE)
Arguments
... |
[vectors] elements to be combined. |
levels |
[logical] should a label for each combination be output as an attribute named levels. |
Value
A matrix, each row being a different combination.
Examples
.combination <- lavaSearch2:::.combination
.combination(1,1)
.combination(1:2,1:2)
.combination(c(1:2,1:2),1:2)
.combination(alpha = 1:2, beta = 3:4)
.combination(alpha = 1:2, beta = 3:4, gamma = 1:4)
.combination(alpha = 1:3, beta = 1:3, gamma = 1:3)
Combine formula
Description
Combine formula by outcome
Usage
combineFormula(ls.formula, as.formula = TRUE, as.unique = FALSE)
Arguments
ls.formula |
a list of formula |
as.formula |
should a list of formula be returned. Otherwise it will be a list of characters. |
as.unique |
should regressors appears at most once in the formula |
Examples
combineFormula(list(Y~X1,Y~X3+X5,Y1~X2))
lava.options(symbols = c("~",","))
combineFormula(list("Y~X1","Y~X3+X5","Y1~X2"))
lava.options(symbols = c("<-","<->"))
combineFormula(list("Y<-X1","Y<-X3+X5","Y1<-X2"))
combineFormula(list(Y~X1,Y~X3+X1,Y1~X2))
combineFormula(list(Y~X1,Y~X3+X1,Y1~X2), as.formula = FALSE)
combineFormula(list(Y~X1,Y~X3+X1,Y1~X2), as.unique = TRUE)
lava.options(symbols = c("~","~~"))
combineFormula(list("Y~X1","Y~X3","Y1~X2"))
Test Linear Hypotheses With Small Sample Correction
Description
Test Linear Hypotheses using Wald statistics in a latent variable model.
Similar to lava::compare
but with small sample correction.
Usage
compare2(
object,
linfct,
rhs,
robust,
cluster,
as.lava,
F.test,
conf.level,
...
)
## S3 method for class 'lvmfit'
compare2(
object,
linfct = NULL,
rhs = NULL,
robust = FALSE,
cluster = NULL,
as.lava = TRUE,
F.test = TRUE,
conf.level = 0.95,
ssc = lava.options()$ssc,
df = lava.options()$df,
...
)
## S3 method for class 'lvmfit2'
compare2(
object,
linfct = NULL,
rhs = NULL,
robust = FALSE,
cluster = NULL,
as.lava = TRUE,
F.test = TRUE,
conf.level = 0.95,
...
)
## S3 method for class 'lvmfit2'
compare(
object,
linfct = NULL,
rhs = NULL,
robust = FALSE,
cluster = NULL,
as.lava = TRUE,
F.test = TRUE,
conf.level = 0.95,
...
)
Arguments
object |
a |
linfct |
[matrix or vector of character] the linear hypotheses to be tested. Same as the argument |
rhs |
[vector] the right hand side of the linear hypotheses to be tested. |
robust |
[logical] should the robust standard errors be used instead of the model based standard errors? |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
as.lava |
[logical] should the output be similar to the one return by |
F.test |
[logical] should a joint test be performed? |
conf.level |
[numeric 0-1] level of the confidence intervals. |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
df |
[character] method used to estimate the degree of freedoms of the Wald statistic: Satterthwaite |
Details
The linfct
argument and rhs
specify the set of linear hypotheses to be tested. They can be written:
linfct * \theta = rhs
where \theta
is the vector of the model coefficients.
The par
argument must contain expression(s) involving the model coefficients.
For example "beta = 0"
or c("-5*beta + alpha = 3","-alpha")
are valid expressions if alpha and beta belong to the set of model coefficients.
A contrast matrix and the right hand side will be generated inside the function.
When directly specified, the contrast matrix must contain as many columns as there are coefficients in the model (mean and variance coefficients).
Each hypothesis correspond to a row in the contrast matrix.
The rhs vector should contain as many elements as there are row in the contrast matrix.
Value
If as.lava=TRUE
an object of class htest
.
Otherwise a data.frame
object.
See Also
createContrast
to create contrast matrices.
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
set.seed(10)
mSim <- lvm(Y~0.1*X1+0.2*X2)
categorical(mSim, labels = c("a","b","c")) <- ~X1
transform(mSim, Id~Y) <- function(x){1:NROW(x)}
df.data <- lava::sim(mSim, 1e2)
#### with lvm ####
m <- lvm(Y~X1+X2)
e.lvm <- estimate(m, df.data)
compare2(e.lvm, linfct = c("Y~X1b","Y~X1c","Y~X2"))
compare2(e.lvm, linfct = c("Y~X1b","Y~X1c","Y~X2"), robust = TRUE)
Confidence Intervals With Small Sample Correction
Description
Extract confidence intervals of the coefficients from a latent variable model.
Similar to lava::confint
but with small sample correction.
Extract estimate, standard error, confidence intervals and p-values associated to each coefficient of a latent variable model.
Similar to lava::confint
but with small sample correction.
Usage
confint2(object, robust, cluster, transform, as.lava, conf.level, ...)
## S3 method for class 'lvmfit'
confint2(
object,
robust = FALSE,
cluster = NULL,
transform = NULL,
as.lava = TRUE,
conf.level = 0.95,
ssc = lava.options()$ssc,
df = lava.options()$df,
...
)
model.tables2(object, robust, cluster, transform, as.lava, conf.level, ...)
Arguments
object |
a |
robust |
[logical] should robust standard errors be used instead of the model based standard errors? Should be |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
transform |
[function] transformation to be applied. |
as.lava |
[logical] when |
conf.level |
[numeric, 0-1] level of the confidence intervals. |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
df |
[character] method used to estimate the degree of freedoms of the Wald statistic: Satterthwaite |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the confidence intervals.
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the confidence intervals.
Value
A data.frame with a row per coefficient.
A data.frame with a row per coefficient.
Examples
#### simulate data ####
set.seed(10)
dW <- sampleRepeated(10, format = "wide")
set.seed(10)
dL <- sampleRepeated(10, format = "long")
dL$time2 <- paste0("visit",dL$time)
#### latent variable models ####
e.lvm <- estimate(lvm(c(Y1,Y2,Y3) ~ 1*eta + X1, eta ~ Z1), data = dW)
confint(e.lvm)
confint2(e.lvm)
confint2(e.lvm, as.lava = FALSE)
Create Rownames for a Contrast Matrix
Description
Create rownames for a contrast matrix using the coefficients and the names of the coefficients. The rownames will be [value * name] == null, e.g. [beta + 4*alpha] = 0.
Usage
.contrast2name(contrast, null = NULL, sep = c("[", "]"))
Arguments
contrast |
[matrix] a contrast matrix defining the left hand side of the linear hypotheses to be tested. |
null |
[vector, optional] the right hand side of the linear hypotheses to be tested. |
sep |
[character of length 2, optional] character used in rownames to wrap the left hand side of the equation. |
Details
When argument NULL
is null then the rownames will not be put into brackets and the right hand side will not be added to the name.
Value
a character vector.
formula character conversion
Description
Conversion of formula into character string or vice versa
Usage
formula2character(f, type = "formula")
Arguments
f |
a formula. |
type |
should the normal formula operator be used ( |
Examples
formula2character(Y1~X1+X2)
formula2character(Y1~X1+X2, type = "symbols")
Create Contrast matrix
Description
Returns a contrast matrix corresponding an object. The contrast matrix will contains the hypotheses in rows and the model coefficients in columns.
Usage
createContrast(object, ...)
## S3 method for class 'character'
createContrast(object, ...)
## S3 method for class 'lvmfit'
createContrast(object, linfct, ...)
## S3 method for class 'lvmfit2'
createContrast(object, linfct, ...)
## S3 method for class 'list'
createContrast(object, linfct = NULL, ...)
## S3 method for class 'mmm'
createContrast(object, linfct = NULL, ...)
Arguments
object |
a |
... |
Argument to be passed to
|
linfct |
[vector of characters] expression defining the linear hypotheses to be tested.
Can also be a regular expression (of length 1) that is used to identify the coefficients to be tested using |
Details
One can initialize an empty contrast matrix setting the argumentlinfct
to character(0)
.
Value
A list containing
contrast [matrix] a contrast matrix corresponding to the left hand side of the linear hypotheses.
null [vector] the right hand side of the linear hypotheses.
Q [integer] the rank of the contrast matrix.
ls.contrast [list, optional] the contrast matrix corresponding to each submodel. Only present when the
argument
object is a list of models.
Examples
## Simulate data
mSim <- lvm(X ~ Age + Treatment,
Y ~ Gender + Treatment,
c(Z1,Z2,Z3) ~ eta, eta ~ treatment,
Age[40:5]~1)
latent(mSim) <- ~eta
categorical(mSim, labels = c("placebo","SSRI")) <- ~Treatment
categorical(mSim, labels = c("male","female")) <- ~Gender
n <- 1e2
set.seed(10)
df.data <- lava::sim(mSim,n)
## Estimate separate models
lmX <- lava::estimate(lvm(X ~ -1 + Age + Treatment), data = df.data)
lmY <- lava::estimate(lvm(Y ~ -1 + Gender + Treatment), data = df.data)
lvmZ <- lava::estimate(lvm(c(Z1,Z2,Z3) ~ -1 + 1*eta, eta ~ -1 + Treatment),
data = df.data)
## Contrast matrix for a given model
createContrast(lmX, linfct = "X~Age")
createContrast(lmX, linfct = c("X~Age=0","X~Age+5*X~TreatmentSSRI=0"))
createContrast(lmX, linfct = c("X~Age=0","X~Age+5*X~TreatmentSSRI=0"), sep = NULL)
createContrast(lmX, linfct = character(0))
## Contrast matrix for the join model
ls.lvm <- list(X = lmX, Y = lmY, Z = lvmZ)
createContrast(ls.lvm, linfct = "TreatmentSSRI=0")
createContrast(ls.lvm, linfct = "TreatmentSSRI=0", rowname.rhs = FALSE)
createContrast(ls.lvm, linfct = character(0))
## Contrast for multigroup models
m <- lava::lvm(Y~Age+Treatment)
e <- lava::estimate(list(m,m), data = split(df.data, df.data$Gender))
print(coef(e))
createContrast(e, linfct = "Y~TreatmentSSRI@1 - Y~TreatmentSSRI@2 = 0")
createContrast(e, linfct = "Y~TreatmentSSRI@2 - Y~TreatmentSSRI@1 = 0")
Create a Mesh for the Integration
Description
Create a mesh for the integration.
Usage
createGrid(n, xmin, xmax, d.y, d.z, zmax, fine, double)
Arguments
n |
[integer >0] the number of points for the mesh in the x direction. |
xmin |
[numeric] the minimal x value. |
xmax |
[numeric] the maximal x value. |
d.y |
[integer >0] the number of dimensions for the triangle. |
d.z |
[integer >0] the number of dimensions. |
zmax |
[numeric >0] the maximal z value (in absolute value). |
fine |
[logical] should the mesh be displayed? |
double |
[logical] should the grid be just outside the region of interest? Otherwise it will be just inside. |
Details
This create a mesh for integrating over a triangular surface using rectangles. The domain is define by constrains on three types of variables:
the x variable: [unidimensional] that varies freely between [-xmax,-xmin] U [xmin,xmax].
the y variables: [dimension d.y] constrained to be lower in absolute value than the x variable. In 2D this corresponds to 2 triangles, and in higher dimension to cones/hypercones.
the z variables: [dimension d.z] constrained to vary between [-zmax,zmax].
The intersection of these three conditions define the domain.
The mesh is obtained slicing the triangles using rectangles.
Value
A list containing:
grid: a
data.frame
where each line corresponds to a point.seqNames.min: a
character
giving the name of the x,y variables (min).seqNames.max: a
character
giving the name of the x,y variables (max).
Examples
createGrid <- lavaSearch2:::createGrid
## no z
gridInt_2d <- createGrid(5, d.y = 1, xmin = 0, xmax = 4,
d.z = 0, fine = FALSE, double = FALSE)
gridExt_2d <- createGrid(5, d.y = 1, xmin = 0, xmax = 4,
d.z = 0, fine = FALSE, double = TRUE)
gridInt_4d <- createGrid(5, d.y = 3, xmin = 0, xmax = 4,
d.z = 0, fine = FALSE, double = FALSE)
gridExt_4d <- createGrid(5, d.y = 3, xmin = 0, xmax = 4,
d.z = 0, fine = FALSE, double = TRUE)
gridInt_2d <- createGrid(5, d.y = 1, xmin = 0, xmax = 4,
d.z = 0, fine = TRUE, double = FALSE)
## z
gridIntZ1_2d <- createGrid(5, d.y = 1, xmin = 0, xmax = 4,
d.z = 1, zmax = 2, fine = FALSE, double = FALSE)
gridExtZ1_2d <- createGrid(5, d.y = 1, xmin = 0, xmax = 4,
d.z = 1, zmax = 2, fine = FALSE, double = TRUE)
gridIntZ2_4d <- createGrid(5, d.y = 3, xmin = 0, xmax = 4,
d.z = 2, zmax = 2, fine = FALSE, double = FALSE)
gridExtZ2_4d <- createGrid(5, d.y = 3, xmin = 0, xmax = 4,
d.z = 2, zmax = 2, fine = FALSE, double = TRUE)
Identify Categorical Links in LVM
Description
Identify categorical links in latent variable models.
Usage
defineCategoricalLink(object, link, data)
## S3 method for class 'lvm'
defineCategoricalLink(object, link = NULL, data = NULL)
Arguments
object |
a |
link |
[character] the links to be analyzed. If |
data |
[data.frame] the dataset that will be used to fit the model. If |
Value
a data.frame
with a description of each link in rows.
The column factitious identify the links that will be replaced with other links
(e.g. "Y1~X1" becomes "Y1~X1b" and "Y1~X1c").
Examples
## Not run:
defineCategoricalLink <- lavaSearch2:::defineCategoricalLink
defineCategoricalLink.lvm <- lavaSearch2:::defineCategoricalLink.lvm
## linear model
m <- lvm(Y1~X1+X2,Y2~X1+X3)
categorical(m, K = 3) <- "X1"
try(defineCategoricalLink(m)) # error
categorical(m, K = 3, labels = 1:3) <- "X1"
defineCategoricalLink(m)
defineCategoricalLink(m, "Y~X1")
defineCategoricalLink(m, "X1:0|1")
defineCategoricalLink(m, "X1:1|2")
defineCategoricalLink(m, c("X1:0|1", "X1:1|2"))
defineCategoricalLink(m, c("Y~X1","Y~X2"))
defineCategoricalLink(m, c("Y~X2","Y~X1"))
## latent variable model
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u~x1+x2
latent(m) <- ~u
covariance(m) <- y1~y2
categorical(m, labels = as.character(1:3)) <- "X1"
defineCategoricalLink(m)
## End(Not run)
Degree of Freedom for the Chi-Square Test
Description
Computation of the degrees of freedom of the chi-squared distribution relative to the model-based variance
Usage
dfSigma(contrast, score, vcov, rvcov, dVcov, dRvcov, keep.param, type)
Arguments
contrast |
[numeric vector] the linear combination of parameters to test |
score |
[numeric matrix] the individual score for each parameter. |
vcov |
[numeric matrix] the model-based variance-covariance matrix of the parameters. |
rvcov |
[numeric matrix] the robust variance-covariance matrix of the parameters. |
dVcov |
[numeric array] the first derivative of the model-based variance-covariance matrix of the parameters. |
dRvcov |
[numeric array] the first derivative of the robust variance-covariance matrix of the parameters. |
keep.param |
[character vector] the name of the parameters with non-zero first derivative of their variance parameter. |
type |
[integer] 1 corresponds to the Satterthwaite approximation of the the degrees of freedom applied to the model-based variance, 2 to the Satterthwaite approximation of the the degrees of freedom applied to the robust variance, 3 to the approximation described in (Pan, 2002) section 2 and 3.1. |
References
Wei Pan and Melanie M. Wall, Small-sample adjustments in using the sandwich variance estiamtor in generalized estimating equations. Statistics in medicine (2002) 21:1429-1441.
Effects Through Pathways With Small Sample Correction
Description
Test whether a path in the latent variable model correspond to a null effect.
Similar to lava::effects
but with small sample correction (if any).
So far it only work for a single path related two variable composed of one or two edges.
Usage
effects2(object, linfct, robust, cluster, conf.level, ...)
## S3 method for class 'lvmfit'
effects2(
object,
linfct,
robust = FALSE,
cluster = NULL,
conf.level = 0.95,
to = NULL,
from = NULL,
df = lava.options()$df,
ssc = lava.options()$ssc,
...
)
## S3 method for class 'lvmfit2'
effects2(
object,
linfct,
robust = FALSE,
cluster = NULL,
conf.level = 0.95,
to = NULL,
from = NULL,
...
)
## S3 method for class 'lvmfit2'
effects(
object,
linfct,
robust = FALSE,
cluster = NULL,
conf.level = 0.95,
to = NULL,
from = NULL,
...
)
Arguments
object |
a |
linfct |
[character vector] The path for which the effect should be assessed (e.g. |
robust |
[logical] should robust standard errors be used instead of the model based standard errors? Should be |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
conf.level |
[numeric, 0-1] level of the confidence intervals. |
... |
additional argument passed to |
from , to |
alternative to argument |
df |
[character] method used to estimate the degree of freedoms of the Wald statistic: Satterthwaite |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the confidence intervals.
Value
A data.frame with a row per path.
Satterthwaite Correction and Small Sample Correction
Description
Correct the bias of the ML estimate of the variance and compute the first derivative of the information matrix.
Compute bias corrected residuals variance covariance matrix and information matrix. Also provides the leverage values and corrected sample size when adjust.n is set to TRUE.
Usage
estimate2(
object,
param,
data,
ssc,
df,
derivative,
hessian,
dVcov.robust,
iter.max,
tol.max,
trace,
...
)
## S3 method for class 'lvm'
estimate2(
object,
param = NULL,
data = NULL,
ssc = lava.options()$ssc,
df = lava.options()$df,
derivative = "analytic",
hessian = FALSE,
dVcov.robust = FALSE,
iter.max = 100,
tol.max = 1e-06,
trace = 0,
...
)
## S3 method for class 'lvmfit'
estimate2(
object,
param = NULL,
data = NULL,
ssc = lava.options()$ssc,
df = lava.options()$df,
derivative = "analytic",
hessian = FALSE,
dVcov.robust = FALSE,
iter.max = 100,
tol.max = 1e-06,
trace = 0,
...
)
## S3 method for class 'list'
estimate2(object, ...)
## S3 method for class 'mmm'
estimate2(object, ...)
.sscResiduals(object, ssc, algorithm = "2")
Arguments
object |
a |
param |
[numeric vector, optional] the values of the parameters at which to perform the correction. |
data |
[data.frame, optional] the dataset relative to which the correction should be performed. |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
df |
[character] method used to estimate the degree of freedoms of the Wald statistic: Satterthwaite |
derivative |
[character] should the first derivative of the information matrix be computed using a formula ( |
hessian |
[logical] should the hessian be stored? Can be |
dVcov.robust |
[logical] should the first derivative of robust variance-covariance matrix be stored? |
iter.max |
[integer >0] the maximum number of iterations used to estimate the bias correction. |
tol.max |
[numeric >0] the largest acceptable absolute difference between two succesive estimates of the bias correction. |
trace |
[logical] should the execution of the function be traced. |
... |
arguments passed to |
Details
The argument value
is equivalent to the argument bias.correct
of the function summary2
.
Examples
#### simulate data ####
set.seed(10)
dW <- sampleRepeated(10, format = "wide")
#### latent variable model ####
m.lvm <- lvm(Y1~X1+X2+Z1)
e2.lvm <- estimate2(m.lvm, data = dW)
summary2(e2.lvm)
Find Object in the Parent Environments
Description
Search an object in the parent environments. For internal use.
Usage
evalInParentEnv(name)
Arguments
name |
[character] the name of the object to get. |
Extract Data From a Latent Variable Model
Description
Extract data from a latent variable model.
Usage
extractData(object, design.matrix, as.data.frame, envir, rm.na)
## S3 method for class 'lvmfit'
extractData(
object,
design.matrix = FALSE,
as.data.frame = TRUE,
envir = environment(),
rm.na = TRUE
)
Arguments
object |
the fitted model. |
design.matrix |
[logical] should the data be extracted after transformation (e.g. conversion of categorical variables to dummy variables)? Otherwise the original data will be returned. |
as.data.frame |
[logical] should the output be converted into a |
envir |
[environment] the environment from which to search the data. |
rm.na |
[logical] should the lines containing missing values in the dataset be removed? |
Value
a dataset.
Examples
#### simulate data ####
set.seed(10)
n <- 101
Y1 <- rnorm(n, mean = 0)
Y2 <- rnorm(n, mean = 0.3)
Id <- findInterval(runif(n), seq(0.1,1,0.1))
data.df <- rbind(data.frame(Y=Y1,G="1",Id = Id),
data.frame(Y=Y2,G="2",Id = Id)
)
#### latent variable model ####
library(lava)
e.lvm <- estimate(lvm(Y ~ G), data = data.df)
extractData(e.lvm)
extractData(e.lvm, design.matrix = TRUE)
Find all New Links Between Variables
Description
Find all new links between variables (adapted from lava::modelsearch).
Usage
findNewLink(object, ...)
## S3 method for class 'lvm'
findNewLink(
object,
data = NULL,
type = "both",
exclude.var = NULL,
rm.latent_latent = FALSE,
rm.endo_endo = FALSE,
rm.latent_endo = FALSE,
output = "names",
...
)
Arguments
object |
a |
... |
[internal] only used by the generic method. |
data |
[optional] a dataset used to identify the categorical variables when not specified in the |
type |
[character vector] the type of links to be considered: |
exclude.var |
[character vector] all links related to these variables will be ignore. |
rm.latent_latent |
[logical] should the links relating two latent variables be ignored? |
rm.endo_endo |
[logical] should the links relating two endogenous variables be ignored? |
rm.latent_endo |
[logical] should the links relating one endogenous variable and one latent variable be ignored? |
output |
[character] Specify |
Value
A list containing:
M.links: a matrix with two columns indicating (by name or position) the exogenous and endogenous variable corresponding to each link.
links: the name of the additional possible links
directional: a logical vector indicating for each link whether the link is unidirectional (
TRUE
, i.e. regression link) or bidirectional (FALSE
, i.e. covariance link).
Examples
library(lava)
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
categorical(m,labels=c("M","F","MF")) <- ~X1
findNewLink(m, rm.endo = FALSE)
findNewLink(m, rm.endo = TRUE)
findNewLink(m, exclude.var = "X1")
regression(m) <- u~x1+x2
latent(m) <- ~u
findNewLink(m, rm.endo = FALSE)
findNewLink(m, rm.endo = TRUE)
findNewLink(m, rm.endo = TRUE, output = "index")
findNewLink(m, type = "covariance")
findNewLink(m, type = "regression")
Estimate LVM With Weights
Description
Estimate LVM with weights.
Usage
gaussian_weight.estimate.hook(x, data, estimator, ...)
gaussian_weight_method.lvm
gaussian_weight_logLik.lvm(object, type = "cond", p, data, weights, ...)
gaussian_weight_objective.lvm(x, ...)
gaussian_weight_score.lvm(
x,
data,
p,
S,
n,
mu = NULL,
weights = NULL,
debug = FALSE,
reindex = FALSE,
mean = TRUE,
constrain = TRUE,
indiv = FALSE,
...
)
gaussian_weight_gradient.lvm(...)
gaussian_weight_hessian.lvm(x, p, n, weights = NULL, ...)
Arguments
x , object |
A latent variable model |
data |
dataset |
estimator |
name of the estimator to be used |
... |
passed to lower level functions. |
type |
must be "cond" |
p |
parameter value |
weights |
weight associated to each iid replicate. |
S |
empirical variance-covariance matrix between variable |
n |
number of iid replicates |
mu |
empirical mean |
debug , reindex , mean , constrain , indiv |
additional arguments not used |
Format
An object of class character
of length 1.
Examples
#### linear regression with weights ####
## data
df <- data.frame(Y = c(1,2,2,1,2),
X = c(1,1,2,2,2),
missing = c(0,0,0,0,1),
weights = c(1,1,2,1,NA))
## using lm
e.lm.GS <- lm(Y~X, data = df)
e.lm.test <- lm(Y~X, data = df[df$missing==0,], weights = df[df$missing==0,"weights"])
## using lvm
m <- lvm(Y~X)
e.GS <- estimate(m, df)
## e.lava.test <- estimate(m, df[df$missing==0,], weights = df[df$missing==0,"weights"])
## warnings!!
e.test <- estimate(m, data = df[df$missing==0,],
weights = df[df$missing==0,"weights"],
estimator = "gaussian_weight")
Identify the Endogenous Variables
Description
Identify the endogenous variables, i.e., returns a vector with length the number of observations, whose values are the index of the repetitions.
Usage
.getIndexOmega(object, data, ...)
## S3 method for class 'lvm'
.getIndexOmega(object, data, ...)
## S3 method for class 'lvmfit'
.getIndexOmega(object, data, ...)
Arguments
object |
a |
data |
dataset. |
... |
[internal] Only used by the generic method. |
Examples
## Not run:
#### simulate data ####
set.seed(10)
dW <- sampleRepeated(10, format = "wide")
set.seed(10)
dL <- sampleRepeated(10, format = "long")
dL$time2 <- paste0("visit",dL$time)
#### lvm model ####
e.lvm <- estimate(lvm(c(Y1,Y2,Y3) ~ 1*eta + X1, eta ~ Z1), data = dW)
## lavaSearch2:::.getIndexOmega(e.lvm, data = dW)
## End(Not run)
Extract the Links that Have Been Found by the modelsearch2.
Description
Extract the links that have been found relevant by modelsearch2.
Usage
getNewLink(object, step)
## S3 method for class 'modelsearch2'
getNewLink(object, step = 1:nStep(object))
Arguments
object |
a |
step |
[logical] which test should be extracted? |
Value
A character vector.
Examples
## Not run:
mSim <- lvm(Y~G+X1+X2)
addvar(mSim) <- ~Z1+Z2+Z3+Z4+Z5+Z6
set.seed(10)
df.data <- lava::sim(mSim, 1e2)
mBase <- lvm(Y~G)
addvar(mBase) <- ~X1+X2+Z1+Z2+Z3+Z4+Z5+Z6
e.lvm <- estimate(mBase, data = df.data)
res <- modelsearch2(e.lvm, method.p.adjust = "holm")
getNewLink(res)
## End(Not run)
Extract the Model that Has Been Retains by the modelsearch2.
Description
Extract the model that has been retained by modelsearch2.
Usage
getNewModel(object, step)
## S3 method for class 'modelsearch2'
getNewModel(object, step = nStep(object))
Arguments
object |
a |
step |
[integer >=0] the step at which the model should be extracted. 0 returns the initial model, i.e. before adding any links. |
Value
A lvmfit
object.
Examples
## Not run:
mSim <- lvm(Y~G+X1+X2)
addvar(mSim) <- ~Z1+Z2+Z3+Z4+Z5+Z6
set.seed(10)
df.data <- lava::sim(mSim, 1e2)
mBase <- lvm(Y~G)
addvar(mBase) <- ~X1+X2+Z1+Z2+Z3+Z4+Z5+Z6
e.lvm <- estimate(mBase, data = df.data)
res <- modelsearch2(e.lvm, method.p.adjust = "holm")
getNewModel(res)
## End(Not run)
Extract one Step From the Sequential Procedure
Description
Extract one step from the sequential procedure.
Usage
getStep(object, step, slot)
## S3 method for class 'modelsearch2'
getStep(object, step = nStep(object), slot = NULL)
Arguments
object |
a |
step |
[integer >0] which test should be extracted? |
slot |
[character] the element from the modelsearch2 object that should be extracted. |
Examples
## Not run:
mSim <- lvm(Y~G+X1+X2)
addvar(mSim) <- ~Z1+Z2+Z3+Z4+Z5+Z6
df.data <- lava::sim(mSim, 1e2)
mBase <- lvm(Y~G)
addvar(mBase) <- ~X1+X2+Z1+Z2+Z3+Z4+Z5+Z6
e.lvm <- estimate(mBase, data = df.data)
res <- modelsearch2(e.lvm, method.p.adjust = "holm")
getStep(res)
getStep(res, slot = "sequenceTest")
getStep(res, step = 1)
## End(Not run)
Residual Variance-Covariance Matrix With Small Sample Correction.
Description
Reconstruct the residual variance-covariance matrix from a latent variable model.
It is similar to nlme::getVarCov
but with small sample correction.
Usage
getVarCov2(object, ...)
## S3 method for class 'lvmfit'
getVarCov2(object, ssc = lava.options()$ssc, ...)
## S3 method for class 'lvmfit2'
getVarCov2(object, ...)
Arguments
object |
a |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the residuals.
Value
A matrix with as many rows and column as the number of endogenous variables
Examples
#### simulate data ####
set.seed(10)
n <- 101
Y1 <- rnorm(n, mean = 0)
Y2 <- rnorm(n, mean = 0.3)
Id <- findInterval(runif(n), seq(0.1,1,0.1))
data.df <- rbind(data.frame(Y=Y1,G="1",Id = Id),
data.frame(Y=Y2,G="2",Id = Id)
)
#### latent variable models ####
library(lava)
e.lvm <- estimate(lvm(Y ~ G), data = data.df)
getVarCov2(e.lvm)
General Linear Hypothesis Testing With Small Sample Correction
Description
Test linear hypotheses on coefficients from a latent variable models with small sample corrections.
Usage
glht2(object, ...)
## S3 method for class 'lvmfit'
glht2(
object,
linfct,
rhs = NULL,
robust = FALSE,
cluster = NULL,
ssc = lava.options()$ssc,
df = lava.options()$df,
...
)
## S3 method for class 'lvmfit2'
glht2(object, linfct, rhs = NULL, robust = FALSE, cluster = NULL, ...)
## S3 method for class 'mmm'
glht2(object, linfct, rhs = 0, robust = FALSE, cluster = NULL, ...)
## S3 method for class 'lvmfit2'
glht(model, linfct, rhs = NULL, robust = FALSE, cluster = NULL, ...)
Arguments
object , model |
a |
... |
[logical] arguments passed to lower level methods. |
linfct |
[matrix or vector of character] the linear hypotheses to be tested. Same as the argument |
rhs |
[vector] the right hand side of the linear hypotheses to be tested. |
robust |
[logical] should robust standard error be used? Otherwise rescale the influence function with the standard error obtained from the information matrix. |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
df |
[character] method used to estimate the degree of freedoms of the Wald statistic: Satterthwaite |
Details
Whenever the argument linfct is not a matrix, it is passed to the function createContrast
to generate the contrast matrix and, if not specified, rhs.
Since only one degree of freedom can be specify in a glht object and it must be an integer, the degree of freedom of the denominator of an F test simultaneously testing all hypotheses is retained, after rounding.
Argument rhs and null are equivalent.
This redondance enable compatibility between lava::compare
, compare2
, multcomp::glht
, and glht2
.
Value
A glht
object.
See Also
createContrast
to create contrast matrices.
estimate2
to pre-compute quantities for the small sample correction.
Examples
library(multcomp)
## Simulate data
mSim <- lvm(c(Y1,Y2,Y3)~ beta * eta, Z1 ~ E, Z2 ~ E, Age[40:5]~1)
latent(mSim) <- "eta"
set.seed(10)
n <- 1e2
df.data <- lava::sim(mSim, n, latent = FALSE, p = c(beta = 1))
#### Inference on a single model ####
e.lvm <- estimate(lvm(Y1~E), data = df.data)
summary(glht2(e.lvm, linfct = c("Y1~E + Y1","Y1")))
#### Inference on separate models ####
## fit separate models
lvmX <- estimate(lvm(Z1 ~ E), data = df.data)
lvmY <- estimate(lvm(Z2 ~ E + Age), data = df.data)
lvmZ <- estimate(lvm(c(Y1,Y2,Y3) ~ eta, eta ~ E),
data = df.data)
#### create mmm object ####
e.mmm <- mmm(X = lvmX, Y = lvmY, Z = lvmZ)
#### create contrast matrix ####
resC <- createContrast(e.mmm, linfct = "E")
#### adjust for multiple comparisons ####
e.glht2 <- glht2(e.mmm, linfct = c(X="E"), df = FALSE)
summary(e.glht2)
Hessian With Small Sample Correction.
Description
Extract the hessian from a latent variable model, with small sample correction
Usage
hessian2(object, indiv, cluster, as.lava, ...)
## S3 method for class 'lvmfit'
hessian2(
object,
indiv = FALSE,
cluster = NULL,
as.lava = TRUE,
ssc = lava.options()$ssc,
...
)
## S3 method for class 'lvmfit2'
hessian2(object, indiv = FALSE, cluster = NULL, as.lava = TRUE, ...)
Arguments
object |
a |
indiv |
[logical] If |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
as.lava |
[logical] if |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the hessian.
Value
An array containing the second derivative of the likelihood relative to each sample (dim 3) and each pair of model coefficients (dim 1,2).
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
n <- 5e1
p <- 3
X.name <- paste0("X",1:p)
link.lvm <- paste0("Y~",X.name)
formula.lvm <- as.formula(paste0("Y~",paste0(X.name,collapse="+")))
m <- lvm(formula.lvm)
distribution(m,~Id) <- Sequence.lvm(0)
set.seed(10)
d <- lava::sim(m,n)
#### latent variable models ####
e.lvm <- estimate(lvm(formula.lvm),data=d)
hessian2(e.lvm)
Compute the Hessian Matrix From the Conditional Moments
Description
Compute the Hessian matrix from the conditional moments.
Usage
.hessian2(
dmu,
dOmega,
d2mu,
d2Omega,
epsilon,
OmegaM1,
missing.pattern,
unique.pattern,
name.pattern,
grid.mean,
grid.var,
grid.hybrid,
name.param,
leverage,
n.cluster,
weights
)
Details
calc_hessian
will perform the computation individually when the
argument index.Omega
is not null.
Influence Function With Small Sample Correction.
Description
Extract the influence function from a latent variable model.
It is similar to lava::iid
but with small sample correction.
Usage
iid2(object, ...)
## S3 method for class 'lvmfit'
iid2(
object,
robust = TRUE,
cluster = NULL,
as.lava = TRUE,
ssc = lava.options()$ssc,
...
)
## S3 method for class 'lvmfit2'
iid2(object, robust = TRUE, cluster = NULL, as.lava = TRUE, ...)
## S3 method for class 'lvmfit2'
iid(x, robust = TRUE, cluster = NULL, as.lava = TRUE, ...)
Arguments
object , x |
a |
... |
additional argument passed to |
robust |
[logical] if |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
as.lava |
[logical] if |
ssc |
[character] method used to correct the small sample bias of the variance coefficients ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the variance-covariance matrix.
Value
A matrix containing the 1st order influence function relative to each sample (in rows) and each model coefficient (in columns).
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
n <- 5e1
p <- 3
X.name <- paste0("X",1:p)
link.lvm <- paste0("Y~",X.name)
formula.lvm <- as.formula(paste0("Y~",paste0(X.name,collapse="+")))
m <- lvm(formula.lvm)
distribution(m,~Id) <- Sequence.lvm(0)
set.seed(10)
d <- sim(m,n)
#### latent variable model ####
e.lvm <- estimate(lvm(formula.lvm),data=d)
iid.tempo <- iid2(e.lvm)
Display the i.i.d. Decomposition
Description
Extract the i.i.d. decomposition and display it along with the corresponding coefficient.
Usage
iid2plot(object, param)
Arguments
object |
a |
param |
[character] name of one of the model parameters. |
Jackknife iid Decomposition from Model Object
Description
Extract iid decomposition (i.e. influence function) from model object.
Usage
iidJack(object, ...)
## Default S3 method:
iidJack(
object,
data = NULL,
grouping = NULL,
cpus = 1,
keep.warnings = TRUE,
keep.error = TRUE,
cl = NULL,
trace = TRUE,
...
)
Arguments
object |
a object containing the model. |
... |
[internal] only used by the generic method. |
data |
[data.frame] dataset used to perform the jackknife. |
grouping |
[vector] variable defining cluster of observations that will be simultaneously removed by the jackknife. |
cpus |
[integer >0] the number of processors to use. If greater than 1, the fit of the model and the computation of the influence function for each jackknife sample is performed in parallel. |
keep.warnings |
[logical] keep warning messages obtained when estimating the model with the jackknife samples. |
keep.error |
[logical]keep error messages obtained when estimating the model with the jackknife samples. |
cl |
[cluster] a parallel socket cluster generated by |
trace |
[logical] should a progress bar be used to trace the execution of the function |
Value
A matrix with in row the samples and in columns the parameters.
Examples
n <- 20
set.seed(10)
mSim <- lvm(c(Y1,Y2,Y3,Y4,Y5) ~ 1*eta)
latent(mSim) <- ~eta
categorical(mSim, K=2) <- ~G
transform(mSim, Id ~ eta) <- function(x){1:NROW(x)}
dW <- lava::sim(mSim, n, latent = FALSE)
#### LVM ####
## Not run:
m1 <- lvm(c(Y1,Y2,Y3,Y4,Y5) ~ 1*eta)
latent(m1) <- ~eta
regression(m1) <- eta ~ G
e <- estimate(m1, data = dW)
iid1 <- iidJack(e)
iid2 <- iid(e)
attr(iid2, "bread") <- NULL
apply(iid1,2,sd)
apply(iid2,2,sd)
quantile(iid2 - iid1)
## End(Not run)
Expected Information With Small Sample Correction.
Description
Extract the expected information matrix from a latent variable model.
Similar to lava::information
but with small sample correction.
Usage
information2(object, as.lava, ssc, ...)
## S3 method for class 'lvmfit'
information2(object, as.lava = TRUE, ssc = lava.options()$ssc, ...)
## S3 method for class 'lvmfit2'
information2(object, as.lava = TRUE, ...)
## S3 method for class 'lvmfit2'
information(x, ...)
Arguments
object , x |
a |
as.lava |
[logical] if |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
... |
additional argument passed to |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the information matrix.
Value
A matrix with as many rows and columns as the number of coefficients.
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
n <- 5e1
p <- 3
X.name <- paste0("X",1:p)
link.lvm <- paste0("Y~",X.name)
formula.lvm <- as.formula(paste0("Y~",paste0(X.name,collapse="+")))
m <- lvm(formula.lvm)
distribution(m,~Id) <- Sequence.lvm(0)
set.seed(10)
d <- lava::sim(m,n)
#### linear models ####
e.lm <- lm(formula.lvm,data=d)
#### latent variable models ####
e.lvm <- estimate(lvm(formula.lvm),data=d)
information(e.lvm)
information2(e.lvm)
information2(e.lvm)[1:4,1:4] - solve(vcov(e.lm))
Compute the Expected Information Matrix From the Conditional Moments
Description
Compute the expected information matrix from the conditional moments.
Usage
.information2(
dmu,
dOmega,
OmegaM1,
missing.pattern,
unique.pattern,
name.pattern,
grid.mean,
grid.var,
name.param,
leverage,
weights = NULL,
n.cluster
)
Details
calc_information
will perform the computation individually when the
argument index.Omega
is not null.
Normalize var1 and var2
Description
Convert var1 and var2 from formula or covariance to character.
Usage
initVarLink(
var1,
var2,
rep.var1 = FALSE,
format = "list",
Slink = c(lava.options()$symbols[1], "~"),
Scov = lava.options()$symbols[2]
)
initVarLinks(var1, format = "list", ...)
Arguments
var1 |
[character or formula] the exogenous variable of the new link or a formula describing the link. |
var2 |
[character] the endogenous variable of the new link.
Disregarded if the argument |
rep.var1 |
[logical] should var1 be duplicated to match var2 length.
Only active if |
format |
[character] should the name of the variable be returned ( |
Slink |
[character] the symbol for regression link. |
Scov |
[character] the symbol for covariance link. |
... |
argument to be passed to |
Value
See argument format
.
Examples
initVarLink(y ~ x1)
initVarLink("y ~ x1")
initVarLink(y ~ x1 + x2)
initVarLink("y ~ x1 + x2")
initVarLink(y ~ x1 + x2, rep.var1 = TRUE)
initVarLink(y ~ x1 + x2, rep.var1 = TRUE, format = "formula")
initVarLink(y ~ x1 + x2, rep.var1 = TRUE, format = "txt.formula")
initVarLink("y", "x1", format = "formula")
initVarLink("y ~ x1:0|1")
initVarLinks(y ~ x1)
initVarLinks("y ~ x1")
initVarLinks(c("y ~ x1","y~ x2"))
initVarLinks(c(y ~ x1,y ~ x2))
initVarLinks(c("y ~ x1","y ~ x2"), format = "formula")
initVarLinks(c(y ~ x1,y ~ x2), format = "formula")
initVarLinks(c("y ~ x1","y~ x2"), format = "txt.formula")
initVarLinks(c(y ~ x1,y ~ x2), format = "txt.formula")
Integrate a Gaussian/Student Density over a Triangle
Description
Consider a univariate random variable X, two multivariate random variables Y and Z, and t1 and t2 two real numbers. This function can compute either P[|X|>t1,|X]>|Y1|,...,|X]>|Yp|] if zmin is not specified, P[|Z1|<t2,...,|Zq|<t2,|X|>t1,|X]>|Y1|,...,|X]>|Yp|] if zmin is specified.
Usage
intDensTri(
mu,
Sigma,
df,
n,
x.min,
z.max = NULL,
type = "double",
proba.min = 1e-06,
prune = NULL,
distribution = "pmvnorm"
)
Arguments
mu |
[numeric vector] the expectation of the joint distribution. |
Sigma |
[matrix] the variance-covariance of the joint distribution. |
df |
[integer > 0] the degree of freedom of the joint Student's t distribution.
Only used when |
n |
[integer > 0] number of points for the numerical integration. |
x.min |
[numeric] the minimum value along the x axis. |
z.max |
[numeric vector, optional] the maximum value along the z axis. Define the dimension of Z. |
type |
[character] the type of mesh to be used.
Can be |
proba.min |
[numeric 0-1] the probability used to find the maximum value along the x axis.
Only used if |
prune |
[integer >0] number of standard deviations after which the domain ends along the x axis. |
distribution |
[character] type of joint distribution.
Can be |
Details
Argument type
:
-
\"raw\"
: mesh with points inside the domain -
\"double\"
: mesh with points outside the domain -
\"fine\"
: mesh with points inside the domain plus additional rectangles trying to fill the missing domain.
Argument Sigma
and mu
:
define the mean and variance-covariance of the random variables X, Y, Z
(in this order). The length of the argument z.max
is used to define the dimension of Z.
The dimension of X is always 1.
Value
A numeric.
Examples
library(mvtnorm)
p <- 2
Sigma <- diag(p)
mu <- rep(0, p)
## bivariate normal distribution
z2 <- qmvt(0.975, mean = mu, sigma = Sigma, df = 1e3)$quantile
# compute integral
intDensTri(mu = mu, Sigma = Sigma, n=5, x.min=0, type = "fine")$value-1/2
intDensTri(mu = mu, Sigma = Sigma, n=30, x.min=0, type = "raw")$value-1/2
intDensTri(mu = mu, Sigma = Sigma, n=50, x.min=0, type = "raw")$value-1/2
intDensTri(mu = mu, Sigma = Sigma, df = 5, n=5, x.min=0, distribution = "pmvt")$value-1/2
res <- intDensTri(mu = mu, Sigma = Sigma, df = 5, n=10, x.min=0, distribution = "pmvt")
res$value-1/2
ggplot2::autoplot(res)
## trivariate normal distribution
## Not run:
p <- 3
Sigma <- diag(p)
mu <- rep(0, p)
res2 <- intDensTri(mu = mu, Sigma = Sigma, n=5, x.min = 0, z.max = 10)
ggplot2::autoplot(res2)
ggplot2::autoplot(res2, coord.plot = c("x","z1"))
res2
## End(Not run)
#### when the distribution is far from 0
## Not run:
eq1 <- intDensTri(mu = c(10,0), Sigma = diag(1,2),
x.min = 2, n=10)
eq1$value-1
ggplot2::autoplot(eq1)
eq2 <- intDensTri(mu = c(10,0,0), Sigma = diag(1,3),
x.min=2, z.max = 10, type = "raw",
n=10)
ggplot2::autoplot(eq2, coord.plot = c("y1","z1"))
eq2$value-1
## more variables
p <- 5
Sigma <- diag(p)
mu <- rep(0, p)
res2 <- intDensTri(mu = mu, Sigma = Sigma, n=5, x.min = 1, z.max = c(2,2))
res2$grid
## End(Not run)
Leverage With Small Sample Correction.
Description
Extract leverage values from a latent variable model, with small sample correction.
Usage
leverage2(object, format, ssc, ...)
## S3 method for class 'lvmfit'
leverage2(object, format = "wide", ssc = lava.options()$ssc, ...)
## S3 method for class 'lvmfit2'
leverage2(object, format = "wide", ...)
Arguments
object |
a |
format |
[character] Use |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
... |
additional argument passed to |
Details
The leverage are defined as the partial derivative of the fitted values with respect to the observations.
leverage_i = \frac{\partial \hat{Y}_i}{\partial Y_i}
See Wei et al. (1998).
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the leverage.
Value
a matrix containing the leverage relative to each sample (in rows) and each endogenous variable (in column).
References
Bo-Cheng Wei et al., Generalized Leverage and its applications (1998), Scandinavian Journal of Statistics 25:1:25-37.
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
set.seed(10)
m <- lvm(Y1~eta,Y2~eta,Y3~eta)
latent(m) <- ~eta
d <- lava::sim(m,20, latent = FALSE)
#### latent variable models ####
e.lvm <- estimate(m, data = d)
leverage2(e.lvm)
Power of a Matrix
Description
Compute the power of a matrix.
Usage
matrixPower(object, power, symmetric, tol = 1e-12, print.warning = TRUE)
Arguments
object |
a matrix. |
power |
[numeric] power to be applied to the matrix. |
symmetric |
[logical] is the matrix symmetric? Argument passed to the function |
tol |
[numeric >0] the threshold under which the eigenvalues are set to 0. |
print.warning |
[logical] should a warning be print when some or the eigenvalues are not strictly positive. |
Value
A matrix.
Examples
## symmetric matrix
set.seed(10)
M <- matrix(rnorm(20*6),20,6)
Sigma <- var(M)
Sigma.half <- matrixPower(Sigma, power = 1/2, symmetric = TRUE)
round(Sigma.half %*% Sigma.half - Sigma,5)
iSigma <- matrixPower(Sigma, power = -1, symmetric = TRUE)
round(iSigma %*% Sigma,5)
iSigma.half <- matrixPower(Sigma, power = -1/2, symmetric = TRUE)
round(iSigma.half %*% iSigma.half - iSigma,5)
## non symmetric matrix
set.seed(10)
M <- matrix(abs(rnorm(9)), 3, 3) + diag(1,3,3)
M-t(M)
iM <- matrixPower(M, power = -1, symmetric = FALSE)
round(iM %*% M,5)
iM.half <- matrixPower(M, power = -1/2, symmetric = FALSE)
round(iM.half %*% iM.half %*% M,5)
Data-driven Extension of a Latent Variable Model
Description
Procedure adding relationship between variables that are supported by the data.
Usage
modelsearch2(
object,
link,
data,
method.p.adjust,
method.maxdist,
n.sample,
na.omit,
alpha,
nStep,
trace,
cpus
)
## S3 method for class 'lvmfit'
modelsearch2(
object,
link = NULL,
data = NULL,
method.p.adjust = "fastmax",
method.maxdist = "approximate",
n.sample = 1e+05,
na.omit = TRUE,
alpha = 0.05,
nStep = NULL,
trace = TRUE,
cpus = 1
)
Arguments
object |
a |
link |
[character, optional for |
data |
[data.frame, optional] the dataset used to identify the model |
method.p.adjust |
[character] the method used to adjust the p.values for multiple comparisons.
Can be any method that is valid for the |
method.maxdist |
[character] the method used to estimate the distribution of the max statistic.
|
n.sample |
[integer, >0] number of samples used in the resampling approach. |
na.omit |
should tests leading to NA for the test statistic be ignored. Otherwise this will stop the selection process. |
alpha |
[numeric 0-1] the significance cutoff for the p-values. When the p-value is below, the corresponding link will be added to the model and the search will continue. Otherwise the search will stop. |
nStep |
the maximum number of links that can be added to the model. |
trace |
[logical] should the execution of the function be traced? |
cpus |
the number of cpus that can be used for the computations. |
Details
method.p.adjust = "max"
computes the p-values based on the distribution of the max statistic.
This max statistic is the max of the square root of the score statistic.
The p-value are computed integrating the multivariate normal distribution.
method.p.adjust = "fastmax"
only compute the p-value for the largest statistic.
It is faster than "max"
and lead to identical results.
method.p.adjust = "gof"
keep adding links until the chi-squared test (of correct specification of the covariance matrix) is no longer significant.
Value
A list containing:
sequenceTest: the sequence of test that has been performed.
sequenceModel: the sequence of models that has been obtained.
sequenceQuantile: the sequence of rejection threshold. Optional.
sequenceIID: the influence functions relative to each test. Optional.
sequenceSigma: the covariance matrix relative to each test. Optional.
initialModel: the model before the sequential search.
statistic: the argument
statistic
.method.p.adjust: the argument
method.p.adjust
.alpha: [numeric 0-1] the significance cutoff for the p-values.
cv: whether the procedure has converged.
Examples
## simulate data
mSim <- lvm()
regression(mSim) <- c(y1,y2,y3,y4)~u
regression(mSim) <- u~x1+x2
categorical(mSim,labels=c("A","B","C")) <- "x2"
latent(mSim) <- ~u
covariance(mSim) <- y1~y2
transform(mSim, Id~u) <- function(x){1:NROW(x)}
set.seed(10)
df.data <- lava::sim(mSim, n = 1e2, latent = FALSE)
## only identifiable extensions
m <- lvm(c(y1,y2,y3,y4)~u)
latent(m) <- ~u
addvar(m) <- ~x1+x2
e <- estimate(m, df.data)
## Not run:
resSearch <- modelsearch(e)
resSearch
resSearch2 <- modelsearch2(e, nStep = 2)
resSearch2
## End(Not run)
## some extensions are not identifiable
m <- lvm(c(y1,y2,y3)~u)
latent(m) <- ~u
addvar(m) <- ~x1+x2
e <- estimate(m, df.data)
## Not run:
resSearch <- modelsearch(e)
resSearch
resSearch2 <- modelsearch2(e)
resSearch2
## End(Not run)
## for instance
mNI <- lvm(c(y1,y2,y3)~u)
latent(mNI) <- ~u
covariance(mNI) <- y1~y2
## estimate(mNI, data = df.data)
## does not converge
Compute Key Quantities of a Latent Variable Model
Description
Compute conditional mean, conditional variance, their first and second derivative regarding model parameters, as well as various derivatives of the log-likelihood.
Usage
moments2(
object,
param,
data,
weights,
Omega,
Psi,
initialize,
usefit,
update.dmoment,
update.d2moment,
score,
information,
hessian,
vcov,
dVcov,
dVcov.robust,
residuals,
leverage,
derivative
)
## S3 method for class 'lvm'
moments2(
object,
param = NULL,
data = NULL,
weights = NULL,
Omega = NULL,
Psi = NULL,
initialize = TRUE,
usefit = TRUE,
update.dmoment = TRUE,
update.d2moment = TRUE,
score = TRUE,
information = TRUE,
hessian = TRUE,
vcov = TRUE,
dVcov = TRUE,
dVcov.robust = TRUE,
residuals = TRUE,
leverage = TRUE,
derivative = "analytic"
)
## S3 method for class 'lvmfit'
moments2(
object,
param = NULL,
data = NULL,
weights = NULL,
Omega = NULL,
Psi = NULL,
initialize = TRUE,
usefit = TRUE,
update.dmoment = TRUE,
update.d2moment = TRUE,
score = TRUE,
information = TRUE,
hessian = TRUE,
vcov = TRUE,
dVcov = TRUE,
dVcov.robust = TRUE,
residuals = TRUE,
leverage = TRUE,
derivative = "analytic"
)
Arguments
object |
a latent variable model. |
param |
[numeric vector] value of the model parameters if different from the estimated ones. |
data |
[data.frame] dataset if different from the one used to fit the model. |
Psi |
[matrix] Average first order bias in the residual variance. Only necessary for computing adjusted residuals. |
initialize |
[logical] Pre-compute quantities dependent on the data but not on the parameters values. |
usefit |
[logical] Compute key quantities based on the parameter values. |
update.dmoment |
[logical] should the first derivative of the moments be computed/updated? |
update.d2moment |
[logical] should the second derivative of the the moments be computed/updated? |
score |
[logical] should the score be output? |
information |
[logical] should the expected information be output? |
hessian |
[logical] should the hessian be output? |
vcov |
[logical] should the variance-covariance matrix based on the expected information be output? |
dVcov |
[logical] should the derivative of the variance-covariance matrix be output? |
dVcov.robust |
[logical] should the derivative of the robust variance-covariance matrix be output? |
Details
For lvmfit objects, there are two levels of pre-computation:
a basic one that do no involve the model coefficient (
conditionalMoment.lvm
).an advanced one that require the model coefficients (
conditionalMoment.lvmfit
).
Examples
m <- lvm(Y1~eta,Y2~eta,Y3~eta)
latent(m) <- ~eta
d <- lava::sim(m,1e2)
e <- estimate(m, d)
## basic pre-computation
res1 <- moments2(e, data = d, initialize = TRUE, usefit = FALSE,
score = TRUE, information = TRUE, hessian = TRUE, vcov = TRUE,
dVcov = TRUE, dVcov.robust = TRUE, residuals = TRUE, leverage = FALSE,
derivative = "analytic")
res1$skeleton$param$Sigma
## full pre-computation
res2 <- moments2(e, param = coef(e), data = d, initialize = TRUE, usefit = TRUE,
score = TRUE, information = TRUE, hessian = TRUE, vcov = TRUE,
dVcov = TRUE, dVcov.robust = TRUE, residuals = TRUE, leverage = FALSE,
derivative = "analytic")
res2$moment$Omega
Find the Number of Steps Performed During the Sequential Testing
Description
Find the number of steps performed during the sequential testing.
Usage
nStep(object)
## S3 method for class 'modelsearch2'
nStep(object)
Arguments
object |
a |
Value
an integer.
Examples
## Not run:
mSim <- lvm(Y~G+X1+X2)
addvar(mSim) <- ~Z1+Z2+Z3+Z4+Z5+Z6
df.data <- lava::sim(mSim, 1e2)
mBase <- lvm(Y~G)
addvar(mBase) <- ~X1+X2+Z1+Z2+Z3+Z4+Z5+Z6
e.lvm <- estimate(mBase, data = df.data)
res <- modelsearch2(e.lvm, method.p.adjust = "holm")
nStep(res)
## End(Not run)
Effective Sample Size.
Description
Extract the effective sample size, i.e. sample size minus the loss in degrees of freedom caused by the estimation of the parameters.
Usage
nobs2(object, ssc, ...)
## S3 method for class 'lvmfit'
nobs2(object, ssc = lava.options()$ssc, ...)
## S3 method for class 'lvmfit2'
nobs2(object, ...)
Arguments
object |
a |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
... |
additional argument passed to |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the leverage.
Value
Numeric vector of length the number of endogenous variables.
See Also
estimate2
to obtain lvmfit2
objects.
Residuals With Small Sample Correction.
Description
Extract residuals from a latent variable model.
Similar to stats::residuals
but with small sample correction.
Usage
residuals2(object, type, format, ssc, ...)
## S3 method for class 'lvmfit'
residuals2(
object,
type = "response",
format = "wide",
ssc = lava.options()$ssc,
...
)
Arguments
object |
a |
type |
[character] the type of residual to extract:
|
format |
[character] Use |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
... |
additional argument passed to |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the residuals.
The raw residuals are defined by observation minus the fitted value:
\varepsilon = (Y_1 - \mu_1, ..., Y_m - \mu_m)
The studentized residuals divided the raw residuals relative to each endogenous variable by the modeled variance of the endogenous variable.
\varepsilon_{stud} =(\frac{Y_1 - \mu_1}{\sigma_1}, ..., \frac{Y_m - \mu_m}{\sigma_m})
The normalized residuals multiply the raw residuals by the inverse of the square root of the modeled residual variance covariance matrix.
\varepsilon_{norm} = \varepsilon \Omega^{-1/2}
Value
a matrix containing the residuals relative to each sample (in rows) and each endogenous variable (in column).
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
set.seed(10)
n <- 101
Y1 <- rnorm(n, mean = 0)
Y2 <- rnorm(n, mean = 0.3)
Id <- findInterval(runif(n), seq(0.1,1,0.1))
data.df <- rbind(data.frame(Y=Y1,G="1",Id = Id),
data.frame(Y=Y2,G="2",Id = Id)
)
#### latent variable models ####
library(lava)
e.lvm <- estimate(lvm(Y ~ G), data = data.df)
residuals(e.lvm)
residuals2(e.lvm)
residuals(e.lvm) - residuals2(e.lvm)
Depreciated Method For Small Sample Correction
Description
Depreciated method for small sample correction, now replaced by the estimate2
method.
Usage
sCorrect(object, ...)
## Default S3 method:
sCorrect(object, ...)
sCorrect(x, ...) <- value
## Default S3 replacement method:
sCorrect(x, ...) <- value
Arguments
object , x |
a |
... |
not used. |
value |
not used. |
Simulate Repeated Measurements over time
Description
Simulate repeated measurements over time (one factor model).
Usage
sampleRepeated(n, n.Xcont = 2, n.Xcat = 2, n.rep = 5, format = "long")
Arguments
n |
[integer] sample size. |
n.Xcont |
[integer] number of continuous covariates acting on the latent variable. |
n.Xcat |
[integer] number of categorical covariates acting on the latent variable. |
n.rep |
[integer] number of measurement of the response variable. |
format |
[character] should the dataset be returned in the |
Value
a data.frame
object.
Examples
sampleRepeated(10, format = "wide")
sampleRepeated(10, format = "long")
Score With Small Sample Correction
Description
Extract the (individual) score a the latent variable model.
Similar to lava::score
but with small sample correction.
Usage
score2(object, indiv, cluster, as.lava, ...)
## S3 method for class 'lvmfit'
score2(
object,
indiv = FALSE,
cluster = NULL,
as.lava = TRUE,
ssc = lava.options()$ssc,
...
)
## S3 method for class 'lvmfit2'
score2(object, indiv = FALSE, cluster = NULL, as.lava = TRUE, ...)
## S3 method for class 'lvmfit2'
score(x, indiv = FALSE, cluster = NULL, as.lava = TRUE, ...)
Arguments
object , x |
a |
indiv |
[logical] If |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
as.lava |
[logical] if |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the confidence intervals.
Value
When argument indiv is TRUE
, a matrix containing the score relative to each sample (in rows)
and each model coefficient (in columns). Otherwise a numeric vector of length the number of model coefficients.
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
n <- 5e1
p <- 3
X.name <- paste0("X",1:p)
link.lvm <- paste0("Y~",X.name)
formula.lvm <- as.formula(paste0("Y~",paste0(X.name,collapse="+")))
m <- lvm(formula.lvm)
distribution(m,~Id) <- Sequence.lvm(0)
set.seed(10)
d <- lava::sim(m,n)
#### linear models ####
e.lm <- lm(Y~X1+X2+X3, data = d)
#### latent variable models ####
m.lvm <- lvm(formula.lvm)
e.lvm <- estimate(m.lvm,data=d)
e2.lvm <- estimate2(m.lvm,data=d)
score.tempo <- score(e2.lvm, indiv = TRUE)
colSums(score.tempo)
Compute the Corrected Score.
Description
Compute the corrected score.
Usage
.score2(
dmu,
dOmega,
epsilon,
OmegaM1,
missing.pattern,
unique.pattern,
name.pattern,
name.param,
name.meanparam,
name.varparam,
n.cluster,
weights
)
Arguments
n.cluster |
[integer >0] the number of observations. |
Regressor of a Formula.
Description
Return the regressor variables contained in the formula
Usage
selectRegressor(object, ...)
## S3 method for class 'formula'
selectRegressor(object, format = "call", ...)
Arguments
object |
a formula |
... |
[internal] Only used by the generic method. |
format |
[character] should an object of format call be returned ( |
Examples
## Not run:
selectRegressor <- lavaSearch2:::selectRegressor
selectRegressor.formula <- lavaSearch2:::selectRegressor.formula
selectRegressor(Y1~X1+X2)
selectRegressor(Y1~X1+X2, format = "vars")
selectRegressor(Y1~X1+Y1)
selectRegressor(Y1+Y2~X1+Y1, format = "vars")
selectRegressor(~X1+X2)
selectRegressor(~X1+X2, format = "vars")
## End(Not run)
Response Variable of a Formula
Description
Return the response variable contained in the formula.
Usage
selectResponse(object, ...)
## S3 method for class 'formula'
selectResponse(object, format = "call", ...)
Arguments
object |
a formula |
... |
[internal] Only used by the generic method. |
format |
[character] should an object of type call be returned ( |
Value
See argument format
.
Examples
## Not run:
selectResponse <- lavaSearch2:::selectResponse
selectResponse.formula <- lavaSearch2:::selectResponse.formula
selectResponse(Y1~X1+X2)
selectResponse(Y1~X1+X2, format = "vars")
selectResponse(Surv(event,time)~X1+X2, format = "vars")
selectResponse(Y1~X1+Y1)
selectResponse(Y1+Y2~X1+Y1, format = "vars")
selectResponse(~X1+X2)
selectResponse(~X1+X2, format = "vars")
## End(Not run)
Set a Link to a Value
Description
Generic interface to set a value to a link in a lvm
object.
Usage
setLink(object, ...)
## S3 method for class 'lvm'
setLink(object, var1, var2, value, warnings = FALSE, ...)
Arguments
object |
a |
... |
[internal] only used by the generic method. |
var1 |
[character or formula] the exogenous variable of the new link or a formula describing the link to be added to the lvm. |
var2 |
[character] the endogenous variable of the new link. Disregarded if the argument |
value |
[numeric] the value at which the link should be set. |
warnings |
[logical] should a warning be displayed if the link is not found in the |
Examples
library(lava)
set.seed(10)
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u~x1+x2
latent(m) <- ~u
covariance(m) <- y1 ~ y2
m1 <- setLink(m, y3 ~ u, value = 1)
estimate(m1, lava::sim(m,1e2))
# m1 <- setLink(m, u ~ y3, value = 1)
m2 <- setLink(m, y1 ~ y2, value = 0.5)
estimate(m2, lava::sim(m,1e2))
Pre-computation for the Score
Description
Pre-compute quantities that are necessary to compute the score of a lvm model.
Usage
skeleton(object, X, endogenous, latent, n.cluster, index.Omega)
skeletonDtheta(
object,
X,
endogenous,
latent,
missing.pattern,
unique.pattern,
name.pattern,
n.cluster,
index.Omega
)
skeletonDtheta2(object)
Arguments
object |
a |
X |
[matrix] design matrix containing the covariates for each endogeneous and latent variable. |
endogenous |
[character vector] the name of the endogeneous variables. |
latent |
[character vector] the name of the latent variables. |
Details
When the user specifies names for the coefficients (e.g. Y1[mu:sigma]) or uses constraints (Y1~beta*X1), as.lava=FALSE
will use the names specified by the user (e.g. mu, sigma, beta)
while as.lava=TRUE
will use the name of the first link defining the coefficient.
Examples
## Not run:
skeleton <- lavaSearch2::skeleton
skeleton.lvm <- lavaSearch2::skeleton.lvm
skeleton.lvmfit <- lavaSearch2::skeleton.lvmfit
## without constrain
m <- lvm(Y1~X1+X2+eta,Y2~X3+eta,Y3~eta)
latent(m) <- ~eta
e <- estimate(m, lava::sim(m,1e2))
M.data <- as.matrix(model.frame(e))
skeleton(e$model, as.lava = TRUE,
name.endogenous = endogenous(e), n.endogenous = 3,
name.latent = latent(e),
update.value = FALSE)
skeleton(e, data = M.data, p = pars(e), as.lava = TRUE,
name.endogenous = endogenous(e), n.endogenous = 3,
name.latent = latent(e),
update.value = TRUE)
## with constrains
m <- lvm(Y[mu:sigma] ~ beta*X1+X2)
e <- estimate(m, lava::sim(m,1e2))
M.data <- as.matrix(model.frame(e))
skeleton(e$model, as.lava = TRUE,
name.endogenous = "Y", n.endogenous = 1,
name.latent = NULL,
update.value = FALSE)$skeleton
skeleton(e, data = M.data, p = pars(e), as.lava = FALSE,
name.endogenous = "Y", n.endogenous = 1,
name.latent = NULL,
update.value = FALSE)$skeleton
## End(Not run)
Display the Type 1 Error Rate
Description
Display the type 1 error rate from the simulation results.
Usage
## S3 method for class 'calibrateType1'
summary(
object,
robust = FALSE,
type = "type1error",
alpha = 0.05,
log.transform = TRUE,
digits = 5,
print = TRUE,
...
)
Arguments
object |
output of the |
robust |
[character] should the results be displayed for both model-based and robust standard errors ( |
type |
[character] should the type 1 error rate be diplayed ( |
alpha |
[numeric, 0-1] the confidence levels. |
log.transform |
[logical] should the confidence intervals be computed on the logit scale. |
digits |
[integer >0] the number of decimal places to use when displaying the summary. |
print |
should the summary be printed in the terminal. |
... |
[internal] only used by the generic method. |
Outcome of Linear Hypothesis Testing
Description
Estimates, p-values, and confidence intevals for linear hypothesis testing, possibly adjusted for multiple comparisons.
Usage
## S3 method for class 'glht2'
summary(
object,
confint = TRUE,
conf.level = 0.95,
transform = NULL,
seed = NULL,
rowname.rhs = TRUE,
...
)
Arguments
object |
a |
confint |
[logical] should confidence intervals be output |
conf.level |
[numeric 0-1] level of the confidence intervals. |
transform |
[function] function to backtransform the estimates, standard errors, null hypothesis, and the associated confidence intervals
(e.g. |
seed |
[integer] value that will be set before adjustment for multiple comparisons to ensure reproducible results.
Can also be |
rowname.rhs |
[logical] when naming the hypotheses, add the right-hand side (i.e. "X1-X2=0" instead of "X1-X2"). |
... |
argument passed to |
summary Method for modelsearch2 Objects
Description
summary method for modelsearch2 objects.
Usage
## S3 method for class 'modelsearch2'
summary(object, print = TRUE, ...)
Arguments
object |
output of the |
print |
should the summary be printed in the terminal. |
... |
[internal] only used by the generic method. |
Details
The column dp.Info
contains the percentage of extended models (i.e. model with one additional link)
for which the information matrix evaluated at the value of the parameters of the initial model is non positive definie.
Latent Variable Model Summary After Small Sample Correction
Description
Summarize a fitted latent variable model.
Similar to stats::summary
with small sample correction.
Usage
summary2(object, robust, cluster, digit, ...)
## S3 method for class 'lvmfit'
summary2(
object,
robust = FALSE,
cluster = NULL,
digit = max(5, getOption("digit")),
ssc = lava.options()$ssc,
df = lava.options()$df,
...
)
## S3 method for class 'lvmfit2'
summary2(
object,
robust = FALSE,
cluster = NULL,
digit = max(5, getOption("digit")),
...
)
## S3 method for class 'lvmfit2'
summary(
object,
robust = FALSE,
cluster = NULL,
digit = max(5, getOption("digit")),
...
)
Arguments
object |
a |
robust |
[logical] should robust standard errors be used instead of the model based standard errors? Should be |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
digit |
[integer > 0] the number of decimal places to use when displaying the summary. |
... |
[logical] arguments passed to lower level methods. |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
df |
[character] method used to estimate the degree of freedoms of the Wald statistic: Satterthwaite |
Details
summary2
is the same as summary
except that it first computes the small sample correction (but does not store it).
So if summary2
is to be called several times,
it is more efficient to pre-compute the quantities for the small sample correction
using sCorrect
and then call summary2
.
summary2
returns an object with an element table2
containing the estimates, standard errors, degrees of freedom,
upper and lower limits of the confidence intervals, test statistics, and p-values.
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
m <- lvm(Y~X1+X2)
set.seed(10)
d <- lava::sim(m, 2e1)
#### latent variable models ####
e.lvm <- estimate(m, data = d)
summary(e.lvm)$coef
summary2(e.lvm)
summary2(e.lvm, ssc = "none")
Symmetrize a Matrix
Description
Complete the upper (or lower) extra-diagonal terms in order to obtain a symmetric matrix.
Usage
symmetrize(M, update.upper = TRUE)
Arguments
M |
a matrix. |
update.upper |
[logical] should the upper extra diagonal terms be updated using the lower extra diagonal terms? |
Examples
symmetrize <- lavaSearch2:::symmetrize
## example
M <- matrix(NA, 4, 4)
M[lower.tri(M)] <- 1:6
symmetrize(M, update.upper = TRUE) # good
M[upper.tri(M, diag = FALSE)] <- M[lower.tri(M, diag = FALSE)]
M # wrong
Apply Transformation to Summary Table
Description
Update summary table according to a transformation, e.g. log-transformtion. P-values are left unchanged but estimates, standard errors, and confidence intervals are updated.
Usage
transformSummaryTable(object, transform = NULL)
Arguments
object |
A data.frame with columns estimate, se, lower, upper. |
transform |
the name of a transformation or a function. |
Value
a data.frame
Run an Expression and Catch Warnings and Errors
Description
Similar to try
but also returns warnings.
Usage
tryWithWarnings(expr)
Arguments
expr |
the line of code to be evaluated |
Details
from https://stackoverflow.com/questions/4948361/how-do-i-save-warnings-and-errors-as-output-from-a-function
Value
A list containing:
value the result of the evaluation of the expression
warnings warning(s) generated during the evaluation of the expression
error error generated during the evaluation of the expression
Examples
FctTest <- function(x){
return(log(x))
}
tryWithWarnings(FctTest(-1))
tryWithWarnings(FctTest(1))
tryWithWarnings(FctTest(xxxx))
Convert Variable Names to Dummy Variables Names.
Description
When dealing with categorical variables, the estimate
function convert the categorical variables into dummy variables.
This function convert a set of variable names to their corresponding name in the model with dummy variables
Usage
var2dummy(object, ...)
## S3 method for class 'list'
var2dummy(object, var, rm.first.factor = TRUE, ...)
## S3 method for class 'lvm'
var2dummy(object, data = NULL, ...)
Arguments
object |
a |
... |
[internal] additional arguments to be passed from |
var |
[character] the variable to be transformed. |
rm.first.factor |
[logical] should the first level of each categorical variable be ignored? |
data |
[data.frame] dataset according to which the model should be updated. |
Examples
## Not run:
var2dummy <- lavaSearch2:::var2dummy
var2dummy.list <- lavaSearch2:::var2dummy.list
var2dummy.lvm <- lavaSearch2:::var2dummy.lvm
m <- lvm()
regression(m) <- c(y1,y2,y3)~u
regression(m) <- u ~ X1+X2
var2dummy(m, var = c("X1","X2"))
categorical(m,labels=c("M","F","MF")) <- ~X1
var2dummy(m, var = c("X1","X2"))
categorical(m,labels=c("1","2","3")) <- ~X2
var2dummy(m, var = c("X1","X2"))
## End(Not run)
Variance-Covariance With Small Sample Correction
Description
Extract the variance-covariance matrix from a latent variable model.
Similar to stats::vcov
but with small sample correction.
Usage
vcov2(object, robust, cluster, as.lava, ...)
## S3 method for class 'lvmfit'
vcov2(
object,
robust = FALSE,
cluster = NULL,
as.lava = TRUE,
ssc = lava.options()$ssc,
...
)
## S3 method for class 'lvmfit2'
vcov2(object, robust = FALSE, cluster = NULL, as.lava = TRUE, ...)
## S3 method for class 'lvmfit2'
vcov(object, robust = FALSE, cluster = NULL, as.lava = TRUE, ...)
Arguments
object |
a |
robust |
[logical] should robust standard errors be used instead of the model based standard errors? Should be |
cluster |
[integer vector] the grouping variable relative to which the observations are iid. |
as.lava |
[logical] if |
... |
additional argument passed to |
ssc |
[character] method used to correct the small sample bias of the variance coefficients: no correction ( |
Details
When argument object is a lvmfit
object, the method first calls estimate2
and then extract the variance-covariance matrix.
Value
A matrix with as many rows and columns as the number of coefficients.
See Also
estimate2
to obtain lvmfit2
objects.
Examples
#### simulate data ####
n <- 5e1
p <- 3
X.name <- paste0("X",1:p)
link.lvm <- paste0("Y~",X.name)
formula.lvm <- as.formula(paste0("Y~",paste0(X.name,collapse="+")))
m <- lvm(formula.lvm)
distribution(m,~Id) <- Sequence.lvm(0)
set.seed(10)
d <- lava::sim(m,n)
#### linear models ####
e.lm <- lm(formula.lvm,data=d)
#### latent variable models ####
e.lvm <- estimate(lvm(formula.lvm),data=d)
vcov0 <- vcov(e.lvm)
vcovSSC <- vcov2(e.lvm)
vcovSSC/vcov0
vcovSSC[1:4,1:4]/vcov(e.lm)
Inverse the Information Matrix
Description
Compute the inverse of the information matrix.
Usage
.info2vcov(information, attr.info = FALSE)
Arguments
information |
[matrix] information matrix to be inverted. |
attr.info |
[logical] should the information matrix be returned as an attribute? |