Repository: CUQI-DTU/CUQIpy

Issue:
Add non-informative/unbounded-uniform distribution
**General description of task**

Add `UnboundedUniform` for non-informative prior.

**Definition of Done (Feature/change)**
- [x] Add `UnboundedUniform` which will be directly subclased from `Distribution`
- [x] Add unit tests for `UnboundedUniform` for evaluating logpdf and gradient, and also for sampling.

**Definition of Done (Mandatory)**
- [x] Documentation added (docstrings on all public methods/classes)
- [x] Unit tests added/updated (and passing!)
- [x] Code reviewed and approved by 2 CUQI-maintainers
- [x] Online documentation is rendered correctly (esp. math) 

**Definition of Done (Optional)**
- [x] HowTo/tutorial added

## The following is the original issue:
## Issue
Currently, a snippet like the following works OK to draw samples based only on the likelihood:
```python
y = cuqi.distribution.Gaussian(A, noise_std**2)
lkh = y(y=y_obs)
np.random.seed(0)
sampler_likelihood = cuqi.experimental.mcmc.NUTS(lkh)
sampler_likelihood.sample(100)
```

However, the code might be difficult to understand as we don't explicitly define an `x` while we are sampling it, as commented by @jeverink.

## Solutions
Possible solutions to reduce/remove the confusion:
1. Explicitly define prior through the existing `Uniform` class:
```python
x = cuqi.distribution.Uniform(low=[-np.inf,-np.inf], high=[np.inf,np.inf])
post = cuqi.distribution.JointDistribution(x,y)(y=y_obs)
```

2. Explicitly define prior through a new `UnboundedUniform`:
```python
x = cuqi.distribution.UnboundedUniform(geometry=A.domain_geometry)
post = cuqi.distribution.JointDistribution(x,y)(y=y_obs)
```
A separate class would help emphasize that 
- the distribution here is an improper distribution, instead of a "normal" `Uniform`;
- and overload "sample()" as we are not supposed to draw samples from an improper distribution.

## Question
I hesitated a bit about which option was better.

Discussion:
@chaozg: @amal-ghamdi @jakobsj, I would like to hear your opinion on this issue.@jakobsj: Thanks for drawing our attention to this @chaozg and @jeverink. 

I think it might make a lot of sense to add the suggested UnboundedUniform as suggested and probably with .sample() disabled as it is an implicit distribution (this should be in line with other implicit distributions).

As discussed this morning, my understanding is in the provided example that we are sampling the likelihood (not posterior). Sampling of the likelihood in CUQIpy is possible through likelihood being subclassed from density, similar to how it works for distributions. As such, the code as it stands ie without defining a prior distribution and not setting up a posterior, reflects that it is indeed sampling from the likelihood.

I think in line 2, the name `data_distb` adds confusion, since it is in fact the result of line 1 which is the data distribution, and at this point ie line 2, after being conditioned on the observed data, it is now the likelihood.

When talking about this, for example in the paper, one should be clear that we are sampling the likelihood, NOT a posterior (even if in some sense with some trivial prior, this would be equivalent).  If instead we want to talk about a posterior, then we would follow the normal workflow and set up the (implicit) distribution for x, form the posterior and sample from it.

It would in fact make for a nice demonstration to show on the example that samples from the resulting posterior are in agreement with samples from the likelihood. @chaozg: @jakobsj Could you please review the DoD here?@jakobsj: Please can you make the "Update unit test for Uniform" more specific (update how and why).  This may relate to the change I commented on in the PR that the logpdf value of Uniform is changed to 1 (unnormalised) from the normalised existing implementation. This need for this change is unclear to me and in any case the DoD should reflect if such a change is to be made. Let us discuss tomorrow.@chaozg: Hi @jakobsj , based on your review comments, I just subclassed `UnboundedUniform` directly from `Distribution`, and the DoDs here are also updated.

Solution patch:
diff --git a/cuqi/distribution/__init__.py b/cuqi/distribution/__init__.py
index 3a51f5472..4842e2cb4 100644
--- a/cuqi/distribution/__init__.py
+++ b/cuqi/distribution/__init__.py
@@ -14,6 +14,6 @@
 from ._normal import Normal
 from ._truncated_normal import TruncatedNormal
 from ._posterior import Posterior
