psychmeta is a comprehensive and open-source R package for computing psychometric meta-analyses, conducting follow-up analyses of meta-analytic results (e.g., cumulative meta-analyses, leave-one-out meta-analyses, bootstrapped meta-analyses, and meta-regressions [with help from metafor]), simulating data impacted by psychometric artifacts, and much more! This tutorial is meant to introduce you to some of the most critical features of psychmeta and familiarize you with workflows within the package. We’ll start with the basics (installation and data import) and work our way up to the more techncial topics.

# Installing and loading psychmeta

Install the official CRAN version of the package:

install.packages("psychmeta")

Or install the GitHub verson (with incremental updates between CRAN releases):

devtools::install_github("psychmeta/psychmeta")

library(psychmeta)

# Data import

Recommended packages for data import:

• readr (many efficient functions for loading delimited datafiles; e.g., CSV)
• haven (imports and exports SPSS, SAS, and Stata files)
• readxl (for Excel files)

Other packages for data import:

• rio (convenient import of many formats; relies on fread() funtion from data.table and can have isues with date variables)
• xlsx (for Excel files; slower than readxl)
• foreign (imports and exports many statistical formats; e.g., SPSS, SAS, Stata; slower than haven)

psychmeta requires data to be in a “long” format. This means that each effect size is coded in its own row of a matrix. For example, a long format looks like this:

# Data conversion and cleaning

Not all data come in a long format, so data manipulation tools can be helpful to convert it to a long format.

Recommended packages for data manipulation:

• dplyr (flexible package for fast, intuitive data manipulation)
• tidyr (functions to convert between wide and long data)

Other packages for data manipulation:

• reshape (outdated package, replaced by tidyr)
• reshape2 (outdated package, replaced by tidyr)

psychmeta can also help with reshaping data!

Sometimes it can be quick to code data in a journal table-like format:

dat_matrix

This matrix-style format can be easily converted to the long-format data psychmeta requires with the reshape_mat2dat() function:

reshape_mat2dat(
var_names = var_names,                # Column of variable names
cor_data = c("X", "Y", "X"),          # Names of correlation columns
common_data = n,                      # Names of columns shared among relationships
unique_data = c("mean", "sd", "rel"), # Names of columns unique to relationships
data = dat_matrix)

Sometimes it can also be quick to code data with multiple correlations in one row (i.e., a wide format):

dat_wide

This wide format can be easily converted to the long format psychmeta requires:

common_vars <- c("sample_id")           # Column names for variables common to all
# relationships
var_names <- c("X", "Y", "Z")
es_design = matrix(NA, 3, 3)            # Matrix containing the column names
es_design[lower.tri(es_design)] <-      # for the intercorrelations among variables
c("rxyi_X_Y", "rxyi_X_Z", "rxyi_Y_Z") # in the lower triangle of the matrix
rownames(es_design) <-
colnames(es_design) <-
var_names
n_design <- "ni"                        # Sample size column name or es_design-like
# matrix
other_design <-                         # Matrix with variable names as row names,
cbind(rel = c("rel_X",                # names of long-format variables as column names,
"rel_Y",               # and column names of dat_wide as elements
"rel_Z"))
rownames(other_design) <- var_names

reshape_wide2long(common_vars = common_vars,
es_design = es_design,
n_design = n_design,
other_design = other_design,
es_name = "rxyi",              # Type of effect size in dat_wide
data = dat_wide)

The convert_es() function can be used to convert a wide variety of statistics to the r or d metric.

convert_es(es = 1,    input_es = "d",       output_es = "r", n1 = 50,  n2 = 50)
convert_es(es = -1.3, input_es = "t",       output_es = "r", n1 = 100, n2 = 140)
convert_es(es = 10.3, input_es = "F",       output_es = "r", n1 = 100, n2 = 150)
convert_es(es = 1.3,  input_es = "chisq",   output_es = "r", n1 = 100, n2 = 100)
convert_es(es = .021, input_es = "p.chisq", output_es = "r", n1 = 100, n2 = 100)
convert_es(es = 4.37, input_es = "or",      output_es = "r", n1 = 100, n2 = 100)
convert_es(es = 1.47, input_es = "lor",     output_es = "r", n1 = 100, n2 = 100)

