PLS-MGA feature for SEMinR (#255)

* Working estimate_pls_mga function

* Test file to check across runs of PLS-MGA

* Faithfully rerun original model for subgroups

- R/estimate_pls.R - rerun.pls_model() to faithfully rerun estimated model
- R/estimate_pls_mga.R - use rerun to re-estimate models for subgroups

* Fixes for CRAN check

* Print function for seminr_pls_mga results

* Fix examples and documentation for CRAN

* Added demo for PLS-MGA

* Update README with re-org of documentation

NAMESPACE

@@ -17,13 +17,15 @@ S3method(plot,summary.predict_pls_model)
 S3method(predict,seminr_model)
 S3method(print,list_output)
 S3method(print,measurement_model_evaluation.seminr_model)
+S3method(print,seminr_pls_mga)
 S3method(print,seminr_theme)
 S3method(print,summary.boot_seminr_model)
 S3method(print,summary.cbsem_model)
 S3method(print,summary.cfa_model)
 S3method(print,summary.predict_pls_model)
 S3method(print,summary.seminr_model)
 S3method(print,table_output)
+S3method(rerun,pls_model)
 S3method(summary,boot_seminr_model)
 S3method(summary,cbsem_model)
 S3method(summary,cfa_model)
@@ -48,6 +50,7 @@ export(estimate_cbsem)
 export(estimate_cfa)
 export(estimate_lavaan_ten_berge)
 export(estimate_pls)
+export(estimate_pls_mga)
 export(fSquared)
 export(get_theme_doc)
 export(higher_composite)
@@ -77,6 +80,7 @@ export(reflective)
 export(regression_weights)
 export(relationships)
 export(report_paths)
+export(rerun)
 export(rho_A)
 export(save_plot)
 export(seminr_theme_create)

R/compute_metrics.R

@@ -10,16 +10,14 @@ compute_AVE <- function(lambdas) {
   sum(lambdas^2) / length(lambdas)
 }
 
-# Returns R-sq of a dv given correlation matrix of ivs, dv
-#
-# @param cor_matrix A correlation matrix that includes ivs and dv
-# @param dv_name Character string of dependent variable
-# @param iv_names Vector of character strings for independent variables
-#
-# @examples
-# cors <- cbsem_summary$descriptives$correlations$constructs
-# cor_rsq(cors, dv_name = "Value", iv_names = c("Image", "Quality"))
-#
+#' Returns R-sq of a dv given correlation matrix of ivs, dv
+#' cors <- cbsem_summary$descriptives$correlations$constructs
+#' cor_rsq(cors, dv_name = "Value", iv_names = c("Image", "Quality"))
+#'
+#' @param cor_matrix A correlation matrix that includes ivs and dv
+#' @param dv_name Character string of dependent variable
+#' @param iv_names Vector of character strings for independent variables
+#'
 cor_rsq <- function(cor_matrix, dv_name, iv_names) {
   iv_cors <- cor_matrix[iv_names, iv_names]
   dv_cors <- cor_matrix[iv_names, dv_name]

R/estimate_pls.R

@@ -104,6 +104,8 @@ estimate_pls <- function(data,
                          missing_value = NA,
                          maxIt=300,
                          stopCriterion=7) {
+  # NOTE: update rerun.pls_model() if parameters change!
+
   message("Generating the seminr model")
   data[data == missing_value] <- NA
   rawdata <- data
@@ -150,9 +152,13 @@ estimate_pls <- function(data,
                                    maxIt=maxIt,
                                    stopCriterion=stopCriterion,
                                    measurement_mode_scheme = measurement_mode_scheme)
+
+  # Store all settings needed for a rerun
   seminr_model$data <- data
   seminr_model$rawdata <- rawdata
   seminr_model$measurement_model <- measurement_model
+  seminr_model$structural_model <- structural_model
+  seminr_model$settings$inner_weights <- inner_weights
   seminr_model$settings$missing_value <- missing_value
   seminr_model$settings$maxIt <- maxIt
   seminr_model$settings$stopCriterion <- stopCriterion
@@ -175,3 +181,76 @@ estimate_pls <- function(data,
   return(seminr_model)
 }
 
+not_null <- function(a, b) {
+  if(!is.null(a)) {
+    a
+  } else {
+    b
+  }
+}
+
+#' Reruns a previously specified seminr model/analysis
+#'
+#' @param x An estimated seminr_model object - refer to specific rerun methods
+#'
+#' @param ... Any parameters to change during the rerun.
+#'
+#' @return A re-estimated model of the same class
+#'
+#' @seealso \code{\link{rerun.pls_model}}
+#'
+#' @export
+rerun <- function (x, ...) {
+  UseMethod("rerun", x)
+}
+
+#' Reruns a previously specified seminr PLS model
+#'
+#' @param x An estimated pls_model object produced by \code{\link{estimate_pls}}
+#'
+#' @param ... Any parameters to change during the re-estimation (e.g., data, measurement_model, etc.)
+#'
+#' @return A re-estimated pls_model object
+#'
+#' @examples
+#'
+#' mobi <- mobi
+#'
+#' mobi_mm <- constructs(
+#'   composite("Image",        multi_items("IMAG", 1:5)),
+#'   composite("Loyalty",      multi_items("CUSL", 1:3))
+#' )
+#'
+#' mobi_sm <- relationships(
+#'   paths(from = "Image",        to = c("Loyalty"))
+#' )
+#'
+#' mobi_pls <- estimate_pls(data = mobi,
+#'                          measurement_model = mobi_mm,
+#'                          structural_model = mobi_sm,
+#'                          missing = mean_replacement,
+#'                          missing_value = NA)
+#'
+#' # Re-estimate model faithfully
+#' mobi_pls2 <- rerun(mobi_pls)
+#'
+#' # Re-estimated model with altered measurement model
+#' mobi_pls3 <- rerun(mobi_pls, measurement_model=as.reflective(mobi_mm))
+#'
+#' @export
+rerun.pls_model <- function(x, ...) {
+  args <- list(...)
+
+  estimate_pls(
+    data              = not_null(args$data,              x$rawdata),
+    measurement_model = not_null(args$measurement_model, x$measurement_model),
+    structural_model  = not_null(args$structural_model,  x$structural_model),
+    model             = not_null(args$model,             x$model),
+    inner_weights     = not_null(args$inner_weights,     x$settings$inner_weights),
+    missing           = not_null(args$missing,           x$settings$missing),
+    missing_value     = not_null(args$missing_value,     x$settings$missing_value),
+    maxIt             = not_null(args$maxIt,             x$settings$maxIt),
+    stopCriterion     = not_null(args$stopCriterion,     x$settings$stopCriterion)
+  )
+}
+

R/estimate_pls_mga.R

@@ -0,0 +1,124 @@
+#' Performs PLS-MGA to report significance of path differences between two subgroups of data
+#'
+#' @param pls_model SEMinR PLS model estimated on the full sample
+#' @param condition logical vector of TRUE/FALSE indicating which rows of sample data are in group 1
+#' @param nboot number of bootstrap resamples to use in PLS-MGA
+#' @param ... any further parameters for bootstrapping (e.g., cores)
+#'
+#' @examples
+#' mobi <- mobi
+#'
+#' #seminr syntax for creating measurement model
+#' mobi_mm <- constructs(
+#'   composite("Image",        multi_items("IMAG", 1:5)),
+#'   composite("Expectation",  multi_items("CUEX", 1:3)),
+#'   composite("Quality",      multi_items("PERQ", 1:7)),
+#'   composite("Value",        multi_items("PERV", 1:2)),
+#'   composite("Satisfaction", multi_items("CUSA", 1:3)),
+#'   composite("Complaints",   single_item("CUSCO")),
+#'   composite("Loyalty",      multi_items("CUSL", 1:3))
+#' )
+#'
+#' #seminr syntax for creating structural model
+#' mobi_sm <- relationships(
+#'   paths(from = "Image",        to = c("Expectation", "Satisfaction", "Loyalty")),
+#'   paths(from = "Expectation",  to = c("Quality", "Value", "Satisfaction")),
+#'   paths(from = "Quality",      to = c("Value", "Satisfaction")),
+#'   paths(from = "Value",        to = c("Satisfaction")),
+#'   paths(from = "Satisfaction", to = c("Complaints", "Loyalty")),
+#'   paths(from = "Complaints",   to = "Loyalty")
+#' )
+#'
+#' mobi_pls <- estimate_pls(data = mobi,
+#'                          measurement_model = mobi_mm,
+#'                          structural_model = mobi_sm,
+#'                          missing = mean_replacement,
+#'                          missing_value = NA)
+#'
+#' # Should usually use nboot ~2000 and don't specify cores for full parallel processing
+#' mobi_mga <- estimate_pls_mga(mobi_pls, mobi$CUEX1 < 8, nboot=100, cores = 2)
+#'
+#' @export
+estimate_pls_mga <- function(pls_model, condition, nboot = 2000, ...) {
+  pls_data <- pls_model$rawdata
+
+  # Given a beta report matrix (paths as rows) get estimates form a path_coef matrix
+  path_estimate <- function(path, path_coef) {
+    path_coef[path["source"], path["target"]]
+  }
+
+  # Get all path estimates of a given beta metrix from a given path_coef matrix
+  # Typically used to apply on 3rd dimension of a 3x3 bootstrap paths array [from,to,boot]
+  boot_paths <- function(path_coef, beta_df) {
+    betas <- apply(beta_df, MARGIN=1, FUN=path_estimate, path_coef = path_coef)
+  }
+
+  # Allocate and Estimate Two Alternative Datasets + Models
+  group1_data <- pls_data[condition, ]
+  group2_data <- pls_data[!condition, ]
+
+  message("Estimating and bootstrapping groups...")
+
+  group1_model <- rerun(pls_model, data = group1_data)
+  group2_model <- rerun(pls_model, data = group2_data)
+
+  group1_boot <- bootstrap_model(seminr_model = group1_model, nboot = nboot, ...)
+  group2_boot <- bootstrap_model(seminr_model = group2_model, nboot = nboot, ...)
+
+  message("Computing similarity of groups")
+  # Produce beta report matrix on all paths (as rows)
+  beta <- as.data.frame(pls_model$smMatrix[,c("source", "target")])
+  path_names <- do.call(paste0, cbind(beta["source"], " -> ", beta["target"]))
+  rownames(beta) <- path_names
+  beta$estimate <- apply(beta, MARGIN = 1, FUN=path_estimate, path_coef = pls_model$path_coef)
+
+  beta$group1_beta <- apply(beta, MARGIN = 1, FUN=path_estimate, path_coef = group1_model$path_coef)
+  beta$group2_beta <- apply(beta, MARGIN = 1, FUN=path_estimate, path_coef = group2_model$path_coef)
+
+  beta_diff <- group1_model$path_coef - group2_model$path_coef
+  beta$diff <- apply(beta, MARGIN = 1, FUN=path_estimate, path_coef = beta_diff)
+
+  # Get bootstrapped paths for both groups
+  boot1_betas <- t(apply(group1_boot$boot_paths, MARGIN=3, FUN=boot_paths, beta_df=beta))
+  colnames(boot1_betas) <- path_names
+
+  boot2_betas <- t(apply(group2_boot$boot_paths, MARGIN=3, FUN=boot_paths, beta_df=beta))
+  colnames(boot2_betas) <- path_names
+
+  # PLSc may not resolve in some bootstrap runs - limit bootstrap paths to resolved number of boots
+  J <- min(dim(boot1_betas)[1], dim(boot2_betas)[1])
+  if (J < nboot) {
+    message(paste("NOTE: Using", J, "bootstrapped results of each group after removing inadmissible runs"))
+  }
+  boot1_betas <- boot1_betas[1:J,]
+  boot2_betas <- boot2_betas[1:J,]
+
+
+  # Insert bootstrap descriptives into beta matrix
+  beta$group1_beta_mean <- apply(boot1_betas, MARGIN=2, FUN=mean)
+  beta$group2_beta_mean <- apply(boot2_betas, MARGIN=2, FUN=mean)
+
+  # beta$group1_beta_sd <- apply(boot1_betas, MARGIN=2, FUN=sd)
+  # beta$group2_beta_sd <- apply(boot2_betas, MARGIN=2, FUN=sd)
+
+  # Compute PLS-MGA p-value
+  # see: Henseler, J., Ringle, C. M., & Sinkovics, R. R. (2009). The use of partial least squares path modeling in international marketing. In New challenges to international marketing. Emerald Group Publishing Limited.
+
+  Theta <- function(s) {
+    ifelse(s > 0, 1, 0)
+  }
+
+  beta_comparison <- function(i, beta, beta1_boots, beta2_boots) {
+    for_all <- expand.grid(beta1_boots[,i], beta2_boots[,i])
+    2*beta$group1_beta_mean[i] - for_all[,1] - 2*beta$group2_beta_mean[i] + for_all[,2]
+  }
+
+  pls_mga_p <- function(i, beta, beta1_boots, beta2_boots) {
+    1 - (sum(Theta(beta_comparison(i, beta, beta1_boots, beta2_boots))) / J^2)
+  }
+
+  beta$pls_mga_p <- sapply(1:nrow(beta), FUN=pls_mga_p, beta=beta, beta1_boots=boot1_betas, beta2_boots=boot2_betas)
+
+  class(beta) <- c("seminr_pls_mga", class(beta))
+  beta
+}

R/report_summary.R

@@ -181,3 +181,19 @@ plot.summary.predict_pls_model <- function(x, indicator, ...) {
   # Grid
   graphics::grid(nx = NULL, ny = NULL, col = "lightgray", lty = "dotted")
 }
+
+#' Summary function for PLS-MGA
+#'
+#' @param x estimated seminr_pls_mga object
+#' @param digits number of digits to print
+#' @param ... any further parameters for printing
+#'
+#' @export
+print.seminr_pls_mga <- function(x, digits=3, ...) {
+  stopifnot(inherits(x, "seminr_pls_mga"))
+  cat("\nPLS-MGA results:\n")
+  mga_report <- data.frame(x$source, .=rep("->", nrow(x)), x$target,
+                           x$group1_beta, x$group2_beta, x$pls_mga_p)
+  colnames(mga_report) <- c("from", "", "to", "group1", "group2","p")
+  print(mga_report, row.names = FALSE, digits = digits)
+}

README.Rmd

@@ -364,14 +364,23 @@ The **vignette** for Seminr can be found on [CRAN](https://cran.r-project.org/pa
 **Demo code** for various use cases with SEMinR can be found in the 
 [seminr/demo/](https://github.com/sem-in-r/seminr/tree/master/demo) folder or by running commands such as `demo("seminr-contained")` after installation.
 
-- [seminr-alternative-models.R](demo/seminr-alternative-models.R): Reuse measurement and structural components to easily create competing models
-- [seminr-cbsem-cfa-ecsi.R](demo/seminr-cbsem-cfa-ecsi.R): Conduct confirmatory model building using CFA and CBSEM
-- [seminr-pls-ecsi.R](demo/seminr-pls-ecsi.R): Conduct PLS path modeling
-- [seminr-pls-higher_order.R](demo/seminr-pls-higher_order.R): Define higher-order composites for PLS models
-- [seminr-pls-interaction.R](demo/seminr-pls-interaction.R): Define interactions between constructs in SEM models
-- [seminr-plsc-ecsi.R](demo/seminr-plsc-ecsi.R): Run PLSc to emulate common factors using in PLSc
-- [seminr-style-contained.R](demo/seminr-style-contained.R): Create and execute a SEM model in one function call
-- [seminr-dot-graph.R](demo/seminr-pls-dot-graph.R): Create a plot from a SEM model
+Model Specification:
+
+-   [seminr-cbsem-cfa-ecsi](demo/seminr-cbsem-cfa-ecsi.R): Conduct confirmatory model building using CFA and CBSEM
+-   [seminr-pls-ecsi](demo/seminr-pls-ecsi.R): Conduct PLS path modeling
+-   [seminr-pls-interaction](demo/seminr-pls-interaction.R): Define interactions between constructs in SEM models
+-   [seminr-pls-higher\_order](demo/seminr-pls-higher_order.R): Define higher-order composites for PLS models
+-   [seminr-pls-mga](demo/seminr-pls-mga.R): Assess structural differences between subgroups using PLS-MGA
+-   [seminr-plsc-ecsi](demo/seminr-plsc-ecsi.R): Run PLSc to emulate common factors using in PLSc
+
+Model Visualization:
+
+-   [seminr-dot-graph](demo/seminr-pls-dot-graph.R): Create a plot from a SEM model
+
+Syntax Style:
+
+-   [seminr-alternative-models](demo/seminr-alternative-models.R): Reuse measurement and structural components in multiple models
+-   [seminr-style-contained](demo/seminr-style-contained.R): Create and execute a SEM model in one function call
 
 ## Sister Projects
 

README.md

@@ -357,23 +357,31 @@ or by running the `vignette("SEMinR")` command after installation.
 folder or by running commands such as `demo("seminr-contained")` after
 installation.
 
--   [seminr-alternative-models.R](demo/seminr-alternative-models.R):
-    Reuse measurement and structural components to easily create
-    competing models
--   [seminr-cbsem-cfa-ecsi.R](demo/seminr-cbsem-cfa-ecsi.R): Conduct
+Model Specification:
+
+-   [seminr-cbsem-cfa-ecsi](demo/seminr-cbsem-cfa-ecsi.R): Conduct
     confirmatory model building using CFA and CBSEM
--   [seminr-pls-ecsi.R](demo/seminr-pls-ecsi.R): Conduct PLS path
-    modeling
--   [seminr-pls-higher\_order.R](demo/seminr-pls-higher_order.R): Define
-    higher-order composites for PLS models
--   [seminr-pls-interaction.R](demo/seminr-pls-interaction.R): Define
+-   [seminr-pls-ecsi](demo/seminr-pls-ecsi.R): Conduct PLS path modeling
+-   [seminr-pls-interaction](demo/seminr-pls-interaction.R): Define
     interactions between constructs in SEM models
--   [seminr-plsc-ecsi.R](demo/seminr-plsc-ecsi.R): Run PLSc to emulate
+-   [seminr-pls-higher\_order](demo/seminr-pls-higher_order.R): Define
+    higher-order composites for PLS models
+-   [seminr-pls-mga](demo/seminr-pls-mga.R): Assess structural
+    differences between subgroups using PLS-MGA
+-   [seminr-plsc-ecsi](demo/seminr-plsc-ecsi.R): Run PLSc to emulate
     common factors using in PLSc
--   [seminr-style-contained.R](demo/seminr-style-contained.R): Create
-    and execute a SEM model in one function call
--   [seminr-dot-graph.R](demo/seminr-pls-dot-graph.R): Create a plot
-    from a SEM model
+
+Model Visualization:
+
+-   [seminr-dot-graph](demo/seminr-pls-dot-graph.R): Create a plot from
+    a SEM model
+
+Syntax Style:
+
+-   [seminr-alternative-models](demo/seminr-alternative-models.R): Reuse
+    measurement and structural components in multiple models
+-   [seminr-style-contained](demo/seminr-style-contained.R): Create and
+    execute a SEM model in one function call
 
 ## Sister Projects
 

demo/00Index

@@ -2,6 +2,7 @@ seminr-style-contained   A demonstration of the contained style of SEMinR syntax
 seminr-plsc-ecsi   A demonstration of Consistent PLS estimation of the ECSI model
 seminr-pls-ecsi   A demonstration of the regular style of SEMinR syntax
 seminr-pls-interaction   A demonstration of modeling an interaction using SEMinR syntax
+seminr-pls-mga   A demonstration of PLS-MGA procedure to assess subgroup differences
 seminr-alternative-models   A demonstration of alternate model specification, quick and easy
 seminr-pls-higher_order   A demonstration of modeling a higher order construct in SEMinR syntax
 seminr-cbsem-cfa-ecsi   A demonstration of the CBSEM estimation of the ECSI model using SEMinR syntax