BlueBrain / diameter-synthesis

Synthesize diameters of neuronal morphologies
Apache License 2.0
3 stars 2 forks source link

How to fix TypeError: expected 1D vector for x #33

Closed riddick-the-furyan closed 1 year ago

riddick-the-furyan commented 1 year ago

The problem seems to be the x value in simple_diamitizer.py. I'm not sure which extra value is required for the vector or maybe a different approach is needed to fix the bug based on the intent of the developers in their design. Any help will be greatly appreciated.

Cell In[1], line 11, in run(output_dir, data_dir)
      9 """Run the example for extracting inputs for synthesis."""
     10 # Generate distribution from directory of neurons
---> 11 distr = extract_input.distributions(
     12     data_dir , "neurons", feature="path_distances", diameter_model="default"
     13 )
     15 # Save distributions in a json file
     16 with open(output_dir / "test_distr.json", "w", encoding="utf-8") as f:

File ~\anaconda3\lib\site-packages\neurots\extract_input\input_distributions.py:102, in distributions(filepath, neurite_types, threshold_sec, diameter_input_morph, feature, diameter_model)
    100     input_distributions["diameter"] = {}
    101     input_distributions["diameter"]["method"] = "default"
--> 102     input_distributions["diameter"]["diameter"] = build_diameter_models(
    103         morphology, config={"models": ["simpler"], "neurite_types": neurite_types}
    104     )
    105 else:
    106     raise NotImplementedError(f"Diameter model {diameter_model} not understood")

File ~\anaconda3\lib\site-packages\diameter_synthesis\build_models.py:82, in build(morphologies, config, with_data)
     71 """Build a model from a set of morphologies and a config file.
     72 
     73 Args:
   (...)
     79     dict: parameter of the models (and possibly the data extracted if `with_data==True`).
     80 """
     81 all_models = _get_model_builder(config)
---> 82 models_params, models_data = all_models[config["models"][0]](morphologies, config)
     83 if with_data:
     84     return models_params, models_data

File ~\anaconda3\lib\site-packages\diameter_synthesis\simpler_diametrizer.py:68, in build_simpler_model(morphologies, config, fit_orders)
     65     all_lengths[neurite_type] = lengths.tolist()
     66     all_diams[neurite_type] += diams.tolist()
---> 68 p, extra = Polynomial.fit(
     69     all_lengths[neurite_type],
     70     all_diams[neurite_type],
     71     fit_orders[neurite_type],
     72     full=True,
     73 )
     74 residues[neurite_type] = float(extra[0][0])
     75 coeffs[neurite_type] = p.convert().coef.tolist()

File ~\anaconda3\lib\site-packages\numpy\polynomial\_polybase.py:980, in ABCPolyBase.fit(cls, x, y, deg, domain, rcond, full, w, window)
    977     window = cls.window
    979 xnew = pu.mapdomain(x, domain, window)
--> 980 res = cls._fit(xnew, y, deg, w=w, rcond=rcond, full=full)
    981 if full:
    982     [coef, status] = res

File ~\anaconda3\lib\site-packages\numpy\polynomial\polynomial.py:1362, in polyfit(x, y, deg, rcond, full, w)
   1214 def polyfit(x, y, deg, rcond=None, full=False, w=None):
   1215     """
   1216     Least-squares fit of a polynomial to data.
   1217 
   (...)
   1360 
   1361     """
-> 1362     return pu._fit(polyvander, x, y, deg, rcond, full, w)

File ~\anaconda3\lib\site-packages\numpy\polynomial\polyutils.py:616, in _fit(vander_f, x, y, deg, rcond, full, w)
    614     raise ValueError("expected deg >= 0")
    615 if x.ndim != 1:
--> 616     raise TypeError("expected 1D vector for x")
    617 if x.size == 0:
    618     raise TypeError("expected non-empty vector for x")

TypeError: expected 1D vector for x
adrien-berchet commented 1 year ago

Could you specify which versions do you use for NeuroTS, TMD, diameter-synthesis and Numpy please?

riddick-the-furyan commented 1 year ago

NeuroTS==3.3.1, TMD==2.2.0, diameter-synthesis==0.5.3 and numpy==1.23.5

Edit: I didn't update them because I was worried that the errors, I fixed in the package would be removed but I updated them, and they weren't. Though I still get the same error mentioned above.

adrien-berchet commented 1 year ago

Did you edit the files? There are weird things in your backtrace... For example, data_dir , "neurons" in the line 12 of

---> 11 distr = extract_input.distributions(
     12     data_dir , "neurons", feature="path_distances", diameter_model="default"
     13 )

should be data_dir / "neurons". And the line 65 of diameter_synthesis\simpler_diametrizer.py should be:

all_lengths[neurite_type] += lengths.tolist()

(+= instead of =)

riddick-the-furyan commented 1 year ago

I got an error when I put "/" and I thought it must have been an error. So I put a comma instead which resolved the error. I forgot what the name of the error but the exception I got indicated that it was that after I had run the code.

I also changed how the path is put into the run function. I'm sorry for being unable to provide the name of the errors but I just put the string fo the path instead of calling the function Path() which resloved that error also. I can't recall from memeory whether the type of error it was, was related to it or not.