convert_es(es = .2,   input_es = "r",       output_es = "d", n1 = 50,  n2 = 50)
convert_es(es = -1.3, input_es = "t",       output_es = "d", n1 = 100, n2 = 140)
convert_es(es = 10.3, input_es = "F",       output_es = "d", n1 = 100, n2 = 150)
convert_es(es = 1.3,  input_es = "chisq",   output_es = "d", n1 = 100, n2 = 100)
convert_es(es = .021, input_es = "p.chisq", output_es = "d", n1 = 100, n2 = 100)
convert_es(es = 4.37, input_es = "or",      output_es = "d", n1 = 100, n2 = 100)
convert_es(es = 1.47, input_es = "lor",     output_es = "d", n1 = 100, n2 = 100)

To convert effect sizes for inclusion in a meta-analysis, pull out the meta_input object:

convert_es(es = c(.4, .3, .25),
input_es = "r", output_es = "d",
n1 = c(50, 110, 65), n2 = c(50, 70, 65)
out_plots$forest$analysis id: 1$barebones # The sensitivity_cumulative() and sensitivity_leave1out() functions also produce plots. out_plots$leave1out$analysis id: 1$barebones$plots out_plots$cumulative$analysis id: 1$barebones$plots ## Funnel plots have been added to 'ma_obj' - use get_plots() to retrieve them. ## Generating reference lists psychmeta can automatically generate a formatted reference list using a BibTeX or Zotero library. First, specify the citekey argument when running ma_r(). Then use the generate_bib() function: ma_obj <- ma_r(ma_method = "ic", rxyi = rxyi, n = n, construct_x = x_name, construct_y = y_name, rxx = rxxi, ryy = ryyi, moderators = moderator, clean_artifacts = FALSE, impute_artifacts = FALSE, citekey = citekey, data = data_r_meas_multi) generate_bib(ma_obj, # The location of your .bib file (normally should be in the same directory # as the analysis script). bib = system.file("templates/sample_bibliography.bib", package="psychmeta"), # Your citation style. Must be the style ID for a style hosted at # http://zotero.org/styles/ style = "apa", # What format to write output as (options are: "word", "html", "pdf", "text", # "Rmd", "biblatex", "citekeys")? output_format = "word", # What filename to output to (use "console" to print to the R console). file = "sources.docx" ) You can also generate a bibliography containing only sources contributing to a specific subset of the meta-analyses: generate_bib(ma_obj, bib = system.file("templates/sample_bibliography.bib", package="psychmeta"), analyses = list(construct_pair = list(c("X", "Y"))) ) # Supplemental tools ## Computing meta-analyses with pre-existing artifact distributions Even if no studies in your meta-analytic database provide artifact information, pre-specified artifact distributions from previous meta-analyses can still be used in psychmeta! There are several ways to use pre-specified distributions, ranging from using them on their on their own to combining them with observed artifact information. Below are four potential methods, but psychmeta allows many other ways to use artifact distributions. 1. You can use the create_ad() function to create artifact distribution objects, which can then be used with the ma_r_ad() and ma_d_ad() functions to correct bare-bones meta-analyses. # Create artifact-distribution objects for X, Y, and Z. ad_x <- create_ad(mean_qxi = 0.8927818, var_qxi = 0.0008095520, k_qxi = 40, mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha") ad_y <- create_ad(mean_qxi = 0.8927818, var_qxi = 0.0008095520, k_qxi = 40, mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha") ad_z <- create_ad(mean_qxi = 0.8962108, var_qxi = 0.0007840593, k_qxi = 40, mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha") # Compute a meta-analysis of X and Y. ma_bb_xy <- ma_r_bb(r = rxyi, n = n, data = filter(data_r_meas_multi, x_name == "X", y_name == "Y")) # Correct the meta-analysis for measurement error. ma_r_ad(ma_obj = ma_bb_xy, ad_obj_x = ad_x, ad_obj_y = ad_y, correct_rr_x = FALSE, correct_rr_y = FALSE) 1. Rather than applying a two-step process (first compute bare-bones, then apply corrections), you can use the ma_r() function to automate your workflow! Just supply the artifact distributions as a list with the supplemental_ads argument. ma_r(ma_method = "ad", rxyi = rxyi, n = n, correct_rr_x = FALSE, correct_rr_y = FALSE, construct_x = x_name, construct_y = y_name, sample_id = sample_id, clean_artifacts = FALSE, impute_artifacts = FALSE, moderators = moderator, data = data_r_meas_multi, # The "supplemental_ads" argument can take a list of artifact-distribution objects. supplemental_ads = list(X = ad_x, Y = ad_y, Z = ad_z)) 1. You can bypass the manual creation of artifact distribution objects and have psychmeta do it for you. Instead of running create_ad(), you can create a list of artifact information to use as the supplemental_ads argument in ma_r() or ma_d(). This list contains artifact information for each construct you want to supply information for. This is a “list of lists” - each construct for which you provide information gets its own list labeled with the construct name. Each construct’s list should contain entries that are named like the arguments to the create_ad() function. ma_r(ma_method = "ad", rxyi = rxyi, n = n, correct_rr_x = FALSE, correct_rr_y = FALSE, construct_x = x_name, construct_y = y_name, sample_id = sample_id, clean_artifacts = FALSE, impute_artifacts = FALSE, moderators = moderator, data = data_r_meas_multi, supplemental_ads = # The "supplemental_ads" argument can also take raw artifact values. # (These values are unrounded so our examples all produce identical results.) list(X = list(mean_qxi = 0.8927818, var_qxi = 0.0008095520, k_qxi = 40, mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha"), Y = list(mean_qxi = 0.8941266, var_qxi = 0.0009367234, k_qxi = 40, mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha"), Z = list(mean_qxi = 0.8962108, var_qxi = 0.0007840593, k_qxi = 40, mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha"))) 1. You can also create a list of artifact distributions using the create_ad_list() function. This list can then be passed to ma_r() as before. ad_list <- create_ad_list(n = n, rxx = rxxi, ryy = ryyi, construct_x = x_name, construct_y = y_name, sample_id = sample_id, data = data_r_meas_multi) ma_r(ma_method = "ad", rxyi = rxyi, n = n, correct_rr_x = FALSE, correct_rr_y = FALSE, construct_x = x_name, construct_y = y_name, sample_id = sample_id, clean_artifacts = FALSE, impute_artifacts = FALSE, moderators = moderator, data = data_r_meas_multi, # The "supplemental_ads" argument can take the output of create_ad_list(). supplemental_ads = ad_list) 1. It’s possible to mix-and-match any of the above methods, as well as use artifacts from your database. # For purposes of illustration, delete all reliability values not associated with "Z": dat_zarts <- data_r_meas_multi dat_zarts[,"rxxi"] <- dat_zarts[dat_zarts$y_name != "Z","ryyi"] <- NA

# Compute a meta-analysis using three different types of artifact formats:
ma_r(ma_method = "ad", rxyi = rxyi, n = n,
# Observed artifacts can be provided along with the "supplemental_ads" argument.
ryy = ryyi,
correct_rr_x = FALSE, correct_rr_y = FALSE,
construct_x = x_name, construct_y = y_name, sample_id = sample_id,
clean_artifacts = FALSE, impute_artifacts = FALSE,
moderators = moderator, data = dat_zarts,
Y = list(mean_qxi = 0.8941266, var_qxi = 0.0009367234, k_qxi = 40,
mean_n_qxi = 11927 / 40, qxi_dist_type = "alpha")))

