Plane Localization in 3-D Fetal Neurosonography for Longitudinal Analysis of the Developing Brain

The parasagittal (PS) plane is a 2-D diagnostic plane used routinely in cranial ultrasonography of the neonatal brain. This paper develops a novel approach to find the PS plane in a 3-D fetal ultrasound scan to allow image-based biomarkers to be tracked from prebirth through the first weeks of postbirth life. We propose an accurate plane-finding solution based on regression forests (RF). The method initially localizes the fetal brain and its midline automatically. The midline on several axial slices is used to detect the midsagittal plane, which is used as a constraint in the proposed RF framework to detect the PS plane. The proposed learning algorithm guides the RF learning method in a novel way by: 1) using informative voxels and voxel informative strength as a weighting within the training stage objective function, and 2) introducing regularization of the RF by proposing a geometrical feature within the training stage. Results on clinical data indicate that the new automated method is more reproducible than manual plane finding obtained by two clinicians.


I. INTRODUCTION
L ONGITUDINAL analysis of the developing brain is an emerging area linking medical image computing and neurology. In the future, this may have impact on clinical management of premature infants and neonates with neurological conditions. Clinicians are exploring tools to track the changes in image-based biomarkers from the womb to the postnatal period [1], [2]. However, promoting and monitoring neurodevelopment remains a significant challenge in the care of premature infants. Recent studies on cohorts of ex-premature infants have shown improved survival but long term neurodevelopmental outcomes have remained unchanged [3]. Volpe postulated that brain injury in premature infants is an amalgam of both destructive and developmental disturbances [4]. "Destructive" pathologies are currently diagnosed and monitored using cranial ultrasonography by identifying conditions such as intraventricular hemorrhages, infarctions, and hydrocephalus. However, developmental disturbances are more difficult to diagnose using cranial sonography. Current practice of assessment of brain ultrasound images for structural development is at best, subjective.
In a recent study [1], prematurity was shown to have a significant effect on thalamocortical development, specifically, lower gestation being associated with thalamic volume reduction. This volume depletion has been shown to be associated with lower cortical volume, decreased hippocampal volume and reduced fractional anisotropy of the corticospinal tracts [2]. This is therefore our clinical motivation for measuring thalamic size serially, which could potentially help identification of disordered brain development in premature infants. Serial MR imaging in premature infants is clinically impractical and expensive. Cranial ultrasonography is a more practical alternative to finding thalamic size with the benefit of being applicable cot-side. Furthermore, the thalamic shadow can be easily identified in the parasagittal plane on ultrasound images. Motivating this paper, a recent study [1] published preliminary findings suggesting a correlation between neonatal gestational age and the thalamic area measured manually on standard 2-D neonatal cranial ultrasound images acquired at the parasagittal plane. Thalamic area measurements were compared with equivalent fetal markers of a normal population as a reference, suggesting that this might provide a clinically useful way to monitor neuro-development of a neonate. Therefore, to derive normative data, PS planes were manually identified from 3-D ultrasound fetal brain volumes of healthy fetuses [5]. Associations with gestation and fetal thalamic cross-sectional area measurements were however weaker. To strengthen this association, automation of PS plane detection was envisaged. Furthermore, in fetal neurosonography, manual assessment of orthogonal planes is important to ensure a proper diagnosis of brain structures.
Recent advances in ultrasound techniques allow accurate measurement of anatomical structures within the fetal and preterm infant brain raising the prospect of identifying and tracking regions of interest correlating to maturational change at the cot-side. While standardized manual image acquisition has been adopted in neonatal cranial ultrasound practice [6], the ability to examine changes over time and correlate them with a reference healthy fetus at the same gestational age is crucially dependent on accurate and reproducible plane localization. Automation of this process would enable reproducible characterization of structural changes within the developing brain and may allow for assessing the effects of therapeutic intervention on brain growth and maturation.
Ultrasound-based longitudinal analysis across the perinatal period to allow correlation between clinical scans taken before and after birth remains a challenge. Antenatal scans (for the basic examination at midpregnancy) are acquired in utero with either a 2-D or a 3-D transducer positioned to obtain axial slices of the fetal brain. However, advanced fetal neurosonography examines the fetal brain in three planes (axial, sagittal, and coronal) to allow better understanding of brain development [7], [8]. Direct visualization of sagittal planes on fetuses is possible. However, it does not only require a lot of experience but also some luck to have the fetus in a perfect position to be able to acquire acceptable sagittal images. Therefore, finding a sagittal plane in a 3-D fetal image that is acquired on the axial view can provide a solution to this issue. Neonatal scans are taken ex-utero, with a different transducer normally limited to 2-D scans acquired through the anterior fontanelle (AF) to provide angled coronal and sagittal planes (cf., Fig. 1). The question we seek to answer in this paper is how fetal and neonatal scans are correlated. The approach outlined in this paper automatically finds, in a 3-D fetal scan, the equivalent (and nonconventional) 2-D plane, called the parasagittal plane (PS), used in neonatal neurosonography, via a machine-learning approach. Detecting such a plane reliably would allow neurological information to be obtained in theory from around 14 weeks gestational age up to at least eight weeks post birth or potentially six months depending on when the AF closes.
In fetal neurosonography, clinicians usually diagnose not only from the axial view of the brain but also the sagittal and coronal views. This assessment is done manually and more importantly, it typically does not assess oblique planes of the brain. In an interesting work by Rizzo et al. [9], the authors investigated the visualization of standard sagittal and coronal planes in 3-D ultrasound of the fetal brain. They studied the agreement rate of two manual reviewers in finding such planes. Therefore, to reduce subjectivity of finding planes manually and to allow oblique plane localization, it is important to develop automatic methods that can potentially detect planes in 3-D ultrasound.
There have been a small number of publications addressing automatic plane or view localization for clinical workflow improvement. On the other hand, [10] proposes a supervised learning method based on marginal space learning to detect several planes in 3-D echocardiography. A boosting framework [11] and regression forest (RF) [12] have been proposed to locate two, three, and four chamber views in 3-D echocardiography. A machine learning method based on the integrated detection network has been developed to detect several fetal brain planes in 3-D ultrasound volumes and estimate the standard clinical 2-D fetal biometric parameters [13]. In our previous work [14], we developed a semiautomatic method to localize the PS plane on 3-D ultrasound scans of the fetal brain. The method builds on a RF framework to estimate the PS plane parameters relying on only "informative" voxels during training and testing. One issue in [14] was the high intraobserver variability on manual localization of the PS plane. To address this issue, we have defined in this paper an objective protocol to find the PS plane.
Notice that, to-date, machine learning methods have been the most popular approach for plane finding in medical imaging, especially ultrasound. Regression forests [12], [14]- [18] has gained a lot of interest in medical imaging recently in applications such as object detection and view localization. RF consists of several trees trained on a set of training examples and possible features. Each split node in a tree optimizes an objective function to create a weak learner, which splits the training examples into two branches. Leaves contain the list of training examples, which traverse all the way from the root. Training examples at a specific leaf node are used to create a distribution from their continuous variable(s). The distribution at a leaf node on each tree is used during testing to vote for the continuous label for an unseen example. More details about RF can be found in [16].
In this paper, we propose a machine learning technique based on RF to provide accurate parasagittal plane localization. The technique constrains the RF method by: 1) forcing the model to use informative voxels and voxel informative strength as a weight when optimizing the objective function used in the training stage, and 2) introducing regularization of the RF by proposing a geometrical feature within the training stage. This new feature is based on the distance between each voxel and the midsagittal (MS) plane. We propose an automatic method to first localize the MS plane, which is then fed to the RF as a constraint feature. We demonstrate that our approach leads to plane detection accuracies, which are as good as manual plane finding and which show promise for their intended use of standardizing the parasagittal plane acquisition.
This study is an extension of our previous semiautomatic approach [14] with a number of important improvements outlined in the following. First, from a technical perspective, the new method is fully automatic whereas [14] is semiautomatic. Second, the evaluation is more thorough: We have doubled our dataset size, increased the range of fetal gestation, and evaluation is performed by two manual observers to investigate intra and interobserver variability. The reproducibility of manually selecting the parasagittal plane has never been studied before, as these planes are normally acquired as 2-D images in clinical practice. Therefore, the quantitative evaluation also includes the manual segmentation of the thalamus on a subset of the planes to investigate the sensitivity of the plane finder for thalamic area estimation.
The rest of the paper is organized as follows. In Section II, we present a framework to detect the midsagittal plane, used as a constraint to automatically find the parasagittal plane. Section III introduces the constrained RF for parasagittal plane finding. Section IV outlines the experimental setup used within this framework and Section V reports the quantitative evaluation. Discussion and conclusions are given in Sections VI and VII, respectively.

