Package 'riskutility'

Description

Provides a comprehensive toolkit for measuring disclosure risk and data utility in anonymized and synthetic datasets. The package supports multivariate Risk-Utility (R-U) evaluation following the framework described in "Beyond the Trade-off Curve" (Thees, Mueller, Templ 2026).

Details

The package covers the following metric families:

Attribution-based risk: CAP-family measures (dcap, tcap, weap, disco) that assess whether adversaries can infer sensitive attributes from quasi-identifiers.

ML-based risk: The rapid metric uses machine learning models to predict sensitive attributes, capturing complex non-linear relationships that CAP methods may miss.

Distance-based risk: Holdout-based metrics (dcr, nndr, ims) detect memorization by comparing distances to training vs. holdout data.

Utility measures: Functions for assessing how well synthetic data preserves the statistical properties of the original, including pMSE, specks, hellinger, energy_distance, ci_proximity, mmd, tstr, copula_fidelity, tail_fidelity, subgroup_utility, regression_fidelity, contingency_fidelity, propscore, gower, and mqs.

Privacy models: Classical privacy model assessments including kanonymity, ldiversity, individual_risk, and suda.

Extended privacy risk: Additional risk models including tcloseness, nnaa, singling_out, linkability, population_uniqueness, epsilon_identifiability, delta_presence, hitting_rate, attacker_risk, drisk, domias, and mia_classifier.

Information-theoretic measures: Entropy and divergence measures (KLDiv, JSDiv, RenyiEntropy, mutualInformation, information_surprisal, max_info_leakage) for quantifying information leakage.

Comparison functions: Direct comparison tools for original and synthetic data including visualization, statistical tests, dimensionality reduction, and predictive model comparisons.

R-U mapping: The rumap function provides comprehensive multivariate Risk-Utility evaluation with Pareto frontier identification.

Visualization Functions

compare_histograms()

Compare distributions using histograms (weighted and unweighted, faceted).

compare_boxplots()

Compare numeric variable distributions using weighted/unweighted boxplots.

Statistical Tests

compare_ks_test()

Kolmogorov-Smirnov test for comparing distributions.

compare_chisq_gof()

Chi-square Goodness-of-Fit test with multi-dimensional frequency tables (weighted/unweighted).

Distributional Comparisons

compare_wasserstein()

Calculate Wasserstein distances between numeric or nominal variables.

compare_multivariate_distribution()

Compare joint distributions using Mahalanobis distance or mutual information.

Summary Statistics

compare_means_frequencies()

Compare means, medians, robust statistics, skewness, kurtosis for numeric data, and frequencies for categorical data.

compare_correlation_matrices()

Compare correlations: Pearson, Spearman, robust correlations for numeric, categorical, and mixed data types.

compare_multivariate_summary_statistics()

Calculate joint summary statistics for numeric and categorical data.

Dimensionality Reduction

compare_pca()

Principal Component Analysis comparison with biplot options.

compare_embedding()

t-SNE, UMAP, or Sammon's Mapping (MDS) visualizations for embedding comparisons.

Machine Learning and Predictive Modeling

compare_model_performance()

Evaluate predictive model performance (accuracy, precision, recall, F1-score, AUC, R-squared, RMSE) using cross-validation.

compare_feature_importance()

Evaluate stability of feature importance measures from Random Forests, Decision Trees, Gradient Boosting, Permutation Importance, and SHAP values.

Dependencies

The package builds on ggplot2 and data.table; uses VIM, MASS, randomForest, and ranger for core computations; and optionally integrates with simPop, synthpop, Rtsne, uwot, caret, and robustbase when installed.

Missing-value handling

Two parameter conventions are used deliberately, according to how much choice a function offers. A logical na.rm (the standard R idiom) appears where the only decision is whether to drop incomplete cases before computing a statistic. A character selector appears where several treatments are supported: na ("remove", "impute", "stop", and in pMSE additionally "indicator") in propscore, pMSE, and mqs; and na_strategy ("constant", "drop", "median") in rapid, which acts specifically on the sensitive target attribute.

Author(s)

Matthias Templ

See Also

ggplot, data.table

Examples

# Create example datasets
set.seed(123)
X <- data.frame(
  income = rnorm(100, mean = 50000, sd = 10000),
  age = rnorm(100, mean = 40, sd = 10),
  gender = sample(c("Male", "Female"), 100, replace = TRUE),
  weight = runif(100, 0.5, 1.5)
)
Y <- data.frame(
  income = rnorm(100, mean = 48000, sd = 12000),
  age = rnorm(100, mean = 42, sd = 11),
  gender = sample(c("Male", "Female"), 100, replace = TRUE),
  weight = runif(100, 0.5, 1.5)
)

# Example for histogram comparison
compare_histograms(X, Y, num_var = "income", cat_vars = c("gender"),
                   weight_X = "weight", weight_Y = "weight")

# Example for PCA comparison
X_pca <- data.frame(
  income = rnorm(100, mean = 50000, sd = 10000),
  age = rnorm(100, mean = 40, sd = 10),
  expenditure = rnorm(100, mean = 2000, sd = 500)
)
Y_pca <- data.frame(
  income = rnorm(100, mean = 48000, sd = 12000),
  age = rnorm(100, mean = 42, sd = 11),
  expenditure = rnorm(100, mean = 2200, sd = 600)
)
res_pca <- compare_pca(X_pca, Y_pca,
                       vars = c("income", "age", "expenditure"),
                       biplot = TRUE)

Description

Computes re-identification risk under three canonical attacker models from Statistical Disclosure Control (SDC). Each model makes different assumptions about the attacker's background knowledge, leading to different risk quantifications based on quasi-identifier (QI) frequencies.

Usage

attacker_risk(X, ...)

## S3 method for class 'synth_pair'
attacker_risk(X, data = c("synthetic", "original"), ...)

## Default S3 method:
attacker_risk(
  X,
  key_vars,
  sampling_fraction = 0.01,
  model = c("all", "prosecutor", "journalist", "marketer"),
  na.rm = TRUE,
  ...
)

Arguments

Details

The three attacker models differ in their assumptions:

Prosecutor Model: The attacker knows a specific individual is in the dataset and attempts to re-identify them. For each record in equivalence class k with frequency f_k, the individual risk is:

$risk_i = 1 / f_k$

The global prosecutor risk is the mean over all records.

Journalist Model: The attacker does not know whether the target is in the dataset. This requires estimating the population frequency F_k from the sample frequency f_k using the sampling fraction:

$F_k \approx f_k / sampling\_fraction$

$risk_i = 1 / F_k = sampling\_fraction / f_k$

This typically yields lower risks than the prosecutor model because population frequencies are larger than sample frequencies.

Marketer Model: The attacker does not target any specific individual but wants to re-identify as many records as possible. The global risk is the expected number of successful re-identifications divided by the total number of records:

$risk = (1/n) \sum_{i=1}^{n} 1/f_k(i)$

Note that for the marketer model, the global risk equals the prosecutor global risk (both are mean 1/f_k).

Interpretation:

• Risk close to 1: High re-identification risk (many unique records)

• Risk > 0.1: Elevated risk, privacy protection may be insufficient

• Risk < 0.05: Generally acceptable

