Type:	Package
Title:	A Collection of R Functions by the Petersen Lab
Version:	1.1.0
Maintainer:	Isaac T. Petersen <isaac-t-petersen@uiowa.edu>
Depends:	R (≥ 4.1.0)
Imports:	stats, graphics, utils, nlme, Hmisc, digest, dplyr, ggplot2, lavaan, mitools, mix, mvtnorm, psych, stringr, xtable, grDevices, plyr, reshape2, RColorBrewer, viridisLite, tidyselect, scales, purrr
Suggests:	testthat (≥ 3.0.0), waldo, withr
Description:	A collection of R functions that are widely used by the Petersen Lab. Included are functions for various purposes, including evaluating the accuracy of judgments and predictions, performing scoring of assessments, generating correlation matrices, conversion of data between various types, data management, psychometric evaluation, extensions related to latent variable modeling, various plotting capabilities, and other miscellaneous useful functions. By making the package available, we hope to make our methods reproducible and replicable by others and to help others perform their data processing and analysis methods more easily and efficiently. The codebase is provided in Petersen (2024) <doi:10.5281/zenodo.7602890> and on CRAN: <doi:10.32614/CRAN.package.petersenlab>. The package is described in "Principles of Psychological Assessment: With Applied Examples in R" (Petersen, 2024, 2025) <doi:10.1201/9781003357421>, <doi:10.25820/work.007199>, <doi:10.5281/zenodo.6466589>.
URL:	https://github.com/DevPsyLab/petersenlab, https://devpsylab.github.io/petersenlab/
BugReports:	https://github.com/DevPsyLab/petersenlab/issues
License:	MIT + file LICENSE
Encoding:	UTF-8
RoxygenNote:	7.3.2
Config/testthat/edition:	3
NeedsCompilation:	no
Packaged:	2025-02-02 23:03:16 UTC; Isaac
Author:	Isaac T. Petersen [aut, cre], Developmental Psychopathology Lab at the University of Iowa [ctb], Angela D. Staples [ctb], Johanna Caskey [ctb], Philipp Doebler [ctb], Loreen Sabel [ctb]
Repository:	CRAN
Date/Publication:	2025-02-03 08:00:02 UTC

NOTIN Operator.

Description

NOTIN operator.

Usage

x %ni% table

Arguments

Details

Determine whether values in one vector are not in another vector.

Value

Vector of TRUE and FALSE, indicating whether values in one vector are not in another vector.

See Also

https://www.r-bloggers.com/2018/07/the-notin-operator/ https://stackoverflow.com/questions/71309487/r-package-documentation-undocumented-arguments-in-documentation-object-for-a?noredirect=1

Examples

# Prepare Data
v1 <- c("Sally","Tom","Barry","Alice")
listToCheckAgainst <- c("Tom","Alice")

v1 %ni% listToCheckAgainst
v1[v1 %ni% listToCheckAgainst]

Statistical Mode.

Description

Calculate statistical mode.

Usage

Mode(x, multipleModes = "all")

Arguments

Details

Calculates statistical mode(s).

Value

Statistical mode(s).

See Also

https://stackoverflow.com/questions/2547402/how-to-find-the-statistical-mode/8189441#8189441

Other computations: kish_ess(), meanSum(), mySum()

Examples

# Prepare Data
v1 <- c(1, 1, 2, 2, 3)

#Calculate Statistical Mode
Mode(v1)
Mode(v1, multipleModes = "mean")
Mode(v1, multipleModes = "first")

Accuracy at a Given Cutoff.

Description

Find the accuracy at a given cutoff. Actuals should be binary, where 1 = present and 0 = absent.

Usage

accuracyAtCutoff(
  predicted,
  actual,
  cutoff,
  UH = NULL,
  UM = NULL,
  UCR = NULL,
  UFA = NULL
)

Arguments

Details

Compute accuracy indices of predicted values in relation to actual values at a given cutoff by specifying the predicted values, actual values, and cutoff value. The target condition is considered present at or above the cutoff value. Optionally, you can also specify the utility of hits, misses, correct rejections, and false alarms to calculate the overall utility of the cutoff. To compute accuracy at each possible cutoff, see accuracyAtEachCutoff.

Value

• cutoff = the cutoff specified

• TP = true positives

• TN = true negatives

• FP = false positives

• FN = false negatives

• SR = selection ratio

• BR = base rate

• percentAccuracy = percent accuracy

• percentAccuracyByChance = percent accuracy by chance

• percentAccuracyPredictingFromBaseRate = percent accuracy from predicting from the base rate

• RIOC = relative improvement over chance

• relativeImprovementOverPredictingFromBaseRate = relative improvement over predicting from the base rate

• SN = sensitivty

• SP = specificity

• TPrate = true positive rate

• TNrate = true negative rate

• FNrate = false negative rate

• FPrate = false positive rate

• HR = hit rate

• FAR = false alarm rate

• PPV = positive predictive value

• NPV = negative predictive value

• FDR = false discovery rate

• FOR = false omission rate

• youdenJ = Youden's J statistic

• balancedAccuracy = balanced accuracy

• f1Score = F1-score

• mcc = Matthews correlation coefficient

• diagnosticOddsRatio = diagnostic odds ratio

• positiveLikelihoodRatio = positive likelihood ratio

• negativeLikelhoodRatio = negative likelihood ratio

• dPrimeSDT = d-Prime index from signal detection theory

• betaSDT = beta index from signal detection theory

• cSDT = c index from signal detection theory

• aSDT = a index from signal detection theory

• bSDT = b index from signal detection theory

• differenceBetweenPredictedAndObserved = difference between predicted and observed values

• informationGain = information gain

• overallUtility = overall utility (if utilities were specified)

See Also

Other accuracy: accuracyAtEachCutoff(), accuracyOverall(), nomogrammer(), optimalCutoff(), posttestOdds()

Examples

# Prepare Data
data("USArrests")
USArrests$highMurderState <- NA
USArrests$highMurderState[which(USArrests$Murder >= 10)] <- 1
USArrests$highMurderState[which(USArrests$Murder < 10)] <- 0

# Calculate Accuracy
accuracyAtCutoff(predicted = USArrests$Assault,
  actual = USArrests$highMurderState, cutoff = 200)
accuracyAtCutoff(predicted = USArrests$Assault,
  actual = USArrests$highMurderState, cutoff = 200,
  UH = 1, UM = 0, UCR = .9, UFA = 0)

Accuracy at Each Cutoff.

Description

Find the accuracy at each possible cutoff. Actuals should be binary, where 1 = present and 0 = absent.

Usage

accuracyAtEachCutoff(
  predicted,
  actual,
  UH = NULL,
  UM = NULL,
  UCR = NULL,
  UFA = NULL
)

Arguments

Details

Compute accuracy indices of predicted values in relation to actual values at each possible cutoff by specifying the predicted values and actual values. The target condition is considered present at or above each cutoff value. Optionally, you can specify the utility of hits, misses, correct rejections, and false alarms to calculate the overall utility of each possible cutoff.

Value

• cutoff = the cutoff specified

• TP = true positives

• TN = true negatives

• FP = false positives

• FN = false negatives

• SR = selection ratio

• BR = base rate

• percentAccuracy = percent accuracy

• percentAccuracyByChance = percent accuracy by chance

• percentAccuracyPredictingFromBaseRate = percent accuracy from predicting from the base rate

• RIOC = relative improvement over chance

• relativeImprovementOverPredictingFromBaseRate = relative improvement over predicting from the base rate

• SN = sensitivty

• SP = specificity

• TPrate = true positive rate

• TNrate = true negative rate

• FNrate = false negative rate

• FPrate = false positive rate

• HR = hit rate

• FAR = false alarm rate

• PPV = positive predictive value

• NPV = negative predictive value

• FDR = false discovery rate

• FOR = false omission rate

• youdenJ = Youden's J statistic

• balancedAccuracy = balanced accuracy

• f1Score = F1-score

• mcc = Matthews correlation coefficient

• diagnosticOddsRatio = diagnostic odds ratio

• positiveLikelihoodRatio = positive likelihood ratio

• negativeLikelhoodRatio = negative likelihood ratio

• dPrimeSDT = d-Prime index from signal detection theory

• betaSDT = beta index from signal detection theory

• cSDT = c index from signal detection theory

• aSDT = a index from signal detection theory

• bSDT = b index from signal detection theory

• differenceBetweenPredictedAndObserved = difference between predicted and observed values

• informationGain = information gain

• overallUtility = overall utility (if utilities were specified)

See Also

Other accuracy: accuracyAtCutoff(), accuracyOverall(), nomogrammer(), optimalCutoff(), posttestOdds()

Examples

# Prepare Data
data("USArrests")
USArrests$highMurderState <- NA
USArrests$highMurderState[which(USArrests$Murder >= 10)] <- 1
USArrests$highMurderState[which(USArrests$Murder < 10)] <- 0

# Calculate Accuracy
accuracyAtEachCutoff(predicted = USArrests$Assault,
  actual = USArrests$highMurderState)
accuracyAtEachCutoff(predicted = USArrests$Assault,
  actual = USArrests$highMurderState,
  UH = 1, UM = 0, UCR = .9, UFA = 0)

Overall Accuracy.

Description

Find overall accuracy.

Usage

accuracyOverall(predicted, actual, dropUndefined = FALSE)

wisdomOfCrowd(predicted, actual, dropUndefined = FALSE)

Arguments

Details

Compute overall accuracy estimates of predicted values in relation to actual values. Estimates of overall accuracy span all cutoffs. Some accuracy estimates can be undefined under various circumstances. Optionally, you can drop undefined values in the calculation of accuracy indices. Note that dropping undefined values changes the meaning of these indices. Use this option at your own risk!

Value

• ME = mean error

• MAE = mean absolute error

• MSE = mean squared error

• RMSE = root mean squared error

• MPE = mean percentage error

• MAPE = mean absolute percentage error

• sMAPE = symmetric mean absolute percentage error

• MASE = mean absolute scaled error

• RMSLE = root mean squared log error

• rsquared = R-squared

• rsquaredAdj = adjusted R-squared

• rsquaredPredictive = predictive R-squared

See Also

Mean absolute scaled error (MASE):
https://stats.stackexchange.com/questions/108734/alternative-to-mape-when-the-data-is-not-a-time-series
https://stats.stackexchange.com/questions/322276/is-mase-specified-only-to-time-series-data
https://stackoverflow.com/questions/31197726/calculate-mase-with-cross-sectional-non-time-series-data-in-r
https://stats.stackexchange.com/questions/401759/how-can-mase-mean-absolute-scaled-error-score-value-be-interpreted-for-non-tim

Predictive R-squared:
https://www.r-bloggers.com/2014/05/can-we-do-better-than-r-squared/

Other accuracy: accuracyAtCutoff(), accuracyAtEachCutoff(), nomogrammer(), optimalCutoff(), posttestOdds()

Examples

# Prepare Data
data("USArrests")

# Calculate Accuracy
accuracyOverall(predicted = USArrests$Assault, actual = USArrests$Murder)
wisdomOfCrowd(predicted = USArrests$Assault, actual = 200)

Add Correlation to Scatterplot.

Description

Add correlation text to scatterplot.

Usage

addText(
  x,
  y,
  xcoord = NULL,
  ycoord = NULL,
  size = 1,
  col = NULL,
  method = "pearson"
)

Arguments

Details

Adds a correlation coefficient and associated p-value to a scatterplot.