II. MIDSAGITTAL PLANE FINDING
One new contribution of the paper is to fully automate the PS plane localization by first localizing the MS plane automatically and using it as a constraint in the PS plane finding learning algorithm. The MS plane is clinically important to assess as it contains several structures which are typically assessed in neurosonographic examination (e.g., corpus callosum and cavum septum pellucidum). In addition, localizing the MS plane is useful itself as it facilitates the study of midline diseases, such as agenesis of corpus callosum or posterior fossa anomalies.
To find the MS plane, we first identify the skull in several 2-D slices of the volume to locate the brain region only (skull stripping). We then localize the midline on these slices since it is the separator between left and right hemispheres. Finally, the MS plane is then constructed as the best plane crossing the correctly detected midlines.

A. Acquisition Protocol of 3-D Fetal Ultrasound Images
According to the acquisition protocol used within this study, all fetal brain volumes were obtained approximately perpendicular to the left or right side of the fetal head. The acquisition aims at locating and measuring the biparietal diameter, found on the transthalamic plane in the middle of the volume. Fig. 1 shows a schematic of the acquisition protocol used, as well as the location of the midsagittal (MS), parasagittal (PS), and transthalamic planes.

B. Preprocessing
The process of skull localization starts by identifying several axial slices in the fetal head. These slices are those forming the middle region of the volume (above and below the middle slice, i.e., the transthalamic plane). We then process each slice independently to find the skull as described in Section II-C. We then check if the detected 2-D skull contours on the selected slices are consistent (i.e., long enough and with similar angles).  Prior to extracting the skull (and hence the midline), we derive the 3-D monogenic signal [19] representation for each volume. Using the monogenic signal is a powerful tool to deal with the appearance inhomogeneities inherent to ultrasound data, which make relying on intensities a hard problem [20]. The monogenic signal of an image I(x, y, z) is defined as where is the band-passedfiltered image by the Gaussian derivative filter b (x, y, z) [21] and ⊗ denotes the convolution operation. From the monogenic signal, local energy (or amplitude) A (x, y, z) and local phase ϕ (x, y, z) can be derived as and (3) The local energy provides a filtered version of the original image, resulting in a visible reduction of the speckle pattern and a clearer edge definition (i.e., skull). The local phase represents structural image information and one of its key properties is that it is invariant to contrast. Hence, local phase determines the type of a feature, whereas local energy determines that there is a feature at one particular pixel. A feature is defined for points of local phase congruency [22]. Fig. 2(c) shows the output of the amplitude on one 2-D slice and Fig. 2(d) shows the local phase corresponding to the same image. The local phase is used to enhance the visibility of the midline inside the brain as a long approximately straight line.

