This vignette introduces the sdcLog package and its main functions, illustrated with various examples. sdcLog provides tools which simplify statistical disclosure control in the context of research data centers (RDC).

The package includes four main functions:

`sdc_descriptives()`

: This function is used for statistical disclosure control of descriptive statistics. It checks the data used for the calculation of descriptive statistics for compliance with the rules and regulations set by the RDC.The calculation of simple extreme values such as minimum and maximum is usually not allowed according to RDC rules, so

`sdc_descriptives()`

cannot be used to check them. Extreme values can only be used if they are calculated with`sdc_min_max()`

.`sdc_min_max()`

: This function is used for the automatic calculation of extreme values according to the rules of the RDC (if possible). It uses the available data and calculates the values for desired variables and groupings in compliance with the rules. The values are calculated as averages of a sufficient number of distinct entities. This helps the researcher to easily follow the rules and simplifies the output control.`sdc_model()`

: This function is used for statistical disclosure control for various types of models such as`lm()`

or`glm()`

. It checks the calculated model and the underlying data for compliance with the rules and regulations of the RDC.`sdc_log()`

: This function is a simple wrapper around`source()`

which makes it easy to run scripts and capture all output (especially output from the other`sdc_*`

functions) in a log file. Usually, this should be used to source R scripts containing one or more of the functions above.

This function performs statistical disclosure control according to two main criteria: On the one hand, it checks for a sufficiently large number of different statistical entities. On the other hand, it checks for dominance, which means that two entities must not account for more than 85 percent of the observed values. How to use `sdc_descriptives()`

is shown below.

To introduce `sdc_descriptives()`

, a simple toy dataset is used. There are 20 observations of 10 distinct entities from two different sectors and values in the years 2019 and 2020 for the variables `val_1`

and `val_2`

.

```
data("sdc_descriptives_DT")
sdc_descriptives_DT#> id id_na sector year val_1 val_2
#> 1: A A S1 2019 NA 9.477642
#> 2: A <NA> S1 2020 94.174449 5.856641
#> 3: B B S1 2019 4.349115 3.697140
#> 4: B B S1 2020 2.589011 6.796527
#> 5: C C S1 2019 6.155680 7.213390
#> 6: C C S1 2020 7.183206 5.948330
#> 7: D D S1 2019 9.173870 3.004441
#> 8: D D S1 2020 4.456933 0.000000
#> 9: E E S1 2019 2.815137 4.137765
#> 10: E E S1 2020 7.928573 0.000000
#> 11: F F S2 2019 9.085507 5.088913
#> 12: F <NA> S2 2020 180.816675 0.000000
#> 13: G G S2 2019 9.502077 2.107123
#> 14: G G S2 2020 7.458567 0.000000
#> 15: H H S2 2019 6.947180 5.059104
#> 16: H H S2 2020 9.927155 3.489741
#> 17: I I S2 2019 6.662026 8.957527
#> 18: I I S2 2020 4.420317 8.618987
#> 19: J J S2 2019 1.556076 4.722792
#> 20: J J S2 2020 7.997007 7.347734
```

Consider the case that the mean for `val_1`

has been calculated and is now to be output as a result:^{1}

```
mean = mean(val_1, na.rm = TRUE))]
sdc_descriptives_DT[, .(#> mean
#> 1: 20.16835
```

Before this result can be released, it must be checked whether all RDC rules for calculating this value have been followed. Thus, the underlying data is checked for compliance with the RDC rules.

This is the simplest case, the descriptive statistic (mean) was calculated for the variable `val_1`

without further specifications. Required arguments of `sdc_descriptives()`

are the data set (`data`

), the ID variable (`id_var`

) and the variable for which the statistics were calculated (`val_var`

):

```
sdc_descriptives(data = sdc_descriptives_DT, id_var = "id", val_var = "val_1")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1 | zero_as_NA: FALSE
#> Output complies to RDC rules.
```

Since there are no problems at this point, the function runs without warnings and returns (invisibly) a list of information containing options, settings and the checked criteria `distinct_ids`

and `dominance`

.

Options and settings are always printed to show that all specifications are set according to RDC rules. From the output above follows that there are at least 5 distinct entities required (`sdc.n_ids: 5`

) and that dominance is defined as 2 entities (`sdc.n_ids_dominance: 2`

) with a value share of more than 85 percent (`sdc.share_dominance: 0.85`

). This reflects the standard values for the options. For details on setting options see the separate vignette on options.

