From 59eb34021f489275f407ce2634ee6340ba366747 Mon Sep 17 00:00:00 2001 From: Alexander Christensen Date: Sat, 9 Mar 2024 09:16:33 -0600 Subject: [PATCH] updates to `net.loads` o correction in `invariance` --- R/invariance.R | 448 +++++++++++++++++++++++----------------------- R/net.loads.R | 42 +---- man/invariance.Rd | 28 +-- 3 files changed, 246 insertions(+), 272 deletions(-) diff --git a/R/invariance.R b/R/invariance.R index 2cee781c..2242b900 100644 --- a/R/invariance.R +++ b/R/invariance.R @@ -4,19 +4,19 @@ #' on all data (across groups) first. After configural variance is established, #' then metric invariance is tested using the community structure that established #' configural invariance (see \strong{Details} for more information on this process) -#' +#' #' @param data Matrix or data frame. #' Should consist only of variables to be used in the analysis -#' +#' #' @param groups Numeric or character vector (length = \code{ncol(data)}). #' Group membership corresponding to each case in data -#' +#' #' @param structure Numeric or character vector (length = \code{ncol(data)}). #' A vector representing the structure (numbers or labels for each item). #' Can be theoretical factors or the structure detected by \code{\link[EGAnet]{EGA}}. #' If supplied, then configural invariance check is skipped (i.e., configural #' invariance is assumed based on the given structure) -#' +#' #' @param iter Numeric (length = 1). #' Number of iterations to perform for the permutation. #' Defaults to \code{500} (recommended) @@ -26,101 +26,101 @@ #' which items should be removed during configural invariance (see \strong{Details} #' for more information). #' Defaults to \code{0.70} (recommended) -#' +#' #' @param configural.type Character (length = 1). #' Type of bootstrap to use for configural invariance in \code{\link[EGAnet]{bootEGA}}. #' Defaults to \code{"parametric"} -#' +#' #' @param corr Character (length = 1). #' Method to compute correlations. #' Defaults to \code{"auto"}. #' Available options: -#' +#' #' \itemize{ -#' +#' #' \item \code{"auto"} --- Automatically computes appropriate correlations for #' the data using Pearson's for continuous, polychoric for ordinal, #' tetrachoric for binary, and polyserial/biserial for ordinal/binary with #' continuous. To change the number of categories that are considered #' ordinal, use \code{ordinal.categories} #' (see \code{\link[EGAnet]{polychoric.matrix}} for more details) -#' -#' \item \code{"cor_auto"} --- Uses \code{\link[qgraph]{cor_auto}} to compute correlations. +#' +#' \item \code{"cor_auto"} --- Uses \code{\link[qgraph]{cor_auto}} to compute correlations. #' Arguments can be passed along to the function -#' -#' \item \code{"pearson"} --- Pearson's correlation is computed for all +#' +#' \item \code{"pearson"} --- Pearson's correlation is computed for all #' variables regardless of categories -#' -#' \item \code{"spearman"} --- Spearman's rank-order correlation is computed +#' +#' \item \code{"spearman"} --- Spearman's rank-order correlation is computed #' for all variables regardless of categories -#' +#' #' } -#' +#' #' For other similarity measures, compute them first and input them #' into \code{data} with the sample size (\code{n}) -#' +#' #' @param na.data Character (length = 1). #' How should missing data be handled? #' Defaults to \code{"pairwise"}. #' Available options: -#' +#' #' \itemize{ -#' +#' #' \item \code{"pairwise"} --- Computes correlation for all available cases between #' two variables -#' +#' #' \item \code{"listwise"} --- Computes correlation for all complete cases in the dataset -#' +#' #' } -#' +#' #' @param model Character (length = 1). #' Defaults to \code{"glasso"}. #' Available options: -#' +#' #' \itemize{ -#' +#' #' \item \code{"BGGM"} --- Computes the Bayesian Gaussian Graphical Model. #' Set argument \code{ordinal.categories} to determine #' levels allowed for a variable to be considered ordinal. #' See \code{?BGGM::estimate} for more details -#' +#' #' \item \code{"glasso"} --- Computes the GLASSO with EBIC model selection. #' See \code{\link[EGAnet]{EBICglasso.qgraph}} for more details -#' +#' #' \item \code{"TMFG"} --- Computes the TMFG method. #' See \code{\link[EGAnet]{TMFG}} for more details -#' +#' #' } -#' -#' @param algorithm Character or +#' +#' @param algorithm Character or #' \code{\link{igraph}} \code{cluster_*} function (length = 1). #' Defaults to \code{"walktrap"}. #' Three options are listed below but all are available #' (see \code{\link[EGAnet]{community.detection}} for other options): -#' +#' #' \itemize{ #' #' \item \code{"leiden"} --- See \code{\link[igraph]{cluster_leiden}} for more details -#' -#' \item \code{"louvain"} --- By default, \code{"louvain"} will implement the Louvain algorithm using -#' the consensus clustering method (see \code{\link[EGAnet]{community.consensus}} +#' +#' \item \code{"louvain"} --- By default, \code{"louvain"} will implement the Louvain algorithm using +#' the consensus clustering method (see \code{\link[EGAnet]{community.consensus}} #' for more information). This function will implement -#' \code{consensus.method = "most_common"} and \code{consensus.iter = 1000} +#' \code{consensus.method = "most_common"} and \code{consensus.iter = 1000} #' unless specified otherwise -#' +#' #' \item \code{"walktrap"} --- See \code{\link[igraph]{cluster_walktrap}} for more details -#' +#' #' } #' #' @param uni.method Character (length = 1). -#' What unidimensionality method should be used? +#' What unidimensionality method should be used? #' Defaults to \code{"louvain"}. #' Available options: -#' +#' #' \itemize{ #' #' \item \code{"expand"} --- Expands the correlation matrix with four variables correlated 0.50. -#' If number of dimension returns 2 or less in check, then the data +#' If number of dimension returns 2 or less in check, then the data #' are unidimensional; otherwise, regular EGA with no matrix #' expansion is used. This method was used in the Golino et al.'s (2020) #' \emph{Psychological Methods} simulation @@ -129,18 +129,18 @@ #' (\code{\link[igraph]{cluster_leading_eigen}}) #' on the empirical correlation matrix. If the number of dimensions is 1, #' then the Leading Eigenvector solution is used; otherwise, regular EGA -#' is used. This method was used in the Christensen et al.'s (2023) +#' is used. This method was used in the Christensen et al.'s (2023) #' \emph{Behavior Research Methods} simulation -#' +#' #' \item \code{"louvain"} --- Applies the Louvain algorithm (\code{\link[igraph]{cluster_louvain}}) -#' on the empirical correlation matrix. If the number of dimensions is 1, -#' then the Louvain solution is used; otherwise, regular EGA is used. +#' on the empirical correlation matrix. If the number of dimensions is 1, +#' then the Louvain solution is used; otherwise, regular EGA is used. #' This method was validated Christensen's (2022) \emph{PsyArXiv} simulation. #' Consensus clustering can be used by specifying either #' \code{"consensus.method"} or \code{"consensus.iter"} -#' +#' #' } -#' +#' #' @param ncores Numeric (length = 1). #' Number of cores to use in computing results. #' Defaults to \code{ceiling(parallel::detectCores() / 2)} or half of your @@ -149,13 +149,13 @@ #' #' If you're unsure how many cores your computer has, #' then type: \code{parallel::detectCores()} -#' +#' #' @param seed Numeric (length = 1). #' Defaults to \code{NULL} or random results. #' Set for reproducible results. #' See \href{https://github.com/hfgolino/EGAnet/wiki/Reproducibility-and-PRNG}{Reproducibility and PRNG} #' for more details on random number generation in \code{\link{EGAnet}} -#' +#' #' @param verbose Boolean (length = 1). #' Should progress be displayed? #' Defaults to \code{TRUE}. @@ -166,17 +166,17 @@ #' \code{\link[EGAnet]{network.estimation}}, #' \code{\link[EGAnet]{community.detection}}, #' \code{\link[EGAnet]{community.consensus}}, -#' \code{\link[EGAnet]{EGA}}, +#' \code{\link[EGAnet]{EGA}}, #' \code{\link[EGAnet]{bootEGA}}, and #' \code{\link[EGAnet]{net.loads}} -#' +#' #' @details In traditional psychometrics, measurement invariance is performed in #' sequential testing from more flexible (more free parameters) to more rigid #' (fewer free parameters) structures. Measurement invariance in network #' psychometrics is no different. -#' +#' #' \strong{Configural Invariance} -#' +#' #' To establish configural invariance, the data are collapsed across groups #' and a common sample structure is identified used \code{\link[EGAnet]{bootEGA}} #' and \code{\link[EGAnet]{itemStability}}. If some variables have a replication @@ -189,11 +189,11 @@ #' variables \emph{are} removed, then configural invariance is not met for the #' original structure. Any removal would suggest only partial configural invariance #' is met. -#' +#' #' \strong{Metric Invariance} -#' +#' #' The variables that remain after configural invariance are submitted to metric -#' invariance. First, each group estimates a network and then network loadings +#' invariance. First, each group estimates a network and then network loadings #' (\code{\link[EGAnet]{net.loads}}) are computed using the assigned #' community memberships (determined during configural invariance). Then, #' the difference between the assigned loadings of the groups is computed. This @@ -207,62 +207,62 @@ #' differences for each variable. Both uncorrected and false discovery rate #' corrected \emph{p}-values are returned in the results. Uncorrected \emph{p}-values #' are flagged for significance along with the direction of group differences. -#' +#' #' \strong{Three or More Groups} -#' +#' #' At this time, only two groups are supported. There is a method proposed to #' test three or more groups in Jamison, Golino, and Christensen (2023) but #' this approach has not been thoroughly vetted and validated. Future versions #' of the package will provide support for three or more groups once there is #' an established consensus for best practice. -#' +#' #' For more details, see Jamison, Golino, and Christensen (2023) #' #' @return Returns a list containing: -#' +#' #' \item{configural.results}{\code{\link[EGAnet]{bootEGA}} results from #' the final run that produced configural invariance. This output will be #' output on the final run of unsuccessful configural invariance runs} -#' +#' #' \item{memberships}{Original memberships provided in \code{structure} #' or from \code{\link[EGAnet]{EGA}} if \code{structure = NULL}} #' #' \item{EGA}{Original \code{\link[EGAnet]{EGA}} results for the full sample} -#' +#' #' \item{groups}{A list containing: -#' +#' #' \itemize{ -#' +#' #' \item \code{\link[EGAnet]{EGA}} --- \code{\link[EGAnet]{EGA}} results for each group -#' +#' #' \item \code{loadings} --- Network loadings (\code{\link[EGAnet]{net.loads}}) for each group -#' +#' #' \item \code{loadingsDifference} --- Difference between the dominant loadings of each group -#' +#' #' } -#' +#' #' } -#' +#' #' \item{permutation}{A list containing: -#' +#' #' \itemize{ -#' +#' #' \item \code{groups} --- Permutated groups acorss iterations -#' +#' #' \item \code{loadings} --- Network loadings (\code{\link[EGAnet]{net.loads}}) for each group for each permutation -#' +#' #' \item \code{loadingsDifference} --- Difference between the dominant loadings of each group for each permutation -#' +#' #' } -#' +#' #' } -#' +#' #' \item{results}{Data frame of the results (which are printed)} #' #' @author Laura Jamison , #' Hudson F. Golino , and #' Alexander P. Christensen , -#' +#' #' @references #' \strong{Original implementation} \cr #' Jamison, L., Golino, H., & Christensen, A. P. (2023). @@ -272,26 +272,26 @@ #' @examples #' # Load data #' wmt <- wmt2[-1,7:24] -#' +#' #' # Groups #' groups <- rep(1:2, each = nrow(wmt) / 2) -#' +#' #' \dontrun{ #' # Measurement invariance #' results <- invariance(wmt, groups, ncores = 2) -#' +#' #' # Plot with uncorrected alpha = 0.05 #' plot(results, p_type = "p", p_value = 0.05) -#' +#' #' # Plot with BH-corrected alpha = 0.10 #' plot(results, p_type = "p_BH", p_value = 0.10)} -#' +#' #' @seealso \code{\link[EGAnet]{plot.EGAnet}} for plot usage in \code{\link{EGAnet}} -#' +#' #' @export #' # Measurement Invariance -# Updated 24.10.2023 +# Updated 09.03.2024 invariance <- function( # `invariance` arguments data, groups, structure = NULL, @@ -300,17 +300,17 @@ invariance <- function( # standard arguments corr = c("auto", "cor_auto", "pearson", "spearman"), na.data = c("pairwise", "listwise"), - model = c("BGGM", "glasso", "TMFG"), + model = c("BGGM", "glasso", "TMFG"), algorithm = c("leiden", "louvain", "walktrap"), uni.method = c("expand", "LE", "louvain"), ncores, seed = NULL, verbose = TRUE, ... ) { - + # Store random state (if there is one) store_state() - + # Check for missing arguments (argument, default, function) configural.type <- set_default(configural.type, "parametric", invariance) corr <- set_default(corr, "auto", invariance) @@ -319,71 +319,71 @@ invariance <- function( algorithm <- set_default(algorithm, "walktrap", community.detection) uni.method <- set_default(uni.method, "louvain", community.unidimensional) if(missing(ncores)){ncores <- ceiling(parallel::detectCores() / 2)} - + # Argument errors (returns 'data' and 'groups') error_return <- invariance_errors( data, groups, iter, configural.threshold, ncores, seed, verbose ) - + # Get data and groups data <- error_return$data; groups <- error_return$groups - + # Ensure data has variable names data <- ensure_dimension_names(data) - + # Get dimensions and dimension names of the data original_dimensions <- dim(data) original_dimension_names <- dimnames(data) - + # Get ellipse arguments ellipse <- list(...) - + # Check for seed if(is.null(seed)){ # Send message about seed message("Argument 'seed' is set to `NULL`. Results will not be reproducible. Set 'seed' for reproducible results\n") } - + # Check for legacy argument 'memberships if("memberships" %in% names(ellipse)){ structure <- ellipse$memberships } - + # Get unique groups unique_groups <- na.omit(unique(groups)) - + # Send warning about only two groups (for now) if(length(unique_groups) > 2){ - + # Send warning warning( "More than two groups is not yet supported. Using the first two groups...", call. = FALSE ) - + # Update unique groups unique_groups <- unique_groups[c(1L, 2L)] - + # Keep indices keep_index <- groups %in% unique_groups - + # Update groups groups <- groups[keep_index] - + # Update data data <- data[keep_index,] - + # Update data dimensions dimensions <- dim(data) - + } - + # If structure is supplied, then skip configural invariance if(is.null(structure)){ - + # Send message message("Testing configural invariance...") - + # Perform configural invariance configural_results <- configural( data = data, iter = iter, structure = structure, @@ -393,113 +393,113 @@ invariance <- function( uni.method = uni.method, ncores = ncores, seed = seed, verbose = verbose, ... ) - + # Check for configural invariance if(!configural_results$configural_flag){ - + # Send message message("\nConfigural invariance was not found. Terminating invariance testing...") - + # Return configural invariance results return(configural_results) - + } - + # Configural invariance was found, continue with metric - + # Send message message(paste( - "\nConfigural invariance was found with", + "\nConfigural invariance was found with", length(configural_results$stable_items), "variables\n" )) - + # Update data data <- configural_results$data - + # Update original EGA original_EGA <- configural_results$boot_object$EGA - + # Set structure based on original `EGA` structure <- original_EGA$wc - + }else{ - + # Remove attributes structure <- remove_attributes(structure) - + # If not `NULL`, then check for error in object type object_error(structure, c("vector", "matrix", "data.frame"), "invariance") - + # Make sure 'structure' is a vector structure <- force_vector(structure) - + # Make sure 'structure' has names names(structure) <- original_dimension_names[[2]] - + } - + # Send message about continuing on with metric invariance message("Testing metric invariance...") - + # Update dimensions and dimension names of the data dimensions <- dim(data) dimension_names <- dimnames(data) # Get community names community_names <- as.character(unique(structure)) - + # Obtain original group EGAs group_ega <- lapply(unique_groups, function(group){ - + # Return `EGA` EGA( - data = data[groups == group,], + data = data[groups == group,], corr = corr, model = model, algorithm = algorithm, uni.method = uni.method, plot.EGA = FALSE, ... ) - + }) - + # Rename list names(group_ega) <- unique_groups - + # Original network loadings group_loadings <- lapply(group_ega, function(x){ - + # Obtain loadings loadings <- as.matrix( net.loads(A = x$network, wc = structure, ...)$std ) - + # Reorder and return loadings return(loadings[dimension_names[[2]], community_names, drop = FALSE]) - + }) - + # Original difference original_difference <- group_loadings[[1]] - group_loadings[[2]] - + # Obtain original assigned difference original_assigned_difference <- ulapply( community_names, function(community){ original_difference[structure == community, community] } ) - + # Ensure same order as original data original_assigned_difference <- original_assigned_difference[dimension_names[[2]]] - + # Get seeds seeds <- reproducible_seeds(iter, seed) - + # Permutate groups perm_groups <- lapply( seeds, function(seed_value){ shuffle(groups, seed = seed_value) } ) - + # Perform permutation estimation of loadings permutated_loadings <- parallel_process( iterations = iter, @@ -509,68 +509,68 @@ invariance <- function( data, corr, model, algorithm, uni.method, ... ){ - + # Estimate loadings return( # By groups lapply(unique_groups, function(group){ - - + + # Get network network <- EGA( - data = data[permutation == group,], + data = data[permutation == group,], corr = corr, model = model, algorithm = algorithm, uni.method = uni.method, plot.EGA = FALSE, ... )$network - + # Obtain loadings loadings <- as.matrix( - net.loads(A = network, wc = structure)$std + net.loads(A = network, wc = structure, ...)$std ) - + # Return loadings return(loadings) - + }) ) - - }, + + }, # Make sure the additional objects get in unique_groups = unique_groups, structure = structure, data = data, corr = corr, model = model, algorithm = algorithm, uni.method = uni.method, ..., ncores = ncores, progress = verbose ) - + # Compute differences (ensure same ordering) difference_list <- lapply(permutated_loadings, function(x){ x[[1]][dimension_names[[2]], community_names, drop = FALSE] - x[[2]][dimension_names[[2]], community_names, drop = FALSE] }) - + # Obtain assigned loadings only assigned_list <- lapply(difference_list, function(one_difference){ - + # Get differences differences <- ulapply( community_names, function(community){ one_difference[structure == community, community] } ) - + # Ensure same order as original data return(differences[dimension_names[[2]]]) - + }) - + # Create results permutation_counts <- lapply(assigned_list, function(x){ abs(x) >= abs(original_assigned_difference) }) - + # Replace the first permutation with all TRUE (original differences) permutation_counts[[1]] <- rep(TRUE, dimensions[2]) - + # Compute p-values p_value <- rowMeans( do.call(cbind, permutation_counts), na.rm = TRUE @@ -583,16 +583,16 @@ invariance <- function( p = round(p_value, 3), p_BH = round(p.adjust(p_value, method = "BH"), 3) ) - + # Order by dimension results_df <- results_df[order(results_df$Membership),] - + # Add significance sig <- swiftelse(results_df$p <= 0.10, ".", "") sig <- swiftelse(results_df$p <= 0.05, "*", sig) sig <- swiftelse(results_df$p <= 0.01, "**", sig) results_df$sig <- swiftelse(results_df$p <= 0.001, "***", sig) - + # Add direction direction <- paste( unique_groups[1], @@ -600,7 +600,7 @@ invariance <- function( unique_groups[2] ) results_df$Direction <- swiftelse(results_df$p <= 0.05, direction, "") - + # Results list results <- list( configural.results = swiftelse( @@ -626,16 +626,16 @@ invariance <- function( ), results = results_df ) - + # Add class class(results) <- "invariance" - + # Restore random state (if there is one) restore_state() - + # Return results return(results) - + } # Bug checking ---- @@ -654,64 +654,64 @@ invariance_errors <- function( ncores, seed, verbose ) { - + # 'data' errors object_error(data, c("matrix", "data.frame", "tibble"), "invariance") - + # Check for tibble if(get_object_type(data) == "tibble"){ data <- as.data.frame(data) } - + # 'groups' errors object_error(groups, c("vector", "matrix", "data.frame"), "invariance") groups <- force_vector(groups) length_error(groups, dim(data)[1], "invariance") - + # 'iter' errors length_error(iter, 1, "invariance") typeof_error(iter, "numeric", "invariance") range_error(iter, c(1, Inf), "invariance") - + # 'configural.threshold' errors length_error(configural.threshold, 1, "invariance") typeof_error(configural.threshold, "numeric", "invariance") range_error(configural.threshold, c(0, 1), "invariance") - + # 'ncores' errors length_error(ncores, 1, "invariance") typeof_error(ncores, "numeric", "invariance") range_error(ncores, c(1, parallel::detectCores()), "invariance") - + # 'seed' errors if(!is.null(seed)){ length_error(seed, 1, "invariance") typeof_error(seed, "numeric", "invariance") range_error(seed, c(0, Inf), "invariance") } - + # 'verbose' errors length_error(verbose, 1, "invariance") typeof_error(verbose, "logical", "invariance") - + # Return usable data and groups return(list(data = usable_data(data, verbose), groups = groups)) - + } -#' @exportS3Method +#' @exportS3Method # S3 Print Method ---- # Updated 10.07.2023 print.invariance <- function(x, ...) { - + # Print results "as-is" print(x$results) cat("----\n") # Add breakspace and significance code cat("Signif. code: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 'n.s.' 1") - + } -#' @exportS3Method +#' @exportS3Method # S3 Summary Method ---- # Updated 10.07.2023 summary.invariance <- function(object, ...) { @@ -727,43 +727,43 @@ group_setup <- function( ... ) { - + # Get ellipse arguments for defaults ellipse <- list(...) - + # Set edge size if(!"edge.size" %in% ellipse){ ellipse$edge.size <- 8 # default in `basic_plot_setup` } - + # Make copy of EGA object EGA_object_copy <- EGA_object - + # Set noninvariant edges to 1 for(i in seq_len(nodes)){ - + # Check for noninvariance if(noninvariant[i]){ - + # Get community index community_index <- EGA_object$wc == EGA_object$wc[i] - + # Target node target_node <- EGA_object$network[i, community_index] - + # Replace non-zero edges with 1 EGA_object_copy$network[i, community_index] <- EGA_object_copy$network[community_index, i] <- swiftelse(target_node != 0, 1, target_node) - + } - + } - + # Set up plot for edges edge_plot <- basic_plot_setup( - network = EGA_object_copy$network, - wc = EGA_object_copy$wc, ..., + network = EGA_object_copy$network, + wc = EGA_object_copy$wc, ..., arguments = TRUE ) @@ -774,36 +774,36 @@ group_setup <- function( plot_ARGS$edge.lty[ edge_plot$ARGS$edge.size == ellipse$edge.size ] <- "dashed" - + # Set up EGA arguments plot_ARGS$net <- NULL plot_ARGS$network <- EGA_object$network plot_ARGS$wc <- EGA_object$wc - + # Return plot arguments return(plot_ARGS) - + } -#' @exportS3Method +#' @exportS3Method # S3 Plot Method ---- # Updated 08.10.2023 plot.invariance <- function(x, p_type = c("p", "p_BH"), p_value = 0.05, ...) { - + # Set default for p-value type p_type <- swiftelse(missing(p_type), "p", match.arg(p_type)) - + # Check for appropriate p-value range range_error(p_value, c(0, 1), "plot.invariance") - + # Ensure same memberships x$groups$EGA[[1]]$wc <- x$memberships x$groups$EGA[[2]]$wc <- x$memberships - + # Obtain noninvariant items noninvariant <- x$results[names(x$memberships), p_type] <= p_value - + # Get number of nodes nodes <- length(noninvariant) @@ -813,24 +813,24 @@ plot.invariance <- function(x, p_type = c("p", "p_BH"), p_value = 0.05, ...) wc = x$groups$EGA[[1]]$wc, ..., arguments = TRUE ) - + # Get plot arguments second_ARGS <- first_group$ARGS - + # Remove some arguments from `first_ARGS` ## Essentially, the same call but allows some freedom second_ARGS[c( "net", "edge.alpha", "edge.color", "edge.lty", "edge.size" )] <- NULL - + # Add network and memberships second_ARGS$network <- x$groups$EGA[[2]]$network second_ARGS$wc <- x$groups$EGA[[2]]$wc second_ARGS$arguments <- TRUE - + # Set up second group plot second_group <- do.call(basic_plot_setup, second_ARGS) - + # Get updated plots for each group ## First group first_group <- do.call( @@ -850,11 +850,11 @@ plot.invariance <- function(x, p_type = c("p", "p_BH"), p_value = 0.05, ...) nodes = nodes, noninvariant = noninvariant ) ) - + # Update `alpha` guide first_group$guides$colour$override.aes$alpha <- 0.25 second_group$guides$colour$override.aes$alpha <- 0.75 - + # Set up p-value title if(p_type == "p") { invariant_title <- bquote( @@ -871,16 +871,16 @@ plot.invariance <- function(x, p_type = c("p", "p_BH"), p_value = 0.05, ...) paste("Noninvariant (", italic(p)[adj.], " < ", .(p_value), ")") ) } - + # Update `title` guide first_group$guides$colour$title <- invariant_title - + second_group$guides$colour$title <- noninvariant_title - + # Update `title.position` guide first_group$guides$colour$title.position <- "top" second_group$guides$colour$title.position <- "top" - + # Adjust size and position first_group <- first_group + ggplot2::theme( @@ -899,7 +899,7 @@ plot.invariance <- function(x, p_type = c("p", "p_BH"), p_value = 0.05, ...) legend = "bottom", common.legend = FALSE ) - + } #' @noRd @@ -908,17 +908,17 @@ plot.invariance <- function(x, p_type = c("p", "p_BH"), p_value = 0.05, ...) configural <- function( data, iter, structure, configural.threshold, configural.type, corr, na.data, model, - algorithm, uni.method, ncores, seed, verbose, + algorithm, uni.method, ncores, seed, verbose, ... ){ - + # Initialize check for items and stability values items <- dim(data)[2] stability <- FALSE - + # Perform `while` loop while(items > 0 & !stability){ - + # Perform bootEGA boot <- bootEGA( data = data, corr = corr, na.data = na.data, @@ -927,28 +927,28 @@ configural <- function( type = configural.type, ncores = ncores, EGA.type = "EGA", typicalStructure = FALSE, plot.typicalStructure = FALSE, - seed = seed, verbose = verbose, + seed = seed, verbose = verbose, clear = TRUE, suppress = TRUE, # additional internal arguments to `bootEGA` ... ) - + # Perform itemStability item_stability <- itemStability(boot, IS.plot = FALSE, structure = structure) - + # Stable items stable_items <- item_stability$item.stability$empirical.dimensions >= configural.threshold - + # Perform checks stability <- all(stable_items) - + # Number of stable items items <- sum(stable_items) - + # Update data data <- data[,stable_items] - + } - + # Send results return( list( @@ -959,6 +959,6 @@ configural <- function( configural_flag = items > 0 ) ) - + } diff --git a/R/net.loads.R b/R/net.loads.R index 1f558105..19273714 100644 --- a/R/net.loads.R +++ b/R/net.loads.R @@ -448,7 +448,7 @@ obtain_signs <- function(target_network) #' @noRd # Experimental loadings ---- -# Updated 25.02.2024 +# Updated 09.03.2024 experimental_loadings <- function( A, wc, nodes, node_names, communities, unique_communities @@ -574,7 +574,7 @@ experimental_loadings <- function( #' @noRd # Standardize loadings ---- -# Updated 08.03.2024 +# Updated 09.03.2024 standardize <- function(unstandardized, loading.method, A, wc) { @@ -583,41 +583,15 @@ standardize <- function(unstandardized, loading.method, A, wc) return(t(t(unstandardized) / sqrt(colSums(abs(unstandardized), na.rm = TRUE)))) }else if(loading.method == "experimental"){ - # Original community order - original_order <- dimnames(unstandardized)[[2]] - - # Set community sequence - community_sequence <- seq_len(dim(unstandardized)[2]) - - # Set diagonal of network to 1 - diag(A) <- 1 - - # Get eigenvectors and eigenvalues - eigens <- eigen(A) - - # Get signs - signs <- sign(unstandardized) - - # Align loadings - alignment <- fungible::faAlign( - F1 = eigens$vectors[,community_sequence], - F2 = unstandardized + # Get sums within-community + sums <- nvapply( + dimnames(unstandardized)[[2]], function(community){ + sum(abs(unstandardized[wc == community, community]), na.rm = TRUE) + } ) - # Pre-compute absolute values for standardization - absolute <- abs(unstandardized) - - # Standardize loadings - standardized <- absolute / (absolute + 1) - - # Get sorted order to pre-multiply by eigenvalues - sorted <- standardized[,alignment$FactorMap["Sorted Order",]] - - # Get scaled loadings - loadings <- t(t(sorted) * sqrt(eigens$values[community_sequence])) - # Return loadings - return(loadings[,original_order] * signs) + return(t(t(unstandardized) / (sums)^(1 / sqrt(fast_table(wc))))) } diff --git a/man/invariance.Rd b/man/invariance.Rd index 6194542c..1b7a7724 100644 --- a/man/invariance.Rd +++ b/man/invariance.Rd @@ -63,13 +63,13 @@ continuous. To change the number of categories that are considered ordinal, use \code{ordinal.categories} (see \code{\link[EGAnet]{polychoric.matrix}} for more details) -\item \code{"cor_auto"} --- Uses \code{\link[qgraph]{cor_auto}} to compute correlations. +\item \code{"cor_auto"} --- Uses \code{\link[qgraph]{cor_auto}} to compute correlations. Arguments can be passed along to the function -\item \code{"pearson"} --- Pearson's correlation is computed for all +\item \code{"pearson"} --- Pearson's correlation is computed for all variables regardless of categories -\item \code{"spearman"} --- Spearman's rank-order correlation is computed +\item \code{"spearman"} --- Spearman's rank-order correlation is computed for all variables regardless of categories } @@ -110,7 +110,7 @@ See \code{\link[EGAnet]{TMFG}} for more details }} -\item{algorithm}{Character or +\item{algorithm}{Character or \code{\link{igraph}} \code{cluster_*} function (length = 1). Defaults to \code{"walktrap"}. Three options are listed below but all are available @@ -120,10 +120,10 @@ Three options are listed below but all are available \item \code{"leiden"} --- See \code{\link[igraph]{cluster_leiden}} for more details -\item \code{"louvain"} --- By default, \code{"louvain"} will implement the Louvain algorithm using -the consensus clustering method (see \code{\link[EGAnet]{community.consensus}} +\item \code{"louvain"} --- By default, \code{"louvain"} will implement the Louvain algorithm using +the consensus clustering method (see \code{\link[EGAnet]{community.consensus}} for more information). This function will implement -\code{consensus.method = "most_common"} and \code{consensus.iter = 1000} +\code{consensus.method = "most_common"} and \code{consensus.iter = 1000} unless specified otherwise \item \code{"walktrap"} --- See \code{\link[igraph]{cluster_walktrap}} for more details @@ -131,14 +131,14 @@ unless specified otherwise }} \item{uni.method}{Character (length = 1). -What unidimensionality method should be used? +What unidimensionality method should be used? Defaults to \code{"louvain"}. Available options: \itemize{ \item \code{"expand"} --- Expands the correlation matrix with four variables correlated 0.50. -If number of dimension returns 2 or less in check, then the data +If number of dimension returns 2 or less in check, then the data are unidimensional; otherwise, regular EGA with no matrix expansion is used. This method was used in the Golino et al.'s (2020) \emph{Psychological Methods} simulation @@ -147,12 +147,12 @@ expansion is used. This method was used in the Golino et al.'s (2020) (\code{\link[igraph]{cluster_leading_eigen}}) on the empirical correlation matrix. If the number of dimensions is 1, then the Leading Eigenvector solution is used; otherwise, regular EGA -is used. This method was used in the Christensen et al.'s (2023) +is used. This method was used in the Christensen et al.'s (2023) \emph{Behavior Research Methods} simulation \item \code{"louvain"} --- Applies the Louvain algorithm (\code{\link[igraph]{cluster_louvain}}) -on the empirical correlation matrix. If the number of dimensions is 1, -then the Louvain solution is used; otherwise, regular EGA is used. +on the empirical correlation matrix. If the number of dimensions is 1, +then the Louvain solution is used; otherwise, regular EGA is used. This method was validated Christensen's (2022) \emph{PsyArXiv} simulation. Consensus clustering can be used by specifying either \code{"consensus.method"} or \code{"consensus.iter"} @@ -184,7 +184,7 @@ Set to \code{FALSE} to not display progress} \code{\link[EGAnet]{network.estimation}}, \code{\link[EGAnet]{community.detection}}, \code{\link[EGAnet]{community.consensus}}, -\code{\link[EGAnet]{EGA}}, +\code{\link[EGAnet]{EGA}}, \code{\link[EGAnet]{bootEGA}}, and \code{\link[EGAnet]{net.loads}}} } @@ -260,7 +260,7 @@ is met. \strong{Metric Invariance} The variables that remain after configural invariance are submitted to metric -invariance. First, each group estimates a network and then network loadings +invariance. First, each group estimates a network and then network loadings (\code{\link[EGAnet]{net.loads}}) are computed using the assigned community memberships (determined during configural invariance). Then, the difference between the assigned loadings of the groups is computed. This