|
6 | 6 | import numpy as np
|
7 | 7 | import getpass
|
8 | 8 | import time
|
9 |
| -from os.path import exists as pexists, realpath |
10 | 9 |
|
11 | 10 | from collections import OrderedDict
|
12 | 11 | from six.moves import xrange
|
@@ -48,122 +47,6 @@ def _fread3_many(fobj, n):
|
48 | 47 | return (b1 << 16) + (b2 << 8) + b3
|
49 | 48 |
|
50 | 49 |
|
51 |
| -def read_aseg_stats(seg_stats_file, |
52 |
| - set_wanted = 'subcortical', |
53 |
| - volumes_only = False): |
54 |
| - """ |
55 |
| - Returns the subcortical stats found in Freesurfer output: subid/stats/aseg.stats |
56 |
| -
|
57 |
| - Tries to match the outputs returned by Freesurfer's Matlab counter part: load_segstats.m |
58 |
| - |
59 |
| - Parameters |
60 |
| - ---------- |
61 |
| - seg_stats_file : str |
62 |
| - Abs path to aseg.stats file. |
63 |
| - |
64 |
| - set_wanted : str |
65 |
| - Which set of volumes to return, among ['subcortical', 'wholebrain', 'etiv_only' ]. |
66 |
| - Default: 'subcortical'. |
67 |
| - The choice 'subcortical' returns the usual subortical segmentations. |
68 |
| - The choice 'wholebrain' returns the volumes in aseg.stats coded as : |
69 |
| - [ 'BrainSegVol', 'BrainSegVolNotVent', 'lhCortexVol', 'rhCortexVol', 'lhCorticalWhiteMatterVol', |
70 |
| - 'rhCorticalWhiteMatterVol', 'SubCortGrayVol', 'TotalGrayVol', 'SupraTentorialVol', |
71 |
| - 'SupraTentorialVolNotVent', 'MaskVol', 'BrainSegVol-to-eTIV', 'MaskVol-to-eTIV', |
72 |
| - 'lhSurfaceHoles', 'rhSurfaceHoles', 'eTIV' ] |
73 |
| - These are noted as 'Measure' in the commented section of stats/aseg.stats file. |
74 |
| - The choice 'etiv_only' returns the value for eTIV (estimated total intra-cranial volume) only. |
75 |
| -
|
76 |
| - volumes_only : bool |
77 |
| - Flag to indicate only the volumes are wanted. |
78 |
| -
|
79 |
| - Default: False, returning all info available, to closely match the outputs returned by Freesurfer's Matlab counter part: |
80 |
| - https://github.com/freesurfer/freesurfer/blob/dev/matlab/load_segstats.m |
81 |
| -
|
82 |
| - Returns |
83 |
| - ------- |
84 |
| - By default (volumes_only=False), three arrays are returned: |
85 |
| -
|
86 |
| - seg_name : numpy array of strings |
87 |
| - Array of segmentation names |
88 |
| - seg_index : numpy array |
89 |
| - Array of indices of segmentations into the Freesurfer color lookup table. |
90 |
| - seg_stats : numpy array |
91 |
| - Matrix of subcortical statistics, with the following 5 columns by default. |
92 |
| - If volumes_only = True, only the volumes in mm^3 are returned. |
93 |
| - Columns in the full output are: |
94 |
| - 1. number of voxels |
95 |
| - 2. volume of voxels (mm^3) -- same as number but scaled by voxvol |
96 |
| - 3. mean intensity over space |
97 |
| - 4. std intensity over space |
98 |
| - 5. min intensity over space |
99 |
| - 6. max intensity over space |
100 |
| - 7. range intensity over space |
101 |
| -
|
102 |
| - When volumes_only=True, only one array is returned containing only volumes. |
103 |
| -
|
104 |
| - """ |
105 |
| - |
106 |
| - seg_stats_file = realpath(seg_stats_file) |
107 |
| - if not pexists(seg_stats_file): |
108 |
| - raise IOError('given path does not exist : {}'.format(seg_stats_file)) |
109 |
| - |
110 |
| - acceptable_choices = ['subcortical', 'wholebrain', 'etiv_only'] |
111 |
| - set_wanted = set_wanted.lower() |
112 |
| - if set_wanted not in acceptable_choices: |
113 |
| - raise ValueError('Invalid choice. Choose one among: {}'.format(acceptable_choices)) |
114 |
| - |
115 |
| - if set_wanted in 'subcortical': |
116 |
| - stats = np.loadtxt(seg_stats_file, dtype="i1,i1,i4,f4,S50,f4,f4,f4,f4,f4") |
117 |
| - if volumes_only: |
118 |
| - out_data = np.array([seg[3] for seg in stats]) |
119 |
| - else: |
120 |
| - # need to ensure both two types return data correspond in seg order |
121 |
| - out_data = stats |
122 |
| - |
123 |
| - elif set_wanted in ['wholebrain', 'etiv_only']: |
124 |
| - wb_regex_pattern = r'# Measure ([\w/+_\- ]+), ([\w/+_\- ]+), ([\w/+_\- ]+), ([\d\.]+), ([\w/+_\-^]+)' |
125 |
| - datatypes = np.dtype('U100,U100,U100,f8,U10') |
126 |
| - stats = np.fromregex(seg_stats_file, wb_regex_pattern, dtype=datatypes) |
127 |
| - if set_wanted in ['etiv_only']: |
128 |
| - out_data = np.array([seg[3] for seg in stats if seg[1] == 'eTIV']) |
129 |
| - else: |
130 |
| - out_data = np.array([seg[3] for seg in stats]) |
131 |
| - |
132 |
| - return out_data |
133 |
| - |
134 |
| - |
135 |
| -def read_aparc_stats(file_path): |
136 |
| - """Read statistics on cortical features (such as thickness, curvature etc) produced by Freesurfer. |
137 |
| -
|
138 |
| - file_path would contain whether it is from the right or left hemisphere. |
139 |
| -
|
140 |
| - """ |
141 |
| - |
142 |
| - # ColHeaders StructName NumVert SurfArea GrayVol ThickAvg ThickStd MeanCurv GausCurv FoldInd CurvInd |
143 |
| - aparc_roi_dtype = [('StructName', 'S50'), ('NumVert', '<i4'), ('SurfArea', '<i4'), ('GrayVol', '<i4'), |
144 |
| - ('ThickAvg', '<f4'), ('ThickStd', '<f4'), ('MeanCurv', '<f4'), ('GausCurv', '<f4'), |
145 |
| - ('FoldInd', '<f4'), ('CurvInd', '<f4')] |
146 |
| - roi_stats = np.genfromtxt(file_path, dtype=aparc_roi_dtype, filling_values=np.NaN) |
147 |
| - subset = ['SurfArea', 'GrayVol', 'ThickAvg', 'ThickStd', 'MeanCurv', 'GausCurv', 'FoldInd', 'CurvInd'] |
148 |
| - roi_stats_values = np.full((len(roi_stats), len(subset)), np.NaN) |
149 |
| - for idx, stat in enumerate(roi_stats): |
150 |
| - roi_stats_values[idx,:] = [ stat[feat] for feat in subset ] |
151 |
| - |
152 |
| - # whole cortex |
153 |
| - # Measure Cortex, NumVert, Number of Vertices, 120233, unitless |
154 |
| - # Measure Cortex, WhiteSurfArea, White Surface Total Area, 85633.5, mm^2 |
155 |
| - # Measure Cortex, MeanThickness, Mean Thickness, 2.59632, mm |
156 |
| - wb_regex_pattern = r'# Measure Cortex, ([\w/+_\- ]+), ([\w/+_\- ]+), ([\d\.]+), ([\w/+_\-^]+)' |
157 |
| - wb_aparc_dtype = np.dtype('U100,U100,f8,U10') |
158 |
| - # wb_aparc_dtype = [('f0', '<U100'), ('f1', '<U100'), ('f2', '<f8'), ('f3', '<U10')] |
159 |
| - wb_stats = np.fromregex(file_path, wb_regex_pattern, dtype=wb_aparc_dtype) |
160 |
| - |
161 |
| - # concatenating while surf total area and global mean thickness |
162 |
| - stats = np.hstack((roi_stats_values.flatten(), (wb_stats[1][2], wb_stats[2][2]))) |
163 |
| - |
164 |
| - return stats |
165 |
| - |
166 |
| - |
167 | 50 | def _read_volume_info(fobj):
|
168 | 51 | """Helper for reading the footer from a surface file."""
|
169 | 52 | volume_info = OrderedDict()
|
|
0 commit comments