Easy way to strip skull form MRI volume?

I'm trying to create a code that separates the skull from the brain for a Serie of MRI Dicom files. I tried to use thresholds combined with selecting the biggest component combined with a mean intensity, but this only works partially. Does anybody has suggestions/knows a method to perform or optimize this kind of operation?

[Walter Roberson]

Have you looked at SPM for this purpose?

[sam CP]

I want to remove skull from brain using morphological operations and thresholding.. i have wrote the following code ,but i didn't get skull stripped MRI.. Please guide me..

clear all; clc [filename pathname] = uigetfile({'*.jpg';'*.bmp'},'Select MRI'); inputimage=strcat(pathname, filename); I = imread(inputimage);

figure,imshow(I), title('input');

se = strel('disk',7);

closeBW = imclose(I,se); figure, imshow(closeBW);title('close operation');

afterOpening = imopen(closeBW,se); figure, imshow(afterOpening,[]);title('open operation');

J=I-afterOpening; figure, imshow(J,[]);title('Skull of MRI');

K=I-J; figure, imshow(K,[]);title('Skull removed MRI');

[Image Analyst]

You forgot to attach the image it failed on. See attached demo.

[sam CP]

i have attached the skull removed MRI ( skull is not removed in this) after running the above code..

[Image Analyst]

Start a new question, and this time, attach the original image, not a screenshot with title, toolbar, etc.

[sam CP]

I want to remove skull from brain using morphological operations and thresholding.. i have wrote the following code ,but the skull didn't extracted properly from the MRI.Hereby i will attach the original MRI input image..Please help me to find a proper output..

clear all;
 clc
 [filename pathname] = uigetfile({'*.jpg';'*.bmp'},'Select MRI'); inputimage=strcat(pathname, filename); I = imread(inputimage);

figure,imshow(I), title('input');

se = strel('disk',7);

closeBW = imclose(I,se); 
figure, imshow(closeBW);title('close operation');

afterOpening = imopen(closeBW,se);
 figure, imshow(afterOpening,[]);title('open operation');

J=I-afterOpening; 
figure, imshow(J,[]);title('Skull of MRI');

K=I-J; 
figure, imshow(K,[]);title('Skull removed MRI');

[Walter Roberson]

Sam, please start a new Question.

[Image Analyst]

See skull stripping code in your question https://www.mathworks.com/matlabcentral/answers/323816-i-want-to-extract-skull-from-brain-mri-using-morphological-operations-but-i-didnt-get-a-useful-out#answer_255279

[sam CP]

Thank you

Answer 1

[Image Analyst]

Without seeing the image, or knowing why you said it didn't work well, that's what I would have suggested. So, other than that, all I can say is to check the published literature: http://iris.usc.edu/Vision-Notes/bibliography/contentsmedical.html#Medical%20Applications,%20CAT,%20MRI,%20Ultrasound,%20Heart%20Models,%20Brain%20Models

[Hugo]

Thanks, the method i tried worked fine on a pre-selected 2D-image, but when i selected a series of images to create a map of the whole brain in some of the images only a part of the brain of appeared, furthermore the skull was still partially visible. The biggest problem is selecting a threshold for the mean intensity, this should rule out any labeled component (in this case the skull) with a higher mean intensity, unfortunately it didn't work, so i'm looking for other methods which aren't depending to much on the intensity of a region.

After applying thresholds i used the following loop:

if true
  % for i=length(fnum):-1:1;
  MRIsubtract = MRI1(:,:,i);
  CC = bwconncomp(MRIsubtract,4);
  STATS = regionprops(CC, MRIsubtract,'PixelIdxList','MeanIntensity');
  C = struct2dataset(STATS);
  C.Intensity = (C.MeanIntensity<4.2e+03);
  numPixels = cellfun(@numel,C.PixelIdxList);
  intensity = C.Intensity;
  C.Biggest = numPixels .* intensity;
  Biggest = C.Biggest;
  [biggest,idx]= max(Biggest);
  MRIsubtract(CC.PixelIdxList{idx}) = 0;
  MRIsubtract2 (:,:,i) = MRIsubtract;
  i=i+1
end

end

[Image Analyst]

First of all you threshold before you label it, so this doesn't make sense: "The biggest problem is selecting a threshold for the mean intensity, this should rule out any labeled component". Next, assuming the threshold is good (which you say is a bad assumption) you can just throw out the outermost object which should be the skull. Next, I don't know what the strcut2dataset is used for - I don't think that's necessary at all, you can just use STATS to do everything. I don't even have that function at all. With your C variable, I don't see you referencing an index - aren't there multiple blobs in the thresholded image? There should be. You say the thresholding method is no good but you don't give the algorithm. Thresholding divides the image into two classes, foreground and background, so it's impossible to do a thresholding without considering the intensity of the image, in at least some part of the image. You can threshold based on some other image that is derived from the image, such as a texture image you mgiht get from stdfilt(), but you're still thresholding based on intensity, you're just doing it on an intermediate image. (By the way, I don't think texture filters will help you find the threshold of the skull).

[Hugo]

Thanks a lot, I used thresholds to filter out the very high and low sequences after which i combined the original volume with the binary threshold/filter images. Thereafter i labeled them and tried to select the preferred label. When i applied labels before filtering some information out the brain and skull were still connected in some of the images. But nevertheless i'm gonna give it another try. I realize the method i described before is creating more areas and makes it more difficult to label and select. So i'm gonna try i the way you described once more

I used the STATS and the C variable to create a dataset in which i was enabled to use a threshold to determine which labels i wanted to use, there were multiple blobs, but the amount of blobs differs between the different images which makes it difficult to just select or cut one out. It worked fine for 1 image, but i think i will try your method for the whole serie.

The code i used to filter (3000 is the lowest intensity i found in the brain and 6000 the highest.)

if true
  % MRI1 = MRI;
MRIsub = MRI>3000 & MRI<6000;
MRIsub2 = bwareaopen(MRIsub,25,4);
MRI1(MRIsub2) = 0;
end

[Image Analyst]

If the skull is disconnected from anything else, and there are no blobs outside the skull, you can just use ismember to create a new labeled image (as shown in my Image Segmentation Tutorial at http://www.mathworks.com/matlabcentral/fileexchange/?term=authorid%3A31862)

newLabeledImage = ismember(labeledImage, 2:numberOfBlobs);

Or, since you know it's blob #1, you can just zero it out directly

newLabeledImage = labeledImage; % Initialize
newLabeledImage(labeledImage==1) = 0;

If it's not blob #1, you can get rid of small blobs outside it with bwareaopen(), or you can try to identify it based on other things like the EquivDiameter, MajorAxisLength, Perimeter, etc. But if the skull connects to the inside, it's a lot more difficult.

[Hugo]

Thanks again, this last one worked fine for me. The skull was connected to the inside, but after applying a threshold and removing some little objects bwlabel only selected the skull as blob #1.