CNAqc uses peak-detection algorithms to QC data; all leverage the idea that VAFs peaks are known for mutations mapped to a segment with given minor/ major allele copies. CNAqc therefore computes expected peaks, and compares them to peaks detected from data. The theory works with minor modifications for both clonal and subclonal segments. Three distinct algorithms are available, each one working with a different type of copy number segment; all analyses are called by this function, which takes care of running all the suitable algorithms based on the input data.
* Simple clonal segments (1:0, 2:0, 1:1, 2:1, 2:2). This QC measures an error for the precision of the current purity estimate, failing a whole sample or a subset of segments the value is over a desired maximum value. The error is determined as a linear combination from the distance between VAF peaks and their theoretical expectation. For this analysis, all mutations mapping across any segment with the same major/minor alleles are pooled. Note that this score can be used to select among alternative copy number solutions, i.e., favouring a solution with lower score. The peaks are determined i) via peak-detection algorithms from the peakPick package, applied to a Gaussian kernel density estimate (gKDE) smooth of the VAF distribution, and ii) via the Bmix Binomial mixture model. Peak-matching (i.e., determining what data peak is closest to the expected peak) has two possible implementations: one matching the closest peaks by euclidean distance, the other ranking peaks from higher to lowr VAFs, and prioritising the former.
* Complex clonal segments. The QC procedure for these “general” segments uses only the KDE and, as for simple segments, pools all mutations mapping across any segment with the same major/minor alleles. In this case,t no segment-level or sample-level scores are produced, and complex segment with many matched peaks is likely to be correct.
* Subclonal simple segments. The QC procedure for these segments uses the KDE and considers 2 subclones with distinct mixing proportions. Differently from clonal CNAs, however, here the analysis is carried out at the level of each segment, i.e., without pooling segments with the same karyotypes. This makes it possible to use subclonal calls from callers that report segment-specific CCF values, e.g., Battenberg. The model in CNAqc ranks the proposed evolutionary alternatives (linear versus branching) based on the number of matched peaks. A subclonal segment with many matched peaks is likely to be correct.
Results from peak-based QC are available via plot_peaks_analysis, and stored
inside the input object.
analyze_peaks(
x,
karyotypes = c("1:0", "1:1", "2:0", "2:1", "2:2"),
min_karyotype_size = 0,
min_absolute_karyotype_mutations = 100,
p_binsize_peaks = 0.005,
matching_epsilon = NULL,
purity_error = 0.05,
VAF_tolerance = 0.015,
n_bootstrap = 1,
kernel_adjust = 1,
matching_strategy = "closest",
KDE = TRUE,
starting_state_subclonal_evolution = "1:1",
cluster_subclonal_CCF = FALSE
)A CNAqc object.
For clonal simple CNAs, the list of segments to test; by default LOH regions (A, AA),
diploid regions (AB), and amplification regions (AAB, AABB) are tested, corresponding to
'1:0', '1:1', '2:1', '2:0', '2:2' in "Major:minor" notation.
For clonal simple CNAs, a filter for the segments to test. The segment size is defined based on the number of mutations mapped, this cut is on the proportion relative to the whole set of segments one wishes to analyse (defined by `karyotypes`). For example, by setting `min_karyotype_size = 0.2` one would QC clonal simple CNAs that contain at least 20 The default of this parameter is `0` (all QCed).
For clonal simple CNAs, as min_karyotype_size but with a cut
measured on absolute mutation counts. For example, by setting `min_absolute_karyotype_mutations = 150` one
would QC clonal simple CNAs that contain at least `150` mutations. The default of this parameter is `100`.
For clonal simple CNAs, peaks detected will be filtered if, in a peak, we map
less than p_binsize_peaks * N mutations. The value N is obtained couting all mutations that map
in all peaks. By default this parameters is `0.005`.
Deprecated parameter.
For clonal simple CNAs, the purity error tolerance to determine QC pass or fail. This can be
set automatically using function auto_tolerance to optimise the analysis based on a desired rate of false
positives matches, as a function of the data coverage and (putative) purity.
For clonal simple CNAs, a tolerance in comparing bands overlaps which is applied to the raw VAF values.
For clonal simple CNAs, the number of times peak detection is bootstrapped (by default 1). This helps sometimes finding peaks that might be visually observable but fail to be detected by the underlying peak-detection heuristics.
For KDE-based matches the adjust density parameter; see density. Note that a
Gaussian kernel is used by setting (kernel = 'gaussian').
For clonal simple CNAs, if "closest" the closest peak will be used to
match the expected peak. If "rightmost" peaks are matched prioritizing
right to left peaks (the higher-VAF gets matched first); this strategy is more correct
in principle but works only if there are no spurious peaks in the estimated
density. By default the "closest" strategy is used.
Deprecated parameter.
For subclonal simple CNAs, the starting state to determine linear versus branching evolutionary models. By default this is an heterozygous diploid `1:1` state.
For subclonal segments, should the tool try to merge segments with similar CCF and the same copy number alteration?
An object of class cnaqc, modified to hold the results from this analysis. For every type
of segment analyzed tables with summary peaks are available in x$peaks_analysis. The most helpful table
is usually the one for simple clonal CNAs `x$peaks_analysis$matches`, which reports several information:
- `mutation_multiplicity`, the number of copies of the mutation (i.e., a phasing information); - `peak`, `x`, `y` the expected peak, and the matched peak (`x` and `y`); - `offset`, `weight` and `score`, the factors of the final score; - `QC`, a pass/fail status for the peak.
The overall sample-level QC result is available in `x$peaks_analysis$QC`.
auto_tolerance, plot_peak_analysis and plot_QC.
data('example_dataset_CNAqc', package = 'CNAqc')
x = init(mutations = example_dataset_CNAqc$mutations, cna = example_dataset_CNAqc$cna, purity = example_dataset_CNAqc$purity)
#>
#> ── CNAqc - CNA Quality Check ───────────────────────────────────────────────────
#>
#> ℹ Using reference genome coordinates for: GRCh38.
#> ✔ Found annotated driver mutations: TTN, CTCF, and TP53.
#> ✔ Fortified calls for 12963 somatic mutations: 12963 SNVs (100%) and 0 indels.
#> ! CNAs have no CCF, assuming clonal CNAs (CCF = 1).
#> ✔ Fortified CNAs for 267 segments: 267 clonal and 0 subclonal.
#> ✔ 12963 mutations mapped to clonal CNAs.
# Note the run outputs
x = analyze_peaks(x)
#>
#> ── Peak analysis: simple CNAs ──────────────────────────────────────────────────
#>
#> ℹ Analysing 9041 mutations mapping to karyotype(s) 2:2 and 2:1.
#> ℹ Mixed type peak detection for karyotype 2:1 (1563 mutations)
#> Warning: replacing previous import ‘cli::num_ansi_colors’ by ‘crayon::num_ansi_colors’ when loading ‘BMix’
#> Warning: replacing previous import ‘crayon::%+%’ by ‘ggplot2::%+%’ when loading ‘BMix’
#> ✔ Loading BMix, 'Binomial and Beta-Binomial univariate mixtures'. Support : <https://caravagnalab.github.io/BMix/>
#> ℹ Mixed type peak detection for karyotype 2:2 (7478 mutations)
#> # A tibble: 4 × 16
#> # Rowwise:
#> mutation_multiplicity karyotype peak delta_vaf x y counts_per_bin
#> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <int>
#> 1 1 2:2 0.235 0.00700 0.24 0.86 65
#> 2 2 2:2 0.471 0.0140 0.461 8.14 625
#> 3 1 2:1 0.308 0.0120 0.307 3.60 59
#> 4 2 2:1 0.616 0.0239 0.602 2.94 60
#> # ℹ 9 more variables: discarded <lgl>, from <chr>, offset_VAF <dbl>,
#> # offset <dbl>, weight <dbl>, epsilon <dbl>, VAF_tolerance <dbl>,
#> # matched <lgl>, QC <chr>
#> ✔ Peak detection PASS with r = 0.00647677487483768 - maximum purity error ε = 0.05.
#> Joining with `by = join_by(Major, minor)`
#> Joining with `by = join_by(karyotype)`
#>
#> ── Peak analysis: complex CNAs ─────────────────────────────────────────────────
#>
#> ℹ Karyotypes 4:2, 3:2, and 3:0 with >100 mutation(s). Using epsilon = 0.05.
#> # A tibble: 3 × 5
#> # Groups: karyotype, matched [3]
#> karyotype n matched mismatched prop
#> <chr> <table[1d]> <int> <dbl> <dbl>
#> 1 3:0 312 3 0 1
#> 2 3:2 1625 3 0 1
#> 3 4:2 1893 4 0 1
#> Adding missing grouping variables: `matched`
#> Joining with `by = join_by(Major, minor, QC_PASS)`
#> Adding missing grouping variables: `matched`
#> Joining with `by = join_by(karyotype, QC_PASS)`
#>
#> ── Peak analysis: subclonal CNAs ───────────────────────────────────────────────
#>
#> ℹ No subclonal CNAs in this sample.
# More precise messages
print(x)
#> ── [ CNAqc ] MySample 12963 mutations in 267 segments (267 clonal, 0 subclonal).
#>
#> ── Clonal CNAs
#>
#> 2:2 [n = 7478, L = 1483 Mb] ■■■■■■■■■■■■■■■■■■■■■■■■■■■ { CTCF }
#> 4:2 [n = 1893, L = 331 Mb] ■■■■■■■
#> 3:2 [n = 1625, L = 357 Mb] ■■■■■■
#> 2:1 [n = 1563, L = 420 Mb] ■■■■■■ { TTN }
#> 3:0 [n = 312, L = 137 Mb] ■
#> 2:0 [n = 81, L = 39 Mb] { TP53 }
#> 16:2 [n = 4, L = 0 Mb]
#> 25:2 [n = 2, L = 1 Mb]
#> 3:1 [n = 2, L = 1 Mb]
#> 106:1 [n = 1, L = 0 Mb]
#>
#> ℹ Sample Purity: 89% ~ Ploidy: 4.
#>
#> ℹ There are 3 annotated driver(s) mapped to clonal CNAs.
#> chr from to ref alt DP NV VAF driver_label is_driver
#> chr2 179431633 179431634 C T 117 77 0.6581197 TTN TRUE
#> chr16 67646006 67646007 C T 120 54 0.4500000 CTCF TRUE
#> chr17 7577106 7577107 G C 84 78 0.9285714 TP53 TRUE
#>
#> ── PASS Peaks QC closest: 200%, λ = 0.0065. Purity correction: 1%. ───────────
#>
#> ℹ 2:2 ~ n = 7478 ( 83%) → PASS -0.034 PASS 0.035
#> ℹ 2:1 ~ n = 1563 ( 17%) → PASS 0.003 PASS 0.028
#>
#> ── General peak QC (3830 mutations): PASS 10 FAIL 0 - epsilon = 0.05. ──────
#>
#> ℹ 3:0 ~ n = 312 ( 8%) → PASS 3 FAIL 0
#> ℹ 3:2 ~ n = 1625 ( 42%) → PASS 3 FAIL 0
#> ℹ 4:2 ~ n = 1893 ( 49%) → PASS 4 FAIL 0
# The tabulars with summary results per peak and segment
print(x$peaks_analysis)
#> $score
#> [1] 0.006476775
#>
#> $fits
#> $fits$`2:1`
#> $fits$`2:1`$xy_peaks
#> # A tibble: 6 × 5
#> x y counts_per_bin discarded from
#> <dbl> <dbl> <int> <lgl> <chr>
#> 1 0.08 0.27 5 FALSE KDE
#> 2 0.3 3.6 78 FALSE KDE
#> 3 0.6 2.94 60 FALSE KDE
#> 4 0.108 0.246 2 FALSE BMix
#> 5 0.602 2.94 60 FALSE BMix
#> 6 0.307 3.60 59 FALSE BMix
#>
#> $fits$`2:1`$density
#>
#> Call:
#> density.default(x = y, adjust = kernel_adjust, kernel = "gaussian", na.rm = T)
#>
#> Data: y (1563 obs.); Bandwidth 'bw' = 0.0341
#>
#> x y
#> Min. :-0.05578 Min. :0.00016
#> 1st Qu.: 0.17384 1st Qu.:0.21293
#> Median : 0.40346 Median :0.51690
#> Mean : 0.40346 Mean :1.08768
#> 3rd Qu.: 0.63308 3rd Qu.:1.99965
#> Max. : 0.86270 Max. :3.60365
#>
#> $fits$`2:1`$matching
#> # A tibble: 2 × 16
#> # Rowwise:
#> mutation_multiplicity karyotype peak delta_vaf x y counts_per_bin
#> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <int>
#> 1 1 2:1 0.308 0.0120 0.307 3.60 59
#> 2 2 2:1 0.616 0.0239 0.602 2.94 60
#> # ℹ 9 more variables: discarded <lgl>, from <chr>, offset_VAF <dbl>,
#> # offset <dbl>, weight <dbl>, epsilon <dbl>, VAF_tolerance <dbl>,
#> # matched <lgl>, QC <chr>
#>
#>
#> $fits$`2:2`
#> $fits$`2:2`$xy_peaks
#> # A tibble: 17 × 5
#> x y counts_per_bin discarded from
#> <dbl> <dbl> <int> <lgl> <chr>
#> 1 0.06 0.94 89 FALSE KDE
#> 2 0.11 0.27 20 TRUE KDE
#> 3 0.22 0.84 69 FALSE KDE
#> 4 0.24 0.86 65 FALSE KDE
#> 5 0.46 8.14 625 FALSE KDE
#> 6 0.67 0.01 0 TRUE KDE
#> 7 0.72 0.02 2 TRUE KDE
#> 8 0.75 0.01 2 TRUE KDE
#> 9 0.79 0.01 1 TRUE KDE
#> 10 0.82 0.01 1 TRUE KDE
#> 11 0.85 0.01 1 TRUE KDE
#> 12 0.91 0.01 1 TRUE KDE
#> 13 0.93 0.01 0 TRUE KDE
#> 14 0.97 0.02 0 TRUE KDE
#> 15 1 0.02 3 TRUE KDE
#> 16 0.176 0.447 33 FALSE BMix
#> 17 0.461 8.14 625 FALSE BMix
#>
#> $fits$`2:2`$density
#>
#> Call:
#> density.default(x = y, adjust = kernel_adjust, kernel = "gaussian", na.rm = T)
#>
#> Data: y (7478 obs.); Bandwidth 'bw' = 0.008033
#>
#> x y
#> Min. :0.02101 Min. :0.000244
#> 1st Qu.:0.27179 1st Qu.:0.007045
#> Median :0.52256 Median :0.203133
#> Mean :0.52256 Mean :0.995946
#> 3rd Qu.:0.77333 3rd Qu.:0.787272
#> Max. :1.02410 Max. :8.138867
#>
#> $fits$`2:2`$matching
#> # A tibble: 2 × 16
#> # Rowwise:
#> mutation_multiplicity karyotype peak delta_vaf x y counts_per_bin
#> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <int>
#> 1 1 2:2 0.235 0.00700 0.24 0.86 65
#> 2 2 2:2 0.471 0.0140 0.461 8.14 625
#> # ℹ 9 more variables: discarded <lgl>, from <chr>, offset_VAF <dbl>,
#> # offset <dbl>, weight <dbl>, epsilon <dbl>, VAF_tolerance <dbl>,
#> # matched <lgl>, QC <chr>
#>
#>
#>
#> $matches
#> # A tibble: 4 × 16
#> # Rowwise:
#> mutation_multiplicity karyotype peak delta_vaf x y counts_per_bin
#> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <int>
#> 1 1 2:2 0.235 0.00700 0.24 0.86 65
#> 2 2 2:2 0.471 0.0140 0.461 8.14 625
#> 3 1 2:1 0.308 0.0120 0.307 3.60 59
#> 4 2 2:1 0.616 0.0239 0.602 2.94 60
#> # ℹ 9 more variables: discarded <lgl>, from <chr>, offset_VAF <dbl>,
#> # offset <dbl>, weight <dbl>, epsilon <dbl>, VAF_tolerance <dbl>,
#> # matched <lgl>, QC <chr>
#>
#> $matching_strategy
#> [1] "closest"
#>
#> $purity_error
#> [1] 0.05
#>
#> $min_karyotype_size
#> [1] 0
#>
#> $p_binsize_peaks
#> [1] 0.005
#>
#> $QC
#> [1] "PASS"
#>
#> $KDE
#> [1] NA
#>
#> $general
#> $general$analysis
#> [1] "3:2" "3:0" "4:2"
#>
#> $general$params
#> $general$params$n_min
#> [1] 100
#>
#> $general$params$epsilon
#> [1] 0.05
#>
#> $general$params$n_bootstrap
#> [1] 1
#>
#>
#> $general$expected_peaks
#> # A tibble: 10 × 9
#> minor Major ploidy multiplicity purity peak karyotype matched n
#> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <chr> <lgl> <table[1d]>
#> 1 2 3 5 1 0.89 0.191 3:2 TRUE 1625
#> 2 2 3 5 2 0.89 0.381 3:2 TRUE 1625
#> 3 2 3 5 3 0.89 0.572 3:2 TRUE 1625
#> 4 0 3 3 1 0.89 0.308 3:0 TRUE 312
#> 5 0 3 3 2 0.89 0.616 3:0 TRUE 312
#> 6 0 3 3 3 0.89 0.924 3:0 TRUE 312
#> 7 2 4 6 1 0.89 0.160 4:2 TRUE 1893
#> 8 2 4 6 2 0.89 0.320 4:2 TRUE 1893
#> 9 2 4 6 3 0.89 0.480 4:2 TRUE 1893
#> 10 2 4 6 4 0.89 0.640 4:2 TRUE 1893
#>
#> $general$data_peaks
#> # A tibble: 20 × 7
#> x y counts_per_bin discarded from karyotype n
#> <dbl> <dbl> <int> <lgl> <chr> <chr> <table[1d]>
#> 1 0.07 1.33 6 FALSE KDE 3:0 312
#> 2 0.17 0.39 1 TRUE KDE 3:0 312
#> 3 0.31 0.72 2 FALSE KDE 3:0 312
#> 4 0.42 0.08 0 TRUE KDE 3:0 312
#> 5 0.48 0.08 0 TRUE KDE 3:0 312
#> 6 0.63 0.42 2 TRUE KDE 3:0 312
#> 7 0.76 0.2 2 TRUE KDE 3:0 312
#> 8 0.92 10.6 34 FALSE KDE 3:0 312
#> 9 0.18 2.13 40 FALSE KDE 3:2 1625
#> 10 0.37 4.39 94 FALSE KDE 3:2 1625
#> 11 0.56 1.61 33 FALSE KDE 3:2 1625
#> 12 0.06 0.89 18 FALSE KDE 4:2 1893
#> 13 0.15 1.62 30 FALSE KDE 4:2 1893
#> 14 0.31 7.93 158 FALSE KDE 4:2 1893
#> 15 0.49 0.03 1 TRUE KDE 4:2 1893
#> 16 0.64 1.65 34 FALSE KDE 4:2 1893
#> 17 0.85 0 0 TRUE KDE 4:2 1893
#> 18 0.86 0 0 TRUE KDE 4:2 1893
#> 19 0.88 0 0 TRUE KDE 4:2 1893
#> 20 1 0.02 1 TRUE KDE 4:2 1893
#>
#> $general$data_densities
#> # A tibble: 1,536 × 4
#> x y karyotype n
#> <dbl> <dbl> <chr> <table[1d]>
#> 1 0.00120 0.00319 3:0 312
#> 2 0.00325 0.00479 3:0 312
#> 3 0.00530 0.00706 3:0 312
#> 4 0.00735 0.0102 3:0 312
#> 5 0.00940 0.0146 3:0 312
#> 6 0.0115 0.0204 3:0 312
#> 7 0.0135 0.0282 3:0 312
#> 8 0.0155 0.0383 3:0 312
#> 9 0.0176 0.0512 3:0 312
#> 10 0.0196 0.0675 3:0 312
#> # ℹ 1,526 more rows
#>
#> $general$summary
#> # A tibble: 3 × 5
#> # Groups: karyotype, matched [3]
#> karyotype n matched mismatched prop
#> <chr> <table[1d]> <int> <dbl> <dbl>
#> 1 3:0 312 3 0 1
#> 2 3:2 1625 3 0 1
#> 3 4:2 1893 4 0 1
#>
#>
# Analysis where simple clonal segments are matched with an alternative algorithm.
x = analyze_peaks(x, matching_strategy = "rightmost")
#>
#> ── Peak analysis: simple CNAs ──────────────────────────────────────────────────
#>
#> ℹ Analysing 9041 mutations mapping to karyotype(s) 2:2 and 2:1.
#> ℹ Mixed type peak detection for karyotype 2:1 (1563 mutations)
#> ℹ Mixed type peak detection for karyotype 2:2 (7478 mutations)
#> # A tibble: 4 × 16
#> # Rowwise:
#> mutation_multiplicity karyotype peak delta_vaf x y counts_per_bin
#> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <int>
#> 1 1 2:2 0.235 0.00700 0.24 0.86 65
#> 2 2 2:2 0.471 0.0140 0.461 8.14 625
#> 3 1 2:1 0.308 0.0120 0.306 3.60 59
#> 4 2 2:1 0.616 0.0239 0.600 2.94 60
#> # ℹ 9 more variables: discarded <lgl>, from <chr>, offset_VAF <dbl>,
#> # offset <dbl>, weight <dbl>, epsilon <dbl>, VAF_tolerance <dbl>,
#> # matched <lgl>, QC <chr>
#> ✔ Peak detection PASS with r = 0.00838466599059719 - maximum purity error ε = 0.05.
#> Joining with `by = join_by(Major, minor, QC_PASS)`
#> Joining with `by = join_by(karyotype, QC_PASS)`
#>
#> ── Peak analysis: complex CNAs ─────────────────────────────────────────────────
#>
#> ℹ Karyotypes 4:2, 3:2, and 3:0 with >100 mutation(s). Using epsilon = 0.05.
#> # A tibble: 3 × 5
#> # Groups: karyotype, matched [3]
#> karyotype n matched mismatched prop
#> <chr> <table[1d]> <int> <dbl> <dbl>
#> 1 3:0 312 3 0 1
#> 2 3:2 1625 3 0 1
#> 3 4:2 1893 4 0 1
#> Adding missing grouping variables: `matched`
#> Joining with `by = join_by(Major, minor, QC_PASS, matched)`
#> Adding missing grouping variables: `matched`
#> Joining with `by = join_by(karyotype, QC_PASS, matched)`
#>
#> ── Peak analysis: subclonal CNAs ───────────────────────────────────────────────
#>
#> ℹ No subclonal CNAs in this sample.
print(x)
#> ── [ CNAqc ] MySample 12963 mutations in 267 segments (267 clonal, 0 subclonal).
#>
#> ── Clonal CNAs
#>
#> 2:2 [n = 7478, L = 1483 Mb] ■■■■■■■■■■■■■■■■■■■■■■■■■■■ { CTCF }
#> 4:2 [n = 1893, L = 331 Mb] ■■■■■■■
#> 3:2 [n = 1625, L = 357 Mb] ■■■■■■
#> 2:1 [n = 1563, L = 420 Mb] ■■■■■■ { TTN }
#> 3:0 [n = 312, L = 137 Mb] ■
#> 2:0 [n = 81, L = 39 Mb] { TP53 }
#> 16:2 [n = 4, L = 0 Mb]
#> 25:2 [n = 2, L = 1 Mb]
#> 3:1 [n = 2, L = 1 Mb]
#> 106:1 [n = 1, L = 0 Mb]
#>
#> ℹ Sample Purity: 89% ~ Ploidy: 4.
#>
#> ℹ There are 3 annotated driver(s) mapped to clonal CNAs.
#> chr from to ref alt DP NV VAF driver_label is_driver
#> chr2 179431633 179431634 C T 117 77 0.6581197 TTN TRUE
#> chr16 67646006 67646007 C T 120 54 0.4500000 CTCF TRUE
#> chr17 7577106 7577107 G C 84 78 0.9285714 TP53 TRUE
#>
#> ── PASS Peaks QC closest: 200%, λ = 0.0084. Purity correction: 1%. ───────────
#>
#> ℹ 2:2 ~ n = 7478 ( 83%) → PASS -0.034 PASS 0.035
#> ℹ 2:1 ~ n = 1563 ( 17%) → PASS 0.01 PASS 0.032
#>
#> ── General peak QC (3830 mutations): PASS 10 FAIL 0 - epsilon = 0.05. ──────
#>
#> ℹ 3:0 ~ n = 312 ( 8%) → PASS 3 FAIL 0
#> ℹ 3:2 ~ n = 1625 ( 42%) → PASS 3 FAIL 0
#> ℹ 4:2 ~ n = 1893 ( 49%) → PASS 4 FAIL 0