• Risk close to 0: Low re-identification risk (large equivalence classes)

Value

An object of class "attacker_risk" containing:

• risk_per_record: data.frame with columns: prosecutor, journalist (if applicable), freq (f_k)

• global_risk: named list with prosecutor (mean 1/f_k), journalist (mean 1/F_k), marketer (1/n * sum(1/f_k))

• n_uniques: count of records with f_k = 1

• pct_uniques: fraction of unique records

• freq_table: table of QI combination frequencies

• model: model(s) used

• key_vars: quasi-identifier variables used

• sampling_fraction: sampling fraction used

• n_records: number of records assessed

• privacy_pass: logical, prosecutor global risk <= 0.1

Author(s)

Matthias Templ

References

Hundepool, A., Domingo-Ferrer, J., Franconi, L., Giessing, S., Schulte Nordholt, E., Spicer, K. & de Wolf, P.-P. (2012). Statistical Disclosure Control. John Wiley & Sons. doi:10.1002/9781118348239

Elliot, M. (2006). A Computational Algorithm for Handling the Special Uniques Problem. International Journal of Uncertainty, Fuzziness and Knowledge-Based Systems, 14(Supp. 1), 45-59.

Templ, M. (2017). Statistical Disclosure Control for Microdata: Methods and Applications in R. Springer. doi:10.1007/978-3-319-50272-4

See Also

kanonymity for k-anonymity assessment, individual_risk for model-based individual risk, suda for special uniques detection

Other privacy-models: delta_presence(), disclosure_report(), domias(), drisk(), epsilon_identifiability(), individual_risk(), inspect_record(), kanonymity(), ldiversity(), linkability(), merge_per_record(), mia_classifier(), population_uniqueness(), recordLinkage(), risk_by_group(), singling_out(), suda(), tcloseness(), top_at_risk()

Examples

# Create example data
set.seed(123)
data <- data.frame(
  age = sample(c("young", "middle", "old"), 200, replace = TRUE),
  gender = sample(c("M", "F"), 200, replace = TRUE),
  region = sample(c("N", "S", "E", "W"), 200, replace = TRUE),
  income = sample(c("low", "medium", "high"), 200, replace = TRUE)
)

# Compute all three attacker models
result <- attacker_risk(data, key_vars = c("age", "gender", "region"))
print(result)
summary(result)


plot(result)

# Prosecutor model only
result_p <- attacker_risk(data, key_vars = c("age", "gender", "region"),
                          model = "prosecutor")
print(result_p)

# Journalist model with different sampling fraction
result_j <- attacker_risk(data, key_vars = c("age", "gender"),
                          model = "journalist", sampling_fraction = 0.05)
print(result_j)

Description

Computes chi-square based utility measures comparing original and synthetic data frequency tables. Implements VW (Voas-Williamson), FT (Freeman-Tukey), G (log-likelihood), and JSD (Jensen-Shannon Divergence) statistics.

Usage

chisq_utility(X, ...)

## S3 method for class 'synth_pair'
chisq_utility(X, vars = NULL, max_cells = 10000, ...)

## Default S3 method:
chisq_utility(
  X,
  Y,
  vars,
  max_cells = 10000,
  weight_X = NULL,
  weight_Y = NULL,
  ...
)

Arguments

Details

These measures compare frequency tables from original and synthetic data.

Chi-Square ($\chi^2$): The standard Pearson chi-square statistic:

$\chi^2 = \sum \frac{(O - E)^2}{E}$

where O is observed (synthetic) and E is expected (original proportions scaled).

Voas-Williamson (VW): A normalized chi-square measure (Voas & Williamson, 2001):

$VW = \frac{\chi^2 - df}{N}$

where df is degrees of freedom and N is sample size. Lower is better; 0 indicates perfect replication.

Freeman-Tukey (FT): Uses a variance-stabilizing transformation:

$FT = 4 \sum (\sqrt{O} - \sqrt{E})^2$

G-Test (Log-Likelihood Ratio):

$G = 2 \sum O \log(O/E)$

Jensen-Shannon Divergence (JSD): A symmetric, bounded divergence measure:

$JSD = \frac{1}{2}[KL(P||M) + KL(Q||M)]$

where M = (P + Q)/2, and KL is Kullback-Leibler divergence.

Interpretation:

• VW close to 0: excellent utility

• VW < 0.01: good utility

• VW > 0.05: potential utility concerns

• JSD ranges 0-1: 0 = identical, 1 = completely different

Value

An object of class "chisq_utility" containing:

• chi2: standard chi-square statistic

• df: degrees of freedom

• p_value: p-value for chi-square test

• VW: Voas-Williamson statistic (normalized chi-square)

• FT: Freeman-Tukey statistic

• G: G-test (log-likelihood ratio) statistic

• JSD: Jensen-Shannon Divergence

• n_cells: number of cells in the table

• n_empty_orig: empty cells in original

• n_empty_synth: empty cells in synthetic

• pct_utility: overall utility percentage (based on VW)

• table_orig: frequency table for original data

• table_synth: frequency table for synthetic data

Author(s)

Matthias Templ

References

Voas, D., & Williamson, P. (2001). Evaluating goodness-of-fit measures for synthetic microdata. Geographical and Environmental Modelling, 5(2), 177-200.

Freeman, M. F., & Tukey, J. W. (1950). Transformations related to the angular and the square root. Annals of Mathematical Statistics, 21, 607-611.

See Also

propscore for propensity score utility, specks for SPECKS utility measure

Other utility: ci_proximity(), contingency_fidelity(), copula_fidelity(), energy_distance(), gower(), hellinger(), mmd(), mqs(), pMSE(), plot.rumap(), propscore(), regression_fidelity(), rumap(), specks(), subgroup_utility(), tail_fidelity(), tstr()

Examples

# Create example data
set.seed(123)
orig <- data.frame(
  age = sample(c("young", "middle", "old"), 500, replace = TRUE),
  gender = sample(c("M", "F"), 500, replace = TRUE),
  region = sample(c("N", "S", "E", "W"), 500, replace = TRUE)
)

# Synthetic data with some differences
synth <- data.frame(
  age = sample(c("young", "middle", "old"), 500, replace = TRUE,
               prob = c(0.35, 0.4, 0.25)),
  gender = sample(c("M", "F"), 500, replace = TRUE),
  region = sample(c("N", "S", "E", "W"), 500, replace = TRUE)
)

# Compute chi-square utility
result <- chisq_utility(orig, synth, vars = c("age", "gender"))
print(result)
summary(result)
plot(result)

Description

Given a sample, how much higher is the change compared to the sampling error.

Usage

ci_overlap(x, y, estimator = median, R = 1000)

Arguments

Value

A list with the following elements:

m_x

The point estimate of the function provided in function argument estimator for the vector x.

m_y

The point estimate of the function provided in function argument estimator for the vector x.

ci_x

The confidence interval of x.

ci_y

The confidence interval of y.

overlap

Overlap in confidence intervals.

overlap_perc

Overlap in confidence intervals in percentage

ratioSE

How much differ is the point estimate of y compared to the half of the confidence interval.

Author(s)

Matthias Templ

See Also

Other comparison: compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