C. Skull and Brain Localization
To find the skull on a 2-D axial brain image, several skull templates were constructed at different sizes to cover a range of fetal gestational ages. Fig. 2(e) illustrates one of the templates used within the framework. We then compute the normalized cross-correlation (NCC) between the output of the amplitude [see Fig. 2(c)] of the 2-D image and each skull template. The template which has the maximum NCC response is used to localize the brain region by masking the 2-D original ultrasound image, amplitude image and local phase image by the bestmatched template of the brain. An extracted brain region is shown in Fig. 2(f) after skull stripping.
Note that the masking does not have to be exact, since we want to retain only the midregion of each slice containing the midline. In further processing, voxels outside the brain regions are neglected and the ones inside the skull are only used for MS and PS plane localization.

D. Midline Localization
Having found the brain region, the next task is to locate the midline. The midline or falx is a structure that appears in the middle of the brain separating the left and right hemispheres. It is typically seen as a bright line in the middle of each brain ultrasound axial slice as in Fig. 3. To find the MS plane, one needs to locate the midline on several slices and then find the plane matching the midline on these slices.
We used the Hough transform to locate the midline on each axial slice of an enhanced amplitude image, called A E , and defined as A E = A × ϕ, where A and ϕ are defined as in (3) and (4), respectively. We enhanced the amplitude image with local phase information because we found that due to the low contrast of ultrasound, speckle, shadows and attenuation of signal, the midline does not always have a strong signal but is locally ridge-like. This corresponds to local phase values of π/2 and 3π/2, for ridges going up or down, respectively. Multiplying the amplitude image by the local phase image enhances the structure of interest and the visualization of the midline is clearer.
The output of this step is the longest line found on every 2-D slice within the brain volume as illustrated in Fig. 3. However, the longest line on some slices does not necessarily represent the midline since some slices contain several structures that appear bright within the region of interest. This might lead to finding incorrect lines on some slices. To handle this, we select lines as follows. We first keep the longest 12 lines from all slices (this number was experimentally chosen). We then find the median angle amongst these lines. The lines whose angles have a small difference with respect to the median angle are retained and the others rejected. The output of this step is a line on a few brain slices corresponding to the brain midlines. An example of the detected midline is shown in Fig. 3(b) as a green line.

