Merge pull request #65 from DelinteNicolas/main

Update

.readthedocs.yalm → docs/.readthedocs.yaml

@@ -1,6 +1,7 @@
 version: 2
 
 build:
+  os: ubuntu-22.04
   tools:
     python: "3.11"
 

docs/conf.py

@@ -12,17 +12,17 @@
 #
 import os
 import sys
-sys.path.insert(0, os.path.abspath('../src/'))
 
+sys.path.insert(0, os.path.abspath('../src/'))
 
 # -- Project information -----------------------------------------------------
 
 project = 'UNRAVEL'
-copyright = '2022, Nicolas Delinte'
+copyright = '2024, Nicolas Delinte'
 author = 'Nicolas Delinte'
 
 # The full version, including alpha/beta/rc tags
-release = '2023'
+release = '2024'
 
 
 # -- General configuration ---------------------------------------------------

docs/index.rst

@@ -4,7 +4,7 @@
    contain the root `toctree` directive.
 
 Welcome to UNRAVEL's documentation!
-================================
+===================================
 
 Welcome to the documentation of UtiliziNg tRActography to uncoVEr muLti-fixel microstructure (UNRAVEL)!
 
@@ -14,6 +14,7 @@ To unravel has two meanings :
 
 - to disentangle the fibers of
 - to resolve the intricacy, complexity, or obscurity of
+
 With the UNRAVEL framework, we utilize tractography to unravel the microstructure of multi-fixel models.
 
 This repository contains the documentation of the code used to combine macroscopic tractography information with microscopic multi-fixel model estimates in order to improve the accuracy in the estimation of the microstructural properties of neural fibers in a specified tract.

docs/requirements.txt

@@ -1,5 +1,6 @@
 numpydoc
-dipy==1.6.0
+dipy
 scikit-learn
 scikit-image
 pyvista
+sphinx-rtd-theme

src/unravel/__init__.py

@@ -5,5 +5,5 @@
 microstructural properties of neural fibers in a specified tract.
 """
 
-__version__ = "1.4.13"
+__version__ = "1.4.20"
 __author__ = 'Nicolas Delinte'

src/unravel/analysis.py

@@ -10,18 +10,21 @@
 from unravel.core import get_microstructure_map, get_weighted_mean
 
 
-def get_metric_along_trajectory(fixel_weights, metric_maps: list, roi_sections):
+def get_metric_along_trajectory(fixel_weights, metric_maps, roi_sections,
+                                weighting: str = 'tsl'):
     '''
 
 
     Parameters
     ----------
     fixel_weights : 4-D array of shape (x,y,z,K)
         Array containing the relative weights of the K fixels in each voxel.
-    metric_maps : list
-        List of K 3-D arrays of shape (x,y,z) containing metric estimations.
+    metric_maps : 4D array of shape (x,y,z,3,k)
+        List of K 4D arrays of shape (x,y,z) containing metric estimations.
     roi_sections : 3D array of size (x,y,z)
         Labeled array containing the volumes of the section of the tract.
+    weighting : str, optional
+        Weighting used for the mean. The default is 'tsl'.
 
     Returns
     -------
@@ -32,24 +35,24 @@ def get_metric_along_trajectory(fixel_weights, metric_maps: list, roi_sections):
 
     '''
 
-    m_array = np.zeros(np.max(roi_sections))
-    std_array = np.zeros(np.max(roi_sections))
+    m_array = np.zeros(np.max(roi_sections)+1)
+    std_array = np.zeros(np.max(roi_sections)+1)
 
-    if len(metric_maps) == 1:
+    if len(metric_maps.shape) <= 3:
         fixel_weights = fixel_weights[..., np.newaxis]
 
     micro_map = get_microstructure_map(fixel_weights, metric_maps)
 
-    for i in range(np.max(roi_sections)):
+    for i in range(np.max(roi_sections)+1):
 
         if i == 0:
             continue
 
         roi = np.where(roi_sections == i, 1, 0)
-        fixel_weights_roi = fixel_weights * np.repeat(roi[..., np.newaxis],
-                                                      1, axis=-1)
+        fixel_weights_roi = fixel_weights * roi[..., np.newaxis]
 
-        mean, std = get_weighted_mean(micro_map, fixel_weights_roi)
+        mean, std = get_weighted_mean(micro_map, fixel_weights_roi,
+                                      weighting=weighting)
         if mean == 0:
             mean = m_array[i-1]
             std = std_array[i-1]

src/unravel/example.py

@@ -4,7 +4,7 @@
 
 @author: Delinte Nicolas
 
-An example pyhton script using the main methods and outputs of UNRAVEL. A
+An example python script using the main methods and outputs of UNRAVEL. A
 tractogram of the middle anterior section of the corpus callosum is used as
 tractography input.
 

src/unravel/stream.py

@@ -61,10 +61,21 @@ def extract_nodes(trk_file: str, level: int = 3, smooth: bool = True):
     m_start = np.mean(start[selec_streamlines], axis=0)
     m_end = np.mean(end[selec_streamlines], axis=0)
 
-    # Re-orders start and end based on main axial direction
-    # !!! does not work on uni-hemisperal left-right tracts
-    main_dir = np.argmax(np.abs(m_start-m_end))
-    if m_start[main_dir] > m_end[main_dir]:
+    # Re-orders start and end based on main axial direction,
+    # first main direction or three-way vote
+    # !!! does not always work
+    diff_abs = np.abs(m_start-m_end)
+    main_dir = np.argmax(diff_abs)
+    small_dir = np.argmin(diff_abs)
+    vote = np.sum(np.where(m_start-m_end > 0, 1, -1))
+    if diff_abs[main_dir] > (np.sum(diff_abs)-diff_abs[main_dir]):
+        if m_start[main_dir] > m_end[main_dir]:
+            m_start, m_end = m_end, m_start
+    elif diff_abs[small_dir] < (np.sum(diff_abs)-diff_abs[small_dir])/4:
+        idx = 1 if small_dir == 0 else 0
+        if m_start[idx] > m_end[idx]:
+            m_start, m_end = m_end, m_start
+    elif vote > 0:
         m_start, m_end = m_end, m_start
 
     # Iterating over specified level ---------------------------------

src/unravel/utils.py

@@ -6,7 +6,7 @@
 """
 
 import numpy as np
