R/helper_test_MGD.R
In M-E-Steiner/cSEM: Composite-Based Structural Equation Modeling
Defines functions
structureTestMGDDecisions
getRelevantParameters
calculatePr
calculateFR
adjustAlpha
calculateParameterDifference
getParameterNames
calculateDistance
Documented in
adjustAlpha
calculateDistance
calculateFR
calculateParameterDifference
calculatePr
getParameterNames
getRelevantParameters
structureTestMGDDecisions
#' Internal: Matrix difference
#'
#' Calculates the average of the differences between all possible pairs of
#' (symmetric) matrices in a list using a given distance measure.
#'
#' `.matrices` must be a list of at least two matrices. If more than two matrices
#' are supplied the arithmetic mean of the differences between all possible pairs of
#' (symmetric) matrices in a list is computed. Mathematically this is
#' n chose 2. Hence, supplying a large number of matrices will
#' become computationally challenging.
#'
#' Currently two distance measures are supported:
#' \describe{
#' \item{`geodesic`}{(Default) The geodesic distance.}
#' \item{`squared_euclidean`}{The squared Euclidean distance}
#' }
#'
#' @usage calculateDistance(
#' .matrices = NULL,
#' .distance = args_default()$.distance
#' )
#'
#' @inheritParams csem_arguments
#'
#' @return A numeric vector of length one containing the (arithmetic) mean of
#' the differences between all possible pairs of matrices supplied via `.matrices`.
#'
#' @keywords internal
calculateDistance <- function(
.matrices = NULL,
.distance = args_default()$.distance
){
### Checks and errors ========================================================
## Check if list and at least of length 2
if(!is.list(.matrices) && length(.matrices) < 2) {
stop2("`.matrices` must be a list of at least length two.")
}
## Check if matrices are all symmetric
if(!all(sapply(.matrices, matrixcalc::is.symmetric.matrix))) {
stop2("All matrices in `.matrices` must be symmetric.")
}
### Calculation ==============================================================
## Combine all matrices into lists of two
temp <- utils::combn(.matrices, 2, simplify = FALSE)
## Compute the distance measure for each group combination
distances <- lapply(temp, function(x) {
switch (.distance,
"geodesic" = {calculateDG(.matrix1 = x[[1]], .matrix2 = x[[2]])},
"squared_euclidean" = {calculateDL(.matrix1 = x[[1]], .matrix2 = x[[2]])}
)
})
## Compute the mean. Since dG and dL are defined as the mean of difference between
## two matrices (1/2 is part of the definition). This must be corrected if
## more than two matrices are to be compared.
if(length(.matrices) > 2) {
out <- mean(2 * unlist(distances))
} else {
out <- unlist(distances)
}
return(out)
}
#' Internal: Parameter names
#'
#' Based on a model in [lavaan model syntax][lavaan::model.syntax], returns the
#' names of the parameters of the structural
#' model, the measurement/composite model and the weight relationship. Used
#' by [testMGD()] to extract the names of the parameters to compare across groups
#' according to the test proposed by \insertCite{Chin2010;textual}{cSEM}.
#'
#' @usage getParameterNames(
#' .object = args_default()$.object,
#' .model = args_default()$.model
#' )
#'
#' @inheritParams csem_arguments
#' @param .model A model in [lavaan model syntax][lavaan::model.syntax] indicating which
#' parameters (i.e, path (`~`), loadings (`=~`), or weights (`<~`)) should be
#' compared across groups. Defaults to `NULL` in which case all parameters of the model
#' are compared.
#'
#' @return A list with elements `names_path`, `names_loadings`, and `names_weights`
#' containing the names of the structural parameters, the loadings,
#' and the weight to compare across groups.
#'
#' @references
#' \insertAllCited{}
#'
#' @keywords internal
getParameterNames <- function(
.object = args_default()$.object,
.model = args_default()$.model
){
## Summarize
x <- summarize(.object)
if(inherits(.object, "cSEMResults_2ndorder")) {
x12 <- x[[1]]$First_stage$Information
x22 <- x[[1]]$Second_stage$Information
# Extract different types of constructs
construct_type <- c(x12$Model$construct_type, x22$Model$construct_type)
construct_type <- construct_type[unique(names(construct_type))]
measurement_org <- x22$Arguments_original$.model$measurement
names_path_org <- x[[1]]$Second_stage$Estimates$Path_estimates$Name
indicators <- x12$Model_original$indicators
cons_exo <- x12$Model_original$cons_exo
names_cor_exo_cons_org <- x[[1]]$Second_stage$Estimates$Exo_construct_correlation$Name
names_cor_indicator_org <- c(x[[1]]$First_stage$Estimates$Indicator_correlation$Name,
x[[1]]$Second_stage$Estimates$Indicator_correlation$Name)
if(nrow(x[[1]]$First_stage$Estimates$Residual_correlation)!=0){
temp <- x[[1]]$First_stage$Estimates$Residual_correlation$Name
}else{
temp <- NULL
}
if(nrow(x[[1]]$Second_stage$Estimates$Residual_correlation)!=0){
temp1 <- x[[1]]$Second_stage$Estimates$Residual_correlation$Name
}else{
temp1 <- NULL
}
names_cor_measurement_error_org <- c(temp,temp1)
} else {#if no second-order model
x22 <- x[[1]]$Information
# Extract different types of constructs
construct_type <- x22$Model$construct_type
measurement_org <- x22$Model$measurement
names_path_org <- x[[1]]$Estimates$Path_estimates$Name
indicators <- x22$Model$indicators
cons_exo <- x22$Model$cons_exo
names_cor_exo_cons_org <- x[[1]]$Estimates$Exo_construct_correlation$Name
names_cor_indicator_org <- x[[1]]$Estimates$Indicator_correlation$Name
if(nrow(x[[1]]$Estimates$Residual_correlation)!=0){
names_cor_measurement_error_org <- x[[1]]$Estimates$Residual_correlation$Name
}else{
names_cor_measurement_error_org<-NULL
}
}
# FOR FUTURE:
# Perhaps it is better to transfer the following code to the testMGD function directly.
# As getParameterNames should just return the names of a model and not internally select a model?
# Parse model that indicates which parameters should be compared.
# If no model indicating the comparisons is provided, all parameters are compared.
# Correlation among the measurement errors and or indicators cannot be compared yet across groups as
# theta is not part of the csem output
# FOR FUTURE: Add all model parameters to the comparison if no model is specified
if(is.null(.model)) {
if(inherits(.object, "cSEMResults_2ndorder")) {
model_compare <- x22$Arguments_original$.model
# add the correlations among the exogenous constructs
} else {
model_compare <- x22$Model
# add the correlations among the exogenous constructs
# add to model_cor_specified and expand
# Attention if the some of the variables are already in that matrix
}
} else {
model_compare <- parseModel(.model, .check_errors = FALSE)
}
# Check whether the constructs specified in the comparison are equal
# to the constructs in the original model
if(!all(rownames(model_compare$structural) %in% rownames(measurement_org))){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one construct appears in the comparison model but",
" not in the original model."
)
}
# Check whether the construct type specified in the comparison are equal
# to the construct types in the original model
construct_type_compare <- model_compare$construct_type[!is.na(model_compare$construct_type)]
if(!all(construct_type_compare == construct_type[names(construct_type_compare)])){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one construct's type in the comparison model differs ",
"from the original model.")
}
# Check wehther the indicators used in the comparison model also appear
# in the original model
if(!all(colnames(model_compare$measurement) %in% colnames(measurement_org))){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one indicator appears in the comparison model and not in ",
" the original model.")
}
# Check whether the indicators used are correctly assigned, i.e.,
# as in the original model
lapply(rownames(model_compare$measurement), function(x) {
names_ind_compare <- colnames(model_compare$measurement[x,,drop=FALSE][,which(model_compare$measurement[x,,drop=FALSE]==1),drop=FALSE])
names_ind_org <- colnames(measurement_org[x,,drop=FALSE][,which(measurement_org[x,,drop=FALSE]==1),drop=FALSE])
if(!all(names_ind_compare %in% names_ind_org)){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one indicator is not correctly assigned in the comparison model.")
}
})
### Extract names ============================================================
# Extract names of the path to be tested
temp <- outer(rownames(model_compare$structural), colnames(model_compare$structural),
FUN = function(x, y) paste(x, y, sep = " ~ "))
names_path_compare <- t(temp)[t(model_compare$structural) != 0]
if(length(names_path_compare) == 0) {
names_path_compare <- NULL
}
# Check whether the path to compare occur also in the original mode
if(!all(names_path_compare %in% names_path_org)){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one of the paths specified for comparison does not ",
" appear in the original model.")
}
# all correlated variables, i.e., that have been specified with ~~
vars_correlated_compare <- rownames(model_compare$cor_specified)
# Select only those that are indicators
ind_correlated_compare <- intersect(vars_correlated_compare, indicators)
# Extract names of the loadings and measurement error to be tested
names_row <- rownames(model_compare$measurement)
## Select only concepts modeled as common factors. Reorder to be have the
## same order as names_row.
i <- intersect(names(which(construct_type == "Common factor")), names_row)
i <- i[match(names_row, i)]
i <- i[!is.na(i)]
if(length(i) != 0) {
temp <- rep(rownames(model_compare$measurement[i, , drop = FALSE]),
times = rowSums(model_compare$measurement[i, , drop = FALSE]))
if(length(temp)!=0){
# Loadings
temp_n <- model_compare$measurement[i, colSums(model_compare$measurement[i, , drop = FALSE]) != 0, drop = FALSE]
names_loadings_compare <- paste0(temp, " =~ ", colnames(temp_n))
# Measurement error
measurement_error_correlated <- intersect(ind_correlated_compare,colnames(temp_n))
cor_measurement_error <- model_compare$cor_specified[ measurement_error_correlated, measurement_error_correlated]
index <- which(cor_measurement_error == 1, arr.ind = TRUE)
names_correlated_measurement_error_compare <- index
if(nrow(names_correlated_measurement_error_compare) ==0 ){
names_correlated_measurement_error_compare <- NULL
}else{
names_correlated_measurement_error_compare <- paste(rownames(cor_measurement_error)[index[,"row"]],
" ~~ ",
colnames(cor_measurement_error)[index[,"col"]],
sep="")
# Select those that have the same order as in the original model
names_correlated_measurement_error_compare <- names_correlated_measurement_error_compare[names_correlated_measurement_error_compare %in% names_cor_measurement_error_org]
}
}else{
names_loadings_compare <- NULL
names_correlated_measurement_error_compare <- NULL
}
} else {
names_loadings_compare <- NULL
names_correlated_measurement_error_compare <- NULL
}
## Extract names of the weights and indicator correlations (composites) to be tested
# Select only concepts modeled as composite. Reorder that it has the
# same order as names_row.
## Indicator correlation, i.e., variables that belong to a composite
i <- intersect(names(which(construct_type == "Composite")), names_row)
i <- i[match(names_row, i)]
i <- i[!is.na(i)]
if(length(i) != 0) {
temp <- rep(rownames(model_compare$measurement[i, , drop = FALSE]),
times = rowSums(model_compare$measurement[i, , drop = FALSE]))
if(length(temp)!=0){
# weights
temp_n <- model_compare$measurement[i, colSums(model_compare$measurement[i, , drop = FALSE]) != 0, drop = FALSE]
names_weights_compare <- paste0(temp, " <~ ", colnames(temp_n))
# indicator correlations
ind_connected_to_composite <- colnames(temp_n)
indicator_correlated <- intersect(ind_correlated_compare,ind_connected_to_composite)
cor_indicator <- model_compare$cor_specified[indicator_correlated, indicator_correlated]
index <- which(cor_indicator == 1, arr.ind = TRUE)
names_correlated_indicator_compare <- index
# In case that no indicator correlations are compared set it to NULL
if(nrow(names_correlated_indicator_compare) ==0 ){
names_correlated_indicator_compare <- NULL
}else{
names_correlated_indicator_compare <- paste(rownames(cor_indicator)[index[,"row"]],
" ~~ ",
colnames(cor_indicator)[index[,"col"]],
sep="")
# Select those that have the same order as in the original model
names_correlated_indicator_compare <- names_correlated_indicator_compare[names_correlated_indicator_compare %in% names_cor_indicator_org]
}
}else{
names_weights_compare <- NULL
names_correlated_indicator_compare <- NULL
}
} else {
names_weights_compare <- NULL
names_correlated_indicator_compare <- NULL
}
## Exogenous construct correlations
cons_exo_correlated_compare <- intersect(vars_correlated_compare, cons_exo)
# Consider only exogenous constructs
cor_cons_exo <- model_compare$cor_specified[cons_exo_correlated_compare,cons_exo_correlated_compare]
# Which are correlated
index <- which(cor_cons_exo == 1,arr.ind = TRUE)
names_correlated_exo_cons_compare <- index
# In case that no measurement error correlations are compared set it to NULL
if(nrow(names_correlated_exo_cons_compare) == 0){
names_correlated_exo_cons_compare <- NULL
}else{
names_correlated_exo_cons_compare <-paste(rownames(cor_cons_exo)[index[,"row"]],
" ~~ ",
colnames(cor_cons_exo)[index[,"col"]], sep="")
# Select those that have the same order as in the original model
names_correlated_exo_cons_compare <- names_correlated_exo_cons_compare[names_correlated_exo_cons_compare %in% names_cor_exo_cons_org]
}
## Return as list
out <- list(
"names_path" = names_path_compare,
"names_weights" = names_weights_compare,
"names_loadings" = names_loadings_compare,
"names_cor_exo_cons" = names_correlated_exo_cons_compare,
"names_cor_indicator" = names_correlated_indicator_compare,
"names_cor_measurement_error" = names_correlated_measurement_error_compare
)
return(out)
}
#' Internal: Parameter differences across groups
#'
#' Calculate the difference between one or more parameter estimates across
#' all possible pairs of groups (data sets) in `.object`.
#'
#' @usage calculateParameterDifference(
#' .object = args_default()$.object,
#' .model = args_default()$.model
#' )
#'
#' @inheritParams csem_arguments
#' @param .model A model in [lavaan model syntax][lavaan::model.syntax] indicating which
#' parameters (i.e., path (`~`), loadings (`=~`), or weights (`<~`)) should be
#' compared across groups. Defaults to `NULL` in which case all parameters of the model
#' are compared.
#'
#' @return A list of length equal to the number of possible pairs of
#' groups in `.object` (mathematically, this is n choose 2, i.e., 3 if there are three
#' groups and 6 if there are 4 groups). Each list elements is itself a list of
#' three. The first list element contains
#' the difference between parameter estimates of the structural model, the second
#' list element the difference between estimated loadings, and the third
#' the difference between estimated weights.
#'
#' @keywords internal
calculateParameterDifference <- function(
.object = args_default()$.object,
.model = args_default()$.model
){
## Summarize
x <- summarize(.object)
## Get names to compare
names <- getParameterNames(.object, .model = .model)
names_path <- names$names_path
names_loadings <- names$names_loadings
names_weights <- names$names_weights
names_cor_exo_cons <- names$names_cor_exo_cons
# names_cor_measurement_error <- names$names_cor_measurement_error
### Compute differences ======================================================
if(inherits(.object, "cSEMResults_2ndorder")) {
path_estimates <- lapply(x, function(y) {y$Second_stage$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Loading_estimates,
y$Second_stage$Estimates$Loading_estimates)
})
weight_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Weight_estimates,
y$Second_stage$Estimates$Weight_estimates)
})
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Second_stage$Estimates$Exo_construct_correlation})
} else { # no second-order model
path_estimates <- lapply(x, function(y) {y$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {y$Estimates$Loading_estimates})
weight_estimates <- lapply(x, function(y) {y$Estimates$Weight_estimates})
# all exogenous construct correlations
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Estimates$Exo_construct_correlation
})
}
## Select
path_estimates <- lapply(path_estimates, function(y) {
y1 <- y[y$Name %in% names_path, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_path
}
y1
})
loading_estimates <- lapply(loading_estimates, function(y) {
y1 <- y[y$Name %in% names_loadings, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_loadings
}
y1
})
weight_estimates <- lapply(weight_estimates, function(y) {
y1 <- y[y$Name %in% names_weights, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_weights
}
y1
})
cor_cons_exo_estimates <- lapply(cor_cons_exo_estimates, function(y) {
y1 <- y[y$Name %in% names_cor_exo_cons, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_cor_exo_cons
}
y1
})
## Path model
temp <- utils::combn(path_estimates, 2, simplify = FALSE)
diff_path <- lapply(temp, function(y) y[[1]] - y[[2]])
names(diff_path) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
## Loadings
temp <- utils::combn(loading_estimates, 2, simplify = FALSE)
diff_loadings <- lapply(temp, function(y) y[[1]] - y[[2]])
names(diff_loadings) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
## Weights
temp <- utils::combn(weight_estimates, 2, simplify = FALSE)
diff_weights <- lapply(temp, function(y) y[[1]] - y[[2]])
names(diff_weights) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
# Exogenous construct correlations
temp <- utils::combn(cor_cons_exo_estimates, 2, simplify = FALSE)
cor_cons_exo <- lapply(temp, function(y) y[[1]] - y[[2]])
names(cor_cons_exo) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
# merge list together
out <- mapply(function(w,x,y,z) c(w,x,y,z),
w = diff_path,
x = diff_loadings,
y = diff_weights,
z = cor_cons_exo,
SIMPLIFY = FALSE)
return(out)
}
#' Internal: Multiple testing correction
#'
#' Adjust a given significance level `.alpha` to accommodate multiple testing.
#' The following corrections are implemented:
#' \describe{
#' \item{`none`}{(Default) No correction is done.}
#' \item{`bonferroni`}{A Bonferroni correction is done, i.e., alpha is divided by the
#' number of comparisons `.nr_comparisons`.}
#' }
#'
#' @usage adjustAlpha(
#' .alpha = args_default()$.alpha,
#' .approach_alpha_adjust = args_default()$.approach_alpha_adjust,
#' .nr_comparisons = args_default()$.nr_comparisons
#' )
#'
#' @inheritParams csem_arguments
#'
#' @return A vector of (possibly adjusted) significance levels.
#'
#' @keywords internal
adjustAlpha <- function(
.alpha = args_default()$.alpha,
.approach_alpha_adjust = args_default()$.approach_alpha_adjust,
.nr_comparisons = args_default()$.nr_comparisons
){
if(.approach_alpha_adjust == 'none'){
return(.alpha)
}
if(.approach_alpha_adjust == 'bonferroni'){
return(.alpha/.nr_comparisons)
}
}
#' Internal: ANOVA F-test statistic
#'
#' Calculate the ANOVA F-test statistic suggested by
#' \insertCite{Sarstedt2011;textual}{cSEM} in the OTG testing procedure.
#'
#' @usage calculateFR(.resample_sarstedt)
#'
#' @inheritParams csem_arguments
#'
#' @return A named scalar, the test statistic of the ANOVA F-test
#'
#' @references
#' \insertAllCited{}
#'
#' @keywords internal
calculateFR <- function(.resample_sarstedt) {
# In Sarstedt et al. (2011), it is assumed that you have the same number of bootstrap estimates per group.
# Consequently, the total number of estimates equals B*G.
# Since, there might be different numbers of bootstrap estimates per group, we adjusted the calculation of the
# F test statistic.
names_param <- colnames(.resample_sarstedt)[-ncol(.resample_sarstedt)]
# Number of grousp
G <- length(unique(.resample_sarstedt[, "group_id"]))
# Total number of estimates
N <- nrow(.resample_sarstedt)
# Order the estimates by the grouping variable
x <- .resample_sarstedt[order(.resample_sarstedt[, "group_id"]), ]
grand_means <- colMeans(x[, names_param, drop = FALSE])
mean_group <- aggregate(x[, names_param, drop = FALSE], by = list(x[, "group_id"]), FUN = mean)
group_mean_matrix <- apply(mean_group[, -1, drop = FALSE], 2, function(y) rep(y, times = table(x[, "group_id"])))
SS_between <- rowSums((t(group_mean_matrix) - grand_means)^2)
MSA_between <- SS_between/ (G -1)
SS_within <- colSums((x[, names_param, drop = FALSE] - group_mean_matrix)^2)
MSA_within <- SS_within/(N - G)
# Note: in the original paper (B -1) was used. This is incorrect. It should
# be B - G (G = the number of groups). See the F-statistic derived from ANOVA
# for comparision.
MSA_between/MSA_within
}
#' Internal: Calculation of the CDF used in Henseler et al. (2009)
#'
#' Calculates the probability that theta^1 is smaller than or equal to theta^2.
#' See Equation (6) in \insertCite{Sarstedt2011;textual}{cSEM}.
#'
#' @usage calculatePr(.resample_centered = NULL, .parameters_to_compare = NULL)
#'
#' @inheritParams csem_arguments
#'
#' @return A named vector
#'
#' @references
#' \insertAllCited{}
#'
#' @keywords internal
calculatePr <- function(
.resample_centered = NULL,
.parameters_to_compare = NULL
){
# Remove names from .parameters_to_compare
names(.parameters_to_compare) <- NULL
group1 <- .resample_centered[[1]][,.parameters_to_compare,drop=FALSE]
group2 <- .resample_centered[[2]][,.parameters_to_compare,drop=FALSE]
ret <- sapply(.parameters_to_compare, function(x){
# Matrix where each element of the first vector is substracted from the complete second vector
# theta^1 <= theta^2
temp12 <- outer(group2[,x],group1[,x],"-")
g1leqg2<-mean((1+sign(temp12))/2)
# theta^2 <= theta^1 is equal to 1 - g1leqg2, therefore we do not need the following calculations
# temp21 <- outer(group1[,x],group2[,x],"-")
# g2leqg1<-mean((1+sign(temp21))/2)
# # The problem is that Henseler's approach is a one-sided test.
# # Therefore, the p-value is flipped if it larger than 0.5
# # Return
# if(g1leqg2<0.5){
# g1leqg2
# }else{
# 1-g1leqg2
# }
})
names(ret) <- .parameters_to_compare
ret
}
#' Internal: Extract relevant parameters from several cSEMResults_multi
#'
#' Extract the relevant parameters from a cSEMResult_multi object in `.object`.
#'
#' @usage getRelevantParameters(
#' .object = args_default()$.object,
#' .model = args_default()$.model
#' )
#'
#' @inheritParams csem_arguments
#' @param .model A model in [lavaan model syntax][lavaan::model.syntax] indicating which
#' parameters (i.e., path (`~`), loadings (`=~`), or weights (`<~`)) should be
#' compared across groups. Defaults to `NULL` in which case all parameters of the model
#' are compared.
#'
#' @return A list of length equal to the number of groups in `.object`.
#' Each list element is itself a list of three. The first list element contains
#' the relevant parameter estimates of the structural model, the second
#' list element the relevant estimated loadings, and the third
#' the relevant estimated weights.
#'
#' @keywords internal
getRelevantParameters <- function(
.object = args_default()$.object,
.model = args_default()$.model
){
## Summarize
x <- summarize(.object)
## Get names to compare
names <- getParameterNames(.object, .model = .model)
names_path <- names$names_path
names_loadings <- names$names_loadings
names_weights <- names$names_weights
names_cor_exo_cons <- names$names_cor_exo_cons
# names_cor_measurement_error <- names$names_cor_measurement_error
### Compute differences ======================================================
if(inherits(.object, "cSEMResults_2ndorder")) {
path_estimates <- lapply(x, function(y) {y$Second_stage$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Loading_estimates,
y$Second_stage$Estimates$Loading_estimates)
})
weight_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Weight_estimates,
y$Second_stage$Estimates$Weight_estimates)
})
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Second_stage$Estimates$Exo_construct_correlation})
} else { # no second-order model
path_estimates <- lapply(x, function(y) {y$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {y$Estimates$Loading_estimates})
weight_estimates <- lapply(x, function(y) {y$Estimates$Weight_estimates})
# all exogenous construct correlations
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Estimates$Exo_construct_correlation
})
}
## Select
path_estimates <- lapply(path_estimates, function(y) {
y1 <- y[y$Name %in% names_path, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_path
}
y1
})
loading_estimates <- lapply(loading_estimates, function(y) {
y1 <- y[y$Name %in% names_loadings, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_loadings
}
y1
})
weight_estimates <- lapply(weight_estimates, function(y) {
y1 <- y[y$Name %in% names_weights, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_weights
}
y1
})
cor_cons_exo_estimates <- lapply(cor_cons_exo_estimates, function(y) {
y1 <- y[y$Name %in% names_cor_exo_cons, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_cor_exo_cons
}
y1
})
# Do not change the order
list(path_estimates = path_estimates,
loading_estimates = loading_estimates,
weight_estimates = weight_estimates,
cor_exo_cons_estimates = cor_cons_exo_estimates)
}
#' Internal: get structured cSEMTestMGD results
#'
#' Convenience function to summarize the results of all tests resulting from a
#' call to [testMGD()] in a user-friendly way.
#'
#' @usage structureTestMGDDecisions(.object)
#'
#' @inheritParams csem_arguments
#'
#' @return A data.frame.
#'
#' @keywords internal
structureTestMGDDecisions <- function(.object){
## Install dplyr, tidyr, and rlang if not already installed
if (!requireNamespace("dplyr", quietly = TRUE)) {
stop2(
"Package `dplyr` required. Use `install.packages(\"dplyr\")` and rerun.")
}
if (!requireNamespace("tidyr", quietly = TRUE)) {
stop2(
"Package `tidyr` (> 1.0.0) required. Use `install.packages(\"tidyr\")` and rerun.")
}
if (!requireNamespace("rlang", quietly = TRUE)) {
stop2(
"Package `rlang` required. Use `install.packages(\"rlang\")` and rerun.")
}
# Check if of class cSEMTestMGD
if(!inherits(.object, "cSEMTestMGD")){
stop2("The following error occured in the `structureTestMGDDecisions()` function:\n",
"Object must be of class `cSEMResults`")
}
# Helper function ------------------------------------------------------------
# Change CI alpha equivlanet (1-CI)
changeCItoAlphaEquivalent <- function(.names){
# .names <- c("99%","95%","90%")
tmpNames <- as.numeric(gsub("%", "", .names))
tmpNames <- 100-tmpNames
tmpNames <- paste0(tmpNames, "%")
return(tmpNames)
}
out <- dplyr::tibble()
# Klesel ---------------------------------------------------------------------
# Test only provides an overall decision
if("Klesel" %in% names(.object)){
klesel <- .object$Klesel$Decision %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "alpha") %>%
# longer format
tidyr::pivot_longer(cols = c(.data$dG, .data$dL),
names_to = "Distance_metric",
values_to = "Decision") %>%
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Klesel", Comparison = "overall")
out <- dplyr::bind_rows(out, klesel)
}
# Chin ----------------------------------------------------------------------
if("Chin" %in% names(.object)){
# overall decision
chin1 <- .object$Chin$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Chin", Comparison = "overall")
chin2 <- .object$Chin$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
# Check all comparisons
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison), all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Chin")
out <- dplyr::bind_rows(out, chin1, chin2)
}
# Sarstedt ------------------------------------------------------------------
if("Sarstedt" %in% names(.object)){
# overall result
sarstedt1 <- .object$Sarstedt$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Sarstedt", Comparison = "overall")
# decision based on single paths
sarstedt2 <- .object$Sarstedt$Decision %>%
purrr::modify_depth(2, dplyr::bind_rows) %>%
purrr::map(dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Sarstedt")
out <- dplyr::bind_rows(out, sarstedt1, sarstedt2)
}
# Keil ----------------------------------------------------------------------
if("Keil" %in% names(.object)){
# overall decision
keil1 <- .object$Keil$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Keil", Comparison = "overall")
# decision based on single paths
keil2 <- .object$Keil$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison),all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Keil")
out <- dplyr::bind_rows(out, keil1, keil2)
}
# Nitzl ---------------------------------------------------------------------
if("Nitzl" %in% names(.object)){
# overall decision
nitzl1 <- .object$Nitzl$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Nitzl", Comparison = "overall")
# decision based on single paths
nitzl2 <- .object$Nitzl$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison),all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Nitzl")
out <- dplyr::bind_rows(out, nitzl1, nitzl2)
}
# Henseler ------------------------------------------------------------------
if("Henseler" %in% names(.object)){
# overall decision
henseler1 <- .object$Henseler$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Henseler", Comparison = "overall")
# decision based on single paths
henseler2 <- .object$Henseler$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison),all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Henseler")
out <- dplyr::bind_rows(out, henseler1, henseler2)
}
# CI para -------------------------------------------------------------------
if("CI_para" %in% names(.object)){
# overall decision
CIpara1 <- .object$CI_para$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_longer(cols = dplyr::contains("CI"),
names_to = "CI_type",
values_to = "Decision") %>%
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivlanet to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_para", Comparison = "overall")
# decision based on single paths
CIpara2 <- .object$CI_para$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(3, ~ dplyr::select(., .data$Name, .data$Decision)) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "CI_type") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "Comparison") %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_wider(names_from = .data$Name, values_from = .data$Decision) %>%
dplyr::group_by(.data$alpha, .data$CI_type) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison), all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$CI_type, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivalent to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_para")
out <- dplyr::bind_rows(out, CIpara1, CIpara2)
}
# CI overlap ----------------------------------------------------------------
if("CI_overlap" %in% names(.object)){
# overall decision
CIoverlap1 <- .object$CI_overlap$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_longer(cols = dplyr::contains("CI"),
names_to = "CI_type",
values_to = "Decision") %>%
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivalent to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_overlap", Comparison = "overall")
# decision based on single paths
CIoverlap2 <- .object$CI_overlap$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(3, ~ dplyr::select(., .data$Name, .data$Decision)) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "CI_type") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "Comparison") %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_wider(names_from = .data$Name, values_from = .data$Decision) %>%
dplyr::group_by(.data$alpha, .data$CI_type) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison), all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$CI_type, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivalent to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_overlap")
out <- dplyr::bind_rows(out, CIoverlap1, CIoverlap2)
}
# Summarize all results ------------------------------------------------------
out <- dplyr::select(out, .data$Test, .data$Comparison, dplyr::contains("%"),
.data$`p-value_correction`, .data$CI_type,
.data$Distance_metric)
return(out)
}
M-E-Steiner/cSEM documentation built on March 18, 2024, 12:18 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions structureTestMGDDecisions getRelevantParameters calculatePr calculateFR adjustAlpha calculateParameterDifference getParameterNames calculateDistance
Documented in adjustAlpha calculateDistance calculateFR calculateParameterDifference calculatePr getParameterNames getRelevantParameters structureTestMGDDecisions
#' Internal: Matrix difference
#'
#' Calculates the average of the differences between all possible pairs of
#' (symmetric) matrices in a list using a given distance measure.
#'
#' `.matrices` must be a list of at least two matrices. If more than two matrices
#' are supplied the arithmetic mean of the differences between all possible pairs of
#' (symmetric) matrices in a list is computed. Mathematically this is
#' n chose 2. Hence, supplying a large number of matrices will
#' become computationally challenging.
#'
#' Currently two distance measures are supported:
#' \describe{
#' \item{`geodesic`}{(Default) The geodesic distance.}
#' \item{`squared_euclidean`}{The squared Euclidean distance}
#' }
#'
#' @usage calculateDistance(
#' .matrices = NULL,
#' .distance = args_default()$.distance
#' )
#'
#' @inheritParams csem_arguments
#'
#' @return A numeric vector of length one containing the (arithmetic) mean of
#' the differences between all possible pairs of matrices supplied via `.matrices`.
#'
#' @keywords internal
calculateDistance <- function(
.matrices = NULL,
.distance = args_default()$.distance
){
### Checks and errors ========================================================
## Check if list and at least of length 2
if(!is.list(.matrices) && length(.matrices) < 2) {
stop2("`.matrices` must be a list of at least length two.")
}
## Check if matrices are all symmetric
if(!all(sapply(.matrices, matrixcalc::is.symmetric.matrix))) {
stop2("All matrices in `.matrices` must be symmetric.")
}
### Calculation ==============================================================
## Combine all matrices into lists of two
temp <- utils::combn(.matrices, 2, simplify = FALSE)
## Compute the distance measure for each group combination
distances <- lapply(temp, function(x) {
switch (.distance,
"geodesic" = {calculateDG(.matrix1 = x[[1]], .matrix2 = x[[2]])},
"squared_euclidean" = {calculateDL(.matrix1 = x[[1]], .matrix2 = x[[2]])}
)
})
## Compute the mean. Since dG and dL are defined as the mean of difference between
## two matrices (1/2 is part of the definition). This must be corrected if
## more than two matrices are to be compared.
if(length(.matrices) > 2) {
out <- mean(2 * unlist(distances))
} else {
out <- unlist(distances)
}
return(out)
}
#' Internal: Parameter names
#'
#' Based on a model in [lavaan model syntax][lavaan::model.syntax], returns the
#' names of the parameters of the structural
#' model, the measurement/composite model and the weight relationship. Used
#' by [testMGD()] to extract the names of the parameters to compare across groups
#' according to the test proposed by \insertCite{Chin2010;textual}{cSEM}.
#'
#' @usage getParameterNames(
#' .object = args_default()$.object,
#' .model = args_default()$.model
#' )
#'
#' @inheritParams csem_arguments
#' @param .model A model in [lavaan model syntax][lavaan::model.syntax] indicating which
#' parameters (i.e, path (`~`), loadings (`=~`), or weights (`<~`)) should be
#' compared across groups. Defaults to `NULL` in which case all parameters of the model
#' are compared.
#'
#' @return A list with elements `names_path`, `names_loadings`, and `names_weights`
#' containing the names of the structural parameters, the loadings,
#' and the weight to compare across groups.
#'
#' @references
#' \insertAllCited{}
#'
#' @keywords internal
getParameterNames <- function(
.object = args_default()$.object,
.model = args_default()$.model
){
## Summarize
x <- summarize(.object)
if(inherits(.object, "cSEMResults_2ndorder")) {
x12 <- x[[1]]$First_stage$Information
x22 <- x[[1]]$Second_stage$Information
# Extract different types of constructs
construct_type <- c(x12$Model$construct_type, x22$Model$construct_type)
construct_type <- construct_type[unique(names(construct_type))]
measurement_org <- x22$Arguments_original$.model$measurement
names_path_org <- x[[1]]$Second_stage$Estimates$Path_estimates$Name
indicators <- x12$Model_original$indicators
cons_exo <- x12$Model_original$cons_exo
names_cor_exo_cons_org <- x[[1]]$Second_stage$Estimates$Exo_construct_correlation$Name
names_cor_indicator_org <- c(x[[1]]$First_stage$Estimates$Indicator_correlation$Name,
x[[1]]$Second_stage$Estimates$Indicator_correlation$Name)
if(nrow(x[[1]]$First_stage$Estimates$Residual_correlation)!=0){
temp <- x[[1]]$First_stage$Estimates$Residual_correlation$Name
}else{
temp <- NULL
}
if(nrow(x[[1]]$Second_stage$Estimates$Residual_correlation)!=0){
temp1 <- x[[1]]$Second_stage$Estimates$Residual_correlation$Name
}else{
temp1 <- NULL
}
names_cor_measurement_error_org <- c(temp,temp1)
} else {#if no second-order model
x22 <- x[[1]]$Information
# Extract different types of constructs
construct_type <- x22$Model$construct_type
measurement_org <- x22$Model$measurement
names_path_org <- x[[1]]$Estimates$Path_estimates$Name
indicators <- x22$Model$indicators
cons_exo <- x22$Model$cons_exo
names_cor_exo_cons_org <- x[[1]]$Estimates$Exo_construct_correlation$Name
names_cor_indicator_org <- x[[1]]$Estimates$Indicator_correlation$Name
if(nrow(x[[1]]$Estimates$Residual_correlation)!=0){
names_cor_measurement_error_org <- x[[1]]$Estimates$Residual_correlation$Name
}else{
names_cor_measurement_error_org<-NULL
}
}
# FOR FUTURE:
# Perhaps it is better to transfer the following code to the testMGD function directly.
# As getParameterNames should just return the names of a model and not internally select a model?
# Parse model that indicates which parameters should be compared.
# If no model indicating the comparisons is provided, all parameters are compared.
# Correlation among the measurement errors and or indicators cannot be compared yet across groups as
# theta is not part of the csem output
# FOR FUTURE: Add all model parameters to the comparison if no model is specified
if(is.null(.model)) {
if(inherits(.object, "cSEMResults_2ndorder")) {
model_compare <- x22$Arguments_original$.model
# add the correlations among the exogenous constructs
} else {
model_compare <- x22$Model
# add the correlations among the exogenous constructs
# add to model_cor_specified and expand
# Attention if the some of the variables are already in that matrix
}
} else {
model_compare <- parseModel(.model, .check_errors = FALSE)
}
# Check whether the constructs specified in the comparison are equal
# to the constructs in the original model
if(!all(rownames(model_compare$structural) %in% rownames(measurement_org))){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one construct appears in the comparison model but",
" not in the original model."
)
}
# Check whether the construct type specified in the comparison are equal
# to the construct types in the original model
construct_type_compare <- model_compare$construct_type[!is.na(model_compare$construct_type)]
if(!all(construct_type_compare == construct_type[names(construct_type_compare)])){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one construct's type in the comparison model differs ",
"from the original model.")
}
# Check wehther the indicators used in the comparison model also appear
# in the original model
if(!all(colnames(model_compare$measurement) %in% colnames(measurement_org))){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one indicator appears in the comparison model and not in ",
" the original model.")
}
# Check whether the indicators used are correctly assigned, i.e.,
# as in the original model
lapply(rownames(model_compare$measurement), function(x) {
names_ind_compare <- colnames(model_compare$measurement[x,,drop=FALSE][,which(model_compare$measurement[x,,drop=FALSE]==1),drop=FALSE])
names_ind_org <- colnames(measurement_org[x,,drop=FALSE][,which(measurement_org[x,,drop=FALSE]==1),drop=FALSE])
if(!all(names_ind_compare %in% names_ind_org)){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one indicator is not correctly assigned in the comparison model.")
}
})
### Extract names ============================================================
# Extract names of the path to be tested
temp <- outer(rownames(model_compare$structural), colnames(model_compare$structural),
FUN = function(x, y) paste(x, y, sep = " ~ "))
names_path_compare <- t(temp)[t(model_compare$structural) != 0]
if(length(names_path_compare) == 0) {
names_path_compare <- NULL
}
# Check whether the path to compare occur also in the original mode
if(!all(names_path_compare %in% names_path_org)){
stop2(
"The following error occured in the `getParameterNames()` function:\n",
"At least one of the paths specified for comparison does not ",
" appear in the original model.")
}
# all correlated variables, i.e., that have been specified with ~~
vars_correlated_compare <- rownames(model_compare$cor_specified)
# Select only those that are indicators
ind_correlated_compare <- intersect(vars_correlated_compare, indicators)
# Extract names of the loadings and measurement error to be tested
names_row <- rownames(model_compare$measurement)
## Select only concepts modeled as common factors. Reorder to be have the
## same order as names_row.
i <- intersect(names(which(construct_type == "Common factor")), names_row)
i <- i[match(names_row, i)]
i <- i[!is.na(i)]
if(length(i) != 0) {
temp <- rep(rownames(model_compare$measurement[i, , drop = FALSE]),
times = rowSums(model_compare$measurement[i, , drop = FALSE]))
if(length(temp)!=0){
# Loadings
temp_n <- model_compare$measurement[i, colSums(model_compare$measurement[i, , drop = FALSE]) != 0, drop = FALSE]
names_loadings_compare <- paste0(temp, " =~ ", colnames(temp_n))
# Measurement error
measurement_error_correlated <- intersect(ind_correlated_compare,colnames(temp_n))
cor_measurement_error <- model_compare$cor_specified[ measurement_error_correlated, measurement_error_correlated]
index <- which(cor_measurement_error == 1, arr.ind = TRUE)
names_correlated_measurement_error_compare <- index
if(nrow(names_correlated_measurement_error_compare) ==0 ){
names_correlated_measurement_error_compare <- NULL
}else{
names_correlated_measurement_error_compare <- paste(rownames(cor_measurement_error)[index[,"row"]],
" ~~ ",
colnames(cor_measurement_error)[index[,"col"]],
sep="")
# Select those that have the same order as in the original model
names_correlated_measurement_error_compare <- names_correlated_measurement_error_compare[names_correlated_measurement_error_compare %in% names_cor_measurement_error_org]
}
}else{
names_loadings_compare <- NULL
names_correlated_measurement_error_compare <- NULL
}
} else {
names_loadings_compare <- NULL
names_correlated_measurement_error_compare <- NULL
}
## Extract names of the weights and indicator correlations (composites) to be tested
# Select only concepts modeled as composite. Reorder that it has the
# same order as names_row.
## Indicator correlation, i.e., variables that belong to a composite
i <- intersect(names(which(construct_type == "Composite")), names_row)
i <- i[match(names_row, i)]
i <- i[!is.na(i)]
if(length(i) != 0) {
temp <- rep(rownames(model_compare$measurement[i, , drop = FALSE]),
times = rowSums(model_compare$measurement[i, , drop = FALSE]))
if(length(temp)!=0){
# weights
temp_n <- model_compare$measurement[i, colSums(model_compare$measurement[i, , drop = FALSE]) != 0, drop = FALSE]
names_weights_compare <- paste0(temp, " <~ ", colnames(temp_n))
# indicator correlations
ind_connected_to_composite <- colnames(temp_n)
indicator_correlated <- intersect(ind_correlated_compare,ind_connected_to_composite)
cor_indicator <- model_compare$cor_specified[indicator_correlated, indicator_correlated]
index <- which(cor_indicator == 1, arr.ind = TRUE)
names_correlated_indicator_compare <- index
# In case that no indicator correlations are compared set it to NULL
if(nrow(names_correlated_indicator_compare) ==0 ){
names_correlated_indicator_compare <- NULL
}else{
names_correlated_indicator_compare <- paste(rownames(cor_indicator)[index[,"row"]],
" ~~ ",
colnames(cor_indicator)[index[,"col"]],
sep="")
# Select those that have the same order as in the original model
names_correlated_indicator_compare <- names_correlated_indicator_compare[names_correlated_indicator_compare %in% names_cor_indicator_org]
}
}else{
names_weights_compare <- NULL
names_correlated_indicator_compare <- NULL
}
} else {
names_weights_compare <- NULL
names_correlated_indicator_compare <- NULL
}
## Exogenous construct correlations
cons_exo_correlated_compare <- intersect(vars_correlated_compare, cons_exo)
# Consider only exogenous constructs
cor_cons_exo <- model_compare$cor_specified[cons_exo_correlated_compare,cons_exo_correlated_compare]
# Which are correlated
index <- which(cor_cons_exo == 1,arr.ind = TRUE)
names_correlated_exo_cons_compare <- index
# In case that no measurement error correlations are compared set it to NULL
if(nrow(names_correlated_exo_cons_compare) == 0){
names_correlated_exo_cons_compare <- NULL
}else{
names_correlated_exo_cons_compare <-paste(rownames(cor_cons_exo)[index[,"row"]],
" ~~ ",
colnames(cor_cons_exo)[index[,"col"]], sep="")
# Select those that have the same order as in the original model
names_correlated_exo_cons_compare <- names_correlated_exo_cons_compare[names_correlated_exo_cons_compare %in% names_cor_exo_cons_org]
}
## Return as list
out <- list(
"names_path" = names_path_compare,
"names_weights" = names_weights_compare,
"names_loadings" = names_loadings_compare,
"names_cor_exo_cons" = names_correlated_exo_cons_compare,
"names_cor_indicator" = names_correlated_indicator_compare,
"names_cor_measurement_error" = names_correlated_measurement_error_compare
)
return(out)
}
#' Internal: Parameter differences across groups
#'
#' Calculate the difference between one or more parameter estimates across
#' all possible pairs of groups (data sets) in `.object`.
#'
#' @usage calculateParameterDifference(
#' .object = args_default()$.object,
#' .model = args_default()$.model
#' )
#'
#' @inheritParams csem_arguments
#' @param .model A model in [lavaan model syntax][lavaan::model.syntax] indicating which
#' parameters (i.e., path (`~`), loadings (`=~`), or weights (`<~`)) should be
#' compared across groups. Defaults to `NULL` in which case all parameters of the model
#' are compared.
#'
#' @return A list of length equal to the number of possible pairs of
#' groups in `.object` (mathematically, this is n choose 2, i.e., 3 if there are three
#' groups and 6 if there are 4 groups). Each list elements is itself a list of
#' three. The first list element contains
#' the difference between parameter estimates of the structural model, the second
#' list element the difference between estimated loadings, and the third
#' the difference between estimated weights.
#'
#' @keywords internal
calculateParameterDifference <- function(
.object = args_default()$.object,
.model = args_default()$.model
){
## Summarize
x <- summarize(.object)
## Get names to compare
names <- getParameterNames(.object, .model = .model)
names_path <- names$names_path
names_loadings <- names$names_loadings
names_weights <- names$names_weights
names_cor_exo_cons <- names$names_cor_exo_cons
# names_cor_measurement_error <- names$names_cor_measurement_error
### Compute differences ======================================================
if(inherits(.object, "cSEMResults_2ndorder")) {
path_estimates <- lapply(x, function(y) {y$Second_stage$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Loading_estimates,
y$Second_stage$Estimates$Loading_estimates)
})
weight_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Weight_estimates,
y$Second_stage$Estimates$Weight_estimates)
})
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Second_stage$Estimates$Exo_construct_correlation})
} else { # no second-order model
path_estimates <- lapply(x, function(y) {y$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {y$Estimates$Loading_estimates})
weight_estimates <- lapply(x, function(y) {y$Estimates$Weight_estimates})
# all exogenous construct correlations
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Estimates$Exo_construct_correlation
})
}
## Select
path_estimates <- lapply(path_estimates, function(y) {
y1 <- y[y$Name %in% names_path, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_path
}
y1
})
loading_estimates <- lapply(loading_estimates, function(y) {
y1 <- y[y$Name %in% names_loadings, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_loadings
}
y1
})
weight_estimates <- lapply(weight_estimates, function(y) {
y1 <- y[y$Name %in% names_weights, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_weights
}
y1
})
cor_cons_exo_estimates <- lapply(cor_cons_exo_estimates, function(y) {
y1 <- y[y$Name %in% names_cor_exo_cons, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_cor_exo_cons
}
y1
})
## Path model
temp <- utils::combn(path_estimates, 2, simplify = FALSE)
diff_path <- lapply(temp, function(y) y[[1]] - y[[2]])
names(diff_path) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
## Loadings
temp <- utils::combn(loading_estimates, 2, simplify = FALSE)
diff_loadings <- lapply(temp, function(y) y[[1]] - y[[2]])
names(diff_loadings) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
## Weights
temp <- utils::combn(weight_estimates, 2, simplify = FALSE)
diff_weights <- lapply(temp, function(y) y[[1]] - y[[2]])
names(diff_weights) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
# Exogenous construct correlations
temp <- utils::combn(cor_cons_exo_estimates, 2, simplify = FALSE)
cor_cons_exo <- lapply(temp, function(y) y[[1]] - y[[2]])
names(cor_cons_exo) <- sapply(temp, function(x) paste0(names(x)[1], '_', names(x)[2]))
# merge list together
out <- mapply(function(w,x,y,z) c(w,x,y,z),
w = diff_path,
x = diff_loadings,
y = diff_weights,
z = cor_cons_exo,
SIMPLIFY = FALSE)
return(out)
}
#' Internal: Multiple testing correction
#'
#' Adjust a given significance level `.alpha` to accommodate multiple testing.
#' The following corrections are implemented:
#' \describe{
#' \item{`none`}{(Default) No correction is done.}
#' \item{`bonferroni`}{A Bonferroni correction is done, i.e., alpha is divided by the
#' number of comparisons `.nr_comparisons`.}
#' }
#'
#' @usage adjustAlpha(
#' .alpha = args_default()$.alpha,
#' .approach_alpha_adjust = args_default()$.approach_alpha_adjust,
#' .nr_comparisons = args_default()$.nr_comparisons
#' )
#'
#' @inheritParams csem_arguments
#'
#' @return A vector of (possibly adjusted) significance levels.
#'
#' @keywords internal
adjustAlpha <- function(
.alpha = args_default()$.alpha,
.approach_alpha_adjust = args_default()$.approach_alpha_adjust,
.nr_comparisons = args_default()$.nr_comparisons
){
if(.approach_alpha_adjust == 'none'){
return(.alpha)
}
if(.approach_alpha_adjust == 'bonferroni'){
return(.alpha/.nr_comparisons)
}
}
#' Internal: ANOVA F-test statistic
#'
#' Calculate the ANOVA F-test statistic suggested by
#' \insertCite{Sarstedt2011;textual}{cSEM} in the OTG testing procedure.
#'
#' @usage calculateFR(.resample_sarstedt)
#'
#' @inheritParams csem_arguments
#'
#' @return A named scalar, the test statistic of the ANOVA F-test
#'
#' @references
#' \insertAllCited{}
#'
#' @keywords internal
calculateFR <- function(.resample_sarstedt) {
# In Sarstedt et al. (2011), it is assumed that you have the same number of bootstrap estimates per group.
# Consequently, the total number of estimates equals B*G.
# Since, there might be different numbers of bootstrap estimates per group, we adjusted the calculation of the
# F test statistic.
names_param <- colnames(.resample_sarstedt)[-ncol(.resample_sarstedt)]
# Number of grousp
G <- length(unique(.resample_sarstedt[, "group_id"]))
# Total number of estimates
N <- nrow(.resample_sarstedt)
# Order the estimates by the grouping variable
x <- .resample_sarstedt[order(.resample_sarstedt[, "group_id"]), ]
grand_means <- colMeans(x[, names_param, drop = FALSE])
mean_group <- aggregate(x[, names_param, drop = FALSE], by = list(x[, "group_id"]), FUN = mean)
group_mean_matrix <- apply(mean_group[, -1, drop = FALSE], 2, function(y) rep(y, times = table(x[, "group_id"])))
SS_between <- rowSums((t(group_mean_matrix) - grand_means)^2)
MSA_between <- SS_between/ (G -1)
SS_within <- colSums((x[, names_param, drop = FALSE] - group_mean_matrix)^2)
MSA_within <- SS_within/(N - G)
# Note: in the original paper (B -1) was used. This is incorrect. It should
# be B - G (G = the number of groups). See the F-statistic derived from ANOVA
# for comparision.
MSA_between/MSA_within
}
#' Internal: Calculation of the CDF used in Henseler et al. (2009)
#'
#' Calculates the probability that theta^1 is smaller than or equal to theta^2.
#' See Equation (6) in \insertCite{Sarstedt2011;textual}{cSEM}.
#'
#' @usage calculatePr(.resample_centered = NULL, .parameters_to_compare = NULL)
#'
#' @inheritParams csem_arguments
#'
#' @return A named vector
#'
#' @references
#' \insertAllCited{}
#'
#' @keywords internal
calculatePr <- function(
.resample_centered = NULL,
.parameters_to_compare = NULL
){
# Remove names from .parameters_to_compare
names(.parameters_to_compare) <- NULL
group1 <- .resample_centered[[1]][,.parameters_to_compare,drop=FALSE]
group2 <- .resample_centered[[2]][,.parameters_to_compare,drop=FALSE]
ret <- sapply(.parameters_to_compare, function(x){
# Matrix where each element of the first vector is substracted from the complete second vector
# theta^1 <= theta^2
temp12 <- outer(group2[,x],group1[,x],"-")
g1leqg2<-mean((1+sign(temp12))/2)
# theta^2 <= theta^1 is equal to 1 - g1leqg2, therefore we do not need the following calculations
# temp21 <- outer(group1[,x],group2[,x],"-")
# g2leqg1<-mean((1+sign(temp21))/2)
# # The problem is that Henseler's approach is a one-sided test.
# # Therefore, the p-value is flipped if it larger than 0.5
# # Return
# if(g1leqg2<0.5){
# g1leqg2
# }else{
# 1-g1leqg2
# }
})
names(ret) <- .parameters_to_compare
ret
}
#' Internal: Extract relevant parameters from several cSEMResults_multi
#'
#' Extract the relevant parameters from a cSEMResult_multi object in `.object`.
#'
#' @usage getRelevantParameters(
#' .object = args_default()$.object,
#' .model = args_default()$.model
#' )
#'
#' @inheritParams csem_arguments
#' @param .model A model in [lavaan model syntax][lavaan::model.syntax] indicating which
#' parameters (i.e., path (`~`), loadings (`=~`), or weights (`<~`)) should be
#' compared across groups. Defaults to `NULL` in which case all parameters of the model
#' are compared.
#'
#' @return A list of length equal to the number of groups in `.object`.
#' Each list element is itself a list of three. The first list element contains
#' the relevant parameter estimates of the structural model, the second
#' list element the relevant estimated loadings, and the third
#' the relevant estimated weights.
#'
#' @keywords internal
getRelevantParameters <- function(
.object = args_default()$.object,
.model = args_default()$.model
){
## Summarize
x <- summarize(.object)
## Get names to compare
names <- getParameterNames(.object, .model = .model)
names_path <- names$names_path
names_loadings <- names$names_loadings
names_weights <- names$names_weights
names_cor_exo_cons <- names$names_cor_exo_cons
# names_cor_measurement_error <- names$names_cor_measurement_error
### Compute differences ======================================================
if(inherits(.object, "cSEMResults_2ndorder")) {
path_estimates <- lapply(x, function(y) {y$Second_stage$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Loading_estimates,
y$Second_stage$Estimates$Loading_estimates)
})
weight_estimates <- lapply(x, function(y) {
rbind(y$First_stage$Estimates$Weight_estimates,
y$Second_stage$Estimates$Weight_estimates)
})
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Second_stage$Estimates$Exo_construct_correlation})
} else { # no second-order model
path_estimates <- lapply(x, function(y) {y$Estimates$Path_estimates})
loading_estimates <- lapply(x, function(y) {y$Estimates$Loading_estimates})
weight_estimates <- lapply(x, function(y) {y$Estimates$Weight_estimates})
# all exogenous construct correlations
cor_cons_exo_estimates <- lapply(x, function(y) {
y$Estimates$Exo_construct_correlation
})
}
## Select
path_estimates <- lapply(path_estimates, function(y) {
y1 <- y[y$Name %in% names_path, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_path
}
y1
})
loading_estimates <- lapply(loading_estimates, function(y) {
y1 <- y[y$Name %in% names_loadings, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_loadings
}
y1
})
weight_estimates <- lapply(weight_estimates, function(y) {
y1 <- y[y$Name %in% names_weights, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_weights
}
y1
})
cor_cons_exo_estimates <- lapply(cor_cons_exo_estimates, function(y) {
y1 <- y[y$Name %in% names_cor_exo_cons, "Estimate"]
if(length(y1) != 0) {
names(y1) <- names_cor_exo_cons
}
y1
})
# Do not change the order
list(path_estimates = path_estimates,
loading_estimates = loading_estimates,
weight_estimates = weight_estimates,
cor_exo_cons_estimates = cor_cons_exo_estimates)
}
#' Internal: get structured cSEMTestMGD results
#'
#' Convenience function to summarize the results of all tests resulting from a
#' call to [testMGD()] in a user-friendly way.
#'
#' @usage structureTestMGDDecisions(.object)
#'
#' @inheritParams csem_arguments
#'
#' @return A data.frame.
#'
#' @keywords internal
structureTestMGDDecisions <- function(.object){
## Install dplyr, tidyr, and rlang if not already installed
if (!requireNamespace("dplyr", quietly = TRUE)) {
stop2(
"Package `dplyr` required. Use `install.packages(\"dplyr\")` and rerun.")
}
if (!requireNamespace("tidyr", quietly = TRUE)) {
stop2(
"Package `tidyr` (> 1.0.0) required. Use `install.packages(\"tidyr\")` and rerun.")
}
if (!requireNamespace("rlang", quietly = TRUE)) {
stop2(
"Package `rlang` required. Use `install.packages(\"rlang\")` and rerun.")
}
# Check if of class cSEMTestMGD
if(!inherits(.object, "cSEMTestMGD")){
stop2("The following error occured in the `structureTestMGDDecisions()` function:\n",
"Object must be of class `cSEMResults`")
}
# Helper function ------------------------------------------------------------
# Change CI alpha equivlanet (1-CI)
changeCItoAlphaEquivalent <- function(.names){
# .names <- c("99%","95%","90%")
tmpNames <- as.numeric(gsub("%", "", .names))
tmpNames <- 100-tmpNames
tmpNames <- paste0(tmpNames, "%")
return(tmpNames)
}
out <- dplyr::tibble()
# Klesel ---------------------------------------------------------------------
# Test only provides an overall decision
if("Klesel" %in% names(.object)){
klesel <- .object$Klesel$Decision %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "alpha") %>%
# longer format
tidyr::pivot_longer(cols = c(.data$dG, .data$dL),
names_to = "Distance_metric",
values_to = "Decision") %>%
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Klesel", Comparison = "overall")
out <- dplyr::bind_rows(out, klesel)
}
# Chin ----------------------------------------------------------------------
if("Chin" %in% names(.object)){
# overall decision
chin1 <- .object$Chin$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Chin", Comparison = "overall")
chin2 <- .object$Chin$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
# Check all comparisons
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison), all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Chin")
out <- dplyr::bind_rows(out, chin1, chin2)
}
# Sarstedt ------------------------------------------------------------------
if("Sarstedt" %in% names(.object)){
# overall result
sarstedt1 <- .object$Sarstedt$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Sarstedt", Comparison = "overall")
# decision based on single paths
sarstedt2 <- .object$Sarstedt$Decision %>%
purrr::modify_depth(2, dplyr::bind_rows) %>%
purrr::map(dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Sarstedt")
out <- dplyr::bind_rows(out, sarstedt1, sarstedt2)
}
# Keil ----------------------------------------------------------------------
if("Keil" %in% names(.object)){
# overall decision
keil1 <- .object$Keil$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Keil", Comparison = "overall")
# decision based on single paths
keil2 <- .object$Keil$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison),all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Keil")
out <- dplyr::bind_rows(out, keil1, keil2)
}
# Nitzl ---------------------------------------------------------------------
if("Nitzl" %in% names(.object)){
# overall decision
nitzl1 <- .object$Nitzl$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Nitzl", Comparison = "overall")
# decision based on single paths
nitzl2 <- .object$Nitzl$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison),all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Nitzl")
out <- dplyr::bind_rows(out, nitzl1, nitzl2)
}
# Henseler ------------------------------------------------------------------
if("Henseler" %in% names(.object)){
# overall decision
henseler1 <- .object$Henseler$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::mutate(Test = "Henseler", Comparison = "overall")
# decision based on single paths
henseler2 <- .object$Henseler$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "Comparison") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "alpha") %>%
dplyr::bind_rows(.id = "p-value_correction") %>%
dplyr::group_by(.data$`p-value_correction`, .data$alpha) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison),all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$`p-value_correction`, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
dplyr::mutate(Test = "Henseler")
out <- dplyr::bind_rows(out, henseler1, henseler2)
}
# CI para -------------------------------------------------------------------
if("CI_para" %in% names(.object)){
# overall decision
CIpara1 <- .object$CI_para$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_longer(cols = dplyr::contains("CI"),
names_to = "CI_type",
values_to = "Decision") %>%
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivlanet to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_para", Comparison = "overall")
# decision based on single paths
CIpara2 <- .object$CI_para$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(3, ~ dplyr::select(., .data$Name, .data$Decision)) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "CI_type") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "Comparison") %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_wider(names_from = .data$Name, values_from = .data$Decision) %>%
dplyr::group_by(.data$alpha, .data$CI_type) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison), all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$CI_type, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivalent to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_para")
out <- dplyr::bind_rows(out, CIpara1, CIpara2)
}
# CI overlap ----------------------------------------------------------------
if("CI_overlap" %in% names(.object)){
# overall decision
CIoverlap1 <- .object$CI_overlap$Decision_overall %>%
purrr::map(dplyr::bind_rows) %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_longer(cols = dplyr::contains("CI"),
names_to = "CI_type",
values_to = "Decision") %>%
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivalent to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_overlap", Comparison = "overall")
# decision based on single paths
CIoverlap2 <- .object$CI_overlap$Decision %>%
purrr::modify_depth(3, dplyr::bind_rows) %>%
purrr::modify_depth(3, ~ dplyr::select(., .data$Name, .data$Decision)) %>%
purrr::modify_depth(2, dplyr::bind_rows, .id = "CI_type") %>%
purrr::modify_depth(1, dplyr::bind_rows, .id = "Comparison") %>%
dplyr::bind_rows(.id = "alpha") %>%
tidyr::pivot_wider(names_from = .data$Name, values_from = .data$Decision) %>%
dplyr::group_by(.data$alpha, .data$CI_type) %>%
dplyr::summarise_at(dplyr::vars(-.data$Comparison), all) %>%
dplyr::ungroup() %>%
# put results into format
tidyr::pivot_longer(cols = (-c(.data$CI_type, .data$alpha)),
names_to = "Comparison",
values_to = "Decision") %>%
# put alphas in cols
tidyr::pivot_wider(names_from = .data$alpha, values_from = .data$Decision) %>%
# Change % as equivalent to alpha
dplyr::rename_at(dplyr::vars(dplyr::contains("%")), changeCItoAlphaEquivalent) %>%
dplyr::mutate(Test = "CI_overlap")
out <- dplyr::bind_rows(out, CIoverlap1, CIoverlap2)
}
# Summarize all results ------------------------------------------------------
out <- dplyr::select(out, .data$Test, .data$Comparison, dplyr::contains("%"),
.data$`p-value_correction`, .data$CI_type,
.data$Distance_metric)
return(out)
}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.