MRI与脑电图

显示一组带有MRI图像的子图，其强度直方图和一些EEG轨迹。

MRI与脑电图示例

  1. import numpy as np
  2. import matplotlib.pyplot as plt
  3. import matplotlib.cbook as cbook
  4. import matplotlib.cm as cm
  6. from matplotlib.collections import LineCollection
  7. from matplotlib.ticker import MultipleLocator
  9. fig = plt.figure("MRI_with_EEG")
  11. # Load the MRI data (256x256 16 bit integers)
  12. with cbook.get_sample_data('s1045.ima.gz') as dfile:
  13.     im = np.fromstring(dfile.read(), np.uint16).reshape((256, 256))
  15. # Plot the MRI image
  16. ax0 = fig.add_subplot(2, 2, 1)
  17. ax0.imshow(im, cmap=cm.gray)
  18. ax0.axis('off')
  20. # Plot the histogram of MRI intensity
  21. ax1 = fig.add_subplot(2, 2, 2)
  22. im = np.ravel(im)
  23. im = im[np.nonzero(im)]  # Ignore the background
  24. im = im / (2**16 - 1)  # Normalize
  25. ax1.hist(im, bins=100)
  26. ax1.xaxis.set_major_locator(MultipleLocator(0.4))
  27. ax1.minorticks_on()
  28. ax1.set_yticks([])
  29. ax1.set_xlabel('Intensity (a.u.)')
  30. ax1.set_ylabel('MRI density')
  32. # Load the EEG data
  33. n_samples, n_rows = 800, 4
  34. with cbook.get_sample_data('eeg.dat') as eegfile:
  35.     data = np.fromfile(eegfile, dtype=float).reshape((n_samples, n_rows))
  36. t = 10 * np.arange(n_samples) / n_samples
  38. # Plot the EEG
  39. ticklocs = []
  40. ax2 = fig.add_subplot(2, 1, 2)
  41. ax2.set_xlim(0, 10)
  42. ax2.set_xticks(np.arange(10))
  43. dmin = data.min()
  44. dmax = data.max()
  45. dr = (dmax - dmin) * 0.7  # Crowd them a bit.
  46. y0 = dmin
  47. y1 = (n_rows - 1) * dr + dmax
  48. ax2.set_ylim(y0, y1)
  50. segs = []
  51. for i in range(n_rows):
  52.     segs.append(np.column_stack((t, data[:, i])))
  53.     ticklocs.append(i * dr)
  55. offsets = np.zeros((n_rows, 2), dtype=float)
  56. offsets[:, 1] = ticklocs
  58. lines = LineCollection(segs, offsets=offsets, transOffset=None)
  59. ax2.add_collection(lines)
  61. # Set the yticks to use axes coordinates on the y axis
  62. ax2.set_yticks(ticklocs)
  63. ax2.set_yticklabels(['PG3', 'PG5', 'PG7', 'PG9'])
  65. ax2.set_xlabel('Time (s)')
  68. plt.tight_layout()
  69. plt.show()

下载这个示例