# Simulations

In addition to its collection of tools for computing meta-analyses, psychmeta also boasts a robust suite of simulation functions. These functions are intended to support high-quality meta-analytic simulations by correctly introducing sampling error, measurement error, and range restriction into a set of parameter values. These functions can be valuable for both research and pedagogical applications.

## Psychometric data sets

You can generate data impacted by psychometric artifacts for a single sample by using the simulate_psych() function. This function reports three datasets per sample: Observed scores, true scores, and error scores. It can also perform selection and indicate (1) which cases would be selected on the basis of observed scores (i.e., actual selection), (2) which cases would be selected on the basis of true scores (i.e., error-free selection), and (3) which cases would be selected on the basis of error scores (i.e., random selection).

sim_dat <- simulate_psych(n = 1000,
rho_mat = reshape_vec2mat(.5),
sigma_vec = c(1, 1),
sr_vec = c(1, .5),
rel_vec = c(.8, .8), var_names = c("Y", "X"))

## Psychometric correlations

You can automate the process of analyzing a simulated sample when you are interested in effect sizes and artifact estimates. For example, you can generate correlations impacted by psychometric artifacts for a single sample of data using the simulate_r_sample() function. This function generates sample statistics when the n argument is a finite value, but it can also analyze parameters when n = Inf!