The settings show again which arguments were specified in the function call and vary depending on the `sdc_function`

. This is important if the result from `sdc_descriptives()`

is not printed right away.

In this and the following section some advanced cases are presented to introduce more arguments and functionalities of `sdc_descriptives()`

.

In this case the descriptive statistics for the variable `val_1`

are grouped by `sector`

:

```
mean = mean(val_1, na.rm = TRUE)), by = "sector"]
sdc_descriptives_DT[, .(#> sector mean
#> 1: S1 15.42511
#> 2: S2 24.43726
```

The mean is computed grouped by sector, so the grouping variable must be specified in `by`

. Checking the results leads to the following:

```
sdc_descriptives(data = sdc_descriptives_DT, id_var = "id", val_var = "val_1", by = "sector")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1 | by: sector | zero_as_NA: FALSE
#> Output complies to RDC rules.
```

The grouped descriptive statistics by sector do not generate a warning and therefore comply with RDC rules. Therefore, the results could be released in this case.

In order to extend this case even further, it is now proposed to group the mean of `val_1`

not only by `sector`

, but also by `year`

:

```
mean = mean(val_1, na.rm = TRUE)), by = c("sector", "year")]
sdc_descriptives_DT[, .(#> sector year mean
#> 1: S1 2019 5.623451
#> 2: S1 2020 23.266434
#> 3: S2 2019 6.750574
#> 4: S2 2020 42.123944
```

To check this result for compliance with RDC rules, use:

```
sdc_descriptives(
data = sdc_descriptives_DT,
id_var = "id",
val_var = "val_1",
by = c("sector", "year")
)#> Warning: DISCLOSURE PROBLEM: Not enough distinct entities.
#> Warning: DISCLOSURE PROBLEM: Dominant entities.
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1 | by: c("sector", "year") | zero_as_NA: FALSE
#> Not enough distinct entities:
#> sector year distinct_ids
#> 1: S1 2019 4
#> 2: S1 2020 5
#> 3: S2 2019 5
#> 4: S2 2020 5
#> Dominant entities:
#> sector year value_share
#> 1: S2 2020 0.9056314
#> 2: S1 2020 0.8776852
#> 3: S1 2019 0.6815011
#> 4: S2 2019 0.5506965
```

Now several warnings appear, as both criteria are violated. For sector `S1`

there are not enough distinct ids in year 2019, as there is a missing value in the data. The dominance criterion for year 2020 is violated in both sectors. As can be seen in the table displayed, the value share of approximately 88 percent for `S1`

and 91 percent for `S2`

are above the 85 percent limit. Therefore, the descriptive statistics for `val_1`

, grouped by `sector`

and `year`

cannot be released.

Now, descriptive statistics are calculated for variable `val_2`

and grouped by sector and year:

```
mean = mean(val_2, na.rm = TRUE)), by = c("sector", "year")]
sdc_descriptives_DT[, .(#> sector year mean
#> 1: S1 2019 5.506076
#> 2: S1 2020 3.720300
#> 3: S2 2019 5.187092
#> 4: S2 2020 3.891292
```

The compliance with the rules can be checked just as in the previous case (only replacing `val_1`

by `val_2`

):

```
sdc_descriptives(
data = sdc_descriptives_DT,
id_var = "id",
val_var = "val_2",
by = c("sector", "year")
)#> A share of 0.2 of 'val_var' are zero. These will be treated as 'NA'.
#> To prevent this behavior and / or avoid this message, set 'zero_as_NA' explicitly.
#> Warning: DISCLOSURE PROBLEM: Not enough distinct entities.
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_2 | by: c("sector", "year") | zero_as_NA: TRUE
#> Not enough distinct entities:
#> sector year distinct_ids
#> 1: S1 2020 3
#> 2: S2 2020 3
#> 3: S1 2019 5
#> 4: S2 2019 5
```

The result indicates that problems exist and the output does not comply to the rules. There are not enough distinct entities and the output cannot be released like this.

An additional message indicates that the value `0`

occurs rather frequently in the data (20 percent of all cases). The message indicates that `0`

is assumed to represent missing values and will be treated as such. Please note that even if `0`

