rv/dist/
students_t.rs

1#[cfg(feature = "serde1")]
2use serde::{Deserialize, Serialize};
3
4use crate::impl_display;
5use crate::misc::ln_gammafn;
6use crate::traits::*;
7use rand::Rng;
8use std::f64::consts::PI;
9use std::fmt;
10
11/// [Student's T distribution](https://en.wikipedia.org/wiki/Student%27s_t-distribution)
12/// over x in (-∞, ∞).
13#[derive(Debug, Clone, PartialEq)]
14#[cfg_attr(feature = "serde1", derive(Serialize, Deserialize))]
15#[cfg_attr(feature = "serde1", serde(rename_all = "snake_case"))]
16pub struct StudentsT {
17    /// Degrees of freedom, ν, in (0, ∞)
18    v: f64,
19}
20
21impl Parameterized for StudentsT {
22    type Parameters = f64;
23
24    fn emit_params(&self) -> Self::Parameters {
25        self.v()
26    }
27
28    fn from_params(v: Self::Parameters) -> Self {
29        Self::new_unchecked(v)
30    }
31}
32
33#[derive(Debug, Clone, PartialEq)]
34#[cfg_attr(feature = "serde1", derive(Serialize, Deserialize))]
35#[cfg_attr(feature = "serde1", serde(rename_all = "snake_case"))]
36pub enum StudentsTError {
37    /// The v parameter is infinite or NaN
38    VNotFinite { v: f64 },
39    /// The v parameter is less than or equal to zero
40    VTooLow { v: f64 },
41}
42
43impl StudentsT {
44    /// Create a new Student's T distribution with degrees of freedom, v.
45    #[inline]
46    pub fn new(v: f64) -> Result<Self, StudentsTError> {
47        if v <= 0.0 {
48            Err(StudentsTError::VTooLow { v })
49        } else if !v.is_finite() {
50            Err(StudentsTError::VNotFinite { v })
51        } else {
52            Ok(StudentsT { v })
53        }
54    }
55
56    /// Creates a new StudentsT without checking whether the parameter is
57    /// valid.
58    #[inline]
59    pub fn new_unchecked(v: f64) -> Self {
60        StudentsT { v }
61    }
62
63    /// Get the degrees of freedom, v
64    #[inline]
65    pub fn v(&self) -> f64 {
66        self.v
67    }
68
69    /// Set the value of v
70    ///
71    /// # Example
72    ///
73    /// ```rust
74    /// use rv::dist::StudentsT;
75    ///
76    /// let mut t = StudentsT::new(1.2).unwrap();
77    /// assert_eq!(t.v(), 1.2);
78    ///
79    /// t.set_v(4.3).unwrap();
80    /// assert_eq!(t.v(), 4.3);
81    /// ```
82    ///
83    /// Will error for invalid values
84    ///
85    /// ```rust
86    /// # use rv::dist::StudentsT;
87    /// # let mut t = StudentsT::new(1.2).unwrap();
88    /// assert!(t.set_v(2.1).is_ok());
89    ///
90    /// // must be greater than zero
91    /// assert!(t.set_v(0.0).is_err());
92    /// assert!(t.set_v(-1.0).is_err());
93    ///
94    ///
95    /// assert!(t.set_v(f64::INFINITY).is_err());
96    /// assert!(t.set_v(f64::NEG_INFINITY).is_err());
97    /// assert!(t.set_v(f64::NAN).is_err());
98    /// ```
99    #[inline]
100    pub fn set_v(&mut self, v: f64) -> Result<(), StudentsTError> {
101        if !v.is_finite() {
102            Err(StudentsTError::VNotFinite { v })
103        } else if v <= 0.0 {
104            Err(StudentsTError::VTooLow { v })
105        } else {
106            self.set_v_unchecked(v);
107            Ok(())
108        }
109    }
110
111    /// Set the value of v without input validation
112    #[inline]
113    pub fn set_v_unchecked(&mut self, v: f64) {
114        self.v = v;
115    }
116}
117
118impl Default for StudentsT {
119    fn default() -> Self {
120        StudentsT { v: 2.0 }
121    }
122}
123
124impl From<&StudentsT> for String {
125    fn from(t: &StudentsT) -> String {
126        format!("Student's({})", t.v)
127    }
128}
129
130impl_display!(StudentsT);
131
132macro_rules! impl_traits {
133    ($kind:ty) => {
134        impl HasDensity<$kind> for StudentsT {
135            fn ln_f(&self, x: &$kind) -> f64 {
136                // TODO: could cache ln(pi*v) and ln_gamma(v/2)
137                let vp1 = (self.v + 1.0) / 2.0;
138                let xf = f64::from(*x);
139                let xterm = -vp1 * (xf * xf / self.v).ln_1p();
140                let zterm = 0.5_f64.mul_add(
141                    -(self.v * PI).ln(),
142                    ln_gammafn(vp1) - ln_gammafn(self.v / 2.0),
143                );
144                zterm + xterm
145            }
146        }
147
148        impl Sampleable<$kind> for StudentsT {
149            fn draw<R: Rng>(&self, rng: &mut R) -> $kind {
150                let t = rand_distr::StudentT::new(self.v).unwrap();
151                rng.sample(t) as $kind
152            }
153
154            fn sample<R: Rng>(&self, n: usize, rng: &mut R) -> Vec<$kind> {
155                let t = rand_distr::StudentT::new(self.v).unwrap();
156                (0..n).map(|_| rng.sample(t) as $kind).collect()
157            }
158        }
159
160        impl Support<$kind> for StudentsT {
161            fn supports(&self, x: &$kind) -> bool {
162                x.is_finite()
163            }
164        }
165
166        impl ContinuousDistr<$kind> for StudentsT {}
167
168        impl Mean<$kind> for StudentsT {
169            fn mean(&self) -> Option<$kind> {
170                if self.v > 1.0 {
171                    Some(0.0)
172                } else {
173                    None
174                }
175            }
176        }
177
178        impl Median<$kind> for StudentsT {
179            fn median(&self) -> Option<$kind> {
180                Some(0.0)
181            }
182        }
183
184        impl Mode<$kind> for StudentsT {
185            fn mode(&self) -> Option<$kind> {
186                Some(0.0)
187            }
188        }
189
190        impl Variance<$kind> for StudentsT {
191            fn variance(&self) -> Option<$kind> {
192                if self.v > 2.0 {
193                    Some((self.v / (self.v - 2.0)) as $kind)
194                } else {
195                    None
196                }
197            }
198        }
199    };
200}
201
202impl Skewness for StudentsT {
203    fn skewness(&self) -> Option<f64> {
204        if self.v > 3.0 {
205            Some(0.0)
206        } else {
207            None
208        }
209    }
210}
211
212impl Kurtosis for StudentsT {
213    fn kurtosis(&self) -> Option<f64> {
214        if self.v > 4.0 {
215            Some(6.0 / (self.v - 4.0))
216        } else if self.v > 2.0 {
217            Some(f64::INFINITY)
218        } else {
219            None
220        }
221    }
222}
223
224impl_traits!(f64);
225impl_traits!(f32);
226
227impl std::error::Error for StudentsTError {}
228
229impl fmt::Display for StudentsTError {
230    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
231        match self {
232            Self::VNotFinite { v } => write!(f, "non-finite v: {}", v),
233            Self::VTooLow { v } => {
234                write!(f, "v ({}) must be greater than zero", v)
235            }
236        }
237    }
238}
239
240#[cfg(test)]
241mod tests {
242    use super::*;
243    use crate::test_basic_impls;
244    use std::f64;
245
246    const TOL: f64 = 1E-12;
247
248    test_basic_impls!(f64, StudentsT);
249
250    #[test]
251    fn new() {
252        let t = StudentsT::new(2.3).unwrap();
253        assert::close(t.v, 2.3, TOL);
254    }
255
256    #[test]
257    fn new_should_reject_v_leq_zero() {
258        assert!(StudentsT::new(f64::MIN_POSITIVE).is_ok());
259        assert!(StudentsT::new(0.0).is_err());
260        assert!(StudentsT::new(-f64::MIN_POSITIVE).is_err());
261        assert!(StudentsT::new(-1.0).is_err());
262    }
263
264    #[test]
265    fn new_should_reject_non_finite_v() {
266        assert!(StudentsT::new(f64::INFINITY).is_err());
267        assert!(StudentsT::new(-f64::NAN).is_err());
268        assert!(StudentsT::new(f64::NEG_INFINITY).is_err());
269    }
270
271    #[test]
272    fn ln_pdf() {
273        let t = StudentsT::new(2.3).unwrap();
274        assert::close(t.ln_pdf(&0.0_f64), -1.024_744_023_893_756_6, TOL);
275        assert::close(t.ln_pdf(&1.0_f64), -1.620_416_044_030_352, TOL);
276        assert::close(t.ln_pdf(&2.5_f64), -3.191_230_587_916_138, TOL);
277        assert::close(t.ln_pdf(&-2.5_f64), -3.191_230_587_916_138, TOL);
278    }
279
280    #[test]
281    fn variance() {
282        let v: f64 = StudentsT::new(2.3).unwrap().variance().unwrap();
283        assert::close(v, 7.666_666_666_666_670_5, TOL);
284    }
285
286    #[test]
287    fn median() {
288        let m: f64 = StudentsT::new(2.3).unwrap().median().unwrap();
289        assert::close(m, 0.0, TOL);
290    }
291
292    #[test]
293    fn mode() {
294        let m: f64 = StudentsT::new(2.3).unwrap().mode().unwrap();
295        assert::close(m, 0.0, TOL);
296    }
297}
rv/dist/students_t.rs

rv/dist/
students_t.rs
