Introduction: Brain Box: Tracking Neural Volume Over Time

By
About: University of California, Irvine - 2nd Year Biomedical Engineering Student

The advance into the frontier of longer human life has brought forth the rise of diseases not seen by civilizations before ours. Among these, Alzheimer's affected approximately 5.3 million living elderly Americans in 2017, or approximately 1 in 10 elderly Americans (https://www.alz.org/facts/) and countless others with dementia. To aid in the fight to understand what is afflicting our elders, this code will equip the future researchers and eager curious with the ability to track brain volume over time.

Step 1: Using Brain Box

To use brain box, one only needs the following:

  • MRI scans of a brain and the name and format of such files (should all have roughly the same dimensions)
  • Length of one scan
  • Distance between each layer (MRI scan)
  • Patient Name (When inputting, do not include spaces and please capitalize first name and last name, like such: FirstnameLastname)

And from this, one has the ability to keep track of an individuals trends in brain volume over time. Thus, figures for Alzheimer's trends can be tracked by this software. The length we used in trial were 180 mm for length of one scan and 5 mm for distance between MRI scans, based off of average figures.

However, the application of brain box need not be constrained to this one task. If the cross sections of a given solid are photograph, like a tumor per se, the trends in changes in volume for these can as well be tracked in the software.

Step 2: Intro: Analyzing Cross Sections

In three dimensional structures, the two dimensional planes that such consist of would be called cross sections. Imagine that a stack of papers make up a rectangular prism, then each piece of paper would be a cross section of the paper. In imagining the brain, we apply the same course of thought. MRI (magnetic resonance imaging) (see information on MRI) captures the cross sections of the brain, and from using the boundaries defined in each "layer" of the brain provided, we can construct a structure to model and find the volume of the brain. We must first build a function to provide information on such boundaries however.

Step 3: Setting Up a Function: LevelCurveTracings.m

First, make sure your computer has MATLAB_R2017b downloaded(download here) and open up MATLAB. In the MATLAB interface, click on the button in the upper left corner of the window that says "New" with a bold yellow plus sign, and select the option "function", to open in the editor window a space that resembles that in the third picture. We will focus on changing the first line to setup the function. Where it says "outputArg1", replace it with "brain", "outputArg2" to say "holes", "untitled2" to "exp2", and "inputArg1" to "image", and delete "inputArg2". You now have a function to be called using "exp2", taking one argument "image" and outputting the boundaries of "brain" and "holes". The first line of the function should resemble the line depicted in the fourth picture. Delete all the code below this initial line.

Step 4: Developing the Bound Algorithm: Finding Bounds

Type in the code as follows below the line. This section of the function does the following line-by-line.

  1. Load in the image "image" into the variable "mri".
  2. Turn "mri" into an image made of values in a range of numbers to ones and zeros (aka binarizing) based on a set threshold value. If the value in a pixel is equal to or greater than 0.1, it is set to one, if not, the value at that pixel is set to zero.
  3. The following four lines turn 10 columns and rows at the edges of the MRI scan layer into zeros, to avoid reading improper values as making bounds (as learned from experimenting with the code).
  4. In the final line, bwboundaries traces the bounds of the binarized image "mri" and sets it equal to "b", an array with the elements whose indices correspond to those of the bounds set to one.

Step 5: Developing the Bound Algorithm: Generating Outer Bound Array

Follow along in the editor window with the following code in the picture. This section of the code does the following line-by-line.

  1. Find the length of each of the rows of the binarized image "b" (cellfun applies the function length to each row).
  2. Set "loc" to store the maximum lengths.
  3. Find the index of the maximum length, set to store in "largestTrace".
  4. Find the size of the image "mri", which consists of the same size as "b", and set to "BWsize".
  5. Find the number of rows in the image's array, set to "ysize".
  6. Find the number of columns in the image's array, set to "xsize".
  7. Generate array "largestTraceMat", a "ysize" by "xsize" matrix of zeros.
  8. Find the equivalent index from the subscripted values corresponding to where the largestTrace x values and y values were, store in vector "lindex".
  9. In the matrix of zeros, "largestTraceMat", turn the elements at the indices that correspond to index values stored as elements in "lindex" into ones.

