RSAmodel {RSAtools}

Estimate polynomial models for response surface analysis

Description

Estimates 37 polynomial families ("STEP1") or user-specific polynomial models ("USER") for use in response surface analysis, based on a 3-step identification strategy (Núñez-Regueiro & Juhel, 2022, 2024).

Usage

RSAmodel(
  formula,
  data = NULL,
  center = "none",
  scale = "none",
  na.rm = FALSE,
  out.rm = TRUE,
  breakline = FALSE,
  models = c("CUBIC", "STEP1"),
  user_model = NULL,
  verbose = TRUE,
  add = "",
  estimator = "MLR",
  se = "robust",
  missing = NA,
  control.variables = NULL,
  center.control.variables = FALSE,
  sampling.weights = NULL,
  group_name = NULL,
  cluster = NULL,
  ...
)

Arguments

formula

A formula in the form z ~ x*y, specifying the variable names used from the data frame, where z is the name of the response variable, and x and y are the names of the predictor variables.

data

A data frame with the variables

center

Method for centering the predictor variables before the analysis. Default option ("none") applies no centering. "pooled" centers the predictor variables on their pooled sample mean. "variablewise" centers the predictor variables on their respective sample mean. You should think carefully before applying the "variablewise" option, as centering or reducing the predictor variables on common values (e.g., their grand means and SDs) can affect the commensurability of the predictor scales.

scale

Method for scaling the predictor variables before the analysis. Default option ("none") applies no scaling. "pooled" scales the predictor variables on their pooled sample SD, which preserves the commensurability of the predictor scales. "variablewise" scales the predictor variables on their respective sample SD. You should think carefully before applying the "variablewise" option, as scaling the predictor variables at different values (e.g., their respective SDs) can affect the commensurability of the predictor scales.

na.rm

Remove missings before proceeding?

out.rm

Should outliers according to Bollen & Jackman (1980) criteria be excluded from the analyses? In large data sets this analysis is the speed bottleneck. If you are sure that no outliers exist, set this option to FALSE for speed improvements.

breakline

Should the breakline in the unconstrained absolute difference model be allowed (the breakline is possible from the model formulation, but empirically rather unrealistic ...). Defaults to FALSE

models

A vector with names of all models that should be computed. Should be any from c("CUBIC","FM1_ONLYX","FM2_ONLYY","FM3_ADDITIVE","FM4_INTER","FM5_QUADX","FM6_QUADY","FM7_CONG","FM8_INCONG","FM9_CURVCONGX","FM10_CURVCONGY","FM11_CURVINCONGX","FM12_CURVINCONGY","FM13_QUADXQUADY","FM14_ROTCONG","FM15_ROTINCONG","FM16_CUBICX","FM17_CUBICY","FM18_LEVDEPQUADX","FM19_LEVDEPQUADY","FM20_ASYMCONG","FM21_ASYMINCONG","FM22_LEVDEPCONG","FM23_LEVDEPINCONG","FM24_PARALLELASYM","FM25_NONPARALLELASYM","FM26_PARALLELASYMWEAK","FM27_PARALLELASYMSTRONG","FM28_NONPARALLELASYMWEAK","FM29_NONPARALLELASYMSTRONG","FM30_ASYMCONGROTY","FM31_ASYMCONGROTX","FM32_ASYMINCONGROTY","FM33_ASYMINCONGROTX","FM34_LEVDEPCONGROTY","FM35_LEVDEPCONGROTX","FM36_LEVDEPINCONGROTY","FM37_LEVDEPINCONGROTX"). For models="STEP1", all polynomial families and the saturated cubic are computed (default), for models="USER" all user-specific models defined in the list "user_model" are computed.

user_model

A list of user-specified polynomial models, defined by setting constraints on the polynomial parameters b1 to b9, using lavaan syntax (e.g., "b1 == 2*b2"). Only parametric constraints are allowed.

verbose

Should additional information during the computation process be printed?

add

Additional syntax that is added to the lavaan model. Can contain, for example, additional constraints, like "p01 == 0; p11 == 0"

estimator

