Threshold and cells

[auesro]

Dear @AlbertPun and @dfriedma ,

These are more questions than issues, I hope it fits here.

-I was wondering how does the network behave in the presence of neuronal bodies present in the images, is it able to remove them? Or should we clean them up before inference?

-The threshold variable is set at a default of 0.01...what does it represent? As far as I understand, if the chunk to be processed in the network has a maximum value higher than threshold then it will be added to the queue. Fine. But what does this value represent? Whatever shade that is not pitch black, right? Is it in float32 or int16?

Thanks!

A

[dfriedma]

(1) It doesn't recognize cell bodies, just fibers. In some cases, when the cell bodies are very bright, it can impact the ability to identify the most local fibers nearby, but I suspect this could be adjusted with a little extra training. Future versions will hopefully have multi-class output for cells and fibers from the same channel.
(2) You're on the right track. The threshold is to exclude chunks that are outside the brain to speed up inference. Since the sample values are normalized between 0-1, a threshold of 0.01 excludes cubes with a max value < 656 (assuming 16-bit raw data).

[auesro]

Thanks for the replies. Very instructive and clear.

I have been testing it on my own dataset and facing some issues which might be related to the acquisition:
-I am looking at the barrel cortex axons coming from the thalamus. Thats a very dense projection of axons where the network is having trouble (to be honest it would be impossible also for me to tell individual axons apart). Do you have any type of solution or future ideas to deal with densely packed fibers? The obvious solution, of course, is to image at higher resolution (2.5X instead of 0.8X for example) but that involves much bigger datasets...
-What function are you using for skeletonization in 3D? Skimage skeletonization_3d? I havent found a perfect solution here. This one works well but...I have the impression I see too many disconnected small segments.

[auesro]

If you prefer I can open a new issue for this but since I mentioned using my dataset, I just came across a weird pattern in the segmented output of the network:

In both test images, the raw data (inverted LUT) on the left, segmented corresponding plane on the right...what is that squared pattern that you can see in the segmented frames? It appears as if the network fails to...properly align back the batches (each "square" is 36by36 pixels)? I can see the checker-board pattern in every single image of the segmented stack. No idea what is going on.

[dfriedma]

-I am looking at the barrel cortex axons coming from the thalamus. Thats a very dense projection of axons where the network is having trouble (to be honest it would be impossible also for me to tell individual axons apart). Do you have any type of solution or future ideas to deal with densely packed fibers? The obvious solution, of course, is to image at higher resolution (2.5X instead of 0.8X for example) but that involves much bigger datasets...

Looking at your examples and hearing your description, I would agree that it’s too dense for the current model. The barrel cortex projections we show in Fig. 4 were imaged at higher optical zoom. As with the layer I cortical axons in Fig. 4, if a “fiber type” is not covered in the training set, it won’t be segmented. While we have good coverage of brain regions and background textures, different cell types will inevitably have different appearances in their terminals, necessitating some transfer learning.
As we expand our user base, we hope to be able to update our model as it improves with new training data.

-What function are you using for skeletonization in 3D? Skimage skeletonization_3d? I havent found a perfect solution here. This one works well but...I have the impression I see too many disconnected small segments.

Yes, that’s the skeletonization function we use. It’s up to you if you want to use our weighted design or not. Small disconnected objects can be identified and later categorized and removed with something like MATLAB’s bwconncomp and regionprops3.

In both test images, the raw data (inverted LUT) on the left, segmented corresponding plane on the right...what is that squared pattern that you can see in the segmented frames? It appears as if the network fails to...properly align back the batches (each "square" is 36by36 pixels)? I can see the checker-board pattern in every single image of the segmented stack. No idea what is going on.

That checkerboard pattern can appear when and object or texture dominates a single cube to the point that the network “gives up” on finding axons in the region. In your examples, the very dense barrels are not being identified as axons and are dominating nearby voxels. Some transfer learning with denser labeling or higher mag imaging might help.

[auesro]

Thanks for your reply @dfriedma
I will play with bwconncomp and regionprops3 to clean up disconnected components.
Hopefully I will not have such dense projections in my experimental data (this was a test brain) but if I find it, I will try to train the network on them. Regarding the training, is the current network available here trained only in datasets with 4.0625 um/px resolution? Or does it include some data also at higher resolution? Does it matter at all actually?