In all_lengths[neurite_type] += lengths.tolist() in the source code on my pc it looked like this all_lengths[neurite_type] /= lengths.tolist() . So I removed the "/" as I thought it was a mistake made. I didn't realize it was "+=" .

There were other things I changed as well but don't know in which .py file. The other were problem with the if statements were if neurite_types = None then it would add a list of two strings in the form of a list of types of dendrites. For some reason it was still None which I think was because of there being two if conditions "on top of each other" in a single if statement. If that makes sense.

Maybe the issue as to why the source code of these packages in regards to errors may be different or why I get errors could be because I installed it via pip? Idk but I hope this helps.

adrien-berchet commented 1 year ago

I suggest we try to solve your issue on the current code because I can't help you on your modified version of the code. 1st step: does the command python extract_synthesis_inputs.py works for you when run from the neurots/examples directory? If it does not, what is the error? If it does, what do you want to achieve exactly and how did you try to do it?

Anyway, the best would be to provide a minimal working example of your issue, so I can reproduce it on my side and try to understand it.

riddick-the-furyan commented 1 year ago

Here's a sample of the code from extract_synthesis_inputs.py:

import json
from pathlib import Path

import neurots
from neurots import extract_input

def run(output_dir, data_dir):
    """Run the example for extracting inputs for synthesis."""
    # Generate distribution from directory of neurons
    distr = extract_input.distributions(
        data_dir , "neurons", feature="path_distances", diameter_model="default"
    )

if __name__ == "__main__":
    result_dir = Path("results_extract_synthesis_inputs")
    result_dir.mkdir(parents=True, exist_ok=True)

    run(result_dir, "C:\\Users\\user\Documents\\Experiment_01\\Data\Human\\Neocortex_interneurons")
    print("DONE!")
riddick-the-furyan commented 1 year ago

Sorry, having a bit of brain fog and didn't read your question properly. I can't run the code from directory. I can only run it from jupyter notebooks. I had a seperate error that occured when running it in any other way besides jupyter notebook. I got an UnboundLocal error each time I did and I had to use jupyter notebook to resolve it.

Edit: Also here is simple_diametrizer.py source code:

"""Simpler diametrizer."""
from functools import partial

import matplotlib.pyplot as plt
import numpy as np
from neurom import NeuriteType
from neurom import iter_sections
from neurom.core.morphology import Morphology
from neurom.core.morphology import Section
from numpy.polynomial import Polynomial
from scipy.stats import pearsonr

def section_path_length(section, cache):
    """Path length from section to root."""
    if section.id not in cache:
        cache[section.id] = sum(s.length for s in section.iupstream())

    return cache[section.id]

def build_simpler_model(morphologies, config, fit_orders=None):
    """Build diameter model."""
    neurite_types = config.get("neurite_types")
    if neurite_types is None:
        neurite_types = ["basal_dendrite", "apical_dendrite"]
    if fit_orders is None:
        fit_orders = {"basal_dendrite": 1, "apical_dendrite": 2, "axon": 1}

    neurite_types = ["basal_dendrite", "apical_dendrite"]

    coeffs = {}
    residues = {}
    all_diams = {n_type: [] for n_type in neurite_types}
    all_lengths = {n_type: [] for n_type in neurite_types}

    for neurite_type in neurite_types:
        for m in morphologies:
            diams = []
            lengths = []
            for neurite in m.neurites:
                cache = {}
                if neurite.type == getattr(NeuriteType, neurite_type):
                    for section in iter_sections(neurite):
                        diams.append(2 * np.mean(section.points[:, 3]))
                        tip_length = max(
                            section_path_length(_section, cache) for _section in section.ipreorder()
                        )
                        lengths.append(
                            tip_length - section_path_length(section, cache) + section.length
                        )
            lengths = np.array(lengths)
            diams = np.array(diams)
            lengths = lengths.max()
            all_lengths[neurite_type] += lengths.tolist()
            all_diams[neurite_type] += diams.tolist()

        p, extra = Polynomial.fit(
            all_lengths[neurite_type],
            all_diams[neurite_type],
            fit_orders[neurite_type],
            full=True
        )
        residues[neurite_type] = float(extra[0][0])
        coeffs[neurite_type] = p.convert().coef.tolist()
    return coeffs, [all_lengths, all_diams, residues]
riddick-the-furyan commented 1 year ago

Here also is the source code of build_models.py

import logging
from functools import partial

from neurom import NeuriteType

import diameter_synthesis.morph_functions as morph_funcs
from diameter_synthesis.distribution_fitting import fit_distribution
from diameter_synthesis.exception import DiameterSynthesisError
from diameter_synthesis.simpler_diametrizer import build_simpler_model

L = logging.getLogger(__name__)

def build(morphologies, config, with_data=False):
    """Build a model from a set of morphologies and a config file.

    Args:
        morphologies (dict): dictionary of morphologies.
        config (dict): configuration dictionary.
        with_data (bool): if set to True, the model data are also returned.

    Returns:
        dict: parameter of the models (and possibly the data extracted if `with_data==True`).
    """
    all_models = _get_model_builder(config)
    models_params, models_data = all_models[config["models"][0]](morphologies, config)
    if with_data:
        return models_params, models_data
    return models_params