Value

Correlation coefficient, degrees of freedom, and p-value printed on scatterplot.

See Also

Other plot: plot2WayInteraction(), ppPlot(), semPlotInteraction(), vwReg()

Other correlations: cor.table(), crossTimeCorrelation(), crossTimeCorrelationDF(), partialcor.table(), vwReg()

Examples

# Prepare Data
data("USArrests")

# Scatterplot
plot(USArrests$Assault, USArrests$Murder)
addText(x = USArrests$Assault, y = USArrests$Murder)

APA Format

Description

Format decimals and leading zeroes. Adapted from the MOTE package.

Usage

apa(value, decimals = 3, leading = TRUE)

Arguments

Details

Formats decimals and leading zeroes for creating reports in scientific style, to be consistent with American Psychological Association (APA) format. This function creates "pretty" character vectors from numeric variables for printing as part of a report. The value can take a single number, matrix, vector, or multiple columns from a data frame, as long as they are numeric. The values will be coerced into numeric if they are characters or logical values, but this process may result in an error if values are truly alphabetical.

Value

Value(s) in the format specified, with the number of decimals places indicated and with or without a leading zero, as indicated.

See Also

https://github.com/doomlab/MOTE

Other formatting: pValue(), specify_decimal(), suppressLeadingZero()

Examples

apa(value = 0.54674, decimals = 3, leading = TRUE)

Attenuation of True Correlation Due to Measurement Error.

Description

Estimate the observed association between the predictor and criterion after accounting for the degree to which a true correlation is attenuated due to measurement error.

Usage

attenuationCorrelation(
  trueAssociation,
  reliabilityOfPredictor,
  reliabilityOfCriterion
)

Arguments

Details

Estimate the association that would be observed between the predictor and criterion after accounting for the degree to which a true correlation is attenuated due to random measurement error (unreliability).

Value

Observed correlation between predictor and criterion.

See Also

Other correlation: disattenuationCorrelation()

Examples

attenuationCorrelation(
  trueAssociation = .7,
  reliabilityOfPredictor = .9,
  reliabilityOfCriterion = .85)

Clean Up Player Names For Merging.

Description

Cleans up names of players for merging.

Usage

cleanUpNames(name)

Arguments

Details

Cleans up names of NFL Football players, including making them all-caps, removing common suffixes, punctuation, spaces, etc. This is helpful for merging multiple datasets.

Value

Vector of cleaned player names.

Examples

oldNames <- c("Peyton Manning","Tom Brady","Marvin Harrison Jr.")
cleanNames <- cleanUpNames(oldNames)
cleanNames

Column Bind and Fill.

Description

Column bind dataframes and fill with NAs.

Usage

columnBindFill(...)

Arguments

Details

Binds columns of two or more dataframes together, and fills in missing rows.

Value

Dataframe with columns binded together.

See Also

https://stackoverflow.com/questions/7962267/cbind-a-dataframe-with-an-empty-dataframe-cbind-fill/7962286#7962286

Other dataManipulation: convert.magic(), dropColsWithAllNA(), dropRowsWithAllNA(), varsDifferentTypes()

Examples

# Prepare Data
df1 <- data.frame(a = rnorm(5), b = rnorm(5))
df2 <- data.frame(c = rnorm(4), d = rnorm(4))

# Column Bind and Fill
columnBindFill(df1, df2)

Simulate Complement Variable.

Description

Simulate data with a specified correlation in relation to an existing variable.

Usage

complement(y, rho, x)

Arguments

Details

Simulates data with a specified correlation in relation to an existing variable.

Value

Vector of a variable that has a specified correlation in relation to a given variable y.

See Also

https://stats.stackexchange.com/a/313138/20338

Other simulation: simulateAUC(), simulateIndirectEffect()

Examples

v1 <- rnorm(100)
complement(y = v1, rho = .5)
complement(y = v1, rho = -.5)

v2 <- complement(y = v1, rho = .85)
plot(v1, v2)

Convert Variable Types.

Description

Converts variable types of multiple columns of a dataframe at once.

Usage

convert.magic(obj, type)

Arguments

Details

Converts variable types of multiple columns of a dataframe at once. Convert variable types to character, numeric, or factor.

Value

Dataframe with columns converted to a particular type.

See Also

https://stackoverflow.com/questions/11261399/function-for-converting-dataframe-column-type/11263399#11263399

Other dataManipulation: columnBindFill(), dropColsWithAllNA(), dropRowsWithAllNA(), varsDifferentTypes()

Other conversion: convertHoursAMPM(), convertToHours(), convertToMinutes(), convertToSeconds(), percentileToTScore(), pom()

Examples

# Prepare Data
data("USArrests")

# Convert variables to character
convert.magic(USArrests, "character")

Convert AM and PM Hours.

Description

Convert hours to 24-hour time.

Usage

convertHoursAMPM(hours, ampm, am = 0, pm = 1, treatMorningAsLate = FALSE)

Arguments

Details

Convert hours to the number of hours in 24-hour time. You can specify whether to treat morning hours (e.g., 1 AM) as late (25 H), e.g., for specifying late bedtimes

Value

Hours in 24-hour-time.

See Also

Other times: convertToHours(), convertToMinutes(), convertToSeconds()

Other conversion: convert.magic(), convertToHours(), convertToMinutes(), convertToSeconds(), percentileToTScore(), pom()

Examples

# Prepare Data
df1 <- data.frame(hours = c(1, 1, 12, 12), ampm = c(0, 0, 1, 1))
df2 <- data.frame(hours = c(1, 1, 12, 12), ampm = c(1, 1, 0, 0))

# Convert AM and PM hours
convertHoursAMPM(hours = df1$hours, ampm = df1$ampm)
convertHoursAMPM(hours = df1$hours, ampm = df1$ampm,
  treatMorningAsLate = TRUE)

convertHoursAMPM(hours = df2$hours, ampm = df2$ampm, am = 1, pm = 0)
convertHoursAMPM(hours = df2$hours, ampm = df2$ampm, am = 1, pm = 0,
  treatMorningAsLate = TRUE)

Convert Time to Hours.

Description

Convert times to hours.

Usage

convertToHours(hours, minutes, seconds, HHMMSS, HHMM)

Arguments

Details

Converts times to hours. To convert times to minutes or seconds, see convertToMinutes or convertToSeconds.

Value

Vector of times in hours.

See Also

Other times: convertHoursAMPM(), convertToMinutes(), convertToSeconds()

Other conversion: convert.magic(), convertHoursAMPM(), convertToMinutes(), convertToSeconds(), percentileToTScore(), pom()

Examples

# Prepare Data
df <- data.frame(hours = c(0,1), minutes = c(15,27), seconds = c(30,13),
  HHMMSS = c("00:15:30","01:27:13"), HHMM = c("00:15","01:27"))

# Convert to Hours
convertToHours(hours = df$hours, minutes = df$minutes, seconds = df$seconds)
convertToHours(HHMMSS = df$HHMMSS)
convertToHours(HHMM = df$HHMM)

Convert Time to Minutes.

Description

Convert times to minutes.

Usage

convertToMinutes(hours, minutes, seconds, HHMMSS, HHMM, MMSS)

Arguments

Details

Converts times to minutes. To convert times to hours or seconds, see convertToHours or convertToSeconds.

Value

Vector of times in minutes.

See Also

Other times: convertHoursAMPM(), convertToHours(), convertToSeconds()

Other conversion: convert.magic(), convertHoursAMPM(), convertToHours(), convertToSeconds(), percentileToTScore(), pom()

Examples

# Prepare Data
df <- data.frame(hours = c(0,1), minutes = c(15,27), seconds = c(30,13),
  HHMMSS = c("00:15:30","01:27:13"), HHMM = c("00:15","01:27"))

# Convert to Minutes
convertToMinutes(hours = df$hours, minutes = df$minutes,
  seconds = df$seconds)
convertToMinutes(HHMMSS = df$HHMMSS)
convertToMinutes(HHMM = df$HHMM)

Convert Time to Seconds.

Description

Convert times to seconds.

Usage

convertToSeconds(hours, minutes, seconds, HHMMSS, HHMM, MMSS)

Arguments

Details

Converts times to seconds. To convert times to hours or minutes, see convertToHours or convertToMinutes.

Value

Vector of times in seconds.

See Also

Other times: convertHoursAMPM(), convertToHours(), convertToMinutes()

Other conversion: convert.magic(), convertHoursAMPM(), convertToHours(), convertToMinutes(), percentileToTScore(), pom()

Examples

# Prepare Data
df <- data.frame(hours = c(0,1), minutes = c(15,27), seconds = c(30,13),
  HHMMSS = c("00:15:30","01:27:13"), HHMM = c("00:15","01:27"),
  MMSS = c("15:30","87:13"))

# Convert to Minutes
convertToSeconds(hours = df$hours, minutes = df$minutes,
  seconds = df$seconds)
convertToSeconds(HHMMSS = df$HHMMSS)
convertToSeconds(HHMM = df$HHMM)
convertToSeconds(MMSS = df$MMSS)

Correlation Matrix.

Description

Function that creates a correlation matrix similar to SPSS output.

Usage

cor.table(x, y, type = "none", dig = 2, correlation = "pearson")

Arguments

Details

Co-created by Angela Staples (astaples@emich.edu) and Isaac Petersen (isaac-t-petersen@uiowa.edu). For a partial correlation matrix, see partialcor.table.

Value

A correlation matrix.

See Also

Other correlations: addText(), crossTimeCorrelation(), crossTimeCorrelationDF(), partialcor.table(), vwReg()

Examples

# Prepare Data
data("mtcars")

# Correlation Matrix
cor.table(mtcars[,c("mpg","cyl","disp")])
cor.table(mtcars[,c("mpg","cyl","disp")])
cor.table(mtcars[,c("mpg","cyl","disp")], dig = 3)
cor.table(mtcars[,c("mpg","cyl","disp")], dig = 3, correlation = "spearman")

cor.table(mtcars[,c("mpg","cyl","disp")], type = "manuscript", dig = 3)
cor.table(mtcars[,c("mpg","cyl","disp")], type = "manuscriptBig")

table1 <- cor.table(mtcars[,c("mpg","cyl","disp")], type = "latex")
table2 <- cor.table(mtcars[,c("mpg","cyl","disp")], type = "latexSPSS")
table3 <- cor.table(mtcars[,c("mpg","cyl","disp")], type = "manuscriptLatex")
table4 <- cor.table(mtcars[,c("mpg","cyl","disp")], type = "manuscriptBigLatex")

cor.table(mtcars[,c("mpg","cyl","disp")], mtcars[,c("drat","qsec")])
cor.table(mtcars[,c("mpg","cyl","disp")], mtcars[,c("drat","qsec")], type = "manuscript", dig = 3)

Cross-Time Correlations.

Description

Calculate the association of a variable across multiple time points.

Usage

crossTimeCorrelation(id = "tcid", time = "age", variable, data)

Arguments

Details

Calculate the association of a variable across multiple time points. It is especially useful when there are three or more time points.

Value

output of cor.test()

See Also

Other correlations: addText(), cor.table(), crossTimeCorrelationDF(), partialcor.table(), vwReg()

Examples

# Prepare Data
df <- expand.grid(ID = 1:100, time = c(1, 2, 3))
df <- df[order(df$ID),]
row.names(df) <- NULL
df$score <- rnorm(nrow(df))

# Cross-Time Correlation
crossTimeCorrelation(id = "ID", time = "time", variable = "score", data = df)