Thus the logical array "largestTraceMat" has the largest bounded region of the given brain scan cross section plotted as ones with a backdrop of zeros

Step 6: Developing the Bound Algorithm: Working With Center Point

Next, we must test to see whether the cross section consists of more than one region(the largest). By testing the alignment of the largest region's centroid, we can see whether there is one contiguous region, which would yield a more centered centroid, or the possibility of multiple regions.

  1. Use "regionProps" to find information on the centroids present, set equal to the structure array "tempStruct"
  2. Form array "centroids" with data from the field "centroid" concatenated vertically
  3. Take the second column values of "centroids" (the horizontal dimension coordinates)
  4. Run a filter to check the alignment of the centroid to the horizontal center

Step 7: Developing the Bound Algorithm: When a Centroid Is Not Centered

In the scenario that the centroid of the largest trace region is not centered, we go through the following steps. As we had observed in the MRI scans, the tendency was to have hemispheres of the brain portrayed in the cross section when not contiguous, so we now continue to plot the second largest trace along with the largest trace in "largestTraceMat"

  1. Set the traced matrix to a new variable "b2"
  2. Initialize empty matrix "b2", with a set indexed by "loc"
  3. Create a conditional, for when a centroid is not centered (i.e. a multiple region layer)
  4. Set a new trace size to be found for each row (traceSize2)
  5. Set "loc2" to find the indices where bounds are present
  6. Let cells specified by "loc2" in "b2" be equal to "largestTrace2"
  7. Convert the subscripts into indices, set to "lindex"
  8. Change elements corresponding to "lindex" in "largestTraceMat" to 1
  9. Initialize empty matrix "b2", with a set indexed by "loc2"

Step 8: Developing the Bound Algorithm: Interference of Holes

In dealing with holes, the values stored in "b2" kept track of structures other than the largest trace, and plotting these onto a filled form of "largestTraceMat" will reveal where there are holes in the brain regions.

  1. Create array "filledMat", which is a filled in form of "largestTraceMat"
  2. Create array "interferenceMat", a "ysize" by "xsize" array of zeros
  3. Create array "interferenceloc", to store the values from "b2", concatenated vertically
  4. Create array "lindex" to store the indices that correspond to "interferenceloc"
  5. For indices in "interferenceMat" that correspond to "lindex", set value to 1, making a different bounded region

Step 9: Developing the Bound Algorithm: Locating Holes, Finalizing Brain Bounds and Hole Bounds

  1. Set array "tempMat" equal to "interferenceMat" plus "filledMat", thus adding each value in the matrix to each other
  2. Set array "holesLoc" equal to the indices where "interferenceMat" and "filledMat" both were equal to one
  3. Setup "holesMat" as a zero matrix of dimensions "ysize" x "xsize"
  4. Set indices in "holesMat" that are equal to "holesLoc" as ones
  5. Set "brain" to "largestTraceMat"
  6. Set "holes" to "holesMat"

With the finding where values of the added matrices were equal to 2, the hole locations were easily secured and plotted onto an empty matrix.

Step 10: Logging Data: Function PatientFiles.m

Much like the setup of the last function, click on the button in the upper left corner of the window that says "New" with a bold yellow plus sign, and select the option "function", to open in the editor window a space that resembles that in the third picture. In the first line, delete the output matrix and replace with merely "output", replace "untitled2" with "patientFiles", delete all the input arguments, and instead follow the formatting specified in the fourth picture of the line of code. The first line of this function should match the formatting of the picture.

Step 11: Logging Data Into Files

To set up a file to log the data found by the main function (yet to be described), we must follow these steps (as prescribed by the code line-by-line).

  1. Check if the input for patientName is a string.
  2. If it is not a string, display that patientName input should be a string.
  3. End the if statement (prevent error).
  4. Set up a string statement "DateandTime" that will give the following format: hour:minutes--month/day/year.
  5. Set variable fileName to the following: patientName.m.

