- Home
- Neuroconductor Tutorials
- Tissue class segmentation
Tissue Class Segmentation
John Muschelli
2021-02-16
All code for this document is located at here.
In this tutorial we will discuss performing tissue class segmentation using atropos
in ANTsR
and it’s wrapper function in extrantsr
, otropos
.
Data Packages
For this analysis, I will use one subject from the Kirby 21 data set. The kirby21.base
and kirby21.fmri
packages are necessary for this analysis and have the data we will be working on. You need devtools to install these. Please refer to installing devtools for additional instructions or troubleshooting.
packages = installed.packages()
packages = packages[, "Package"]
if (!"kirby21.base" %in% packages) {
source("https://neuroconductor.org/neurocLite.R")
neuroc_install("kirby21.base")
}
if (!"kirby21.fmri" %in% packages) {
source("https://neuroconductor.org/neurocLite.R")
neuroc_install("kirby21.fmri")
}
Loading Data
We will use the get_image_filenames_df
function to extract the filenames on our hard disk for the T1 image and the fMRI images (4D).
library(kirby21.t1)
library(kirby21.base)
fnames = get_image_filenames_df(ids = 113,
modalities = c("T1"),
visits = c(1),
long = FALSE)
t1_fname = fnames$T1[1]
Using information from the T1 image
Brain extracted image
Please visit the brain extraction tutorial on how to extract a brain from this image. We will use the output from fslbet_robust
from that tutorial.
outfile = nii.stub(t1_fname, bn = TRUE)
outfile = file.path("..", "brain_extraction", outfile)
outfile = paste0(outfile, "_SS.nii.gz")
ss = readnii(outfile)
ss_red = dropEmptyImageDimensions(ss)
ortho2(ss_red)
Again, we can see the zoomed-in view of the image now.
Tissue-Class Segmentation with Atropos
outfile = nii.stub(t1_fname, bn = TRUE)
prob_files = paste0(outfile,
"_prob_", 1:3,
".nii.gz")
seg_outfile = paste0(outfile, "_Seg.nii.gz")
if (!all(file.exists(
c(seg_outfile, prob_files)
))) {
seg = extrantsr::otropos(
ss_red,
x = ss_red > 0,
v = 1)
hard_seg = seg$segmentation
writenii(hard_seg, seg_outfile)
for (i in seq_along(seg$probabilityimages)) {
writenii(seg$probabilityimages[[i]], prob_files[i])
}
# writenii(seg, )
} else {
hard_seg = readnii(seg_outfile)
seg = vector(mode = "list", length = 2)
names(seg) = c("segmentation", "probabilityimages")
seg$segmentation = hard_seg
seg$probabilityimages = vector(mode = "list", length = 3)
for (i in 1:3) {
seg$probabilityimages[[i]] = readnii(prob_files[i])
}
}
Atropos results
Now we have a hard segmentation, which assigns a class with the maximum probability to that voxel. We also have a separate probability image for each tissue class.
double_ortho(ss_red, hard_seg)
We see that much of the structures have been segmented well, but there may be errors.
Atropos intensity histograms
We can also look at the distribution of intensities (marginally) for each tissue class. In atropos
, the classes are ordered by mean intensity, so we can re-assign them to the corresponding tissue class
df = data.frame(value = ss_red[ss_red > 0],
class = hard_seg[ss_red > 0])
df$class = c("CSF", "GM", "WM")[df$class]
ggplot(df, aes(x = value, colour = factor(class))) + geom_line(stat = "density")