simulate_r_sample(n = 1000,
# Correlation parameter matrix
rho_mat = reshape_vec2mat(.5, order = 5),
# Reliability parameter vector
rel_vec = rep(.8, 5),
# Selection ratio vector
sr_vec = c(1, 1, 1, 1, .5),
# Number of items in each scale
k_items_vec = 1:5,
# Matrix of weights to use in creating composites
wt_mat = cbind(c(0, 0, 0, .3, 1),
c(1, .3, 0, 0, 0)),
# Selection ratios for composites
sr_composites = c(.7, .5))

simulate_r_sample(n = Inf,
rho_mat = reshape_vec2mat(.5, order = 5),
rel_vec = rep(.8, 5),
sr_vec = c(1, 1, 1, 1, .5),
k_items_vec = 1:5,
wt_mat = cbind(c(0, 0, 0, .3, 1),
c(1, .3, 0, 0, 0)),
sr_composites = c(.7, .5))

One of the most powerful simulation functions in psychmeta is simulate_r_database(), which can generate entire databases of correlations for use in meta-analysis simulations. This function is a wrapper for simulate_r_sample() and it can sample parameter values from whichever distributions you specify.

# Note the varying methods for defining parameters:
simulate_r_database(k = 10,

# Sample-size parameters
# Parameters can be defined as functions.
n_params = function(n) rgamma(n, shape = 100),

# Correlation parameters (one parameter distribution per correlation)
# Parameters can also be defined as vectors...
rho_params = list(c(.1, .3, .5),
# ...as a mean + a standard deviation...
c(mean = .3, sd = .05),
# ...or as a matrix of values and weights.
rbind(value = c(.1, .3, .5),
weight = c(1, 2, 1))),

# Reliability parameters
rel_params = list(c(.7, .8, .9),
c(mean = .8, sd = .05),
rbind(value = c(.7, .8, .9),
weight = c(1, 2, 1))),

# Selection-ratio parameters
sr_params = list(1, 1, c(.5, .7)),

# Measure-length parameters
k_items_params = list(5, 8, 10),

# Composite weight parameters
wt_params = list(list(1, 1, 0)),

# Selection-ratio parameters for composites
sr_composite_params = list(1),

# Variable names
var_names = c("X", "Y", "Z"))

## Psychometric d values

psychmeta includes d-value simulation functions that parallel the correlation simulation functions described above. These functions automate the simulation of subgroup comparisons and can introduce sampling error, measurement error, and range-restriction artifacts into d-value parameters.

To simulate a single sample of results, use the simulate_d_sample() function. This function can generate data for any number of variables with any number of groups and it can form composite variables and perform selection.

## Simulate statistics by providing integers as "n_vec":
simulate_d_sample(n_vec = c(200, 100),

# List of rho matrices - one per group
rho_mat_list = list(reshape_vec2mat(.5),
reshape_vec2mat(.4)),

# Matrix of group means (groups on rows, variables on columns)
mu_mat = rbind(c(1, .5),
c(0, 0)),

# Matrix of group SDs (groups on rows, variables on columns)
sigma_mat = rbind(c(1, 1),
c(1, 1)),

# Matrix of group reliabilities (groups on rows, variables on columns)
rel_mat = rbind(c(.8, .7),
c(.7, .7)),

# Vector of selection ratios
sr_vec = c(1, .5),

# Number of items in each scale
k_items_vec = c(5, 10),

# Group names
group_names = c("A", "B"))

Like simulate_r_sample(), simulate_d_sample() can also analyze parameters. To get parameter estimates, simply define the sample sizes as proportions (with the n_vec argument) rather than as integers.