E. Midsagittal Plane Detection
Having found one line on each of the axial slices in a subset of the volume, we compute all points on every line as p_line(x, y, z). We then fit a plane to pass through all the line points on the retained axial slices. We use least-squares-fitting to compute the plane normal to these points. The output of this step is an MS plane normal and a location on the plane. This is used later within the RF technique to compute distance features from voxels to the MS plane in each volume.

A. Problem Description
Given a 3-D fetal ultrasound scan, the midsagittal plane is detected as described in Section II. The goal now is to use the fetal brain volume and the midsagittal plane to locate the parasagittal plane. We propose a machine learning method to solve this problem. The method is based on the RF technique. In the training stage, a set of 3-D training images with their corresponding manual planes are used to build the forest. In our formulation, we can generate a training example for each voxel. A training example for voxel v i has the following form: , d y , d z ) .
Here, f j is a feature computed for the voxel v i . n is the total number of possible features. d = (dx, dy, dz) is a distance vector which represents the continuous class label for the training example. d is the perpendicular distance between v i and the manual PS plane. The PS plane is represented as three parameters for the normal vector and three parameters to represent a point on the plane. We compute d as in [12].
The list of training examples from several ultrasound volumes alongside their continuous class labels are used to train an RF model. The forest is then used to test unseen ultrasound volumes.