x <- rnorm(100, 10)
y <- x + runif(100, 0, 1)

ci_overlap(x, y, R = 1000)

x <- factor(sample(1:2, size = 100, replace = TRUE))
y <- factor(sample(1:2, size = 100, replace = TRUE))

my_func <- function(x){
  table(x)[1]
}

ci_overlap(x, y, estimator = my_func)

Description

Computes the proximity of confidence intervals for numeric variables between original and synthetic datasets. This utility measure assesses how well the synthetic data preserves point estimates and their uncertainty.

Usage

ci_proximity(X, ...)

## S3 method for class 'synth_pair'
ci_proximity(X, ...)

## Default S3 method:
ci_proximity(
  X,
  Y,
  vars = NULL,
  conf.level = 0.95,
  weight_X = NULL,
  weight_Y = NULL,
  na.rm = TRUE,
  ...
)

Arguments

Details

For each numeric variable, this function computes:

• CI Overlap: Proportion of the CIs that overlap, ranging from 0 (no overlap) to 1 (complete overlap or one contains the other)

• Relative Error: |mean_Y - mean_X| / |mean_X|, measuring how close the synthetic mean is to the original mean

• Proximity Score: Combined measure accounting for both mean accuracy and CI overlap

The CI overlap coefficient is computed as:

$Overlap = \frac{\max(0, \min(U_X, U_Y) - \max(L_X, L_Y))}{\max(U_X, U_Y) - \min(L_X, L_Y)}$

where L and U are lower and upper CI bounds.

For utility assessment, higher proximity scores indicate better preservation of statistical properties in synthetic data.

Value

An object of class "ci_proximity" containing:

• per_variable: data.frame with CI comparison for each variable

• proximity_mean: mean proximity score across all variables

• overlap_mean: mean overlap coefficient across all variables

• relative_error_mean: mean relative error of synthetic means

• n_vars: number of variables compared

• vars: variable names used

• conf.level: confidence level used

Author(s)

Matthias Templ

See Also

ci_overlap for basic CI overlap calculation, compare_ks_test for distribution comparison

Other utility: chisq_utility(), contingency_fidelity(), copula_fidelity(), energy_distance(), gower(), hellinger(), mmd(), mqs(), pMSE(), plot.rumap(), propscore(), regression_fidelity(), rumap(), specks(), subgroup_utility(), tail_fidelity(), tstr()

Examples

set.seed(123)
# Original data
X <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  age = rnorm(500, mean = 40, sd = 10),
  score = rnorm(500, mean = 100, sd = 15)
)

# Good synthetic data (similar means)
Y_good <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  age = rnorm(500, mean = 40, sd = 10),
  score = rnorm(500, mean = 100, sd = 15)
)

# Poor synthetic data (shifted means)
Y_poor <- data.frame(
  income = rnorm(500, mean = 55000, sd = 10000),
  age = rnorm(500, mean = 45, sd = 10),
  score = rnorm(500, mean = 90, sd = 15)
)

result_good <- ci_proximity(X, Y_good)
print(result_good)
summary(result_good)

result_poor <- ci_proximity(X, Y_poor)
print(result_poor)

Description

This function compares the distribution of a numeric variable between two datasets (e.g., an original dataset and an anonymized/synthetic version) using boxplots. It allows conditional faceting on categorical variables and supports different visualization styles, including side-by-side, overlapping, and separate boxplots. Sampling weights can be considered to compute weighted medians and quartiles.

Usage

compare_boxplots(X, ...)

## S3 method for class 'synth_pair'
compare_boxplots(X, ...)

## Default S3 method:
compare_boxplots(
  X,
  Y,
  num_var,
  cat_vars,
  weight_X = NULL,
  weight_Y = NULL,
  facet_type = "wrap",
  plot_type = "side",
  ...
)

Arguments

Value

A ggplot2 boxplot visualization comparing the distributions.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  age = sample(20:80, 500, replace = TRUE),
  gender = sample(c("Male", "Female"), 500, replace = TRUE),
  income = rnorm(500, mean = 50000, sd = 10000),
  weight = runif(500, 0.5, 1.5)
)

Y <- data.frame(
  age = sample(20:80, 1000, replace = TRUE),
  gender = sample(c("Male", "Female"), 1000, replace = TRUE),
  income = rnorm(1000, mean = 48000, sd = 12000),
  weight = runif(1000, 0.5, 1.5)
)

compare_boxplots(X, Y, num_var = "income", cat_vars = c("gender"),
                 weight_X = "weight", weight_Y = "weight",
                 facet_type = "wrap", plot_type = "side")

compare_boxplots(X, Y, num_var = "income", cat_vars = c("gender"),
                 weight_X = "weight", weight_Y = "weight",
                 facet_type = "wrap", plot_type = "separate")

Description

This function performs a Chi-square goodness-of-fit test to compare the observed frequencies in an anonymized/synthetic dataset (Y) with the expected frequencies from an original dataset (X), based on one or more categorical variables (cat_vars). If sampling weights are provided via weight_X and weight_Y, weighted frequency tables are computed using xtabs with a weight formula; otherwise, unweighted counts are computed. Additionally, if group_vars is provided, the Chi-square test is computed separately for each unique group defined by the grouping variables. For each group (or overall if no grouping is provided), the expected frequencies from X are scaled to match the total count in Y before performing the test.

Usage

compare_chisq_gof(X, ...)

## S3 method for class 'synth_pair'
compare_chisq_gof(X, ...)

## Default S3 method:
compare_chisq_gof(
  X,
  Y,
  cat_vars,
  group_vars = NULL,
  weight_X = NULL,
  weight_Y = NULL,
  simulate_p = TRUE,
  B = 2000,
  ...
)

Arguments

Details

The function computes frequency tables for the variables in cat_vars for both datasets X and Y. If sampling weights are provided, xtabs is used with a formula of the form "weight ~ var1 + var2 + ..."; otherwise, unweighted counts are computed using xtabs. A complete grid of all possible combinations of factor levels is created to ensure that missing cells are filled with zeros. The expected frequencies from X are scaled to match the total count in Y before performing the Chi-square goodness-of-fit test. When group_vars is provided, the function iterates over each unique combination of grouping variables and performs the test on each subset. The simulation of p-values (when simulate_p is TRUE) helps to handle cases with small expected frequencies.

Value

A data.table containing the grouping variables (if provided), the Chi-square statistic, degrees of freedom, and p-value for each group (or one overall row if group_vars is NULL).

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  gender = factor(sample(c("Male", "Female"), 500, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 500, replace = TRUE)),
  occupation = factor(sample(c("Engineer", "Doctor", "Artist", "Teacher"), 500, replace = TRUE)),
  weight = runif(500, 0.5, 1.5)
)

Y <- data.frame(
  gender = factor(sample(c("Male", "Female"), 1000, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 1000, replace = TRUE)),
  occupation = factor(sample(c("Engineer", "Doctor", "Artist", "Teacher"), 500, replace = TRUE)),
  weight = runif(1000, 0.5, 1.5)
)

# Overall test (weighted) with simulated p-values
result_overall <- compare_chisq_gof(X, Y, cat_vars = c("gender", "region"),
                                    weight_X = "weight", weight_Y = "weight",
                                    simulate_p = FALSE, B = 2000)