Cross-Time Correlations Dataframe.

Description

Dataframe used to compute cross-time correlations.

Usage

crossTimeCorrelationDF(id = "tcid", time = "age", variable, data)

Arguments

Details

Dataframe used to calculate the association of a variable across multiple time points. It is especially useful when there are three or more time points.

Value

dataframe with three columns in the form of: ID, time1, time2

See Also

Other correlations: addText(), cor.table(), crossTimeCorrelation(), partialcor.table(), vwReg()

Examples

# Prepare Data
df <- expand.grid(ID = 1:100, time = c(1, 2, 3))
df <- df[order(df$ID),]
row.names(df) <- NULL
df$score <- rnorm(nrow(df))

# Cross-Time Correlation
crossTimeCorrelationDF(id = "ID", time = "time", variable = "score", data = df)

Item and Test Information from Zero-Inflated Negative Binomial Model.

Description

Estimate item and test information from Bayesian zero-inflated negative binomial model that was fit using the brms package.

Usage

deriv_d_negBinom(n, alpha, beta, theta, phi)

d_negBinom(n, alpha, beta, theta, phi)

log_gen_binom(n, phi)

deriv_logd_negBinom(n, alpha, beta, theta, phi)

info_neg_binom_analytical(
  theta = seq(-2.5, 2.5, length.out = 101),
  alpha,
  beta,
  phi,
  varpi
)

item_info_NB_zero_analytical(theta, alpha, beta, phi, varpi)

item_info_NB_analytical(theta, alpha, beta, phi, varpi = NULL, zero = FALSE)

test_info_NB(theta, alpha, beta, phi, varpi = NULL, zero = FALSE)

error_variance_NB(
  lower = -Inf,
  upper = Inf,
  alpha,
  beta,
  phi,
  varpi = NULL,
  zero = FALSE
)

reliability_NB(alpha, beta, phi, varpi = NULL, zero = FALSE)

Arguments

Details

Created by Philipp Doebler (doebler@statistik.tu-dortmund.de) and Loreen Sabel (loreen.sabel@tu-dortmund.de).

Value

The amount of information for a given item (or the test as a whole) at each of the values of theta specified. Based on test information, one can estimate error variance and marginal reliability using error_variance_NB() and reliability_NB(), respectively.

See Also

Other bayesian: pA()

Other IRT: discriminationToFactorLoading(), fourPL(), itemInformation(), reliabilityIRT(), standardErrorIRT(), test_info_4PL()

Examples

## Not run: 
library("brms")
library("rstan")

coef_bayesianMixedEffectsGRM_gam <- coef(bayesianMixedEffectsGRM_gam)
str(coef_bayesianMixedEffectsGRM_gam)
itempars <- coef_bayesianMixedEffectsGRM_gam$item[,1,1:4]

# define a grid of thetas for the computations:
theta_seq <- seq(-4, 4, length.out = 201)

# item information for all items
# The resulting matrix has length(theta_seq) columns and a row per item.
# We use a loop for the calculations
item_info <- matrix(NA, nrow = nrow(itempars), ncol = length(theta_seq))
for(i in 1:nrow(itempars)){
  item_info[i, ] <- item_info_NB_zero_analytical(
    theta = theta_seq,
    alpha = itempars[i, "alpha_Intercept"],
    beta = itempars[i, "beta_Intercept"],
    phi = exp(itempars[i, "shape_Intercept"]),
    varpi = plogis(itempars[i, "zi_Intercept"]))
}

test_info <- data.frame(
  theta = theta_seq,
  testInformation = colSums(item_info)
)

# Or, alternatively:
test_info_NB(
  theta = compareTestInfo$theta,
  alpha = itempars[,"alpha_Intercept"],
  beta = itempars[,"beta_Intercept"],
  phi = exp(itempars[,"shape_Intercept"]),
  varpi = plogis(itempars[,"zi_Intercept"]),
  zero = TRUE)

# Test Standard Error of Measurement in Different Theta Ranges
error_variance_NB(
  lower = -4,
  upper = 4,
  alpha = itempars[,"alpha_Intercept"],
  beta = itempars[,"beta_Intercept"],
  phi = exp(itempars[,"shape_Intercept"]),
  varpi = plogis(itempars[,"zi_Intercept"]),
  zero = TRUE
)

error_variance_NB(
  lower = -4,
  upper = 0,
  alpha = itempars[,"alpha_Intercept"],
  beta = itempars[,"beta_Intercept"],
  phi = exp(itempars[,"shape_Intercept"]),
  varpi = plogis(itempars[,"zi_Intercept"]),
  zero = TRUE
)

error_variance_NB(
  lower = 0,
  upper = 1.5,
  alpha = itempars[,"alpha_Intercept"],
  beta = itempars[,"beta_Intercept"],
  phi = exp(itempars[,"shape_Intercept"]),
  varpi = plogis(itempars[,"zi_Intercept"]),
  zero = TRUE
)

error_variance_NB(
  lower = 1.5,
  upper = 4,
  alpha = itempars[,"alpha_Intercept"],
  beta = itempars[,"beta_Intercept"],
  phi = exp(itempars[,"shape_Intercept"]),
  varpi = plogis(itempars[,"zi_Intercept"]),
  zero = TRUE
)

# One-Number Summary of Test Reliability
reliability_NB(
  alpha = itempars[,"alpha_Intercept"],
  beta = itempars[,"beta_Intercept"],
  phi = exp(itempars[,"shape_Intercept"]),
  varpi = plogis(itempars[,"zi_Intercept"]),
  zero = TRUE)

## End(Not run)

Disattenuation of Observed Correlation Due to Measurement Error.

Description

Estimate the true association between the predictor and criterion after accounting for the degree to which a true correlation is attenuated due to measurement error.

Usage

disattenuationCorrelation(
  observedAssociation,
  reliabilityOfPredictor,
  reliabilityOfCriterion
)

Arguments

Details

Estimate the true association between the predictor and criterion after accounting for the degree to which a true correlation is attenuated due to random measurement error (unreliability).

Value

True association between predictor and criterion.

See Also

Other correlation: attenuationCorrelation()

Examples

disattenuationCorrelation(
  observedAssociation = .7,
  reliabilityOfPredictor = .9,
  reliabilityOfCriterion = .85)

Discrimination (IRT) to Standardized Factor Loading.

Description

Convert a discrimination parameter in item response theory to a standardized factor loading.

Usage

discriminationToFactorLoading(a, model = "probit")

Arguments

Details

Convert a discrimination parameter in item response theory to a standardized factor loading

Value

Standardized factor loading.

See Also

https://aidenloe.github.io/introToIRT.html https://stats.stackexchange.com/questions/228629/conversion-of-irt-logit-discrimination-parameter-to-factor-loading-metric

Other IRT: deriv_d_negBinom(), fourPL(), itemInformation(), reliabilityIRT(), standardErrorIRT(), test_info_4PL()

Examples

discriminationToFactorLoading(0.5)
discriminationToFactorLoading(1.3)
discriminationToFactorLoading(1.3, model = "logit")

Drop NA columns.

Description

Drop columns with all missing (NA) values.

Usage

dropColsWithAllNA(data, ignore = NULL)

Arguments

Details

Drop columns that have no observed values, i.e., all values in the column are missing (NA), excluding the ignored columns.

Value

A dataframe with columns removed that had all missing values in non-ignored columns.

See Also

Other dataManipulation: columnBindFill(), convert.magic(), dropRowsWithAllNA(), varsDifferentTypes()

Other dataEvaluations: dropRowsWithAllNA(), is.nan.data.frame(), not_all_na(), not_any_na()

Examples

# Prepare Data
df <- expand.grid(ID = 1:100, time = c(1, 2, 3), rater = c(1, 2),
naCol1 = NA, naCol2 = NA)
df <- df[order(df$ID),]
row.names(df) <- NULL
df$score1 <- rnorm(nrow(df))
df$score2 <- rnorm(nrow(df))
df$score3 <- rnorm(nrow(df))
df[sample(1:nrow(df), size = 100), c("score1","score2","score3")] <- NA

# Drop Rows with All NA in Non-Ignored Columns
dropColsWithAllNA(df)
dropColsWithAllNA(df, ignore = c("naCol2"))

Drop NA rows.

Description

Drop rows with all missing (NA) values.

Usage

dropRowsWithAllNA(data, ignore = NULL)

Arguments

Details

Drop rows that have no observed values, i.e., all values in the row are missing (NA), excluding the ignored columns.

Value

A dataframe with rows removed that had all missing values in non-ignored columns.

See Also

Other dataManipulation: columnBindFill(), convert.magic(), dropColsWithAllNA(), varsDifferentTypes()

Other dataEvaluations: dropColsWithAllNA(), is.nan.data.frame(), not_all_na(), not_any_na()

Examples

# Prepare Data
df <- expand.grid(ID = 1:100, time = c(1, 2, 3))
df <- df[order(df$ID),]
row.names(df) <- NULL
df$score1 <- rnorm(nrow(df))
df$score2 <- rnorm(nrow(df))
df$score3 <- rnorm(nrow(df))
df[sample(1:nrow(df), size = 100), c("score1","score2","score3")] <- NA

# Drop Rows with All NA in Non-Ignored Columns
dropRowsWithAllNA(df, ignore = c("ID","time"))

Chi-Square Equivalence Test for Structural Equation Models.

Description

Function that performs a chi-square equivalence test for structural equation models.

Usage

equiv_chi(alpha = 0.05, chi, df, m, N_sample, popRMSEA = 0.08)

Arguments

Details

Created by Counsell et al. (2020): Counsell, A., Cribbie, R. A., & Flora, D. B. (2020). Evaluating equivalence testing methods for measurement invariance. Multivariate Behavioral Research, 55(2), 312-328. https://doi.org/10.1080/00273171.2019.1633617

Value

p-value indicating whether to reject the null hypothesis that the model is a poor fit to the data.

See Also

Other structural equation modeling: make_esem_model(), puc(), satorraBentlerScaledChiSquareDifferenceTestStatistic(), semPlotInteraction()

Examples

# Prepare Data
data("mtcars")

# Fit structural equation model

# Extract statistics
N1 <- 1222
m <- 1
Tml1 <- 408.793
df1 <- 80

# Equivalence test
equiv_chi(alpha = .05, chi = Tml1, df = df1, m = 1, N_sample = N1, popRMSEA = .08)

4-Parameter Logistic Curve.

Description

4-parameter logistic curve for item response theory.

Usage

fourPL(a = 1, b, c = 0, d = 1, theta)

Arguments

Details

Estimates the probability of item endorsement as function of the four-parameter logistic (4PL) curve and the person's level on the construct (theta).

Value

Probability of item endorsement (or expected value on the item).

See Also

doi:10.1177/0146621613475471

Other IRT: deriv_d_negBinom(), discriminationToFactorLoading(), itemInformation(), reliabilityIRT(), standardErrorIRT(), test_info_4PL()

Examples

fourPL(b = 2, theta = -4:4) #1PL
fourPL(b = 2, a = 1.5, theta = -4:4) #2PL
fourPL(b = 2, a = 1.5, c = 0.10, theta = -4:4) #3PL
fourPL(b = 2, a = 1.5, c = 0.10, d = 0.95, theta = -4:4) #4PL

Package Dependencies.

Description

Determine package dependencies.

Usage

getDependencies(packs)

Arguments

Details

Determine which packages depend on a target package (or packages).

