How to extract features from fMOST data with neuron_morphology¶
To obtain the package, please see https://neuron-morphology.readthedocs.io/en/readthedocs/.
This step-by-step guide will walk you through the process of extracting features from fMOST data.
- Load fMOST data from file
- Clean up swc to prevent errors
- Extract the neuron_morphology default_features
- Add a custom feature to calculate as well
import sys
sys.path.insert(0, "../")
from io import StringIO
import copy
import matplotlib.pyplot as plt
import neuron_morphology.swc_io as swcio
from neuron_morphology.morphology import Morphology
from neuron_morphology.swc_io import morphology_from_swc
from neuron_morphology.feature_extractor.data import Data
from neuron_morphology.feature_extractor.feature_extractor import FeatureExtractor
from neuron_morphology.features.default_features import default_features
from neuron_morphology.constants import (
SOMA, AXON, BASAL_DENDRITE, APICAL_DENDRITE
)
import json
import numpy as np
import neuron_morphology.feature_extractor.feature_writer as fw
Pick your neuron (.swc/.csv)¶
Here select one SWC file as an example.
fmost_swc_file = "../../tests/data/17545-6151-X24259-Y36270.swc"
Ensure one SWC with only one soma point¶
This is developed for current registered fMOST SWC file. It will adjusted and added as a module into neuron_morphology later.
def prepare_neuron_tree(swc_data):
nodes = swc_data.to_dict('record')
replace_type = 2 # default node type
for node in nodes:
node['parent'] = int(node['parent'])
node['id'] = int(node['id'])
node['type'] = int(node['type'])
if node['parent'] == -1 and node['type'] != 1:
replace_type = node['type']
if node['type'] == 1 and node['parent'] != -1:
node['type'] = replace_type
soma_list = []
for node in nodes:
if node['type'] == 1:
soma_list.append(node)
# create a new soma point
if len(soma_list) > 1:
x = 0
y = 0
z = 0
n = len(soma_list)
for node in soma_list:
x += node['x']
y += node['y']
z += node['z']
soma = copy.deepcopy(soma_list[0])
soma['id'] = nodes[-1]['id']
soma['x'] = x/n
soma['y'] = y/n
soma['z'] = z/n
nodes.append(soma)
for node in soma_list:
node['parent'] = soma['id']
node['type'] = replace_type
return nodes
Load SWC file as morphology data¶
We can load the SWC file into our morphology data object and then calculate features on it.
# load
swc_data = swcio.read_swc(fmost_swc_file)
nodes = prepare_neuron_tree(swc_data)
test_data = Data(Morphology(nodes, node_id_cb=lambda node: node['id'], parent_id_cb=lambda node: node['parent']))
# visualize
x = [node['x'] for node in nodes]
y = [node['y'] for node in nodes]
z = [node['z'] for node in nodes]
fig, ax = plt.subplots(1, 2)
ax[0].scatter(x, y, s=0.1)
ax[0].set_title('x-y view')
ax[1].scatter(z, y, s=0.1)
ax[1].set_title('z-y view')
Define interested features set¶
First, we instantiate a FeatureExtractor Second, we will use register the default_features
fe = FeatureExtractor()
fe.register_features(default_features)
Then, let's add a custom feature. We will use a feature included in the neurom_morphology package which calculates the number of stems exiting from the soma, by neurite type.
For reference, here is the entire function:
@marked(RequiresSoma)
@marked(RequiresRoot)
def calculate_number_of_stems(data: Data, node_types: Optional[List[int]]):
"""
Calculate the number of soma stems.
This is defined as the total number of non-soma child nodes on soma nodes.
Parameters
----------
data: Data Object containing a morphology
node_types: list (AXON, BASAL_DENDRITE, APICAL_DENDRITE)
Type to restrict search to
Returns
-------
Scalar value
"""
soma = data.morphology.get_soma()
return len(data.morphology.children_of(soma))This function has been marked to indicate that it needs a soma in order to be valid, and is also specializable based on node types. We will use the specialize() function and the NEURITE_SPECIALIZATION to tell the feature extractor to calculate this feature for each neurite_type (AXON, DENDRITE, BASAL_DENDRITE, APICAL_DENDRITE, ALL_NEURITES)
from neuron_morphology.features.soma import calculate_number_of_stems
from neuron_morphology.feature_extractor.marked_feature import specialize
from neuron_morphology.feature_extractor.feature_specialization import NEURITE_SPECIALIZATIONS
fe.register_features([specialize(calculate_number_of_stems, NEURITE_SPECIALIZATIONS)])
Extracting features¶
Now that we have registered the features that we are interested in, we can call feature extractor on our swc. This will create a dictionary of results, which we can unnest and print below.
Calling fe.extract() on another test_data will calculate the same feature set, so it can be used to process many swcs.
feature_extraction_run = fe.extract(test_data)
results = feature_extraction_run.results
from neuron_morphology.feature_extractor.utilities import unnest
unnest(results)
Save the features as h5 or csv file¶
We can also write these features to a file using the FeatureWriter. The FeatureWriter can write the features as a nested json or as a flat csv. If a feature returns a large amount of data, it can be output to a h5 file instead.
heavy_path = "test_features.h5"
table_path = "test_features.csv"
features_writer = fw.FeatureWriter(heavy_path, table_path)
features_writer.add_run("test", feature_extraction_run.serialize())
features_writer.write_table()