s are actual `0`

s in the data, this assumption might be correct in the context of statistical disclosure control. For example, if most of the cases are `0`

, it might be known publicly which entities do not have a value of `0`

for this specific variable. So treating those `0`

as `NA`

is correct in this context. Since this is the more defensive interpretation of `0`

s, it’s the default.

However, it might be the case that it is accurate according to the data basis to treat values of `0`

as zero (instead of `NA`

). Then, specifying the argument `zero_as_NA = FALSE`

circumvents the default behavior and treats `0`

like other numeric values:

```
sdc_descriptives(
data = sdc_descriptives_DT,
id_var = "id",
val_var = "val_2",
by = c("sector", "year"),
zero_as_NA = FALSE
)#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_2 | by: c("sector", "year") | zero_as_NA: FALSE
#> Output complies to RDC rules.
```

Now `0`

is not recognized as `NA`

anymore. In this case the criterion of distinct entities is not longer violated. Therefore, the output could be released (assuming it is actually correct to treat `0`

s as usual numeric values).

This function automatically calculates extreme values that comply with the rules of the RDC. It checks the criteria of distinct entities and dominance. The values are calculated as averages of a sufficiently large number of observations. It is based on an iterative procedure that aggregates data until there are enough distinct entities to calculate the extreme values and no problems with dominance occur.

The function always starts the iteration process with the lowest possible number of observations for each extreme value (here `5`

, since at least five distinct statistical units must be included in the calculation according to the rules of the RDC). Furthermore, the function checks that the subsets of data for minimum and maximum do not overlap.

If there are no problems with the calculation, the function returns a `data.table`

with the extreme values. Maximum and minimum are always output together, none of the two can be calculated separately. If it is not possible to calculate extreme values under these criteria, a corresponding message is printed and the result is filled with `NA`

.

To introduce `sdc_min_max()`

, another simple dataset is used. We have 20 observations of 10 different entities, for which the corresponding sector is given and values for the variables `val_1`

, `val_2`

, `val_3`

in the years 2019 and 2020, respectively.

```
data("sdc_min_max_DT")
sdc_min_max_DT#> id sector year val_1 val_2 val_3
#> 1: A S1 2019 20 200 NA
#> 2: B S1 2020 19 190 19
#> 3: C S1 2019 18 18 18
#> 4: D S1 2020 17 17 17
#> 5: E S1 2019 16 16 16
#> 6: F S1 2020 15 15 15
#> 7: G S1 2019 14 14 14
#> 8: H S1 2020 13 13 13
#> 9: I S1 2019 12 12 12
#> 10: J S1 2020 11 11 11
#> 11: A S2 2019 10 10 10
#> 12: B S2 2020 9 9 9
#> 13: C S2 2019 8 8 8
#> 14: D S2 2020 7 7 7
#> 15: E S2 2019 6 6 6
#> 16: F S2 2020 5 5 5
#> 17: G S2 2019 4 4 4
#> 18: H S2 2020 3 3 3
#> 19: I S2 2019 2 2 2
#> 20: J S2 2020 1 1 1
```

In this simple case, extreme values should be calculated for variable `val_1`

. This can be done with `sdc_min_max()`

by specifying the dataset (`data`

), the id variable (`id_var`

) and the variable for which extreme values are to be calculated (`val_var`

).

```
sdc_min_max(data = sdc_min_max_DT, id_var = "id", val_var = "val_1")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1
#> val_var min distinct_ids_min max distinct_ids_max
#> 1: val_1 3 5 18 5
```

Since no problems occur, the function (invisibly) returns a list with the options, settings and extreme values and prints the calculated extreme values. As shown in the output, the extreme values could be calculated and 5 distinct entities were used for each value. Thus, no additional entities had to be included in the calculation.

In this case minimum and maximum values are to be calculated for variable `val_2`

:

```
sdc_min_max(data = sdc_min_max_DT, id_var = "id", val_var = "val_2")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_2
#> val_var min distinct_ids_min max distinct_ids_max
#> 1: val_2 3 5 67.14286 7
```

When we look at the output, we see that values from 5 distinct entities were used to calculate the minimum and 7 distinct entities to calculate the maximum. This is because the dominance criterion would be violated if only 5 distinct entities were considered for the maximum. Thus, 7 distinct entities are automatically taken into account.

