Skip to main content

sklears_datasets/validation/
distributions.rs

1//! Distribution testing and validation functions
2//!
3//! This module provides functions for testing whether data follows
4//! specific statistical distributions and performing goodness-of-fit tests.
5
6use super::types::{ValidationConfig, ValidationReport, ValidationResult};
7
8/// Distribution types for testing
9#[derive(Debug, Clone)]
10pub enum DistributionType {
11    /// Normal
12    Normal(f64, f64), // mean, std
13    /// Uniform
14    Uniform(f64, f64), // min, max
15    /// Exponential
16    Exponential(f64), // rate
17    /// Custom
18    Custom(Vec<f64>), // reference samples
19}
20
21/// Kolmogorov-Smirnov test for distribution comparison
22pub fn kolmogorov_smirnov_test(
23    data: &[f64],
24
25    expected_dist: &DistributionType,
26    config: &ValidationConfig,
27) -> ValidationResult {
28    let mut sorted_data = data.to_vec();
29    sorted_data.sort_by(|a, b| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal));
30
31    let n = sorted_data.len() as f64;
32    let mut max_diff: f64 = 0.0;
33
34    for (i, &value) in sorted_data.iter().enumerate() {
35        let empirical_cdf = (i + 1) as f64 / n;
36        let theoretical_cdf = match expected_dist {
37            DistributionType::Normal(mean, std) => {
38                // Approximate normal CDF using error function approximation
39                let z = (value - mean) / std;
40                0.5 * (1.0
41                    + z / (1.0
42                        + 0.278393 * z.abs()
43                        + 0.230389 * z.abs().powi(2)
44                        + 0.000972 * z.abs().powi(3)
45                        + 0.078108 * z.abs().powi(4))
46                    .powi(4))
47            }
48            DistributionType::Uniform(min, max) => {
49                if value < *min {
50                    0.0
51                } else if value > *max {
52                    1.0
53                } else {
54                    (value - min) / (max - min)
55                }
56            }
57            DistributionType::Exponential(rate) => {
58                if value < 0.0 {
59                    0.0
60                } else {
61                    1.0 - (-rate * value).exp()
62                }
63            }
64            DistributionType::Custom(ref samples) => {
65                let pos = samples.iter().filter(|&&x| x <= value).count();
66                pos as f64 / samples.len() as f64
67            }
68        };
69
70        let diff = (empirical_cdf - theoretical_cdf).abs();
71        max_diff = max_diff.max(diff);
72    }
73
74    // Critical value for α = 0.05
75    let critical_value = 1.36 / (n.sqrt());
76    let passed = max_diff < critical_value;
77
78    ValidationResult {
79        property: "Kolmogorov-Smirnov Test".to_string(),
80        passed,
81        expected: critical_value,
82        actual: max_diff,
83        tolerance: config.tolerance,
84        message: if passed {
85            "Distribution matches expected distribution".to_string()
86        } else {
87            "Distribution does not match expected distribution".to_string()
88        },
89    }
90}
91
92/// Chi-square goodness of fit test
93pub fn chi_square_goodness_of_fit_test(
94    data: &[f64],
95    expected_dist: &DistributionType,
96    num_bins: usize,
97    config: &ValidationConfig,
98) -> ValidationResult {
99    let min_val = data.iter().cloned().fold(f64::INFINITY, f64::min);
100    let max_val = data.iter().cloned().fold(f64::NEG_INFINITY, f64::max);
101    let bin_width = (max_val - min_val) / num_bins as f64;
102
103    let mut observed_counts = vec![0; num_bins];
104    let mut expected_counts = vec![0.0; num_bins];
105
106    // Count observed frequencies
107    for &value in data {
108        let bin_idx = if value == max_val {
109            num_bins - 1
110        } else {
111            ((value - min_val) / bin_width).floor() as usize
112        };
113        observed_counts[bin_idx] += 1;
114    }
115
116    // Calculate expected frequencies
117    // i is used both for arithmetic (bin_start/bin_end) and for indexing
118    #[allow(clippy::needless_range_loop)]
119    for i in 0..num_bins {
120        let bin_start = min_val + i as f64 * bin_width;
121        let bin_end = min_val + (i + 1) as f64 * bin_width;
122
123        let prob = match expected_dist {
124            DistributionType::Normal(mean, std) => {
125                // Simple approximation for normal distribution probability in bin
126                let z1 = (bin_start - mean) / std;
127                let z2 = (bin_end - mean) / std;
128                let cdf1 =
129                    0.5 * (1.0 + z1 / (1.0 + 0.278393 * z1.abs() + 0.230389 * z1.abs().powi(2)));
130                let cdf2 =
131                    0.5 * (1.0 + z2 / (1.0 + 0.278393 * z2.abs() + 0.230389 * z2.abs().powi(2)));
132                cdf2 - cdf1
133            }
134            DistributionType::Uniform(min, max) => {
135                let overlap_start = bin_start.max(*min);
136                let overlap_end = bin_end.min(*max);
137                if overlap_start < overlap_end {
138                    (overlap_end - overlap_start) / (max - min)
139                } else {
140                    0.0
141                }
142            }
143            DistributionType::Exponential(rate) => {
144                let cdf1 = if bin_start < 0.0 {
145                    0.0
146                } else {
147                    1.0 - (-rate * bin_start).exp()
148                };
149                let cdf2 = if bin_end < 0.0 {
150                    0.0
151                } else {
152                    1.0 - (-rate * bin_end).exp()
153                };
154                cdf2 - cdf1
155            }
156            DistributionType::Custom(_) => {
157                // For custom distributions, use uniform as fallback
158                1.0 / num_bins as f64
159            }
160        };
161
162        expected_counts[i] = prob * data.len() as f64;
163    }
164
165    // Calculate chi-square statistic
166    let mut chi_square = 0.0;
167    for i in 0..num_bins {
168        if expected_counts[i] > 0.0 {
169            let diff = observed_counts[i] as f64 - expected_counts[i];
170            chi_square += diff * diff / expected_counts[i];
171        }
172    }
173
174    // Critical value for α = 0.05, df = num_bins - 1
175    let degrees_of_freedom = num_bins - 1;
176    let critical_value = match degrees_of_freedom {
177        1 => 3.841,
178        2 => 5.991,
179        3 => 7.815,
180        4 => 9.488,
181        5 => 11.070,
182        _ => 3.841 + degrees_of_freedom as f64 * 2.0, // Rough approximation
183    };
184
185    let passed = chi_square < critical_value;
186
187    ValidationResult {
188        property: "Chi-Square Goodness of Fit".to_string(),
189        passed,
190        expected: critical_value,
191        actual: chi_square,
192        tolerance: config.tolerance,
193        message: if passed {
194            "Distribution passes chi-square goodness of fit test".to_string()
195        } else {
196            "Distribution fails chi-square goodness of fit test".to_string()
197        },
198    }
199}
200
201/// Test for uniform distribution
202pub fn validate_uniform_distribution(
203    data: &[f64],
204    expected_min: f64,
205    expected_max: f64,
206    config: &ValidationConfig,
207) -> ValidationReport {
208    let mut report = ValidationReport::new();
209
210    let min_val = data.iter().cloned().fold(f64::INFINITY, f64::min);
211    let max_val = data.iter().cloned().fold(f64::NEG_INFINITY, f64::max);
212
213    // Check if data range is within expected bounds
214    let min_passed = (min_val - expected_min).abs() <= config.tolerance;
215    let max_passed = (max_val - expected_max).abs() <= config.tolerance;
216
217    report.add_result(ValidationResult {
218        property: "Uniform Distribution Min".to_string(),
219        passed: min_passed,
220        expected: expected_min,
221        actual: min_val,
222        tolerance: config.tolerance,
223        message: if min_passed {
224            "Minimum value matches expected".to_string()
225        } else {
226            "Minimum value does not match expected".to_string()
227        },
228    });
229
230    report.add_result(ValidationResult {
231        property: "Uniform Distribution Max".to_string(),
232        passed: max_passed,
233        expected: expected_max,
234        actual: max_val,
235        tolerance: config.tolerance,
236        message: if max_passed {
237            "Maximum value matches expected".to_string()
238        } else {
239            "Maximum value does not match expected".to_string()
240        },
241    });
242
243    // Use Kolmogorov-Smirnov test
244    let ks_result = kolmogorov_smirnov_test(
245        data,
246        &DistributionType::Uniform(expected_min, expected_max),
247        config,
248    );
249    report.add_result(ks_result);
250
251    report
252}
253
254/// Test for normal distribution with enhanced checks
255pub fn validate_normal_distribution(
256    data: &[f64],
257    expected_mean: f64,
258    expected_std: f64,
259    config: &ValidationConfig,
260) -> ValidationReport {
261    let mut report = ValidationReport::new();
262
263    let n = data.len() as f64;
264    let mean = data.iter().sum::<f64>() / n;
265    let variance = data.iter().map(|&x| (x - mean).powi(2)).sum::<f64>() / n;
266    let std_dev = variance.sqrt();
267
268    // Check mean and standard deviation
269    let mean_passed = (mean - expected_mean).abs() <= config.tolerance;
270    let std_passed = (std_dev - expected_std).abs() <= config.tolerance;
271
272    report.add_result(ValidationResult {
273        property: "Normal Distribution Mean".to_string(),
274        passed: mean_passed,
275        expected: expected_mean,
276        actual: mean,
277        tolerance: config.tolerance,
278        message: if mean_passed {
279            "Mean matches expected value".to_string()
280        } else {
281            "Mean does not match expected value".to_string()
282        },
283    });
284
285    report.add_result(ValidationResult {
286        property: "Normal Distribution Std Dev".to_string(),
287        passed: std_passed,
288        expected: expected_std,
289        actual: std_dev,
290        tolerance: config.tolerance,
291        message: if std_passed {
292            "Standard deviation matches expected value".to_string()
293        } else {
294            "Standard deviation does not match expected value".to_string()
295        },
296    });
297
298    // Use Kolmogorov-Smirnov test
299    let ks_result = kolmogorov_smirnov_test(
300        data,
301        &DistributionType::Normal(expected_mean, expected_std),
302        config,
303    );
304    report.add_result(ks_result);
305
306    // Use Chi-square test
307    let chi_square_result = chi_square_goodness_of_fit_test(
308        data,
309        &DistributionType::Normal(expected_mean, expected_std),
310        10,
311        config,
312    );
313    report.add_result(chi_square_result);
314
315    report
316}
317
318/// Test for exponential distribution
319pub fn validate_exponential_distribution(
320    data: &[f64],
321    expected_rate: f64,
322    config: &ValidationConfig,
323) -> ValidationReport {
324    let mut report = ValidationReport::new();
325
326    let n = data.len() as f64;
327    let mean = data.iter().sum::<f64>() / n;
328    let expected_mean = 1.0 / expected_rate;
329
330    // Check mean (should be 1/rate for exponential distribution)
331    let mean_passed = (mean - expected_mean).abs() <= config.tolerance;
332
333    report.add_result(ValidationResult {
334        property: "Exponential Distribution Mean".to_string(),
335        passed: mean_passed,
336        expected: expected_mean,
337        actual: mean,
338        tolerance: config.tolerance,
339        message: if mean_passed {
340            "Mean matches expected value for exponential distribution".to_string()
341        } else {
342            "Mean does not match expected value for exponential distribution".to_string()
343        },
344    });
345
346    // Use Kolmogorov-Smirnov test
347    let ks_result =
348        kolmogorov_smirnov_test(data, &DistributionType::Exponential(expected_rate), config);
349    report.add_result(ks_result);
350
351    report
352}