Value

Vector of packages that depend on the target package(s).

See Also

https://stackoverflow.com/questions/52929114/install-packages-in-r-without-internet-connection-with-all-dependencies/52935020#52935020

Other packages: load_or_install()

Examples

old <- options("repos")
options(repos = "https://cran.r-project.org")
getDependencies("tidyverse")
options(old)

Combine Results from Mixed Effect Imputation Models.

Description

Function that combines lme results across multiple imputation runs.

Usage

imputationCombine(model, dig = 3)

Arguments

Details

[INSERT].

Value

Summary of model fit and information for mixed effect imputation models.

See Also

Other multipleImputation: imputationModelCompare(), imputationPRV(), lmCombine()

Examples

#INSERT

Compare Mixed Effect Imputation Models.

Description

Function that compares two nested lme() models from multiple imputation using likelihood ratio test.

Usage

imputationModelCompare(model1, model2)

Arguments

Details

[INSERT].

Value

Likelihood ratio test for model comparison of two mixed effect imputation models.

See Also

Other multipleImputation: imputationCombine(), imputationPRV(), lmCombine()

Examples

#INSERT

Proportional Reduction of Variance from Imputation Models.

Description

Calculate the proportional reduction of variance in imputation models.

Usage

imputationPRV(baseline, full, baselineTime = 1, fullTime = 1)

Arguments

Details

[INSERT].

Value

The proportional reduction of variance from a baseline mixed-effects model to a full mixed effects model.

See Also

Other multipleImputation: imputationCombine(), imputationModelCompare(), lmCombine()

Examples

#INSERT

NaN (Not a Number).

Description

Check whether a value is "Not A Number" (NaN) in a dataframe.

Usage

## S3 method for class 'data.frame'
is.nan(x)

Arguments

Details

[INSERT].

Value

TRUE or FALSE, indicating whether values in a dataframe are Not a Number (NA).

See Also

https://stackoverflow.com/questions/18142117/how-to-replace-nan-value-with-zero-in-a-huge-data-frame/18143097#18143097

Other dataEvaluations: dropColsWithAllNA(), dropRowsWithAllNA(), not_all_na(), not_any_na()

Examples

# Prepare Data
df <- data.frame(item1 = rnorm(1000), item2 = rnorm(1000), item3 = rnorm(1000))
df[sample(1:nrow(df), size = 100), c("item1","item2","item3")] <- NaN

# Calculate Missingness-Adjusted Row Sum
is.nan(df)

Item Information.

Description

Item information in item response theory.

Usage

itemInformation(a = 1, b, c = 0, d = 1, theta)

Arguments

Details

Estimates the amount of information provided by a given item as function of the item parameters and the person's level on the construct (theta).

Value

Amount of item information.

See Also

doi:10.1177/0146621613475471

Other IRT: deriv_d_negBinom(), discriminationToFactorLoading(), fourPL(), reliabilityIRT(), standardErrorIRT(), test_info_4PL()

Examples

itemInformation(b = 2, theta = -4:4) #1PL
itemInformation(b = 2, a = 1.5, theta = -4:4) #2PL
itemInformation(b = 2, a = 1.5, c = 0.10, theta = -4:4) #3PL
itemInformation(b = 2, a = 1.5, c = 0.10, d = 0.95, theta = -4:4) #4PL

Weighted Quantiles.

Description

Computes weighted quantiles. whdquantile() uses a weighted Harrell-Davis quantile estimator. wthdquantile() uses a weighted trimmed Harrell-Davis quantile estimator. wquantile() uses a weighted traditional quantile estimator.

Usage

kish_ess(w)

wquantile_generic(x, w, probs, cdf)

whdquantile(x, w, probs)

wthdquantile(x, w, probs, width = 1/sqrt(kish_ess(w)))

wquantile(x, w, probs, type = 7)

Arguments

Details

Computes weighted quantiles according to Akinshin (2023).

Value

Numeric vector of specified quantiles.

See Also

doi:10.48550/arXiv.2304.07265

Other computations: Mode(), meanSum(), mySum()

Examples

mydata <- c(1:100, 1000)
mydataWithNAs <- mydata
mydataWithNAs[c(1,5,7)] <- NA
weights <- rep(1, length(mydata))
quantiles <- c(0.01, 0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95, 0.99)

whdquantile(
  x = mydata,
  w = weights,
  probs = quantiles)

wthdquantile(
  x = mydata,
  w = weights,
  probs = quantiles)

wquantile(
  x = mydata,
  w = weights,
  probs = quantiles)

whdquantile(
  x = mydataWithNAs,
  w = weights,
  probs = quantiles)

wthdquantile(
  x = mydataWithNAs,
  w = weights,
  probs = quantiles)

wquantile(
  x = mydataWithNAs,
  w = weights,
  probs = quantiles)

Combine Results from Multiple Regression Imputation Models.

Description

Function that combines lm() results across multiple imputation runs.

Usage

lmCombine(model, dig = 3)

Arguments

Details

[INSERT].

Value

Summary of multiple regression imputation models.

See Also

Other multipleImputation: imputationCombine(), imputationModelCompare(), imputationPRV()

Other multipleRegression: plot2WayInteraction(), ppPlot(), semPlotInteraction(), update_nested()

Examples

#INSERT

Summarize mixed effects model.

Description

Summarizes the results of a model fit by the lme() function of the nlme package.

Usage

lmeSummary(model, dig = 3)

Arguments

Details

Summarizes the results of a model fit by the lme() function of the nlme package. Includes summary of parameters, pseudo-r-squared, and whether model is positive definite.

Value

Output summary of lme() model object.

Examples

# Fit Model
library("nlme")
model <- lme(distance ~ age + Sex, data = Orthodont, random = ~ 1 + age)

# Model Summary
summary(model)
lmeSummary(model)

Load or Install Packages.

Description

Loads packages or, if not already installed, installs and loads packages.

Usage

load_or_install(package_names, ...)

Arguments

Details

Loads packages that are already installed, and if the packages are not already installed, it installs and then loads them.

Value

Loaded packages.

See Also

https://www.r-bloggers.com/2012/05/loading-andor-installing-packages-programmatically/ https://stackoverflow.com/questions/4090169/elegant-way-to-check-for-missing-packages-and-install-them

Other packages: getDependencies()

Examples

## Not run: 
old <- options("repos")
options(repos = "https://cran.r-project.org")
# Warning: the command below installs packages that are not already installed
load_or_install(c("tidyverse","nlme"))
options(old)

## End(Not run)

Make ESEM Model.

Description

Make lavaan syntax for exploratory structural equation model (ESEM).

Usage

make_esem_model(loadings, anchors)

Arguments

Details

Makes syntax for exploratory structural equation model (ESEM) to be fit in lavaan.

Value

lavaan model syntax.

See Also

https://msilvestrin.me/post/esem/

Other structural equation modeling: equiv_chi(), puc(), satorraBentlerScaledChiSquareDifferenceTestStatistic(), semPlotInteraction()

Examples

# Prepare Data
data("HolzingerSwineford1939", package = "lavaan")

# Specify EFA Syntax
efa_syntax <- '
  # EFA Factor Loadings
  efa("efa1")*f1 +
  efa("efa1")*f2 +
  efa("efa1")*f3 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9
'

# Fit EFA Model
mplusRotationArgs <- list(rstarts = 30,
  row.weights = "none",
  algorithm = "gpa",
  orthogonal = FALSE,
  jac.init.rot = TRUE,
  std.ov = TRUE, # row standard = correlation
  geomin.epsilon = 0.0001)

efa_fit <- lavaan::sem(
  efa_syntax,
  data = HolzingerSwineford1939,
  information = "observed",
  missing = "ML",
  estimator = "MLR",
  rotation = "geomin",
  # mimic Mplus
  meanstructure = TRUE,
  rotation.args = mplusRotationArgs)

# Extract Factor Loadings
esem_loadings <- lavaan::parameterEstimates(
  efa_fit,
  standardized = TRUE
) |>
  dplyr::filter(efa == "efa1") |>
  dplyr::select(lhs, rhs, est) |>
  dplyr::rename(item = rhs, latent = lhs, loading = est)

# Specify Anchor Item for Each Latent Factor
anchors <- c(f1 = "x3", f2 = "x5", f3 = "x7")

# Generate ESEM Syntax
esemModel_syntax <- make_esem_model(esem_loadings, anchors)

# Fit ESEM Model
lavaan::sem(
  esemModel_syntax,
  data = HolzingerSwineford1939,
  missing = "ML",
  estimator = "MLR")

Mean Sum.

Description

Compute a missingness-adjusted row sum.

Usage

meanSum(x)

Arguments

Details

Take row mean across columns (items) and then multiply by number of items to account for missing (NA) values.

Value

Missingness-adjusted row sum.

See Also

Other computations: Mode(), kish_ess(), mySum()

Examples

# Prepare Data
df <- data.frame(item1 = rnorm(1000), item2 = rnorm(1000), item3 = rnorm(1000))

# Calculate Missingness-Adjusted Row Sum
df$missingnessAdjustedSum <- meanSum(df)

Mortgage Principal and Interest.

Description

Amount of principal and interest payments on a mortgage.

Usage

mortgage(balance, interest, term = 30, n = 12)

Arguments

Details

Calculates the amount of principal and interest payments on a mortgage.

Value

Amount of principal and interest payments.

Examples

mortgage(balance = 300000, interest = .05)
mortgage(balance = 300000, interest = .04)
mortgage(balance = 300000, interest = .06)
mortgage(balance = 300000, interest = .05, term = 15)

My Sum.

Description

Compute a row sum and retain NAs when all values in the row are NA.

Usage

mySum(data)

Arguments

Details

Compute a row sum and set the row sum to be missing (not zero) when all values in the row are missing (NA).

Value

Modified row sum to set row sum to be missing when all values in the row are missing (NA).

See Also

Other computations: Mode(), kish_ess(), meanSum()

Examples

# Prepare Data
df <- data.frame(item1 = rnorm(1000), item2 = rnorm(1000), item3 = rnorm(1000))
df[sample(1:nrow(df), size = 100), c("item1","item2","item3")] <- NA

# Calculate Missingness-Adjusted Row Sum
df$sum <- mySum(df)

Sorts loadings from exploratory factor analysis.

Description

Sorts items' loadings based on their loadings from exploratory factor analysis fit with the psych::fa() function.

Usage

my_loadings_sorter(
  fit,
  sort_type = "largest_loading",
  nchar = 40,
  return_blocks = FALSE,
  showlatentcor = TRUE,
  itemLabels = NULL
)

Arguments

Details

Adapted from code by Philipp Doebler (doebler@statistik.tu-dortmund.de).

Value

Sorted loadings from exploratory factor analysis model.

Create Nomogram.

Description

Create nomogram plot.

Usage

nomogrammer(
  TP = NULL,
  TN = NULL,
  FP = NULL,
  FN = NULL,
  pretestProb = NULL,
  selectionRate = NULL,
  SN = NULL,
  SP = NULL,
  FPR = NULL,
  PLR = NULL,
  NLR = NULL,
  Detail = FALSE,
  NullLine = FALSE,
  LabelSize = (14/5),
  Verbose = FALSE
)

Arguments

Details

Create nomogram plot from the following at a given cut point:

• 1) true positives (TP), true negatives (TN), false positives (FP), and false negatives (FN)

• 2) pretest probability (pretestProb), sensitivity (SN), and specificity (SP), OR