-from ._uniform import Uniform
+from ._uniform import Uniform, UnboundedUniform
 from ._custom import UserDefinedDistribution, DistributionGallery
 from ._joint_distribution import JointDistribution, _StackedJointDistribution, MultipleLikelihoodPosterior
diff --git a/cuqi/distribution/_uniform.py b/cuqi/distribution/_uniform.py
index a61ae71e6..b1f3b60fc 100644
--- a/cuqi/distribution/_uniform.py
+++ b/cuqi/distribution/_uniform.py
@@ -1,5 +1,6 @@
 import numpy as np
 from cuqi.distribution import Distribution
+from cuqi.geometry import Geometry
 
 class Uniform(Distribution):
 
@@ -57,4 +58,37 @@ def _sample(self,N=1, rng=None):
         else:
             s = np.random.uniform(self.low, self.high, (N,self.dim)).T
 
-        return s
\ No newline at end of file
+        return s
+
+class UnboundedUniform(Distribution):
+    """
+    Unbounded uniform distribution. This is a special case of the
+    Uniform distribution, where the lower and upper bounds are set to
+    -inf and inf, respectively. This distribution is not normalizable,
+    and therefore cannot be sampled from. It is mainly used for
+    initializing non-informative priors.
+    Parameters
+    ----------
+    geometry : int or Geometry
+        The geometry of the distribution. If an integer is given, it is
+        interpreted as the dimension of the distribution. If a
+        Geometry object is given, its par_dim attribute is used.
+    """
+    def __init__(self, geometry, is_symmetric=True, **kwargs):
+        super().__init__(geometry=geometry, is_symmetric=is_symmetric, **kwargs) 
+
+    def logpdf(self, x):
+        """
+        Evaluate the logarithm of the unnormalized PDF at the given values of x.
+        """
+        # Always return 1.0 (the unnormalized log PDF)
+        return 1.0
+
+    def gradient(self, x):
+        """
+        Computes the gradient of logpdf at the given values of x.
+        """
+        return np.zeros_like(x)
+
+    def _sample(self, N=1, rng=None):
+        raise NotImplementedError("Cannot sample from UnboundedUniform distribution")
\ No newline at end of file


Test patch:
diff --git a/tests/test_distribution.py b/tests/test_distribution.py
index 358d4a265..961afbad3 100644
--- a/tests/test_distribution.py
+++ b/tests/test_distribution.py
@@ -203,7 +203,6 @@ def test_Uniform_logpdf(low, high, toeval, expected):
     UD = cuqi.distribution.Uniform(low, high)
     assert np.allclose(UD.logpdf(toeval), expected) 
 
-
 @pytest.mark.parametrize("low,high,expected",[(np.array([1, .5]), 
                                                np.array([2, 2.5]),
                                                np.array([[1.5507979 , 1.29090474, 1.89294695],
@@ -217,6 +216,19 @@ def test_Uniform_sample(low, high, expected):
     print(cuqi_samples)
     assert np.allclose(cuqi_samples.samples, expected) 
 
+def test_UnboundedUniform_logpdf_and_gradient():
+    x1 = cuqi.distribution.UnboundedUniform(geometry=1)
+    assert np.allclose(x1.logpdf(np.array([1])), 1.0)
+    assert np.allclose(x1.gradient(np.array([0])), np.array([0]))
+    x2 = cuqi.distribution.UnboundedUniform(geometry=cuqi.geometry.Discrete(2))
+    assert np.allclose(x2.logpdf(np.array([1, 2])), 1.0)
+    assert np.allclose(x2.gradient(np.array([0, 0])), np.array([0, 0]))
+
+def test_UnboundedUniform_sample():
+    x = cuqi.distribution.UnboundedUniform(geometry=2)
+    with pytest.raises(NotImplementedError, match="Cannot sample from UnboundedUniform distribution"):
+        x.sample()
+
 @pytest.mark.parametrize("distribution, kwargs",
                          [(cuqi.distribution.Uniform, 
                           {'low':np.array([2, 2.5, 3,5]),
diff --git a/tests/test_distributions_shape.py b/tests/test_distributions_shape.py
index 7d0f69d1c..7a5b13780 100644
--- a/tests/test_distributions_shape.py
+++ b/tests/test_distributions_shape.py
@@ -13,6 +13,7 @@
     "Posterior",
     "Distribution",
     "JointGaussianSqrtPrec",
+    "UnboundedUniform",
 ]
 
 # Define cases to skip (these are TODO)