riddick-the-furyan commented 1 year ago

Lastly here's the code for input_distributions.py:

"""Input distributions."""

# Copyright (C) 2021  Blue Brain Project, EPFL
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <https://www.gnu.org/licenses/>.

import logging

import tmd
from diameter_synthesis.build_models import build as build_diameter_models
from neurom import NeuriteType
from neurom import load_morphologies

from neurots.extract_input import from_diameter
from neurots.extract_input.from_neurom import number_neurites
from neurots.extract_input.from_neurom import soma_data
from neurots.extract_input.from_neurom import trunk_neurite
from neurots.extract_input.from_TMD import persistent_homology_angles
from neurots.utils import format_values
from neurots.utils import neurite_type_warning

L = logging.getLogger(__name__)

def _append_dicts(*args):
    """Merge all dicts into the first one."""
    ret = args[0]
    for other_dict in args[1:]:
        ret.update(other_dict)
    return ret

def distributions(
    filepath,
    neurite_types=None,
    threshold_sec=2,
    diameter_input_morph=None,
    feature="path_distances",
    diameter_model=None,
):
    """Extracts the input distributions from an input population.

    The population is defined by a directory of swc or h5 files.

    Args:
        filepath (str): the morphology file.
        neurite_types (list[str]): the neurite types to consider.
        threshold_sec (int): defines the minimum accepted number of terminations.
        diameter_input_morph (str): if input set of morphologies is provided it will be used for the
            generation of diameter model, if no input is provided no diameter model will be
            generated.
        feature (str): defines the TMD feature that will be used to extract the persistence barcode
            (can be `radial_distances`, `path_distances` or `trunk_length`). It is also possible to
            define one different feature per neurite type using a dict like
            ``{<neurite type 1>: <feature 1>, ...}``.
        diameter_model (str): model for diameters, internal models are `M1`, `M2`, `M3`, `M4` and
            `M5`. Can be set to `external` for external model.

    Returns:
        dict: The input distributions.
    """
    if neurite_types is None:
        neurite_types = ["basal_dendrite", "apical_dendrite", "axon"]

    for i, neurite_type in enumerate(neurite_types):
        if neurite_type in ("basal", "apical"):
            neurite_type_warning(neurite_type)
            neurite_types[i] = neurite_type + "_dendrite"

    pop_tmd = tmd.io.load_population(filepath, use_morphio=True)
    pop_nm = load_morphologies(filepath)

    input_distributions = {"soma": {}, "basal_dendrite": {}, "apical_dendrite": {}, "axon": {}}
    input_distributions["soma"] = soma_data(pop_nm)

    if diameter_input_morph is None:
        diameter_input_morph = filepath
    morphology = load_morphologies(diameter_input_morph)
    if isinstance(diameter_model, str) and diameter_model.startswith("M"):
        input_distributions["diameter"] = from_diameter.model(morphology)
        input_distributions["diameter"]["method"] = diameter_model

    elif hasattr(diameter_model, "__call__"):
        input_distributions["diameter"] = diameter_model(morphology)
        input_distributions["diameter"]["method"] = "external"
    elif (
        isinstance(diameter_model, str) and diameter_model == "default"
    ) or diameter_model is None:
        input_distributions["diameter"] = {}
        input_distributions["diameter"]["method"] = "default"
        input_distributions["diameter"]["diameter"] = build_diameter_models(
            morphology, config={"models": ["simpler"], "neurite_types": neurite_types}
        )
    else:
        raise NotImplementedError(f"Diameter model {diameter_model} not understood")

    for neurite_type in neurite_types:
        if isinstance(feature, str):
            type_feature = feature
        else:
            type_feature = feature.get(neurite_type, "path_distances")
        nm_type = getattr(NeuriteType, neurite_type)

        input_distributions[neurite_type] = _append_dicts(
            trunk_neurite(pop_nm, nm_type), number_neurites(pop_nm, nm_type)
        )
        if type_feature in ["path_distances", "radial_distances"]:
            _append_dicts(
                input_distributions[neurite_type],
                persistent_homology_angles(
                    pop_tmd,
                    threshold=threshold_sec,
                    neurite_type=neurite_type,
                    feature=type_feature,
                ),
                {"filtration_metric": type_feature},
            )
    return format_values(input_distributions)
arnaudon commented 1 year ago

Hello! Thanks for trying to use this code! Here is a bit of advice before going to 'wild'.

Did you modify the code from a git clone? If yes, you could do git reset --hard origin/main to revert your changes, and try again to send us the original error. if you modified the files directly in your install dir, you can re-install with pip to revert your changes. In general, if you report a bug on a code on github, it should be the code that is on github, otherwise, we cannot help you. If you want to fix a but you found, you can open a Pull Request on, and we can then look at you have modified and see if it can be included in the code.

Alternatively, we can setup a call to discuss what you want to achieve, and if we can be of any help. Good luck with your research!