B. Sampling Informative Voxels
During the RF training on image volumes, millions of features need to be computed on typically millions of examples (voxels in our case). RF is responsible for finding the best set of features to build each tree. Although many of the training examples are noninformative, RF treats all training examples equally during training. During a node split in the forest, an objective function is optimized over all training examples which reached that node. If several of these examples are noninformative, the objective function could fail to find a good split. As a result, this could create weaker trees and hence a weaker forest which means lower regression accuracy. Additionally during testing of an unseen volume, relying on all testing examples (voxels) to vote for plane parameters could potentially yield erroneous results because many voxels are noninformative.
We only sample voxels from inside the brain since we extract the brain region (see Section II). This allows us to exclude all voxels from outside the brain that do not necessarily contain meaningful information. In fact, due to the existence of maternal tissues and the appearance heterogeneity outside the fetal brain region, voxels from outside the brain are misleading to the learning algorithm. Therefore, we restrict the learning algorithm to utilize brain voxels only to locate the PS plane.
We also propose a technique which finds the informative voxels within the fetal brain and uses their strength to weigh their influence during tree training and testing. We present a solution to overcome a common problem observed in machine learning in ultrasound of otherwise generating a suboptimal RF because of the existence of many noninformative voxels.
An informative voxel is a voxel that is discriminative within its neighborhood. Here, we use feature asymmetry (FA) to highlight informative voxels motivated by its successful application to 3-D echo feature enhancement and detection, such as [23], [24]. FA is computed as where , . sets to zero the negative values, s represents the scale, N is the total number of scales considered, ε is a constant used to avoid division by zero (typically ε = 0.01), and T s is an orientation independent threshold that controls the spurious responses to noise at scale s [22]. FA computes the step edge strength at a specific voxel and provides a value in the range [0, 1], thus providing a direct measure of the strength at a specific voxel. The region-of-interest (ROI) highlighted in red in Fig. 4 shows an area with speckle texture while the ROI highlighted in green shows an area containing fetal brain structures (part of the falx). Most importantly, even within areas of apparent weak image appearance (red ROI in Fig. 4); FA identifies the voxels which Fig. 4. Localization of informative voxels. Blue ellipse represents the brain region which we sample from. Red is an ROI with less informative voxels while green is an ROI with highly informative voxels. Ultrasound slice of a fetal brain. Corresponding FA image. Corresponding magnified regions from the 2-D slices. may be discriminative. Fig. 4 motivates the choice of using FA measure to represent voxels informative strength.
We use FA to create stronger trees during training by: 1) masking out all voxels which have FA measure below a specific threshold and 2) using the FA value to weight the importance of each training example when finding a splitting score in a tree node.
The traditional function that is optimized is where z is the variable that defines the optimal split of a training set. Instead, we normalize the traditional variable by the FA measure, i.e., z = z FA . Note that we divide rather than multiply since we are minimizing an objective function. We find the best split of a training set S (v is a training example) to left (L) and right (R), denoted (S L and S R ).z is the mean value of z for all training points reaching a specific tree node.

C. Feature Sets
Careful selection of features is an important part in any machine learning algorithm. In the medical image analysis literature, several types of features have been used [11], [17], [18], [25]- [29]. In this study, we use traditional appearance features which have shown promising results in other medical imaging domains [17], [29]. However, because of the complexity of our problem, we propose a constrained distance feature set to aid appearance features when localizing the PS plane.
1) Appearance Features: Appearance features are computed within a neighborhood of each voxel. We deploy features like voxel intensity, mean intensity within a cuboid and difference of mean intensity of two cuboids. These types of features have been successfully used in related classification and RF applications and showed promising results [17], [29].
2) Constrained Distance Features: The use of appearance features within RF creates a learner which works well if image appearance is clear and not ambiguous. However, in ultrasound imaging this is not always the case as previously discussed. To address this issue we introduce a geometric constraint expressed as a distance feature to help constrain the RF training. This type of feature embeds the perpendicular distance between the voxel of interest and the midsagittal plane. Given the midsagittal plane normal n = (n x , n y , n z ), a point on the plane p = (p x ,p y ,p z ) and a point in space v = (v x , v y , v z ), the perpendicular distance D between v and the plane parameterized by n and p can be expressed as This feature clusters training examples depending on their distance from the midsagittal plane. However, this distance feature alone is not sufficient to capture the variability of plane appearance in different fetuses especially when dealing with different gestational ages. Therefore, it is important to incorporate both appearance and distance features inside the RF and to let the optimization of the objective function choose a feature at a specific depth depending on its ability to split the training points properly.

D. Localizing the Parasagittal Plane
In the testing stage, the informative voxels in a test image traverse the trees until each voxel ends up in a leaf. Only voxels that reach confident leaves are allowed to vote for the output plane parameters. The confidence level is measured using the variance of the training examples that reached a specific leaf during training. The mean parameters of all informative voxels reaching confident leaves are used to output the final plane parameters by fitting a mean plane which represents the estimated parasagittal plane position.