simulate_d_sample(n_vec = c(2/3, 1/3),
rho_mat_list = list(reshape_vec2mat(.5),
reshape_vec2mat(.4)),
mu_mat = rbind(c(1, .5),
c(0, 0)),
sigma_mat = rbind(c(1, 1),
c(1, 1)),
rel_mat = rbind(c(.8, .7),
c(.7, .7)),
sr_vec = c(1, .5),
k_items_vec = c(5, 10),
group_names = c("A", "B"))

The simulate_d_database() function is the d-value counterpart to simulate_r_database(). However, each parameter argument now requires a set of parameters for each group (e.g., a distribution of correlations, means, reliabilities, etc. must be specified for each group).

simulate_d_database(k = 5,
# "Group1" and "Group2" labels are not required for parameter arguments:
# They are shown only for enhanced interpretability.

# Sample-size parameter distributions
n_params = list(Group1 = c(mean = 200, sd = 20),
Group2 = c(mean = 100, sd = 20)),

# Correlation parameter distributions
rho_params = list(Group1 = list(c(.3, .4, .5)),
Group2 = list(c(.3, .4, .5))),

# Mean parameter distributions
mu_params = list(Group1 = list(c(mean = .5, sd = .5),
c(-.5, 0, .5)),
Group2 = list(c(mean = 0, sd = .5),
c(-.2, 0, .2))),

# Reliability parameter distributions
rel_params = list(Group1 = list(.8, .8),
Group2 = list(c(.7, .8, .8),
c(.7, .8, .8))),

# Group names
group_names = c("Group1", "Group2"),

# Variable names
var_names = c("Y", "Z"))

# Multivariate range-restriction corrections

In addition to its array of univariate and bivariate range-restriction corrections, psychmeta includes functions to correct data for multivariate range restriction. The correct_matrix_mvrr() corrects covariance matrices for range restriction and the correct_means_mvrr() function corrects vectors of means for range restriction.

# Create example matrices with different variances and covariances:
Sigma_i <- reshape_vec2mat(cov = .2, var = .8, order = 4)
Sigma_xx_a <- reshape_vec2mat(cov = .5, order = 2)

# Correct "Sigma_i" for range restriction in variables 1 & 2:
correct_matrix_mvrr(Sigma_i = Sigma_i,
Sigma_xx_a = Sigma_xx_a,
x_col = 1:2)

# Correct the means of variables 1 & 2 for range restriction:
correct_means_mvrr(Sigma = Sigma_i,
means_i = c(2, 2, 1, 1),
means_x_a = c(1, 1),
x_col = 1:2)

# Composites

## Spearman-Brown reliability estimates

The Spearman-Brown formula can be used to estimate the reliability of a measure that has been lengthened (or shortened) k times. This formula is implemented in the estimate_rel_sb() function and its inverse (which estimates the number of times a measure needs to be lengthened to achieve a desired reliability) is implemented in the estimate_length_sb() function.

estimate_rel_sb(rel_initial = .5, k = 4)

estimate_length_sb(rel_initial = .5, rel_desired = .8)

## Correlations between a variable and a composite or between two composites

The composite_r_matrix() and composite_r_scalar() functions can compute composite correlations using matrix and scalar inputs, respectively. Both functions are capable of computing correlations between a composite and a single variable or between two composites.

# Correlation between a variable and a composite
composite_r_matrix(r_mat = reshape_vec2mat(.4, order = 5), x_col = 2:5, y_col = 1)

composite_r_scalar(mean_rxy = .4, k_vars_x = 4, mean_intercor_x = .4)

# Correlation between two composites
composite_r_matrix(r_mat = reshape_vec2mat(.3, order = 5), x_col = 1:3, y_col = 4:5)

composite_r_scalar(mean_rxy = .3,
k_vars_x = 3, mean_intercor_x = .3, k_vars_y = 2, mean_intercor_y = .3)

## d values of composites

The composite_d_matrix() and composite_d_scalar() functions can compute composite d values using matrix and scalar inputs, respectively.

composite_d_matrix(d_vec = c(1, 1),
r_mat = reshape_vec2mat(.7),
wt_vec = c(1, 1),
p = .5)

composite_d_scalar(mean_d = 1, mean_intercor = .7, k_vars = 2, p = .5)

## Reliability of composites

The composite_rel_matrix() and composite_rel_scalar() functions can compute composite reliabilities (i.e., Mosier reliabilities) using matrix and scalar inputs, respectively.