print(result_overall)

# Test conditionally by grouping variables (e.g., "region" and "gender")
result_by_group <- compare_chisq_gof(X, Y,
  cat_vars = c("gender"),
  group_vars = c("region", "occupation"),
  weight_X = "weight", weight_Y = "weight",
  simulate_p = TRUE, B = 2000)
print(result_by_group)

Description

This function computes and compares the correlation matrices for a specified set of variables between two datasets (X and Y), such as an original dataset and an anonymized/synthetic version. For continuous variables, the user can choose Pearson, Spearman, or robust correlation methods (using either MM‐estimation or the MCD estimator). For categorical variables, the function computes Cramér’s V. For variable pairs with one continuous and one nominal variable, if the nominal variable is binary, the point‐biserial correlation is computed; otherwise, the correlation ratio (eta squared) is used and transformed into a correlation‐like measure.

Usage

compare_correlation_matrices(X, ...)

## S3 method for class 'synth_pair'
compare_correlation_matrices(X, ...)

## Default S3 method:
compare_correlation_matrices(
  X,
  Y,
  vars,
  method = "pearson",
  mixed_method = "point_biserial",
  ...
)

Arguments

Details

The function determines the type of each variable (continuous or categorical) based on its class. Continuous variables are correlated using the specified method. For robust methods, the MCD estimator from the robustbase package or MM-estimation via MASS::cov.rob is used. For categorical variables, Cramér’s V is computed from the chi-square statistic. For mixed pairs, if one variable is continuous and the other nominal, point-biserial correlation is used when the nominal variable has two levels, and the correlation ratio (eta squared) is computed otherwise.

Value

A list with three elements:

corr_X

The correlation matrix computed from dataset X.

corr_Y

The correlation matrix computed from dataset Y.

diff

The absolute difference matrix between the two correlation matrices.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  age = rnorm(500, mean = 40, sd = 10),
  gender = factor(sample(c("Male", "Female"), 500, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 500, replace = TRUE))
)

Y <- data.frame(
  income = rnorm(1000, mean = 48000, sd = 12000),
  age = rnorm(1000, mean = 42, sd = 11),
  gender = factor(sample(c("Male", "Female"), 1000, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 1000, replace = TRUE))
)

# Continuous correlation using Pearson
res1 <- compare_correlation_matrices(X, Y, vars = c("income", "age"), method = "pearson")
print(res1$corr_X)
print(res1$corr_Y)
print(res1$diff)

# Categorical correlation using Cramér's V (computed automatically)
res2 <- compare_correlation_matrices(X, Y, vars = c("gender", "region"), method = "pearson")
print(res2$corr_X)
print(res2$corr_Y)

# Mixed example: income (continuous) vs gender (nominal); point-biserial or eta-squared is computed.
res3 <- compare_correlation_matrices(X, Y,
  vars = c("income", "gender"),
  method = "pearson",
  mixed_method = "point_biserial")
print(res3$corr_X)
print(res3$corr_Y)

Description

This function calculates and compares the empirical cumulative distribution functions (ECDF) of a sample (X) and a population (Y) for specified variables. It can handle weighted samples and provides options for conditional CDFs, and approximation.

Usage

compare_distributions_cont(X, ...)

## S3 method for class 'synth_pair'
compare_distributions_cont(X, ...)

## Default S3 method:
compare_distributions_cont(
  X,
  Y,
  vars = NULL,
  kind = NULL,
  conditional = NULL,
  weights = NULL,
  approx = c(FALSE, TRUE),
  n_approx = 10000,
  bounds = TRUE,
  ...
)

Arguments

Details

The following default calculations are taken:

If vars links to one or two numeric variables: the (weighted) ECDF is calculated for both data sets. With kind one can change to "density", "density_bayes", "density_ratio" and "density_ratio_bayes".

If vars links to one, two or three categorical variables: barcharts and mosaic plots are made.

If vars links to one numeric and one categorical variable: boxplots statistics are calculated.

If weights are provided, the weighted versions are computed. The function supports conditional estimates, which are computed by conditioning on the specified variable.

If the approx argument is set to TRUE, the ECDFs are approximated using the specified number of points (n_approx).

Value

A list with class "compare" containing the following components:

formula

A formula representing the CDF comparison.

ecdf

A data frame with the ECDF values.

kind

A character string indicating the type of comparison ("numeric").

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

# Example data 1: Complex sample and population data
S <- data.frame(
  age = sample(20:80, 100, replace = TRUE),
  income = rnorm(100, mean = 50000, sd = 10000),
  weights = runif(100, 0.5, 2)
)
U <- data.frame(
  age = sample(20:80, 10000, replace = TRUE),
  income = rnorm(10000, mean = 50000, sd = 10000)
)

# Compare ECDFs for age and income
result <- compare_distributions_cont(S, U, vars = c("age", "income"), weights = "weights")

# View the result
print(result$formula)
head(result$ecdf)
plot(result)
plot(result, which = "density")
plot(result, which = "density_bayes")

# Example data 2: Complex sample and complex sample data
X <- data.frame(
  age = sample(20:80, 100, replace = TRUE),
  income = rnorm(100, mean = 50000, sd = 10000),
  weights = runif(100, 0.5, 2)
)
Y <- data.frame(
  age = sample(20:80, 100, replace = TRUE),
  income = rnorm(10000, mean = 50000, sd = 10000),
  weights = runif(100, 0.5, 2)
)

# Compare ECDFs for age and income
result <- compare_distributions_cont(X, Y, vars = c("age", "income"), weights = "weights")

# View the result
print(result$formula)
head(result$ecdf)
plot(result)
plot(result, which = "density")
plot(result, which = "density_bayes")
plot(result, which = "density_ratio")

# Example data 3: Sample and sample data
X <- data.frame(
  age = sample(20:80, 100, replace = TRUE),
  income = rnorm(100, mean = 50000, sd = 10000)
)
Y <- data.frame(
  age = sample(20:80, 100, replace = TRUE),
  income = rnorm(100, mean = 50000, sd = 10000)
)

# Compare ECDFs for age and income
result <- compare_distributions_cont(X, Y, vars = c("age", "income"), n_approx = 100)

# View the result
print(result$formula)
head(result$ecdf)
plot(result)
plot(result, which = "density")
plot(result, which = "density_bayes")

Description

This function compares two datasets using dimensionality reduction methods: t-distributed Stochastic Neighbor Embedding (t-SNE), Uniform Manifold Approximation and Projection (UMAP), and Multidimensional Scaling (MDS) with Sammon's mapping. It allows visualizing complex data structures, such as clusters or outliers, in two-dimensional space.

Usage

compare_embedding(X, ...)

## S3 method for class 'synth_pair'
compare_embedding(X, ...)

## Default S3 method:
compare_embedding(
  X,
  Y,
  vars,
  method = "tsne",
  perplexity = 30,
  n_neighbors = 15,
  side_by_side = FALSE,
  seed = NULL,
  ...
)

Arguments

Value

A list containing the embeddings and a ggplot2 visualization.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

