In-vivo magnetic resonance imaging (MRI) of laminae in the human cortex.

The human neocortex is organized radially into six layers which differ in their myelination and the density and arrangement of neuronal cells. This cortical cyto- and myeloarchitecture plays a central role in the anatomical and functional neuroanatomy but is primarily accessible through invasive histology only. To overcome this limitation, several non-invasive MRI approaches have been, and are being, developed to resolve the anatomical cortical layers. As a result, recent studies on large populations and structure-function relationships at the laminar level became possible. Early proof-of-concept studies targeted conspicuous laminar structures such as the stria of Gennari in the primary visual cortex. Recent work characterized the laminar structure outside the visual cortex, investigated the relationship between laminar structure and function, and demonstrated layer-specific maturation effects. This paper reviews the methods and in-vivo MRI studies on the anatomical layers in the human cortex based on conventional and quantitative MRI (excluding diffusion imaging). A focus is on the related challenges, promises and potential future developments. The rapid development of MRI scanners, motion correction techniques, analysis methods and biophysical modeling promise to overcome the challenges of spatial resolution, precision and specificity of systematic imaging of cortical laminae.


Need for non-invasive mapping of cortical lamination
The neocortex of the human brain consists of six layers which differ in their pattern of myelination (Nieuwenhuys, 2013;Vogt and Vogt, 1919) and the density and arrangement of neuronal cells (Brodmann, 1909;Zilles and Amunts, 2010). Those differences gave rise to the broad research field of myelo-and cyto-architecture, which is based on invasive histology. To facilitate studies of the human brain in large populations and longitudinally, non-invasive methods are much needed for visualizing intracortical anatomy in-vivo, since histological methods can only be applied ex-vivo. A promising approach is magnetic resonance imaging (MRI) as it is non-invasive and highly sensitive to the presence of myelin (and iron) in the cortex (Stüber et al., 2014). It is expected to reflect the myeloarchitecture rather than the cyto-architecture, although there is a strong correspondence between the two .

Scope of the review
This paper reviews the recent in-vivo MRI studies on the anatomical layers in the human neocortex. An additional focus is on the related challenges, promises and potential future developments. Studies using diffusion weighted MRI are excluded from this review, since they are covered by Assaf et al. (2017) in the same special issue on "Prospects for cortical laminar MRI: functional and anatomical imaging of cerebral cortical layers".