Now to the next section of the function: Does a file of this name already exist?

1)Suppose the file of this name already exists:

  1. Run the file to get the values from the past queued up
  2. Add the "DateandTime" data of the current iteration as a new cell in the cell array of x values (index end+1)
  3. Add the current "brainVolume" value as a new cell in the cell array of y values (index end+1)
  4. Save the current variables loaded in the file.

2)Suppose the file of this name does not exist:

  1. Create a new file with the name stored in the variable "patientName"
  2. Add the current "DateandTime" data as a cell into the empty cell array of x values
  3. Add the current "brainVolume" data as a cell into the empty cell array of y values
  4. Save the current variables loaded in the file.

Step 12: Logging Data: Displaying a Plot of Brain Volume Over Time

  1. Convert the x values array (xVals) into a categorical array (xValsCategorical), to allow plotting
  2. Generate figure window 5
  3. Plot the points designated by "xValsCategorical" and "yVals"(containing brain volume), using hollow circles to indicate points and to be connected by dashed lines
  4. Title the plot as: patientName Brain Volume Data
  5. Label the x axis as shown in the picture
  6. Label the y axis as shown in the picture
  7. Let figure 5 equal to output

From this, the function patientName being called will yield a file with edited data keeping track of brain volume over time, and a plot displaying trends.

Step 13: Closing Gaps in Subplots: Subplotclose.m

The function, adapted from code from http://www.briandalessandro.com, functions to close the gaps in between the subplot figures of the main code, when the figures displaying the MRI images and the brain layers are created. The subplot function used within subplotclose.m adjusts the position of the given subplots to fit snugly agains one another in the aspect of the longer dimension. For example, if the code intends a 7 x 3 matrix, the rows will fit snug as the row dimension is longer. If the code intends a 3 x 7 matrix, the columns will fit snug, with gaps in the rows, as displayed in the figures of our main code.

Step 14: The Main Code: Clearing All and Prompting for Inputs

To start the main code, click on the same button that says "New" in the upper left corner of the window, and select "Script" instead of "Function" from the earlier sections. Type the code as shown in the picture within the editor window. The lines of code do the following tasks in order:

  1. Close all files open except for 0, 1, and 2.
  2. Close all figure windows.
  3. Clear all variables in the workspace.
  4. Clear the Command Window.
  5. Display in the Command Window: Please input the following dimensions for the MRI scans:
  6. On a new line in the Command Window, ask: Length of one scan in milimeters: . The response put in by the user will be set to the variable "lengthMM".
  7. On a new line, ask: Distance between MRI scans in milimeters: . The response put in by the user will be set to the variable "ZStacks".

Step 15: The Main Code: Batch Processing the Images

In this section, the code will load up the images (consisting of the MRI scans of the cross sections of the brain) and store the names of each image file in the variable "Base" and display each of the MRI scans. Please follow with the code in the picture, which does the following:

  1. Create structure array "BrainImages" that contains information on all files within the current folder that fit the name format of MRI_().png
  2. Set variable "NumberofImages" equal to the number of elements in the structure array "BrainImages"
  3. Open figure window 1
  4. Set a for loop to cycle through for the number of images counted in the file
  5. For each loop, "CurrentImage" is the respective name of each file MRI_i.png, with the iteration number as 'i'
  6. Generate the a 3 x 7 subplot to display the 19 images to be loaded by "imshow"
  7. Display each image as another element in the subplot figure window
  8. Title each subplot element as Level_, where blank is the iteration number of the for loop.
  9. End the for loop (avoiding error)

This will display in figure window 1 all of the MRI scans in raw form in a 3 x 7 configuration with no gaps in the x orientation.

Step 16: The Main Code: Padding

