File: multi_subject_parcellation.py

package info (click to toggle)
nipy 0.6.1-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 7,352 kB
  • sloc: python: 39,115; ansic: 30,931; makefile: 210; sh: 93
file content (63 lines) | stat: -rwxr-xr-x 1,926 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
#!/usr/bin/env python3
# emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*-
# vi: set ft=python sts=4 ts=4 sw=4 et:
__doc__ = """
This script contains a quick demo on  a multi-subject parcellation on a toy 2D
example.

Note how the middle parcels adapt to the individual configuration.

Needs matplotlib
"""
print(__doc__)

import numpy as np

try:
    import matplotlib.pyplot as plt
except ImportError:
    raise RuntimeError("This script needs the matplotlib library")

import nipy.labs.spatial_models.discrete_domain as dom
import nipy.labs.spatial_models.hierarchical_parcellation as hp
import nipy.labs.utils.simul_multisubject_fmri_dataset as simul

# step 1:  generate some synthetic data
n_subj = 10
shape = (60, 60)
pos = 3 * np.array([[6, 7],
                  [10, 10],
                  [15, 10]])
ampli = np.array([5, 7, 6])
sjitter = 6.0
dataset = simul.surrogate_2d_dataset(n_subj=n_subj, shape=shape, pos=pos,
                                     ampli=ampli, width=10.0)
# dataset represents 2D activation images from n_subj subjects,

# step 2 : prepare all the information for the parcellation
nbparcel = 10
ldata = np.reshape(dataset, (n_subj, np.prod(shape), 1))
domain = dom.grid_domain_from_shape(shape)

# step 3 : run the algorithm
Pa = hp.hparcel(domain, ldata, nbparcel, mu=3.0)
# note: play with mu to change the 'stiffness of the parcellation'

# step 4:  look at the results
Label = np.array([np.reshape(Pa.individual_labels[:, s], shape)
                  for s in range(n_subj)])

plt.figure(figsize=(8, 4))
plt.title('Input data')
for s in range(n_subj):
    plt.subplot(2, 5, s + 1)
    plt.imshow(dataset[s], interpolation='nearest')
    plt.axis('off')

plt.figure(figsize=(8, 4))
plt.title('Resulting parcels')
for s in range(n_subj):
    plt.subplot(2, 5, s+1)
    plt.imshow(Label[s], interpolation='nearest', vmin=-1, vmax=nbparcel)
    plt.axis('off')
plt.show()