Trying to understand the code

[smauermann]

Dear Bertrand, I am trying to understand the parcellation code published here, especially the spectral clustering part. I plan to implement something closely related but using structural MRI instead. This might be a long shot, but I would be super happy if you could comment on a couple of lines of code, e.g. from group_parcellation.py:

if method == 'spectral':
        i, j = connectivity.nonzero()
        sigma = np.sum((Xv[i] - Xv[j]) ** 2, 1).mean()
        connectivity.data = np.exp( - np.sum((Xv[i] - Xv[j]) ** 2, 1) /
                                      (2 * sigma))
        connectivity = connectivity.copy() + dia_matrix(
            (1.e-3 * np.ones(n_voxels), [0]), 
            shape=(n_voxels, n_voxels)).tocsr()

Especially the last operation bugs me as I cannot seem to understand why you're adding 1.e-3 to the diagonal elements of connectivity.

Many thanks in advance!

Best, Stephan

[bthirion]

----- Mail original -----

De: "neurotronix" notifications@github.com À: "bthirion/frontiers_2014" frontiers_2014@noreply.github.com Envoyé: Mardi 9 Août 2016 20:46:36 Objet: [bthirion/frontiers_2014] Trying to understand the code (#1)

Dear Bertrand, I am trying to understand the parcellation code published here, especially the spectral clustering part. I plan to implement something closely related but using structural MRI instead. This might be a long shot, but I would be super happy if you could comment on a couple of lines of code, e.g. from group_parcellation.py : if method == 'spectral':

define the support of connectivity

i, j = connectivity.nonzero()

define the scale parameter as mean distance between neighbors

sigma = np.sum((Xv[i] - Xv[j]) \ 2, 1).mean()

RBF kernel-based similarity on the features

connectivity.data = np.exp( - np.sum((Xv[i] - Xv[j]) * 2, 1) / (2 \ sigma))

add a diagonal matrix to prevent numerical underflows (linear algebra solvers). Note that this is unrelated to the model and should not have impact on the results (it ossets the eigenvalues, but leaves eigenvectors unchanged).

connectivity = connectivity.copy() + dia_matrix( (1.e-3 * np.ones(n_voxels), [0]), shape=(n_voxels, n_voxels)).tocsr()

Especially the last operation bugs me as I cannot seem to understand why you're adding 1.e-3 to the diagonal elements of connectivity. Please let me know if this is still unclear. Best,

Bertrand

[smauermann]

Many thanks for the quick reply! I was mainly confused by the step which adds the diagonal matrix. I'll have a closer look in the office tomorrow and in case there are still questions I get back to you.

Have a good evening!