Studies of the stria of Gennari
As the cortex is only 2-4 mm thick and convoluted (Fischl and Dale, 2000), mapping cortical layers requires sub-millimeter isotropic resolution. Another challenge is that the differences in myelo-architecture between cortical layers and between different brain areas are rather subtle (Braitenberg, 1962) and significantly smaller than between cortex and white matter or other brain structures. Unsurprisingly, the first target of intracortical MRI was the hallmark of the primary visual cortex V1 or Brodmann area (BA) 17 (Clark et al., 1992) the strongly myelinated stria of Gennari (Gennari 1782, Vicq d'Azyr 1786. V1 was also among the first areas for which a structure-function relationship was successfully established (Bolton, 1900), inspiring MRI studies combining structural and functional MRI (Hinds et al., 2009). Although the stria of Gennari makes up only around 0.3 mm of V1's 2 mm thickness (von Economo and Koskinas, 1929), it is even visible to the naked eye in untreated, freshly cut cadaver brain sections as a white strip running through the cortical gray matter (Funkhouser, 1915;Gennari, 1782, Vicq d'Azyr 1786. Consequently, in-vivo MRI studies successfully visualizing this layer were performed at 1.5 T (Clark et al., 1992;Walters et al., 2003) and 3 T (Barbier et al., 2002;Bridge et al., 2005;Clare and Bridge, 2005;Turner et al., 2008). Different MR sequences were employed in the studies. Whereas Clark et al. (1992) used proton density (PD) weighting for visualizing the Gennari stripe, in the other studies a preparation pulse was used to establish a high T1-weighted contrast, taking advantage of the high myelination of the stria of Gennari and the resulting decreased T1 relaxation time Koenig et al. 1990). The contrast between the Gennari stripe and the surrounding layers of gray matter can be further increased by using a Turbo Spin Echo (TSE) readout after the inversion pulse, to add magnetization transfer (MT) contrast to the T1 weighting (Turner et al., 2008).
In studies at 1.5 T and 3 T, the required Signal-to-Noise Ratio (SNR) could only be achieved by averaging, resulting in long scan times of 40 min and more. The increased field strength of 4.7 T and sophisticated radio-frequency (RF) hardware increased the SNR substantially and thus reduced the scan time for visualizing the stria of Gennari to only 6 min for 17 slices covering the occipital lobe (Carmichael et al., 2006). However, even with such a setup the use of 2 mm thick slices, requiring a perfect perpendicular alignment of the imaging slice to the cortex minimizing partial volume effects, was necessary to achieve sufficient SNR. Even higher fields of 7 T enable the robust visualization of the stria of Gennari in about 10 min using isotropic voxels of 0.5 mm resolution (Trampel et al., 2011). Going beyond proof-of-concept studies, the high resolution and short scan times now allow for anatomical examination such as the investigation of the Gennari stripe in congenitally blind people (Trampel et al., 2011), which was found to form and persist independently of visual sensation. The studies performed at 4.7 and 7 T employed TSE or GRASE  acquisition techniques, and therefore yielded a mixed contrast of MT, PD and T2 relaxation weighting resulting in a strong contrast between the myelin rich Gennari stripe and the adjacent layers. Interestingly, identifying the stria of Gennari by a combination of magnetization transfer contrast and nuclear Overhauser enhancement (NOE) suggested that the NOE signal is more sensitive to myelination than MT (Mougin et al., 2013).
A high isotropic resolution of 0.4 mm, obtained in a 15-min FLASH acquisition, was achieved by optimizing for the shorter T2* relaxation time of the stria of Gennari compared to the surrounding layers of gray matter (S anchez-Panchuelo et al., 2012). However, exceptional high intracortical differentiation can be achieved by analyzing not only the magnitude but also the phase of the MRI signal acquired by a T2*-weighted gradient echo sequence (Deistung et al., 2013;Duyn et al., 2007). These images are not only sensitive to the increased myelination in the Gennari stripe but also its higher iron content, since both contribute to its susceptibility difference to the surrounding gray matter (Duyn et al., 2007). The influence of iron on the MR signal was elegantly investigated by removing iron from post mortem samples of human brain tissue (Fukunaga et al., 2010). MR images of the stria of Gennari obtained at different field strengths using different MR sequences and contrast weightings are depicted in Fig. 1.
The emergence of functional MRI (fMRI) in 1992 greatly facilitated the systematic study of human functional neuroanatomy and structurefunction relationships. The use of ultra-high field strengths of 7 T and above and the resulting high SNR help linking of structure and function even on a single subject basis (S anchez-Panchuelo et al., 2012;Turner and Geyer, 2014). For the reasons mentioned above it is not surprising that the primary visual cortex is an area where a correlation between cortical lamination and executed function was studied in humans at 3 T Clare and Bridge, 2005;Koopmans et al., 2010) and at 7 T (S anchez-Panchuelo et al., 2012) in-vivo. By using the stria of Gennari as a landmark, it was possible to determine the border between primary visual cortex V1 and higher visual area V2 anatomically S anchez-Panchuelo et al., 2012). Whereas, Bridge at al.   Barbier et al. (2002) with permission from the publishers. C) Field strength: 7 T, 24-channel phased array coil, Main contrast: MT þ T2 þ PD, Voxel size: 0.5 mm Â 0.5 mm x 0.5 mm, Acquisition time:~10 min. Taken from Trampel et al. (2011) with permission from the publishers. D) Field strength: 7 T, 32-channel phased array coil, Main contrast: Signal phase, Voxel size: 0.24 mm Â 0.24 mm x 1 mm, Acquisition time: 6.5-13 min (area outlined by yellow box is used for later analysis in the respective study). Taken from Duyn et al. (2007) with permission from the publishers, Copyright (2007) National Academy of Sciences.