composite_rel_matrix(rel_vec = c(.8, .8), r_mat = reshape_vec2mat(.4), sd_vec = c(1, 1))

composite_rel_scalar(mean_rel = .8, mean_intercor = .4, k_vars = 2)

# Matrix regression

The lm_mat() function mimics the structure and appearance of the lm() function from the stats package, but uses correlation/covariance matrices to compute regression models instead of raw data sets.

# For purposes of demonstration, generate a dataset from the following matrix:
S <- reshape_vec2mat(cov = c(.3 * 2 * 3,
.4 * 2 * 4,
.5 * 3 * 4),
var = c(2, 3, 4)^2,
var_names = c("X", "Y", "Z"))
mean_vec <- setNames(c(1, 2, 3), colnames(S))
dat <- data.frame(MASS::mvrnorm(n = 100, mu = mean_vec, Sigma = S, empirical = TRUE))

The same formulas can be used with lm_mat() as with lm(), you just have to provide cov_mat and n arguments (and an optional mean_vec argument) instead of data. However, lm_mat() can only compute moderated regressions if the product term is already in cov_mat.

# Compute regression models using one predictor:
lm_out1 <- lm(formula = Y ~ X, data = dat)
lm_mat_out1 <- lm_mat(formula = Y ~ X, cov_mat = S, mean_vec = mean_vec, n = nrow(dat))

# Compute regression models using two predictors:
lm_out2 <- lm(formula = Y ~ X + Z, data = dat)
lm_mat_out2 <- lm_mat(formula = Y ~ X + Z, cov_mat = S, mean_vec = mean_vec, n = nrow(dat))

The summary() method works for both lm() and lm_mat() results and the printouts look identical.

summary(lm_out1)
summary(lm_mat_out1)

lm_mat() results can also be used with functions such as anova(), predict(), and confint() to yield output comparable to what one would get with lm() results.

anova(lm_out1, lm_out2)
anova(lm_mat_out1, lm_mat_out2)

# Matrix regression applied to meta-analytic correlations

The get_matrix() function can extract correlation matrices from meta-analysis objects for use in subsequent multivariate analyses.

mat_array <- get_matrix(ma_obj = ma_obj)
R <- mat_array$individual_correction$moderator_comb: 1$true_score$mean_rho
n <- mat_array$individual_correction$moderator_comb: 1$true_score$N[1,3]

The lm_mat() function can be used to compute regressions from meta-analytic matrices, using the mean N or harmonic mean N as a sample size.

lm_x <- lm_mat(formula = Y ~ X, cov_mat = R, n = n)
lm_xz <- lm_mat(formula = Y ~ X + Z, cov_mat = R, n = n)

Meta-analytic incremental validity analyses are simple with the anova() method applied to lm_mat() results.

summary(lm_x)
summary(lm_xz)
anova(lm_x, lm_xz)

# Best practices for publishing with psychmeta:

1. Include your R script as an appendix or online supplement.
2. Include your data table as a CSV in an online supplement.
3. Cite psychmeta.
• To cite psychmeta, please reference our in-press software announcement in Applied Psychological Measurement (Dahlke & Wiernik, in press, see the bibliography below for reference information).

# Bibliography

Dahlke, J. A., & Wiernik, B. M. (2018). One of these artifacts is not like the others: Accounting for indirect range restriction in organizational and psychological research. Manuscript submitted for publication.

Dahlke, J. A., & Wiernik, B. M. (in press). psychmeta: An R package for psychometric meta-analysis. Applied Psychological Measurement.

Gonzalez-Mulé, E., Mount, M. K., & Oh, I.-S. (2014). A meta-analysis of the relationship between general mental ability and nontask performance. Journal of Applied Psychology, 99(6), 1222–1243. https://doi.org/10.1037/a0037547

Schmidt, F. L. (2017). Statistical and measurement pitfalls in the use of meta-regression in meta-analysis. Career Development International, 22(5), 469–476. https://doi.org/10/gcnjf4

Wiernik, B. M., Kostal, J. W., Wilmot, M. P., Dilchert, S., & Ones, D. S. (2017). Empirical benchmarks for interpreting effect size variability in meta-analysis. Industrial and Organizational Psychology, 10(3), 472–479. https://doi.org/10/ccnv