..

dweindl · dweindl · commit 0b704ba7d9e0 · 2024-12-10T18:20:10.000+01:00
diff --git a/doc/example/distributions.ipynb b/doc/example/distributions.ipynb
@@ -177,20 +177,11 @@
    "source": [
     "plot(Prior(NORMAL, (10, 1), bounds=(6, 14), transformation=\"log10\"))\n",
     "plot(Prior(PARAMETER_SCALE_NORMAL, (10, 1), bounds=(10**6, 10**14), transformation=\"log10\"))\n",
-    "plot(Prior(LAPLACE, (10, 2), bounds=(6, 14)))\n",
-    "\n"
+    "plot(Prior(LAPLACE, (10, 2), bounds=(6, 14)))"
    ],
    "id": "581e1ac431860419",
    "outputs": [],
    "execution_count": null
-  },
-  {
-   "metadata": {},
-   "cell_type": "code",
-   "source": "",
-   "id": "633733651bbc3ef0",
-   "outputs": [],
-   "execution_count": null
   }
  ],
  "metadata": {
diff --git a/petab/v1/distributions.py b/petab/v1/distributions.py
@@ -4,35 +4,76 @@
 import abc
 
 import numpy as np
-from scipy.stats import laplace, lognorm, loguniform, norm, uniform
+from scipy.stats import laplace, norm, uniform
 
 __all__ = [
     "Distribution",
     "Normal",
-    "LogNormal",
     "Uniform",
-    "LogUniform",
     "Laplace",
-    "LogLaplace",
 ]
 
 
 class Distribution(abc.ABC):
     """A univariate probability distribution."""
 
-    @abc.abstractmethod
+    def __init__(self, log: bool | float = False):
+        if log is True:
+            log = np.exp(1)
+        self._log = log
+
+    def _undo_log(self, x: np.ndarray | float) -> np.ndarray | float:
+        """Undo the log transformation.
+
+        :param x: The sample to transform.
+        :return: The transformed sample
+        """
+        if self._log is False:
+            return x
+        return self._log**x
+
+    def _apply_log(self, x: np.ndarray | float) -> np.ndarray | float:
+        """Apply the log transformation.
+
+        :param x: The value to transform.
+        :return: The transformed value.
+        """
+        if self._log is False:
+            return x
+        return np.log(x) / np.log(self._log)
+
     def sample(self, shape=None) -> np.ndarray:
         """Sample from the distribution.
 
         :param shape: The shape of the sample.
         :return: A sample from the distribution.
         """
-        ...
+        sample = self._sample(shape)
+        return self._undo_log(sample)
 
     @abc.abstractmethod
+    def _sample(self, shape=None) -> np.ndarray:
+        """Sample from the underlying distribution.
+
+        :param shape: The shape of the sample.
+        :return: A sample from the underlying distribution,
+            before applying, e.g., the log transformation.
+        """
+        ...
+
     def pdf(self, x):
         """Probability density function at x.
 
+        :param x: The value at which to evaluate the PDF.
+        :return: The value of the PDF at ``x``.
+        """
+        chain_rule_factor = (1 / (x * np.log(self._log))) if self._log else 1
+        return self._pdf(self._apply_log(x)) * chain_rule_factor
+
+    @abc.abstractmethod
+    def _pdf(self, x):
+        """Probability density function of the underlying distribution at x.
+
         :param x: The value at which to evaluate the PDF.
         :return: The value of the PDF at ``x``.
         """
@@ -44,71 +85,39 @@ class Normal(Distribution):
 
     def __init__(
         self,
-        mean: float,
-        std: float,
+        loc: float,
+        scale: float,
         truncation: tuple[float, float] | None = None,
+        log: bool | float = False,
     ):
-        super().__init__()
-        self._mean = mean
-        self._std = std
+        super().__init__(log=log)
+        self._loc = loc
+        self._scale = scale
         self._truncation = truncation
 
         if truncation is not None:
             raise NotImplementedError("Truncation is not yet implemented.")
 
     def __repr__(self):
-        return (
-            f"Normal(mean={self._mean}, std={self._std}, "
-            f"truncation={self._truncation})"
-        )
-
-    def sample(self, shape=None):
-        return np.random.normal(loc=self._mean, scale=self._std, size=shape)
-
-    def pdf(self, x):
-        return norm.pdf(x, loc=self._mean, scale=self._std)
-
-
-class LogNormal(Distribution):
-    """A log-normal distribution.
-
-    :param mean: The mean of the underlying normal distribution.
-    :param std: The standard deviation of the underlying normal distribution.
-
-    """
-
-    def __init__(
-        self,
-        mean: float,
-        std: float,
-        truncation: tuple[float, float] | None = None,
-        base: float = np.exp(1),
-    ):
-        super().__init__()
-        self._mean = mean
-        self._std = std
-        self._truncation = truncation
-        self._base = base
+        trunc = f", truncation={self._truncation}" if self._truncation else ""
+        log = f", log={self._log}" if self._log else ""
+        return f"Normal(loc={self._loc}, scale={self._scale}{trunc}{log})"
 
-        if truncation is not None:
-            raise NotImplementedError("Truncation is not yet implemented.")
+    def _sample(self, shape=None):
+        return np.random.normal(loc=self._loc, scale=self._scale, size=shape)
 
-        if base != np.exp(1):
-            raise NotImplementedError("Only base e is supported.")
+    def _pdf(self, x):
+        return norm.pdf(x, loc=self._loc, scale=self._scale)
 
-    def __repr__(self):
-        return (
-            f"LogNormal(mean={self._mean}, std={self._std}, "
-            f"base={self._base}, truncation={self._truncation})"
-        )
-
-    def sample(self, shape=None):
-        return np.random.lognormal(
-            mean=self._mean, sigma=self._std, size=shape
-        )
+    @property
+    def loc(self):
+        """The location parameter of the underlying distribution."""
+        return self._loc
 
-    def pdf(self, x):
-        return lognorm.pdf(x, scale=np.exp(self._mean), s=self._std)
+    @property
+    def scale(self):
+        """The scale parameter of the underlying distribution."""
+        return self._scale
 
 
 class Uniform(Distribution):
@@ -118,58 +127,24 @@ def __init__(
         self,
         low: float,
         high: float,
+        *,
+        log: bool | float = False,
     ):
-        super().__init__()
+        super().__init__(log=log)
         self._low = low
         self._high = high
 
     def __repr__(self):
-        return f"Uniform(low={self._low}, high={self._high})"
+        log = f", log={self._log}" if self._log else ""
+        return f"Uniform(low={self._low}, high={self._high}{log})"
 
-    def sample(self, shape=None):
+    def _sample(self, shape=None):
         return np.random.uniform(low=self._low, high=self._high, size=shape)
 
-    def pdf(self, x):
+    def _pdf(self, x):
         return uniform.pdf(x, loc=self._low, scale=self._high - self._low)
 
 
-class LogUniform(Distribution):
-    """A log-uniform distribution.
-
-    :param low: The lower bound of the underlying normal distribution.
-    :param high: The upper bound of the underlying normal distribution.
-    """
-
-    def __init__(
-        self,
-        low: float,
-        high: float,
-        base: float = np.exp(1),
-    ):
-        super().__init__()
-        self._low = low
-        self._high = high
-        self._base = base
-        # re-scaled distribution parameters as required by
-        #  scipy.stats.loguniform
-        self._low_internal = np.exp(np.log(base) * low)
-        self._high_internal = np.exp(np.log(base) * high)
-
-    def __repr__(self):
-        return (
-            f"LogUniform(low={self._low}, high={self._high}, "
-            f"base={self._base})"
-        )
-
-    def sample(self, shape=None):
-        return loguniform.rvs(
-            self._low_internal, self._high_internal, size=shape
-        )
-
-    def pdf(self, x):
-        return loguniform.pdf(x, self._low_internal, self._high_internal)
-
-
 class Laplace(Distribution):
     """A Laplace distribution."""
 
@@ -178,65 +153,32 @@ def __init__(
         loc: float,
         scale: float,
         truncation: tuple[float, float] | None = None,
+        log: bool | float = False,
     ):
-        super().__init__()
+        super().__init__(log=log)
         self._loc = loc
         self._scale = scale
         self._truncation = truncation
         if truncation is not None:
             raise NotImplementedError("Truncation is not yet implemented.")
 
-    def sample(self, shape=None):
+    def __repr__(self):
+        trunc = f", truncation={self._truncation}" if self._truncation else ""
+        log = f", log={self._log}" if self._log else ""
+        return f"Laplace(loc={self._loc}, scale={self._scale}{trunc}{log})"
+
+    def _sample(self, shape=None):
         return np.random.laplace(loc=self._loc, scale=self._scale, size=shape)
 
-    def pdf(self, x):
+    def _pdf(self, x):
         return laplace.pdf(x, loc=self._loc, scale=self._scale)
 
-
-class LogLaplace(Distribution):
-    """A log-Laplace distribution."""
-
-    def __init__(
-        self,
-        loc: float,
-        scale: float,
-        truncation: tuple[float, float] | None = None,
-        base: float = np.exp(1),
-    ):
-        super().__init__()
-        self._loc = loc
-        self._scale = scale
-        self._truncation = truncation
-        self._base = base
-        if truncation is not None:
-            raise NotImplementedError("Truncation is not yet implemented.")
-        if base != np.exp(1):
-            raise NotImplementedError("Only base e is supported.")
-
-    def __repr__(self):
-        return (
-            f"LogLaplace(loc={self._loc}, scale={self._scale}, "
-            f"base={self._base}, truncation={self._truncation})"
-        )
-
     @property
     def loc(self):
-        """The mean of the underlying Laplace distribution."""
+        """The location parameter of the underlying distribution."""
         return self._loc
 
     @property
     def scale(self):
-        """The scale of the underlying Laplace distribution."""
+        """The scale parameter of the underlying distribution."""
         return self._scale
-
-    def sample(self, shape=None):
-        return np.exp(
-            np.random.laplace(loc=self._loc, scale=self._scale, size=shape)
-        )
-
-    def pdf(self, x):
-        return (
-            1
-            / (2 * self.scale * x)
-            * np.exp(-np.abs(np.log(x) - self._loc) / self._scale)
-        )
diff --git a/petab/v1/priors.py b/petab/v1/priors.py
@@ -97,21 +97,17 @@ def __init__(
             case (C.LAPLACE, _) | (C.PARAMETER_SCALE_LAPLACE, C.LIN):
                 self.distribution = Laplace(*parameters)
             case (C.PARAMETER_SCALE_UNIFORM, C.LOG):
-                self.distribution = LogUniform(*parameters)
+                self.distribution = Uniform(*parameters, log=True)
             case (C.LOG_NORMAL, _) | (C.PARAMETER_SCALE_NORMAL, C.LOG):
-                self.distribution = LogNormal(*parameters)
+                self.distribution = Normal(*parameters, log=True)
             case (C.LOG_LAPLACE, _) | (C.PARAMETER_SCALE_LAPLACE, C.LOG):
-                self.distribution = LogLaplace(*parameters)
+                self.distribution = Laplace(*parameters, log=True)
             case (C.PARAMETER_SCALE_UNIFORM, C.LOG10):
-                self.distribution = LogUniform(*parameters, base=10)
+                self.distribution = Uniform(*parameters, log=10)
             case (C.PARAMETER_SCALE_NORMAL, C.LOG10):
-                self.distribution = LogNormal(
-                    np.log(10) * parameters[0], np.log(10) * parameters[1]
-                )
+                self.distribution = Normal(*parameters, log=10)
             case (C.PARAMETER_SCALE_LAPLACE, C.LOG10):
-                self.distribution = LogLaplace(
-                    np.log(10) * parameters[0], np.log(10) * parameters[1]
-                )
+                self.distribution = Laplace(*parameters, log=10)
             case _:
                 raise ValueError(
                     "Unsupported distribution type / transformation: "
diff --git a/tests/v1/test_distributions.py b/tests/v1/test_distributions.py
@@ -14,12 +14,14 @@
     list(
         product(
             [
-                Normal(1, 1),
-                LogNormal(2, 1),
+                Normal(2, 1),
+                Normal(2, 1, log=True),
+                Normal(2, 1, log=10),
                 Uniform(2, 4),
-                LogUniform(1, 2),
+                Uniform(-2, 4, log=True),
+                Uniform(2, 4, log=10),
                 Laplace(1, 2),
-                LogLaplace(1, 0.5),
+                Laplace(1, 0.5, log=True),
             ],
             [LIN, LOG, LOG10],
         )