# MDS (Sammon) embedding - doesn't require optional packages
set.seed(123)
X <- data.frame(
  var1 = rnorm(100, 50, 10),
  var2 = rnorm(100, 40, 8),
  var3 = rnorm(100, 30, 5)
)
Y <- data.frame(
  var1 = rnorm(100, 48, 12),
  var2 = rnorm(100, 42, 9),
  var3 = rnorm(100, 28, 6)
)
res_mds <- compare_embedding(X, Y, vars = c("var1", "var2", "var3"),
                             method = 'mds')


# t-SNE (requires Rtsne package)
if (requireNamespace("Rtsne", quietly = TRUE)) {
  res_tsne <- compare_embedding(X, Y, vars = c("var1", "var2", "var3"),
                                method = 'tsne')
}

# UMAP (requires uwot package)
if (requireNamespace("uwot", quietly = TRUE)) {
  res_umap <- compare_embedding(X, Y, vars = c("var1", "var2", "var3"),
                                method = 'umap')
}

Description

Compares feature importance rankings from models trained on real (X) and synthetic/anonymized (Y) datasets.

Usage

compare_feature_importance(X, ...)

## S3 method for class 'synth_pair'
compare_feature_importance(X, ...)

## Default S3 method:
compare_feature_importance(
  X,
  Y,
  formula,
  method,
  importance_type = c("model", "permutation", "shap"),
  pred_wrapper = NULL,
  metric = "Accuracy",
  nsim = 5,
  ...
)

Arguments

Details

This function supports three types of feature importance measures:

• Model-based: Uses caret::varImp.

• Permutation-based: Uses vip::vi_permute. A prediction wrapper is required.

• SHAP-based: Computes mean absolute Shapley values using the kernelshap package.

• If the target variable in X is a factor, the target variable in Y is coerced to a factor with the same levels.

• For SHAP importance, the function uses the kernelshap package. The returned importance is the mean absolute SHAP value per feature.

• For classification tasks in the SHAP branch, if no pred_wrapper is provided, predictions are taken as the probability of the first class.

Value

A list containing:

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

# Model-based importance (Decision Tree) - basic example
set.seed(123)
X <- data.frame(
  age = rnorm(100, 40, 10),
  income = 50000 + rnorm(100, 0, 10000),
  outcome = factor(sample(c("Yes", "No"), 100, replace = TRUE))
)
Y <- data.frame(
  age = rnorm(100, 42, 12),
  income = 48000 + rnorm(100, 0, 12000),
  outcome = factor(sample(c("Yes", "No"), 100, replace = TRUE))
)
res_tree <- compare_feature_importance(X, Y, outcome ~ .,
                                       method = "rpart",
                                       importance_type = "model")


# Permutation importance (requires vip package)
if (requireNamespace("vip", quietly = TRUE)) {
  pred_fn <- function(object, newdata) predict(object, newdata, type = "raw")
  res_perm <- compare_feature_importance(
    X, Y, outcome ~ ., method = "rf",
    importance_type = "permutation",
    pred_wrapper = pred_fn, metric = "Accuracy", nsim = 5
  )
}

Description

This function compares the distribution of a numeric variable between two datasets (e.g., an original dataset and an anonymized/synthetic version) using histograms. It allows conditional faceting on categorical variables and supports different visualization styles, including overlapping, side-by-side, stacked, and separate histograms.

Usage

compare_histograms(X, ...)

## S3 method for class 'synth_pair'
compare_histograms(X, ...)

## Default S3 method:
compare_histograms(
  X,
  Y,
  num_var,
  cat_vars,
  weight_X = NULL,
  weight_Y = NULL,
  facet_type = "wrap",
  bins = 30,
  transparency = 0.4,
  plot_type = "overlap",
  ...
)

Arguments

Value

A ggplot2 histogram visualization comparing the distributions.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  age = sample(20:80, 500, replace = TRUE),
  gender = sample(c("Male", "Female"), 500, replace = TRUE),
  income = rnorm(500, mean = 50000, sd = 10000),
  weight = runif(500, 0.5, 1.5)
)

Y <- data.frame(
  age = sample(20:80, 1000, replace = TRUE),
  gender = sample(c("Male", "Female"), 1000, replace = TRUE),
  income = rnorm(1000, mean = 48000, sd = 12000),
  weight = runif(1000, 0.5, 1.5)
)

# Overlapping histograms with transparency
compare_histograms(X, Y, num_var = "income", cat_vars = c("gender"),
                   weight_X = "weight", weight_Y = "weight",
                   facet_type = "wrap", plot_type = "overlap")
# Side-by-side histograms
compare_histograms(X, Y, num_var = "income", cat_vars = c("gender"),
                   weight_X = "weight", weight_Y = "weight",
                   facet_type = "wrap", plot_type = "side")

# Stacked histograms
compare_histograms(X, Y, num_var = "income", cat_vars = c("gender"),
                   weight_X = "weight", weight_Y = "weight",
                   facet_type = "wrap", plot_type = "stack")

# Separate histograms
compare_histograms(X, Y, num_var = "income", cat_vars = c("gender"),
                   weight_X = "weight", weight_Y = "weight",
                   facet_type = "wrap", plot_type = "separate")

Description

This function performs a non-parametric Kolmogorov-Smirnov (KS) test to compare the cumulative distribution functions of a numeric variable between an original dataset (X) and an anonymized/synthetic dataset (Y). The test measures the maximum deviation between the two empirical cumulative distribution functions. The function can also perform the test separately for groups defined by one or more categorical variables.

Usage

compare_ks_test(X, ...)

## S3 method for class 'synth_pair'
compare_ks_test(X, ...)

## Default S3 method:
compare_ks_test(X, Y, num_var, cat_vars = NULL, ...)

Arguments

Value

A data.table with columns for the grouping variables (if provided), the KS test statistic, and the corresponding p-value. If no grouping is provided, a data.table with one row is returned.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  age = sample(20:80, 500, replace = TRUE),
  gender = sample(c("Male", "Female"), 500, replace = TRUE),
  region = sample(c("North", "South", "East", "West"), 500, replace = TRUE),
  income = rnorm(500, mean = 50000, sd = 10000)
)

Y <- data.frame(
  age = sample(20:80, 1000, replace = TRUE),
  gender = sample(c("Male", "Female"), 1000, replace = TRUE),
  region = sample(c("North", "South", "East", "West"), 500, replace = TRUE),
  income = rnorm(1000, mean = 48000, sd = 12000)
)

# Perform KS test on the entire datasets
compare_ks_test(X, Y, num_var = "income")

# Perform KS test within groups defined by gender
compare_ks_test(X, Y, num_var = "income", cat_vars = c("gender"))
compare_ks_test(X, Y, num_var = "income", cat_vars = c("gender","region"))

Description

This function compares the central tendencies of continuous variables and the frequency distributions of categorical variables between two datasets (e.g., an original dataset and an anonymized/synthetic dataset). For continuous variables, the function computes a set of summary statistics for each variable, including arithmetic mean, median, Huber mean, standard deviation, interquartile range (IQR), median absolute deviation (MAD), skewness, and kurtosis. Weighted estimates are computed if sampling weights are provided. In addition, conditional summaries can be produced if grouping variables are specified. For categorical variables, frequency and relative frequency tables are computed, using weights if available.

Usage

compare_means_frequencies(X, ...)