With padding, we avoid the issue of slight discrepancies in the image sizes that may yield error for dimension mismatch in the case that one picture is slightly larger than another.

  1. Open figure window 2
  2. Load the image matrix from MRI_1.png into the variable "padBase"
  3. Find the size of the image's matrix and set to "OriginalXPixels" (for number of rows) and "OriginalYPixels" (for number of columns)
  4. Setup the matrix "BrainMat" to consist of all zeros with 20 more rows and 20 more columns for each plane, and 19 total cross sections, one per plane.
  5. Set up "HolesMat" to consist of the same three dimensional array of zeros to input hole coordinates later
  6. Create "zeroMat" to be the size of pad plus twenty rows and twenty columns, a two dimensional array of zeros.

Step 17: The Main Code: Determining Boundaries

  1. Set a for loop to go through the data from each image loaded earlier
  2. In the same manner of batch processing earlier, "CurrentImage" loads up files with "MRI_i.png", where i is the iteration number
  3. Run each image through the processing function "LevelCurveTracings2.m" that you made earlier
  4. Find the size of the output "Brain", set the number of rows to "Currentrow" and the number of columns to "Currentcolumns"
  5. Set "CurrentMat" to a matrix of zeros with the dimensions specified by "Currentrow" and "Currentcolumns"
  6. Center the data from "Brain" in "CurrentMat", with a margin of 10 rows on all sides
  7. Generate a subplot of dimensions 3 x 7, to display the boundaries of the images
  8. Title each of the subplot elements in the figure window
  9. Generate three dimensional matrix "BrainMat" comprised of each bounds layer "CurrentMat"
  10. End the for loop (for avoinding errors)

The subsection following fills in the holes left at the top and bottom of the proposed three dimensional shape

  1. Set "LevelCurve1" equal to the first layer of "BrainMat" (bottom of solid)
  2. Set "LevelCurveEnd" equal to the final layer of "BrainMat" (top of solid)
  3. Overwrite "LevelCurve1" with a filled in layer
  4. Overwrite "LevelCurveEnd" with a filled in layer
  5. Set the filled in layer as the bottom layer of "BrainMat"
  6. Set the filled in layer as the top layer of "BrainMat"

Step 18: The Main Code: Determining the Proper Z Dimension

The first three lines consist of setting up an empty array "z", and doing simple conversion operations (divide pixels by length) to get a proper reading of volume in mm^3

  1. Create a for loop to cycle through each layer
  2. Find the number of ones in a given layer
  3. Convert the z coordinates for the ones into values that are scaled to the proper ratio, set to "tempz", a column vector
  4. Add z value for the level curve to the vector z

With this z coordinates are adjusted properly.

Step 19: The Main Code: Determining X and Y Coordinates

Now to determine the x and y positions of each of the points in the boundaries.

  1. Initialize "xBrain" as an empty array
  2. Initialize "yBrain" as an empty array
  3. Set up a for loop to loop through each image loaded
  4. Compile a two column matrix to store the planar coordinates of each point in the bound, represented by the column vectors "RowBrain" and "ColumnBrain"
  5. Append "xBrain" with the currently found "RowBrain" coordinates
  6. Append "yBrain" with the currently found "ColumnBrain" coordinates
  7. End the for loop (for avoiding error)

Step 20: The Main Code: Plotting a Three Dimensional Structure, Finding Volume, and Logging Data

Using the function alphaShape, we will create a three dimensional structure from which we may calculate the volume of the brain.

  1. Use function alphaShape, plug in the vectors "xBrain", "yBrain" and "z" for the x, y, and z coordinates, and set equal to "BrainPolyhedron"
  2. Open figure window 3
  3. Plot the calculated alpha shape "BrainPolyhedron", display in the figure window
  4. Calculate the volume of the alpha shape, using a function "volume" that works for alpha shapes
  5. Convert the volume into mm^3
  6. Print the volume of the solid in the Command Window
  7. Prompt for a patient name to be set as an input
  8. Get current date and time with clock, and set to "DateandTime"
  9. Call function "patientFiles" to log and plot the data calculated

From here the the second and third pictures should show the figures that pop up, and the fourth picture that which should be displayed in the Command Window.