Type of estimator that should be used by lavaan. Defaults to "MLR", which provides robust standard errors, a robust scaled test statistic, and can handle missing values. If you want to reproduce standard OLS estimates, use estimator="ML" and se="standard"

se

Type of standard errors. This parameter gets passed through to the sem function of the lavaan package. See options there. By default, robust SEs are computed. If you use se="boot", lavaan provides CIs and p-values based on the bootstrapped standard error. If you use confint(..., method="boot"), in contrast, you get CIs and p-values based on percentile bootstrap.

missing

Handling of missing values (this parameter is passed to the lavaan sem function). By default (missing=NA), Full Information Maximum Likelihood (FIML) is employed in case of missing values. If cases with missing values should be excluded, use missing = "listwise".

control.variables

A string vector with variable names from data. These variables are added as linear predictors to the model (in order "to control for them"). No interactions with the other variables are modeled.

center.control.variables

Should the control variables be centered before analyses? This can improve interpretability of the intercept, which will then reflect the predicted outcome value at the point (X,Y)=(0,0) when all control variables take their respective average values.

sampling.weights

Name of variable containing sampling weights. Needs to be added here (not in ...) to be included in the analysis dataset.

group_name

Name of variable defining groups, for multigroup modeling.

cluster

Name of variable for clusters, for cluster-level variance correction.

...

Additional parameters passed to the sem function.

Details

This function implements a comparative framework for identifying best-fitting RSA solutions (Núñez-Regueiro & Juhel, 2022, 2024). The default feature ("STEP1") involves the comparison of 37 polynomial families predefined by parametric constraints, to identify likely candidates for best-fitting solution. Step 2 involves probing variants within the retained best-fitting family, by testing user-specific constraints ("USER") on lower-order polynomials that do not define the family. Step 3 (optional) can be conducted on the final variant by bootstrapping validation (using bootstrapLavaan(RSA_object$models$name_final, FUN="coef"), see examples) or cross-validation data.

Value

A list of objects containing polynomials models and names of variables

References

Núñez-Regueiro, F., Juhel, J. (2022). Model-Building Strategies in Response Surface Analysis. Manuscript submitted for publication.

Núñez-Regueiro, F., Juhel, J. (2024). Response Surface Analysis for the Social Sciences I: Identifying Best-Fitting Polynomial Solutions. Manuscript submitted for publication.

See Also

plotting.rsa, plotting.ext

Examples

###Simulation data from RSAtools
summary(sim_NSfit)

##### STEP 1: Estimate and compare polynomial families
##NB: to estimate all families, use "STEP1" in the "model" option
RSA_step1 <-  RSAmodel(engagement ~ needs*supplies,
data= sim_NSfit, model= c("CUBIC","FM20_ASYMCONG",
"FM21_ASYMINCONG","FM26_PARALLELASYMWEAK"))
names(RSA_step1$models)
##Compare solutions according to parsimony (wAIC=Akaike weight)
RSA_step1_fit <- best.rsa(RSA_step1,order=c("wAIC"))
head(RSA_step1_fit)
##Inspect best-fitting family model
summary(RSA_step1$models$FM26_PARALLELASYMWEAK)

##### STEP 2: Test variant within best-fitting family
##Define variant as constraints
list_variant <- list()
list_variant[["variant1"]] <- c('
####First-order polynomials: variant-specific			
			b1 == -1/4*b2
####Second-order polynomials: variant-specific			
			b5 == b3/2
			b4 == 0
####Third-order polynomials: define family FM26			
			b6 == 0
			b7 == 0
			b9 == b8/-3		
')
RSA_NSfit_Step2  <- RSAmodel(formula= engagement ~ needs*supplies,
data= sim_NSfit, model= c("FM26_PARALLELASYMWEAK","USER"),
user_model= list_variant)
##Compare variant to best-fitting family (e.g., LRT_pvalue p > .05)
best.rsa2(RSA_NSfit_Step2,m1="variant1",m2="FM26_PARALLELASYMWEAK")[2,1:3]

######PROBE RESPONSE SURFACE EXTREMA: see ident.ext()

######PLOT RESPONSE SURFACE: see plotting.rsa, plotting.ext