• 3) pretest probability (pretestProb), sensitivity (SN), and false positive rate (FPR), OR

• 4) pretest probability (pretestProb), sensitivity (SN), and selection rate (selectionRate), OR

• 5) pretest probability (pretestProb), positive likelihood ratio (PLR), and negative likelihood ratio (NLR)

Value

ggplot object of nomogram plot.

See Also

https://github.com/achekroud/nomogrammer

Other accuracy: accuracyAtCutoff(), accuracyAtEachCutoff(), accuracyOverall(), optimalCutoff(), posttestOdds()

Examples

nomogrammer(
  TP = 253,
  TN = 386,
  FP = 14,
  FN = 347)

nomogrammer(
  pretestProb = .60,
  SN = 0.421,
  SP = 0.965)

nomogrammer(
  pretestProb = .60,
  SN = 0.421,
  FPR = 0.035)

nomogrammer(
  pretestProb = .60,
  SN = 0.421,
  selectionRate = 0.267)

nomogrammer(
  pretestProb = .60,
  PLR = 12,
  NLR = 0.6)

Any Rows Not NA.

Description

Check if any rows for a column are not NA.

Usage

not_all_na(x)

Arguments

Details

Determine whether any rows for a column (or vector) are not missing (NA).

Value

TRUE or FALSE

See Also

Other dataEvaluations: dropColsWithAllNA(), dropRowsWithAllNA(), is.nan.data.frame(), not_any_na()

Examples

# Prepare Data
data("USArrests")

# Check if any rows are not NA
not_all_na(USArrests$Murder)

Not Any NA.

Description

Check if all rows for a column are NA.

Usage

not_any_na(x)

Arguments

Details

[INSERT].

Value

TRUE or FALSE, indicating whether the whole column does not have any missing values (NA).

See Also

Other dataEvaluations: dropColsWithAllNA(), dropRowsWithAllNA(), is.nan.data.frame(), not_all_na()

Examples

# Prepare Data
df <- data.frame(item1 = rnorm(1000), item2 = rnorm(1000), item3 = rnorm(1000))
df[sample(1:nrow(df), size = 100), "item2"] <- NA
df[,"item3"] <- NA

# Check if Not Any NA
not_any_na(df$item1)
not_any_na(df$item2)
not_any_na(df$item3)

Optimal Cutoff.

Description

Find the optimal cutoff for different aspects of accuracy. Actuals should be binary, where 1 = present and 0 = absent.

Usage

optimalCutoff(predicted, actual, UH = NULL, UM = NULL, UCR = NULL, UFA = NULL)

Arguments

Details

Identify the optimal cutoff for different aspects of accuracy of predicted values in relation to actual values by specifying the predicted values and actual values. Optionally, you can specify the utility of hits, misses, correct rejections, and false alarms to calculate the overall utility of each possible cutoff.

Value

The optimal cutoff and optimal accuracy index at that cutoff based on:

• percentAccuracy = percent accuracy

• percentAccuracyByChance = percent accuracy by chance

• RIOC = relative improvement over chance

• relativeImprovementOverPredictingFromBaseRate = relative improvement over predicting from the base rate

• PPV = positive predictive value

• NPV = negative predictive value

• youdenJ = Youden's J statistic

• balancedAccuracy = balanced accuracy

• f1Score = F1-score

• mcc = Matthews correlation coefficient

• diagnosticOddsRatio = diagnostic odds ratio

• positiveLikelihoodRatio = positive likelihood ratio

• negativeLikelhoodRatio = negative likelihood ratio

• dPrimeSDT = d-Prime index from signal detection theory

• betaSDT = beta index from signal detection theory

• cSDT = c index from signal detection theory

• aSDT = a index from signal detection theory

• bSDT = b index from signal detection theory

• differenceBetweenPredictedAndObserved = difference between predicted and observed values

• informationGain = information gain

• overallUtility = overall utility (if utilities were specified)

See Also

Other accuracy: accuracyAtCutoff(), accuracyAtEachCutoff(), accuracyOverall(), nomogrammer(), posttestOdds()

Examples

# Prepare Data
data("USArrests")
USArrests$highMurderState <- NA
USArrests$highMurderState[which(USArrests$Murder >= 10)] <- 1
USArrests$highMurderState[which(USArrests$Murder < 10)] <- 0

# Determine Optimal Cutoff
optimalCutoff(predicted = USArrests$Assault,
  actual = USArrests$highMurderState)
optimalCutoff(predicted = USArrests$Assault,
  actual = USArrests$highMurderState,
  UH = 1, UM = 0, UCR = .9, UFA = 0)

Bayes' Theorem.

Description

Estimate marginal and conditional probabilities using Bayes theorem.

Usage

pA(pAgivenB, pB, pAgivenNotB)

pB(pBgivenA, pA, pBgivenNotA)

pAgivenB(pBgivenA, pA, pB = NULL, pBgivenNotA = NULL)

pBgivenA(pAgivenB, pB, pA = NULL, pAgivenNotB = NULL)

pAgivenNotB(pAgivenB, pA, pB)

pBgivenNotA(pBgivenA, pA, pB)

Arguments

Details

Estimates marginal or conditional probabilities using Bayes theorem.

Value

The requested marginal or conditional probability. One of:

• the marginal probability of A

• the marginal probability of B

• the conditional probability of A given B

• the conditional probability of B given A

• the conditional probability of A given NOT B

• the conditional probability of B given NOT A

See Also

Other bayesian: deriv_d_negBinom()

Examples

pA(pAgivenB = .95, pB = .285, pAgivenNotB = .007171515)

pB(pBgivenA = .95, pA = .285, pBgivenNotA = .007171515)

pAgivenB(pBgivenA = .95, pA = .285, pB = .2758776)
pAgivenB(pBgivenA = .95, pA = .285, pBgivenNotA = .007171515)
pAgivenB(pBgivenA = .95, pA = .003, pBgivenNotA = .007171515)

pBgivenA(pAgivenB = .95, pB = .285, pA = .2758776)
pBgivenA(pAgivenB = .95, pB = .285, pAgivenNotB = .007171515)
pBgivenA(pAgivenB = .95, pB = .003, pAgivenNotB = .007171515)

pAgivenNotB(pAgivenB = .95, pB = .003, pA = .01)

pBgivenNotA(pBgivenA = .95, pA = .003, pB = .01)

p-values.

Description

Suppress the leading zero when printing p-values.

Usage

pValue(value, digits = 3)

Arguments

Details

[INSERT].

Value

p-value.

See Also

Other formatting: apa(), specify_decimal(), suppressLeadingZero()

Examples

pValue(0.70)
pValue(0.04)
pValue(0.00002)

Partial Correlation Matrix.

Description

Function that creates a partial correlation matrix similar to SPSS output.

Usage

partialcor.table(
  x,
  y,
  z = NULL,
  type = "none",
  dig = 2,
  correlation = "pearson"
)

Arguments

Details

Co-created by Angela Staples (astaples@emich.edu) and Isaac Petersen (isaac-t-petersen@uiowa.edu). Creates a partial correlation matrix, controlling for one or more covariates. For a standard correlation matrix, see cor.table.

Value

A partial correlation matrix.

See Also

Other correlations: addText(), cor.table(), crossTimeCorrelation(), crossTimeCorrelationDF(), vwReg()

Examples

# Prepare Data
data("mtcars")

#Correlation Matrix
partialcor.table(mtcars[,c("mpg","cyl","disp")], z = mtcars$hp)
partialcor.table(mtcars[,c("mpg","cyl","disp")], z = mtcars[,c("hp","wt")])
partialcor.table(mtcars[,c("mpg","cyl","disp")], z = mtcars[,c("hp","wt")],
  dig = 3)
partialcor.table(mtcars[,c("mpg","cyl","disp")], z = mtcars[,c("hp","wt")],
  dig = 3, correlation = "spearman")

partialcor.table(mtcars[,c("mpg","cyl","disp")], z = mtcars[,c("hp","wt")],
  type = "manuscript", dig = 3)
partialcor.table(mtcars[,c("mpg","cyl","disp")], z = mtcars[,c("hp","wt")],
  type = "manuscriptBig")

table1 <- partialcor.table(mtcars[,c("mpg","cyl","disp")],
  z = mtcars[,c("hp","wt")], type = "latex")
table2 <- partialcor.table(mtcars[,c("mpg","cyl","disp")],
  z = mtcars[,c("hp","wt")], type = "latexSPSS")
table3 <- partialcor.table(mtcars[,c("mpg","cyl","disp")],
  z = mtcars[,c("hp","wt")], type = "manuscriptLatex")
table4 <- partialcor.table(mtcars[,c("mpg","cyl","disp")],
  z = mtcars[,c("hp","wt")], type = "manuscriptBigLatex")

partialcor.table(mtcars[,c("mpg","cyl","disp")], mtcars[,c("drat","qsec")],
  mtcars[,c("hp","wt")])
partialcor.table(mtcars[,c("mpg","cyl","disp")], mtcars[,c("drat","qsec")],
  mtcars[,c("hp","wt")], type = "manuscript", dig = 3)

Person Months.

Description

Calculate perons months for personnel effort in grants.

Usage

percentEffort(
  academicMonths = NULL,
  calendarMonths = NULL,
  summerMonths = NULL,
  appointment = 9
)

personMonths(
  academicMonths = NULL,
  calendarMonths = NULL,
  summerMonths = NULL,
  effortAcademic = NULL,
  effortCalendar = NULL,
  effortSummer = NULL,
  appointment = 9
)

Arguments

Details

Calculate person months for personnel effort in grant proposals from academic months, calendar months, and summer months.

Value

The person months of effort.

See Also

https://nexus.od.nih.gov/all/2015/05/27/how-do-you-convert-percent-effort-into-person-months/

Examples

# Specify Values
appointmentDuration <- 9 #(in months)

# Specify either Set 1 (months) or Set 2 (percent effort) below:

#Set 1: Months
academicMonths <- 1.3 #AY (academic year) months (should be between 0 to appointmentDuration)
calendarMonths <- 0 #CY (calendar year) months (should be between 0-12)
summerMonths <- 0.5 #SM (summer) months (should be between 0 to [12-appointmentDuration])

# Set 2: Percent Effort
percentEffortAcademic <- 0.1444444 #(a proportion; should be between 0-1)
percentEffortCalendar <- 0 #(a proportion; should be between 0-1)
percentEffortSummer <- 0.1666667 #(a proportion; should be between 0-1)

# Calculations
summerDuration <- 12 - appointmentDuration

# Percent effort (in proportion)
percentEffort(academicMonths = academicMonths)
percentEffort(calendarMonths = calendarMonths)
percentEffort(summerMonths = summerMonths)

# Person-Months From NIH Website
(percentEffort(academicMonths = academicMonths) * appointmentDuration) +
 (percentEffort(calendarMonths = calendarMonths) * 12) +
 (percentEffort(summerMonths = summerMonths) * summerDuration)

# Person-Months from Academic/Calendar/Summer Months
personMonths(academicMonths = academicMonths,
             calendarMonths = calendarMonths,
             summerMonths = summerMonths)

# Person-Months from Percent Effort
personMonths(effortAcademic = percentEffortAcademic,
             effortCalendar = percentEffortCalendar,
             effortSummer = percentEffortSummer)

Percentile to T-Score Conversion.

Description

Conversion of percentile ranks to T-scores.

Usage

percentileToTScore(percentileRank)

Arguments

Details