[dfriedma]

It all is at 4.0625 um/vx, but in practice, augmentation allows for axons of different apparent size. It's likely that there just has to be a little separation between the fibers in order to identify them as fibers and not splotches in background. Your clearing quality looks nice! It'll work out. I'm working on getting a new model for higher density regions now as well, but haven't finalized things yet.

[auesro]

Sounds good, looking forward!

[auesro]

Reopening this, because I came up with another question:
How do you deal with the variability in fluorescence intensity of axons and background among different brains? I am thinking at 2 levels:
-For the network, I guess it will take into account the difference between background and signal to determine whether a row of pixels is an axon or not, right? I mean that it might be difficult for the network to identify axons if the brain has a high level of background compared to another with the same axonal intensity but lower background.
-For quantification purposes comparing several brains, I saw that in Figure 7E from "Serotonin Neurons in the Dorsal and Medial Raphe Nuclei: from Single-Cell Transcriptomes to Whole-Brain Projections" for example, according to the figure legend I gather you take the cleaned up skeleton and normalize each region to the region volume and to "total axon content per brain"...what is exactly this measure and how do you calculate it?

[saskra]

-What function are you using for skeletonization in 3D? Skimage skeletonization_3d? I havent found a perfect solution here. This one works well but...I have the impression I see too many disconnected small segments.

Yes, that’s the skeletonization function we use.

We use this function, too, but are currently looking for an alternative due to this bug: scikit-image/scikit-image#3757 Do you think this could influence the results here?

[albert597]

Reopening this, because I came up with another question:
How do you deal with the variability in fluorescence intensity of axons and background among different brains? I am thinking at 2 levels:
-For the network, I guess it will take into account the difference between background and signal to determine whether a row of pixels is an axon or not, right? I mean that it might be difficult for the network to identify axons if the brain has a high level of background compared to another with the same axonal intensity but lower background.
-For quantification purposes comparing several brains, I saw that in Figure 7E from "Serotonin Neurons in the Dorsal and Medial Raphe Nuclei: from Single-Cell Transcriptomes to Whole-Brain Projections" for example, according to the figure legend I gather you take the cleaned up skeleton and normalize each region to the region volume and to "total axon content per brain"...what is exactly this measure and how do you calculate it?

Very sorry for such a late reply, I only saw you reopened this now, but I thought it's worth answering if anyone has the same question! We did not normalize the volumes because we lost too much information (raw intensity is important for axon identification). For example, a noisy image of very low intensity would begin to look like axons if you normalized it.

Instead, we took our training set and augmented with random scaling and constants additions. This created a training set that was diverse enough in axon intensity to perform well on our validation set.

For quantification, if you see the skeletonization_3d is leading to disconnected components, you can lower the confidence from 0.5 to say 0.3. If you are still unhappy with the results, you can also use some classical post-processing filters such as removing small components or filtering out shapes.

Let me know if you have anymore questions!

[auesro]

Hi Albert, thanks for the reply!

Only if you could reply to the second question, that would be great:

-For quantification purposes comparing several brains, I saw that in Figure 7E from "Serotonin Neurons in the Dorsal and Medial Raphe Nuclei: from Single-Cell Transcriptomes to Whole-Brain Projections" for example, according to the figure legend I gather you take the cleaned up skeleton and normalize each region to the region volume and to "total axon content per brain"...what is exactly this measure and how do you calculate it?

Thanks!

[dfriedma]

I can speak to the quantification strategy--
We have two separate normalizations in the analysis. Step (1) is to the volume of the annotated region. e.g. the CP is huge compared to the SCN and we want to start with a "density" metric rather than a "total axon length" metric. Step (2) is to the total sum of cleaned-up skeleton. In practice this levels out variability across animals due to injection volume for example. It does introduce a bias in the interpretation of the results (which we viewed as favorable, but depending on your data it may be better to leave this step out). Namely, a projection pattern with very specific projections only to very small nuclei will have inflated numbers. This reduces the ability to compare "densities" across animals, but emphasizes instead a value related to the "importance" of that target region in this specific set of projections. Hope that helps!