A. Dataset and Manual Annotation
161 3-D fetal brain ultrasounds between 23 and 33 weeks of gestation from normal fetuses were acquired with a Philips HD9 ultrasound machine [5]. These fetal brain volumes were not acquired specifically for our study but for standard fetal biometry, and hence were not optimized to our application. The number of volumes used in each gestational week is presented in Fig. 5. Volumes vary in size but they are approximately (250 × 250 × 200) voxels. There were no specific settings used to acquire the 3-D images because these images were acquired for a clinical study [5] and before the start of our analysis study. We did not use volumes for fetuses older than 33 weeks because this study was set up to correlate to neonatal scans of preterm babies, and the age range of interest was until 33 weeks of gestation, as these are defined to be more critical than later on.
The quantitative evaluation reported in this paper is more extensive than in our earlier work [14]. Two clinicians participated in the evaluation by manually selecting the PS and MS planes. Clinician 1 manually localized the PS plane twice to derive intraobserver variability. Clinician 2 manually localized the PS plane once to investigate interobserver variability. Our model was trained by the second manual localization of the PS planes performed by the clinician 1. A strict protocol to find the PS plane was defined prior to manual selection of planes. Manual localization of the MS plane was carried out once to allow us to investigate the accuracy of detecting this plane automatically.

B. Evaluation Measures
To evaluate the accuracy of plane-finding, two evaluation metrics were used. The first is the angular distance between the manual and the detected planes defined as where n m is the normal to the manually defined plane and n d is the normal to the automatically found plane. The angular distance provides insight into how close two planes are in the 3-D space but does not take into consideration their line of intersection -two planes may have a small angle between them but their intersection can be far outside the brain which means the error in diagnostic plane finding is large. To circumvent this problem, we include, as a second evaluation metric, the Euclidian distance between the centers of the manual and the automatic PS planes. Since the fetal brain appears approximately in the middle of the volume in most cases (given this requirement is part of the scanning protocol), the distance metric expresses how far the center of the two planes are from each other. The smaller the distance and angle, the higher the accuracy. In addition, we measured the angle between the manual and automatic MS planes as well as the distance between the center of the MS plane to the center of the volume.
To evaluate the PS plane localization, we report six comparisons: 1) manual-to-manual plane agreement, which includes intra and interexpert plane comparison; 2) manual-to-automatic plane agreement with the RF trained on appearance features only; 3) manual-to-automatic plane agreement using the RF trained on appearance features only but constrained by the use of informative voxels (CRF-FA); 4) manual-to-automatic plane agreement using RF trained on all examples but constrained by automatic MS plane distance features (CRF-Dist); 5) manual-to-automatic plane agreement using RF constrained by the use of informative voxels and automatic MS plane distance features (CRF-FA-Dist-A); and finally, 6) manual-tosemiautomatic plane agreement using RF constrained by the use of informative voxels and manual MS plane distance features (CRF-FA-Dist-M). Fig. 6. Distribution of distances for the manual MS planes (x-axis) and the automatic MS planes. The distance is measured as the distance between the center of the volume to the center of the plane.

C. Implementation Details
We performed three-fold cross-validation to ensure that every volume is included in testing. In all experiments we optimized the method parameters on the classic RF and fixed them for the proposed RFs to allow a fair comparison. We used 15 trees with a maximum depth of 20. We optimized 300 candidate features during node creation. We set the FA threshold to 0.2 to allow, as much as possible, informative voxels during the training stage. In addition, when comparing multiple FA thresholds, visual inspection showed that the chosen threshold provides a suitable spread of voxels in different images. This implies informative voxels from inside the brain were selected during the training and testing stages. We used the 3-D Gaussian derivative [21] of scaleσ, as a filter to compute FA [see (4)] and we experimentally set the scales to σ = 4 Finally, we used the C# language to build the methods in a parallel fashion.

