Skip to main content

causal_hub/random/datasets/table/gaussian/
incomplete.rs

1use ndarray::prelude::*;
2use ndarray_stats::{Quantile1dExt, interpolate::Nearest};
3use noisy_float::prelude::*;
4use rand::{Rng, RngExt};
5use rand_distr::{Uniform, num_traits::ToPrimitive};
6
7use crate::{
8    datasets::{Dataset, GaussIncTable, GaussTable, GaussType, IncDataset, MissingMechanism},
9    models::HasLabels,
10    random::Random,
11    types::{Error, Result},
12};
13
14/// A struct representing a random incomplete gaussian table dataset generator.
15pub struct RngGaussIncTable<'a, R> {
16    rng: &'a mut R,
17    dataset: &'a GaussTable,
18    missing_mechanism: &'a MissingMechanism,
19    p_min: f64,
20    p_max: f64,
21}
22
23impl<'a, R: Rng> RngGaussIncTable<'a, R> {
24    /// Creates a new `RngGaussIncTable` instance.
25    ///
26    /// # Arguments
27    ///
28    /// * `rng` - A mutable reference to a random number generator.
29    /// * `dataset` - A reference to the complete gaussian table dataset.
30    /// * `missing_mechanism` - A reference to the missingness mechanism.
31    /// * `p_min` - The minimum probability of missingness.
32    /// * `p_max` - The maximum probability of missingness.
33    ///
34    /// # Returns
35    ///
36    /// A new `RngGaussIncTable` instance.
37    pub fn new(
38        rng: &'a mut R,
39        dataset: &'a GaussTable,
40        missing_mechanism: &'a MissingMechanism,
41        p_min: f64,
42        p_max: f64,
43    ) -> Result<Self> {
44        // Check that dataset labels are equals to missing mechanism labels.
45        if dataset.labels() != missing_mechanism.labels() {
46            return Err(Error::InvalidParameter(
47                "missing_mechanism",
48                "labels do not match dataset labels",
49            ));
50        }
51        // Check that p_min and p_max are in [0, 1].
52        if !(0.0..=1.0).contains(&p_min) {
53            return Err(Error::InvalidParameter("p_min", "must be in [0, 1]"));
54        }
55        if !(0.0..=1.0).contains(&p_max) {
56            return Err(Error::InvalidParameter("p_max", "must be in [0, 1]"));
57        }
58        // Check that p_min is less than or equal to p_max.
59        if p_min > p_max {
60            return Err(Error::InvalidParameter(
61                "p_min",
62                "must be less than or equal to p_max",
63            ));
64        }
65
66        Ok(Self {
67            rng,
68            dataset,
69            missing_mechanism,
70            p_min,
71            p_max,
72        })
73    }
74}
75
76impl<R: Rng> Random for RngGaussIncTable<'_, R> {
77    type Output = Result<GaussIncTable>;
78
79    fn random(&mut self) -> Self::Output {
80        // Get the missing indicator.
81        const M: GaussType = GaussIncTable::MISSING;
82        // Get dataset labels.
83        let labels = self.dataset.labels().clone();
84        // Get dataset values.
85        let mut values = self.dataset.values().clone();
86
87        // Define the uniform distributions for sampling.
88        let p_s = Uniform::new_inclusive(0., 1.)?;
89        let p_u = Uniform::new_inclusive(self.p_min, self.p_max)?;
90
91        // Iterate over the missing mechanism.
92        for (&x, pa_x) in self.missing_mechanism {
93            // Get mutable reference to the column X.
94            let mut c_x = values.column_mut(x);
95
96            // Check if the variable has no parents.
97            if pa_x.is_empty() {
98                // Sample a missingness probability.
99                let p_x = self.rng.sample(p_u);
100                // Modify the corresponding column.
101                c_x.iter_mut().for_each(|x| {
102                    if self.rng.sample(p_s) < p_x {
103                        *x = M;
104                    }
105                });
106                // Continue to the next variable.
107                continue;
108            }
109
110            // For each parent ...
111            for &z in pa_x {
112                // Get reference to the column Z.
113                let c_z = self.dataset.values().column(z);
114                // Get the threshold of the parent via quantile.
115                let s_z = match c_z
116                    .mapv(n64)
117                    .quantile_mut(n64(0.2), &Nearest)
118                    .map(|q| q.to_f64())
119                {
120                    Ok(Some(q)) => q,
121                    // If quantile computation fails, skip this parent.
122                    _ => continue,
123                };
124                // Modify the corresponding column.
125                azip!((x in &mut c_x, z in &c_z) {
126                    // Sample a missingness probability for X given parent Z.
127                    if
128                        (self.rng.sample(p_s) < self.p_max && *z <  s_z) ||
129                        (self.rng.sample(p_s) < self.p_min && *z >= s_z)
130                    {
131                        *x = M;
132                    }
133                });
134            }
135        }
136
137        // Return the incomplete dataset.
138        GaussIncTable::new(labels, values)
139    }
140}