If you specify `max_obs = 5`

, there is no feasible solution:

```
sdc_min_max(data = sdc_min_max_DT, id_var = "id", val_var = "val_2", max_obs = 5)
#> It is impossible to compute extreme values for variable 'val_2' that comply to RDC rules.
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_2
#> val_var min distinct_ids_min max distinct_ids_max
#> 1: val_2 NA NA NA NA
```

Note that `max_obs`

controls the maximum number of observations, not distinct entities.

It is also possible to calculate minimum and maximum values by groups. In the following, these are calculated by `year`

and `sector`

, separately.

```
sdc_min_max(data = sdc_min_max_DT, id_var = "id", val_var = "val_1", by = "year")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1 | by: year
#> val_var year min distinct_ids_min max distinct_ids_max
#> 1: val_1 2019 6 5 16 5
#> 2: val_1 2020 5 5 15 5
sdc_min_max(data = sdc_min_max_DT, id_var = "id", val_var = "val_1", by = "sector")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1 | by: sector
#> val_var sector min distinct_ids_min max distinct_ids_max
#> 1: val_1 S1 13 5 18 5
#> 2: val_1 S2 3 5 8 5
```

No problems occur, so minimum and maximum values are calculated and shown for each group.

This can also be done for several grouping variables. In the following, extreme values for variable `val_1`

are to be calculated by `year`

and `sector`

.

```
<- sdc_min_max(
res data = sdc_min_max_DT,
id_var = "id",
val_var = "val_1",
by = c("sector", "year")
)#> It is impossible to compute extreme values for variable 'val_1' that comply to RDC rules.
```

Now a message occurs, explaining that RDC rules would be violated for the calculation of these values. For programming purposes, please note that the structure of the resulting `data.table`