-from dipy.io.streamline import load_tractogram
+from dipy.io.streamline import load_tractogram, save_tractogram
 
 
 def tract_to_ROI(trk_file: str):
@@ -76,10 +76,14 @@ def peaks_to_RGB(peaks, frac=None, fvf=None, order: str = 'rgb'):
 
     if frac is None:
         frac = np.ones(peaks.shape[:3]+(K,))/K
+    elif len(frac.shape) <= 3:
+        frac = frac[..., np.newaxis]
     frac = np.stack((frac,)*3, axis=3)
 
     if fvf is None:
         fvf = np.ones(peaks.shape[:3]+(K,))
+    elif len(fvf.shape) <= 3:
+        fvf = fvf[..., np.newaxis]
     fvf = np.stack((fvf,)*3, axis=3)
 
     rgb = np.sum(abs(peaks)*frac*fvf, axis=-1)
@@ -209,32 +213,6 @@ def tensor_to_DTI(t):
     return FA, AD, RD, MD
 
 
-def xyz_to_spherical(xyz):
-    '''
-    X,y,z coordinates to spherical coordinates.
-
-    Parameters
-    ----------
-    xyz : array of size (n,3)
-        X,y,z coordinates of n points
-
-    Returns
-    -------
-    r : array of size (n)
-        DESCRIPTION.
-    theta : array of size (n)
-        DESCRIPTION.
-    phi : array of size (n)
-        DESCRIPTION.
-
-    '''
-    xy = xyz[:, 0]**2 + xyz[:, 1]**2
-    r = np.sqrt(xy + xyz[:, 2]**2)
-    theta = np.arctan2(np.sqrt(xy), xyz[:, 2])
-    phi = np.arctan2(xyz[:, 1], xyz[:, 0])
-    return r, theta, phi
-
-
 def get_streamline_count(trk) -> int:
     '''
     Returns the number of streamlines in a tractogram.
@@ -504,20 +482,20 @@ def plot_streamline_trajectory(trk, resolution_increase: int = 1,
 
 def xyz_to_spherical(xyz):
     '''
-
+    X,y,z coordinates to spherical coordinates.
 
     Parameters
     ----------
     xyz : array of size (n,3)
-        DESCRIPTION.
+        X,y,z coordinates of n points
 
     Returns
     -------
-    r : TYPE
+    r : array of size (n)
         DESCRIPTION.
-    theta : TYPE
+    theta : array of size (n)
         DESCRIPTION.
-    phi : TYPE
+    phi : array of size (n)
         DESCRIPTION.
 
     '''
@@ -535,3 +513,39 @@ def spherical_to_xyz(theta, phi):
     z = np.cos(theta)
 
     return x, y, z
+
+
+def fuse_trk(trk_file_1: str, trk_file_2: str, output_file: str):
+    '''
+    Creates a new .trk file with all streamlines contained in the two input .trk
+    files. The input files must be in the same space.
+
+    Parameters
+    ----------
+    trk_file_1 : str
+        Filename of first input .trk file.
+    trk_file_2 : str
+        Filename of second input .trk file.
+    output_file : str
+        Filename of output .trk file.
+
+    Returns
+    -------
+    None.
+
+    '''
+
+    trk = load_tractogram(trk_file_1, 'same')
+    trk.to_vox()
+    trk.to_corner()
+    streams_1 = trk.streamlines
+
+    trk = load_tractogram(trk_file_2, 'same')
+    trk.to_vox()
+    trk.to_corner()
+    streams_2 = trk.streamlines
+
+    streams_1.extend(streams_2)
+
+    trk_new = trk.from_sft(streams_1, trk)
+    save_tractogram(trk_new, output_file)