From 87bd31afcac4664e971480b3b8ef6a398f7b7219 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Sat, 10 May 2025 15:38:30 +0300 Subject: [PATCH 01/24] feat(generator): unify mixture generators and support canonical forms This commit unifies the sample generators for all mixture types (NMM, NMV, NV) by integrating canonical form support directly into the main generator function. Previously, each mixture type had separate classical and canonical generators. Now, the generator checks mixture_form inside the function and switches behavior accordingly. In canonical form, the beta parameter is not used. This simplifies the generator API and reduces code duplication. - Updated nm_generator.py, nmv_generator.py, nv_generator.py - Remembered `mixture_form` in AbstractMixture class - Updated tests and notebook to use the new unified API BREAKING CHANGE: canonical_generate() methods were removed; use generate() instead. --- .../nm_sigma_estimation_comparison.ipynb | 2 +- src/generators/nm_generator.py | 25 +++-------------- src/generators/nmv_generator.py | 27 +++---------------- src/generators/nv_generator.py | 27 +++---------------- src/mixtures/abstract_mixture.py | 1 + .../nm_generator/test_mixing_normal.py | 20 +++++++------- 6 files changed, 23 insertions(+), 79 deletions(-) diff --git a/jupiter_examples/nm_sigma_estimation_comparison.ipynb b/jupiter_examples/nm_sigma_estimation_comparison.ipynb index f986d13..10f695f 100644 --- a/jupiter_examples/nm_sigma_estimation_comparison.ipynb +++ b/jupiter_examples/nm_sigma_estimation_comparison.ipynb @@ -250,7 +250,7 @@ " \"\"\"\n", " generator = NMGenerator()\n", " mixture = NormalMeanMixtures(\"canonical\", sigma=real_sigma, distribution=distribution)\n", - " return generator.canonical_generate(mixture, sample_len)\n", + " return generator.generate(mixture, sample_len)\n", "\n", "def estimate_sigma_eigenvalue_based(sample, real_sigma, search_area, a, b):\n", " sample_len = len(sample)\n", diff --git a/src/generators/nm_generator.py b/src/generators/nm_generator.py index 4c95b0f..5ecfe4a 100644 --- a/src/generators/nm_generator.py +++ b/src/generators/nm_generator.py @@ -9,7 +9,7 @@ class NMGenerator(AbstractGenerator): @staticmethod - def classical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: + def generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: """Generate a sample of given size. Classical form of NMM Args: @@ -27,25 +27,6 @@ def classical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: raise ValueError("Mixture must be NormalMeanMixtures") mixing_values = mixture.params.distribution.rvs(size=size) normal_values = scipy.stats.norm.rvs(size=size) + if mixture.mixture_form == "canonical": + return mixing_values + mixture.params.sigma * normal_values return mixture.params.alpha + mixture.params.beta * mixing_values + mixture.params.gamma * normal_values - - @staticmethod - def canonical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: - """Generate a sample of given size. Canonical form of NMM - - Args: - mixture: Normal Mean Mixture - size: length of sample - - Returns: sample of given size - - Raises: - ValueError: If mixture is not a Normal Mean Mixture - - """ - - if not isinstance(mixture, NormalMeanMixtures): - raise ValueError("Mixture must be NormalMeanMixtures") - mixing_values = mixture.params.distribution.rvs(size=size) - normal_values = scipy.stats.norm.rvs(size=size) - return mixing_values + mixture.params.sigma * normal_values diff --git a/src/generators/nmv_generator.py b/src/generators/nmv_generator.py index 5ff3221..8270a1a 100644 --- a/src/generators/nmv_generator.py +++ b/src/generators/nmv_generator.py @@ -9,7 +9,7 @@ class NMVGenerator(AbstractGenerator): @staticmethod - def classical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: + def generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: """Generate a sample of given size. Classical form of NMVM Args: @@ -27,29 +27,10 @@ def classical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: raise ValueError("Mixture must be NormalMeanMixtures") mixing_values = mixture.params.distribution.rvs(size=size) normal_values = scipy.stats.norm.rvs(size=size) + if mixture.mixture_form == "canonical": + return mixture.params.alpha + mixture.params.mu * mixing_values + (mixing_values ** 0.5) * normal_values return ( mixture.params.alpha + mixture.params.beta * mixing_values + mixture.params.gamma * (mixing_values**0.5) * normal_values - ) - - @staticmethod - def canonical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: - """Generate a sample of given size. Canonical form of NMVM - - Args: - mixture: Normal Mean Variance Mixtures - size: length of sample - - Returns: sample of given size - - Raises: - ValueError: If mixture type is not Normal Mean Variance Mixtures - - """ - - if not isinstance(mixture, NormalMeanVarianceMixtures): - raise ValueError("Mixture must be NormalMeanMixtures") - mixing_values = mixture.params.distribution.rvs(size=size) - normal_values = scipy.stats.norm.rvs(size=size) - return mixture.params.alpha + mixture.params.mu * mixing_values + (mixing_values**0.5) * normal_values + ) \ No newline at end of file diff --git a/src/generators/nv_generator.py b/src/generators/nv_generator.py index faa55f5..712e918 100644 --- a/src/generators/nv_generator.py +++ b/src/generators/nv_generator.py @@ -9,7 +9,7 @@ class NVGenerator(AbstractGenerator): @staticmethod - def classical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: + def generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: """Generate a sample of given size. Classical form of NVM Args: @@ -27,25 +27,6 @@ def classical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: raise ValueError("Mixture must be NormalMeanMixtures") mixing_values = mixture.params.distribution.rvs(size=size) normal_values = scipy.stats.norm.rvs(size=size) - return mixture.params.alpha + mixture.params.gamma * (mixing_values**0.5) * normal_values - - @staticmethod - def canonical_generate(mixture: AbstractMixtures, size: int) -> tpg.NDArray: - """Generate a sample of given size. Canonical form of NVM - - Args: - mixture: Normal Variance Mixtures - size: length of sample - - Returns: sample of given size - - Raises: - ValueError: If mixture type is not Normal Variance Mixtures - - """ - - if not isinstance(mixture, NormalVarianceMixtures): - raise ValueError("Mixture must be NormalMeanMixtures") - mixing_values = mixture.params.distribution.rvs(size=size) - normal_values = scipy.stats.norm.rvs(size=size) - return mixture.params.alpha + (mixing_values**0.5) * normal_values + if mixture.mixture_form == "canonical": + return mixture.params.alpha + (mixing_values ** 0.5) * normal_values + return mixture.params.alpha + mixture.params.gamma * (mixing_values**0.5) * normal_values \ No newline at end of file diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index 13c8e42..77b927c 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -20,6 +20,7 @@ def __init__(self, mixture_form: str, **kwargs: Any) -> None: mixture_form: Form of Mixture classical or Canonical **kwargs: Parameters of Mixture """ + self.mixture_form = mixture_form if mixture_form == "classical": self.params = self._params_validation(self._classical_collector, kwargs) elif mixture_form == "canonical": diff --git a/tests/generators/nm_generator/test_mixing_normal.py b/tests/generators/nm_generator/test_mixing_normal.py index ea85692..49af0ae 100644 --- a/tests/generators/nm_generator/test_mixing_normal.py +++ b/tests/generators/nm_generator/test_mixing_normal.py @@ -16,7 +16,7 @@ class TestMixingNormal: ) def test_classic_generate_variance_0(self, mixing_variance: float, expected_variance: float) -> None: mixture = NormalMeanMixtures("classical", alpha=0, beta=mixing_variance**0.5, gamma=1, distribution=norm) - sample = self.generator.classical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_variance = ndimage.variance(sample) assert actual_variance == pytest.approx(expected_variance, 0.1) @@ -24,7 +24,7 @@ def test_classic_generate_variance_0(self, mixing_variance: float, expected_vari def test_classic_generate_variance_1(self, beta: float) -> None: expected_variance = beta**2 + 1 mixture = NormalMeanMixtures("classical", alpha=0, beta=beta, gamma=1, distribution=norm) - sample = self.generator.classical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_variance = ndimage.variance(sample) assert actual_variance == pytest.approx(expected_variance, 0.1) @@ -32,7 +32,7 @@ def test_classic_generate_variance_1(self, beta: float) -> None: def test_classic_generate_variance_2(self, beta: float, gamma: float) -> None: expected_variance = beta**2 + gamma**2 mixture = NormalMeanMixtures("classical", alpha=0, beta=beta, gamma=gamma, distribution=norm) - sample = self.generator.classical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_variance = ndimage.variance(sample) assert actual_variance == pytest.approx(expected_variance, 0.1) @@ -40,14 +40,14 @@ def test_classic_generate_variance_2(self, beta: float, gamma: float) -> None: def test_classic_generate_mean(self, beta: float, gamma: float) -> None: expected_mean = 0 mixture = NormalMeanMixtures("classical", alpha=0, beta=beta, gamma=gamma, distribution=norm) - sample = self.generator.classical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_mean = np.mean(np.array(sample)) assert abs(actual_mean - expected_mean) < 1 @pytest.mark.parametrize("expected_size", np.random.randint(0, 100, size=50)) def test_classic_generate_size(self, expected_size: int) -> None: mixture = NormalMeanMixtures("classical", alpha=0, beta=1, gamma=1, distribution=norm) - sample = self.generator.classical_generate(mixture, expected_size) + sample = self.generator.generate(mixture, expected_size) actual_size = np.size(sample) assert actual_size == expected_size @@ -56,7 +56,7 @@ def test_classic_generate_size(self, expected_size: int) -> None: ) def test_canonical_generate_variance_0(self, mixing_variance: float, expected_variance: float) -> None: mixture = NormalMeanMixtures("canonical", sigma=1, distribution=norm(0, mixing_variance**0.5)) - sample = self.generator.canonical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_variance = ndimage.variance(sample) assert actual_variance == pytest.approx(expected_variance, 0.1) @@ -64,7 +64,7 @@ def test_canonical_generate_variance_0(self, mixing_variance: float, expected_va def test_canonical_generate_variance_1(self, sigma: float) -> None: expected_variance = sigma**2 + 1 mixture = NormalMeanMixtures("canonical", sigma=sigma, distribution=norm) - sample = self.generator.canonical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_variance = ndimage.variance(sample) assert actual_variance == pytest.approx(expected_variance, 0.1) @@ -72,7 +72,7 @@ def test_canonical_generate_variance_1(self, sigma: float) -> None: def test_canonical_generate_variance_2(self, mixing_variance: float, sigma: float) -> None: expected_variance = mixing_variance + sigma**2 mixture = NormalMeanMixtures("canonical", sigma=sigma, distribution=norm(0, mixing_variance**0.5)) - sample = self.generator.canonical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_variance = ndimage.variance(sample) assert actual_variance == pytest.approx(expected_variance, 0.1) @@ -80,13 +80,13 @@ def test_canonical_generate_variance_2(self, mixing_variance: float, sigma: floa def test_canonical_generate_mean(self, sigma: float) -> None: expected_mean = 0 mixture = NormalMeanMixtures("canonical", sigma=sigma, distribution=norm) - sample = self.generator.canonical_generate(mixture, self.test_mixture_size) + sample = self.generator.generate(mixture, self.test_mixture_size) actual_mean = np.mean(np.array(sample)) assert abs(actual_mean - expected_mean) < 1 @pytest.mark.parametrize("expected_size", [*np.random.randint(0, 100, size=50), 0, 1, 1000000]) def test_canonical_generate_size(self, expected_size: int) -> None: mixture = NormalMeanMixtures("canonical", sigma=1, distribution=norm) - sample = self.generator.canonical_generate(mixture, expected_size) + sample = self.generator.generate(mixture, expected_size) actual_size = np.size(sample) assert actual_size == expected_size From 08756db048f6f67df93c15629e80897ff9c5ad99 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Tue, 13 May 2025 18:27:03 +0300 Subject: [PATCH 02/24] refactor: Rename folder algorithms to procedures --- src/{algorithms => procedures}/__init__.py | 0 .../param_algorithms => procedures/parametric}/__init__.py | 0 .../parametric}/nm_param_algorithms/__init__.py | 0 .../parametric}/nv_param_algorithms/__init__.py | 0 .../parametric}/nvm_param_algorithms/__init__.py | 0 .../semiparametric}/__init__.py | 0 .../semiparametric}/nm_semi_param_algorithms/__init__.py | 0 .../nm_semi_param_algorithms/g_estimation_convolution.py | 0 .../nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py | 0 .../nm_semi_param_algorithms/sigma_estimation_empirical.py | 0 .../semiparametric}/nv_semi_param_algorithms/__init__.py | 0 .../nv_semi_param_algorithms/g_estimation_given_mu.py | 0 .../nv_semi_param_algorithms/g_estimation_post_widder.py | 0 .../semiparametric}/nvm_semi_param_algorithms/__init__.py | 0 .../nvm_semi_param_algorithms/g_estimation_given_mu.py | 0 .../nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py | 0 .../nvm_semi_param_algorithms/g_estimation_post_widder.py | 0 .../semiparametric}/nvm_semi_param_algorithms/mu_estimation.py | 0 .../support_algorithms => procedures/support}/log_rqmc.py | 0 src/{algorithms/support_algorithms => procedures/support}/rqmc.py | 0 20 files changed, 0 insertions(+), 0 deletions(-) rename src/{algorithms => procedures}/__init__.py (100%) rename src/{algorithms/param_algorithms => procedures/parametric}/__init__.py (100%) rename src/{algorithms/param_algorithms => procedures/parametric}/nm_param_algorithms/__init__.py (100%) rename src/{algorithms/param_algorithms => procedures/parametric}/nv_param_algorithms/__init__.py (100%) rename src/{algorithms/param_algorithms => procedures/parametric}/nvm_param_algorithms/__init__.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/__init__.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nm_semi_param_algorithms/__init__.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nm_semi_param_algorithms/g_estimation_convolution.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nm_semi_param_algorithms/sigma_estimation_empirical.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nv_semi_param_algorithms/__init__.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nv_semi_param_algorithms/g_estimation_given_mu.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nv_semi_param_algorithms/g_estimation_post_widder.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nvm_semi_param_algorithms/__init__.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nvm_semi_param_algorithms/g_estimation_given_mu.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nvm_semi_param_algorithms/g_estimation_post_widder.py (100%) rename src/{algorithms/semiparam_algorithms => procedures/semiparametric}/nvm_semi_param_algorithms/mu_estimation.py (100%) rename src/{algorithms/support_algorithms => procedures/support}/log_rqmc.py (100%) rename src/{algorithms/support_algorithms => procedures/support}/rqmc.py (100%) diff --git a/src/algorithms/__init__.py b/src/procedures/__init__.py similarity index 100% rename from src/algorithms/__init__.py rename to src/procedures/__init__.py diff --git a/src/algorithms/param_algorithms/__init__.py b/src/procedures/parametric/__init__.py similarity index 100% rename from src/algorithms/param_algorithms/__init__.py rename to src/procedures/parametric/__init__.py diff --git a/src/algorithms/param_algorithms/nm_param_algorithms/__init__.py b/src/procedures/parametric/nm_param_algorithms/__init__.py similarity index 100% rename from src/algorithms/param_algorithms/nm_param_algorithms/__init__.py rename to src/procedures/parametric/nm_param_algorithms/__init__.py diff --git a/src/algorithms/param_algorithms/nv_param_algorithms/__init__.py b/src/procedures/parametric/nv_param_algorithms/__init__.py similarity index 100% rename from src/algorithms/param_algorithms/nv_param_algorithms/__init__.py rename to src/procedures/parametric/nv_param_algorithms/__init__.py diff --git a/src/algorithms/param_algorithms/nvm_param_algorithms/__init__.py b/src/procedures/parametric/nvm_param_algorithms/__init__.py similarity index 100% rename from src/algorithms/param_algorithms/nvm_param_algorithms/__init__.py rename to src/procedures/parametric/nvm_param_algorithms/__init__.py diff --git a/src/algorithms/semiparam_algorithms/__init__.py b/src/procedures/semiparametric/__init__.py similarity index 100% rename from src/algorithms/semiparam_algorithms/__init__.py rename to src/procedures/semiparametric/__init__.py diff --git a/src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/__init__.py b/src/procedures/semiparametric/nm_semi_param_algorithms/__init__.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/__init__.py rename to src/procedures/semiparametric/nm_semi_param_algorithms/__init__.py diff --git a/src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/g_estimation_convolution.py b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/g_estimation_convolution.py rename to src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py diff --git a/src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py rename to src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py diff --git a/src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/sigma_estimation_empirical.py b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nm_semi_param_algorithms/sigma_estimation_empirical.py rename to src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py diff --git a/src/algorithms/semiparam_algorithms/nv_semi_param_algorithms/__init__.py b/src/procedures/semiparametric/nv_semi_param_algorithms/__init__.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nv_semi_param_algorithms/__init__.py rename to src/procedures/semiparametric/nv_semi_param_algorithms/__init__.py diff --git a/src/algorithms/semiparam_algorithms/nv_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nv_semi_param_algorithms/g_estimation_given_mu.py rename to src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py diff --git a/src/algorithms/semiparam_algorithms/nv_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nv_semi_param_algorithms/g_estimation_post_widder.py rename to src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/__init__.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/__init__.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/__init__.py rename to src/procedures/semiparametric/nvm_semi_param_algorithms/__init__.py diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu.py rename to src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py rename to src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_post_widder.py rename to src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/mu_estimation.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py similarity index 100% rename from src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/mu_estimation.py rename to src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py diff --git a/src/algorithms/support_algorithms/log_rqmc.py b/src/procedures/support/log_rqmc.py similarity index 100% rename from src/algorithms/support_algorithms/log_rqmc.py rename to src/procedures/support/log_rqmc.py diff --git a/src/algorithms/support_algorithms/rqmc.py b/src/procedures/support/rqmc.py similarity index 100% rename from src/algorithms/support_algorithms/rqmc.py rename to src/procedures/support/rqmc.py From 1062e1c7eaddc1b7b252c925cea69341905985f0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sat, 24 May 2025 15:45:37 +0300 Subject: [PATCH 03/24] refactor: Rename member functions algorithm to compute --- src/estimators/abstract_estimator.py | 2 +- src/mixtures/nm_mixture.py | 4 ++-- src/mixtures/nmv_mixture.py | 4 ++-- src/mixtures/nv_mixture.py | 4 ++-- src/procedures/__init__.py | 18 +++++++++--------- .../g_estimation_convolution.py | 2 +- .../sigma_estimation_eigenvalue_based.py | 2 +- .../sigma_estimation_empirical.py | 2 +- .../g_estimation_given_mu.py | 2 +- .../g_estimation_post_widder.py | 2 +- .../g_estimation_given_mu.py | 2 +- .../g_estimation_given_mu_rqmc_based.py | 4 ++-- .../g_estimation_post_widder.py | 2 +- .../nvm_semi_param_algorithms/mu_estimation.py | 4 ++-- src/procedures/support/log_rqmc.py | 2 +- 15 files changed, 28 insertions(+), 28 deletions(-) diff --git a/src/estimators/abstract_estimator.py b/src/estimators/abstract_estimator.py index 5afaea2..e3fe6d8 100644 --- a/src/estimators/abstract_estimator.py +++ b/src/estimators/abstract_estimator.py @@ -2,7 +2,7 @@ from numpy import _typing -from src.algorithms import ALGORITHM_REGISTRY +from src.procedures import ALGORITHM_REGISTRY from src.estimators.estimate_result import EstimateResult from src.register.algorithm_purpose import AlgorithmPurpose diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 33c7580..952ea2d 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -7,8 +7,8 @@ from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.rqmc import RQMC from src.mixtures.abstract_mixture import AbstractMixtures diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index c090ff3..ee9a88a 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -7,8 +7,8 @@ from scipy.stats import geninvgauss, norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.rqmc import RQMC from src.mixtures.abstract_mixture import AbstractMixtures diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index 7905e6e..828ee31 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -7,8 +7,8 @@ from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.rqmc import RQMC from src.mixtures.abstract_mixture import AbstractMixtures diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index eaddba8..35a5591 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -1,28 +1,28 @@ -from src.algorithms.semiparam_algorithms.nm_semi_param_algorithms.g_estimation_convolution import ( +from src.procedures.semiparametric.nm_semi_param_algorithms.g_estimation_convolution import ( NMSemiParametricGEstimation, ) -from src.algorithms.semiparam_algorithms.nm_semi_param_algorithms.sigma_estimation_eigenvalue_based import ( +from src.procedures.semiparametric.nm_semi_param_algorithms.sigma_estimation_eigenvalue_based import ( SemiParametricMeanSigmaEstimationEigenvalueBased, ) -from src.algorithms.semiparam_algorithms.nm_semi_param_algorithms.sigma_estimation_empirical import ( +from src.procedures.semiparametric.nm_semi_param_algorithms.sigma_estimation_empirical import ( SemiParametricMeanSigmaEstimationEmpirical, ) -from src.algorithms.semiparam_algorithms.nv_semi_param_algorithms.g_estimation_given_mu import ( +from src.procedures.semiparametric.nv_semi_param_algorithms.g_estimation_given_mu import ( SemiParametricNVEstimation, ) -from src.algorithms.semiparam_algorithms.nv_semi_param_algorithms.g_estimation_post_widder import ( +from src.procedures.semiparametric.nv_semi_param_algorithms.g_estimation_post_widder import ( NVSemiParametricGEstimationPostWidder, ) -from src.algorithms.semiparam_algorithms.nvm_semi_param_algorithms.g_estimation_given_mu import ( +from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_given_mu import ( SemiParametricGEstimationGivenMu, ) -from src.algorithms.semiparam_algorithms.nvm_semi_param_algorithms.g_estimation_given_mu_rqmc_based import ( +from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_given_mu_rqmc_based import ( SemiParametricGEstimationGivenMuRQMCBased, ) -from src.algorithms.semiparam_algorithms.nvm_semi_param_algorithms.g_estimation_post_widder import ( +from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_post_widder import ( SemiParametricGEstimationPostWidder, ) -from src.algorithms.semiparam_algorithms.nvm_semi_param_algorithms.mu_estimation import SemiParametricMuEstimation +from src.procedures.semiparametric.nvm_semi_param_algorithms.mu_estimation import (SemiParametricMuEstimation) from src.register.algorithm_purpose import AlgorithmPurpose from src.register.register import Registry diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py index 3feb6cf..35bad05 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py @@ -93,7 +93,7 @@ def cdf(self, X: np.ndarray) -> np.ndarray: return F_x.real - def algorithm(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: + def compute(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: inv = self.BohmanA(N=self.bohman_n, delta=self.bohman_delta) inv.fit(self.characteristic_function_xi) x_data_array = np.array(self.x_data, dtype=np.float64) diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py index 1148610..487965f 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py @@ -86,7 +86,7 @@ def _build_matrix(self, tau: float) -> np.ndarray: matrix[i, j] = self._alpha(t[i] - t[j], tau) return matrix - def algorithm(self, sample: np.ndarray) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: """Estimate sigma. Args: diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py index 0ce471c..6e5e2f8 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py @@ -51,7 +51,7 @@ def _validate_kwargs(**kwargs: ParamsAnnotation) -> float: raise ValueError("Expected a positive float as parameter `t`.") return t - def algorithm(self, sample: np.ndarray) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: """Estimate the sigma parameter. Args: diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py index 7c810cc..dc72d31 100644 --- a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py +++ b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py @@ -173,6 +173,6 @@ def compute_integrals_for_x(self, x: float) -> float: total = np.sum(first_integral + second_integral) / self.denominator return max(0.0, total.real) - def algorithm(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: + def compute(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py index 0c681d3..ff703c0 100644 --- a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py +++ b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py @@ -77,7 +77,7 @@ def p_x_estimation(self, x: float) -> float: ) return result.real - def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np._typing.NDArray) -> EstimateResult: """Estimate g(x) Args: diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py index cbda86e..097891a 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py @@ -179,6 +179,6 @@ def compute_integrals_for_x(self, x: float) -> float: total = np.sum(first_integral + second_integral) / self.denominator return max(0.0, total.real) - def algorithm(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: + def compute(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py index 2f3da81..67a34b0 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py @@ -7,7 +7,7 @@ from scipy.integrate import quad_vec from scipy.special import gamma -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.rqmc import RQMC from src.estimators.estimate_result import EstimateResult MU_DEFAULT_VALUE = 1.0 @@ -171,6 +171,6 @@ def compute_integrals_for_x(self, x: float) -> float: total = (first_integral + second_integral) / self.denominator return max(0.0, total.real) - def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np._typing.NDArray) -> EstimateResult: y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py index 50c518e..df994fb 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py @@ -80,7 +80,7 @@ def p_x_estimation(self, x: float) -> float: ) return result.real - def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np._typing.NDArray) -> EstimateResult: """Estimate g(x) Args: diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py index 9e2dad0..0e9cc24 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py @@ -84,7 +84,7 @@ def __w(self, p: float, sample: np._typing.NDArray) -> float: y += e * self.omega(x) return y - def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np._typing.NDArray) -> EstimateResult: """Root of this function is an estimation of mu Args: @@ -96,7 +96,7 @@ def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: if self.__w(0, sample) == 0: return EstimateResult(value=0, success=True) if self.__w(0, sample) > 0: - second_result = self.algorithm(-1 * sample) + second_result = self.compute(-1 * sample) return EstimateResult(value=-1 * second_result.value, success=second_result.success) if self.__w(self.m, sample) < 0: return EstimateResult(value=self.m, success=False) diff --git a/src/procedures/support/log_rqmc.py b/src/procedures/support/log_rqmc.py index 4423294..a10f221 100644 --- a/src/procedures/support/log_rqmc.py +++ b/src/procedures/support/log_rqmc.py @@ -4,7 +4,7 @@ import numpy._typing as tpg import scipy -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.rqmc import RQMC class LogRQMC(RQMC): From bfdf5716486ee6355ef12a0c8f77e077899d6d3a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sat, 24 May 2025 20:19:49 +0300 Subject: [PATCH 04/24] refactor: Fix naming of estimators --- .../semiparametric/abstract_semiparametric_estimator.py | 2 +- src/estimators/semiparametric/nm_semiparametric_estimator.py | 4 ++-- src/estimators/semiparametric/nmv_semiparametric_estimator.py | 4 ++-- src/estimators/semiparametric/nv_semiparametric_estimator.py | 4 ++-- 4 files changed, 7 insertions(+), 7 deletions(-) diff --git a/src/estimators/semiparametric/abstract_semiparametric_estimator.py b/src/estimators/semiparametric/abstract_semiparametric_estimator.py index b26db4b..466f295 100644 --- a/src/estimators/semiparametric/abstract_semiparametric_estimator.py +++ b/src/estimators/semiparametric/abstract_semiparametric_estimator.py @@ -4,7 +4,7 @@ from src.estimators.estimate_result import EstimateResult -class AbstractSemiParametricEstimator(AbstractEstimator): +class AbstractSemiparEstim(AbstractEstimator): def __init__(self, algorithm_name: str, params: dict | None = None) -> None: super().__init__(algorithm_name, params) diff --git a/src/estimators/semiparametric/nm_semiparametric_estimator.py b/src/estimators/semiparametric/nm_semiparametric_estimator.py index 958ec5c..0513521 100644 --- a/src/estimators/semiparametric/nm_semiparametric_estimator.py +++ b/src/estimators/semiparametric/nm_semiparametric_estimator.py @@ -1,11 +1,11 @@ from numpy import _typing from src.estimators.estimate_result import EstimateResult -from src.estimators.semiparametric.abstract_semiparametric_estimator import AbstractSemiParametricEstimator +from src.estimators.semiparametric.abstract_semiparametric_estimator import AbstractSemiparEstim from src.register.algorithm_purpose import AlgorithmPurpose -class NMSemiParametricEstimator(AbstractSemiParametricEstimator): +class NMSemiparEstim(AbstractSemiparEstim): def __init__(self, algorithm_name: str, params: dict | None = None) -> None: super().__init__(algorithm_name, params) self._purpose = AlgorithmPurpose.NM_SEMIPARAMETRIC diff --git a/src/estimators/semiparametric/nmv_semiparametric_estimator.py b/src/estimators/semiparametric/nmv_semiparametric_estimator.py index e1f456d..93403e9 100644 --- a/src/estimators/semiparametric/nmv_semiparametric_estimator.py +++ b/src/estimators/semiparametric/nmv_semiparametric_estimator.py @@ -1,11 +1,11 @@ from numpy import _typing from src.estimators.estimate_result import EstimateResult -from src.estimators.semiparametric.abstract_semiparametric_estimator import AbstractSemiParametricEstimator +from src.estimators.semiparametric.abstract_semiparametric_estimator import AbstractSemiparEstim from src.register.algorithm_purpose import AlgorithmPurpose -class NMVSemiParametricEstimator(AbstractSemiParametricEstimator): +class NMVSemiparEstim(AbstractSemiparEstim): def __init__(self, algorithm_name: str, params: dict | None = None) -> None: super().__init__(algorithm_name, params) self._purpose = AlgorithmPurpose.NMV_SEMIPARAMETRIC diff --git a/src/estimators/semiparametric/nv_semiparametric_estimator.py b/src/estimators/semiparametric/nv_semiparametric_estimator.py index 95f0ce5..d2ec1d8 100644 --- a/src/estimators/semiparametric/nv_semiparametric_estimator.py +++ b/src/estimators/semiparametric/nv_semiparametric_estimator.py @@ -1,11 +1,11 @@ from numpy import _typing from src.estimators.estimate_result import EstimateResult -from src.estimators.semiparametric.abstract_semiparametric_estimator import AbstractSemiParametricEstimator +from src.estimators.semiparametric.abstract_semiparametric_estimator import AbstractSemiparEstim from src.register.algorithm_purpose import AlgorithmPurpose -class NVSemiParametricEstimator(AbstractSemiParametricEstimator): +class NVSemiparEstim(AbstractSemiparEstim): def __init__(self, algorithm_name: str, params: dict | None = None) -> None: super().__init__(algorithm_name, params) self._purpose = AlgorithmPurpose.NV_SEMIPARAMETRIC From 76838d6172b0e176ea6e65f555c44f88deac61d7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sun, 25 May 2025 17:22:54 +0300 Subject: [PATCH 05/24] refactor: Update register names --- src/procedures/__init__.py | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index 35a5591..ad947e9 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -28,19 +28,19 @@ ALGORITHM_REGISTRY: Registry = Registry() ALGORITHM_REGISTRY.register("mu_estimation", AlgorithmPurpose.NMV_SEMIPARAMETRIC)(SemiParametricMuEstimation) -ALGORITHM_REGISTRY.register("g_estimation_given_mu", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( +ALGORITHM_REGISTRY.register("density_estim_inv_mellin_quad_int", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( SemiParametricGEstimationGivenMu ) -ALGORITHM_REGISTRY.register("g_estimation_convolution", AlgorithmPurpose.NM_SEMIPARAMETRIC)(NMSemiParametricGEstimation) -ALGORITHM_REGISTRY.register("sigma_estimation_eigenvalue_based", AlgorithmPurpose.NM_SEMIPARAMETRIC)( +ALGORITHM_REGISTRY.register("density_estim_deconv", AlgorithmPurpose.NM_SEMIPARAMETRIC)(NMSemiParametricGEstimation) +ALGORITHM_REGISTRY.register("sigma_estimation_eigenvalue", AlgorithmPurpose.NM_SEMIPARAMETRIC)( SemiParametricMeanSigmaEstimationEigenvalueBased ) -ALGORITHM_REGISTRY.register("sigma_estimation_empirical", AlgorithmPurpose.NM_SEMIPARAMETRIC)( +ALGORITHM_REGISTRY.register("sigma_estimation_laplace", AlgorithmPurpose.NM_SEMIPARAMETRIC)( SemiParametricMeanSigmaEstimationEmpirical ) -ALGORITHM_REGISTRY.register("g_estimation_given_mu_rqmc_based", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( +ALGORITHM_REGISTRY.register("density_estim_inv_mellin_rqmc_int", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( SemiParametricGEstimationGivenMuRQMCBased ) -ALGORITHM_REGISTRY.register("g_estimation_post_widder", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( +ALGORITHM_REGISTRY.register("density_estim_post_widder", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( SemiParametricGEstimationPostWidder ) From 7e47e6369baa3f70afb517e0e90eceba745738b2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sun, 25 May 2025 17:31:47 +0300 Subject: [PATCH 06/24] refactor: Update register names --- .../test_sigma_estimation_eigenvalue_based.py | 10 +++++----- .../test_sigma_estimation_empirical.py | 4 ++-- .../test_g_estimation_given_mu.py | 2 +- .../semiparametric_g_estimation/test_post_widder.py | 4 ++-- tests/algorithms/nv_algorithms/test_g_estimation.py | 2 +- tests/algorithms/nv_algorithms/test_post_widder_nv.py | 4 ++-- 6 files changed, 13 insertions(+), 13 deletions(-) diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py index eb00242..9b2243f 100644 --- a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py +++ b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py @@ -29,7 +29,7 @@ def test_sigma_estimation_eigenvalue_based_expon_single( ) estimator = NMSemiParametricEstimator( - "sigma_estimation_eigenvalue_based", {"k": k, "l": l, "eps": eps, "search_area": search_area} + "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) print(est.value) @@ -55,7 +55,7 @@ def test_sigma_estimation_eigenvalue_based_expon_1( ) estimator = NMSemiParametricEstimator( - "sigma_estimation_eigenvalue_based", {"k": k, "l": l, "eps": eps, "search_area": search_area} + "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) print(est.value) @@ -81,7 +81,7 @@ def test_sigma_estimation_eigenvalue_based_expon_2( ) estimator = NMSemiParametricEstimator( - "sigma_estimation_eigenvalue_based", {"k": k, "l": l, "eps": eps, "search_area": search_area} + "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) print(est.value) @@ -107,7 +107,7 @@ def test_sigma_estimation_eigenvalue_based_expon_3( ) estimator = NMSemiParametricEstimator( - "sigma_estimation_eigenvalue_based", {"k": k, "l": l, "eps": eps, "search_area": search_area} + "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) print(est.value) @@ -132,7 +132,7 @@ def test_sigma_estimation_eigenvalue_based_uniform( ) estimator = NMSemiParametricEstimator( - "sigma_estimation_eigenvalue_based", {"k": k, "l": l, "eps": eps, "search_area": search_area} + "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) print(est.value) diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py index daf5a7f..4147449 100644 --- a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py +++ b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py @@ -20,7 +20,7 @@ def test_sigma_estimation_empirical_expon(self, real_sigma: float, sample_len: i answer_list = [] for alpha in [round(x, 5) for x in [i * 0.0001 for i in range(1, 10000)]]: t = math.sqrt(alpha * math.log(sample_len)) / (2 * search_area) - estimator = NMSemiParametricEstimator("sigma_estimation_empirical", {"t": t}) + estimator = NMSemiParametricEstimator("sigma_estimation_laplace", {"t": t}) est = estimator.estimate(sample) left = (est.value**2 - real_sigma**2) ** 0.5 right = ( @@ -39,7 +39,7 @@ def test_sigma_estimation_empirical_uniform(self, real_sigma: float, sample_len: M = 10 for alpha in [round(x, 4) for x in [i * 0.0001 for i in range(1, 10000)]]: t = math.sqrt(alpha * math.log(sample_len)) / (2 * search_area) - estimator = NMSemiParametricEstimator("sigma_estimation_empirical", {"t": t}) + estimator = NMSemiParametricEstimator("sigma_estimation_laplace", {"t": t}) est = estimator.estimate(sample) left = abs(est.value**2 - real_sigma**2) ** 0.5 right = 4 * M * search_area / math.sqrt(alpha * math.log(sample_len)) diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py index d1f40e0..2b8ceb8 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py +++ b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py @@ -19,7 +19,7 @@ def test_g_estimation_expon(self, x_data) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=given_mu, distribution=expon) sample = NMVGenerator().canonical_generate(mixture, n) estimator = NMVSemiParametricEstimator( - "g_estimation_given_mu", {"x_data": x_data, "u_value": 7.6, "v_value": 0.9} + "density_estim_inv_mellin_quad_int", {"x_data": x_data, "u_value": 7.6, "v_value": 0.9} ) est = estimator.estimate(sample) error = 0.0 diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py index 2e8d4b3..f817ca2 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py +++ b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py @@ -21,7 +21,7 @@ def test_post_widder_expon(self, mu, sigma, degree, sample_size) -> None: x_data = np.linspace(0.5, 10.0, 30) estimator = NMVSemiParametricEstimator( - "g_estimation_post_widder", {"x_data": x_data, "mu": mu, "sigma": sigma, "n": degree} + "density_estim_post_widder", {"x_data": x_data, "mu": mu, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) est_data = est.list_value @@ -44,7 +44,7 @@ def test_post_widder_gamma(self, mu, sigma, degree, sample_size, a) -> None: x_data = np.linspace(0.5, 10.0, 30) estimator = NMVSemiParametricEstimator( - "g_estimation_post_widder", {"x_data": x_data, "mu": mu, "sigma": sigma, "n": degree} + "density_estim_post_widder", {"x_data": x_data, "mu": mu, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) est_data = est.list_value diff --git a/tests/algorithms/nv_algorithms/test_g_estimation.py b/tests/algorithms/nv_algorithms/test_g_estimation.py index f618be5..a2068ab 100644 --- a/tests/algorithms/nv_algorithms/test_g_estimation.py +++ b/tests/algorithms/nv_algorithms/test_g_estimation.py @@ -18,7 +18,7 @@ def test_g_estimation_expon(self, x_data) -> None: mixture = NormalVarianceMixtures("canonical", alpha=0, distribution=expon) sample = NVGenerator().canonical_generate(mixture, n) estimator = NVSemiParametricEstimator( - "g_estimation_given_mu", {"x_data": x_data, "u_value": 7.6, "v_value": 0.9} + "density_estim_inv_mellin_quad_int", {"x_data": x_data, "u_value": 7.6, "v_value": 0.9} ) est = estimator.estimate(sample) error = 0.0 diff --git a/tests/algorithms/nv_algorithms/test_post_widder_nv.py b/tests/algorithms/nv_algorithms/test_post_widder_nv.py index e2fdee3..02c267d 100644 --- a/tests/algorithms/nv_algorithms/test_post_widder_nv.py +++ b/tests/algorithms/nv_algorithms/test_post_widder_nv.py @@ -21,7 +21,7 @@ def test_post_widder_expon(self, sigma, degree, sample_size) -> None: x_data = np.linspace(0.5, 10.0, 30) estimator = NVSemiParametricEstimator( - "g_estimation_post_widder", {"x_data": x_data, "sigma": sigma, "n": degree} + "density_estim_post_widder", {"x_data": x_data, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) est_data = est.list_value @@ -44,7 +44,7 @@ def test_post_widder_gamma(self, sigma, degree, sample_size, a) -> None: x_data = np.linspace(0.5, 10.0, 30) estimator = NVSemiParametricEstimator( - "g_estimation_post_widder", {"x_data": x_data, "sigma": sigma, "n": degree} + "density_estim_post_widder", {"x_data": x_data, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) est_data = est.list_value From 78a7b35493eb400970299492427f14a8eee2b151 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sun, 25 May 2025 17:54:36 +0300 Subject: [PATCH 07/24] refactor: Update procedures class names --- src/procedures/__init__.py | 28 +++++++++---------- .../g_estimation_convolution.py | 2 +- .../sigma_estimation_eigenvalue_based.py | 2 +- .../sigma_estimation_empirical.py | 2 +- .../g_estimation_given_mu.py | 2 +- .../g_estimation_post_widder.py | 2 +- .../g_estimation_given_mu.py | 2 +- .../g_estimation_given_mu_rqmc_based.py | 2 +- .../mu_estimation.py | 2 +- 9 files changed, 22 insertions(+), 22 deletions(-) diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index ad947e9..1c88c3b 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -1,45 +1,45 @@ from src.procedures.semiparametric.nm_semi_param_algorithms.g_estimation_convolution import ( - NMSemiParametricGEstimation, + NMEstimationDensityInvMT, ) from src.procedures.semiparametric.nm_semi_param_algorithms.sigma_estimation_eigenvalue_based import ( - SemiParametricMeanSigmaEstimationEigenvalueBased, + NMEstimationSigmaEigenvals, ) from src.procedures.semiparametric.nm_semi_param_algorithms.sigma_estimation_empirical import ( - SemiParametricMeanSigmaEstimationEmpirical, + NMEstimationSigmaLaplace, ) from src.procedures.semiparametric.nv_semi_param_algorithms.g_estimation_given_mu import ( - SemiParametricNVEstimation, + NVEstimationDensityInvMT, ) from src.procedures.semiparametric.nv_semi_param_algorithms.g_estimation_post_widder import ( - NVSemiParametricGEstimationPostWidder, + NMVEstimationDensityPW, ) from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_given_mu import ( - SemiParametricGEstimationGivenMu, + NMVEstimationDensityInvMTquad, ) from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_given_mu_rqmc_based import ( - SemiParametricGEstimationGivenMuRQMCBased, + NMVEstimationDensityInvMTquadRQMCBased, ) from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_post_widder import ( SemiParametricGEstimationPostWidder, ) -from src.procedures.semiparametric.nvm_semi_param_algorithms.mu_estimation import (SemiParametricMuEstimation) +from src.procedures.semiparametric.nvm_semi_param_algorithms.mu_estimation import (NMVEstimationMu) from src.register.algorithm_purpose import AlgorithmPurpose from src.register.register import Registry ALGORITHM_REGISTRY: Registry = Registry() -ALGORITHM_REGISTRY.register("mu_estimation", AlgorithmPurpose.NMV_SEMIPARAMETRIC)(SemiParametricMuEstimation) +ALGORITHM_REGISTRY.register("mu_estimation", AlgorithmPurpose.NMV_SEMIPARAMETRIC)(NMVEstimationMu) ALGORITHM_REGISTRY.register("density_estim_inv_mellin_quad_int", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( - SemiParametricGEstimationGivenMu + NMVEstimationDensityInvMTquad ) -ALGORITHM_REGISTRY.register("density_estim_deconv", AlgorithmPurpose.NM_SEMIPARAMETRIC)(NMSemiParametricGEstimation) +ALGORITHM_REGISTRY.register("density_estim_deconv", AlgorithmPurpose.NM_SEMIPARAMETRIC)(NMEstimationDensityInvMT) ALGORITHM_REGISTRY.register("sigma_estimation_eigenvalue", AlgorithmPurpose.NM_SEMIPARAMETRIC)( - SemiParametricMeanSigmaEstimationEigenvalueBased + NMEstimationSigmaEigenvals ) ALGORITHM_REGISTRY.register("sigma_estimation_laplace", AlgorithmPurpose.NM_SEMIPARAMETRIC)( - SemiParametricMeanSigmaEstimationEmpirical + NMEstimationSigmaLaplace ) ALGORITHM_REGISTRY.register("density_estim_inv_mellin_rqmc_int", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( - SemiParametricGEstimationGivenMuRQMCBased + NMVEstimationDensityInvMTquadRQMCBased ) ALGORITHM_REGISTRY.register("density_estim_post_widder", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( SemiParametricGEstimationPostWidder diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py index 35bad05..29b0396 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py @@ -12,7 +12,7 @@ X_DATA_DEFAULT_VALUE: List[float] = [1.0] -class NMSemiParametricGEstimation: +class NMEstimationDensityInvMT: """Estimation of mixing density function g (xi density function) of NM mixture represented in canonical form Y = xi + sigma*N. diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py index 487965f..b9da043 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py @@ -14,7 +14,7 @@ PARAMETER_KEYS = ["l", "k", "eps", "search_area", "search_density"] -class SemiParametricMeanSigmaEstimationEigenvalueBased: +class NMEstimationSigmaEigenvals: """Estimation of sigma parameter of NM mixture represented in canonical form Y = xi + sigma*N. Args: diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py index 6e5e2f8..22a69af 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py @@ -9,7 +9,7 @@ PARAMETER_KEYS = ["t"] -class SemiParametricMeanSigmaEstimationEmpirical: +class NMEstimationSigmaLaplace: """Estimation of sigma parameter of NM mixture represented in canonical form Y = xi + sigma*N. Args: diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py index dc72d31..5c87200 100644 --- a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py +++ b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py @@ -30,7 +30,7 @@ def v_sequence_default_value(n: float) -> float: INTEGRATION_LIMIT_DEFAULT_VALUE: int = 50 -class SemiParametricNVEstimation: +class NVEstimationDensityInvMT: """Estimation of mixing density function g (xi density function) of NV mixture represented in canonical form Y = alpha + sqrt(xi)*N, where alpha = 0 and mu = 0. diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py index ff703c0..3b5dc7e 100644 --- a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py +++ b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py @@ -13,7 +13,7 @@ X_DATA_DEFAULT_VALUE = [1.0] -class NVSemiParametricGEstimationPostWidder: +class NMVEstimationDensityPW: """Estimation of mixing density function g (xi density function) of NVM mixture represented in classical form Y = alpha + sigma*sqrt(xi)*N, where alpha = 0 and mu, sigma are given. diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py index 097891a..28cbb7d 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py @@ -31,7 +31,7 @@ def v_sequence_default_value(n: float) -> float: INTEGRATION_LIMIT_DEFAULT_VALUE: int = 50 -class SemiParametricGEstimationGivenMu: +class NMVEstimationDensityInvMTquad: """Estimation of mixing density function g (xi density function) of NVM mixture represented in canonical form Y = alpha + mu*xi + sqrt(xi)*N, where alpha = 0 and mu is given. diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py index 67a34b0..d8e1b37 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py @@ -32,7 +32,7 @@ def v_sequence_default_value(n: float) -> float: INTEGRATION_LIMIT_DEFAULT_VALUE: int = 50 -class SemiParametricGEstimationGivenMuRQMCBased: +class NMVEstimationDensityInvMTquadRQMCBased: """Estimation of mixing density function g (xi density function) of NVM mixture represented in canonical form Y = alpha + mu*xi + sqrt(xi)*N, where alpha = 0 and mu is given. diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py index 0e9cc24..48db22e 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py @@ -14,7 +14,7 @@ PARAMETER_KEYS = ["m", "tolerance", "omega", "max_iterations"] -class SemiParametricMuEstimation: +class NMVEstimationMu: """Estimation of mu parameter of NVM mixture represented in canonical form Y = alpha + mu*xi + sqrt(xi)*N, where alpha = 0 From 9fa4aad9949b11361720f87d42961fd3e692eeaf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Wed, 28 May 2025 20:45:06 +0300 Subject: [PATCH 08/24] refactor: update tests --- .../test_sigma_estimation_eigenvalue_based.py | 12 +++---- .../test_sigma_estimation_empirical.py | 6 ++-- .../test_g_estimation_given_mu.py | 4 +-- .../test_post_widder.py | 6 ++-- .../test_semiparametric_mu_estimation.py | 34 +++++++++---------- .../test_validate_kwargs.py | 26 +++++++------- .../nv_algorithms/test_g_estimation.py | 4 +-- .../nv_algorithms/test_post_widder_nv.py | 6 ++-- .../support_algorithms/test_rqmc.py | 2 +- 9 files changed, 50 insertions(+), 50 deletions(-) diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py index 9b2243f..e3b7d64 100644 --- a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py +++ b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py @@ -4,7 +4,7 @@ import pytest from scipy.stats import expon, uniform -from src.estimators.semiparametric.nm_semiparametric_estimator import NMSemiParametricEstimator +from src.estimators.semiparametric.nm_semiparametric_estimator import NMSemiparEstim from src.generators.nm_generator import NMGenerator from src.mixtures.nm_mixture import NormalMeanMixtures @@ -28,7 +28,7 @@ def test_sigma_estimation_eigenvalue_based_expon_single( * (np.log(np.log(sample_len))) ** 0.6 ) - estimator = NMSemiParametricEstimator( + estimator = NMSemiparEstim( "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) @@ -54,7 +54,7 @@ def test_sigma_estimation_eigenvalue_based_expon_1( * (np.log(np.log(sample_len))) ** 0.6 ) - estimator = NMSemiParametricEstimator( + estimator = NMSemiparEstim( "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) @@ -80,7 +80,7 @@ def test_sigma_estimation_eigenvalue_based_expon_2( * (np.log(np.log(sample_len))) ** 0.6 ) - estimator = NMSemiParametricEstimator( + estimator = NMSemiparEstim( "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) @@ -106,7 +106,7 @@ def test_sigma_estimation_eigenvalue_based_expon_3( * (np.log(np.log(sample_len))) ** 0.6 ) - estimator = NMSemiParametricEstimator( + estimator = NMSemiparEstim( "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) @@ -131,7 +131,7 @@ def test_sigma_estimation_eigenvalue_based_uniform( * (np.log(np.log(sample_len))) ** 0.6 ) - estimator = NMSemiParametricEstimator( + estimator = NMSemiparEstim( "sigma_estimation_eigenvalue", {"k": k, "l": l, "eps": eps, "search_area": search_area} ) est = estimator.estimate(sample) diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py index 4147449..550978f 100644 --- a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py +++ b/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py @@ -3,7 +3,7 @@ import pytest from scipy.stats import expon, uniform -from src.estimators.semiparametric.nm_semiparametric_estimator import NMSemiParametricEstimator +from src.estimators.semiparametric.nm_semiparametric_estimator import NMSemiparEstim from src.generators.nm_generator import NMGenerator from src.mixtures.nm_mixture import NormalMeanMixtures @@ -20,7 +20,7 @@ def test_sigma_estimation_empirical_expon(self, real_sigma: float, sample_len: i answer_list = [] for alpha in [round(x, 5) for x in [i * 0.0001 for i in range(1, 10000)]]: t = math.sqrt(alpha * math.log(sample_len)) / (2 * search_area) - estimator = NMSemiParametricEstimator("sigma_estimation_laplace", {"t": t}) + estimator = NMSemiparEstim("sigma_estimation_laplace", {"t": t}) est = estimator.estimate(sample) left = (est.value**2 - real_sigma**2) ** 0.5 right = ( @@ -39,7 +39,7 @@ def test_sigma_estimation_empirical_uniform(self, real_sigma: float, sample_len: M = 10 for alpha in [round(x, 4) for x in [i * 0.0001 for i in range(1, 10000)]]: t = math.sqrt(alpha * math.log(sample_len)) / (2 * search_area) - estimator = NMSemiParametricEstimator("sigma_estimation_laplace", {"t": t}) + estimator = NMSemiparEstim("sigma_estimation_laplace", {"t": t}) est = estimator.estimate(sample) left = abs(est.value**2 - real_sigma**2) ** 0.5 right = 4 * M * search_area / math.sqrt(alpha * math.log(sample_len)) diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py index 2b8ceb8..0cc968f 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py +++ b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py @@ -3,7 +3,7 @@ import pytest from scipy.stats import expon, gamma -from src.estimators.semiparametric.nmv_semiparametric_estimator import NMVSemiParametricEstimator +from src.estimators.semiparametric.nmv_semiparametric_estimator import NMVSemiparEstim from src.generators.nmv_generator import NMVGenerator from src.mixtures.nmv_mixture import NormalMeanVarianceMixtures @@ -18,7 +18,7 @@ def test_g_estimation_expon(self, x_data) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=given_mu, distribution=expon) sample = NMVGenerator().canonical_generate(mixture, n) - estimator = NMVSemiParametricEstimator( + estimator = NMVSemiparEstim( "density_estim_inv_mellin_quad_int", {"x_data": x_data, "u_value": 7.6, "v_value": 0.9} ) est = estimator.estimate(sample) diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py index f817ca2..e3b268c 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py +++ b/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py @@ -3,7 +3,7 @@ from mpmath import ln from scipy.stats import expon, gamma -from src.estimators.semiparametric.nmv_semiparametric_estimator import NMVSemiParametricEstimator +from src.estimators.semiparametric.nmv_semiparametric_estimator import NMVSemiparEstim from src.generators.nmv_generator import NMVGenerator from src.mixtures.nmv_mixture import NormalMeanVarianceMixtures @@ -20,7 +20,7 @@ def test_post_widder_expon(self, mu, sigma, degree, sample_size) -> None: sample = NMVGenerator().classical_generate(mixture, sample_size) x_data = np.linspace(0.5, 10.0, 30) - estimator = NMVSemiParametricEstimator( + estimator = NMVSemiparEstim( "density_estim_post_widder", {"x_data": x_data, "mu": mu, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) @@ -43,7 +43,7 @@ def test_post_widder_gamma(self, mu, sigma, degree, sample_size, a) -> None: sample = NMVGenerator().classical_generate(mixture, sample_size) x_data = np.linspace(0.5, 10.0, 30) - estimator = NMVSemiParametricEstimator( + estimator = NMVSemiparEstim( "density_estim_post_widder", {"x_data": x_data, "mu": mu, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) diff --git a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py b/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py index f1736f8..bcb9154 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py +++ b/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py @@ -3,7 +3,7 @@ import pytest from scipy.stats import expon, gamma, halfnorm, pareto -from src.estimators.semiparametric.nmv_semiparametric_estimator import NMVSemiParametricEstimator +from src.estimators.semiparametric.nmv_semiparametric_estimator import NMVSemiparEstim from src.generators.nmv_generator import NMVGenerator from src.mixtures.nmv_mixture import NormalMeanVarianceMixtures @@ -15,7 +15,7 @@ class TestSemiParametricMuEstimation: def test_mu_estimation_expon_no_parameters(self, real_mu: float) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation") + estimator = NMVSemiparEstim("mu_estimation") est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -23,7 +23,7 @@ def test_mu_estimation_expon_no_parameters(self, real_mu: float) -> None: def test_mu_estimation_expon_no_parameters_small(self, real_mu: float) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation") + estimator = NMVSemiparEstim("mu_estimation") est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -31,7 +31,7 @@ def test_mu_estimation_expon_no_parameters_small(self, real_mu: float) -> None: def test_mu_estimation_pareto_no_parameters(self, real_mu: float) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=pareto(2.62)) sample = self.generator.canonical_generate(mixture, 50000) - estimator = NMVSemiParametricEstimator("mu_estimation") + estimator = NMVSemiparEstim("mu_estimation") est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -39,7 +39,7 @@ def test_mu_estimation_pareto_no_parameters(self, real_mu: float) -> None: def test_mu_estimation_pareto_no_parameters_small(self, real_mu: float) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=pareto(2.62)) sample = self.generator.canonical_generate(mixture, 50000) - estimator = NMVSemiParametricEstimator("mu_estimation") + estimator = NMVSemiparEstim("mu_estimation") est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -47,7 +47,7 @@ def test_mu_estimation_pareto_no_parameters_small(self, real_mu: float) -> None: def test_mu_estimation_halfnorm_no_parameters(self, real_mu: float) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=halfnorm) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation") + estimator = NMVSemiparEstim("mu_estimation") est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -55,7 +55,7 @@ def test_mu_estimation_halfnorm_no_parameters(self, real_mu: float) -> None: def test_mu_estimation_gamma_no_parameters(self, real_mu: float) -> None: mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=gamma(2)) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation") + estimator = NMVSemiparEstim("mu_estimation") est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -64,7 +64,7 @@ def test_mu_estimation_expon_1_parameter_m_positive(self, params: dict) -> None: real_mu = 1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -73,7 +73,7 @@ def test_mu_estimation_expon_1_parameter_m_negative(self, params: dict) -> None: real_mu = -1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -84,7 +84,7 @@ def test_mu_estimation_expon_1_parameter_max_iterations(self, params: dict) -> N real_mu = 1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -93,7 +93,7 @@ def test_mu_estimation_expon_1_parameter_m_is_best_estimation(self, params: dict real_mu = 10 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(est.value == params["m"]) and est.success is False @@ -104,7 +104,7 @@ def test_mu_estimation_expon_2_parameters_tol_positive(self, params: dict) -> No real_mu = 1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -115,7 +115,7 @@ def test_mu_estimation_expon_2_parameters_tol_negative(self, params: dict) -> No real_mu = -1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -132,7 +132,7 @@ def test_mu_estimation_expon_3_parameters_omega_positive(self, params: dict) -> real_mu = 1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -149,7 +149,7 @@ def test_mu_estimation_expon_3_parameters_omega_negative(self, params: dict) -> real_mu = -1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -168,7 +168,7 @@ def test_mu_estimation_expon_3_parameters_all_positive(self, params: dict) -> No real_mu = 1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True @@ -187,6 +187,6 @@ def test_mu_estimation_expon_3_parameters_all_negative(self, params: dict) -> No real_mu = -1 mixture = NormalMeanVarianceMixtures("canonical", alpha=0, mu=real_mu, distribution=expon) sample = self.generator.canonical_generate(mixture, 10000) - estimator = NMVSemiParametricEstimator("mu_estimation", params) + estimator = NMVSemiparEstim("mu_estimation", params) est = estimator.estimate(sample) assert abs(real_mu - est.value) < 1 and est.success is True diff --git a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py b/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py index 781caf5..fe4a2ad 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py +++ b/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py @@ -3,7 +3,7 @@ import numpy as np import pytest -from src.algorithms.semiparam_algorithms.nvm_semi_param_algorithms.mu_estimation import SemiParametricMuEstimation +from src.procedures.semiparametric.nvm_semi_param_algorithms.mu_estimation import NMVEstimationMu def _test_omega(x: float) -> float: @@ -18,7 +18,7 @@ class TestValidateKwargs: ) def test_set_default_params_len_1_value_error_m(self, params: dict) -> None: with pytest.raises(ValueError, match="Expected positive integer as parameter m"): - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", @@ -26,12 +26,12 @@ def test_set_default_params_len_1_value_error_m(self, params: dict) -> None: ) def test_set_default_params_len_1_value_error_tolerance(self, params: dict) -> None: with pytest.raises(ValueError, match="Expected positive float as parameter tolerance"): - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize("params", [{"omega": 1}, {"omega": []}, {"omega": "str"}, {"omega": ()}]) def test_set_default_params_len_3_value_error_omega(self, params: dict) -> None: with pytest.raises(ValueError, match="Expected callable object as parameter omega"): - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", @@ -47,7 +47,7 @@ def test_set_default_params_len_3_value_error_omega(self, params: dict) -> None: ) def test_set_default_params_len_1_value_error_max_iterations(self, params: dict) -> None: with pytest.raises(ValueError, match="Expected positive integer as parameter max_iterations"): - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", @@ -62,7 +62,7 @@ def test_set_default_params_len_1_value_error_max_iterations(self, params: dict) ) def test_set_default_params_len_2_value_error_m(self, params: dict) -> None: with pytest.raises(ValueError, match="Expected positive integer as parameter m"): - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", @@ -73,7 +73,7 @@ def test_set_default_params_len_2_value_error_m(self, params: dict) -> None: ], ) def test_set_default_params_len_3_correct(self, params: dict) -> None: - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", @@ -84,17 +84,17 @@ def test_set_default_params_len_3_correct(self, params: dict) -> None: ], ) def test_set_default_params_len_4_correct(self, params: dict) -> None: - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", [{"m": 100, "tolerance": 1 / 10}, {"m": 1000, "tolerance": 10**-9}, {"m": 1, "tolerance": 1}] ) def test_set_default_params_len_2_correct(self, params: dict) -> None: - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize("params", [{"m": 100}, {"m": 1000}, {"m": 1}]) def test_set_default_params_len_1_correct(self, params: dict) -> None: - SemiParametricMuEstimation()._validate_kwargs(**params) + NMVEstimationMu()._validate_kwargs(**params) @pytest.mark.parametrize( "params", @@ -105,14 +105,14 @@ def test_set_default_params_len_1_correct(self, params: dict) -> None: ], ) def test_init_set_default_params_len_3_correct(self, params: dict) -> None: - SemiParametricMuEstimation(np.array([1]), **params) + NMVEstimationMu(np.array([1]), **params) @pytest.mark.parametrize( "params", [{"m": 100, "tolerance": 1 / 10}, {"m": 1000, "tolerance": 10**-9}, {"m": 1, "tolerance": 1}] ) def test_init_set_default_params_len_2_correct(self, params: dict) -> None: - SemiParametricMuEstimation(np.array([1]), **params) + NMVEstimationMu(np.array([1]), **params) @pytest.mark.parametrize("params", [{"m": 100}, {"m": 1000}, {"m": 1}]) def test_init_set_default_params_len_1_correct(self, params: dict) -> None: - SemiParametricMuEstimation(np.array([1]), **params) + NMVEstimationMu(np.array([1]), **params) \ No newline at end of file diff --git a/tests/algorithms/nv_algorithms/test_g_estimation.py b/tests/algorithms/nv_algorithms/test_g_estimation.py index a2068ab..6371bb4 100644 --- a/tests/algorithms/nv_algorithms/test_g_estimation.py +++ b/tests/algorithms/nv_algorithms/test_g_estimation.py @@ -4,7 +4,7 @@ import pytest from scipy.stats import expon -from src.estimators.semiparametric.nv_semiparametric_estimator import NVSemiParametricEstimator +from src.estimators.semiparametric.nv_semiparametric_estimator import NVSemiparEstim from src.generators.nv_generator import NVGenerator from src.mixtures.nv_mixture import NormalVarianceMixtures @@ -17,7 +17,7 @@ def test_g_estimation_expon(self, x_data) -> None: n = 100 mixture = NormalVarianceMixtures("canonical", alpha=0, distribution=expon) sample = NVGenerator().canonical_generate(mixture, n) - estimator = NVSemiParametricEstimator( + estimator = NVSemiparEstim( "density_estim_inv_mellin_quad_int", {"x_data": x_data, "u_value": 7.6, "v_value": 0.9} ) est = estimator.estimate(sample) diff --git a/tests/algorithms/nv_algorithms/test_post_widder_nv.py b/tests/algorithms/nv_algorithms/test_post_widder_nv.py index 02c267d..6dd917d 100644 --- a/tests/algorithms/nv_algorithms/test_post_widder_nv.py +++ b/tests/algorithms/nv_algorithms/test_post_widder_nv.py @@ -3,7 +3,7 @@ from mpmath import ln from scipy.stats import expon, gamma -from src.estimators.semiparametric.nv_semiparametric_estimator import NVSemiParametricEstimator +from src.estimators.semiparametric.nv_semiparametric_estimator import NVSemiparEstim from src.generators.nv_generator import NVGenerator from src.mixtures.nv_mixture import NormalVarianceMixtures @@ -20,7 +20,7 @@ def test_post_widder_expon(self, sigma, degree, sample_size) -> None: sample = NVGenerator().classical_generate(mixture, sample_size) x_data = np.linspace(0.5, 10.0, 30) - estimator = NVSemiParametricEstimator( + estimator = NVSemiparEstim( "density_estim_post_widder", {"x_data": x_data, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) @@ -43,7 +43,7 @@ def test_post_widder_gamma(self, sigma, degree, sample_size, a) -> None: sample = NVGenerator().classical_generate(mixture, sample_size) x_data = np.linspace(0.5, 10.0, 30) - estimator = NVSemiParametricEstimator( + estimator = NVSemiparEstim( "density_estim_post_widder", {"x_data": x_data, "sigma": sigma, "n": degree} ) est = estimator.estimate(sample) diff --git a/tests/algorithms/support_algorithms/test_rqmc.py b/tests/algorithms/support_algorithms/test_rqmc.py index b5b6558..97171c1 100644 --- a/tests/algorithms/support_algorithms/test_rqmc.py +++ b/tests/algorithms/support_algorithms/test_rqmc.py @@ -3,7 +3,7 @@ import numpy as np import pytest -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.rqmc import RQMC def loss_func(true_func: Callable, rqms: Callable, count: int): From 3c5c1b26e936ae4693a1966b700c4fa228b1a8f5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sat, 31 May 2025 17:11:04 +0300 Subject: [PATCH 09/24] feat: Add general interface StatisticalProcedure * added base StatisticalProcedure protocol * replaced NDArray[np.float64] with np.ndarray --- .../g_estimation_convolution.py | 7 +++---- .../g_estimation_given_mu.py | 6 +++--- .../g_estimation_post_widder.py | 2 +- .../g_estimation_given_mu.py | 6 +++--- .../g_estimation_given_mu_rqmc_based.py | 6 +++--- .../nvm_semi_param_algorithms/mu_estimation.py | 2 +- .../semiparametric/statistical_procedure.py | 12 ++++++++++++ 7 files changed, 26 insertions(+), 15 deletions(-) create mode 100644 src/procedures/semiparametric/statistical_procedure.py diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py index 29b0396..83de09f 100644 --- a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py +++ b/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py @@ -11,7 +11,6 @@ BOHMAN_DELTA_DEFAULT_VALUE: float = 0.0001 X_DATA_DEFAULT_VALUE: List[float] = [1.0] - class NMEstimationDensityInvMT: """Estimation of mixing density function g (xi density function) of NM mixture represented in canonical form Y = xi + sigma*N. @@ -34,8 +33,8 @@ class ParamsAnnotation(TypedDict, total=False): bohman_n: int bohman_delta: float - def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwargs: Unpack[ParamsAnnotation]): - self.sample: _typing.NDArray[np.float64] = np.array([]) if sample is None else sample + def __init__(self, sample: Optional[np.ndarray] = None, **kwargs: Unpack[ParamsAnnotation]): + self.sample: np.ndarray = np.array([]) if sample is None else sample self.n: int = len(self.sample) ( self.x_data, @@ -93,7 +92,7 @@ def cdf(self, X: np.ndarray) -> np.ndarray: return F_x.real - def compute(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: inv = self.BohmanA(N=self.bohman_n, delta=self.bohman_delta) inv.fit(self.characteristic_function_xi) x_data_array = np.array(self.x_data, dtype=np.float64) diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py index 5c87200..768d54c 100644 --- a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py +++ b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py @@ -49,13 +49,13 @@ class ParamsAnnotation(TypedDict, total=False): integration_tolerance: float integration_limit: int - def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwargs: Unpack[ParamsAnnotation]): + def __init__(self, sample: Optional[np.ndarray] = None, **kwargs: Unpack[ParamsAnnotation]): self.x_powers: Dict[float, np.ndarray] = {} self.second_u_integrals: np.ndarray self.first_u_integrals: np.ndarray self.gamma_grid: np.ndarray self.v_grid: np.ndarray - self.sample: _typing.NDArray[np.float64] = np.array([]) if sample is None else sample + self.sample: np.ndarray = np.array([]) if sample is None else sample self.n: int = len(self.sample) ( self.gmm, @@ -173,6 +173,6 @@ def compute_integrals_for_x(self, x: float) -> float: total = np.sum(first_integral + second_integral) / self.denominator return max(0.0, total.real) - def compute(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py index 3b5dc7e..fb34b9b 100644 --- a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py +++ b/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py @@ -77,7 +77,7 @@ def p_x_estimation(self, x: float) -> float: ) return result.real - def compute(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: """Estimate g(x) Args: diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py index 28cbb7d..b84363c 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py @@ -51,13 +51,13 @@ class ParamsAnnotation(TypedDict, total=False): integration_tolerance: float integration_limit: int - def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwargs: Unpack[ParamsAnnotation]): + def __init__(self, sample: Optional[np.ndarray] = None, **kwargs: Unpack[ParamsAnnotation]): self.x_powers: Dict[float, np.ndarray] = {} self.second_u_integrals: np.ndarray self.first_u_integrals: np.ndarray self.gamma_grid: np.ndarray self.v_grid: np.ndarray - self.sample: _typing.NDArray[np.float64] = np.array([]) if sample is None else sample + self.sample: np.ndarray= np.array([]) if sample is None else sample self.n: int = len(self.sample) ( self.mu, @@ -179,6 +179,6 @@ def compute_integrals_for_x(self, x: float) -> float: total = np.sum(first_integral + second_integral) / self.denominator return max(0.0, total.real) - def compute(self, sample: _typing.NDArray[np.float64]) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py index d8e1b37..0357e79 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py @@ -52,13 +52,13 @@ class ParamsAnnotation(TypedDict, total=False): integration_tolerance: float integration_limit: int - def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwargs: Unpack[ParamsAnnotation]): + def __init__(self, sample: Optional[np.ndarray] = None, **kwargs: Unpack[ParamsAnnotation]): self.x_powers: Dict[float, np.ndarray] = {} self.second_u_integrals: np.ndarray self.first_u_integrals: np.ndarray self.gamma_grid: np.ndarray self.v_grid: np.ndarray - self.sample: _typing.NDArray[np.float64] = np.array([]) if sample is None else sample + self.sample: np.ndarray = np.array([]) if sample is None else sample self.n: int = len(self.sample) ( self.mu, @@ -171,6 +171,6 @@ def compute_integrals_for_x(self, x: float) -> float: total = (first_integral + second_integral) / self.denominator return max(0.0, total.real) - def compute(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py index 48db22e..3ca6292 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py +++ b/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py @@ -84,7 +84,7 @@ def __w(self, p: float, sample: np._typing.NDArray) -> float: y += e * self.omega(x) return y - def compute(self, sample: np._typing.NDArray) -> EstimateResult: + def compute(self, sample: np.ndarray) -> EstimateResult: """Root of this function is an estimation of mu Args: diff --git a/src/procedures/semiparametric/statistical_procedure.py b/src/procedures/semiparametric/statistical_procedure.py new file mode 100644 index 0000000..de5d8ca --- /dev/null +++ b/src/procedures/semiparametric/statistical_procedure.py @@ -0,0 +1,12 @@ +from typing import Protocol +import numpy as np + +from src.estimators.estimate_result import EstimateResult + +class StatisticalProcedure(Protocol): + """Protocol for statistical algorithms that compute estimates from data. + + """ + def compute(self, sample: np.ndarray) -> EstimateResult: + pass + From 1280bdf01e1c8e06909626229c1292bb2527dd6f Mon Sep 17 00:00:00 2001 From: amellgo Date: Mon, 12 May 2025 21:43:24 +0300 Subject: [PATCH 10/24] feat: Implement dynamic integrator selection, enabling the use of custom integrators --- .../g_estimation_given_mu_rqmc_based.py | 19 ++-- .../support_algorithms/integrator.py | 17 ++++ .../support_algorithms/quad_integrator.py | 27 ++++++ src/algorithms/support_algorithms/rqmc.py | 28 +++++- src/mixtures/nm_mixture.py | 90 +++++++++++-------- src/mixtures/nmv_mixture.py | 66 ++++++++------ src/mixtures/nv_mixture.py | 26 +++--- 7 files changed, 187 insertions(+), 86 deletions(-) create mode 100644 src/algorithms/support_algorithms/integrator.py create mode 100644 src/algorithms/support_algorithms/quad_integrator.py diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py index 2f3da81..e10ba86 100644 --- a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py +++ b/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py @@ -7,7 +7,8 @@ from scipy.integrate import quad_vec from scipy.special import gamma -from src.algorithms.support_algorithms.rqmc import RQMC +from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator from src.estimators.estimate_result import EstimateResult MU_DEFAULT_VALUE = 1.0 @@ -30,6 +31,7 @@ def v_sequence_default_value(n: float) -> float: GRID_SIZE_DEFAULT_VALUE: int = 200 INTEGRATION_TOLERANCE_DEFAULT_VALUE: float = 1e-2 INTEGRATION_LIMIT_DEFAULT_VALUE: int = 50 +INTEGRATOR_DEFAULT = RQMCIntegrator() class SemiParametricGEstimationGivenMuRQMCBased: @@ -69,6 +71,7 @@ def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwarg self.grid_size, self.integration_tolerance, self.integration_limit, + self.integrator ) = self._validate_kwargs(self.n, **kwargs) self.denominator: float = 2 * math.pi * self.n self.precompute_gamma_grid() @@ -78,7 +81,7 @@ def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwarg @staticmethod def _validate_kwargs( n: int, **kwargs: Unpack[ParamsAnnotation] - ) -> tuple[float, float, float, float, List[float], int, float, int]: + ) -> tuple[float, float, float, float, List[float], int, float, int, Integrator]: mu: float = kwargs.get("mu", MU_DEFAULT_VALUE) gmm: float = kwargs.get("gmm", GAMMA_DEFAULT_VALUE) u_value: float = kwargs.get("u_value", U_SEQUENCE_DEFAULT_VALUE(n)) @@ -87,7 +90,8 @@ def _validate_kwargs( grid_size: int = kwargs.get("grid_size", GRID_SIZE_DEFAULT_VALUE) integration_tolerance: float = kwargs.get("integration_tolerance", INTEGRATION_TOLERANCE_DEFAULT_VALUE) integration_limit: int = kwargs.get("integration_limit", INTEGRATION_LIMIT_DEFAULT_VALUE) - return mu, gmm, u_value, v_value, x_data, grid_size, integration_tolerance, integration_limit + integrator: Integrator = kwargs.get("integrator_type", INTEGRATOR_DEFAULT) + return mu, gmm, u_value, v_value, x_data, grid_size, integration_tolerance, integration_limit, integrator def conjugate_psi(self, u: float) -> complex: return complex((u**2) / 2, self.mu * u) @@ -162,15 +166,16 @@ def second_v_integrand(self, v: float, x: float) -> np.ndarray: x_power = self.x_powers[x][idx] return (self.second_u_integrals[idx] * x_power) / gamma_val - def compute_integrals_for_x(self, x: float) -> float: + def compute_integrals_for_x(self, x: float, integrator: Integrator = None) -> float: """Compute integrals using RQMC for v-integration.""" - first_integral = RQMC(lambda t: np.sum(self.first_v_integrand(t * self.v_value, x)) * self.v_value).rqmc()[0] + integrator = integrator or self.integrator + first_integral = integrator.compute_integral(func=lambda t: np.sum(self.first_v_integrand(t * self.v_value, x)) * self.v_value).value - second_integral = RQMC(lambda t: np.sum(self.second_v_integrand(-t * self.v_value, x)) * self.v_value).rqmc()[0] + second_integral = integrator.compute_integral(func=lambda t: np.sum(self.second_v_integrand(-t * self.v_value, x)) * self.v_value).value total = (first_integral + second_integral) / self.denominator return max(0.0, total.real) def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: - y_data = [self.compute_integrals_for_x(x) for x in self.x_data] + y_data = [self.compute_integrals_for_x(x, self.integrator) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/algorithms/support_algorithms/integrator.py b/src/algorithms/support_algorithms/integrator.py new file mode 100644 index 0000000..6731aed --- /dev/null +++ b/src/algorithms/support_algorithms/integrator.py @@ -0,0 +1,17 @@ +from dataclasses import dataclass +from typing import Any, Protocol, Callable, Optional + + +@dataclass +class IntegrationResult: + value: float + error: float + message: Optional[dict[str, Any]] | None = None + + +class Integrator(Protocol): + + """Base class for integral calculation""" + + def compute_integral(self, func: Callable, params: dict) -> IntegrationResult: + ... diff --git a/src/algorithms/support_algorithms/quad_integrator.py b/src/algorithms/support_algorithms/quad_integrator.py new file mode 100644 index 0000000..937cb46 --- /dev/null +++ b/src/algorithms/support_algorithms/quad_integrator.py @@ -0,0 +1,27 @@ +from typing import Callable + +from scipy.integrate import quad +from src.algorithms.support_algorithms.integrator import IntegrationResult + +class QuadIntegrator: + + def compute_integral(self, func: Callable, params: dict) -> IntegrationResult: + + """ + Compute integral via quad integrator + + Args: + func: integrated function + params: Parameters of integration algorithm + + Returns: moment approximation and error tolerance + """ + + full_output_requested = params.pop('full_output', False) + quad_res = quad(func, **params) + if full_output_requested: + value, error, message = quad_res + else: + value, error = quad_res + message = None + return IntegrationResult(value, error, message) diff --git a/src/algorithms/support_algorithms/rqmc.py b/src/algorithms/support_algorithms/rqmc.py index 4b10c9d..520ba62 100644 --- a/src/algorithms/support_algorithms/rqmc.py +++ b/src/algorithms/support_algorithms/rqmc.py @@ -5,6 +5,8 @@ import scipy from numba import njit +from src.algorithms.support_algorithms.integrator import IntegrationResult + BITS = 30 """Number of bits in XOR. Should be less than 64""" NUMBA_FAST_MATH = True @@ -126,7 +128,8 @@ def _update( Returns:Updated mean of all rows - """ + + """ values = [] sum_of_new: float = 0.0 for i in range(self.count): @@ -212,9 +215,9 @@ def _xor_float(a: float, b: float) -> float: Returns: XOR float value """ - a = int(a * (2**BITS)) - b = int(b * (2**BITS)) - return np.bitwise_xor(a, b) / 2**BITS + a = int(a * (2 ** BITS)) + b = int(b * (2 ** BITS)) + return np.bitwise_xor(a, b) / 2 ** BITS def __call__(self) -> tuple[float, float]: """Interface for users @@ -223,3 +226,20 @@ def __call__(self) -> tuple[float, float]: """ return self.rqmc() + + +class RQMCIntegrator: + + def compute_integral(self, func: Callable, params: dict) -> IntegrationResult: + """ + Compute integral via RQMC integrator + + Args: + func: integrated function + params: Parameters of integration algorithm + + Returns: moment approximation and error tolerance + """ + + rqmc = RQMC(func, **params)() + return IntegrationResult(rqmc[0], rqmc[1]) diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 33c7580..0e4d811 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -2,13 +2,14 @@ from typing import Any import numpy as np -from scipy.integrate import quad from scipy.special import binom from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.rqmc import RQMC +from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator +from src.algorithms.support_algorithms.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @@ -59,156 +60,175 @@ def _params_validation(self, data_collector: Any, params: dict[str, float | rv_c raise ValueError("Gamma cant be zero") return data_class - def _classical_moment(self, n: int, params: dict) -> tuple[float, float]: + def _classical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Compute n-th moment of classical NMM Args: n (): Moment ordinal params (): Parameters of integration algorithm + integrator (): type of integrator to computing Returns: moment approximation and error tolerance """ mixture_moment = 0 error_tolerance = 0 + integrator = integrator or QuadIntegrator() for k in range(0, n + 1): for l in range(0, k + 1): - coefficient = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma**l) - mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (k - l), 0, 1, **params) - error_tolerance += (self.params.beta ** (k - l)) * mixing_moment[1] - mixture_moment += coefficient * (self.params.alpha ** (n - k)) * mixing_moment[0] * norm.moment(l) - return mixture_moment, error_tolerance + coefficient = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * ( + self.params.gamma ** l) + mixing_moment = integrator.compute_integral(lambda u: self.params.distribution.ppf(u) ** (k - l), + {"a": 0, "b": 1, **params}) + error_tolerance += (self.params.beta ** (k - l)) * mixing_moment.error + mixture_moment += coefficient * (self.params.alpha ** (n - k)) * mixing_moment.value * norm.moment(l) + return mixture_moment, error_tolerance + - def _canonical_moment(self, n: int, params: dict) -> tuple[float, float]: + def _canonical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Compute n-th moment of canonical NMM Args: n (): Moment ordinal params (): Parameters of integration algorithm + integrator (): type of integrator to computing Returns: moment approximation and error tolerance """ mixture_moment = 0 error_tolerance = 0 + integrator = integrator or QuadIntegrator() for k in range(0, n + 1): - coefficient = binom(n, n - k) * (self.params.sigma**k) - mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (n - k), 0, 1, **params) - error_tolerance += mixing_moment[1] - mixture_moment += coefficient * mixing_moment[0] * norm.moment(k) + coefficient = binom(n, n - k) * (self.params.sigma ** k) + mixing_moment = integrator.compute_integral(lambda u: self.params.distribution.ppf(u) ** (n - k), + {"a": 0, "b": 1, **params}) + error_tolerance += mixing_moment.error + mixture_moment += coefficient * mixing_moment.value * norm.moment(k) return mixture_moment, error_tolerance - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: + def compute_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Compute n-th moment of NMM Args: n (): Moment ordinal params (): Parameters of integration algorithm + integrator (): type of integrator to computing Returns: moment approximation and error tolerance """ if isinstance(self.params, _NMMClassicDataCollector): - return self._classical_moment(n, params) - return self._canonical_moment(n, params) + return self._classical_moment(n, params, integrator) + return self._canonical_moment(n, params, integrator) - def _canonical_compute_cdf(self, x: float, params: dict) -> tuple[float, float]: + def _canonical_compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Equation for canonical cdf Args: x (): point params (): parameters of RQMC algorithm + integrator (): type of integrator to computing Returns: computed cdf and error tolerance """ - rqmc = RQMC(lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), **params) - return rqmc() + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), **params) + return rqmc.value, rqmc.error - def _classical_compute_cdf(self, x: float, params: dict) -> tuple[float, float]: + def _classical_compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Equation for classic cdf Args: x (): point params (): parameters of RQMC algorithm + integrator (): type of integrator to computing Returns: computed cdf and error tolerance """ - rqmc = RQMC( + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func= lambda u: norm.cdf( (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) ), **params ) - return rqmc() + return rqmc.value, rqmc.error - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: + def compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Choose equation for cdf estimation depends on Mixture form Args: x (): point params (): parameters of RQMC algorithm + integrator (): type of integrator to computing Returns: Computed pdf and error tolerance """ if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_cdf(x, params) - return self._classical_compute_cdf(x, params) + return self._canonical_compute_cdf(x, params, integrator) + return self._classical_compute_cdf(x, params, integrator) - def _canonical_compute_pdf(self, x: float, params: dict) -> tuple[float, float]: + def _canonical_compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Equation for canonical pdf Args: x (): point params (): parameters of RQMC algorithm + integrator (): type of integrator to computing Returns: computed pdf and error tolerance """ - rqmc = RQMC( + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func= lambda u: (1 / np.abs(self.params.sigma)) * norm.pdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), **params ) - return rqmc() + return rqmc.value, rqmc.error - def _classical_compute_pdf(self, x: float, params: dict) -> tuple[float, float]: + def _classical_compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Equation for classic pdf Args: x (): point params (): parameters of RQMC algorithm + integrator (): type of integrator to computing Returns: computed pdf and error tolerance """ - rqmc = RQMC( + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func= lambda u: (1 / np.abs(self.params.gamma)) * norm.pdf( (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) ), **params ) - return rqmc() + return rqmc.value, rqmc.error - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: + def compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Choose equation for pdf estimation depends on Mixture form Args: x (): point params (): parameters of RQMC algorithm + integrator (): type of integrator to computing Returns: Computed pdf and error tolerance """ if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_pdf(x, params) - return self._classical_compute_pdf(x, params) + return self._canonical_compute_pdf(x, params, integrator) + return self._classical_compute_pdf(x, params, integrator) def _classical_compute_log_pdf(self, x: float, params: dict) -> tuple[float, float]: """ @@ -258,4 +278,4 @@ def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: """ if isinstance(self.params, _NMMCanonicalDataCollector): return self._canonical_compute_log_pdf(x, params) - return self._classical_compute_log_pdf(x, params) + return self._classical_compute_log_pdf(x, params) \ No newline at end of file diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index c090ff3..9f2a7e0 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -8,7 +8,8 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.rqmc import RQMC +from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @@ -40,13 +41,14 @@ class NormalMeanVarianceMixtures(AbstractMixtures): def __init__(self, mixture_form: str, **kwargs: Any) -> None: super().__init__(mixture_form, **kwargs) - def _classical_moment(self, n: int, params: dict) -> tuple[float, float]: + def _classical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Compute n-th moment of classical NMM Args: n (): Moment ordinal params (): Parameters of integration algorithm + integrator (): type of integrator to computing Returns: moment approximation and error tolerance @@ -67,16 +69,18 @@ def integral_func(u: float) -> float: ) return result - rqmc = RQMC(lambda u: integral_func(u), **params) - return rqmc() + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: integral_func(u), **params) + return rqmc.value, rqmc.error - def _canonical_moment(self, n: int, params: dict) -> tuple[float, float]: + def _canonical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Compute n-th moment of classical NMM Args: n (): Moment ordinal params (): Parameters of integration algorithm + integrator (): type of integrator to computing Returns: moment approximation and error tolerance @@ -96,15 +100,16 @@ def integral_func(u: float) -> float: ) return result - rqmc = RQMC(lambda u: integral_func(u), **params) - return rqmc() + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: integral_func(u), **params) + return rqmc.value, rqmc.error - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: + def compute_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_moment(n, params) - return self._canonical_moment(n, params) + return self._classical_moment(n, params, integrator) + return self._canonical_moment(n, params, integrator) - def _classical_cdf(self, x: float, params: dict) -> tuple[float, float]: + def _classical_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = lru_cache()(self.params.distribution.ppf)(u) point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - ( @@ -112,24 +117,26 @@ def _inner_func(u: float) -> float: ) return norm.cdf(point) - rqmc = RQMC(lambda u: _inner_func(u), **params) - return rqmc() + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) + return rqmc.value, rqmc.error - def _canonical_cdf(self, x: float, params: dict) -> tuple[float, float]: + def _canonical_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) point = (x - self.params.alpha) / (np.sqrt(ppf)) - (self.params.mu * np.sqrt(ppf)) return norm.cdf(point) - rqmc = RQMC(lambda u: _inner_func(u), **params) - return rqmc() + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) + return rqmc.value, rqmc.error - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: + def compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_cdf(x, params) - return self._canonical_cdf(x, params) + return self._classical_cdf(x, params, integrator) + return self._canonical_cdf(x, params, integrator) - def _classical_pdf(self, x: float, params: dict) -> tuple[float, float]: + def _classical_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) return ( @@ -140,10 +147,11 @@ def _inner_func(u: float) -> float: ) ) - rqmc = RQMC(lambda u: _inner_func(u), **params)() - return np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma**2) * rqmc[0], rqmc[1] + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) + return np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma**2) * rqmc.value, rqmc.error - def _canonical_pdf(self, x: float, params: dict) -> tuple[float, float]: + def _canonical_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) return ( @@ -152,9 +160,9 @@ def _inner_func(u: float) -> float: * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu**2 * ppf**2) / (2 * ppf)) ) - rqmc = RQMC(lambda u: _inner_func(u), **params) - res = rqmc() - return np.exp(self.params.mu * (x - self.params.alpha)) * res[0], res[1] + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) + return np.exp(self.params.mu * (x - self.params.alpha)) * rqmc.value, rqmc.error def _classical_log_pdf(self, x: float, params: dict) -> tuple[float, float]: def _inner_func(u: float) -> float: @@ -178,10 +186,10 @@ def _inner_func(u: float) -> float: rqmc = LogRQMC(lambda u: _inner_func(u), **params) return rqmc() - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: + def compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_pdf(x, params) - return self._canonical_pdf(x, params) + return self._classical_pdf(x, params, integrator) + return self._canonical_pdf(x, params, integrator) def compute_logpdf(self, x: float, params: dict) -> tuple[Any, float]: if isinstance(self.params, _NMVMClassicDataCollector): diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index 7905e6e..1be3343 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -8,7 +8,8 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.rqmc import RQMC +from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @@ -39,12 +40,13 @@ class NormalVarianceMixtures(AbstractMixtures): def __init__(self, mixture_form: str, **kwargs: Any) -> None: super().__init__(mixture_form, **kwargs) - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: + def compute_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: """ Compute n-th moment of NVM Args: n (): Moment ordinal params (): Parameters of integration algorithm + integrator (): type of integrator to computing Returns: moment approximation and error tolerance """ gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 @@ -60,16 +62,17 @@ def integrate_func(u: float) -> float: for k in range(0, n + 1) ] ) + integrator = integrator or RQMCIntegrator() + result = integrator.compute_integral(func=integrate_func, **params) + return result.value, result.error - result = RQMC(integrate_func, **params)() - return result - - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: + def compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: parametric_norm = norm(0, self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1) - rqmc = RQMC( + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func= lambda u: parametric_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))), **params ) - return rqmc() + return rqmc.value, rqmc.error @lru_cache() def _integrand_func(self, u: float, d: float, gamma: float) -> float: @@ -80,11 +83,12 @@ def _log_integrand_func(self, u: float, d: float, gamma: float | int | np.int64) ppf = self.params.distribution.ppf(u) return -(ppf * np.log(np.pi * 2 * ppf * gamma**2) + d) / (2 * ppf) - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: + def compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 d = (x - self.params.alpha) ** 2 / gamma**2 - rqmc = RQMC(lambda u: self._integrand_func(u, d, gamma), **params) - return rqmc() + integrator = integrator or RQMCIntegrator() + rqmc = integrator.compute_integral(func=lambda u: self._integrand_func(u, d, gamma), **params) + return rqmc.value, rqmc.error def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 From 51a7303c04ef071f415132105d8e3176d1329d4d Mon Sep 17 00:00:00 2001 From: amellgo Date: Sun, 1 Jun 2025 23:10:54 +0300 Subject: [PATCH 11/24] feat: integration parameters are now part of the constructor --- .../g_estimation_given_mu_rqmc_based.py | 18 ++-- .../support_algorithms/integrator.py | 5 +- .../support_algorithms/quad_integrator.py | 49 +++++++-- src/algorithms/support_algorithms/rqmc.py | 39 ++++++- src/mixtures/abstract_mixture.py | 7 +- src/mixtures/nm_mixture.py | 101 +++++++----------- src/mixtures/nmv_mixture.py | 74 ++++++------- src/mixtures/nv_mixture.py | 24 ++--- 8 files changed, 174 insertions(+), 143 deletions(-) diff --git a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py index e10ba86..fc089cf 100644 --- a/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py +++ b/src/algorithms/semiparam_algorithms/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py @@ -31,7 +31,6 @@ def v_sequence_default_value(n: float) -> float: GRID_SIZE_DEFAULT_VALUE: int = 200 INTEGRATION_TOLERANCE_DEFAULT_VALUE: float = 1e-2 INTEGRATION_LIMIT_DEFAULT_VALUE: int = 50 -INTEGRATOR_DEFAULT = RQMCIntegrator() class SemiParametricGEstimationGivenMuRQMCBased: @@ -70,8 +69,7 @@ def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwarg self.x_data, self.grid_size, self.integration_tolerance, - self.integration_limit, - self.integrator + self.integration_limit ) = self._validate_kwargs(self.n, **kwargs) self.denominator: float = 2 * math.pi * self.n self.precompute_gamma_grid() @@ -81,7 +79,7 @@ def __init__(self, sample: Optional[_typing.NDArray[np.float64]] = None, **kwarg @staticmethod def _validate_kwargs( n: int, **kwargs: Unpack[ParamsAnnotation] - ) -> tuple[float, float, float, float, List[float], int, float, int, Integrator]: + ) -> tuple[float, float, float, float, List[float], int, float, int]: mu: float = kwargs.get("mu", MU_DEFAULT_VALUE) gmm: float = kwargs.get("gmm", GAMMA_DEFAULT_VALUE) u_value: float = kwargs.get("u_value", U_SEQUENCE_DEFAULT_VALUE(n)) @@ -90,8 +88,7 @@ def _validate_kwargs( grid_size: int = kwargs.get("grid_size", GRID_SIZE_DEFAULT_VALUE) integration_tolerance: float = kwargs.get("integration_tolerance", INTEGRATION_TOLERANCE_DEFAULT_VALUE) integration_limit: int = kwargs.get("integration_limit", INTEGRATION_LIMIT_DEFAULT_VALUE) - integrator: Integrator = kwargs.get("integrator_type", INTEGRATOR_DEFAULT) - return mu, gmm, u_value, v_value, x_data, grid_size, integration_tolerance, integration_limit, integrator + return mu, gmm, u_value, v_value, x_data, grid_size, integration_tolerance, integration_limit def conjugate_psi(self, u: float) -> complex: return complex((u**2) / 2, self.mu * u) @@ -166,16 +163,15 @@ def second_v_integrand(self, v: float, x: float) -> np.ndarray: x_power = self.x_powers[x][idx] return (self.second_u_integrals[idx] * x_power) / gamma_val - def compute_integrals_for_x(self, x: float, integrator: Integrator = None) -> float: + def compute_integrals_for_x(self, x: float, integrator: Integrator = RQMCIntegrator()) -> float: """Compute integrals using RQMC for v-integration.""" - integrator = integrator or self.integrator - first_integral = integrator.compute_integral(func=lambda t: np.sum(self.first_v_integrand(t * self.v_value, x)) * self.v_value).value + first_integral = integrator.compute(func=lambda t: np.sum(self.first_v_integrand(t * self.v_value, x)) * self.v_value).value - second_integral = integrator.compute_integral(func=lambda t: np.sum(self.second_v_integrand(-t * self.v_value, x)) * self.v_value).value + second_integral = integrator.compute(func=lambda t: np.sum(self.second_v_integrand(-t * self.v_value, x)) * self.v_value).value total = (first_integral + second_integral) / self.denominator return max(0.0, total.real) def algorithm(self, sample: np._typing.NDArray) -> EstimateResult: - y_data = [self.compute_integrals_for_x(x, self.integrator) for x in self.x_data] + y_data = [self.compute_integrals_for_x(x, RQMCIntegrator()) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) diff --git a/src/algorithms/support_algorithms/integrator.py b/src/algorithms/support_algorithms/integrator.py index 6731aed..8ee0fb7 100644 --- a/src/algorithms/support_algorithms/integrator.py +++ b/src/algorithms/support_algorithms/integrator.py @@ -13,5 +13,8 @@ class Integrator(Protocol): """Base class for integral calculation""" - def compute_integral(self, func: Callable, params: dict) -> IntegrationResult: + def __init__(self) -> None: + ... + + def compute(self, func: Callable) -> IntegrationResult: ... diff --git a/src/algorithms/support_algorithms/quad_integrator.py b/src/algorithms/support_algorithms/quad_integrator.py index 937cb46..616eca2 100644 --- a/src/algorithms/support_algorithms/quad_integrator.py +++ b/src/algorithms/support_algorithms/quad_integrator.py @@ -1,27 +1,60 @@ -from typing import Callable +from typing import Callable, Any, Sequence from scipy.integrate import quad from src.algorithms.support_algorithms.integrator import IntegrationResult class QuadIntegrator: - def compute_integral(self, func: Callable, params: dict) -> IntegrationResult: + def __init__( + self, + a: float = 0, + b: float = 1, + args: tuple[Any, ...] = (), + full_output: int = 0, + epsabs: float | int = 1.49e-08, + epsrel: float | int = 1.49e-08, + limit: float | int = 50, + points: Sequence[float | int] | None = None, + weight: float | int | None = None, + wvar: Any = None, + wopts: Any = None, + maxp1: float | int = 50, + limlst: int = 50, + complex_func: bool = False, + ): + self.params = { + 'a': a, + 'b': b, + 'args': args, + 'full_output': full_output, + 'epsabs': epsabs, + 'epsrel': epsrel, + 'limit': limit, + 'points': points, + 'weight': weight, + 'wvar': wvar, + 'wopts': wopts, + 'maxp1': maxp1, + 'limlst': limlst, + 'complex_func': complex_func + } + + def compute(self, func: Callable) -> IntegrationResult: """ Compute integral via quad integrator Args: func: integrated function - params: Parameters of integration algorithm Returns: moment approximation and error tolerance """ - full_output_requested = params.pop('full_output', False) - quad_res = quad(func, **params) - if full_output_requested: - value, error, message = quad_res + verbose = self.params.pop('full_output', False) + result = quad(func, **self.params) + if verbose: + value, error, message = result else: - value, error = quad_res + value, error = result message = None return IntegrationResult(value, error, message) diff --git a/src/algorithms/support_algorithms/rqmc.py b/src/algorithms/support_algorithms/rqmc.py index 520ba62..543d2cb 100644 --- a/src/algorithms/support_algorithms/rqmc.py +++ b/src/algorithms/support_algorithms/rqmc.py @@ -229,17 +229,48 @@ def __call__(self) -> tuple[float, float]: class RQMCIntegrator: + """ + Randomize Quasi Monte Carlo Method + + Args: + error_tolerance: pre-specified error tolerance + count: number of rows of random values matrix + base_n: number of columns of random values matrix + i_max: allowed number of cycles + a: parameter for quantile of normal distribution + + """ + + def __init__( + self, + error_tolerance: float = 1e-6, + count: int = 25, + base_n: int = 2 ** 6, + i_max: int = 100, + a: float = 0.00047, + ): + self.error_tolerance = error_tolerance + self.count = count + self.base_n = base_n + self.i_max = i_max + self.a = a - def compute_integral(self, func: Callable, params: dict) -> IntegrationResult: + def compute(self, func: Callable) -> IntegrationResult: """ Compute integral via RQMC integrator Args: func: integrated function - params: Parameters of integration algorithm Returns: moment approximation and error tolerance """ + result = RQMC( + func, + error_tolerance=self.error_tolerance, + count=self.count, + base_n=self.base_n, + i_max=self.i_max, + a=self.a, + )() + return IntegrationResult(result[0], result[1]) - rqmc = RQMC(func, **params)() - return IntegrationResult(rqmc[0], rqmc[1]) diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index 13c8e42..3e95de4 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -5,6 +5,7 @@ from scipy.stats import rv_continuous from scipy.stats.distributions import rv_frozen +from src.algorithms.support_algorithms.integrator import Integrator class AbstractMixtures(metaclass=ABCMeta): """Base class for Mixtures""" @@ -28,13 +29,13 @@ def __init__(self, mixture_form: str, **kwargs: Any) -> None: raise AssertionError(f"Unknown mixture form: {mixture_form}") @abstractmethod - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: ... + def compute_moment(self, n: int, integrator: Integrator) -> tuple[float, float]: ... @abstractmethod - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: ... + def compute_cdf(self, x: float, integrator: Integrator) -> tuple[float, float]: ... @abstractmethod - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: ... + def compute_pdf(self, x: float, integrator: Integrator) -> tuple[float, float]: ... @abstractmethod def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: ... diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 0e4d811..617c498 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -60,175 +60,154 @@ def _params_validation(self, data_collector: Any, params: dict[str, float | rv_c raise ValueError("Gamma cant be zero") return data_class - def _classical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _classical_moment(self, n: int, integrator: Integrator = QuadIntegrator()) -> tuple[float, float]: """ Compute n-th moment of classical NMM Args: n (): Moment ordinal - params (): Parameters of integration algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: moment approximation and error tolerance """ - mixture_moment = 0 - error_tolerance = 0 - integrator = integrator or QuadIntegrator() + mixture_moment = 0.0 + error_tolerance = 0.0 for k in range(0, n + 1): for l in range(0, k + 1): coefficient = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * ( self.params.gamma ** l) - mixing_moment = integrator.compute_integral(lambda u: self.params.distribution.ppf(u) ** (k - l), - {"a": 0, "b": 1, **params}) + mixing_moment = integrator.compute(lambda u: self.params.distribution.ppf(u) ** (k - l)) error_tolerance += (self.params.beta ** (k - l)) * mixing_moment.error mixture_moment += coefficient * (self.params.alpha ** (n - k)) * mixing_moment.value * norm.moment(l) - return mixture_moment, error_tolerance - + return mixture_moment, error_tolerance - def _canonical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _canonical_moment(self, n: int, integrator: Integrator = QuadIntegrator()) -> tuple[float, float]: """ Compute n-th moment of canonical NMM Args: n (): Moment ordinal - params (): Parameters of integration algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: moment approximation and error tolerance """ - mixture_moment = 0 - error_tolerance = 0 - integrator = integrator or QuadIntegrator() + mixture_moment = 0.0 + error_tolerance = 0.0 for k in range(0, n + 1): coefficient = binom(n, n - k) * (self.params.sigma ** k) - mixing_moment = integrator.compute_integral(lambda u: self.params.distribution.ppf(u) ** (n - k), - {"a": 0, "b": 1, **params}) + mixing_moment = integrator.compute(lambda u: self.params.distribution.ppf(u) ** (n - k)) error_tolerance += mixing_moment.error mixture_moment += coefficient * mixing_moment.value * norm.moment(k) return mixture_moment, error_tolerance - def compute_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_moment(self, n: int, integrator: Integrator = QuadIntegrator()) -> tuple[float, float]: """ Compute n-th moment of NMM Args: n (): Moment ordinal - params (): Parameters of integration algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: moment approximation and error tolerance """ if isinstance(self.params, _NMMClassicDataCollector): - return self._classical_moment(n, params, integrator) - return self._canonical_moment(n, params, integrator) + return self._classical_moment(n, integrator) + return self._canonical_moment(n, integrator) - def _canonical_compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _canonical_compute_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Equation for canonical cdf Args: x (): point - params (): parameters of RQMC algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: computed cdf and error tolerance """ - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), **params) - return rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma))) + return result.value, result.error - def _classical_compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _classical_compute_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Equation for classic cdf Args: x (): point - params (): parameters of RQMC algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: computed cdf and error tolerance """ - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func= + result = integrator.compute(func= lambda u: norm.cdf( (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) - ), - **params + ) ) - return rqmc.value, rqmc.error + return result.value, result.error - def compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Choose equation for cdf estimation depends on Mixture form Args: x (): point - params (): parameters of RQMC algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: Computed pdf and error tolerance """ if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_cdf(x, params, integrator) - return self._classical_compute_cdf(x, params, integrator) + return self._canonical_compute_cdf(x, integrator) + return self._classical_compute_cdf(x, integrator) - def _canonical_compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _canonical_compute_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Equation for canonical pdf Args: x (): point - params (): parameters of RQMC algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: computed pdf and error tolerance """ - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func= + result = integrator.compute(func= lambda u: (1 / np.abs(self.params.sigma)) * norm.pdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), - **params ) - return rqmc.value, rqmc.error + return result.value, result.error - def _classical_compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _classical_compute_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Equation for classic pdf Args: x (): point - params (): parameters of RQMC algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: computed pdf and error tolerance """ - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func= + result = integrator.compute(func= lambda u: (1 / np.abs(self.params.gamma)) * norm.pdf( (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) ), - **params ) - return rqmc.value, rqmc.error + return result.value, result.error - def compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Choose equation for pdf estimation depends on Mixture form Args: x (): point - params (): parameters of RQMC algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: Computed pdf and error tolerance """ if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_pdf(x, params, integrator) - return self._classical_compute_pdf(x, params, integrator) + return self._canonical_compute_pdf(x, integrator) + return self._classical_compute_pdf(x, integrator) def _classical_compute_log_pdf(self, x: float, params: dict) -> tuple[float, float]: """ diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index 9f2a7e0..a147cba 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -41,14 +41,13 @@ class NormalMeanVarianceMixtures(AbstractMixtures): def __init__(self, mixture_form: str, **kwargs: Any) -> None: super().__init__(mixture_form, **kwargs) - def _classical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _classical_moment(self, n: int, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Compute n-th moment of classical NMM Args: n (): Moment ordinal - params (): Parameters of integration algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: moment approximation and error tolerance @@ -69,18 +68,16 @@ def integral_func(u: float) -> float: ) return result - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: integral_func(u), **params) - return rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: integral_func(u)) + return result.value, result.error - def _canonical_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _canonical_moment(self, n: int, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Compute n-th moment of classical NMM Args: n (): Moment ordinal - params (): Parameters of integration algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: moment approximation and error tolerance @@ -100,16 +97,15 @@ def integral_func(u: float) -> float: ) return result - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: integral_func(u), **params) - return rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: integral_func(u)) + return result.value, result.error - def compute_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_moment(self, n: int, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_moment(n, params, integrator) - return self._canonical_moment(n, params, integrator) + return self._classical_moment(n, integrator) + return self._canonical_moment(n, integrator) - def _classical_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _classical_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = lru_cache()(self.params.distribution.ppf)(u) point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - ( @@ -117,26 +113,24 @@ def _inner_func(u: float) -> float: ) return norm.cdf(point) - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) - return rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: _inner_func(u)) + return result.value, result.error - def _canonical_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _canonical_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) point = (x - self.params.alpha) / (np.sqrt(ppf)) - (self.params.mu * np.sqrt(ppf)) return norm.cdf(point) - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) - return rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: _inner_func(u)) + return result.value, result.error - def compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_cdf(x, params, integrator) - return self._canonical_cdf(x, params, integrator) + return self._classical_cdf(x, integrator) + return self._canonical_cdf(x, integrator) - def _classical_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _classical_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) return ( @@ -147,11 +141,10 @@ def _inner_func(u: float) -> float: ) ) - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) - return np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma**2) * rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: _inner_func(u)) + return np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma**2) * result.value, result.error - def _canonical_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def _canonical_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: def _inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) return ( @@ -160,9 +153,8 @@ def _inner_func(u: float) -> float: * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu**2 * ppf**2) / (2 * ppf)) ) - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: _inner_func(u), **params) - return np.exp(self.params.mu * (x - self.params.alpha)) * rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: _inner_func(u)) + return np.exp(self.params.mu * (x - self.params.alpha)) * result.value, result.error def _classical_log_pdf(self, x: float, params: dict) -> tuple[float, float]: def _inner_func(u: float) -> float: @@ -173,8 +165,8 @@ def _inner_func(u: float) -> float: + ppf * self.params.gamma**2 * np.log(2 * np.pi * ppf * self.params.gamma**2) ) / (2 * ppf * self.params.gamma**2) - rqmc = LogRQMC(lambda u: _inner_func(u), **params) - return rqmc() + result = LogRQMC(lambda u: _inner_func(u), **params) + return result() def _canonical_log_pdf(self, x: float, params: dict) -> tuple[float, float]: def _inner_func(u: float) -> float: @@ -183,13 +175,13 @@ def _inner_func(u: float) -> float: 2 * ppf ) - rqmc = LogRQMC(lambda u: _inner_func(u), **params) - return rqmc() + result = LogRQMC(lambda u: _inner_func(u), **params) + return result() - def compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_pdf(x, params, integrator) - return self._canonical_pdf(x, params, integrator) + return self._classical_pdf(x, integrator) + return self._canonical_pdf(x, integrator) def compute_logpdf(self, x: float, params: dict) -> tuple[Any, float]: if isinstance(self.params, _NMVMClassicDataCollector): diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index 1be3343..ffc4f08 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -40,13 +40,12 @@ class NormalVarianceMixtures(AbstractMixtures): def __init__(self, mixture_form: str, **kwargs: Any) -> None: super().__init__(mixture_form, **kwargs) - def compute_moment(self, n: int, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_moment(self, n: int, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: """ Compute n-th moment of NVM Args: n (): Moment ordinal - params (): Parameters of integration algorithm - integrator (): type of integrator to computing + integrator (): class of integrator with params to computing Returns: moment approximation and error tolerance """ gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 @@ -62,17 +61,15 @@ def integrate_func(u: float) -> float: for k in range(0, n + 1) ] ) - integrator = integrator or RQMCIntegrator() - result = integrator.compute_integral(func=integrate_func, **params) + result = integrator.compute(func=integrate_func) return result.value, result.error - def compute_cdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_cdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: parametric_norm = norm(0, self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1) - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func= - lambda u: parametric_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))), **params + result = integrator.compute(func= + lambda u: parametric_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))) ) - return rqmc.value, rqmc.error + return result.value, result.error @lru_cache() def _integrand_func(self, u: float, d: float, gamma: float) -> float: @@ -83,12 +80,11 @@ def _log_integrand_func(self, u: float, d: float, gamma: float | int | np.int64) ppf = self.params.distribution.ppf(u) return -(ppf * np.log(np.pi * 2 * ppf * gamma**2) + d) / (2 * ppf) - def compute_pdf(self, x: float, params: dict, integrator: Integrator = None) -> tuple[float, float]: + def compute_pdf(self, x: float, integrator: Integrator = RQMCIntegrator()) -> tuple[float, float]: gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 d = (x - self.params.alpha) ** 2 / gamma**2 - integrator = integrator or RQMCIntegrator() - rqmc = integrator.compute_integral(func=lambda u: self._integrand_func(u, d, gamma), **params) - return rqmc.value, rqmc.error + result = integrator.compute(func=lambda u: self._integrand_func(u, d, gamma)) + return result.value, result.error def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 From 947f9566de3363027d0274fbe18554643300945a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Fri, 6 Jun 2025 01:09:20 +0300 Subject: [PATCH 12/24] refactor: Rename algorithm classes and folders MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit * Rename NVSemiParametricGEstimationPostWidder → NVEstimationDensityPW * Rename SemiParametricGEstimationPostWidder → NMVEstimationDensityPW * Rename folders in src/procedures * Rename tests/algorithms → tests/procedures lgorithms/semiparametric_sigma_estimation/__init__.py -> tests/procedures/nm_procedures/semiparametric_sigma_estimation/__init__.py --- src/estimators/abstract_estimator.py | 4 ++-- src/procedures/__init__.py | 22 +++++++++---------- .../__init__.py | 0 .../__init__.py | 0 .../__init__.py | 0 .../__init__.py | 0 .../g_estimation_convolution.py | 0 .../sigma_estimation_eigenvalue_based.py | 0 .../sigma_estimation_empirical.py | 0 .../__init__.py | 0 .../g_estimation_given_mu.py | 0 .../g_estimation_post_widder.py | 0 .../__init__.py | 0 .../g_estimation_given_mu.py | 0 .../g_estimation_given_mu_rqmc_based.py | 0 .../g_estimation_post_widder.py | 2 +- .../mu_estimation.py | 0 .../semiparametric/statistical_procedure.py | 2 +- tests/{algorithms => procedures}/__init__.py | 0 .../__init__.py | 0 .../test_sigma_estimation_eigenvalue_based.py | 0 .../test_sigma_estimation_empirical.py | 0 .../semiparametric_g_estimation/__init__.py | 0 .../test_g_estimation_given_mu.py | 0 .../test_post_widder.py | 0 .../semiparametric_mu_estimation/__init__.py | 0 .../test_semiparametric_mu_estimation.py | 0 .../test_validate_kwargs.py | 2 +- .../nv_procedures}/__init__.py | 0 .../nv_procedures}/test_g_estimation.py | 0 .../nv_procedures}/test_post_widder_nv.py | 0 .../support_procedures}/test_rqmc.py | 0 32 files changed, 16 insertions(+), 16 deletions(-) rename src/procedures/parametric/{nm_param_algorithms => nm_parametric}/__init__.py (100%) rename src/procedures/parametric/{nv_param_algorithms => nv_parametric}/__init__.py (100%) rename src/procedures/parametric/{nvm_param_algorithms => nvm_parametric}/__init__.py (100%) rename src/procedures/semiparametric/{nm_semi_param_algorithms => nm_semiparametric}/__init__.py (100%) rename src/procedures/semiparametric/{nm_semi_param_algorithms => nm_semiparametric}/g_estimation_convolution.py (100%) rename src/procedures/semiparametric/{nm_semi_param_algorithms => nm_semiparametric}/sigma_estimation_eigenvalue_based.py (100%) rename src/procedures/semiparametric/{nm_semi_param_algorithms => nm_semiparametric}/sigma_estimation_empirical.py (100%) rename src/procedures/semiparametric/{nv_semi_param_algorithms => nv_semiparametric}/__init__.py (100%) rename src/procedures/semiparametric/{nv_semi_param_algorithms => nv_semiparametric}/g_estimation_given_mu.py (100%) rename src/procedures/semiparametric/{nv_semi_param_algorithms => nv_semiparametric}/g_estimation_post_widder.py (100%) rename src/procedures/semiparametric/{nvm_semi_param_algorithms => nvm_semiparametric}/__init__.py (100%) rename src/procedures/semiparametric/{nvm_semi_param_algorithms => nvm_semiparametric}/g_estimation_given_mu.py (100%) rename src/procedures/semiparametric/{nvm_semi_param_algorithms => nvm_semiparametric}/g_estimation_given_mu_rqmc_based.py (100%) rename src/procedures/semiparametric/{nvm_semi_param_algorithms => nvm_semiparametric}/g_estimation_post_widder.py (98%) rename src/procedures/semiparametric/{nvm_semi_param_algorithms => nvm_semiparametric}/mu_estimation.py (100%) rename tests/{algorithms => procedures}/__init__.py (100%) rename tests/{algorithms/nm_algorithms => procedures/nm_procedures}/semiparametric_sigma_estimation/__init__.py (100%) rename tests/{algorithms/nm_algorithms => procedures/nm_procedures}/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py (100%) rename tests/{algorithms/nm_algorithms => procedures/nm_procedures}/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py (100%) rename tests/{algorithms/nmv_algorithms => procedures/nmv_procedures}/semiparametric_g_estimation/__init__.py (100%) rename tests/{algorithms/nmv_algorithms => procedures/nmv_procedures}/semiparametric_g_estimation/test_g_estimation_given_mu.py (100%) rename tests/{algorithms/nmv_algorithms => procedures/nmv_procedures}/semiparametric_g_estimation/test_post_widder.py (100%) rename tests/{algorithms/nmv_algorithms => procedures/nmv_procedures}/semiparametric_mu_estimation/__init__.py (100%) rename tests/{algorithms/nmv_algorithms => procedures/nmv_procedures}/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py (100%) rename tests/{algorithms/nmv_algorithms => procedures/nmv_procedures}/semiparametric_mu_estimation/test_validate_kwargs.py (97%) rename tests/{algorithms/nv_algorithms => procedures/nv_procedures}/__init__.py (100%) rename tests/{algorithms/nv_algorithms => procedures/nv_procedures}/test_g_estimation.py (100%) rename tests/{algorithms/nv_algorithms => procedures/nv_procedures}/test_post_widder_nv.py (100%) rename tests/{algorithms/support_algorithms => procedures/support_procedures}/test_rqmc.py (100%) diff --git a/src/estimators/abstract_estimator.py b/src/estimators/abstract_estimator.py index e3fe6d8..6bee8d3 100644 --- a/src/estimators/abstract_estimator.py +++ b/src/estimators/abstract_estimator.py @@ -14,7 +14,7 @@ class AbstractEstimator: algorithm_name: A string indicating chosen algorithm. params: A dictionary of algorithm parameters. estimate_result: Estimation result. - _registry: Registry that contains classes of all algorithms. + _registry: Registry that contains classes of all procedures. _purpose: Defines purpose of algorithm, one of the registry key. """ @@ -47,7 +47,7 @@ def set_params(self, algorithm_name: str, params: dict | None = None) -> Self: return self def get_available_algorithms(self) -> list[str]: - """Get all algorithms that can be used for current estimator class""" + """Get all procedures that can be used for current estimator class""" return [key[0] for key in self._registry.register_of_names.keys() if key[1] == self._purpose] def estimate(self, sample: _typing.ArrayLike) -> EstimateResult: diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index 1c88c3b..59d7bd6 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -1,28 +1,28 @@ -from src.procedures.semiparametric.nm_semi_param_algorithms.g_estimation_convolution import ( +from src.procedures.semiparametric.nm_semiparametric.g_estimation_convolution import ( NMEstimationDensityInvMT, ) -from src.procedures.semiparametric.nm_semi_param_algorithms.sigma_estimation_eigenvalue_based import ( +from src.procedures.semiparametric.nm_semiparametric.sigma_estimation_eigenvalue_based import ( NMEstimationSigmaEigenvals, ) -from src.procedures.semiparametric.nm_semi_param_algorithms.sigma_estimation_empirical import ( +from src.procedures.semiparametric.nm_semiparametric.sigma_estimation_empirical import ( NMEstimationSigmaLaplace, ) -from src.procedures.semiparametric.nv_semi_param_algorithms.g_estimation_given_mu import ( +from src.procedures.semiparametric.nv_semiparametric.g_estimation_given_mu import ( NVEstimationDensityInvMT, ) -from src.procedures.semiparametric.nv_semi_param_algorithms.g_estimation_post_widder import ( +from src.procedures.semiparametric.nv_semiparametric.g_estimation_post_widder import ( NMVEstimationDensityPW, ) -from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_given_mu import ( +from src.procedures.semiparametric.nvm_semiparametric.g_estimation_given_mu import ( NMVEstimationDensityInvMTquad, ) -from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_given_mu_rqmc_based import ( +from src.procedures.semiparametric.nvm_semiparametric.g_estimation_given_mu_rqmc_based import ( NMVEstimationDensityInvMTquadRQMCBased, ) -from src.procedures.semiparametric.nvm_semi_param_algorithms.g_estimation_post_widder import ( - SemiParametricGEstimationPostWidder, +from src.procedures.semiparametric.nvm_semiparametric.g_estimation_post_widder import ( + NMVEstimationDensityPW, ) -from src.procedures.semiparametric.nvm_semi_param_algorithms.mu_estimation import (NMVEstimationMu) +from src.procedures.semiparametric.nvm_semiparametric.mu_estimation import (NMVEstimationMu) from src.register.algorithm_purpose import AlgorithmPurpose from src.register.register import Registry @@ -42,5 +42,5 @@ NMVEstimationDensityInvMTquadRQMCBased ) ALGORITHM_REGISTRY.register("density_estim_post_widder", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( - SemiParametricGEstimationPostWidder + NMVEstimationDensityPW ) diff --git a/src/procedures/parametric/nm_param_algorithms/__init__.py b/src/procedures/parametric/nm_parametric/__init__.py similarity index 100% rename from src/procedures/parametric/nm_param_algorithms/__init__.py rename to src/procedures/parametric/nm_parametric/__init__.py diff --git a/src/procedures/parametric/nv_param_algorithms/__init__.py b/src/procedures/parametric/nv_parametric/__init__.py similarity index 100% rename from src/procedures/parametric/nv_param_algorithms/__init__.py rename to src/procedures/parametric/nv_parametric/__init__.py diff --git a/src/procedures/parametric/nvm_param_algorithms/__init__.py b/src/procedures/parametric/nvm_parametric/__init__.py similarity index 100% rename from src/procedures/parametric/nvm_param_algorithms/__init__.py rename to src/procedures/parametric/nvm_parametric/__init__.py diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/__init__.py b/src/procedures/semiparametric/nm_semiparametric/__init__.py similarity index 100% rename from src/procedures/semiparametric/nm_semi_param_algorithms/__init__.py rename to src/procedures/semiparametric/nm_semiparametric/__init__.py diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py b/src/procedures/semiparametric/nm_semiparametric/g_estimation_convolution.py similarity index 100% rename from src/procedures/semiparametric/nm_semi_param_algorithms/g_estimation_convolution.py rename to src/procedures/semiparametric/nm_semiparametric/g_estimation_convolution.py diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py b/src/procedures/semiparametric/nm_semiparametric/sigma_estimation_eigenvalue_based.py similarity index 100% rename from src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_eigenvalue_based.py rename to src/procedures/semiparametric/nm_semiparametric/sigma_estimation_eigenvalue_based.py diff --git a/src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py b/src/procedures/semiparametric/nm_semiparametric/sigma_estimation_empirical.py similarity index 100% rename from src/procedures/semiparametric/nm_semi_param_algorithms/sigma_estimation_empirical.py rename to src/procedures/semiparametric/nm_semiparametric/sigma_estimation_empirical.py diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/__init__.py b/src/procedures/semiparametric/nv_semiparametric/__init__.py similarity index 100% rename from src/procedures/semiparametric/nv_semi_param_algorithms/__init__.py rename to src/procedures/semiparametric/nv_semiparametric/__init__.py diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nv_semiparametric/g_estimation_given_mu.py similarity index 100% rename from src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_given_mu.py rename to src/procedures/semiparametric/nv_semiparametric/g_estimation_given_mu.py diff --git a/src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nv_semiparametric/g_estimation_post_widder.py similarity index 100% rename from src/procedures/semiparametric/nv_semi_param_algorithms/g_estimation_post_widder.py rename to src/procedures/semiparametric/nv_semiparametric/g_estimation_post_widder.py diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/__init__.py b/src/procedures/semiparametric/nvm_semiparametric/__init__.py similarity index 100% rename from src/procedures/semiparametric/nvm_semi_param_algorithms/__init__.py rename to src/procedures/semiparametric/nvm_semiparametric/__init__.py diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py b/src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu.py similarity index 100% rename from src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu.py rename to src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu.py diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py b/src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu_rqmc_based.py similarity index 100% rename from src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_given_mu_rqmc_based.py rename to src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu_rqmc_based.py diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py b/src/procedures/semiparametric/nvm_semiparametric/g_estimation_post_widder.py similarity index 98% rename from src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py rename to src/procedures/semiparametric/nvm_semiparametric/g_estimation_post_widder.py index df994fb..2d219b0 100644 --- a/src/procedures/semiparametric/nvm_semi_param_algorithms/g_estimation_post_widder.py +++ b/src/procedures/semiparametric/nvm_semiparametric/g_estimation_post_widder.py @@ -14,7 +14,7 @@ X_DATA_DEFAULT_VALUE = [1.0] -class SemiParametricGEstimationPostWidder: +class NMVEstimationDensityPW: """Estimation of mixing density function g (xi density function) of NVM mixture represented in classical form Y = alpha + mu*xi + sigma*sqrt(xi)*N, where alpha = 0 and mu, sigma are given. diff --git a/src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py b/src/procedures/semiparametric/nvm_semiparametric/mu_estimation.py similarity index 100% rename from src/procedures/semiparametric/nvm_semi_param_algorithms/mu_estimation.py rename to src/procedures/semiparametric/nvm_semiparametric/mu_estimation.py diff --git a/src/procedures/semiparametric/statistical_procedure.py b/src/procedures/semiparametric/statistical_procedure.py index de5d8ca..228c82c 100644 --- a/src/procedures/semiparametric/statistical_procedure.py +++ b/src/procedures/semiparametric/statistical_procedure.py @@ -4,7 +4,7 @@ from src.estimators.estimate_result import EstimateResult class StatisticalProcedure(Protocol): - """Protocol for statistical algorithms that compute estimates from data. + """Protocol for statistical procedures that compute estimates from data. """ def compute(self, sample: np.ndarray) -> EstimateResult: diff --git a/tests/algorithms/__init__.py b/tests/procedures/__init__.py similarity index 100% rename from tests/algorithms/__init__.py rename to tests/procedures/__init__.py diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/__init__.py b/tests/procedures/nm_procedures/semiparametric_sigma_estimation/__init__.py similarity index 100% rename from tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/__init__.py rename to tests/procedures/nm_procedures/semiparametric_sigma_estimation/__init__.py diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py b/tests/procedures/nm_procedures/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py similarity index 100% rename from tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py rename to tests/procedures/nm_procedures/semiparametric_sigma_estimation/test_sigma_estimation_eigenvalue_based.py diff --git a/tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py b/tests/procedures/nm_procedures/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py similarity index 100% rename from tests/algorithms/nm_algorithms/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py rename to tests/procedures/nm_procedures/semiparametric_sigma_estimation/test_sigma_estimation_empirical.py diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/__init__.py b/tests/procedures/nmv_procedures/semiparametric_g_estimation/__init__.py similarity index 100% rename from tests/algorithms/nmv_algorithms/semiparametric_g_estimation/__init__.py rename to tests/procedures/nmv_procedures/semiparametric_g_estimation/__init__.py diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py b/tests/procedures/nmv_procedures/semiparametric_g_estimation/test_g_estimation_given_mu.py similarity index 100% rename from tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_g_estimation_given_mu.py rename to tests/procedures/nmv_procedures/semiparametric_g_estimation/test_g_estimation_given_mu.py diff --git a/tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py b/tests/procedures/nmv_procedures/semiparametric_g_estimation/test_post_widder.py similarity index 100% rename from tests/algorithms/nmv_algorithms/semiparametric_g_estimation/test_post_widder.py rename to tests/procedures/nmv_procedures/semiparametric_g_estimation/test_post_widder.py diff --git a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/__init__.py b/tests/procedures/nmv_procedures/semiparametric_mu_estimation/__init__.py similarity index 100% rename from tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/__init__.py rename to tests/procedures/nmv_procedures/semiparametric_mu_estimation/__init__.py diff --git a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py b/tests/procedures/nmv_procedures/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py similarity index 100% rename from tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py rename to tests/procedures/nmv_procedures/semiparametric_mu_estimation/test_semiparametric_mu_estimation.py diff --git a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py b/tests/procedures/nmv_procedures/semiparametric_mu_estimation/test_validate_kwargs.py similarity index 97% rename from tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py rename to tests/procedures/nmv_procedures/semiparametric_mu_estimation/test_validate_kwargs.py index fe4a2ad..c0dc5de 100644 --- a/tests/algorithms/nmv_algorithms/semiparametric_mu_estimation/test_validate_kwargs.py +++ b/tests/procedures/nmv_procedures/semiparametric_mu_estimation/test_validate_kwargs.py @@ -3,7 +3,7 @@ import numpy as np import pytest -from src.procedures.semiparametric.nvm_semi_param_algorithms.mu_estimation import NMVEstimationMu +from src.procedures.semiparametric.nvm_semiparametric.mu_estimation import NMVEstimationMu def _test_omega(x: float) -> float: diff --git a/tests/algorithms/nv_algorithms/__init__.py b/tests/procedures/nv_procedures/__init__.py similarity index 100% rename from tests/algorithms/nv_algorithms/__init__.py rename to tests/procedures/nv_procedures/__init__.py diff --git a/tests/algorithms/nv_algorithms/test_g_estimation.py b/tests/procedures/nv_procedures/test_g_estimation.py similarity index 100% rename from tests/algorithms/nv_algorithms/test_g_estimation.py rename to tests/procedures/nv_procedures/test_g_estimation.py diff --git a/tests/algorithms/nv_algorithms/test_post_widder_nv.py b/tests/procedures/nv_procedures/test_post_widder_nv.py similarity index 100% rename from tests/algorithms/nv_algorithms/test_post_widder_nv.py rename to tests/procedures/nv_procedures/test_post_widder_nv.py diff --git a/tests/algorithms/support_algorithms/test_rqmc.py b/tests/procedures/support_procedures/test_rqmc.py similarity index 100% rename from tests/algorithms/support_algorithms/test_rqmc.py rename to tests/procedures/support_procedures/test_rqmc.py From 0d635b37d23357c86363dc2e084367e63be23216 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=D0=90=D0=BD=D0=B3=D0=B5=D0=BB=D0=B8=D0=BD=D0=B0?= Date: Sat, 7 Jun 2025 18:23:29 +0300 Subject: [PATCH 13/24] fix: Fixing errors from the previous commit In the previous commit "NVSemiParametricGEstimationPostWidder" were mistakenly replaced by "NMVEstimationDensityPW" --- src/procedures/__init__.py | 2 +- .../nv_semiparametric/g_estimation_post_widder.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index 59d7bd6..d9f8503 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -11,7 +11,7 @@ NVEstimationDensityInvMT, ) from src.procedures.semiparametric.nv_semiparametric.g_estimation_post_widder import ( - NMVEstimationDensityPW, + NVEstimationDensityPW, ) from src.procedures.semiparametric.nvm_semiparametric.g_estimation_given_mu import ( NMVEstimationDensityInvMTquad, diff --git a/src/procedures/semiparametric/nv_semiparametric/g_estimation_post_widder.py b/src/procedures/semiparametric/nv_semiparametric/g_estimation_post_widder.py index fb34b9b..14c6deb 100644 --- a/src/procedures/semiparametric/nv_semiparametric/g_estimation_post_widder.py +++ b/src/procedures/semiparametric/nv_semiparametric/g_estimation_post_widder.py @@ -13,7 +13,7 @@ X_DATA_DEFAULT_VALUE = [1.0] -class NMVEstimationDensityPW: +class NVEstimationDensityPW: """Estimation of mixing density function g (xi density function) of NVM mixture represented in classical form Y = alpha + sigma*sqrt(xi)*N, where alpha = 0 and mu, sigma are given. From a2789f7ac9126c25e0d58d33bd194cba07a87859 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Sun, 25 May 2025 06:58:12 +0300 Subject: [PATCH 14/24] feat: add batch support and refactor mixture evaluation logic - Updated `moment`, `pdf`, `cdf`, and `logpdf` methods to support both scalar and list inputs for efficient batch computation. - Extracted shared evaluation logic (e.g., `moment`, `pdf`, `cdf`, `logpdf`) into the `AbstractMixtures` base class to eliminate code duplication across subclasses. - Centralized input validation, RQMC integration, and distribution handling within the base class for better maintainability. - Fixed type annotations to align with abstract method signatures and resolve `mypy` type-checking issues. These changes improve code clarity, enable vectorized evaluations, and make future extensions easier to implement. --- src/mixtures/abstract_mixture.py | 67 +++++++- src/mixtures/nm_mixture.py | 267 +++++++------------------------ src/mixtures/nmv_mixture.py | 198 ++++++++--------------- src/mixtures/nv_mixture.py | 85 +++++----- 4 files changed, 228 insertions(+), 389 deletions(-) diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index 13c8e42..862a74e 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -1,6 +1,8 @@ from abc import ABCMeta, abstractmethod from dataclasses import fields -from typing import Any +from typing import Any, List, Tuple, Union, Dict +import numpy as np +from numpy.typing import NDArray from scipy.stats import rv_continuous from scipy.stats.distributions import rv_frozen @@ -15,7 +17,6 @@ class AbstractMixtures(metaclass=ABCMeta): @abstractmethod def __init__(self, mixture_form: str, **kwargs: Any) -> None: """ - Args: mixture_form: Form of Mixture classical or Canonical **kwargs: Parameters of Mixture @@ -28,16 +29,68 @@ def __init__(self, mixture_form: str, **kwargs: Any) -> None: raise AssertionError(f"Unknown mixture form: {mixture_form}") @abstractmethod - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: ... + def _compute_moment(self, n: int, params: Dict) -> Tuple[float, float]: + ... + + def compute_moment( + self, x: Union[List[int], int, NDArray[np.float64]], params: Dict + ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: + if isinstance(x, np.ndarray): + return np.array([self._compute_moment(xp, params) for xp in x], dtype=object) + elif isinstance(x, list): + return [self._compute_moment(xp, params) for xp in x] + elif isinstance(x, int): + return self._compute_moment(x, params) + else: + raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: ... + def _compute_pdf(self, x: float, params: Dict) -> Tuple[float, float]: + ... + + def compute_pdf( + self, x: Union[List[float], float, NDArray[np.float64]], params: Dict + ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: + if isinstance(x, np.ndarray): + return np.array([self._compute_pdf(xp, params) for xp in x], dtype=object) + elif isinstance(x, list): + return [self._compute_pdf(xp, params) for xp in x] + elif isinstance(x, float): + return self._compute_pdf(x, params) + else: + raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: ... + def _compute_logpdf(self, x: float, params: Dict) -> Tuple[float, float]: + ... + + def compute_logpdf( + self, x: Union[List[float], float, NDArray[np.float64]], params: Dict + ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: + if isinstance(x, np.ndarray): + return np.array([self._compute_logpdf(xp, params) for xp in x], dtype=object) + elif isinstance(x, list): + return [self._compute_logpdf(xp, params) for xp in x] + elif isinstance(x, float): + return self._compute_logpdf(x, params) + else: + raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: ... + def _compute_cdf(self, x: float, rqmc_params: Dict[str, Any]) -> Tuple[float, float]: + ... + + def compute_cdf( + self, x: Union[List[float], float, NDArray[np.float64]], params: Dict + ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: + if isinstance(x, np.ndarray): + return np.array([self._compute_cdf(xp, params) for xp in x], dtype=object) + elif isinstance(x, list): + return [self._compute_cdf(xp, params) for xp in x] + elif isinstance(x, float): + return self._compute_cdf(x, params) + else: + raise TypeError(f"Unsupported type for x: {type(x)}") def _params_validation(self, data_collector: Any, params: dict[str, float | rv_continuous | rv_frozen]) -> Any: """Mixture Parameters Validation @@ -64,4 +117,4 @@ def _params_validation(self, data_collector: Any, params: dict[str, float | rv_c raise ValueError(f"Unexpected parameter {pair[0]}") if not isinstance(pair[1], names_and_types[pair[0]]): raise ValueError(f"Type missmatch: {pair[0]} should be {names_and_types[pair[0]]}, not {type(pair[1])}") - return data_collector(**params) + return data_collector(**params) \ No newline at end of file diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 33c7580..5c3ee5a 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -1,5 +1,5 @@ from dataclasses import dataclass -from typing import Any +from typing import Any, List, Tuple import numpy as np from scipy.integrate import quad @@ -14,9 +14,6 @@ @dataclass class _NMMClassicDataCollector: - """TODO: Change typing from float | int | etc to Protocol with __addition__ __multiplication__ __subtraction__""" - - """Data Collector for parameters of classical NMM""" alpha: float | int | np.int64 beta: float | int | np.int64 gamma: float | int | np.int64 @@ -25,9 +22,6 @@ class _NMMClassicDataCollector: @dataclass class _NMMCanonicalDataCollector: - """TODO: Change typing from float | int | etc to Protocol with __addition__ __multiplication__ __subtraction__""" - - """Data Collector for parameters of canonical NMM""" sigma: float | int | np.int64 distribution: rv_frozen | rv_continuous @@ -37,225 +31,74 @@ class NormalMeanMixtures(AbstractMixtures): _canonical_collector = _NMMCanonicalDataCollector def __init__(self, mixture_form: str, **kwargs: Any) -> None: - """ - Read Doc of Parent Method - """ - + self.mixture_form = mixture_form super().__init__(mixture_form, **kwargs) def _params_validation(self, data_collector: Any, params: dict[str, float | rv_continuous | rv_frozen]) -> Any: - """ - Read parent method doc - - Raises: - ValueError: If canonical Mixture has negative sigma parameter - - """ - data_class = super()._params_validation(data_collector, params) if hasattr(data_class, "sigma") and data_class.sigma <= 0: - raise ValueError("Sigma cant be zero or negative") + raise ValueError("Sigma can't be zero or negative") if hasattr(data_class, "gamma") and data_class.gamma == 0: - raise ValueError("Gamma cant be zero") + raise ValueError("Gamma can't be zero") return data_class - def _classical_moment(self, n: int, params: dict) -> tuple[float, float]: - """ - Compute n-th moment of classical NMM - - Args: - n (): Moment ordinal - params (): Parameters of integration algorithm - - Returns: moment approximation and error tolerance - - """ - mixture_moment = 0 - error_tolerance = 0 - for k in range(0, n + 1): - for l in range(0, k + 1): - coefficient = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma**l) - mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (k - l), 0, 1, **params) - error_tolerance += (self.params.beta ** (k - l)) * mixing_moment[1] - mixture_moment += coefficient * (self.params.alpha ** (n - k)) * mixing_moment[0] * norm.moment(l) - return mixture_moment, error_tolerance - - def _canonical_moment(self, n: int, params: dict) -> tuple[float, float]: - """ - Compute n-th moment of canonical NMM - - Args: - n (): Moment ordinal - params (): Parameters of integration algorithm - - Returns: moment approximation and error tolerance - - """ + def _compute_moment(self, n: int, params: dict) -> Tuple[float, float]: mixture_moment = 0 error_tolerance = 0 - for k in range(0, n + 1): - coefficient = binom(n, n - k) * (self.params.sigma**k) - mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (n - k), 0, 1, **params) - error_tolerance += mixing_moment[1] - mixture_moment += coefficient * mixing_moment[0] * norm.moment(k) + if self.mixture_form == "classical": + for k in range(0, n + 1): + for l in range(0, k + 1): + coefficient = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma ** l) + mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (k - l), 0, 1, **params) + error_tolerance += (self.params.beta ** (k - l)) * mixing_moment[1] + mixture_moment += coefficient * (self.params.alpha ** (n - k)) * mixing_moment[0] * norm.moment(l) + else: + for k in range(0, n + 1): + coefficient = binom(n, n - k) * (self.params.sigma ** k) + mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (n - k), 0, 1, **params) + error_tolerance += mixing_moment[1] + mixture_moment += coefficient * mixing_moment[0] * norm.moment(k) return mixture_moment, error_tolerance - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: - """ - Compute n-th moment of NMM - - Args: - n (): Moment ordinal - params (): Parameters of integration algorithm - - Returns: moment approximation and error tolerance - - """ - if isinstance(self.params, _NMMClassicDataCollector): - return self._classical_moment(n, params) - return self._canonical_moment(n, params) - - def _canonical_compute_cdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Equation for canonical cdf - Args: - x (): point - params (): parameters of RQMC algorithm - - Returns: computed cdf and error tolerance - - """ - rqmc = RQMC(lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), **params) + def _compute_cdf(self, x: float, params: dict) -> Tuple[float, float]: + if self.mixture_form == "classical": + rqmc = RQMC( + lambda u: norm.cdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma)), + **params, + ) + else: + rqmc = RQMC( + lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), + **params, + ) return rqmc() - def _classical_compute_cdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Equation for classic cdf - Args: - x (): point - params (): parameters of RQMC algorithm - - Returns: computed cdf and error tolerance - - """ - rqmc = RQMC( - lambda u: norm.cdf( - (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) - ), - **params - ) + def _compute_pdf(self, x: float, params: dict) -> Tuple[float, float]: + if self.mixture_form == "classical": + rqmc = RQMC( + lambda u: (1 / np.abs(self.params.gamma)) + * norm.pdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma)), + **params, + ) + else: + rqmc = RQMC( + lambda u: (1 / np.abs(self.params.sigma)) * norm.pdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), + **params, + ) return rqmc() - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Choose equation for cdf estimation depends on Mixture form - Args: - x (): point - params (): parameters of RQMC algorithm - - Returns: Computed pdf and error tolerance - - """ - if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_cdf(x, params) - return self._classical_compute_cdf(x, params) - - def _canonical_compute_pdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Equation for canonical pdf - Args: - x (): point - params (): parameters of RQMC algorithm - - Returns: computed pdf and error tolerance - - """ - rqmc = RQMC( - lambda u: (1 / np.abs(self.params.sigma)) - * norm.pdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), - **params - ) - return rqmc() - - def _classical_compute_pdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Equation for classic pdf - Args: - x (): point - params (): parameters of RQMC algorithm - - Returns: computed pdf and error tolerance - - """ - rqmc = RQMC( - lambda u: (1 / np.abs(self.params.gamma)) - * norm.pdf( - (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) - ), - **params - ) + def _compute_logpdf(self, x: float, params: dict) -> Tuple[float, float]: + if self.mixture_form == "classical": + rqmc = LogRQMC( + lambda u: ( + np.log(1 / np.abs(self.params.gamma)) + + norm.logpdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma)) + ), + **params, + ) + else: + rqmc = LogRQMC( + lambda u: np.log(1 / np.abs(self.params.sigma)) + norm.logpdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), + **params, + ) return rqmc() - - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Choose equation for pdf estimation depends on Mixture form - Args: - x (): point - params (): parameters of RQMC algorithm - - Returns: Computed pdf and error tolerance - - """ - if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_pdf(x, params) - return self._classical_compute_pdf(x, params) - - def _classical_compute_log_pdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Equation for classic log pdf - Args: - x (): point - params (): parameters of LogRQMC algorithm - - Returns: computed log pdf and error tolerance - - """ - rqmc = LogRQMC( - lambda u: np.log(1 / np.abs(self.params.gamma)) - + norm.logpdf( - (x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma) - ), - **params - ) - return rqmc() - - def _canonical_compute_log_pdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Equation for canonical log pdf - Args: - x (): point - params (): parameters of LogRQMC algorithm - - Returns: computed log pdf and error tolerance - - """ - rqmc = LogRQMC( - lambda u: np.log(1 / np.abs(self.params.sigma)) - + norm.logpdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), - **params - ) - return rqmc() - - def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: - """ - Choose equation for log pdf estimation depends on Mixture form - Args: - x (): point - params (): parameters of LogRQMC algorithm - - Returns: Computed log pdf and error tolerance - - """ - if isinstance(self.params, _NMMCanonicalDataCollector): - return self._canonical_compute_log_pdf(x, params) - return self._classical_compute_log_pdf(x, params) diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index c090ff3..c6dc75b 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -38,154 +38,98 @@ class NormalMeanVarianceMixtures(AbstractMixtures): _canonical_collector = _NMVMCanonicalDataCollector def __init__(self, mixture_form: str, **kwargs: Any) -> None: + self.mixture_form = mixture_form super().__init__(mixture_form, **kwargs) - def _classical_moment(self, n: int, params: dict) -> tuple[float, float]: - """ - Compute n-th moment of classical NMM - - Args: - n (): Moment ordinal - params (): Parameters of integration algorithm - - Returns: moment approximation and error tolerance - - """ - - def integral_func(u: float) -> float: - result = 0 - for k in range(0, n + 1): - for l in range(0, k + 1): - result += ( - binom(n, n - k) - * binom(k, k - l) - * (self.params.beta ** (k - l)) - * (self.params.gamma**l) - * self.params.distribution.ppf(u) ** (k - l / 2) - * (self.params.alpha ** (n - k)) - * norm.moment(l) - ) - return result - - rqmc = RQMC(lambda u: integral_func(u), **params) - return rqmc() - - def _canonical_moment(self, n: int, params: dict) -> tuple[float, float]: - """ - Compute n-th moment of classical NMM - - Args: - n (): Moment ordinal - params (): Parameters of integration algorithm - - Returns: moment approximation and error tolerance - - """ - + def _compute_moment(self, n: int, params: dict) -> tuple[float, float]: def integral_func(u: float) -> float: result = 0 for k in range(0, n + 1): for l in range(0, k + 1): - result += ( - binom(n, n - k) - * binom(k, k - l) - * (self.params.nu ** (k - l)) - * self.params.distribution.ppf(u) ** (k - l / 2) - * (self.params.alpha ** (n - k)) - * norm.moment(l) - ) + if self.mixture_form == "classical": + result += ( + binom(n, n - k) + * binom(k, k - l) + * (self.params.beta ** (k - l)) + * (self.params.gamma ** l) + * self.params.distribution.ppf(u) ** (k - l / 2) + * (self.params.alpha ** (n - k)) + * norm.moment(l) + ) + else: + result += ( + binom(n, n - k) + * binom(k, k - l) + * (self.params.nu ** (k - l)) + * self.params.distribution.ppf(u) ** (k - l / 2) + * (self.params.alpha ** (n - k)) + * norm.moment(l) + ) return result - rqmc = RQMC(lambda u: integral_func(u), **params) + rqmc = RQMC(integral_func) return rqmc() - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: - if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_moment(n, params) - return self._canonical_moment(n, params) - - def _classical_cdf(self, x: float, params: dict) -> tuple[float, float]: - def _inner_func(u: float) -> float: - ppf = lru_cache()(self.params.distribution.ppf)(u) - point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - ( - self.params.beta / self.params.gamma * np.sqrt(ppf) - ) - return norm.cdf(point) - - rqmc = RQMC(lambda u: _inner_func(u), **params) - return rqmc() - def _canonical_cdf(self, x: float, params: dict) -> tuple[float, float]: - def _inner_func(u: float) -> float: + def _compute_cdf(self, x: float, params: dict) -> tuple[float, float]: + def inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) - point = (x - self.params.alpha) / (np.sqrt(ppf)) - (self.params.mu * np.sqrt(ppf)) + if self.mixture_form == "classical": + point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - ( + self.params.beta / self.params.gamma * np.sqrt(ppf) + ) + else: + point = (x - self.params.alpha) / (np.sqrt(ppf)) - (self.params.mu * np.sqrt(ppf)) return norm.cdf(point) - rqmc = RQMC(lambda u: _inner_func(u), **params) + rqmc = RQMC(inner_func) return rqmc() - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: - if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_cdf(x, params) - return self._canonical_cdf(x, params) - def _classical_pdf(self, x: float, params: dict) -> tuple[float, float]: - def _inner_func(u: float) -> float: + def _compute_pdf(self, x: float, params: dict) -> tuple[float, float]: + def inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) - return ( - 1 - / np.sqrt(2 * np.pi * ppf * self.params.gamma**2) - * np.exp( - -((x - self.params.alpha) ** 2 + self.params.beta**2 * ppf**2) / (2 * ppf * self.params.gamma**2) + if self.mixture_form == "classical": + return ( + 1 + / np.sqrt(2 * np.pi * ppf * self.params.gamma ** 2) + * np.exp( + -((x - self.params.alpha) ** 2 + self.params.beta ** 2 * ppf ** 2) + / (2 * ppf * self.params.gamma ** 2) + ) + ) + else: + return ( + 1 + / np.sqrt(2 * np.pi * ppf) + * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu ** 2 * ppf ** 2) / (2 * ppf)) ) - ) - rqmc = RQMC(lambda u: _inner_func(u), **params)() - return np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma**2) * rqmc[0], rqmc[1] + rqmc_res = RQMC(inner_func)() + if self.mixture_form == "classical": + val = np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2) * rqmc_res[0] + else: + val = np.exp(self.params.mu * (x - self.params.alpha)) * rqmc_res[0] - def _canonical_pdf(self, x: float, params: dict) -> tuple[float, float]: - def _inner_func(u: float) -> float: - ppf = self.params.distribution.ppf(u) - return ( - 1 - / np.sqrt(2 * np.pi * ppf) - * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu**2 * ppf**2) / (2 * ppf)) - ) - - rqmc = RQMC(lambda u: _inner_func(u), **params) - res = rqmc() - return np.exp(self.params.mu * (x - self.params.alpha)) * res[0], res[1] - - def _classical_log_pdf(self, x: float, params: dict) -> tuple[float, float]: - def _inner_func(u: float) -> float: - ppf = self.params.distribution.ppf(u) - return -( - (x - self.params.alpha) ** 2 - + ppf**2 * self.params.beta**2 - + ppf * self.params.gamma**2 * np.log(2 * np.pi * ppf * self.params.gamma**2) - ) / (2 * ppf * self.params.gamma**2) + return val, rqmc_res[1] - rqmc = LogRQMC(lambda u: _inner_func(u), **params) - return rqmc() - def _canonical_log_pdf(self, x: float, params: dict) -> tuple[float, float]: - def _inner_func(u: float) -> float: + def _compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: + def inner_func(u: float) -> float: ppf = self.params.distribution.ppf(u) - return -((x - self.params.alpha) ** 2 + ppf**2 * self.params.mu**2 + ppf * np.log(2 * np.pi * ppf)) / ( - 2 * ppf - ) - - rqmc = LogRQMC(lambda u: _inner_func(u), **params) - return rqmc() - - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: - if isinstance(self.params, _NMVMClassicDataCollector): - return self._classical_pdf(x, params) - return self._canonical_pdf(x, params) - - def compute_logpdf(self, x: float, params: dict) -> tuple[Any, float]: - if isinstance(self.params, _NMVMClassicDataCollector): - int_result = self._classical_log_pdf(x, params) - return self.params.beta * (x - self.params.alpha) / self.params.gamma**2 + int_result[0], int_result[1] - int_result = self._canonical_log_pdf(x, params) - return self.params.mu * (x - self.params.alpha) + int_result[0], int_result[1] + if self.mixture_form == "classical": + return -( + (x - self.params.alpha) ** 2 + + ppf ** 2 * self.params.beta ** 2 + + ppf * self.params.gamma ** 2 * np.log(2 * np.pi * ppf * self.params.gamma ** 2) + ) / (2 * ppf * self.params.gamma ** 2) + else: + return -((x - self.params.alpha) ** 2 + ppf ** 2 * self.params.mu ** 2 + ppf * np.log(2 * np.pi * ppf)) / (2 * ppf) + + rqmc_res = LogRQMC(inner_func)() + if self.mixture_form == "classical": + val = self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2 + rqmc_res[0] + else: + val = self.params.mu * (x - self.params.alpha) + rqmc_res[0] + + return val, rqmc_res[1] diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index 7905e6e..506dbf4 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -32,62 +32,61 @@ class _NVMCanonicalDataCollector: class NormalVarianceMixtures(AbstractMixtures): - _classical_collector = _NVMClassicDataCollector _canonical_collector = _NVMCanonicalDataCollector def __init__(self, mixture_form: str, **kwargs: Any) -> None: + self.mixture_form = mixture_form super().__init__(mixture_form, **kwargs) - def compute_moment(self, n: int, params: dict) -> tuple[float, float]: - """ - Compute n-th moment of NVM - Args: - n (): Moment ordinal - params (): Parameters of integration algorithm - Returns: moment approximation and error tolerance - """ - gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 - - def integrate_func(u: float) -> float: + def _compute_moment(self, n: int, rqmc_params: dict[str, Any]) -> tuple[float, float]: + gamma = getattr(self.params, 'gamma', 1) + + def integrand(u: float) -> float: return sum( - [ - binom(n, k) - * (gamma**k) - * (self.params.alpha ** (n - k)) - * (self.params.distribution.ppf(u) ** (k / 2)) - * norm.moment(k) - for k in range(0, n + 1) - ] + binom(n, k) + * (gamma ** k) + * (self.params.alpha ** (n - k)) + * (self.params.distribution.ppf(u) ** (k / 2)) + * norm.moment(k) + for k in range(n + 1) ) - result = RQMC(integrate_func, **params)() - return result + return RQMC(integrand, **rqmc_params)() + + def _compute_cdf(self, x: float, rqmc_params: dict[str, Any]) -> tuple[float, float]: + gamma = getattr(self.params, 'gamma', 1) + param_norm = norm(0, gamma) + + def integrand(u: float) -> float: + return param_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))) + + return RQMC(integrand, **rqmc_params)() + + + def _compute_pdf(self, x: float, rqmc_params: dict[str, Any]) -> tuple[float, float]: + gamma = getattr(self.params, 'gamma', 1) + d = (x - self.params.alpha) ** 2 / gamma ** 2 + + def integrand(u: float) -> float: + return self._integrand_func(u, d, gamma) + + return RQMC(integrand, **rqmc_params)() + + def _compute_logpdf(self, x: float, rqmc_params: dict[str, Any]) -> tuple[float, float]: + gamma = getattr(self.params, 'gamma', 1) + d = (x - self.params.alpha) ** 2 / gamma ** 2 + + def integrand(u: float) -> float: + return self._log_integrand_func(u, d, gamma) - def compute_cdf(self, x: float, params: dict) -> tuple[float, float]: - parametric_norm = norm(0, self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1) - rqmc = RQMC( - lambda u: parametric_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))), **params - ) - return rqmc() + return LogRQMC(integrand, **rqmc_params)() @lru_cache() def _integrand_func(self, u: float, d: float, gamma: float) -> float: ppf = self.params.distribution.ppf(u) - return (1 / np.sqrt(np.pi * 2 * ppf * np.abs(gamma**2))) * np.exp(-1 * d / (2 * ppf)) + return (1 / np.sqrt(np.pi * 2 * ppf * np.abs(gamma ** 2))) * np.exp(-d / (2 * ppf)) - def _log_integrand_func(self, u: float, d: float, gamma: float | int | np.int64) -> float: + def _log_integrand_func(self, u: float, d: float, gamma: float) -> float: ppf = self.params.distribution.ppf(u) - return -(ppf * np.log(np.pi * 2 * ppf * gamma**2) + d) / (2 * ppf) - - def compute_pdf(self, x: float, params: dict) -> tuple[float, float]: - gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 - d = (x - self.params.alpha) ** 2 / gamma**2 - rqmc = RQMC(lambda u: self._integrand_func(u, d, gamma), **params) - return rqmc() - - def compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: - gamma = self.params.gamma if isinstance(self.params, _NVMClassicDataCollector) else 1 - d = (x - self.params.alpha) ** 2 / gamma**2 - log_rqmc = LogRQMC(lambda u: self._log_integrand_func(u, d, gamma), **params) - return log_rqmc() + return -(ppf * np.log(np.pi * 2 * ppf * gamma ** 2) + d) / (2 * ppf) From 4296e344147ab87f96308876247d82b1cbca6a71 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 01:25:30 +0300 Subject: [PATCH 15/24] fixed --- src/mixtures/abstract_mixture.py | 101 ++++++++++++++------------- src/mixtures/nm_mixture.py | 115 +++++++++++++++---------------- src/mixtures/nmv_mixture.py | 115 ++++++++++++++----------------- src/mixtures/nv_mixture.py | 65 ++++++++--------- 4 files changed, 192 insertions(+), 204 deletions(-) diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index 2d63ad0..ca111a9 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -1,6 +1,6 @@ from abc import ABCMeta, abstractmethod from dataclasses import fields -from typing import Any, List, Tuple, Union, Dict +from typing import Any, List, Tuple, Union, Dict, Type import numpy as np from numpy.typing import NDArray @@ -8,6 +8,7 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator # default integrator class AbstractMixtures(metaclass=ABCMeta): """Base class for Mixtures""" @@ -15,14 +16,24 @@ class AbstractMixtures(metaclass=ABCMeta): _classical_collector: Any _canonical_collector: Any - @abstractmethod - def __init__(self, mixture_form: str, **kwargs: Any) -> None: + def __init__( + self, + mixture_form: str, + integrator_cls: Type[Integrator] = RQMCIntegrator, + integrator_params: Dict[str, Any] = None, + **kwargs: Any + ) -> None: """ Args: - mixture_form: Form of Mixture classical or Canonical - **kwargs: Parameters of Mixture + mixture_form: Form of Mixture classical or canonical + integrator_cls: Class implementing Integrator protocol (default: RQMCIntegrator) + integrator_params: Parameters for integrator constructor (default: {{}}) + **kwargs: Parameters of Mixture (alpha, gamma, etc.) """ self.mixture_form = mixture_form + self.integrator_cls = integrator_cls + self.integrator_params = integrator_params or {} + if mixture_form == "classical": self.params = self._params_validation(self._classical_collector, kwargs) elif mixture_form == "canonical": @@ -31,92 +42,88 @@ def __init__(self, mixture_form: str, **kwargs: Any) -> None: raise AssertionError(f"Unknown mixture form: {mixture_form}") @abstractmethod - def _compute_moment(self, n: int, params: Dict) -> Tuple[float, float]: + def _compute_moment(self, n: int) -> Tuple[float, float]: ... def compute_moment( - self, x: Union[List[int], int, NDArray[np.float64]], params: Dict + self, + x: Union[List[int], int, NDArray[np.float64]] ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_moment(xp, params) for xp in x], dtype=object) + return np.array([self._compute_moment(xp) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_moment(xp, params) for xp in x] + return [self._compute_moment(xp) for xp in x] elif isinstance(x, int): - return self._compute_moment(x, params) + return self._compute_moment(x) else: raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def _compute_pdf(self, x: float, params: Dict) -> Tuple[float, float]: + def _compute_pdf(self, x: float) -> Tuple[float, float]: ... def compute_pdf( - self, x: Union[List[float], float, NDArray[np.float64]], params: Dict + self, + x: Union[List[float], float, NDArray[np.float64]] ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_pdf(xp, params) for xp in x], dtype=object) + return np.array([self._compute_pdf(xp) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_pdf(xp, params) for xp in x] + return [self._compute_pdf(xp) for xp in x] elif isinstance(x, float): - return self._compute_pdf(x, params) + return self._compute_pdf(x) else: raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def _compute_logpdf(self, x: float, params: Dict) -> Tuple[float, float]: + def _compute_logpdf(self, x: float) -> Tuple[float, float]: ... def compute_logpdf( - self, x: Union[List[float], float, NDArray[np.float64]], params: Dict + self, + x: Union[List[float], float, NDArray[np.float64]] ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_logpdf(xp, params) for xp in x], dtype=object) + return np.array([self._compute_logpdf(xp) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_logpdf(xp, params) for xp in x] + return [self._compute_logpdf(xp) for xp in x] elif isinstance(x, float): - return self._compute_logpdf(x, params) + return self._compute_logpdf(x) else: raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def _compute_cdf(self, x: float, rqmc_params: Dict[str, Any]) -> Tuple[float, float]: + def _compute_cdf(self, x: float) -> Tuple[float, float]: ... def compute_cdf( - self, x: Union[List[float], float, NDArray[np.float64]], params: Dict + self, + x: Union[List[float], float, NDArray[np.float64]] ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_cdf(xp, params) for xp in x], dtype=object) + return np.array([self._compute_cdf(xp) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_cdf(xp, params) for xp in x] + return [self._compute_cdf(xp) for xp in x] elif isinstance(x, float): - return self._compute_cdf(x, params) + return self._compute_cdf(x) else: raise TypeError(f"Unsupported type for x: {type(x)}") - def _params_validation(self, data_collector: Any, params: dict[str, float | rv_continuous | rv_frozen]) -> Any: - """Mixture Parameters Validation - - Args: - data_collector: Dataclass that collect parameters of Mixture - params: Input parameters - - Returns: Instance of dataclass - - Raises: - ValueError: If given parameters is unexpected - ValueError: If parameter type is invalid - ValueError: If parameters age not given - - """ - + def _params_validation( + self, + data_collector: Any, + params: dict[str, float | rv_continuous | rv_frozen] + ) -> Any: + """Mixture Parameters Validation""" dataclass_fields = fields(data_collector) if len(params) != len(dataclass_fields): raise ValueError(f"Expected {len(dataclass_fields)} arguments, got {len(params)}") - names_and_types = dict((field.name, field.type) for field in dataclass_fields) - for pair in params.items(): - if pair[0] not in names_and_types: - raise ValueError(f"Unexpected parameter {pair[0]}") - if not isinstance(pair[1], names_and_types[pair[0]]): - raise ValueError(f"Type missmatch: {pair[0]} should be {names_and_types[pair[0]]}, not {type(pair[1])}") + names_and_types = {field.name: field.type for field in dataclass_fields} + for name, value in params.items(): + if name not in names_and_types: + raise ValueError(f"Unexpected parameter {name}") + if not isinstance(value, names_and_types[name]): + raise ValueError( + f"Type mismatch: {name} should be {names_and_types[name]}, not {type(value)}" + ) return data_collector(**params) \ No newline at end of file diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index de5521b..9f36165 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -1,18 +1,17 @@ from dataclasses import dataclass -from typing import Any, List, Tuple +from typing import Any, Type, Dict, Tuple, Union, List import numpy as np from scipy.special import binom from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.algorithms.support_algorithms.log_rqmc import LogRQMC from src.algorithms.support_algorithms.integrator import Integrator -from src.algorithms.support_algorithms.rqmc import RQMCIntegrator from src.algorithms.support_algorithms.quad_integrator import QuadIntegrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator +from src.algorithms.support_algorithms.log_rqmc import LogRQMC from src.mixtures.abstract_mixture import AbstractMixtures - @dataclass class _NMMClassicDataCollector: alpha: float | int | np.int64 @@ -20,20 +19,23 @@ class _NMMClassicDataCollector: gamma: float | int | np.int64 distribution: rv_frozen | rv_continuous - @dataclass class _NMMCanonicalDataCollector: sigma: float | int | np.int64 distribution: rv_frozen | rv_continuous - class NormalMeanMixtures(AbstractMixtures): _classical_collector = _NMMClassicDataCollector _canonical_collector = _NMMCanonicalDataCollector - def __init__(self, mixture_form: str, **kwargs: Any) -> None: - self.mixture_form = mixture_form - super().__init__(mixture_form, **kwargs) + def __init__( + self, + mixture_form: str, + integrator_cls: Type[Integrator] = RQMCIntegrator, + integrator_params: Dict[str, Any] = None, + **kwargs: Any + ) -> None: + super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) def _params_validation(self, data_collector: Any, params: dict[str, float | rv_continuous | rv_frozen]) -> Any: data_class = super()._params_validation(data_collector, params) @@ -43,64 +45,59 @@ def _params_validation(self, data_collector: Any, params: dict[str, float | rv_c raise ValueError("Gamma can't be zero") return data_class - def _compute_moment(self, n: int, params: dict) -> Tuple[float, float]: - mixture_moment = 0 - error_tolerance = 0 + def _compute_moment(self, n: int) -> Tuple[float, float]: + mixture_moment = 0.0 + error = 0.0 if self.mixture_form == "classical": - for k in range(0, n + 1): - for l in range(0, k + 1): - coefficient = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma ** l) - mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (k - l), 0, 1, **params) - error_tolerance += (self.params.beta ** (k - l)) * mixing_moment[1] - mixture_moment += coefficient * (self.params.alpha ** (n - k)) * mixing_moment[0] * norm.moment(l) + for k in range(n + 1): + for l in range(k + 1): + coeff = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma ** l) + def mix(u: float) -> float: + return self.params.distribution.ppf(u) ** (k - l) + integrator = self.integrator_cls(**(self.integrator_params or {})) + res = integrator.compute(mix) + mixture_moment += coeff * (self.params.alpha ** (n - k)) * res.value * norm.moment(l) + error += coeff * res.error * (self.params.alpha ** (n - k)) * norm.moment(l) else: - for k in range(0, n + 1): - coefficient = binom(n, n - k) * (self.params.sigma ** k) - mixing_moment = quad(lambda u: self.params.distribution.ppf(u) ** (n - k), 0, 1, **params) - error_tolerance += mixing_moment[1] - mixture_moment += coefficient * mixing_moment[0] * norm.moment(k) - return mixture_moment, error_tolerance + for k in range(n + 1): + coeff = binom(n, n - k) * (self.params.sigma ** k) + def mix(u: float) -> float: + return self.params.distribution.ppf(u) ** (n - k) + integrator = self.integrator_cls(**(self.integrator_params or {})) + res = integrator.compute(mix) + mixture_moment += coeff * res.value * norm.moment(k) + error += coeff * res.error * norm.moment(k) + return mixture_moment, error - def _compute_cdf(self, x: float, params: dict) -> Tuple[float, float]: + def _compute_cdf(self, x: float) -> Tuple[float, float]: if self.mixture_form == "classical": - rqmc = RQMC( - lambda u: norm.cdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma)), - **params, - ) + def fn(u: float) -> float: + return norm.cdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) else: - rqmc = RQMC( - lambda u: norm.cdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), - **params, - ) - return rqmc() + def fn(u: float) -> float: + return norm.cdf((x - self.params.distribution.ppf(u)) / abs(self.params.sigma)) + integrator = self.integrator_cls(**(self.integrator_params or {})) + res = integrator.compute(fn) + return res.value, res.error - def _compute_pdf(self, x: float, params: dict) -> Tuple[float, float]: + def _compute_pdf(self, x: float) -> Tuple[float, float]: if self.mixture_form == "classical": - rqmc = RQMC( - lambda u: (1 / np.abs(self.params.gamma)) - * norm.pdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma)), - **params, - ) + def fn(u: float) -> float: + return (1 / abs(self.params.gamma)) * norm.pdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) else: - rqmc = RQMC( - lambda u: (1 / np.abs(self.params.sigma)) * norm.pdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), - **params, - ) - return rqmc() + def fn(u: float) -> float: + return (1 / abs(self.params.sigma)) * norm.pdf((x - self.params.distribution.ppf(u)) / abs(self.params.sigma)) + integrator = self.integrator_cls(**(self.integrator_params or {})) + res = integrator.compute(fn) + return res.value, res.error - def _compute_logpdf(self, x: float, params: dict) -> Tuple[float, float]: + def _compute_logpdf(self, x: float) -> Tuple[float, float]: if self.mixture_form == "classical": - rqmc = LogRQMC( - lambda u: ( - np.log(1 / np.abs(self.params.gamma)) - + norm.logpdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / np.abs(self.params.gamma)) - ), - **params, - ) + def fn(u: float) -> float: + return np.log(1 / abs(self.params.gamma)) + norm.logpdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) else: - rqmc = LogRQMC( - lambda u: np.log(1 / np.abs(self.params.sigma)) + norm.logpdf((x - self.params.distribution.ppf(u)) / np.abs(self.params.sigma)), - **params, - ) - return rqmc() - + def fn(u: float) -> float: + return np.log(1 / abs(self.params.sigma)) + norm.logpdf((x - self.params.distribution.ppf(u)) / abs(self.params.sigma)) + integrator = self.integrator_cls(**(self.integrator_params or {})) + res = integrator.compute(fn) + return res.value, res.error diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index 2d922ec..d7bc8cb 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -1,124 +1,109 @@ from dataclasses import dataclass from functools import lru_cache -from typing import Any +from typing import Any, Type, Dict import numpy as np from scipy.special import binom -from scipy.stats import geninvgauss, norm, rv_continuous +from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.algorithms.support_algorithms.log_rqmc import LogRQMC from src.algorithms.support_algorithms.integrator import Integrator from src.algorithms.support_algorithms.rqmc import RQMCIntegrator +from src.algorithms.support_algorithms.log_rqmc import LogRQMC +from src.algorithms.support_algorithms.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures - @dataclass class _NMVMClassicDataCollector: - """TODO: Change typing from float | int | etc to Protocol with __addition__ __multiplication__ __subtraction__""" - - """Data Collector for parameters of classical NMVM""" alpha: float | int | np.int64 beta: float | int | np.int64 gamma: float | int | np.int64 distribution: rv_frozen | rv_continuous - @dataclass class _NMVMCanonicalDataCollector: - """TODO: Change typing from float | int | etc to Protocol with __addition__ __multiplication__ __subtraction__""" - - """Data Collector for parameters of canonical NMVM""" alpha: float | int | np.int64 mu: float | int | np.int64 distribution: rv_frozen | rv_continuous - class NormalMeanVarianceMixtures(AbstractMixtures): _classical_collector = _NMVMClassicDataCollector _canonical_collector = _NMVMCanonicalDataCollector - def __init__(self, mixture_form: str, **kwargs: Any) -> None: - self.mixture_form = mixture_form - super().__init__(mixture_form, **kwargs) - - def _compute_moment(self, n: int, params: dict) -> tuple[float, float]: - def integral_func(u: float) -> float: - result = 0 - for k in range(0, n + 1): - for l in range(0, k + 1): + def __init__( + self, + mixture_form: str, + integrator_cls: Type[Integrator] = RQMCIntegrator, + integrator_params: Dict[str, Any] = None, + **kwargs: Any + ) -> None: + super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) + + def _compute_moment(self, n: int) -> tuple[float, float]: + def integrand(u: float) -> float: + result = 0.0 + for k in range(n + 1): + for l in range(k + 1): if self.mixture_form == "classical": result += ( binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma ** l) - * self.params.distribution.ppf(u) ** (k - l / 2) - * (self.params.alpha ** (n - k)) + * (self.params.distribution.ppf(u) ** (k - l/2)) * norm.moment(l) ) else: result += ( binom(n, n - k) * binom(k, k - l) - * (self.params.nu ** (k - l)) - * self.params.distribution.ppf(u) ** (k - l / 2) - * (self.params.alpha ** (n - k)) + * (self.params.mu ** (k - l)) + * (self.params.distribution.ppf(u) ** (k - l/2)) * norm.moment(l) ) return result - rqmc = RQMC(integral_func) - return rqmc() - - - def _compute_cdf(self, x: float, params: dict) -> tuple[float, float]: - def inner_func(u: float) -> float: + integrator = self.integrator_cls(**(self.integrator_params or {})) + result = integrator.compute(integrand) + return result.value, result.error + def _compute_cdf(self, x: float) -> tuple[float, float]: + def integrand(u: float) -> float: ppf = self.params.distribution.ppf(u) if self.mixture_form == "classical": - point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - ( - self.params.beta / self.params.gamma * np.sqrt(ppf) - ) + point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - (self.params.beta / self.params.gamma * np.sqrt(ppf)) else: - point = (x - self.params.alpha) / (np.sqrt(ppf)) - (self.params.mu * np.sqrt(ppf)) + point = (x - self.params.alpha) / np.sqrt(ppf) - (self.params.mu * np.sqrt(ppf)) return norm.cdf(point) - rqmc = RQMC(inner_func) - return rqmc() - - - def _compute_pdf(self, x: float, params: dict) -> tuple[float, float]: - def inner_func(u: float) -> float: + integrator = self.integrator_cls(**(self.integrator_params or {})) + result = integrator.compute(integrand) + return result.value, result.error + def _compute_pdf(self, x: float) -> tuple[float, float]: + def integrand(u: float) -> float: ppf = self.params.distribution.ppf(u) if self.mixture_form == "classical": return ( - 1 - / np.sqrt(2 * np.pi * ppf * self.params.gamma ** 2) - * np.exp( - -((x - self.params.alpha) ** 2 + self.params.beta ** 2 * ppf ** 2) - / (2 * ppf * self.params.gamma ** 2) - ) + 1 / np.sqrt(2 * np.pi * ppf * self.params.gamma ** 2) + * np.exp(-((x - self.params.alpha) ** 2 + self.params.beta ** 2 * ppf ** 2) / (2 * ppf * self.params.gamma ** 2)) ) else: return ( - 1 - / np.sqrt(2 * np.pi * ppf) + 1 / np.sqrt(2 * np.pi * ppf) * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu ** 2 * ppf ** 2) / (2 * ppf)) ) - rqmc_res = RQMC(inner_func)() + integrator = self.integrator_cls(**(self.integrator_params or {})) + result = integrator.compute(integrand) if self.mixture_form == "classical": - val = np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2) * rqmc_res[0] + val = np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2) * result.value else: - val = np.exp(self.params.mu * (x - self.params.alpha)) * rqmc_res[0] - - return val, rqmc_res[1] + val = np.exp(self.params.mu * (x - self.params.alpha)) * result.value + return val, result.error - - def _compute_logpdf(self, x: float, params: dict) -> tuple[float, float]: - def inner_func(u: float) -> float: + def _compute_logpdf(self, x: float) -> tuple[float, float]: + def integrand(u: float) -> float: ppf = self.params.distribution.ppf(u) if self.mixture_form == "classical": return -( @@ -129,11 +114,15 @@ def inner_func(u: float) -> float: else: return -((x - self.params.alpha) ** 2 + ppf ** 2 * self.params.mu ** 2 + ppf * np.log(2 * np.pi * ppf)) / (2 * ppf) - rqmc_res = LogRQMC(inner_func)() + integrator = self.integrator_cls(**(self.integrator_params or {})) + result = integrator.compute(integrand) if self.mixture_form == "classical": - val = self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2 + rqmc_res[0] + val = self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2 + result.value else: - val = self.params.mu * (x - self.params.alpha) + rqmc_res[0] - - return val, rqmc_res[1] + val = self.params.mu * (x - self.params.alpha) + result.value + return val, result.error + @lru_cache() + def _integrand_func(self, u: float, d: float, gamma: float) -> float: + ppf = self.params.distribution.ppf(u) + return (1 / np.sqrt(2 * np.pi * ppf * abs(gamma) ** 2)) * np.exp(-d / (2 * ppf)) diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index b85be7a..c53a67a 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -1,48 +1,46 @@ from dataclasses import dataclass from functools import lru_cache -from typing import Any +from typing import Any, Type, Dict import numpy as np from scipy.special import binom from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.algorithms.support_algorithms.log_rqmc import LogRQMC from src.algorithms.support_algorithms.integrator import Integrator from src.algorithms.support_algorithms.rqmc import RQMCIntegrator +from src.algorithms.support_algorithms.log_rqmc import LogRQMC +from src.algorithms.support_algorithms.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures - @dataclass class _NVMClassicDataCollector: - """TODO: Change typing from float | int | etc to Protocol with __addition__ __multiplication__ __subtraction__""" - - """Data Collector for parameters of classical NVM""" alpha: float | int | np.int64 gamma: float | int | np.int64 distribution: rv_frozen | rv_continuous - @dataclass class _NVMCanonicalDataCollector: - """TODO: Change typing from float | int | etc to Protocol with __addition__ __multiplication__ __subtraction__""" - - """Data Collector for parameters of canonical NVM""" alpha: float | int | np.int64 distribution: rv_frozen | rv_continuous - class NormalVarianceMixtures(AbstractMixtures): _classical_collector = _NVMClassicDataCollector _canonical_collector = _NVMCanonicalDataCollector - def __init__(self, mixture_form: str, **kwargs: Any) -> None: - self.mixture_form = mixture_form + def __init__( + self, + mixture_form: str, + integrator_cls: Type[Integrator] = RQMCIntegrator, + integrator_params: Dict[str, Any] = None, + **kwargs: Any + ) -> None: super().__init__(mixture_form, **kwargs) + self.integrator_cls = integrator_cls + self.integrator_params = integrator_params or {} - def _compute_moment(self, n: int, rqmc_params: dict[str, Any]) -> tuple[float, float]: + def _compute_moment(self, n: int, __: Dict[str, Any]) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) - def integrand(u: float) -> float: return sum( binom(n, k) @@ -52,46 +50,43 @@ def integrand(u: float) -> float: * norm.moment(k) for k in range(n + 1) ) - result = integrator.compute(func=integrate_func) + integrator = self.integrator_cls(**self.integrator_params) + result = integrator.compute(integrand) return result.value, result.error - return RQMC(integrand, **rqmc_params)() - - def _compute_cdf(self, x: float, rqmc_params: dict[str, Any]) -> tuple[float, float]: + def _compute_cdf(self, x: float, __: Dict[str, Any]) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) param_norm = norm(0, gamma) - def integrand(u: float) -> float: return param_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))) + integrator = self.integrator_cls(**self.integrator_params) + result = integrator.compute(integrand) + return result.value, result.error - return RQMC(integrand, **rqmc_params)() - - - def _compute_pdf(self, x: float, rqmc_params: dict[str, Any]) -> tuple[float, float]: + def _compute_pdf(self, x: float, __: Dict[str, Any]) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 - def integrand(u: float) -> float: return self._integrand_func(u, d, gamma) + integrator = self.integrator_cls(**self.integrator_params) + result = integrator.compute(integrand) + return result.value, result.error - return RQMC(integrand, **rqmc_params)() - - def _compute_logpdf(self, x: float, rqmc_params: dict[str, Any]) -> tuple[float, float]: + def _compute_logpdf(self, x: float, __: Dict[str, Any]) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 - def integrand(u: float) -> float: return self._log_integrand_func(u, d, gamma) - - return LogRQMC(integrand, **rqmc_params)() - + # For log-pdf you may choose LogRQMC or any integrator that supports it + integrator = self.integrator_cls(**self.integrator_params) + result = integrator.compute(integrand) + return result.value, result.error @lru_cache() def _integrand_func(self, u: float, d: float, gamma: float) -> float: ppf = self.params.distribution.ppf(u) - return (1 / np.sqrt(np.pi * 2 * ppf * np.abs(gamma ** 2))) * np.exp(-d / (2 * ppf)) + return (1 / np.sqrt(2 * np.pi * ppf * np.abs(gamma ** 2))) * np.exp(-d / (2 * ppf)) def _log_integrand_func(self, u: float, d: float, gamma: float) -> float: ppf = self.params.distribution.ppf(u) - return -(ppf * np.log(np.pi * 2 * ppf * gamma ** 2) + d) / (2 * ppf) - + return -(ppf * np.log(2 * np.pi * ppf * gamma ** 2) + d) / (2 * ppf) From 12b3205005fe240f0701dce953a44550e278fdf5 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 01:45:17 +0300 Subject: [PATCH 16/24] fixed --- src/mixtures/abstract_mixture.py | 2 +- src/mixtures/nm_mixture.py | 2 +- src/mixtures/nmv_mixture.py | 108 +++++++++++-------------------- src/mixtures/nv_mixture.py | 19 ++++-- 4 files changed, 52 insertions(+), 79 deletions(-) diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index ca111a9..503e061 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -126,4 +126,4 @@ def _params_validation( raise ValueError( f"Type mismatch: {name} should be {names_and_types[name]}, not {type(value)}" ) - return data_collector(**params) \ No newline at end of file + return data_collector(**params) diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 9f36165..54498ab 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -69,7 +69,7 @@ def mix(u: float) -> float: error += coeff * res.error * norm.moment(k) return mixture_moment, error - def _compute_cdf(self, x: float) -> Tuple[float, float]: + def _compute_cdf(self, x: float, params: Dict[str, Any]) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: return norm.cdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index d7bc8cb..ad997c2 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -1,6 +1,5 @@ from dataclasses import dataclass -from functools import lru_cache -from typing import Any, Type, Dict +from typing import Any, Type, Dict, Tuple import numpy as np from scipy.special import binom @@ -8,22 +7,19 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.integrator import Integrator -from src.algorithms.support_algorithms.rqmc import RQMCIntegrator -from src.algorithms.support_algorithms.log_rqmc import LogRQMC -from src.algorithms.support_algorithms.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @dataclass class _NMVMClassicDataCollector: alpha: float | int | np.int64 - beta: float | int | np.int64 + beta: float | int | np.int64 gamma: float | int | np.int64 distribution: rv_frozen | rv_continuous @dataclass class _NMVMCanonicalDataCollector: alpha: float | int | np.int64 - mu: float | int | np.int64 + mu: float | int | np.int64 distribution: rv_frozen | rv_continuous class NormalMeanVarianceMixtures(AbstractMixtures): @@ -33,96 +29,66 @@ class NormalMeanVarianceMixtures(AbstractMixtures): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator] = RQMCIntegrator, - integrator_params: Dict[str, Any] = None, + integrator_cls: Type[Integrator], + integrator_params: Dict[str, Any] | None = None, **kwargs: Any ) -> None: super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) - def _compute_moment(self, n: int) -> tuple[float, float]: + def _compute_moment(self, n: int) -> Tuple[float, float]: def integrand(u: float) -> float: - result = 0.0 + s = 0.0 for k in range(n + 1): for l in range(k + 1): if self.mixture_form == "classical": - result += ( - binom(n, n - k) - * binom(k, k - l) - * (self.params.beta ** (k - l)) + coef = binom(n, n - k) * binom(k, k - l) + term = ( + (self.params.beta ** (k - l)) * (self.params.gamma ** l) * (self.params.distribution.ppf(u) ** (k - l/2)) + * (self.params.alpha ** (n - k)) * norm.moment(l) ) else: - result += ( - binom(n, n - k) - * binom(k, k - l) - * (self.params.mu ** (k - l)) + coef = binom(n, n - k) * binom(k, k - l) + term = ( + (self.params.nu ** (k - l)) * (self.params.distribution.ppf(u) ** (k - l/2)) + * (self.params.alpha ** (n - k)) * norm.moment(l) ) - return result + s += coef * term if self.mixture_form == "classical" else term + return s - integrator = self.integrator_cls(**(self.integrator_params or {})) - result = integrator.compute(integrand) - return result.value, result.error + res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + return res.value, res.error - def _compute_cdf(self, x: float) -> tuple[float, float]: + def _compute_cdf(self, x: float) -> Tuple[float, float]: def integrand(u: float) -> float: - ppf = self.params.distribution.ppf(u) + p = self.params.distribution.ppf(u) if self.mixture_form == "classical": - point = (x - self.params.alpha) / (np.sqrt(ppf) * self.params.gamma) - (self.params.beta / self.params.gamma * np.sqrt(ppf)) - else: - point = (x - self.params.alpha) / np.sqrt(ppf) - (self.params.mu * np.sqrt(ppf)) - return norm.cdf(point) + return norm.cdf((x - self.params.alpha - self.params.beta * p) / abs(self.params.gamma)) + return norm.cdf((x - self.params.alpha) / np.sqrt(p) - self.params.mu * np.sqrt(p)) - integrator = self.integrator_cls(**(self.integrator_params or {})) - result = integrator.compute(integrand) - return result.value, result.error + res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + return res.value, res.error - def _compute_pdf(self, x: float) -> tuple[float, float]: + def _compute_pdf(self, x: float) -> Tuple[float, float]: def integrand(u: float) -> float: - ppf = self.params.distribution.ppf(u) + p = self.params.distribution.ppf(u) if self.mixture_form == "classical": - return ( - 1 / np.sqrt(2 * np.pi * ppf * self.params.gamma ** 2) - * np.exp(-((x - self.params.alpha) ** 2 + self.params.beta ** 2 * ppf ** 2) / (2 * ppf * self.params.gamma ** 2)) - ) - else: - return ( - 1 / np.sqrt(2 * np.pi * ppf) - * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu ** 2 * ppf ** 2) / (2 * ppf)) - ) + return (1/abs(self.params.gamma)) * norm.pdf((x - self.params.alpha - self.params.beta * p)/abs(self.params.gamma)) + return (1/abs(self.params.mu)) * norm.pdf((x - p)/abs(self.params.mu)) - integrator = self.integrator_cls(**(self.integrator_params or {})) - result = integrator.compute(integrand) - if self.mixture_form == "classical": - val = np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2) * result.value - else: - val = np.exp(self.params.mu * (x - self.params.alpha)) * result.value - return val, result.error + res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + return res.value, res.error - def _compute_logpdf(self, x: float) -> tuple[float, float]: + def _compute_logpdf(self, x: float) -> Tuple[float, float]: def integrand(u: float) -> float: - ppf = self.params.distribution.ppf(u) + p = self.params.distribution.ppf(u) if self.mixture_form == "classical": - return -( - (x - self.params.alpha) ** 2 - + ppf ** 2 * self.params.beta ** 2 - + ppf * self.params.gamma ** 2 * np.log(2 * np.pi * ppf * self.params.gamma ** 2) - ) / (2 * ppf * self.params.gamma ** 2) - else: - return -((x - self.params.alpha) ** 2 + ppf ** 2 * self.params.mu ** 2 + ppf * np.log(2 * np.pi * ppf)) / (2 * ppf) + return np.log(1/abs(self.params.gamma)) + norm.logpdf((x - self.params.alpha - self.params.beta * p)/abs(self.params.gamma)) + return np.log(1/abs(self.params.mu)) + norm.logpdf((x - p)/abs(self.params.mu)) - integrator = self.integrator_cls(**(self.integrator_params or {})) - result = integrator.compute(integrand) - if self.mixture_form == "classical": - val = self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2 + result.value - else: - val = self.params.mu * (x - self.params.alpha) + result.value - return val, result.error - - @lru_cache() - def _integrand_func(self, u: float, d: float, gamma: float) -> float: - ppf = self.params.distribution.ppf(u) - return (1 / np.sqrt(2 * np.pi * ppf * abs(gamma) ** 2)) * np.exp(-d / (2 * ppf)) + res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + return res.value, res.error diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index c53a67a..2d7cbad 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -35,12 +35,13 @@ def __init__( integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: - super().__init__(mixture_form, **kwargs) + super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) self.integrator_cls = integrator_cls self.integrator_params = integrator_params or {} - def _compute_moment(self, n: int, __: Dict[str, Any]) -> tuple[float, float]: + def _compute_moment(self, n: int) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) + def integrand(u: float) -> float: return sum( binom(n, k) @@ -50,34 +51,40 @@ def integrand(u: float) -> float: * norm.moment(k) for k in range(n + 1) ) + integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error - def _compute_cdf(self, x: float, __: Dict[str, Any]) -> tuple[float, float]: + def _compute_cdf(self, x: float) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) param_norm = norm(0, gamma) + def integrand(u: float) -> float: return param_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))) + integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error - def _compute_pdf(self, x: float, __: Dict[str, Any]) -> tuple[float, float]: + def _compute_pdf(self, x: float) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 + def integrand(u: float) -> float: return self._integrand_func(u, d, gamma) + integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error - def _compute_logpdf(self, x: float, __: Dict[str, Any]) -> tuple[float, float]: + def _compute_logpdf(self, x: float) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 + def integrand(u: float) -> float: return self._log_integrand_func(u, d, gamma) - # For log-pdf you may choose LogRQMC or any integrator that supports it + integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error From 1b5138da62526ec664e45abcefc0b29997eb195b Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 01:48:05 +0300 Subject: [PATCH 17/24] fixed --- src/mixtures/nm_mixture.py | 60 +++++++++++++++++++------------------- 1 file changed, 30 insertions(+), 30 deletions(-) diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 54498ab..22d9ac5 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -1,5 +1,5 @@ from dataclasses import dataclass -from typing import Any, Type, Dict, Tuple, Union, List +from typing import Any, Type, Dict, Tuple import numpy as np from scipy.special import binom @@ -7,15 +7,12 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.integrator import Integrator -from src.algorithms.support_algorithms.quad_integrator import QuadIntegrator -from src.algorithms.support_algorithms.rqmc import RQMCIntegrator -from src.algorithms.support_algorithms.log_rqmc import LogRQMC from src.mixtures.abstract_mixture import AbstractMixtures @dataclass class _NMMClassicDataCollector: alpha: float | int | np.int64 - beta: float | int | np.int64 + beta: float | int | np.int64 gamma: float | int | np.int64 distribution: rv_frozen | rv_continuous @@ -31,7 +28,7 @@ class NormalMeanMixtures(AbstractMixtures): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator] = RQMCIntegrator, + integrator_cls: Type[Integrator] = Integrator, integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: @@ -51,53 +48,56 @@ def _compute_moment(self, n: int) -> Tuple[float, float]: if self.mixture_form == "classical": for k in range(n + 1): for l in range(k + 1): - coeff = binom(n, n - k) * binom(k, k - l) * (self.params.beta ** (k - l)) * (self.params.gamma ** l) + coeff = binom(n, n - k) * binom(k, k - l) def mix(u: float) -> float: - return self.params.distribution.ppf(u) ** (k - l) - integrator = self.integrator_cls(**(self.integrator_params or {})) - res = integrator.compute(mix) - mixture_moment += coeff * (self.params.alpha ** (n - k)) * res.value * norm.moment(l) - error += coeff * res.error * (self.params.alpha ** (n - k)) * norm.moment(l) + return ( + self.params.distribution.ppf(u) ** (k - l) + ) + res = self.integrator_cls(**(self.integrator_params or {})).compute(mix) + mixture_moment += coeff * (self.params.beta ** (k - l)) * (self.params.gamma ** l) * (self.params.alpha ** (n - k)) * res.value * norm.moment(l) + error += coeff * (self.params.beta ** (k - l)) * (self.params.gamma ** l) * (self.params.alpha ** (n - k)) * res.error * norm.moment(l) else: for k in range(n + 1): - coeff = binom(n, n - k) * (self.params.sigma ** k) + coeff = binom(n, n - k) def mix(u: float) -> float: return self.params.distribution.ppf(u) ** (n - k) - integrator = self.integrator_cls(**(self.integrator_params or {})) - res = integrator.compute(mix) - mixture_moment += coeff * res.value * norm.moment(k) - error += coeff * res.error * norm.moment(k) + res = self.integrator_cls(**(self.integrator_params or {})).compute(mix) + mixture_moment += coeff * (self.params.sigma ** k) * res.value * norm.moment(k) + error += coeff * (self.params.sigma ** k) * res.error * norm.moment(k) return mixture_moment, error - def _compute_cdf(self, x: float, params: Dict[str, Any]) -> Tuple[float, float]: + def _compute_cdf(self, x: float) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: - return norm.cdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) + p = self.params.distribution.ppf(u) + return norm.cdf((x - self.params.alpha - self.params.beta * p) / abs(self.params.gamma)) else: def fn(u: float) -> float: - return norm.cdf((x - self.params.distribution.ppf(u)) / abs(self.params.sigma)) - integrator = self.integrator_cls(**(self.integrator_params or {})) - res = integrator.compute(fn) + p = self.params.distribution.ppf(u) + return norm.cdf((x - p) / abs(self.params.sigma)) + res = self.integrator_cls(**(self.integrator_params or {})).compute(fn) return res.value, res.error def _compute_pdf(self, x: float) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: - return (1 / abs(self.params.gamma)) * norm.pdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) + p = self.params.distribution.ppf(u) + return (1 / abs(self.params.gamma)) * norm.pdf((x - self.params.alpha - self.params.beta * p) / abs(self.params.gamma)) else: def fn(u: float) -> float: - return (1 / abs(self.params.sigma)) * norm.pdf((x - self.params.distribution.ppf(u)) / abs(self.params.sigma)) - integrator = self.integrator_cls(**(self.integrator_params or {})) - res = integrator.compute(fn) + p = self.params.distribution.ppf(u) + return (1 / abs(self.params.sigma)) * norm.pdf((x - p) / abs(self.params.sigma)) + res = self.integrator_cls(**(self.integrator_params or {})).compute(fn) return res.value, res.error def _compute_logpdf(self, x: float) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: - return np.log(1 / abs(self.params.gamma)) + norm.logpdf((x - self.params.alpha - self.params.beta * self.params.distribution.ppf(u)) / abs(self.params.gamma)) + p = self.params.distribution.ppf(u) + return np.log(1 / abs(self.params.gamma)) + norm.logpdf((x - self.params.alpha - self.params.beta * p) / abs(self.params.gamma)) else: def fn(u: float) -> float: - return np.log(1 / abs(self.params.sigma)) + norm.logpdf((x - self.params.distribution.ppf(u)) / abs(self.params.sigma)) - integrator = self.integrator_cls(**(self.integrator_params or {})) - res = integrator.compute(fn) + p = self.params.distribution.ppf(u) + return np.log(1 / abs(self.params.sigma)) + norm.logpdf((x - p) / abs(self.params.sigma)) + res = self.integrator_cls(**(self.integrator_params or {})).compute(fn) return res.value, res.error From 5f7295c772351b256116ea020b23cb0be1c05b23 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 01:55:42 +0300 Subject: [PATCH 18/24] fixed --- src/mixtures/nm_mixture.py | 3 +- src/mixtures/nmv_mixture.py | 59 ++++++++++++++++++++++++------------- 2 files changed, 41 insertions(+), 21 deletions(-) diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 22d9ac5..a2eaf76 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -7,6 +7,7 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @dataclass @@ -28,7 +29,7 @@ class NormalMeanMixtures(AbstractMixtures): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator] = Integrator, + integrator_cls: Type[Integrator] = RQMCIntegrator, integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index ad997c2..bc0f575 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -1,4 +1,5 @@ from dataclasses import dataclass +from functools import lru_cache from typing import Any, Type, Dict, Tuple import numpy as np @@ -7,6 +8,8 @@ from scipy.stats.distributions import rv_frozen from src.algorithms.support_algorithms.integrator import Integrator +from src.algorithms.support_algorithms.rqmc import RQMCIntegrator +from src.algorithms.support_algorithms.log_rqmc import LogRQMC from src.mixtures.abstract_mixture import AbstractMixtures @dataclass @@ -29,35 +32,37 @@ class NormalMeanVarianceMixtures(AbstractMixtures): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator], - integrator_params: Dict[str, Any] | None = None, + integrator_cls: Type[Integrator] = RQMCIntegrator, + integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) def _compute_moment(self, n: int) -> Tuple[float, float]: + gamma = getattr(self.params, 'gamma', None) + def integrand(u: float) -> float: s = 0.0 for k in range(n + 1): for l in range(k + 1): - if self.mixture_form == "classical": - coef = binom(n, n - k) * binom(k, k - l) + if self.mixture_form == 'classical': term = ( - (self.params.beta ** (k - l)) + binom(n, n - k) + * binom(k, k - l) + * (self.params.beta ** (k - l)) * (self.params.gamma ** l) * (self.params.distribution.ppf(u) ** (k - l/2)) - * (self.params.alpha ** (n - k)) * norm.moment(l) ) else: - coef = binom(n, n - k) * binom(k, k - l) term = ( - (self.params.nu ** (k - l)) + binom(n, n - k) + * binom(k, k - l) + * (self.params.mu ** (k - l)) * (self.params.distribution.ppf(u) ** (k - l/2)) - * (self.params.alpha ** (n - k)) * norm.moment(l) ) - s += coef * term if self.mixture_form == "classical" else term + s += term return s res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) @@ -66,8 +71,8 @@ def integrand(u: float) -> float: def _compute_cdf(self, x: float) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) - if self.mixture_form == "classical": - return norm.cdf((x - self.params.alpha - self.params.beta * p) / abs(self.params.gamma)) + if self.mixture_form == 'classical': + return norm.cdf((x - self.params.alpha) / (np.sqrt(p) * self.params.gamma)) return norm.cdf((x - self.params.alpha) / np.sqrt(p) - self.params.mu * np.sqrt(p)) res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) @@ -76,19 +81,33 @@ def integrand(u: float) -> float: def _compute_pdf(self, x: float) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) - if self.mixture_form == "classical": - return (1/abs(self.params.gamma)) * norm.pdf((x - self.params.alpha - self.params.beta * p)/abs(self.params.gamma)) - return (1/abs(self.params.mu)) * norm.pdf((x - p)/abs(self.params.mu)) + if self.mixture_form == 'classical': + return ( + 1 / np.sqrt(2 * np.pi * p * self.params.gamma ** 2) + * np.exp(-((x - self.params.alpha) ** 2 + self.params.beta ** 2 * p ** 2) / (2 * p * self.params.gamma ** 2)) + ) + return ( + 1 / np.sqrt(2 * np.pi * p) + * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu ** 2 * p ** 2) / (2 * p)) + ) res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) - return res.value, res.error + if self.mixture_form == 'classical': + val = np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2) * res.value + else: + val = np.exp(self.params.mu * (x - self.params.alpha)) * res.value + return val, res.error def _compute_logpdf(self, x: float) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) - if self.mixture_form == "classical": - return np.log(1/abs(self.params.gamma)) + norm.logpdf((x - self.params.alpha - self.params.beta * p)/abs(self.params.gamma)) - return np.log(1/abs(self.params.mu)) + norm.logpdf((x - p)/abs(self.params.mu)) + if self.mixture_form == 'classical': + return -((x - self.params.alpha) ** 2 + p ** 2 * self.params.beta ** 2 + p * self.params.gamma ** 2 * np.log(2 * np.pi * p * self.params.gamma ** 2)) / (2 * p * self.params.gamma ** 2) + return -((x - self.params.alpha) ** 2 + p ** 2 * self.params.mu ** 2 + p * np.log(2 * np.pi * p)) / (2 * p) res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) - return res.value, res.error + if self.mixture_form == 'classical': + val = self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2 + res.value + else: + val = self.params.mu * (x - self.params.alpha) + res.value + return val, res.error From 8d481af33354d8d9f465dab5cec7f294501e0538 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 03:22:26 +0300 Subject: [PATCH 19/24] fixed --- .../nvm_semiparametric/g_estimation_given_mu_rqmc_based.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu_rqmc_based.py b/src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu_rqmc_based.py index d166970..2bf1687 100644 --- a/src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu_rqmc_based.py +++ b/src/procedures/semiparametric/nvm_semiparametric/g_estimation_given_mu_rqmc_based.py @@ -32,7 +32,7 @@ def v_sequence_default_value(n: float) -> float: INTEGRATION_LIMIT_DEFAULT_VALUE: int = 50 -class NMVEstimationDensityInvMTquadRQMCBased: +class SemiParametricGEstimationGivenMuRQMCBased: """Estimation of mixing density function g (xi density function) of NVM mixture represented in canonical form Y = alpha + mu*xi + sqrt(xi)*N, where alpha = 0 and mu is given. @@ -162,8 +162,9 @@ def second_v_integrand(self, v: float, x: float) -> np.ndarray: x_power = self.x_powers[x][idx] return (self.second_u_integrals[idx] * x_power) / gamma_val - def compute_integrals_for_x(self, x: float, integrator: Integrator = RQMCIntegrator()) -> float: + def compute_integrals_for_x(self, x: float) -> float: """Compute integrals using RQMC for v-integration.""" + integrator = RQMCIntegrator() first_integral = integrator.compute(func=lambda t: np.sum(self.first_v_integrand(t * self.v_value, x)) * self.v_value).value second_integral = integrator.compute(func=lambda t: np.sum(self.second_v_integrand(-t * self.v_value, x)) * self.v_value).value @@ -172,5 +173,5 @@ def compute_integrals_for_x(self, x: float, integrator: Integrator = RQMCIntegra return max(0.0, total.real) def algorithm(self, sample: np.ndarray) -> EstimateResult: - y_data = [self.compute_integrals_for_x(x, Integrator()) for x in self.x_data] + y_data = [self.compute_integrals_for_x(x) for x in self.x_data] return EstimateResult(list_value=y_data, success=True) From b1c944c063d331b960781b5186fe7132450573bf Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 03:26:38 +0300 Subject: [PATCH 20/24] fixed --- src/procedures/__init__.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index d9f8503..18c1af9 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -17,7 +17,7 @@ NMVEstimationDensityInvMTquad, ) from src.procedures.semiparametric.nvm_semiparametric.g_estimation_given_mu_rqmc_based import ( - NMVEstimationDensityInvMTquadRQMCBased, + SemiParametricGEstimationGivenMuRQMCBased, ) from src.procedures.semiparametric.nvm_semiparametric.g_estimation_post_widder import ( NMVEstimationDensityPW, From d4e7c0631bed0a9343e20b8ebd98f5be9b1ab1bd Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 03:52:04 +0300 Subject: [PATCH 21/24] fixed --- src/procedures/__init__.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/procedures/__init__.py b/src/procedures/__init__.py index 18c1af9..00d3e1b 100644 --- a/src/procedures/__init__.py +++ b/src/procedures/__init__.py @@ -39,7 +39,7 @@ NMEstimationSigmaLaplace ) ALGORITHM_REGISTRY.register("density_estim_inv_mellin_rqmc_int", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( - NMVEstimationDensityInvMTquadRQMCBased + SemiParametricGEstimationGivenMuRQMCBased ) ALGORITHM_REGISTRY.register("density_estim_post_widder", AlgorithmPurpose.NMV_SEMIPARAMETRIC)( NMVEstimationDensityPW From cb90a631c6443838fcbc44e3cea4b63aaad720c1 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Thu, 26 Jun 2025 04:03:12 +0300 Subject: [PATCH 22/24] fixed --- README.md | 2 +- src/procedures/support/rqmc.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index 3b71d37..de58717 100644 --- a/README.md +++ b/README.md @@ -120,7 +120,7 @@ Computes given integral with the limits of integration being 0 and 1 using Rando ### Usage: ```Python -from src.algorithms.support_algorithms.rqmc import RQMC +from src.procedures.support.rqmc import RQMC rqmc = RQMC(lambda x: x**3 - x**2 + 1, error_tolerance=1e-5) ``` So `rqmc()[0]` is estimated integral value and `rqmc()[1]` is current error tolerance. diff --git a/src/procedures/support/rqmc.py b/src/procedures/support/rqmc.py index 543d2cb..4254615 100644 --- a/src/procedures/support/rqmc.py +++ b/src/procedures/support/rqmc.py @@ -5,7 +5,7 @@ import scipy from numba import njit -from src.algorithms.support_algorithms.integrator import IntegrationResult +from src.procedures.support.integrator import IntegrationResult BITS = 30 """Number of bits in XOR. Should be less than 64""" From de0d635cb270119837a6da27dbfcbab126dca1f8 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Fri, 27 Jun 2025 04:45:54 +0300 Subject: [PATCH 23/24] fixed --- src/mixtures/abstract_mixture.py | 49 ++++++++++++++++---------------- src/mixtures/nm_mixture.py | 24 ++++++++-------- src/mixtures/nmv_mixture.py | 20 ++++++------- src/mixtures/nv_mixture.py | 12 +++----- 4 files changed, 50 insertions(+), 55 deletions(-) diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index fea8edb..abf5797 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -19,8 +19,6 @@ class AbstractMixtures(metaclass=ABCMeta): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator] = RQMCIntegrator, - integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: """ @@ -31,8 +29,7 @@ def __init__( **kwargs: Parameters of Mixture (alpha, gamma, etc.) """ self.mixture_form = mixture_form - self.integrator_cls = integrator_cls - self.integrator_params = integrator_params or {} + if mixture_form == "classical": self.params = self._params_validation(self._classical_collector, kwargs) @@ -42,70 +39,74 @@ def __init__( raise AssertionError(f"Unknown mixture form: {mixture_form}") @abstractmethod - def _compute_moment(self, n: int) -> Tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator) -> Tuple[float, float]: ... def compute_moment( self, - x: Union[List[int], int, NDArray[np.float64]] + x: Union[List[int], int, NDArray[np.float64]], + integrator: Integrator ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_moment(xp) for xp in x], dtype=object) + return np.array([self._compute_moment(xp, integrator) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_moment(xp) for xp in x] + return [self._compute_moment(xp, integrator) for xp in x] elif isinstance(x, int): - return self._compute_moment(x) + return self._compute_moment(x, integrator) else: raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def _compute_pdf(self, x: float) -> Tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator) -> Tuple[float, float]: ... def compute_pdf( self, - x: Union[List[float], float, NDArray[np.float64]] + x: Union[List[float], float, NDArray[np.float64]], + integrator: Integrator ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_pdf(xp) for xp in x], dtype=object) + return np.array([self._compute_pdf(xp, integrator) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_pdf(xp) for xp in x] + return [self._compute_pdf(xp, integrator) for xp in x] elif isinstance(x, float): - return self._compute_pdf(x) + return self._compute_pdf(x, integrator) else: raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def _compute_logpdf(self, x: float) -> Tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator) -> Tuple[float, float]: ... def compute_logpdf( self, - x: Union[List[float], float, NDArray[np.float64]] + x: Union[List[float], float, NDArray[np.float64]], + integrator: Integrator ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_logpdf(xp) for xp in x], dtype=object) + return np.array([self._compute_logpdf(xp, integrator) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_logpdf(xp) for xp in x] + return [self._compute_logpdf(xp, integrator) for xp in x] elif isinstance(x, float): - return self._compute_logpdf(x) + return self._compute_logpdf(x, integrator) else: raise TypeError(f"Unsupported type for x: {type(x)}") @abstractmethod - def _compute_cdf(self, x: float) -> Tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator) -> Tuple[float, float]: ... def compute_cdf( self, - x: Union[List[float], float, NDArray[np.float64]] + x: Union[List[float], float, NDArray[np.float64]], + integrator: Integrator ) -> Union[List[Tuple[float, float]], Tuple[float, float], NDArray[Any]]: if isinstance(x, np.ndarray): - return np.array([self._compute_cdf(xp) for xp in x], dtype=object) + return np.array([self._compute_cdf(xp, integrator) for xp in x], dtype=object) elif isinstance(x, list): - return [self._compute_cdf(xp) for xp in x] + return [self._compute_cdf(xp, integrator) for xp in x] elif isinstance(x, float): - return self._compute_cdf(x) + return self._compute_cdf(x, integrator) else: raise TypeError(f"Unsupported type for x: {type(x)}") diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index 54e7e98..db84538 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -9,6 +9,7 @@ from src.procedures.support.integrator import Integrator from src.procedures.support.rqmc import RQMCIntegrator from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @dataclass @@ -30,11 +31,9 @@ class NormalMeanMixtures(AbstractMixtures): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator] = RQMCIntegrator, - integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: - super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) + super().__init__(mixture_form, **kwargs) def _params_validation(self, data_collector: Any, params: dict[str, float | rv_continuous | rv_frozen]) -> Any: data_class = super()._params_validation(data_collector, params) @@ -44,7 +43,7 @@ def _params_validation(self, data_collector: Any, params: dict[str, float | rv_c raise ValueError("Gamma can't be zero") return data_class - def _compute_moment(self, n: int) -> Tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator=QuadIntegrator) -> Tuple[float, float]: mixture_moment = 0.0 error = 0.0 if self.mixture_form == "classical": @@ -55,7 +54,8 @@ def mix(u: float) -> float: return ( self.params.distribution.ppf(u) ** (k - l) ) - res = self.integrator_cls(**(self.integrator_params or {})).compute(mix) + + res = integrator.compute(mix) mixture_moment += coeff * (self.params.beta ** (k - l)) * (self.params.gamma ** l) * (self.params.alpha ** (n - k)) * res.value * norm.moment(l) error += coeff * (self.params.beta ** (k - l)) * (self.params.gamma ** l) * (self.params.alpha ** (n - k)) * res.error * norm.moment(l) else: @@ -63,12 +63,12 @@ def mix(u: float) -> float: coeff = binom(n, n - k) def mix(u: float) -> float: return self.params.distribution.ppf(u) ** (n - k) - res = self.integrator_cls(**(self.integrator_params or {})).compute(mix) + res = integrator.compute(mix) mixture_moment += coeff * (self.params.sigma ** k) * res.value * norm.moment(k) error += coeff * (self.params.sigma ** k) * res.error * norm.moment(k) return mixture_moment, error - def _compute_cdf(self, x: float) -> Tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: p = self.params.distribution.ppf(u) @@ -77,10 +77,10 @@ def fn(u: float) -> float: def fn(u: float) -> float: p = self.params.distribution.ppf(u) return norm.cdf((x - p) / abs(self.params.sigma)) - res = self.integrator_cls(**(self.integrator_params or {})).compute(fn) + res = integrator.compute(fn) return res.value, res.error - def _compute_pdf(self, x: float) -> Tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: p = self.params.distribution.ppf(u) @@ -89,10 +89,10 @@ def fn(u: float) -> float: def fn(u: float) -> float: p = self.params.distribution.ppf(u) return (1 / abs(self.params.sigma)) * norm.pdf((x - p) / abs(self.params.sigma)) - res = self.integrator_cls(**(self.integrator_params or {})).compute(fn) + res = integrator.compute(fn) return res.value, res.error - def _compute_logpdf(self, x: float) -> Tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMC) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: p = self.params.distribution.ppf(u) @@ -101,5 +101,5 @@ def fn(u: float) -> float: def fn(u: float) -> float: p = self.params.distribution.ppf(u) return np.log(1 / abs(self.params.sigma)) + norm.logpdf((x - p) / abs(self.params.sigma)) - res = self.integrator_cls(**(self.integrator_params or {})).compute(fn) + res = integrator.compute(fn) return res.value, res.error diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index 0422b08..41aca43 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -32,13 +32,11 @@ class NormalMeanVarianceMixtures(AbstractMixtures): def __init__( self, mixture_form: str, - integrator_cls: Type[Integrator] = RQMCIntegrator, - integrator_params: Dict[str, Any] = None, **kwargs: Any ) -> None: - super().__init__(mixture_form, integrator_cls=integrator_cls, integrator_params=integrator_params, **kwargs) + super().__init__(mixture_form, **kwargs) - def _compute_moment(self, n: int) -> Tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: gamma = getattr(self.params, 'gamma', None) def integrand(u: float) -> float: @@ -65,20 +63,20 @@ def integrand(u: float) -> float: s += term return s - res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + res = integrator.compute(integrand) return res.value, res.error - def _compute_cdf(self, x: float) -> Tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) if self.mixture_form == 'classical': return norm.cdf((x - self.params.alpha) / (np.sqrt(p) * self.params.gamma)) return norm.cdf((x - self.params.alpha) / np.sqrt(p) - self.params.mu * np.sqrt(p)) - res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + res = integrator.compute(integrand) return res.value, res.error - def _compute_pdf(self, x: float) -> Tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) if self.mixture_form == 'classical': @@ -91,21 +89,21 @@ def integrand(u: float) -> float: * np.exp(-((x - self.params.alpha) ** 2 + self.params.mu ** 2 * p ** 2) / (2 * p)) ) - res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + res = integrator.compute(integrand) if self.mixture_form == 'classical': val = np.exp(self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2) * res.value else: val = np.exp(self.params.mu * (x - self.params.alpha)) * res.value return val, res.error - def _compute_logpdf(self, x: float) -> Tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMC) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) if self.mixture_form == 'classical': return -((x - self.params.alpha) ** 2 + p ** 2 * self.params.beta ** 2 + p * self.params.gamma ** 2 * np.log(2 * np.pi * p * self.params.gamma ** 2)) / (2 * p * self.params.gamma ** 2) return -((x - self.params.alpha) ** 2 + p ** 2 * self.params.mu ** 2 + p * np.log(2 * np.pi * p)) / (2 * p) - res = self.integrator_cls(**(self.integrator_params or {})).compute(integrand) + res = integrator.compute(integrand) if self.mixture_form == 'classical': val = self.params.beta * (x - self.params.alpha) / self.params.gamma ** 2 + res.value else: diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index 97a4eab..3d991e2 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -40,7 +40,7 @@ def __init__( self.integrator_cls = integrator_cls self.integrator_params = integrator_params or {} - def _compute_moment(self, n: int) -> tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator=QuadIntegrator) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) def integrand(u: float) -> float: @@ -53,40 +53,36 @@ def integrand(u: float) -> float: for k in range(n + 1) ) - integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error - def _compute_cdf(self, x: float) -> tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator=QuadIntegrator) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) param_norm = norm(0, gamma) def integrand(u: float) -> float: return param_norm.cdf((x - self.params.alpha) / np.sqrt(self.params.distribution.ppf(u))) - integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error - def _compute_pdf(self, x: float) -> tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator=QuadIntegrator) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 def integrand(u: float) -> float: return self._integrand_func(u, d, gamma) - integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error - def _compute_logpdf(self, x: float) -> tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMC) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 def integrand(u: float) -> float: return self._log_integrand_func(u, d, gamma) - integrator = self.integrator_cls(**self.integrator_params) result = integrator.compute(integrand) return result.value, result.error From 34a76ea166776684a096349827152d55275dd534 Mon Sep 17 00:00:00 2001 From: plidan123 Date: Fri, 27 Jun 2025 05:09:14 +0300 Subject: [PATCH 24/24] fixed --- src/mixtures/abstract_mixture.py | 1 - src/mixtures/nm_mixture.py | 11 +-- src/mixtures/nmv_mixture.py | 11 +-- src/mixtures/nv_mixture.py | 11 +-- src/procedures/support/log_rqmc.py | 122 ++++++++++------------------- 5 files changed, 59 insertions(+), 97 deletions(-) diff --git a/src/mixtures/abstract_mixture.py b/src/mixtures/abstract_mixture.py index abf5797..de31713 100644 --- a/src/mixtures/abstract_mixture.py +++ b/src/mixtures/abstract_mixture.py @@ -8,7 +8,6 @@ from scipy.stats.distributions import rv_frozen from src.procedures.support.integrator import Integrator -from src.procedures.support.rqmc import RQMCIntegrator # default integrator class AbstractMixtures(metaclass=ABCMeta): """Base class for Mixtures""" diff --git a/src/mixtures/nm_mixture.py b/src/mixtures/nm_mixture.py index db84538..1f1e15f 100644 --- a/src/mixtures/nm_mixture.py +++ b/src/mixtures/nm_mixture.py @@ -8,7 +8,8 @@ from src.procedures.support.integrator import Integrator from src.procedures.support.rqmc import RQMCIntegrator -from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.rqmc import RQMC +from src.procedures.support.log_rqmc import LogRQMCIntegrator from src.procedures.support.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @@ -43,7 +44,7 @@ def _params_validation(self, data_collector: Any, params: dict[str, float | rv_c raise ValueError("Gamma can't be zero") return data_class - def _compute_moment(self, n: int, integrator: Integrator=QuadIntegrator) -> Tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator=QuadIntegrator()) -> Tuple[float, float]: mixture_moment = 0.0 error = 0.0 if self.mixture_form == "classical": @@ -68,7 +69,7 @@ def mix(u: float) -> float: error += coeff * (self.params.sigma ** k) * res.error * norm.moment(k) return mixture_moment, error - def _compute_cdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator=RQMCIntegrator()) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: p = self.params.distribution.ppf(u) @@ -80,7 +81,7 @@ def fn(u: float) -> float: res = integrator.compute(fn) return res.value, res.error - def _compute_pdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator=RQMCIntegrator()) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: p = self.params.distribution.ppf(u) @@ -92,7 +93,7 @@ def fn(u: float) -> float: res = integrator.compute(fn) return res.value, res.error - def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMC) -> Tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMCIntegrator()) -> Tuple[float, float]: if self.mixture_form == "classical": def fn(u: float) -> float: p = self.params.distribution.ppf(u) diff --git a/src/mixtures/nmv_mixture.py b/src/mixtures/nmv_mixture.py index 41aca43..6bd4dcb 100644 --- a/src/mixtures/nmv_mixture.py +++ b/src/mixtures/nmv_mixture.py @@ -9,7 +9,8 @@ from src.procedures.support.integrator import Integrator from src.procedures.support.rqmc import RQMCIntegrator -from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.rqmc import RQMC +from src.procedures.support.log_rqmc import LogRQMCIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @dataclass @@ -36,7 +37,7 @@ def __init__( ) -> None: super().__init__(mixture_form, **kwargs) - def _compute_moment(self, n: int, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator=RQMCIntegrator()) -> Tuple[float, float]: gamma = getattr(self.params, 'gamma', None) def integrand(u: float) -> float: @@ -66,7 +67,7 @@ def integrand(u: float) -> float: res = integrator.compute(integrand) return res.value, res.error - def _compute_cdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator=RQMCIntegrator()) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) if self.mixture_form == 'classical': @@ -76,7 +77,7 @@ def integrand(u: float) -> float: res = integrator.compute(integrand) return res.value, res.error - def _compute_pdf(self, x: float, integrator: Integrator=RQMCIntegrator) -> Tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator=RQMCIntegrator()) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) if self.mixture_form == 'classical': @@ -96,7 +97,7 @@ def integrand(u: float) -> float: val = np.exp(self.params.mu * (x - self.params.alpha)) * res.value return val, res.error - def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMC) -> Tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMCIntegrator()) -> Tuple[float, float]: def integrand(u: float) -> float: p = self.params.distribution.ppf(u) if self.mixture_form == 'classical': diff --git a/src/mixtures/nv_mixture.py b/src/mixtures/nv_mixture.py index 3d991e2..20f75ae 100644 --- a/src/mixtures/nv_mixture.py +++ b/src/mixtures/nv_mixture.py @@ -7,10 +7,11 @@ from scipy.stats import norm, rv_continuous from scipy.stats.distributions import rv_frozen -from src.procedures.support.log_rqmc import LogRQMC +from src.procedures.support.log_rqmc import LogRQMCIntegrator from src.procedures.support.integrator import Integrator from src.procedures.support.rqmc import RQMCIntegrator +from src.procedures.support.rqmc import RQMC from src.procedures.support.quad_integrator import QuadIntegrator from src.mixtures.abstract_mixture import AbstractMixtures @@ -40,7 +41,7 @@ def __init__( self.integrator_cls = integrator_cls self.integrator_params = integrator_params or {} - def _compute_moment(self, n: int, integrator: Integrator=QuadIntegrator) -> tuple[float, float]: + def _compute_moment(self, n: int, integrator: Integrator=QuadIntegrator()) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) def integrand(u: float) -> float: @@ -56,7 +57,7 @@ def integrand(u: float) -> float: result = integrator.compute(integrand) return result.value, result.error - def _compute_cdf(self, x: float, integrator: Integrator=QuadIntegrator) -> tuple[float, float]: + def _compute_cdf(self, x: float, integrator: Integrator=QuadIntegrator()) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) param_norm = norm(0, gamma) @@ -66,7 +67,7 @@ def integrand(u: float) -> float: result = integrator.compute(integrand) return result.value, result.error - def _compute_pdf(self, x: float, integrator: Integrator=QuadIntegrator) -> tuple[float, float]: + def _compute_pdf(self, x: float, integrator: Integrator=QuadIntegrator()) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 @@ -76,7 +77,7 @@ def integrand(u: float) -> float: result = integrator.compute(integrand) return result.value, result.error - def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMC) -> tuple[float, float]: + def _compute_logpdf(self, x: float, integrator: Integrator=LogRQMCIntegrator()) -> tuple[float, float]: gamma = getattr(self.params, 'gamma', 1) d = (x - self.params.alpha) ** 2 / gamma ** 2 diff --git a/src/procedures/support/log_rqmc.py b/src/procedures/support/log_rqmc.py index a10f221..af267b0 100644 --- a/src/procedures/support/log_rqmc.py +++ b/src/procedures/support/log_rqmc.py @@ -1,121 +1,81 @@ from typing import Callable - import numpy as np import numpy._typing as tpg import scipy - +from src.procedures.support.integrator import IntegrationResult from src.procedures.support.rqmc import RQMC -class LogRQMC(RQMC): +class LogRQMCIntegrator: + """Log-sum-exp stabilized randomized Quasi-Monte Carlo integrator.""" + def __init__( self, - func: Callable, error_tolerance: float = 1e-6, count: int = 25, - base_n: int = 2**6, + base_n: int = 2 ** 6, i_max: int = 100, a: float = 0.00047, ): - super().__init__(func, error_tolerance, count, base_n, i_max, a) + self.error_tolerance = error_tolerance + self.count = count + self.base_n = base_n + self.i_max = i_max + self.a = a @staticmethod def lse(args: tpg.NDArray) -> float: """ - Compute LSE - Args: - args (): Values - - Returns: LSE result - + Log-Sum-Exp stabilization helper. """ max_value = max(args) return max_value + np.log(np.sum(np.exp(args - max_value))) - def _independent_estimator(self, values: np._typing.NDArray) -> float: - """Apply function to row of matrix and find mean of row - - Args: - values: row of random values matrix - - Returns: mean of row - - """ + def _independent_estimator(self, values: tpg.NDArray) -> float: vfunc = np.vectorize(self.func) return -np.log(len(values)) + self.lse(vfunc(values)) - def _estimator(self, random_matrix: np._typing.NDArray) -> tuple[float, np._typing.NDArray]: - """Find mean of all rows - - Args: - random_matrix: matrix of random values - - Returns: mean of all rows - - """ + def _estimator(self, random_matrix: tpg.NDArray) -> tuple[float, tpg.NDArray]: values = np.array(list(map(self._independent_estimator, random_matrix))) return -np.log(self.count) + self.lse(values), values - def _update_independent_estimator(self, i: int, old_value: float, new_values: np._typing.NDArray) -> float: - """Update mean of row - - Args: - i: step count - old_value: previous value of row on i-1 step - new_values: new generated row of random values - - Returns: Updated mean of row - - """ - + def _update_independent_estimator( + self, i: int, old_value: float, new_values: tpg.NDArray + ) -> float: return -np.log(i + 1) + self.lse( np.array(i * [old_value] + [self._independent_estimator(new_values[i * self.base_n :])]) ) def _update( - self, j: int, old_values: np._typing.NDArray, random_matrix: np._typing.NDArray - ) -> tuple[float, np._typing.NDArray]: - """Update mean of all rows - - Args: - j: step count - old_values: previous values of row on i-1 step - random_matrix: new generated matrix of random values - - Returns:Updated mean of all rows - - """ + self, j: int, old_values: tpg.NDArray, random_matrix: tpg.NDArray + ) -> tuple[float, tpg.NDArray]: values = [] - for i in range(self.count): - old_value, new_values = old_values[i], random_matrix[i] - value = self._update_independent_estimator(j, old_value, new_values) - values.append(value) + for idx in range(self.count): + old_val, new_vals = old_values[idx], random_matrix[idx] + values.append(self._update_independent_estimator(j, old_val, new_vals)) np_values = np.array(values) return -np.log(self.count) + self.lse(np_values), np_values - def rqmc(self) -> tuple[float, float]: - """Main function of algorithm - - Returns: approximation for integral of function from 0 to 1 - - """ - sobol_sampler = scipy.stats.qmc.Sobol(d=1, scramble=False) - sobol_sample = np.repeat(sobol_sampler.random(self.base_n).transpose(), self.count, axis=0) - xor_sample = np.array(np.random.rand(1, self.count)[0]) - sample = self._digital_shift(sobol_sample, xor_sample) - approximation, values = self._estimator(sample) - current_error_tolerance = self._sigma(values, approximation) * self.z - for i in range(1, self.i_max): - if current_error_tolerance < self.error_tolerance: - return approximation, current_error_tolerance - sobol_sampler.reset() - sobol_sample = np.repeat(sobol_sampler.random(self.base_n * (i + 1)).transpose(), self.count, axis=0) - sample = self._digital_shift(sobol_sample, xor_sample) - approximation, values = self._update(i, values, sample) - current_error_tolerance = self._sigma(values, approximation) * self.z - return approximation, current_error_tolerance + def compute(self, func: Callable) -> IntegrationResult: + """Compute integral of `func` over [0,1] with log-RQMC.""" + # Assign the integrand + self.func = func + # Instantiate underlying RQMC with same parameters + rqmc_inst = RQMC( + func, + error_tolerance=self.error_tolerance, + count=self.count, + base_n=self.base_n, + i_max=self.i_max, + a=self.a, + ) + # __call__ returns (approximation, error) + approximation, error = rqmc_inst() + return IntegrationResult(approximation, error, None) if __name__ == "__main__": - log_rqmc = LogRQMC(lambda x: x**3 - x**2 + 1, i_max=100) - print(log_rqmc()) + # Sample usage + integrator = LogRQMCIntegrator(i_max=100) + result = integrator.compute(lambda x: x**3 - x**2 + 1) + print(result)