Converts percentile ranks to the equivalent T-scores.

Value

Vector of T-scores.

See Also

Other conversion: convert.magic(), convertHoursAMPM(), convertToHours(), convertToMinutes(), convertToSeconds(), pom()

Examples

percentileRanks <- c(1, 10, 20, 30, 40, 50, 60, 70, 80, 90, 99)

percentileToTScore(percentileRanks)

Plot 2-way interaction.

Description

Generates a plot of a 2-way interaction.

Usage

plot2WayInteraction(
  predictor,
  outcome,
  moderator,
  predictorLabel = "predictor",
  outcomeLabel = "outcome",
  moderatorLabel = "moderator",
  varList,
  varTypes,
  values = NA,
  interaction = "normal",
  legendLabels = NA,
  legendLocation = "topright",
  ylim = NA,
  pvalues = TRUE,
  data
)

Arguments

Details

Generates a plot of a 2-way interaction: the association between a predictor and an outcome at two levels of the moderator.

Value

Plot of two-way interaction.

See Also

Other plot: addText(), ppPlot(), semPlotInteraction(), vwReg()

Other multipleRegression: lmCombine(), ppPlot(), semPlotInteraction(), update_nested()

Examples

# Prepare Data
predictor <- rnorm(1000, 10, 3)
moderator <- rnorm(1000, 50, 10)
outcome <- (1.7 * predictor) + (1.3 * moderator) +
  (1.5 * predictor * moderator) + rnorm(1000, sd = 3)
covariate <- rnorm(1000)
df <- data.frame(predictor, moderator, outcome, covariate)

# Linear Regression
lmModel <- lm(outcome ~ predictor + moderator + predictor:moderator,
  data = df, na.action = "na.exclude")
summary(lmModel)

# 1. Plot 2-Way Interaction
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    data = df)

# 2. Specify y-axis Range
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),                                  #new
                    data = df)

# 3. Add Variable Labels
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),
                    predictorLabel = "Stress",                        #new
                    outcomeLabel = "Aggression",                      #new
                    moderatorLabel = "Gender",                        #new
                    data = df)

# 4. Change Legend Labels
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),
                    predictorLabel = "Stress",
                    outcomeLabel = "Aggression",
                    moderatorLabel = "Gender",
                    legendLabels = c("Boys","Girls"),                 #new
                    data = df)

# 5. Move Legend Location
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),
                    predictorLabel = "Stress",
                    outcomeLabel = "Aggression",
                    moderatorLabel = "Gender",
                    legendLabels = c("Boys","Girls"),
                    legendLocation = "topleft",                       #new
                    data = df)

#6. Turn Off p-Values
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),
                    predictorLabel = "Stress",
                    outcomeLabel = "Aggression",
                    moderatorLabel = "Gender",
                    legendLabels = c("Boys","Girls"),
                    legendLocation = "topleft",
                    pvalues = FALSE,                                  #new
                    data = df)

#7. Get Regression Output from Mean-Centered Predictor and Moderator
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),
                    predictorLabel = "Stress",
                    outcomeLabel = "Aggression",
                    moderatorLabel = "Gender",
                    legendLabels = c("Boys","Girls"),
                    legendLocation = "topleft",
                    interaction = "meancenter",                       #new
                    data = df)

#8. Get Regression Output from Orthogonalized Interaction Term
plot2WayInteraction(predictor = "predictor",
                    outcome = "outcome",
                    moderator = "moderator",
                    varList = c("predictor","moderator","covariate"),
                    varTypes = c("sd","binary","mean"),
                    ylim = c(10,50),
                    predictorLabel = "Stress",
                    outcomeLabel = "Aggression",
                    moderatorLabel = "Gender",
                    legendLabels = c("Boys","Girls"),
                    legendLocation = "topleft",
                    interaction = "orthogonalize",                    #new
                    data = df)

Proportion of Maximum (POM).

Description

Calculate the proportion of maximum (POM) score given a minimum and maximum score.

Usage

pom(data, min = NULL, max = NULL)

Arguments

Details

The minimum and maximum score for calculating the proportion of maximum could be the possible or observed minimum and maximum, respectively. Using the possible minimum and maximum would yield the proportion of maximum possible score. Using the observed minimum and maximum would yield the proportion of minimum and maximum observed score. If the minimum and maximum possible scores are not specified, the observed minimum and maximum are used.

Value

Proportion of maximum possible or observed values.

See Also

Other conversion: convert.magic(), convertHoursAMPM(), convertToHours(), convertToMinutes(), convertToSeconds(), percentileToTScore()

Examples

# Prepare Data
v1 <- sample(1:9, size = 1000, replace = TRUE)

# Calculate Proportion of Maximum Possible (by specifying the minimum and maximum possible)
pom(v1, min = 0, max = 10)

# Calculate Proportion of Maximum Observed
pom(v1)

Posttest Odds & Probability.

Description

Estimate posttest odds and posttest probability.

Usage

posttestOdds(
  TP,
  TN,
  FP,
  FN,
  pretestProb = NULL,
  SN = NULL,
  SP = NULL,
  likelihoodRatio = NULL
)

posttestProbability(
  TP,
  TN,
  FP,
  FN,
  pretestProb = NULL,
  SN = NULL,
  SP = NULL,
  likelihoodRatio = NULL
)

Arguments

Details

Estimates posttest odds or posttest probability.

Value

The requested posttest odds or pottest probability.

See Also

Other accuracy: accuracyAtCutoff(), accuracyAtEachCutoff(), accuracyOverall(), nomogrammer(), optimalCutoff()

Examples

posttestOdds(
  TP = 26,
  TN = 56,
  FP = 14,
  FN = 14)

posttestOdds(
  pretestProb = 0.3636364,
  SN = 0.65,
  SP = 0.80)

posttestOdds(
  pretestProb = 0.3636364,
  likelihoodRatio = 3.25)

posttestProbability(
  TP = 26,
  TN = 56,
  FP = 14,
  FN = 14)

posttestProbability(
  pretestProb = 0.3636364,
  SN = 0.65,
  SP = 0.80)

posttestProbability(
  pretestProb = 0.3636364,
  likelihoodRatio = 3.25)

PP Plot.

Description

Normal Probability (P-P) Plot.

Usage

ppPlot(model)

Arguments

Details

A normal probability (P-P) plot compares the empirical cumulative distribution to the theoretical cumulative distribution.

Value

Normal probability (P-P) plot.

See Also

https://www.r-bloggers.com/2009/12/r-tutorial-series-graphic-analysis-of-regression-assumptions/

Other plot: addText(), plot2WayInteraction(), semPlotInteraction(), vwReg()

Other multipleRegression: lmCombine(), plot2WayInteraction(), semPlotInteraction(), update_nested()

Examples

# Prepare Data
predictor1 <- rnorm(100)
predictor2 <- rnorm(100)
outcome <- rnorm(100)

# Fit Model
lmModel <- lm(outcome ~ predictor1 + predictor2)

# P-P Plot
ppPlot(lmModel)

Percent of Uncontaminated Correlations (PUC).

Description

Percent of uncontaminated correlations (PUC) from bifactor model.

Usage

puc(numItems, numSpecificFactors)

Arguments

Details

Estimates the percent of uncontaminated correlations (PUC) from a bifactor model. The PUC represents the percentage of correlations (i.e., covariance terms) that reflect variance from only the general factor (i.e., not variance from a specific factor). Correlations that are explained by the specific factors are considered "contaminated" by multidimensionality.

Value

Percent of Uncontaminated Correlations (PUC).

See Also

doi:10.31234/osf.io/6tf7j doi:10.1177/0013164412449831 doi:10.1037/met0000045

Other structural equation modeling: equiv_chi(), make_esem_model(), satorraBentlerScaledChiSquareDifferenceTestStatistic(), semPlotInteraction()

Examples

puc(
  numItems = 9,
  numSpecificFactors = 3
)

mydata <- data.frame(
  numItems = c(9,18,18,36,36,36),
  numSpecificFactors = c(3,3,6,3,6,12)
)

puc(
  numItems = mydata$numItems,
  numSpecificFactors = mydata$numSpecificFactors
)

Read Encrypted Data.

Description

Read data from encrypted file.

Usage

read.aes(filename, key)

Arguments

Details

Reads data from an encrypted file. To write an data to an encrypted file, see write.aes.

Value

Unencrypted data.

See Also

https://stackoverflow.com/questions/25318800/how-do-i-read-an-encrypted-file-from-disk-with-r/25321586#25321586

Other encrypted: write.aes()

Examples

# Location of Encryption Key on Local Computer (where only you should have access to it)
#encryptionKeyLocation <- file.path(getwd(), "/encryptionKey.RData",
#  fsep = "") #Can change to a different path, e.g.: "C:/Users/[USERNAME]/"

# Generate a Temporary File Path for Encryption Key
encryptionKeyLocation <- tempfile(fileext = ".RData")

# Generate Encryption Key
key <- as.raw(sample(1:16, 16))

# Save Encryption Key
save(key, file = encryptionKeyLocation)

# Specify Credentials
credentials <- "Insert My Credentials Here"

# Generate a Temporary File Path for Encrypted Credentials
encryptedCredentialsLocation <- tempfile(fileext = ".txt")

# Save Encrypted Credentials
#write.aes(
#  df = credentials,
#  filename = file.path(getwd(), "/encrypytedCredentials.txt", fsep = ""),
#  key = key) # Change the file location to save this on the lab drive

write.aes(
  df = credentials,
  filename = encryptedCredentialsLocation,
  key = key)

rm(credentials)
rm(key)

# Read and Unencrypt the Credentials Using the Encryption Key
load(encryptionKeyLocation)

#credentials <- read.aes(
#  filename = file.path(getwd(), "/encrypytedCredentials.txt", fsep = ""),
#  key = key)

credentials <- read.aes(
  filename = encryptedCredentialsLocation,
  key = key)

Recode Intensity.

Description

Recode intensity of behavior based on frequency of behavior.

Usage

recode_intensity(intensity, did_not_occur = NULL, frequency = NULL)

mark_intensity_as_zero(
  item_names,
  data,
  did_not_occur_vars = NULL,
  frequency_vars = NULL
)

Arguments

Details

Recodes the intensity of behavior to zero if the frequency of the behavior is zero (i.e., if the behavior has not occurred).

Value

The intensity of the behavior.

Progress Bar for REDCap.

Description

Function that identifies the values for a progress bar in REDCap.

Usage

redcapProgressBar(numSurveys, beginning = 2, end = 99)

Arguments

Details

A progress bar in REDCap can be created using the following code:

  Progress:
  <div style="width:100%;border:0;margin:0;padding:0;background-color:
  #A9BAD1;text-align:center;"><div style="width:2%;border: 0;margin:0;
  padding:0;background-color:#8491A2"><span style="color:#8491A2">.
  </span></div></div>

where width:2% specifies the progress (out of 100%).

Value

sequence of numbers for the progress bar in REDCap.

Examples

redcapProgressBar(numSurveys = 6)
redcapProgressBar(6)
redcapProgressBar(4)
redcapProgressBar(numSurveys = 7, beginning = 1, end = 99)

Reliability (IRT).

Description

Estimate the reliability in item response theory.

Usage

reliabilityIRT(information, varTheta = 1)

Arguments

Details