remains the same (but is filled with `NA`

:

```
# extreme_vals
res#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id | val_var: val_1 | by: c("sector", "year")
#> val_var sector year min distinct_ids_min max distinct_ids_max
#> 1: val_1 S1 2019 NA NA NA NA
#> 2: val_1 S1 2020 NA NA NA NA
#> 3: val_1 S2 2019 NA NA NA NA
#> 4: val_1 S2 2020 NA NA NA NA
```

This function checks if your model complies to RDC rules. The criterion of distinct entities is also checked here. In addition, it is checked whether there are enough different entities for each attribute or value level. For continuous variables, `sdc_model()`

distinguishes between `<zero>`

and `<non-zero>`

values. The function can be used to check a broad range of models like `lm`

, `glm`

and various others. In fact, anything which can be handled by `broom::augment()`

can also be handled by `sdc_model()`

. For a list of supported models see `?generics::augment`

.

To introduce `sdc_model()`

, another dataset with different variables is used, which includes dummy-variables.

We have 80 observations of 10 different entities for the variables `y`

, `x_1`

, `x_2`

, `x_3`

, `x_4`

and additional information on sector, year and country (dummy variables). A summary of the data set is given below.

```
data("sdc_model_DT")
print(skim(sdc_model_DT))
#> ── Data Summary ────────────────────────
#> Values
#> Name sdc_model_DT
#> Number of rows 80
#> Number of columns 9
#> Key id
#> _______________________
#> Column type frequency:
#> factor 4
#> numeric 5
#> ________________________
#> Group variables None
#>
#> ── Variable type: factor ──────────────────────────────────────────────────────────────────────
#> skim_variable n_missing complete_rate ordered n_unique top_counts
#> 1 id 0 1 FALSE 10 A: 8, B: 8, C: 8, D: 8
#> 2 dummy_1 0 1 FALSE 2 S1: 40, S2: 40
#> 3 dummy_2 0 1 FALSE 8 Y1: 10, Y2: 10, Y3: 10, Y4: 10
#> 4 dummy_3 0 1 FALSE 4 ES: 36, BE: 20, DE: 20, FR: 4
#>
#> ── Variable type: numeric ─────────────────────────────────────────────────────────────────────
#> skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
#> 1 y 0 1 121. 7.21 102. 116. 121. 125. 139. ▁▅▇▆▁
#> 2 x_1 0 1 122. 13.3 95.9 111. 122. 133. 152. ▅▇▇▇▂
#> 3 x_2 0 1 116. 19.4 86.1 99.1 113. 128. 152. ▇▇▆▃▆
#> 4 x_3 48 0.4 126. 59.0 36.0 81.3 120. 173. 230. ▆▆▅▇▃
#> 5 x_4 0 1 2544. 5944. 5.98 91.4 164. 227. 25059. ▇▁▁▁▁
```

Various simple linear models are specified from this dataset for illustration purposes.

```
<- lm(y ~ x_1 + x_2, data = sdc_model_DT)
model_1 <- lm(y ~ x_1 + x_2 + x_3, data = sdc_model_DT)
model_2 <- lm(y ~ x_1 + x_2 + dummy_1 * dummy_2, data = sdc_model_DT)
model_3 <- lm(y ~ x_1 + x_2 + dummy_1 * dummy_3, data = sdc_model_DT) model_4
```

These models are now checked for compliance with the rules of the RDC. It is checked if there are enough distinct entities in the whole model and if every level of each variable is checked for compliance with the rules.

A selection of problematic and unproblematic models has been made to better explain the differences. To check for compliance, the model object (`model`

), the data used (`data`

) and the ID variable (`id_var`

) must be specified in `sdc_model()`

.

A check of `model_1`

and `model_3`

is shown below.

```
sdc_model(data = sdc_model_DT, model = model_1, id_var = "id")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id
#> Output complies to RDC rules.
sdc_model(data = sdc_model_DT, model = model_3, id_var = "id")
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id
#> Output complies to RDC rules.
```

As we see, there are no problems and the models could be released as output. Note that `sdc_log()`

supports the interaction term in `model_3`

.

Now we turn to the problematic cases. We are checking the models `model_2`

and `model_4`

:

```
sdc_model(data = sdc_model_DT, model = model_2, id_var = "id")
#> Warning: DISCLOSURE PROBLEM: Not enough distinct entities.
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id
#> Not enough distinct entities:
#> distinct_ids
#> 1: 4
#> $x_1
#> Not enough distinct entities:
#> x_1 distinct_ids
#> 1: <non-zero> 4
#>
#> $x_2
#> Not enough distinct entities:
#> x_2 distinct_ids
#> 1: <non-zero> 4
#>
#> $x_3
#> Not enough distinct entities:
#> x_3 distinct_ids
#> 1: <non-zero> 4
```

Some difficulties occur with these models, but which?

`model_2`

leads to problems with the number of distinct entities. This problem arises with the inclusion of variable `x_3`

due to a high number of `NA`

s.

```
sdc_model(data = sdc_model_DT, model = model_4, id_var = "id")
#> Warning: DISCLOSURE PROBLEM: Not enough distinct entities.
#> OPTIONS: sdc.n_ids: 5 | sdc.n_ids_dominance: 2 | sdc.share_dominance: 0.85
#> SETTINGS: id_var: id
#> $dummy_3
#> Not enough distinct entities:
#> dummy_3 distinct_ids
#> 1: FR 4
#> 2: BE 10
#> 3: DE 10
#> 4: ES 10
#>
#> $`dummy_1:dummy_3`
#> Not enough distinct entities:
#> dummy_1:dummy_3 distinct_ids
#> 1: S2:FR 4
#> 2: S1:BE 10
#> 3: S1:DE 10
#> 4: S2:ES 10
```

For `model_4`

the problem stems from a small number of distinct entities for the value level `FR`

of `dummy_3`

. This also leads to a problem in the interaction term. Therefore the respective coefficients cannot be released either. Please note that this last case is probably the most common problem to occur when checking models.

This function serves to create Stata-like log files from R Scripts. The function is called to wrap an R script containing your analysis to write the corresponding code and console output into a log file. It can handle single files or a list of files at once.

A character vector containing the path(s) of the R script(s) which should be run must be specified as well as a character vector containing the path(s) of the text file(s) where the log(s) should be stored. To replace existing log files, one can specify the argument `replace = TRUE`

.

A simple call of this function could look as follows:

```
sdc_log(
r_scripts = "/home/my_project/R/my_script.R",
log_files = "/home/my_project/log/my_script.txt"
)
```

Even though this seems trivial, creating logs for scripts is essential because a log file bundles all information needed by the RDC for output control.

Since sdcLog heavily relies on

`data.table`

, all examples will use`data.table`

syntax as well.↩︎