Introduction: Moyamoya Image Processing
Moyamoya, "puff of smoke," is a rare disease that is caused by the blockage of arteries at the basal ganglia, which is an area at the base of the brain. The disease is a progressive cerebrovascular disease that mostly affects children. Symptoms of moyamoya include an initial stroke, constant mini-strokes, muscular weakness, paralysis, or seizures as a result of the progressive narrowing of arteries. Without treatment, moyamoya will cause problems with speech, sensory impairment, and damaged consciousness. In our project, we will utilize MATLAB in order to preprocess an MRI or MRA image using various filters to reduce noise within the image to locate the affected area. In addition, we will use a feature enhancement to locate the affected areas more precisely. Moreover, we will then run an independent samples t-test to determine whether there is a significant difference between the amount of blood vessels in a normal brain compared to a brain affected by moyamoya.
Step 1: Find MRI and MRA Scans of a Normal Brain and a Brain Affected With Moyamoya
These images are the scans we used for the project that we found online. The two images with the blood vessels located in the middle are MRA scans, while the other two images are MRI scans.
The following links are where these images are found:
Step 2: Load Images Onto MATLAB and Assign Images to a Variable to Display Images
To begin the process, start by clearing the command window, close all the possible figures and graphs that may be already open and clear variables already assigned in your workspace.
After, create a for loop from 1 to 2 using command i = [1:2]
After, load the MRA images using the command imread(sprintf('filename%.filetype', i)) to read the images from the files specified by the file name followed by the number used for the loop using sprintf for batch loading and assign it to a variable.
Then to display the image in a figure, use the command imshow(I).
To assign a gray colormap, use the command colormap(gray).
To completely eliminate color and convert the 3D matrix for the images into 2D, use command rgb2gray(I) and assign it to a separate variable.
Then load the MRI images by using the command previously stated or imread(sprintf('filename%.filetype', i)) and assign it to a new variable
Repeat the rgb2gray command with the new variable used for the MRI images.
If necessary, you can resize an image using the command imresize(A, scale) and assign to a separate variable.
Step 3: Enhance Elongated Structures in Intensity Image With Multiscale Filtering
Using a new variable, use the command fibermetric(A) to enhance the tubular structures in the images
With the previous variable, use the command histeq(B) to enhance the histogram equalizations by transforming the intensity of the images and assign it to a new variable.
Display the histogram using the command imhist(B)
Create a new variable to create a threshold for the filter. In this case assign the previous variable> 0.875, filtering out pixel intensity under the value of 0.875
After, create a new figure and use the command imshow(A) to display the new filtered image.
Step 4: Run a 2D Median Filter
Using command medfilt2(A,[m n]), run a 2D median filter, where each output pixel contains the median value in the mxn boundary around the respective pixel in the input image.
Create a new figure and use imshow(A) to display the median filtered image.
Step 5: Mask the Image
Using the median filtered image, use the command [labeledImage, numberOfBlots] = bwlabel(A) to count the number of white blots in the image
Then, use the region props function states = regionprops(labeledImage, 'Area') to calculate the areas of each blot or blood vessel
Assign all the areas into one variable
Then using another variable, count the number of blots that exceed 50 pixels
After, sort out any blots that are below 50 pixels in descending order using command [sortedAreas, sortedIndicies] = sort(Areas, 'descend')
Then, using another variable, use the command ismember(labeledImage, sortedIndicies(1:numberToExtract)) to return an array with elements of labeledImage are found in sortedIndicies from number 1 to the number of blood vessels to return a logical 1 (true) or a logical 0 (false).
With the variable in the previous step, find the points that are true (values > 0) and create a logical array to make a binary image and assign it to a new variable.
Create a new figure and use imshow(A) the new binary image.
Then, invert the image using command imcomplement(A) and assign it to a different variable.
To create a masked image, use a new variable with command resizedimage.*uint8(invertedimage)
Create a new figure and use imshow(A) to display masked image.
To end the entire code, be sure to use the command 'end' to end the entire for loop
Step 6: Select the MRA Scans for Statistical Testing
To prepare for statistical testing, select the MRA scans that are to be used for the independent samples t-test. Because our two samples will be brains affected Moyamoya, and normal brains, select a decent amount of MRA scans of each group.
Step 7: Calculate the Area of Blood Vessels in Preparation for Statistical Testing
The statistical test will focus on the length or amount of blood vessels shown on the MRA scans. Thus, we must calculate the area of the blood vessels before comparison.
Start by filtering the MRAs of normal brains and calculating the amount of blood vessels. To do this, run a for loop. Because there are three images, the condition will be i = [1:3].
Open the image with the imread command and assign it to a variable.
Next, create an if/else statement with the if, else command. For the if statement, use the command size(A,3)==3, where A is the variable used to open to the image, to create an if statement for when the third dimension of the array is 3. Then convert the image to 2D and get rid of color using the command rgb2gray(A) and assign it to a new variable. Use the command imresize(A, [m n]) to resize the image. In this case, we resized the images to 1024 x 1024 matrix. To enhance the tubular structures of the image, use the fibermetric command again, and assign it to a new variable.
The following is for the else statement. If the image is not a 3D matrix, we want to skip conversion. Do the same as the if statement, but without the rgb2gray(A) command.
Create a new variable, setting it equal to the variable from the fibermetric step greater than 0.15. This thresholds the image for intensities greater than 0.15.
We will repeat the lines of codes from steps 4 and 5 of the instructable from the median filter line until the imshow(I) line. After, use the command sum(I(:)) to add up all pixels that make up blood vessels, and assign it to a separate variable. Name a new variable NormalBloodVessels(i) and set it equal to the variable from the sum(I(:)) command. This adds the data to the matrix.
End the loop and repeat but for the MRAs of brains affected by Moyamoya. Name the variable in the end MoyaMoyaBloodVessels(i) to not confuse it with the normal brain MRAs.
Step 8: Run an Independent Samples T-test
Since there are two independent samples, and a small population, run an independent samples t-test.
Create and name a function that runs an independent samples t-test to determine whether the amount of blood vessels in the MRAs of normal brains are significantly equal or not to that of MRAs of brains affected by Moyamoya.
Display the hypothesis set for the test by using the command disp('X'). On the first line, display, "Hypotheses for two sample t test." On the second line, display, "H0 = The amount of blood vessels of a normal brain equals the amount of blood vessels of a brain with Moyamoya disease," to state the null hypothesis. On the third line, display, "HA = The amount of blood vessels of a normal brain does not equal the amount of blood vessels of a brain with Moyamoya disease." to state the alternative hypothesis.
Using a 95% confidence interval and a sample size of 3, calculate the t score by using the command tinv([0.025 0.975], 2) and assign to variable t. Use command var(NormalBloodVessels) and var(MoyaMoyaBloodVessels) and assign them to variables to calculate the variances of both data sets.
Check if the variances are close to equal or not. Do this by creating an if/else statement with the if, else command. For the condition in the if statement, write A / B == [0.25:4], where A is the variable that accounts for the variance of normal blood vessels and B is the variable that accounts for the variance of Moyamoya blood vessels. 0.25 and 4 come from a general estimate for determining whether variances are equal or not. Then run the two sample t test with [h,p] = ttest2(A, B, 0.05, 'both', 'equal'), with A and B being the same variables as mentioned before. For the else statement, use [h,p] = ttest2(A, B, 0.05, 'both','unequal') to run a two sample t test in the case where the variances are not equal. End the if/else statement. This will calculate p.
Create an if/else statement that will display a conclusion based on the value of p. The condition for the if statement will be p > 0.05. Since we typically fail to reject the null hypothesis when the value of p is greater than 0.05, use the disp('X') command to display "Because the p value is greater than 0.05 we fail to reject the null hypothesis," and "Therefore we fail to reject that the amount of blood vessels of a normal brain equals that of a brain with Moyamoya disease." In the else statement, since we typically reject the null hypothesis when the value of p is under 0.05, use the disp('X') command to display "Because the p value is less than 0.05 we reject the null hypothesis," and "Therefore we fail to reject that the amount of blood vessels of a normal brain not equal to that of a brain with Moyamoya disease." End the if/else statement.
