.. _smat_intro: Custom score matrices ********************* NBLAST calculates the similarity between neurons' morphology. For more information on NBLAST in ``navis``, see the other tutorials. The core of the algorithm is the function which converts point matches (a distance and a dot product) into a score expressing how likely they are to have come from the same point cloud. This is generally a 2D lookup table, referred to as the score matrix, generated by "training" it on neurons known to be matching or non-matching. ``navis`` provides (and uses by default) the score matrix used in the original publication (Costa et al., 2016) which is based on FlyCircuit, a light-level database of Drosophila neurons. This works quite well in many cases. However, how appropriate it is for your data depends on a number of factors: * How big your neurons are (commonly addressed by scaling the distance axis of the built-in score matrix) * How you have pre-processed your neurons (pruning dendrites, resampling etc.) * Your actual task (matching left-right pairs, finding lineages etc.) * How distinct you expect your matches and non-matches to be (e.g. how large a body volume you're drawing neurons from) Utilities in ``navis.nbl.smat`` allow you to train your own score matrix. .. code:: ipython3 import navis from navis.nbl.smat import Lookup2d, LookupDistDotBuilder, Digitizer :class:`~navis.nbl.smat.Lookup2d` is the lookup table, and can be given to any nblast-related class or function. Each axis represented by a :class:`~navis.nbl.smat.Digitizer`, which converts continuous values into discrete indices. These indices are used to look up score values in an array. The :class:`~navis.nbl.smat.LookupDistDotBuilder` is a class which generates :class:`~navis.nbl.smat.Lookup2d` instances from training data. First, we need some training data. Let's augment our set of example neurons by randomly mutating them: .. code:: ipython3 import numpy as np # use a replicable random number generator RNG = np.random.default_rng(2021) def augment_neuron(nrn: navis.TreeNeuron, scale_sigma=0.1, translation_sigma=50, jitter_sigma=10): nrn = nrn.copy(deepcopy=True) nrn.name += "_aug" dims = list("xyz") coords = nrn.nodes[dims].to_numpy() # translate whole neuron coords += RNG.normal(scale=translation_sigma, size=coords.shape[-1]) # jitter individual coordinates coords += RNG.normal(scale=jitter_sigma, size=coords.shape) # rescale mean = np.mean(coords, axis=0) coords -= mean coords *= RNG.normal(loc=1.0, scale=scale_sigma) coords += mean nrn.nodes[dims] = coords return nrn original = list(navis.example_neurons()) jittered = [augment_neuron(n) for n in original] dotprops = [navis.make_dotprops(n, k=5, resample=False) for n in original + jittered] matching_pairs = [[idx, idx + len(original)] for idx in range(len(original))] The score matrix builder needs some neurons as a list of :class:`navis.Dotprops` objects, and then to know which neurons should match with each other as indices into that list. It's assumed that matches are relatively rare among the total set of possible pairings, so non-matching pairs are drawn randomly (although a non-matching list can be given explicitly). Then it needs to know where to draw the boundaries between bins in the output lookup table. These can be given explicitly as a list of 2 `Digitizer`s, or can be inferred from the data: bins will be drawn to evenly partition the matching neuron scores. The resulting ``Lookup2d`` can be imported/exported as a ``pandas.DataFrame`` for ease of viewing and storing. .. code:: ipython3 builder = LookupDistDotBuilder( dotprops, matching_pairs, use_alpha=True, seed=2021 ).with_bin_counts([8, 5]) smat = builder.build() as_table = smat.to_dataframe() as_table .. parsed-literal:: [(7, 5), (7, 4), (4, 8), (9, 3), (6, 4), (5, 7), (4, 9), (2, 9), (5, 9), (6, 1), (5, 4)] .. raw:: html
[0.0,0.1456628620624542) [0.1456628620624542,0.29889603853225716) [0.29889603853225716,0.4806328892707825) [0.4806328892707825,0.7351819753646852) [0.7351819753646852,0.9988406300544739)
[2.014085054397583,57.849366188049316) 0.490599 0.379599 0.365107 0.466307 0.860021
[57.849366188049316,81.31283569335938) 0.408948 0.335051 0.405051 0.419681 0.807174
[81.31283569335938,104.08576202392578) 0.262346 0.305339 0.216944 0.254417 1.012022
[104.08576202392578,128.14104461669922) 0.263291 0.096855 0.164301 0.104572 1.304963
[128.14104461669922,155.36119651794434) -0.140083 -0.100454 -0.129114 0.250419 1.989136
[155.36119651794434,202.6728515625) -0.577435 -0.447474 -0.476505 -0.221284 1.413861
[202.6728515625,395.9569091796875) -1.027309 -0.941450 -0.863832 -0.623969 0.051871
[395.9569091796875,4709.61474609375) -0.618927 -0.731949 -0.679048 -0.456665 -0.197261
Now that we can have this score matrix, we can use it for a problem which can be solved by NBLAST: we've mixed up a bag of neurons which look very similar to some of our examples, and need to know which they match with. .. code:: ipython3 original_dps = dotprops[:len(original)] new_dps = [ navis.make_dotprops( augment_neuron(n), k=5, resample=False ) for n in original ] RNG.shuffle(new_dps) result = navis.nblast( original_dps, new_dps, use_alpha=True, scores="mean", normalized=True, smat=smat, n_cores=1, ) result.index = [dp.name for dp in original_dps] result.columns = [dp.name for dp in new_dps] result .. parsed-literal:: NBLAST is optimized for data in microns and it looks like your queries are not in microns. NBLAST is optimized for data in microns and it looks like your targets are not in microns. .. parsed-literal:: .. raw:: html
754538881_aug 722817260_aug 1734350908_aug 1734350788_aug 754534424_aug
1734350788 -0.271504 -0.361381 -0.272018 0.159048 -0.325920
1734350908 -0.400432 -0.491376 0.858478 -0.437819 -0.216206
722817260 -0.157433 0.127931 -0.365195 -0.407436 -0.393940
754534424 -0.237769 -0.413532 -0.193794 -0.401148 0.857898
754538881 -0.021153 -0.213349 -0.200214 -0.303390 -0.115017
You can see that while there aren't any particularly good scores (this is a very small amount of training data, and one would normally preprocess the neurons), in each case the original's best match is its augmented partner.