## S3 method for class 'synth_pair'
compare_means_frequencies(X, ...)

## Default S3 method:
compare_means_frequencies(
  X,
  Y,
  cont_vars,
  cat_vars,
  group_vars = NULL,
  weight_X = NULL,
  weight_Y = NULL,
  stats = c("mean", "median", "sd", "IQR", "mad", "huber", "skewness", "kurtosis"),
  ...
)

Arguments

Details

For continuous variables, if sampling weights are provided the function computes weighted estimates. The arithmetic mean and standard deviation are computed using weighted.mean and a weighted variance formula. The median is computed using Hmisc::wtd.quantile when weights are provided. The Huber mean is computed via MASS::huber (weighted Huber is not implemented). The IQR and MAD are computed using base functions. For skewness and kurtosis, if weights are provided, the function computes them using custom formulas; otherwise, it computes the unweighted versions.

For categorical variables, frequencies are computed as the sum of weights (if provided) or as raw counts, with relative frequencies calculated accordingly.

Value

A list containing three elements:

continuous_summary

A data.table with overall summary statistics for each continuous variable in X and Y.

conditional_summary

A data.table with conditional summary statistics computed by the grouping variables (if provided); otherwise, NULL.

categorical_summary

A data.table with frequency and relative frequency tables for each categorical variable in X and Y.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  age = sample(20:80, 500, replace = TRUE),
  income = rnorm(500, mean = 50000, sd = 10000),
  gender = sample(c("Male", "Female"), 500, replace = TRUE),
  region = sample(c("North", "South", "East", "West"), 500, replace = TRUE),
  weight = runif(500, 0.5, 1.5)
)

Y <- data.frame(
  age = sample(20:80, 1000, replace = TRUE),
  income = rnorm(1000, mean = 48000, sd = 12000),
  gender = sample(c("Male", "Female"), 1000, replace = TRUE),
  region = sample(c("North", "South", "East", "West"), 1000, replace = TRUE),
  weight = runif(1000, 0.5, 1.5)
)

result <- compare_means_frequencies(X, Y,
             cont_vars = c("age", "income"),
             cat_vars = c("gender", "region"),
             group_vars = c("region"),
             weight_X = "weight", weight_Y = "weight",
             stats = c("mean", "median", "sd", "IQR", "mad", "huber", "skewness", "kurtosis"))

print(result$continuous_summary)
print(result$conditional_summary)
print(result$categorical_summary)

Description

This function compares missing value patterns in an original dataset (X) and a synthetic/anonymized dataset (Y). It assesses whether the synthetic data mimics the same percentage and distribution of missingness across variables.

Usage

compare_missing_values(X, ...)

## S3 method for class 'synth_pair'
compare_missing_values(X, ...)

## Default S3 method:
compare_missing_values(X, Y, method = "percentage", ...)

Arguments

Details

Supported methods include:

• "percentage": Computes and compares missing value percentages for each variable.

• "pattern": Summarizes the frequency of missingness patterns across rows.

Both X and Y should be numeric data frames or matrices with the same structure.

Value

An object of class "missingCompare":

• For "percentage": a list with a data frame (summary) reporting per-variable missing counts and percentages.

• For "pattern": a list with two tables (patterns_X and patterns_Y) reporting the frequency of missingness patterns.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(a = c(rnorm(95), rep(NA, 5)), b = c(rnorm(90), rep(NA, 10)))
Y <- data.frame(a = c(rnorm(95), rep(NA, 5)), b = c(rnorm(90), rep(NA, 10)))

# Compare missing value percentages
res_pct <- compare_missing_values(X, Y, method = "percentage")
print(res_pct)
summary(res_pct)
plot(res_pct)

# Compare missing value patterns
res_pat <- compare_missing_values(X, Y, method = "pattern")
print(res_pat)
summary(res_pat)
plot(res_pat)

Description

This function evaluates and compares predictive model performance between an original dataset (X) and a synthetic/anonymized dataset (Y). The model is trained separately within each dataset using cross-validation, and various performance metrics are computed.

Usage

compare_model_performance(X, ...)

## S3 method for class 'synth_pair'
compare_model_performance(X, ...)

## Default S3 method:
compare_model_performance(
  X,
  Y,
  formula,
  method = "rpart",
  metric = ifelse(is.factor(X[[as.character(formula[[2]])]]), "Accuracy", "RMSE"),
  trControl = NULL,
  ...
)

Arguments

Value

List containing trained models, performance metrics for X and Y, and comparison results.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

if (requireNamespace("caret", quietly = TRUE)) {
library(caret)
set.seed(123)
X <- data.frame(income = rnorm(500, 50000, 10000),
                age = rnorm(500, 40, 10),
                gender = factor(sample(c("M", "F"), 500, replace = TRUE)),
                target = factor(sample(c("Yes", "No"), 500, replace = TRUE)))

Y <- data.frame(income = rnorm(500, 48000, 12000),
                age = rnorm(500, 42, 11),
                gender = factor(sample(c("M", "F"), 500, replace = TRUE)),
                target = factor(sample(c("Yes", "No"), 500, replace = TRUE)))

result <- compare_model_performance(X, Y,
           formula = target ~ income + age + gender,
           method = 'rpart',
           metric = 'Accuracy')

print(result$comparison)
}

Description

This function compares the joint distribution of a set of variables between two datasets (X and Y) using one of two methods. For continuous variables, it computes the weighted (if weights are provided) Mahalanobis distance between the mean vectors, using a pooled covariance matrix (or its generalized inverse). For nominal variables, it discretizes the data (if necessary) and computes a divergence based on the Jensen-Shannon divergence of the joint frequency distributions, which serves as a symmetric measure of distributional difference.

Usage

compare_multivariate_distribution(X, ...)

## S3 method for class 'synth_pair'
compare_multivariate_distribution(X, ...)

## Default S3 method:
compare_multivariate_distribution(
  X,
  Y,
  vars,
  method = "mahalanobis",
  weight_X = NULL,
  weight_Y = NULL,
  n_bins = 10,
  ...
)

Arguments

Details

For the "mahalanobis" method, if sampling weights are provided, the function computes the weighted mean vector for each dataset and a weighted pooled covariance matrix. The Mahalanobis distance is then computed as

$D_M = \sqrt{(\mu_X - \mu_Y)^T S^{-1} (\mu_X - \mu_Y)}$

where $S$ is the pooled covariance matrix (or its generalized inverse if singular). For the "mutual_information" method, numeric variables are discretized into n_bins bins using cut(), and joint frequency tables are constructed. The tables are normalized into probability vectors, and the Jensen-Shannon divergence is computed as

$JS(p,q) = \frac{1}{2} KL(p||m) + \frac{1}{2} KL(q||m)$

where $m = \frac{1}{2}(p+q)$ and $KL$ is the Kullback-Leibler divergence.

Value

A list containing:

distance

The computed distance measure (Mahalanobis distance or Jensen-Shannon divergence).

method

The method used ("mahalanobis" or "mutual_information").

additional

A list of additional computed quantities (e.g., weighted means and pooled covariance for Mahalanobis, or probability vectors for mutual_information).

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

