Skip to content

Reading freesurfer stats files correctly #1152

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 7 commits into
base: master
Choose a base branch
from
Open
Changes from 3 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
54 changes: 53 additions & 1 deletion nibabel/freesurfer/io.py
Original file line number Diff line number Diff line change
@@ -1,10 +1,12 @@
""" Read / write FreeSurfer geometry, morphometry, label, annotation formats
""" Read / write FreeSurfer geometry, morphometry, label, stats, annotation formats
"""
from __future__ import annotations

import warnings
import numpy as np
import getpass
import time
import re

from collections import OrderedDict
from ..openers import Opener
Expand Down Expand Up @@ -622,3 +624,53 @@ def _serialize_volume_info(volume_info):
strings.append(
f'{key:6s} = {val[0]:.10g} {val[1]:.10g} {val[2]:.10g}\n'.encode('utf-8'))
return b''.join(strings)

def read_stats_file(file_path):
"""Extracts stats from stats files except '*curv.stats' files

Parameters
----------
file_path: str, required

Examples
--------
>>> stats_a2009, column_names = read_stats_file(r'lh.aparc.a2009s.stats')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If we're going to have example code, it needs to run. We could drop an example file in nibabel/tests/data/ and access it with

from nibabel.testing import test_data
stats_a2009 = read_stats_file(test_data(fname='lh.aparc.a2009s.stats'))


"""
with open(file_path, 'r') as f:
for line in f:
if re.findall(r'ColHeaders .*', line):
parameters = line.split()
break
f.close()
stats = np.loadtxt(file_path, comments='#', dtype=str)
column_names = parameters[2:]
return stats, column_names


def read_stats_file_both_hemispheres(file_path: str):
"""Extracts stats data of both hemispheres and merges them

Parameters
----------
file_path: str, required
Path of the stats file belong to left (lh) or right(rh) hemisphere

Returns
-------
stats_both_hemispheres: ndarray
Stats data of both hemisphers
column_naems: ndarray
Name of columns

Examples
--------
>>> stats_a2009, column_names = read_stats_file_both_hemispheres(r'lh.aparc.a2009s.stats')

"""
stats_left, columns_left = read_stats_file(file_path.replace('rh', 'lh'))
stats_right, columns_right = read_stats_file(file_path.replace('lh', 'rh'))
stats_both_hemispheres = np.concatenate((stats_left, stats_right[:, 1:]), axis=1)
column_names = [col_name + '_left' for col_name in columns_left] + [col_name + '_right' for col_name in
columns_right[1:]]
Comment on lines +662 to +666
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You're no longer returning column names. I think you can use numpy.lib.recfunctions.join_by to do what you want.

I would also use the following pattern for doctests:

    .. testsetup::

        >>> import os
        >>> from nibabel.testing import get_test_data
        >>> cwd = os.getcwd()
        >>> os.chdir(get_test_data())

    >>> stats_a2009 = read_stats_file_both_hemispheres(r'lh.aparc.a2009s.stats')

    .. testcleanup::

        >>> os.chdir(cwd)

return stats_both_hemispheres, column_names