A. MS Plane Localization
We compared the manual MS planes to the automatically generated ones. The mean angle between the two was found to be 5.6 ± 2.6°. Although the angle between two planes provides an insight on how close they are to each other, it does not provide a physical distance between them. Therefore, we measured the distance between the center of the volume and the center of the manual and automatic MS planes. Fig. 6 shows the correlation between the manual and automatic distances. The mean of the absolute difference between the manual and automatic distances is 2.1 ± 1.9 mm. The correlation coefficient between the distances is 0.83 (p-value < 0.005).

B. Manual to Manual PS Plane Localization
The first two rows in Table I show the mean angle and distance between the manual planes localized by the same clinician (M-M intra) and by two clinicians (M-M inter), respectively. M-M intra reflects the intraobserver variability in detecting the PS plane while M-M inter evaluates the interobserver variability. Table I presents the mean and standard deviation of angles and distances between different PS plane localization methods.  Fig. 7(a) shows the distribution of angles between manual and automatic methods at all gestational ages we have. We also present in Fig. 7(b) the distribution of distances for the proposed methods. Finally, Fig. 8 shows the PS planes localized manually and by the proposed methods for visual comparison.

A. Importance of the Preprocessing
Since the quality of the ultrasound images suffers from several issues like shadowing, weak edges, speckle, we have employed the monogenic signal to create suitable features for analysis. This is not a major new contribution to this paper as this form of preprocessing has been done before. However, it is an important step which contributes to the success of our method.

B. MS Plane Localization
We presented an automatic method to detect the MS plane in 3-D ultrasound images of the fetal brain. We first locate the skull and extract the brain region. We then localize the midline at several slices. Midlines at several brain slices are then used to fit the best plane which represents the estimated MS plane. Evaluation results confirm that we can detect the MS plane with high accuracy.

C. Automatic PS Plane Localization
As shown in Table I and Fig. 7, the classic RF struggles to perform well especially at later gestations. The use of informa-tive voxels and voxels within the brain (CRF-FA) improves the PS plane localization results compared to the conventional RF. Adding the newly proposed MS plane distance feature (CRF-Dist) also improves localization over the classic RF. However, both CRF-FA and CRF-Dist are still less accurate than the interobserver manual planes agreement. On the other hand, combining both ideas (CRF-FA-Dist-A) not only improves the localization accuracy but also provides better localization than the interobserver variability obtained by two clinicians.
The best proposed method is CRF-FA-Dist-M; however, it is semiautomatic. In addition, it is marginally better than its automatic equivalent (CRF-FA-Dist-A). Therefore, taking into account its accuracy and no need for user interaction, CRF-FA-Dist-A is the most promising method of all.
Note that in Fig. 7, the accuracy of the localization of the PS plane decreases with increased gestational age. This is mainly because fetal brains of later gestation are bigger and therefore the thalamus and its surrounding structures are larger in size and hence present more room for localization errors. In fact, if we compare our previous work in [14] with our current study, we find that the localization accuracy in [14] using CRF-FA-Dist-M is better than here since that paper used a smaller gestational age range (23-27 weeks) and smaller number of fetal brain volumes (87 volumes).

VII. CONCLUSION AND FUTURE WORK
This paper introduces for the first time an automatic method to localize the PS plane in 3-D ultrasound volumes of the fetal brain. This is a very difficult task for a sonographer to do well. Results show that the automatic method provides more precise and reproducible PS planes than the interexpert variability between two clinicians.
One way to possibly improve the MS plane localization would be to detect the midline on axial and coronal slices and then fit the best plane between these. This could provide a more accurate orientation of the MS plane since it is detected using different views. The improvement on the MS plane detection could also improve the PS plane detection as suggested by Table I where best results were obtained with manual MS plane detection in the main algorithm.
Although we showed that finding the PS plane is important to facilitate the measurement of the thalamus, the PS plane contains several other structures which would be interesting to assess in the future, such as the lateral ventricle wall and choroid plexus. In addition, we applied our method to volumes from one ultrasound machine, but the method is generalizable to images from other ultrasound machines and potentially other imaging modalities (e.g., MRI). Finally, to validate the utility as a tool for biomarker assessment, a thalamic area comparison on manual and automatic planes needs to be investigated. This is the subject of future clinical work.