# Continuous example using Mahalanobis distance (weighted)
set.seed(123)
X <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  age = rnorm(500, mean = 40, sd = 10),
  weight = runif(500, 0.5, 1.5)
)
Y <- data.frame(
  income = rnorm(1000, mean = 48000, sd = 12000),
  age = rnorm(1000, mean = 42, sd = 11),
  weight = runif(1000, 0.5, 1.5)
)
result_mahal <- compare_multivariate_distribution(X, Y, vars = c("income", "age"),
                                                  method = "mahalanobis",
                                                  weight_X = "weight", weight_Y = "weight")
print(result_mahal)

# Nominal example using mutual information (Jensen-Shannon divergence)
set.seed(456)
X_nom <- data.frame(
  gender = factor(sample(c("Male", "Female"), 500, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 500, replace = TRUE))
)
Y_nom <- data.frame(
  gender = factor(sample(c("Male", "Female"), 1000, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 1000, replace = TRUE))
)
result_js <- compare_multivariate_distribution(X_nom, Y_nom, vars = c("gender", "region"),
                                               method = "mutual_information", n_bins = 5)
print(result_js)

Description

This function computes higher-order summary statistics for multiple variables in two datasets (e.g., an original dataset and an anonymized/synthetic dataset) to assess whether the synthetic data replicate the joint properties of the original data. For continuous variables, the function calculates the multivariate mean vector, variance–covariance matrix, and correlation matrix. When sampling weights are provided, weighted estimates are computed. For categorical variables, joint frequency tables (and relative frequencies) are constructed. The function returns the statistics for each dataset along with the differences between them.

Usage

compare_multivariate_summary_statistics(X, ...)

## S3 method for class 'synth_pair'
compare_multivariate_summary_statistics(X, ...)

## Default S3 method:
compare_multivariate_summary_statistics(
  X,
  Y,
  cont_vars,
  cat_vars,
  weight_X = NULL,
  weight_Y = NULL,
  ...
)

Arguments

Details

For continuous variables, if sampling weights are provided the function computes weighted means using weighted.mean and weighted covariance matrices using a custom function. The correlation matrices are derived from the covariance matrices. For categorical variables, if weights are provided the frequencies are computed as the sum of weights; otherwise, raw counts are used. Joint frequency tables are normalized to obtain relative frequencies.

Value

A list with two elements:

continuous

A list containing:

mean_X

The multivariate mean vector for X.

mean_Y

The multivariate mean vector for Y.

mean_diff

The difference between the mean vectors (X minus Y).

cov_X

The variance–covariance matrix for X.

cov_Y

The variance–covariance matrix for Y.

cov_diff

The difference between the covariance matrices (X minus Y).

cor_X

The correlation matrix for X.

cor_Y

The correlation matrix for Y.

cor_diff

The absolute difference between the correlation matrices.

categorical

A list containing:

freq_X

The joint frequency table for the categorical variables in X.

freq_Y

The joint frequency table for the categorical variables in Y.

rel_freq_X

The relative frequency table for X.

rel_freq_Y

The relative frequency table for Y.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_outliers(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  age = rnorm(500, mean = 40, sd = 10),
  gender = factor(sample(c("Male", "Female"), 500, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 500, replace = TRUE)),
  weight = runif(500, 0.5, 1.5)
)

Y <- data.frame(
  income = rnorm(1000, mean = 48000, sd = 12000),
  age = rnorm(1000, mean = 42, sd = 11),
  gender = factor(sample(c("Male", "Female"), 1000, replace = TRUE)),
  region = factor(sample(c("North", "South", "East", "West"), 1000, replace = TRUE)),
  weight = runif(1000, 0.5, 1.5)
)

result <- compare_multivariate_summary_statistics(X, Y,
              cont_vars = c("income", "age"),
              cat_vars = c("gender", "region"),
              weight_X = "weight", weight_Y = "weight")

print(result$continuous)
print(result$categorical)

Description

This function compares the prevalence of outliers in an original dataset (X) and an anonymized/synthetic dataset (Y) using various outlier detection methods. Supported methods include:

• "zscore": Uses z-scores (default threshold = 3) to flag outliers.

• "iqr": Flags values outside 1.5 times the interquartile range.

• "dbscan": Uses DBSCAN clustering to identify noise points as outliers.

• "robust": Uses robust Mahalanobis distances based on the Minimum Covariance Determinant.

Usage

compare_outliers(X, ...)

## S3 method for class 'synth_pair'
compare_outliers(X, ...)

## Default S3 method:
compare_outliers(X, Y, method = "zscore", threshold = NULL, ...)

Arguments

Details

Input datasets must be numeric and of the same structure. For multivariate methods, the comparison is performed on the entire dataset.

Value

A list with:

• summary: A data frame summarizing outlier counts and proportions in X and Y.

