Score diagnostics for pepVet scoring models

score_diagnostics() runs three analyses on a batch_evaluate() result to quantify multicollinearity, dimensionality, and component contributions in the scoring model. Use it to validate that the five (or six) component scores are measuring orthogonal dimensions and that each contributes meaningfully to the composite.

Usage

score_diagnostics(batch_result, weights = NULL)

Arguments

batch_result: A tibble returned by batch_evaluate(). Must contain at least two rows and two component-score columns (prefixed S_). If NULL, too few rows, or missing required columns, raises an error.
weights: Optional named numeric vector of component weights. When NULL (default), weights are inferred from the detected component columns using pepVet's default scoring weights (protein-only or proteome-aware, depending on presence of S_unique).

Value

A named list with six elements:

vif: A named numeric vector of VIF values, one per component. NA when too few proteins for reliable estimation and Inf for perfect collinearity.
pca: A list with var_explained (numeric vector), loadings (rotation matrix), sdev (standard deviations), and x (PCA scores matrix).
ablation: A data.frame with columns component, weight, mean_drop, sd_drop, max_drop, and n_verdict_flipped.
n_proteins: Number of proteins in the input.
n_components: Number of component scores detected.
weights: The resolved weight vector used for ablation.

Details

Three analyses are performed:

VIF (Variance Inflation Factor): For each component, a linear regression is fit using all other components as predictors. VIF = 1 / (1 - R^2). VIF < 5 is acceptable. VIF > 10 indicates problematic collinearity that may warrant component removal or merger. VIF requires more than n_components + 1 observations; when this condition is not met for an input with at least two rows, NA values are returned with a warning. Perfect collinearity is represented by Inf rather than an unclassed regression warning.

PCA (Principal Component Analysis): Centered and scaled PCA on the component-score matrix. The proportion of variance explained by each principal component reveals the effective dimensionality of the scoring model. If PC1 + PC2 explain > 80% of variance, the model has low effective dimensionality (expected for a well-designed multi-criteria score).

Ablation: Each component is set to 0 (its minimum possible value) for all proteins while others remain at their actual values, and the composite score is recomputed. The mean and maximum drop in composite score, plus the number of verdicts that flip, quantify each component's marginal contribution. Components with mean drops near zero are candidates for removal or down-weighting.

Limitations

Diagnostics are meaningful only on batch results with enough proteins (at least 10-20 recommended; VIF requires more than n_components + 1). Results are specific to the enzyme, missed-cleavage setting, and protein set used. Running diagnostics on a different batch may give different VIF and ablation values.

Examples

small_fasta <- system.file(
  "extdata", "small_proteome_50_proteins.fasta",
  package = "pepVet"
)
batch <- batch_evaluate(small_fasta, enzyme = "trypsin")
diag <- score_diagnostics(batch)
diag$vif
#>   S_length S_coverage    S_count    S_hydro   S_charge 
#>   5.776149   2.972499   4.640217   1.688987   2.203141 
diag$ablation
#>    component weight  mean_drop    sd_drop   max_drop n_verdict_flipped
#> 1   S_length  0.200 0.09700884 0.01976096 0.12432432                26
#> 2 S_coverage  0.348 0.23973064 0.06094034 0.34124272                45
#> 3    S_count  0.226 0.20766686 0.04160020 0.22600000                43
#> 4    S_hydro  0.138 0.08726940 0.01952765 0.12075000                20
#> 5   S_charge  0.088 0.06177470 0.01058011 0.07372973                14

# Proteome-aware mode
prot_digest <- digest_protein(small_fasta, enzyme = "trypsin")
batch_pa <- batch_evaluate(small_fasta, enzyme = "trypsin",
  proteome = prot_digest)
diag_pa <- score_diagnostics(batch_pa)
diag_pa$vif
#>   S_length S_coverage    S_count    S_hydro   S_charge   S_unique 
#>   5.870346   3.106780   4.675022   2.100193   3.555974   3.907455