# Passing normal distributions as mean for

**URL:** <https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539>\
**Category:** version agnostic\
**Created:** [December 27, 2023, 10:09pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539 "2023-12-27T22:09:27Z")\
**Posts on this page:** 12\
**Page:** 1

<div class="post-metadata">

**Author:** ![zepater](https://avatars.discourse-cdn.com/v4/letter/z/d07c76/32.png) [@zepater](https://discourse.pymc.io/u/zepater)\
**Post date:** [December 27, 2023, 10:09pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/1 "2023-12-27T22:09:27Z")

</div>

Is it possible to use a vector of normal distributions as mean for a multivariate normal?

```auto
#Defined data
x_data
y_data
z_data
with pm.Model() as basic_model:
    x = pm.Normal('xval', mu=x_data.mean(), sigma=x_data.std())
    y = pm.Normal('yval', mu=y_data.mean(), sigma=y_data.std())
    z = pm.Normal('zval', mu=z_data.mean(), sigma=z_data.std())
    
    likelihood = pm.MvNormal(
        "observations", mu=np.array([x, y, z]), cov=np.eye(3,3) * 0.1, 
        observed=np.array([x_data, y_data, z_data])

    trace.sample()

```

---

<div class="post-metadata">

**Author:** ![jessegrabowski](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/jessegrabowski/32/5010_2.png) [@jessegrabowski](https://discourse.pymc.io/u/jessegrabowski)\
**Post date:** [December 28, 2023, 1:38am UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/2 "2023-12-28T01:38:31Z")

</div>

Hi Zepater,

You can do exactly this, but you need to use pytensor functions instead of numpy. You can see some examples and explanation here.

In your case, you could also simplify things greatly by using vectorized distributions. Basically this:

```auto
    x = pm.Normal('xval', mu=0, sigma=1)
    y = pm.Normal('yval', mu=1, sigma=2)
    z = pm.Normal('zval', mu=2, sigma=3)
    d = pt.stack([x, y, z])

```

Is the same as this:

```auto
d = pm.Normal('d', mu=[0,1,2], sigma=[1,2,3])

```

So in your case, you could write:

```auto
# Matrix of data:
X = np.c_[x_data, y_data, z_data]

with pm.Model() as basic_model:
    mu = pm.Normal('mu', mu=X.mean(axis=0), sigma=X.std(axis=0))
    y_hat = pm.MvNormal('y_hat', mu=mu, cov=pt.eye(3) * 0.1, observed=X)
    idata = pm.sample()

```

---

<div class="post-metadata">

**Author:** ![zepater](https://avatars.discourse-cdn.com/v4/letter/z/d07c76/32.png) [@zepater](https://discourse.pymc.io/u/zepater)\
**Post date:** [December 28, 2023, 6:10am UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/3 "2023-12-28T06:10:36Z")

</div>

Thank you for the help!

Before I pass x, y, z to the data, I perform some operations on it.

```auto
x = pm.Normal('xval', mu=0, sigma=1)
y = pm.Normal('yval', mu=1, sigma=2)
z = pm.Normal('zval', mu=2, sigma=3)

#Perform some rotations of the x,y,z here

d = pt.stack([x, y, z])
# Matrix of data:
X = np.c_[x_data, y_data, z_data]

with pm.Model() as basic_model:
    y_hat = pm.MvNormal('y_hat', mu=d, cov=pt.eye(3) * 0.1, observed=X)
    idata = pm.sample()

```

An error I got during the process

```auto
ValueError: Incompatible Elemwise input shapes [(length_of_xyz_data, 3), (3,3)]

```

---

<div class="post-metadata">

**Author:** ![jessegrabowski](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/jessegrabowski/32/5010_2.png) [@jessegrabowski](https://discourse.pymc.io/u/jessegrabowski)\
**Post date:** [December 28, 2023, 12:51pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/4 "2023-12-28T12:51:09Z")

</div>

This snippet runs for me without problems:

```auto
X = np.random.normal(size=(100, 3))

with pm.Model() as basic_model:
    x = pm.Normal('xval', mu=0, sigma=1)
    y = pm.Normal('yval', mu=1, sigma=2)
    z = pm.Normal('zval', mu=2, sigma=3)

    #Perform some rotations of the x,y,z here

    mu = pm.Deterministic('mu', pt.stack([x, y, z]))
    cov = pt.eye(3) * 0.1

    y_hat = pm.MvNormal('y_hat', mu=mu, cov=cov, observed=X)
    idata = pm.sample()

```

If that runs for you, we could take a look at the transformations of x, y, z you are doing?

---

<div class="post-metadata">

**Author:** ![zepater](https://avatars.discourse-cdn.com/v4/letter/z/d07c76/32.png) [@zepater](https://discourse.pymc.io/u/zepater)\
**Post date:** [December 28, 2023, 6:57pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/5 "2023-12-28T18:57:35Z")

</div>

Yes, that ran without problems!

Here’s the transformations I am doing where I define three angles.

```auto
with pm.Model() as basic_model: 
    x = pm.Normal('xval', mu=0, sigma=1)
    y = pm.Normal('yval', mu=1, sigma=2)
    z = pm.Normal('zval', mu=2, sigma=3)
    angle1 = pm.Norma('angle1', mu=0, sigma=np.pi/4)
    angle2 = pm.Norma('angle2', mu=0, sigma=np.pi/4)
    angle3 = pm.Norma('angle2', mu=0, sigma=np.pi/4)

    #Perform some rotations of the x,y,z here

    phi = np.linspace(0,2 * np.pi, 100)
    theta = np.linspace(0, np.pi, 50)
    phi, theta = np.meshgrid(phi, theta)

    xc = x*np.sin(theta) * np.cos(phi)
    yc = y*np.sin(theta) * np.sin(phi)
    zc = z*np.cos(theta)

    R1 = np.array([np.cos(angle1), -np.sin(angle1), 0], [np.sin(angle1), np.cos(angle1), 0], [0,0,1]])
    R2 = np.array([1,0,0], [0, np.cos(angle2), -np.sin(angle2)], [0, np.sin(angle2), np.cos(angle2)]])
    R3 = np.array([[np.cos(angle3), -np.sin(angle3), 0], [np.sin(angle3), np.cos(angle3), 0], [0,0,1]])

   rotated_coords = np.dot(R3, np.dot(R2, np.dot(R1, np.array([xc.flatten(), yc.flatten(), zc.flatten()]))))

   rotated_coords = rotated_coords.reshape(3, -1)

    mu = pm.Deterministic('mu', pt.stack([rotated_coords[0][0], rotated_coords[1][0], rotated_coords[2][0]]))
    cov = pt.eye(3) * 0.1

    y_hat = pm.MvNormal('y_hat', mu=mu, cov=cov, observed=X)
    idata = pm.sample()

```

---

<div class="post-metadata">

**Author:** ![jessegrabowski](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/jessegrabowski/32/5010_2.png) [@jessegrabowski](https://discourse.pymc.io/u/jessegrabowski)\
**Post date:** [December 28, 2023, 9:01pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/6 "2023-12-28T21:01:42Z")

</div>

You want to get into the habit of using the `pytensor` functions instead of numpy inside pymc models. Sometimes things work as expected, but it’s not technically supported. In particular, all of those `np.array` calls with pymc variables will not work as expected.

I guess the only thing missing is `meshgrid`, which just looks like this:

```auto
def pt_meshgrid(*xi):
    ndim = len(xi)
    s0 = (1,) * ndim
    outputs = [x.reshape(s0[:i] + (-1,) + s0[i + 1:]) for i, x in enumerate(xi)]
    return pt.broadcast_arrays(*outputs)

```

So you should rewrite you code as:

```auto
with pm.Model() as basic_model: 
    x = pm.Normal('xval', mu=0, sigma=1)
    y = pm.Normal('yval', mu=1, sigma=2)
    z = pm.Normal('zval', mu=2, sigma=3)
    angle1 = pm.Normal('angle1', mu=0, sigma=np.pi/4)
    angle2 = pm.Normal('angle2', mu=0, sigma=np.pi/4)
    angle3 = pm.Normal('angle3', mu=0, sigma=np.pi/4)

    #Perform some rotations of the x,y,z here

    phi = pt.linspace(0,2 * np.pi, 100)
    theta = pt.linspace(0, np.pi, 50)
    phi, theta = pt_meshgrid(phi, theta)

    xc = x*pt.sin(theta) * pt.cos(phi)
    yc = y*pt.sin(theta) * pt.sin(phi)
    zc = z*pt.cos(theta)

    R1 = pt.stack([[pt.cos(angle1), -pt.sin(angle1), 0], [pt.sin(angle1), pt.cos(angle1), 0], [0,0,1]])
    R2 = pt.stack([[1,0,0], [0, pt.cos(angle2), -pt.sin(angle2)], [0, pt.sin(angle2), pt.cos(angle2)]])
    R3 = pt.stack([[pt.cos(angle3), -pt.sin(angle3), 0], [pt.sin(angle3), pt.cos(angle3), 0], [0,0,1]])

    rotated_coords = pt.dot(R3, pt.dot(R2, pt.dot(R1, pt.stack([xc.flatten(), yc.flatten(), zc.flatten()]))))

    rotated_coords = rotated_coords.reshape((3, -1))

    mu = pm.Deterministic('mu', pt.stack([rotated_coords[0][0], rotated_coords[1][0], rotated_coords[2][0]]))
    cov = pt.eye(3) * 0.1

    y_hat = pm.MvNormal('y_hat', mu=mu, cov=cov, observed=X)
    idata = pm.sample()

```

That meshgrid function would make a good first PR, you’d just have to write a couple tests for it 🙂

---

<div class="post-metadata">

**Author:** ![zepater](https://avatars.discourse-cdn.com/v4/letter/z/d07c76/32.png) [@zepater](https://discourse.pymc.io/u/zepater)\
**Post date:** [January 3, 2024, 5:20am UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/7 "2024-01-03T05:20:51Z")

</div>

Okay, I will look into that! Thank you! When I evaluated rotated\_coords, I obtain the following error.

```auto
Value Error: Cannot reshape array of size 100 into shape ()

```

---

<div class="post-metadata">

**Author:** ![jessegrabowski](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/jessegrabowski/32/5010_2.png) [@jessegrabowski](https://discourse.pymc.io/u/jessegrabowski)\
**Post date:** [January 3, 2024, 12:50pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/8 "2024-01-03T12:50:01Z")

</div>

The code I posted runs for me without problems. What versions of pymc/pytensor are you running?

---

<div class="post-metadata">

**Author:** ![zepater](https://avatars.discourse-cdn.com/v4/letter/z/d07c76/32.png) [@zepater](https://discourse.pymc.io/u/zepater)\
**Post date:** [January 3, 2024, 8:04pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/9 "2024-01-03T20:04:18Z")

</div>

Version 5.10.3. I think the issue that I found was that the size of the observed data (as opposed to X in this minimum working example) is (2600,3). When I only use a smaller subset, say (50,3), it works. Is there a size limit?

---

<div class="post-metadata">

**Author:** ![jessegrabowski](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/jessegrabowski/32/5010_2.png) [@jessegrabowski](https://discourse.pymc.io/u/jessegrabowski)\
**Post date:** [January 3, 2024, 9:17pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/10 "2024-01-03T21:17:12Z")

</div>

No, there shouldn’t be

---

<div class="post-metadata">

**Author:** ![zepater](https://avatars.discourse-cdn.com/v4/letter/z/d07c76/32.png) [@zepater](https://discourse.pymc.io/u/zepater)\
**Post date:** [January 3, 2024, 11:19pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/11 "2024-01-03T23:19:49Z")

</div>

What I noticed was that the size of phi, theta in the linspace definitions seem to have an effect of the number of observed values allowed. Could you increase the size of the random X to 200 and see if that causes any issues on your end?

---

<div class="post-metadata">

**Author:** ![jessegrabowski](https://yyz2.discourse-cdn.com/flex036/user_avatar/discourse.pymc.io/jessegrabowski/32/5010_2.png) [@jessegrabowski](https://discourse.pymc.io/u/jessegrabowski)\
**Post date:** [January 3, 2024, 11:26pm UTC](https://discourse.pymc.io/t/passing-normal-distributions-as-mean-for/13539/12 "2024-01-03T23:26:53Z")

</div>

In the code I posted all the computation is independent of X (it’s just used for observations). Could you share exactly the code you are running?