• details: Method-specific details (e.g., per-variable counts or clustering output).

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_pca(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(a = rnorm(100), b = rnorm(100))
Y <- data.frame(a = rnorm(100, mean = 0.1), b = rnorm(100, mean = -0.1))

# Using z-score method (no extra packages required)
res_z <- compare_outliers(X, Y, method = "zscore", threshold = 3)

# Using IQR method
res_iqr <- compare_outliers(X, Y, method = "iqr")

# Using DBSCAN (requires the 'dbscan' package)
if (requireNamespace("dbscan", quietly = TRUE)) {
  res_db <- compare_outliers(X, Y, method = "dbscan", eps = 0.5, minPts = 5)
}

# Using robust Mahalanobis (requires the 'robustbase' package)
if (requireNamespace("robustbase", quietly = TRUE)) {
  res_robust <- compare_outliers(X, Y, method = "robust", threshold = 0.975)
}

Description

This function performs Principal Component Analysis (PCA) on two datasets (X and Y) for a specified set of numeric variables. When side_by_side is FALSE, a combined PCA is performed on the union of the data for visualization of the point scores, while separate PCA analyses are computed for X and Y to obtain loadings. The function then overlays loadings arrows from X (blue) and Y (red) onto the combined PCA biplot. When side_by_side is TRUE, separate biplots for X and Y are produced and arranged side by side.

Usage

compare_pca(X, ...)

## S3 method for class 'synth_pair'
compare_pca(X, ...)

## Default S3 method:
compare_pca(
  X,
  Y,
  vars,
  center = TRUE,
  scale = TRUE,
  biplot = TRUE,
  side_by_side = FALSE,
  ...
)

Arguments

Details

When side_by_side is FALSE, the function first adds a dataset indicator to X and Y, combines them, and performs PCA on the combined dataset for the point scores. Then, PCA is computed separately for X and Y. Loadings from these separate analyses are scaled using the range of the combined PCA scores and overlaid on the combined biplot with blue arrows for X and red arrows for Y.

Value

A list with two elements:

pca

If side_by_side = FALSE, a list with three PCA objects: combined, X, and Y. If side_by_side = TRUE, a list with two PCA objects: X and Y.

plot

A ggplot2 object: either a combined biplot with separate loadings or separate plots arranged side by side.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_wasserstein(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  age = rnorm(500, mean = 40, sd = 10)
)
# Add a positively correlated variable
X$expenses <- X$income * 0.5 + rnorm(n = 500, mean = 1000, sd = 500)
Y <- data.frame(
  income = rnorm(1000, mean = 48000, sd = 12000),
  age = rnorm(1000, mean = 42, sd = 11)
)
Y$expenses <- Y$income * 0.5 + rnorm(n = 1000, mean = 1000, sd = 500)
# Combined PCA biplot with separate loadings for X (blue) and Y (red)
res_combined <- compare_pca(X, Y,
  vars = c("income", "age", "expenses"),
  biplot = TRUE, side_by_side = FALSE)
print(res_combined$pca)
print(res_combined$plot)

# Separate PCA biplots for X and Y arranged side by side
res_separate <- compare_pca(X, Y, vars = c("income", "age"), biplot = TRUE, side_by_side = TRUE)
print(res_separate$pca)
print(res_separate$plot)

Description

This function computes the one-dimensional Wasserstein distance between a numeric variable in two datasets, X (original) and Y (anonymized or synthetic). For continuous variables, the Wasserstein distance is approximated by integrating the absolute differences between the quantile functions computed on a grid of probabilities. If sampling weights are provided, weighted quantiles are computed using Hmisc::wtd.quantile. For nominal (categorical) variables, the function computes the total variation distance (half the L1 distance between the probability distributions) assuming a unit cost for mismatches.

Usage

compare_wasserstein(X, ...)

## S3 method for class 'synth_pair'
compare_wasserstein(X, ...)

## Default S3 method:
compare_wasserstein(
  X,
  Y,
  num_var,
  cat_vars = NULL,
  weight_X = NULL,
  weight_Y = NULL,
  var_type = "auto",
  n_grid = 1000,
  ...
)

Arguments

Details

When grouping variables are provided (via cat_vars), the Wasserstein distance is computed within each group.

Value

A data.table with the computed Wasserstein distance for each group (or overall if no grouping is provided). For continuous variables, the Wasserstein distance is approximated by the mean absolute difference between the quantile functions. For nominal variables, the result is the total variation distance.

Author(s)

Matthias Templ

See Also

Other comparison: ci_overlap(), compare_boxplots(), compare_chisq_gof(), compare_correlation_matrices(), compare_distributions_cont(), compare_embedding(), compare_feature_importance(), compare_histograms(), compare_ks_test(), compare_means_frequencies(), compare_missing_values(), compare_model_performance(), compare_multivariate_distribution(), compare_multivariate_summary_statistics(), compare_outliers(), compare_pca(), densitydiff_1d_num(), densitydiff_kl_num(), densitydiff_pca()

Examples

set.seed(123)
X <- data.frame(
  income = rnorm(500, mean = 50000, sd = 10000),
  gender = sample(c("Male", "Female"), 500, replace = TRUE),
  region = sample(c("North", "South", "East", "West"), 500, replace = TRUE),
  weight = runif(500, 0.5, 1.5)
)

Y <- data.frame(
  income = rnorm(1000, mean = 48000, sd = 12000),
  gender = sample(c("Male", "Female"), 1000, replace = TRUE),
  region = sample(c("North", "South", "East", "West"), 1000, replace = TRUE),
  weight = runif(1000, 0.5, 1.5)
)

# Compare overall continuous distributions (weighted)
compare_wasserstein(X, Y, num_var = "income", var_type = "continuous",
                    weight_X = "weight", weight_Y = "weight")

# Compare distributions within groups defined by gender
compare_wasserstein(X, Y, num_var = "income", cat_vars = c("gender"),
                    var_type = "continuous", weight_X = "weight", weight_Y = "weight")
compare_wasserstein(X, Y, num_var = "income", cat_vars = c("gender","region"),
                    var_type = "continuous", weight_X = "weight", weight_Y = "weight")

# Example for nominal variables

# Create a synthetic dataset X with a nominal variable "gender"
set.seed(123)
X_nom <- data.frame(
  gender = sample(c("Male", "Female"), 500, replace = TRUE),
  weight = runif(500, 0.5, 1.5)
)

# Create a synthetic dataset Y with a nominal variable "gender"
Y_nom <- data.frame(
  gender = sample(c("Male", "Female"), 1000, replace = TRUE),
  weight = runif(1000, 0.5, 1.5)
)

# Compare the distributions of the nominal variable using the Wasserstein distance.
# Here, var_type is explicitly set to "nominal" so that the function uses frequency tables.
result_nominal <- compare_wasserstein(X_nom, Y_nom, num_var = "gender",
                                      var_type = "nominal",
                                      weight_X = "weight", weight_Y = "weight")

print(result_nominal)

Description

Computes confidence intervals for the RAPID disclosure risk measure (proportion of records at risk). Three methods are available: Wald (normal approximation), Wilson (score interval, recommended), and bootstrap (percentile).

Usage

## S3 method for class 'rapid'
confint(
  object,
  parm = NULL,
  level = 0.95,
  method = c("wilson", "wald", "bootstrap"),
  n_bootstrap = 1000,
  seed = NULL,
  ...
)

Arguments

Details

The RAPID score is a sample proportion (number at risk / total). Because the model is trained on separate synthetic data and applied to the original records, the per-record risk indicators are conditionally i.i.d. Bernoulli, making binomial confidence intervals valid.

wald

Standard normal approximation $\hat p \pm z_{\alpha/2}\sqrt{\hat p(1-\hat p)/n}$. Can produce intervals outside $[0,1]$.

wilson

Score interval (Wilson, 1927). Has better coverage near 0 and 1; recommended as the default.

bootstrap

Resamples the at_risk indicators from object$records and returns percentile confidence limits. No model refitting is needed.

Value

A $1 \times 2$ matrix with columns for the lower and upper confidence limits, following the standard confint convention.

Author(s)

Oscar Thees, Matthias Templ

References

Wilson, E.B. (1927). Probable Inference, the Law of Succession, and Statistical Inference. Journal of the American Statistical Association, 22(158), 209–212.

See Also

rapid, rapid_test

Other rapid: rapid(), rapid_synthesizer_cv(), rapid_test(), rapid_threshold_select()

Examples

# Small runnable example
set.seed(42)
X <- data.frame(
  age = sample(20:60, 80, replace = TRUE),
  sex = sample(c("M", "F"), 80, replace = TRUE),
  income = rnorm(80, 50000, 10000)
)
Y <- X
Y$income <- Y$income + rnorm(80, 0, 5000)
r <- rapid(X, Y, key_vars = c("age", "sex"),
           target_var = "income", model_type = "lm")
confint(r)
confint(r, method = "wald")


# With random forest model and bootstrap CI
r2 <- rapid(X, Y, key_vars = c("age", "sex"),
            target_var = "income", model_type = "rf")
confint(r2)
confint(r2, method = "bootstrap", n_bootstrap = 500)

Description

Compares bivariate contingency tables (joint distributions) for all pairs of categorical variables between original and synthetic data. For each pair, the total variation (TV) distance between the proportion tables is computed. The mean TV distance across all pairs summarizes how well the synthetic data preserves the bivariate categorical dependence structure.

Usage

contingency_fidelity(X, ...)

## S3 method for class 'synth_pair'
contingency_fidelity(X, ...)

## Default S3 method:
contingency_fidelity(X, Y, vars = NULL, na.rm = TRUE, ...)