riccardo ha revisionato questo gist 9 months ago. Vai alla revisione
1 file changed, 2 insertions, 1 deletion
mp_distr.py
| @@ -10,7 +10,8 @@ class MarchenkoPastur: | |||
| 10 | 10 | Parameters | |
| 11 | 11 | ---------- | |
| 12 | 12 | ratio : float | |
| 13 | - | The ratio between the number of variables (columns) and the size of the sample (rows) contained in the data matrix. For numerical stability, it should be less than 1. | |
| 13 | + | The ratio between the number of variables (columns) and the size of the sample (rows) contained in the data matrix. | |
| 14 | + | For numerical stability, it should be less than 1. | |
| 14 | 15 | sigma : float | |
| 15 | 16 | The standard deviation of the distribution of values, by default, 1.0. | |
| 16 | 17 | ||
riccardo ha revisionato questo gist 9 months ago. Vai alla revisione
Nessuna modifica
riccardo ha revisionato questo gist 11 months ago. Vai alla revisione
1 file changed, 2 insertions, 2 deletions
mp_distr.py
| @@ -27,8 +27,8 @@ class MarchenkoPastur: | |||
| 27 | 27 | self.sigma = sigma | |
| 28 | 28 | ||
| 29 | 29 | # Compute the limits of the distribution | |
| 30 | - | self.l_bottom = sigma**2 * (1.0 - np.sqrt(self.ratio))**2 | |
| 31 | - | self.l_upper = sigma**2 * (1.0 + np.sqrt(self.ratio))**2 | |
| 30 | + | self.l_bottom = self.sigma**2 * (1.0 - np.sqrt(self.ratio))**2 | |
| 31 | + | self.l_upper = self.sigma**2 * (1.0 + np.sqrt(self.ratio))**2 | |
| 32 | 32 | ||
| 33 | 33 | def pdf(self, x: float | ArrayLike) -> float | ArrayLike: | |
| 34 | 34 | """ | |
riccardo ha revisionato questo gist 11 months ago. Vai alla revisione
1 file changed, 4 insertions, 1 deletion
mp_distr.py
| @@ -47,7 +47,10 @@ class MarchenkoPastur: | |||
| 47 | 47 | if not np.isscalar(x): | |
| 48 | 48 | return np.vectorize(self.pdf, otypes=[float])(x) | |
| 49 | 49 | ||
| 50 | + | if x == 0.0: | |
| 51 | + | return 0.0 | |
| 52 | + | ||
| 50 | 53 | num = np.sqrt(max(self.l_upper - x, 0.0) * max(x - self.l_bottom, 0.0)) | |
| 51 | 54 | den = 2.0 * np.pi * self.sigma**2 * self.ratio * x | |
| 52 | 55 | ||
| 53 | - | return float(num / den) if den != 0.0 else 0.0 | |
| 56 | + | return float(num / den) | |
riccardo ha revisionato questo gist 11 months ago. Vai alla revisione
1 file changed, 53 insertions
mp_distr.py(file creato)
| @@ -0,0 +1,53 @@ | |||
| 1 | + | import numpy as np | |
| 2 | + | ||
| 3 | + | from numpy.typing import ArrayLike | |
| 4 | + | ||
| 5 | + | class MarchenkoPastur: | |
| 6 | + | """Definition of a Marchenko-Pastur distribution""" | |
| 7 | + | ||
| 8 | + | def __init__(self, ratio: float, sigma: float = 1.0): | |
| 9 | + | """ | |
| 10 | + | Parameters | |
| 11 | + | ---------- | |
| 12 | + | ratio : float | |
| 13 | + | The ratio between the number of variables (columns) and the size of the sample (rows) contained in the data matrix. For numerical stability, it should be less than 1. | |
| 14 | + | sigma : float | |
| 15 | + | The standard deviation of the distribution of values, by default, 1.0. | |
| 16 | + | ||
| 17 | + | Raises | |
| 18 | + | ------ | |
| 19 | + | ValueError | |
| 20 | + | If ratio or sigma are not strictly positive. | |
| 21 | + | """ | |
| 22 | + | if ratio <= 0.0: | |
| 23 | + | raise ValueError("The ratio must be strictly positive, but found %s <= 0.0!" % ratio) | |
| 24 | + | self.ratio = ratio | |
| 25 | + | if sigma <= 0.0: | |
| 26 | + | raise ValueError("The standard deviation must be strictly positive, but found %s <= 0.0!" % sigma) | |
| 27 | + | self.sigma = sigma | |
| 28 | + | ||
| 29 | + | # Compute the limits of the distribution | |
| 30 | + | self.l_bottom = sigma**2 * (1.0 - np.sqrt(self.ratio))**2 | |
| 31 | + | self.l_upper = sigma**2 * (1.0 + np.sqrt(self.ratio))**2 | |
| 32 | + | ||
| 33 | + | def pdf(self, x: float | ArrayLike) -> float | ArrayLike: | |
| 34 | + | """ | |
| 35 | + | Return the value of the probability distribution function. | |
| 36 | + | ||
| 37 | + | Parameters | |
| 38 | + | ---------- | |
| 39 | + | x : float | ArrayLike | |
| 40 | + | The value(s) at which to compute the PDF. | |
| 41 | + | ||
| 42 | + | Returns | |
| 43 | + | ------- | |
| 44 | + | float | ArrayLike | |
| 45 | + | The value(s) of the PDF | |
| 46 | + | """ | |
| 47 | + | if not np.isscalar(x): | |
| 48 | + | return np.vectorize(self.pdf, otypes=[float])(x) | |
| 49 | + | ||
| 50 | + | num = np.sqrt(max(self.l_upper - x, 0.0) * max(x - self.l_bottom, 0.0)) | |
| 51 | + | den = 2.0 * np.pi * self.sigma**2 * self.ratio * x | |
| 52 | + | ||
| 53 | + | return float(num / den) if den != 0.0 else 0.0 | |