How to sample a symmetric matrix with parameters at particular locations?

Hi everyone,
I want to generate a symmetric matrix called K. At some given locations of this matrix, the elements are parameters kij drawn from normal distributions; at other locations, the elements are 0. The uploaded file “adjacent_half_binary.csv” is an upper triangular matrix. It has the same size as matrix K. This upper triangular matrix tells us the locations where the symmetric matrix K contains parameters. Specifically, the “adjacent_half_binary.csv” contains 0 and 1. The locations with 1 are the locations that contain parameters and the locations with 0 contain 0. Below is part of my code used to generate such kind of symmetric matrix with the use of the given upper triangular matrix. However, when I used the obtained K to do further calculations, I obtain the following error:

TypeError: Tensor type field must be a TensorType; found <class ‘theano.graph.type.Generic’>.

Therefore, I think my way to generate the symmetric matrix K is not correct. Could you please provide some advice? Thank you so much!

adjacent_half_binary.csv (14.3 KB)

with pm.Model() as model:

    BoundedNormal_k = pm.Bound(pm.Normal, lower = 0)
    n_ks = np.sum(np.sum(connectome_hcp_half_binary))
    ks = []
    for i in range(n_ks):
        ks.append(BoundedNormal_k('k'+str(i+1),0,0.5))
    k_loc = np.array(adjacent_half_binary).reshape(-1)
    K_half_vec = []
    j = 0
    for i in range(len(k_loc)):
        if k_loc[i] == 1:
            K_half_vec.append(ks[j]) 
            j = j + 1
        else:
            K_half_vec.append(0)

    K_half = np.array(K_half_vec).reshape((84,84))
    K = K_half.T + K_half- np.diag(np.diag(K_half))
    K = theano.shared(K)

This is the simplest way I can think of to generate a symmetric matrix with desired sparsity and independent, normally-distributed entries. The basic idea is that we take the upper triangular matrix, add it to its transpose to create a binary mask, and use this mask to zero out elements of a symmetric matrix. We make sure it is symmetric by summing a matrix with its transpose. Then, since the sum of independent normals is also normal, we are guaranteed to have normally-distributed elements in the entries of K. If you want to tweak the standard deviation of the prior for those entires, you’ll just have to remember that the variance of the sum is the sum of the variances. One downside of this approach is that it samples extra variables in K_raw which are unused since we immediately zero them out in the operation mask * K_raw.

import pymc3 as pm
import numpy as np

d = 4

# This is the upper triangular matrix like in the CSV file
tri  = np.triu(np.ones([d,d]), k=1)

# We convert this into a mask for the symmetric matrix
mask = tri + tri.T + np.eye(d)
zero_locations = [(0,1), (1,2), (2,3)]

for row, col in zero_locations:
    mask[row, col] = 0
    mask[col, row] = 0

with pm.Model():
    K_raw = pm.HalfNormal('K_raw', shape=[d,d])
    K     = pm.Deterministic('K', mask * (K_raw + K_raw.T))
    trace = pm.sample()
1 Like