Studies of lamination outside V1
The stria of Gennari is a particularly prominent variant of the outer band of Baillarger, (1840), which is a myelin sheet that can be found almost throughout the entire isocortex at a microscopic level. Using sufficiently sensitive acquisition and analysis techniques, the layered structure, particularly the inner band of Baillarger, can be observed in MRI data outside V1. Successful visualization of intracortical anatomy was achieved in motion sensitive visual area V5 (S anchez-Panchuelo et al., 2012;Walters et al., 2003Walters et al., , 2007 and auditory cortex (De Martino et al., 2015;Dick et al., 2012). The MR contrasts used for structurally identifying cortical layers with different myelination included T1 at 1.5 T (Walters et al., 2003(Walters et al., , 2007, T2* at 7 T (S anchez-Panchuelo et al., 2012), and a combination of T1, T2* and PD at 7 T (De Martino et al., 2015). In contrast to the conventional "weighting" towards certain MR parameters, Dick et al. used quantitative maps of the relaxation rate R1 (¼ 1/T1) to obtain intracortical profiles of the auditory core (Dick et al., 2012), which improves comparability across brain areas.
In particular, the bands of Baillarger were identified also in the parietal and frontal lobes (Fracasso et al., 2016). In this study MRI data were compared to histology of post mortem specimens and hypo-intense bands in gray matter, consistent with the bands of Baillarger, were observed. Here again, the sensitivity of T1 relaxation to myelin content was leveraged using a bias field corrected T1-weighted MPRAGE sequence. Barazany and Assaf (2012) and Lifshits et al. (2017) used multi inversion recovery to map T1 distributions, in order to reveal sub-voxel layer components, which show different T1 values, and achieve a super-resolution effect. Using a FLAIR sequence optimized for ultra-high field strengths, a multi-layered appearance throughout the cerebral cortex with different regional profiles could be visualized at 7 T (Zwanenburg et al., 2012).
Being able to detect the myelo-architecture and cortical profiles invivo throughout the entire brain offers new opportunities for parcellating the human cortex (Dick et al., 2012;Geyer et al., 2011;Marques et al., 2017) dubbed "in-vivo Brodmann mapping" , a particular aspect of "in-vivo histology" using MRI (hMRI, Clare, 2006, Deistung et al., 2013;Weiskopf et al., 2015). Cortical profiles of quantitative values, such as T1, T2* or susceptibility, may thereby help to detect variations of myelin and iron on an intracortical level. Using area-specific models, estimating cyto-architecture in addition to myelo-architecture may be feasible . The potential advantage of individual parcellation of the cortex over standardized atlases is the more accurate definition of brain areas (Lutti et al., 2014;Tardif et al., 2015;Waehnert et al., 2016). The approach can also help to indirectly determine specific functional areas, making lengthy and complicated functional localizer scans unnecessary, e.g., retinotopic mapping of V1 being complemented or replaced by mapping the striate cortex S anchez-Panchuelo et al., 2012;Sereno et al., 2013).

Challenges and recent developments
The reliable imaging of cortical layers poses substantial challenges for the entire imaging chain, i.e., data acquisition, image reconstruction, image processing and biophysical modeling. It is still an open question which MR contrast is best suited for visualizing the layered-structure of the cortex. Candidates include approaches based on the magnitude and phase of gradient echo acquisitions (Deistung et al., 2013;Duyn et al., 2007), and turbo spin echo sequences (Trampel et al., 2011). Exceptional intracortical contrast also arises from those TSE sequences which are additionally equipped with a T1-preparation (Turner et al., 2008) or adjusted to achieve a FLAIR-like weighting (Zwanenburg et al., 2012). Recently, a study showed an improved contrast between the Gennari stripe and the adjacent cortex layers for TSE acquisitions compared to T1-weighted images and T1 maps obtained with the MP2RAGE ( However, for comparing the layering in different cortical areas or across subjects quantitative techniques are advantageous over conventional "weighted" imaging (see section 4.2).
Obviously, imaging cortical layers requires imaging at the maximal spatial resolution possible, since the layers are only 200-300 μm thick. Since in-vivo MRI is primarily SNR and scan time limited, it requires methods maximizing the SNR per time unit. Therefore, techniques for increasing SNR will help future studies to gain even more insight into the human cerebral cortex.
Conceptually, layer-specific imaging requires handling a broad spatial scale ranging from the microstructure over the mesoscopic to the macroscopic scale. This requires novel approaches to image processing and biophysical modeling. Characterizing the layer-specific myelo-architecture requires information about the intracortical myelin fibers whose size is in the range of a few micrometers, i.e., information about the microstructure. The definition of the layers is based on delineating these 200-300 μm thick structures, i.e., information about the mesoscopic structure. In order to account for the extent of typical brain areas a field of view of at least a few centimeters is optimal, i.e., a macroscopic view. We note that the range of spatial scales introduces a transition from discrete structures such as cell bodies to continuous metrics such as average MRI signal intensities or parameter values within an imaging voxel.

MRI acquisition methods and hardware
A straightforward approach to increase the SNR is the averaging of several acquisitions and/or extending scan times. For long scan times this approach is only successful if data acquisition is corrected for subject motion. A promising motion correction technique providing the necessary precision in the range of a 1/10th voxel size (i.e., approximately 50 μm) and speed (approximate sampling rate 100 Hz) is prospective motion correction based on optical tracking (see Fig. 2; Callaghan et al., 2015;Maclaren et al., 2012;Schulz et al., 2012). However, sophisticated hardware, like MR compatible video cameras, are necessary. Techniques based on MR data acquisition modifications, e.g., using the fat signal as a navigator , do not need additional hardware but require MR sequence modifications which prolong acquisition times for sequences without sufficiently long "dead times" or require additional RF pulses, increasing the specific absorption rate which may be especially disadvantageous for high field imaging. However, the improvement in data quality promises the robust delineation of intracortical structures (see Fig. 3; Federau and Gallichan, 2016). Higher resolutions can be attained by averaging together the data acquired across multiple sessions (Lüsebrink et al., 2017;Pine et al., 2017) but precise non-linear image registration is required to account for short-term changes in brain morphometry. The SNR performance of MR pulse sequences can also be improved by increasing the data acquisition duty cycle using multi-echo readouts (Helms and Dechent, 2009), tailoring the readout bandwidth to the application, or using different steady state magnetization readouts, e.g., capitalizing on the high steady-state magnetization of SSFP acquisitions (Scheffler and Lehnhardt, 2003).
The SNR increases supralinearly with the static magnetic field strength, making field strengths of 7 T and beyond attractive for high resolution imaging (Pohmann et al., 2016). Susceptibility-related contrast also increases at higher fields, further improving the contrast-to-noise ratio for the related contrasts. Consequently, T2*-, phase-, and susceptibility-weighted contrast obtained at 9.4 T already proved that the higher spatial resolution improves the visualization of intracortical details (Budde et al., 2011(Budde et al., , 2014. The recent and upcoming installations of human MRI scanners with a field strength of 10.5 T or more promise further improvements in SNR and spatial resolution (Eryaman et al., 2017). However, it should be noted that the increased specific absorption rate and pronounced transmit field inhomogeneities are still an issue for ultra-high field imaging (Collins et al., 2004;Van der Moortele et al., 2005;van der Zwaag et al., 2016).
The use of high performance RF coils can also significantly improve the SNR, particularly in the cortex. The advent of multi-channel RF receive coils with 32 or more channels has led to improvements in SNR by more than a factor 2-3 in the cortex compared to single channel birdcage RF coils (Keil and Wald, 2013), while still maintaining a whole head coverage and ease of use. At the same time, the increase in channel count improves the parallel imaging performance used to shorten the lengthy ultra high resolution scans (Keil and Wald, 2013). Unsurprisingly, the multi-channel coils are now routinely used and became a mainstay for high resolution imaging. Despite this success special RF coil developments are ongoing to increase the local SNR even further. An example is the construction of flexible RF surface coils arrays that facilitate the optimal close fit of the coil to the head and thus maximize the SNR (Corea et al., 2016;Kriegl et al., 2015).

Quantitative MRI and in-vivo histology using MRI (hMRI)
The specificity of the cortical MRI microstructure estimates plays an equally important role as the spatial resolution for mapping of cortical layers. The majority of studies on intracortical myelination relied on MRI whose image contrast is weighted towards multiple relaxation time constants such as T1, T2* or similar. This conventional approach has the disadvantage that the image intensity values are not quantitative and the images yield a mixed contrast sensitive to different underlying microstructural components. For example, an increased iron concentration and myelin density may appear the same in T1-weighted images and cannot be well distinguished.
Quantitative anatomical MRI yields parameter maps of relaxation times, susceptibility or (MT) saturation effects in standardized physical units that are comparable across time points and sites (Deistung et al., 2013;Marques et al., 2017;Weiskopf et al., 2013Weiskopf et al., , 2015. This can significantly reduce the measurement variability in multi-center studies . Frequently, they provide information about multiple parameters at the same time, which offer different views on the underlying microstructure . Quantitative MRI can determine the different MR parameters (e.g., relaxation times) underlying the image contrast. This results in a more direct relationship with specific microstructural components, although not a trivial one-to-one relationship. Studies on the contribution of myelin and iron to different relaxation time constants showed that maps of R1 (¼1/T1) values are mainly dominated by myelin concentration whereas R2* (¼1/T2*) or quantitative susceptibility maps are particularly sensitive to iron (Wharton and Bowtell, 2012) in deep gray matter structures but sensitive to both myelin and iron in cortical gray matter (Stüber et al., 2014;Wallace et al., 2016). Proton density maps, however, can be used to visualize macromolecular tissue volume (Mezer et al., 2013). The contrast in weighted images, however, is generated by different MR parameters (e.g., the contrast of a T1-weighted image is a mixture of proton density, longitudinal and apparent transverse relaxation) hampering the inference upon the underlying microstructure. An alternative to quantitative imaging used for parcellation of the cortex is the use of ratios between T1-and T2-weighted images producing maps reflecting myeloarchitecture. However, it has not been extended to layer-specific studies (Glasser and Van Essen, 2011). An example of T1 maps obtained at different cortical depths across the whole brain and acquired using the bias-field corrected MP2RAGE sequence (Marques et al., 2010) is given in Fig. 4.
The more biologically specific information provided by quantitative MRI was used to access the intracortical myelination. In particular, the Fig. 3. Quantitative maps of PD and R1 with a 400 μm isotropic voxel size from a healthy, 28 year old female volunteer (Pine et al., 2017). Whole brain maps were created from two multi-echo 3D gradient echo volumes with PD-and T1-weightings and calibration data to correct RF transmit field non-uniformity (field strength: 7 T, total acquisition time: 70 min), using bespoke MATLAB tools within SPM12 (http://www.fil.ion.ucl.ac.uk/spm). The high resolution allows excellent delineation of intracortical features, such as the stria of Gennari (insets; hypointense in PD and hyperintense in R1 maps). longitudinal relaxation time/rate (T1 ¼ 1/R1) demonstrated clear layerspecific patterns (see Fig. 5B; Dick et al., 2012;Lutti et al., 2014;Sereno et al., 2013) such as the line of Gennari or brain area specific profiles. Also MT saturation parameters were used to track intracortical maturation and myelination in adolescence (see Fig. 5C; Whitaker et al., 2016). Cortical profiles obtained at different brain areas using different MR contrasts are shown in Fig. 5.
Different research groups work on extending MRI and quantitative MRI into a tool for directly mapping the brain microstructure, e.g., biologically relevant parameters such as fiber densities, axon diameters, relative axonal myelination (g-ratio) or myelin concentration (for a review, see Weiskopf et al., 2015). Such a tool may provide similar information as invasive histological approaches but allow for in-vivo investigations due to its non-invasive nature, i.e., enable in-vivo histology using MRI (hMRI). Most approaches focused on estimating white matter microstructure so far. For example, methods based on diffusion weighted images were developed for estimating axonal diameters and densities in white matter (e.g., CHARMED, ActiveAx, AxCaliber; Assaf et al., 2013). Various contrasts were explored and also calibrated for measuring the myelin volume fraction (e.g., T1, T2*, MT and T2 based metrics; Stüber et al., 2014) or the related macromolecular tissue volume (Mezer et al., 2013). A recent study on myeloarchitecture of the human cortex at different cortical depths showed, however, that a extrapolation from ex-vivo results obtained for a single specimen (Stüber et al., 2014) to in-vivo whole brain data cannot easily be performed (Marques et al., 2017). Whereas a linear relationship between R2*/susceptibility and myelin volume could be confirmed throughout the cortex the results were less clear for R1 which is thought to be modulated by the existing myelin surface fraction rather than its volume (Marques et al., 2017). This suggests that future work is still needed in finessing the models.
Combining the information gleaned from different contrasts widens the possibilities for in-vivo histology further. Recently, myelin volume fraction measurements based on MT were combined with fiber density estimates from diffusion imaging to provide estimates of the relative axonal myelination Stikov et al., 2015). Despite the success and more or less widespread use of these approaches, the careful and systematic validation against gold standard data remains a major challenge. For example, comparisons against ex-vivo histology methods are complicated by the vagaries of histological staining, limited field of views and the difficulty to create comprehensive 3D sets of histological information. This requires the development and application of novel ex-vivo histology approaches such as quantitative PIXE microscopy (Stüber et al., 2014) or 3D microscopy enabled by clearing techniques .

Defining the laminar depth coordinate system
To combine the information across the spatial scales ranging from microstructure to gross anatomy, a robust definition of cortical layers, i.e., a laminar depth coordinate system, is essential. Given anatomical images of adequate resolution, automated or manual segmentations of the cortical boundaries provide a natural surface-based coordinate system to map and explore cortical signals (Dale et al., 1999). Generating such segmentations can be a challenge at higher resolutions, and new MR contrasts usually require re-tuning of the segmentation algorithms Goebel, 2012).
Assuming we have such segmentations, an important question is how to define cortical depth between them in order to locate laminar features and define cortical profiles. Early studies simply built straight lines between surfaces, even sometimes matching vertices between surface Fig. 4. Whole brain quantitative T1 maps at four different cortical depths acquired with 500 μm isotropic resolution in-vivo (group average; N ¼ 10). The T1 maps were registered using multi-contrast multi-scale surface registration with T1 contrast only. The T1 maps in the first two rows are represented using the same colormap, highlighting the T1 gradient with cortical depth. The colormaps of the T1 maps in the bottom two rows have been adjusted for optimal contrast at each depth individually. The contrast between cortical areas varies with cortical depth, and reflects myeloarchitectonic descriptions of the cortex. Taken from Tardif et al. (2015) with permission from the publishers, Copyright Elsevier (2015). meshes, but the high level of curvature in the cerebral cortex makes this approach unrealistic beyond single profiles at individual locations. Instead, one would require the cortical depth to follow as closely as possible the laminar organization of the cortex within a given myeloarchitectonic area. Early systematic works by Bok (1929) established that folding of the cortex should preserve local volume of the layers, indicating that cortex maintains the local density of microstructure by modulating the relative thickness of layers with their curvature. This view was later discarded to favor an angle-preserving formulation, which guarantees that cortical profiles are strictly orthogonal to the lamination (Jones et al., 2000;Schleicher et al., 2005;Yezzi and Prince, 2003), with an elegant mathematical framework based on the Laplace equation.
While those differences were inconsequential at coarser resolutions, they impact the ability to locate laminar features and derive average profiles across regions with varying curvature (see Fig. 6; Waehnert et al., 2014Waehnert et al., , 2016. In order to mitigate this problem, two main approaches have been proposed: regressing mean curvature linearly from the data sampled at depths measured as the relative distance to the cortical boundaries (Sereno et al., 2013), and modeling the volume-preserving principle (Kok et al., 2016;Leprince et al., 2015;Waehnert et al., 2014). These methods have shown improved accuracy in following layer boundaries visible with MRI in ex-vivo samples at ultra high resolution, and removed most of the influence of curvature on sampled values. Noting that the relationship modeled by Bok is not linear, the volume-preserving approaches are preferable, though more complex to estimate in practice (see Fig. 6; open software implementations have been released 1 , 2 , 3 ). The biomechanical properties of cortical folding may also be more complex and require further refinement of the cortical depth models as resolution increases further: with both a physical gel-based brain model and numerical simulations, Tallinen et al., 2016 replicated the gyrification patterns occurring during development. While their study did not include layers or measure the evolution of cortical thickness with gyrification, the underlying neo-Hookerian hyper-elastic model could be used to derive a more elaborate model of thickness variations under various levels of cortical folding. Fig. 5. Dependence of myelin-sensitive MR contrast parameters (R1, MT) and T2w signal intensity on cortical depth, and its use in studying microscopic neuroanatomy. Note that layers with high myelination appear as local maxima in R1 and MT maps but as local minima in T2w images. A) T2-weighted signal intensity profiles obtained in primary visual cortex (averaged across N ¼ 11 sighted and blind subjects, respectively, and across primary visual cortex) show reduced signal in the stria of Gennari due to the higher myelination. Effects are similar in congenitally blind and sighted controls, indicating that this anatomical feature is not a developmental result of visual input, and it does not degenerate in the absence of visual input. Taken from Trampel et al. (2011) with permission from the publishers. B) R1 (¼1/T1) values as a function of cortical depth obtained in the auditory cortex (averaged across N ¼ 6 subjects and across the auditory cortex) reflect the myelination profile expected from ex-vivo myelin stains (bottom left inset) and distinguish the auditory core (TE1.0) from other auditory areas (TE1.1, TE1.2). Taken from Dick et al. (2012) with permission from the publishers. C) MT values at different cortical depths (estimated from N ¼ 297 14-24 y old volunteers and across entire cortex) increase monotonically from the pial surface to the GM/WM boundary (left). The adolescent change in MT (ΔMT over the age range of 14-24 y old) peaks at 70% cortical depth, indicating specific maturation processes (right). For comparison a sketch of the cytoarchitectonic laminae of human neocortex is provided (middle). Taken from Whitaker et al. (2016) with permission from the publishers, Copyright (2016) National Academy of Sciences.

Deriving characteristic anatomical profiles
The definition of cortical depth and layering goes hand in hand with that of cortical profiles, which would intuitively follow the gradient of depth from white matter to pial surface. Another reason for the interest in Laplacian depth was the fact that the related profiles could be defined as exactly orthogonal to the depth and not intersecting, while it is not possible for any other formulation. In order to define non-intersecting profiles for the more accurate volume-preserving depth model, the constraint of orthogonality needs to be relaxed, but current methods only provide an approximation based on following the gradient of depth. Given a profile definition, there is a second problem due to the fact that profile sample points are defined in continuous space, requiring some form of interpolation when sampling data, which inherently smooths the signal of interest across the depth of the cortex.
An alternative advocated in (Leprince et al., 2015) is to define small regions of interest that include enough voxels to sample the entire depth of the cortex but small enough to be homogeneous, although the resulting profiles have varying numbers of samples at different depths depending on the local cortical geometry. A second improvement over basic sampling is to define a "spatial GLM", i.e., measure the partial contribution of a set number of depth layers to each voxel, and invert the corresponding design matrix to obtain unmixed values. This method has been used so far for region-of-interest analysis over larger regions (Kok et al., 2016).
Independently of the method used to measure cortical profiles, there are important confounds to take into account. First, the definition of the white matter and pial boundaries are always imperfect, even in highresolution data and with careful manual delineation. In the best cases, the uncertainty in boundary location is about half the size of a voxel, which is still around 10-20% of the cortical depth for data acquired with the best in-vivo resolutions to date. Second, especially at very high resolutions, the penetrating arteries and draining veins of the cortical vasculature become visible in the MR images, depending on the contrast of choice. As these structures are oriented like cortical profiles, they generate entire profile outliers that need to be taken into account and cannot be readily removed by anisotropic spatial filters. At lower resolutions where they are not visible, they are still contributing to the signal and thus biasing the profile estimates in their vicinity. These problems will need to be addressed with more elaborate geometrical and statistical modeling as well as dedicated vascular filters (e.g., Bazin et al., 2014Bazin et al., , 2016. Alternatively, a combination of simplified measures such as the overall amount of myelination or the myelinated thickness ratio (Rowley et al., 2015) may provide a coarse but robust characterization for the cortical microarchitecture profiles.

Conclusion and outlook
The feasibility of characterizing certain cortical laminae was well established by proof-of-concept studies, which mainly focused on the stria of Gennari in the primary visual cortex (e.g., Barbier et al., 2002;Clark et al., 1992;Duyn et al., 2007;Turner et al., 2008). Recently, imaging was extended into other brain areas and the scope of the studies was extended to assessing neuroscientific questions such as layer specific plasticity. For example, Whitaker et al. (2016) demonstrated intracortical myelination changes due to maturation during adolescence. The changes, which affect the intracortical layers differently, peaked approximately at the internal layer of projection neurons in line with a genetically patterned process of consolidating anatomical network hubs.
The robust mapping of the detailed myeloarchitecture and the cortical layers across the entire cortex is not achieved yet. Several conceptual and technical challenges need to be overcome including even higher resolution scanning, understanding and modeling of the MRI contrast dependence on the cortical microstructure and advanced image processing methods handling these new types of data . A conceptual challenge lies in the large spatial resolution ranging from micrometers to centimeters when pursuing myeloarchitectural studies. The validation of the novel methods also calls for intensive post-mortem histology and creative in-vivo validation. Recent approaches for ex-vivo histology that enable quantification of microstructure (e.g., PIXE; Stüber et al., 2014) and 3D histology  facilitate the comparison to quantitative 3D MRI data. Taken together, the new developments in the areas of MRI, histology and neuroscientific validations promise a rapid development of laminar mapping and its validation.
Reliable laminar mapping would open up completely novel ways to study the human brain. Currently, investigations are mostly limited to post-mortem histology. An in-vivo approach would allow for studying subtle changes in the cortical microstructure due to pathological processes, such as in neurodegenerative diseases, and physiological changes, such as in training induced plasticity . Together with layer-specific functional MRI, a new era of studying structure-function relationships at the laminar level would become feasible. This would complement animal research and allow for studying the relationship in human cognitive functions such as distinguishing the different sub-processes involved in language production (Friederici, 2011).
The high availability of MRI would allow for large scale population studies on individual differences of structures and structure-function relationships. The approach also holds promise for clinical research and diagnostics by potentially providing sensitive biomarkers that indicate and determine the progress of pathological changes earlier than current markers.