Estimate the reliability in item response theory using the test information (i.e., the sum of all items' information).

Value

Reliability for that amount of test information.

See Also

https://groups.google.com/g/mirt-package/c/ZAgpt6nq5V8/m/R3OEeEqdAQAJ

Other IRT: deriv_d_negBinom(), discriminationToFactorLoading(), fourPL(), itemInformation(), standardErrorIRT(), test_info_4PL()

Examples

# Calculate information for 4 items
item1 <- itemInformation(b = -2, a = 0.6, theta = -4:4)
item2 <- itemInformation(b = -1, a = 1.2, theta = -4:4)
item3 <- itemInformation(b = 1, a = 1.5, theta = -4:4)
item4 <- itemInformation(b = 2, a = 2, theta = -4:4)

items <- data.frame(item1, item2, item3, item4)

# Calculate test information
items$testInformation <- rowSums(items)

# Estimate reliability
reliabilityIRT(items$testInformation)

Reliability of Difference Score.

Description

Estimate the reliability of a difference score.

Usage

reliabilityOfDifferenceScore(x, y, reliabilityX, reliabilityY)

Arguments

Details

Estimates the reliability of a difference score.

Value

Reliability of the difference score that is computed from the difference of x and y.

See Also

Other reliability: repeatability()

Examples

v1 <- rnorm(1000, mean = 100, sd = 15)
v2 <- v1 + rnorm(1000, mean = 1, sd = 15)
reliabilityOfDifferenceScore(x = v1, y = v2,
 reliabilityX = .7, reliabilityY = .8)

Repeatability.

Description

Estimate the repeatability of a measure's scores across two time points.

Usage

repeatability(measure1, measure2)

Arguments

Details

Estimates the coefficient of repeatability (CR), bias, and the lower and upper limits of agreement (LOA).

Value

Dataframe with the coefficient of repeatability (CR), bias, the lower limit of agreement (lowerLOA), and the upper limit of agreement (upperLOA). Also generates a Bland-Altman plot with a solid black reference line (indicating a difference of zero), a dashed red line indicating the bias, and dashed blue lines indicating the limits of agreement.

See Also

Other reliability: reliabilityOfDifferenceScore()

Examples

v1 <- rnorm(1000, mean = 100, sd = 15)
v2 <- v1 + rnorm(1000, mean = 1, sd = 3)
repeatability(v1, v2)

Reverse Score Variables.

Description

Reverse score variables using either the theoretical min and max, or the observed max.

Usage

reverse_score(
  data,
  variables,
  theoretical_max = NULL,
  theoretical_min = NULL,
  append_string = NULL
)

Arguments

Details

Reverse scores variables using either the theoretical min and max (by subtracting the theoretical maximum from each score and adding the theoretical minimum to each score) or by subtracting each score from the maximum score for that variable.

Value

Dataframe with reverse-scored variables.

Examples

mydata <- data.frame(
  var1 = c(1, 2, NA, 4, 5),
  var2 = c(NA, 4, 3, 2, 1)
)

variables_to_reverse_score <- c("var1", "var2")

reverse_score(
  mydata,
  variables = variables_to_reverse_score)

reverse_score(
  mydata,
  variables = variables_to_reverse_score,
  append_string = ".R")

reverse_score(
  mydata,
  variables = variables_to_reverse_score,
  theoretical_max = 7)

reverse_score(
  mydata,
  variables = variables_to_reverse_score,
  theoretical_max = 7,
  theoretical_min = 1)

Satorra-Bentler Scaled Chi-Square Difference Test Statistic.

Description

Function that computes the Satorra-Bentler Scaled Chi-Square Difference Test statistic.

Usage

satorraBentlerScaledChiSquareDifferenceTestStatistic(T0, c0, d0, T1, c1, d1)

Arguments

Details

Computes the Satorra-Bentler Scaled Chi-Square Difference Test statistic between two structural equation models.

Value

Satorra-Bentler Scaled Chi-Square Difference Test statistic.

See Also

Other structural equation modeling: equiv_chi(), make_esem_model(), puc(), semPlotInteraction()

Examples

# Fit structural equation model
HS.model <- '
 visual =~ x1 + x2 + x3
 textual =~ x4 + x5 + x6
 speed =~ x7 + x8 + x9
'

fit1 <- lavaan::cfa(HS.model, data = lavaan::HolzingerSwineford1939,
 estimator = "MLR")
fit0 <- lavaan::cfa(HS.model, data = lavaan::HolzingerSwineford1939,
 orthogonal = TRUE, estimator = "MLR")

# Chi-square difference test
# lavaan::anova(fit1, fit0)
satorraBentlerScaledChiSquareDifferenceTestStatistic(
 T0 = lavaan::fitMeasures(fit0)["chisq.scaled"],
 c0 = lavaan::fitMeasures(fit0)["chisq.scaling.factor"],
 d0 = lavaan::fitMeasures(fit0)["df.scaled"],
 T1 = lavaan::fitMeasures(fit1)["chisq.scaled"],
 c1 = lavaan::fitMeasures(fit1)["chisq.scaling.factor"],
 d1 = lavaan::fitMeasures(fit1)["df.scaled"])

Plot interaction from SEM model.

Description

Generates a plot of a 2-way interaction from a structural equation model (SEM) that was estimated using the lavaan package.

Usage

semPlotInteraction(
  data,
  fit,
  predictor,
  centered_predictor,
  moderator,
  centered_moderator,
  interaction,
  outcome,
  covariates = NULL,
  predStr = NULL,
  modStr = NULL,
  outStr = NULL
)

Arguments

Details

Created by Johanna Caskey (johanna-caskey@uiowa.edu).

Value

Plot of two-way interaction from structural equation model.

See Also

Other plot: addText(), plot2WayInteraction(), ppPlot(), vwReg()

Other multipleRegression: lmCombine(), plot2WayInteraction(), ppPlot(), update_nested()

Other structural equation modeling: equiv_chi(), make_esem_model(), puc(), satorraBentlerScaledChiSquareDifferenceTestStatistic()

Examples

states <- as.data.frame(state.x77)
names(states)[which(names(states) == "HS Grad")] <- "HS.Grad"
states$Income_rescaled <- states$Income/100

# Mean Center Predictors
states$Illiteracy_centered <- scale(states$Illiteracy, scale = FALSE)
states$Murder_centered <- scale(states$Murder, scale = FALSE)

# Compute Interaction Term
states$interaction <- states$Illiteracy_centered * states$Murder_centered

# Specify model syntax
moderationModel <- '
  Income_rescaled ~ Illiteracy_centered + Murder_centered + interaction +
  HS.Grad
'

# Fit the model
moderationFit <- lavaan::sem(
  moderationModel,
  data = states,
  missing = "ML",
  estimator = "MLR",
  fixed.x = FALSE)

# Pass model to function (unlabeled plot)
semPlotInteraction(
  data = states,
  fit = moderationFit,
  predictor = "Illiteracy",
  centered_predictor = "Illiteracy_centered",
  moderator = "Murder",
  centered_moderator = "Murder_centered",
  interaction = "interaction",
  outcome = "Income_rescaled",
  covariates = "HS.Grad")

# Pass model to function (labeled plot)
semPlotInteraction(
  data = states,
  fit = moderationFit,
  predictor = "Illiteracy",
  centered_predictor = "Illiteracy_centered",
  moderator = "Murder",
  centered_moderator = "Murder_centered",
  interaction = "interaction",
  outcome = "Income_rescaled",
  covariates = "HS.Grad",
  predStr = "Illiteracy Level",
  modStr = "Murder Rate",
  outStr = "Income")

Set Lab Path.

Description

Sets the path directory to the lab drive.

Usage

setLabPath()

Details

Sets the path directory to the lab drive, and saves it in the object petersenLab.

Value

The object petersenLab with containing the path directory to the lab drive.

Examples

petersenLabPath <- setLabPath()

Simulate Area Under the ROC Curve (AUC).

Description

Simulate data with a specified area under the receiver operating characteristic curve—i.e., the AUC of an ROC curve.

Usage

simulateAUC(auc, n)

Arguments

Details

Simulates data with a specified area under the receiver operating characteristic curve—i.e., the AUC of an ROC curve.

Value

Dataframe with two columns:

• x is the predictor variable.

• y is the dichotomous criterion variable.

See Also

https://stats.stackexchange.com/questions/422926/generate-synthetic-data-given-auc/424213

Other simulation: complement(), simulateIndirectEffect()

Examples

simulateAUC(.60, 50000)
simulateAUC(.70, 50000)
simulateAUC(.80, 50000)
simulateAUC(.90, 50000)
simulateAUC(.95, 50000)
simulateAUC(.99, 50000)

Simulate Indirect Effect.

Description

Simulate indirect effect from mediation analyses.

Usage

simulateIndirectEffect(
  N = NA,
  x = NA,
  m = NA,
  XcorM = NA,
  McorY = NA,
  corTotal = NA,
  proportionMediated = NA,
  seed = NA
)

Arguments

Details

Co-created by Robert G. Moulder Jr. and Isaac T. Petersen

Value

• the correlation between the predictor variable (x) and the mediating variable (m).

• the correlation between the mediating variable (m) and the outcome variable (Y).

• the correlation between the predictor variable (x) and the outcome variable (Y).

• the direct correlation between the predictor variable (x) and the outcome variable (Y), while controlling for the mediating variable (m).

• the indirect correlation between the predictor variable (x) and the outcome variable (Y) through the mediating variable (m).

• the total correlation between the predictor variable (x) and the outcome variable (Y): i.e., the sum of the direct correlation and the indirect correlation.

• the proportion of the correlation between the predictor variable (x) and the outcome variable (Y) that is mediated through the mediating variable (m).

See Also

Other simulation: complement(), simulateAUC()

Examples

#INSERT

Specify Decimals.

Description

Specify the number of decimals to print.

Usage

specify_decimal(x, k)

Arguments

Details

[INSERT].

Value

Character vector of numbers with the specified number of decimal places.

See Also

Other formatting: apa(), pValue(), suppressLeadingZero()

Examples

# Prepare Data
v1 <- rnorm(1000)

# Specify Decimals
specify_decimal(v1, 2)

Standard Error of Measurement (IRT).

Description

Estimate the standard error of measurement in item response theory.

Usage

standardErrorIRT(information)

Arguments

Details

Estimate the standard error of measurement in item response theory using the test information (i.e., the sum of all items' information).

Value

Standard error of measurement for that amount of test information.

See Also

doi:10.1177/0146621613475471

Other IRT: deriv_d_negBinom(), discriminationToFactorLoading(), fourPL(), itemInformation(), reliabilityIRT(), test_info_4PL()

Examples

# Calculate information for 4 items
item1 <- itemInformation(b = -2, a = 0.6, theta = -4:4)
item2 <- itemInformation(b = -1, a = 1.2, theta = -4:4)
item3 <- itemInformation(b = 1, a = 1.5, theta = -4:4)
item4 <- itemInformation(b = 2, a = 2, theta = -4:4)

items <- data.frame(item1, item2, item3, item4)

# Calculate test information
items$testInformation <- rowSums(items)

# Calculate standard error of measurement
standardErrorIRT(items$testInformation)

Suppress Leading Zero.

Description

Suppress leading zero of numbers.

Usage

suppressLeadingZero(value)

Arguments

Details

[INSERT].

Value

Character vector of numbers without leading zeros.

See Also

Other formatting: apa(), pValue(), specify_decimal()

Examples

# Prepare Data
v1 <- rnorm(1000)

# Suppress Leading Zero
suppressLeadingZero(v1)

Test Information from Logistic IRT Model.

Description

Estimate test information from logistic item response theory model.

Usage

test_info_4PL(
  theta,
  alpha,
  beta,
  gamma = rep(0, length(alpha)),
  delta = rep(1, length(alpha))
)

error_variance_4PL(
  lower = -Inf,
  upper = Inf,
  alpha,
  beta,
  gamma = rep(0, length(alpha)),
  delta = rep(1, length(alpha)),
  mean = 0,
  sd = 1,
  density_cutoff = 1e-10
)

reliability_4PL(
  alpha,
  beta,
  gamma = rep(0, length(alpha)),
  delta = rep(1, length(alpha))
)

Arguments

Details

Created by Philipp Doebler (doebler@statistik.tu-dortmund.de) and Loreen Sabel (loreen.sabel@tu-dortmund.de).

Value

The amount of information for a given the test as a whole at each of the values of theta specified. Based on test information, one can estimate error variance and marginal reliability using error_variance_4PL() and reliability_4PL(), respectively.

See Also

Other IRT: deriv_d_negBinom(), discriminationToFactorLoading(), fourPL(), itemInformation(), reliabilityIRT(), standardErrorIRT()

Examples

test_info_4PL(0,1,0,0,1) # 0.25
test_info_4PL(-0.849, 1.1, -1, 0.2, 0.95) # Magis, 2013, Fig. 2
optimize(function(x)- test_info_4PL(x, 1.1, -1, 0.2, 0.95), c(-3, 3))

# test
set.seed(23)
# parameters (some are totally unrealistic)
alpha <- runif(20,0.5,2.5)
beta <- runif(20,-2,2)
gamma <- runif(20,0,0.3)
delta <- runif(20,0.8,1)
error_variance_4PL(
  lower = -Inf, upper = Inf,
  alpha, beta, gamma, delta)

error_variance_4PL(
  lower= -Inf, upper= Inf,
  alpha, beta, gamma, delta,
  density_cutoff = 1e-9)

error_variance_4PL(
  lower= -Inf, upper= Inf,
  alpha, beta, gamma, delta,
  density_cutoff = 1e-8)

error_variance_4PL(
  lower = -Inf, upper= Inf,
  alpha, beta, gamma, delta,
  density_cutoff = 1e-7)

reliability_4PL(alpha, beta, gamma, delta)

theta <- seq(-4, 4, length.out = 101)

plot(theta, test_info_4PL(theta, alpha, beta, gamma, delta))

Frequency Per Duration.

Description

Estimate frequency of a behavior for a particular duration.

Usage

timesPerInterval(
  num_occurrences = NULL,
  interval = NULL,
  duration = "month",
  not_occurred_past_year = NULL
)

timesPerLifetime(num_occurrences = NULL, never_occurred = NULL)

computeItemFrequencies(
  item_names,
  data,
  duration = "month",
  frequency_vars,
  interval_vars,
  not_in_past_year_vars
)

computeLifetimeFrequencies(
  item_names,
  data,
  frequency_vars,
  never_occurred_vars
)

Arguments

Details

Estimates the frequency of a given behavior for a particular duration, given a specified number of times it occurred during a specified interval.

Value

The frequency of the behavior for the specified duration.

Examples

timesPerInterval(
  num_occurrences = 2,
  interval = 3,
  duration = "month",
  not_occurred_past_year = 0
)

timesPerInterval(
  duration = "month",
  not_occurred_past_year = 1
)

timesPerLifetime(
  num_occurrences = 2,
  never_occurred = 0
)

timesPerLifetime(
  never_occurred = 1
)

Update Nested Models in Hierarchical Regression.

Description

Wrapper function to ensure the same observations are used for each updated model as were used in the first model.

Usage

update_nested(object, formula., ..., evaluate = TRUE)

Arguments

Details

Convenience wrapper function to ensure the same observations are used for each updated model as were used in the first model, to ensure comparability of models.

Value

lm model

See Also

https://stackoverflow.com/a/37341927

https://stackoverflow.com/a/37416336

https://stackoverflow.com/a/47195348

Other multipleRegression: lmCombine(), plot2WayInteraction(), ppPlot(), semPlotInteraction()

Examples

# Prepare Data
data("mtcars")

dat <- mtcars

# Create some missing values in mtcars
dat[1, "wt"] <- NA
dat[5, "cyl"] <- NA
dat[7, "hp"] <- NA

m1 <- lm(mpg ~ wt + cyl + hp, data = dat)
m2 <- update_nested(m1, . ~ . - wt)  # Remove wt
m3 <- update_nested(m1, . ~ . - cyl) # Remove cyl
m4 <- update_nested(m1, . ~ . - wt - cyl) # Remove wt and cyl
m5 <- update_nested(m1, . ~ . - wt - cyl - hp) # Remove all three variables
# (i.e., model with intercept only)

anova(m1, m2, m3, m4, m5)

Identify Variables of Different Types.

Description

Identifies the variables in common across two dataframes that have different types.

Usage

varsDifferentTypes(df1, df2)

Arguments

Details

Identifies the variables that have the same name across two dataframes that have different types, which can pose challenges for merging two dataframes.

Value

Dataframe with columns for the variable name, the variable type in df1 and the variable type in df2.

See Also

Other dataManipulation: columnBindFill(), convert.magic(), dropColsWithAllNA(), dropRowsWithAllNA()

Examples

# Prepare Data
df1 <- data.frame(
  A = 1:3,
  B = 2:4,
  C = 3:5
)

df2 <- data.frame(
  A = as.character(1:3),
  B = 2:4,
  C = as.factor(3:5)
)

# Check if any rows are not NA
varsDifferentTypes(df1, df2)

Visually Weighted Regression.

Description

Create watercolor plot to visualize weighted regression.

Usage

vwReg(
  formula,
  data,
  title = "",
  B = 1000,
  shade = TRUE,
  shade.alpha = 0.1,
  spag = FALSE,
  spag.color = "darkblue",
  mweight = TRUE,
  show.lm = FALSE,
  show.median = TRUE,
  median.col = "white",
  shape = 21,
  show.CI = FALSE,
  method = loess,
  bw = FALSE,
  slices = 200,
  palette = colorRampPalette(c("#FFEDA0", "#DD0000"), bias = 2)(20),
  ylim = NULL,
  quantize = "continuous",
  add = FALSE,
  ...
)

Arguments

Details

Creates a watercolor plot to visualize weighted regression.

Value

plot

See Also

https://www.nicebread.de/visually-weighted-regression-in-r-a-la-solomon-hsiang/

https://www.nicebread.de/visually-weighted-watercolor-plots-new-variants-please-vote/

http://www.fight-entropy.com/2012/07/visually-weighted-regression.html

http://www.fight-entropy.com/2012/08/visually-weighted-confidence-intervals.html

http://www.fight-entropy.com/2012/08/watercolor-regression.html

https://papers.ssrn.com/sol3/papers.cfm?abstract_id=2265501

Other plot: addText(), plot2WayInteraction(), ppPlot(), semPlotInteraction()

Other correlations: addText(), cor.table(), crossTimeCorrelation(), crossTimeCorrelationDF(), partialcor.table()

Examples

# Prepare Data
data("mtcars")
df <- data.frame(x = mtcars$hp, y = mtcars$mpg)

## Visually Weighted Regression

# Default
vwReg(y ~ x, df)

# Shade
vwReg(y ~ x, df, shade = TRUE, show.lm = TRUE, show.CI = TRUE,
quantize = "continuous")
vwReg(y ~ x, df, shade = TRUE, show.lm = TRUE, show.CI = TRUE,
quantize = "SD")

# Spaghetti
vwReg(y ~ x, df, shade = FALSE, spag = TRUE, show.lm = TRUE, show.CI = TRUE)
vwReg(y ~ x, df, shade = FALSE, spag = TRUE)

# Black/white
vwReg(y ~ x, df, shade = TRUE, spag = FALSE, show.lm = TRUE, show.CI = TRUE,
bw = TRUE, quantize = "continuous")
vwReg(y ~ x, df, shade = TRUE, spag = FALSE, show.lm = TRUE, show.CI = TRUE,
bw = TRUE, quantize = "SD")
vwReg(y ~ x, df, shade = FALSE, spag = TRUE, show.lm = TRUE, show.CI = TRUE,
bw = TRUE, quantize = "SD")

# Change the bootstrap smoothing
vwReg(y ~ x, df, family = "symmetric") # use an M-estimator for
# bootstrap smoothers. Usually yields wider confidence intervals
vwReg(y ~ x, df, span = 1.7) # increase the span of the smoothers
vwReg(y ~ x, df, span = 0.5) # decrease the span of the smoothers

# Change the color scheme
vwReg(y ~ x, df, palette = viridisLite::viridis(4)) # viridis
vwReg(y ~ x, df, palette = viridisLite::magma(4)) # magma
vwReg(y ~ x, df, palette = RColorBrewer::brewer.pal(9, "YlGnBu")) # change the
# color scheme, using a predefined ColorBrewer palette. You can see all
# available palettes by using this command:
# `library(RColorBrewer); display.brewer.all()`
vwReg(y ~ x, df, palette = grDevices::colorRampPalette(c("white","yellow",
"green","red"))(20)) # use a custom-made palette
vwReg(y ~ x, df, palette = grDevices::colorRampPalette(c("white","yellow",
"green","red"), bias = 3)(20)) # use a custom-made palette, with the
# parameter bias you can shift the color ramp to the “higher” colors
vwReg(y ~ x, df, bw = TRUE) # black and white version
vwReg(y ~ x, df, shade.alpha = 0, palette = grDevices::colorRampPalette(
c("black","grey30","white"), bias = 4)(20)) # Milky-Way Plot
vwReg(y ~ x, df, shade.alpha = 0, slices = 400, palette =
grDevices::colorRampPalette(c("black","green","yellow","red"),
bias = 5)(20), family = "symmetric") # Northern Light Plot/ fMRI plot
vwReg(y ~ x, df, quantize = "SD") # 1-2-3-SD plot

Write Encrypted Data.

Description

Write data to encrypted file.

Usage

write.aes(df, filename, key)

Arguments

Details

Writes data to an encrypted file. To read data from an encrypted file, see read.aes.

Value

A file with encrypted data.

See Also

https://stackoverflow.com/questions/25318800/how-do-i-read-an-encrypted-file-from-disk-with-r/25321586#25321586

Other encrypted: read.aes()

Examples

# Location Where to Save Encryption Key on Local Computer
  #(where only you should have access to it)
#encryptionKeyLocation <- file.path(getwd(), "/encryptionKey.RData",
#  fsep = "") #Can change to a different path, e.g.: "C:/Users/[USERNAME]/"

# Generate a Temporary File Path for Encryption Key
encryptionKeyLocation <- tempfile(fileext = ".RData")

# Generate Encryption Key
key <- as.raw(sample(1:16, 16))

# Save Encryption Key
save(key, file = encryptionKeyLocation)

# Specify Credentials
credentials <- "Insert My Credentials Here"

# Generate a Temporary File Path for Encrypted Credentials
encryptedCredentialsLocation <- tempfile(fileext = ".txt")

# Save Encrypted Credentials
#write.aes(
#  df = credentials,
#  filename = file.path(getwd(), "/encrypytedCredentials.txt", fsep = ""),
#  key = key) #Change the file location to save this on the lab drive

write.aes(
  df = credentials,
  filename = encryptedCredentialsLocation,
  key = key)

rm(credentials)
rm(key)