**Optimized Signal Separation for 3D-Polarized Light Imaging**

Jürgen Dammers, Lukas Breuer, Giuseppe Tabbì and Markus Axer

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/55246

### **1. Introduction**

#### **1.1. Anatomical connectivity mapping**

"Why should we bother about connectivity in the age of functional imaging, at a time when magnets of ever increasing strength promise to detect the location of even the faintest thought? Isn't it enough to locate cortical areas engaged in deception, introspection, empathy? Do we really have to worry about their connections? The answer is "yes". In the case of the nervous system, the unit of relational architecture that allows the whole to exceed the sum of the parts is known as large-scale network. Its elucidation requires an elaborate understanding of connectivity patterns" [1]. Despite considerable advances in experimental techniques and in our understanding of animal anatomy over the last decades, the real connectivity of the human brain has essentially remained a mystery. It is the human brain's multiscale topology that poses a particular challenge to any neuroimaging technique and prevented the neuroscientists from unraveling the connectome so far.

However, it is also the brain's architecture that allows different morphological entities to be defined at different scales depending on the spatial resolution provided by the available neuroimaging techniques and the scientific objectives. Consequently, a comprehensive description of neuronal networks and their intricate fiber connections requires a multimodal approach based on complementary imaging techniques targeting different levels of organiza‐ tion (microscale, mesoscale, and macroscale) [2,3].

MR-based diffusion imaging is the most frequently used method to visualize fiber pathways in both the living and the postmortem human brain (for a comprehensive introduction to the field cf. [4,5]). Diffusion imaging contributes to the understanding of the macroscopic connec‐

© 2013 Dammers et al.; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. © 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

tivity (i.e., at the millimeter scale) in the living brain, while postmortem studies already explore the upper mesoscale (i.e., the sub-millimeter scale). Hence, not surprisingly, diffusion imaging is of great appeal to neuroscientists as a method for the visualization of connectivity patterns in both clinical and basic research. However, complex fiber networks and small fiber tracts are difficult to be disentangled reliably at present. Furthermore, the termination of fields of pathways emanating from cortical areas no larger than a few millimeters in size cannot be demonstrated with the required precision.

of-section inclination angle *α*. Basic principles of optics (Snell's law and Huygens-Fresnel principle) and the Jones calculus [21] mathematically link the measured light intensity profile

The amplitude of the profile quantifies the phase retardation *δ* induced to the light wave by the myelin. This phase retardation is a function of the light wavelength *λ*, the section thickness *d*, the birefringence *Δn* of the myelin, and the inclination angle *α*. The transmittance *I*<sup>0</sup> denotes the intensity of the incident light modified by local extinction effects. While *ρ* is the azimuth angle of the transmission axis of the first polarizer, *φ* represents the projection of the fiber axis onto the section plane with respect to the starting angle of the polarimeter. As a consequence, the fundamental data structure gained by 3D-PLI is a 3D vector field description of fibers and

For estimating the sinusoidal profile of the 3D-PLI signal Discrete Fourier Analysis (DFA) can

In 3D-PLI, the light intensity profile (cf. Figure 1) represents a crucial data set in terms of fiber orientation determination, since peak position and signal amplitude of the profile are directly related to the in-section direction and out-of-section inclination of the fiber orientation, respectively. A precise and undisturbed recording of the light intensity passing through a microtome section is therefore mandatory for the reliable reconstruction of high-resolution fiber tracts. The signal quality is however influenced by several conditions. Thermal effects and electrical noise in both the light source and the CCD-electronics deteriorate the PLI signal characteristics at each image pixel. Filter inhomogeneities of the polarizers or retarder may also manipulate the intensity profile. Therefore, a standardized image calibration technique using flat fields is usually applied to all raw images prior to analysis to compensate pixel-wise for inhomogeneities across the field-of-view [22]. Depending on the section thickness the pixelby-pixel intensity is also influenced by the scatter properties of the investigated object. Another possible source of artifacts is dust on the polarizer. Although the polarimeter should be operated in a shielded construction to prevent for external light and dust particles, dust cannot completely be avoided. As a result of the rotation of the polarizer dust particles will deteriorate the measured light intensity only once within a half circle, if and only if they are located on

fiber tract orientations – the basis for subsequent tractography and fiber modeling.

<sup>2</sup> ⋅ 1 + *sin*(2*ρ* - 2*φ*) ⋅ *sinδ* , (1)

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

357

*I* =*a*<sup>0</sup> + *a*<sup>1</sup> ⋅ sin (2*ρ*) + *b*<sup>1</sup> ⋅ cos (2*ρ*) (2)

<sup>2</sup> ⋅ sin (*δ*) ⋅ sin (2*φ*).

*I* to the fiber orientation via

*<sup>λ</sup>* <sup>⋅</sup> *cos*<sup>2</sup>

be used to deduce *I* from Equation (1):

<sup>2</sup> , *<sup>a</sup>*<sup>1</sup> <sup>=</sup> *<sup>I</sup>*<sup>0</sup>

*α*.

with *<sup>δ</sup>* <sup>≈</sup>2*<sup>π</sup>* <sup>⋅</sup> *<sup>d</sup>* ⋅Δ*<sup>n</sup>*

with *a*<sup>0</sup> <sup>=</sup> *<sup>I</sup>*<sup>0</sup>

*<sup>I</sup>* <sup>=</sup> *<sup>I</sup>*<sup>0</sup>

<sup>2</sup> <sup>⋅</sup> sin (*δ*) <sup>⋅</sup> cos (2*φ*), *<sup>b</sup>*<sup>1</sup> <sup>=</sup> *<sup>I</sup>*<sup>0</sup>

**1.3. Relevance and restoration of the light intensity profile**

Conversely, microscopic techniques generate data sets of impressing neuroanatomical detail, but they are limited to small sample sizes (i.e., small areas of interest in a small number of subjects) of postmortem tissue. This substantially restricts their predictive power. In the recent years, anatomical connections in the human postmortem brains were studied with dissection techniques [6,7], in myelin stained sections of adult human brains [8], or of immature brains taking advantage of heterochronic myelination of different fiber tracts during pre- and early postnatal development [9], in lesioned brains using various techniques for staining degener‐ ating fibers [10,11], and using tract-tracing methods for discovering local connections [12,13]. These studies have enriched our knowledge about human brain fiber tracts, but all of them suffer from severe restrictions if the 3D courses of fiber pathways are to be mapped in the adult human brain. In contrast to studies in animals, the tight packing of different fiber tracts in the white substance, and the lack of specific tracers for in vitro tracking of long distance fibers made comprehensive fiber tract mapping impossible in the adult human brain [14].

### **1.2. 3D-polarized light imaging (3D-PLI)**

Recently, a neuroimaging technique referred to as *3D-polarized light imaging (3D-PLI)* has been introduced that has opened up new avenues to study human brain regions with complex fiber architecture as well as small cortical fibers at the micrometer level [15–17]. This technique is applicable to microtome sections of postmortem brains and utilizes the optically birefringence of nerve fibers, which is induced by the optical anisotropy of the myelin sheaths surrounding axons [18–20]. 3D-PLI provides a three-dimensional description of the anatomical connectivity in form of a vector field indicating the prevailing fiber orientation per voxel. Depending on the used imaging setup and the section thickness the acquirable spatial resolution is 1.6 – 100 µm. Hence, the method bridges the microscopic and the macroscopic description of the anatomical connectivity of the human brain.

The birefringence of brain tissue is measured by passing linearly polarized light through histological brain sections and by detecting local changes in the polarization state of light using a polarimeter setup (Figure 1a-b). The polarimeter is equipped with a pair of crossed polarizers, a tilting specimen stage, a quarter-wave retarder, an LED light source (with a narrow-band green wavelength spectrum), and a charge-coupled device (CCD) camera. By rotating the optical devices simultaneously around the stationary brain section and by imaging the sample at discrete rotation angles *ρ*, a sinusoidal variation of the measured light intensity (i.e., the light intensity profile) is observed for each image pixel (Figure 1c). The individual course of a light intensity profile essentially depends on the locally prevailing 3D fiber orientation described by an orientation unit vector which is defined by the in-section direction angle *φ* and the outof-section inclination angle *α*. Basic principles of optics (Snell's law and Huygens-Fresnel principle) and the Jones calculus [21] mathematically link the measured light intensity profile *I* to the fiber orientation via

$$I = \frac{I\_0}{2} \cdot \left[1 + \sin(2\rho \cdot 2\rho) \cdot \sin \delta\right] \tag{1}$$

with *<sup>δ</sup>* <sup>≈</sup>2*<sup>π</sup>* <sup>⋅</sup> *<sup>d</sup>* ⋅Δ*<sup>n</sup> <sup>λ</sup>* <sup>⋅</sup> *cos*<sup>2</sup> *α*.

tivity (i.e., at the millimeter scale) in the living brain, while postmortem studies already explore the upper mesoscale (i.e., the sub-millimeter scale). Hence, not surprisingly, diffusion imaging is of great appeal to neuroscientists as a method for the visualization of connectivity patterns in both clinical and basic research. However, complex fiber networks and small fiber tracts are difficult to be disentangled reliably at present. Furthermore, the termination of fields of pathways emanating from cortical areas no larger than a few millimeters in size cannot be

Conversely, microscopic techniques generate data sets of impressing neuroanatomical detail, but they are limited to small sample sizes (i.e., small areas of interest in a small number of subjects) of postmortem tissue. This substantially restricts their predictive power. In the recent years, anatomical connections in the human postmortem brains were studied with dissection techniques [6,7], in myelin stained sections of adult human brains [8], or of immature brains taking advantage of heterochronic myelination of different fiber tracts during pre- and early postnatal development [9], in lesioned brains using various techniques for staining degener‐ ating fibers [10,11], and using tract-tracing methods for discovering local connections [12,13]. These studies have enriched our knowledge about human brain fiber tracts, but all of them suffer from severe restrictions if the 3D courses of fiber pathways are to be mapped in the adult human brain. In contrast to studies in animals, the tight packing of different fiber tracts in the white substance, and the lack of specific tracers for in vitro tracking of long distance fibers

made comprehensive fiber tract mapping impossible in the adult human brain [14].

Recently, a neuroimaging technique referred to as *3D-polarized light imaging (3D-PLI)* has been introduced that has opened up new avenues to study human brain regions with complex fiber architecture as well as small cortical fibers at the micrometer level [15–17]. This technique is applicable to microtome sections of postmortem brains and utilizes the optically birefringence of nerve fibers, which is induced by the optical anisotropy of the myelin sheaths surrounding axons [18–20]. 3D-PLI provides a three-dimensional description of the anatomical connectivity in form of a vector field indicating the prevailing fiber orientation per voxel. Depending on the used imaging setup and the section thickness the acquirable spatial resolution is 1.6 – 100 µm. Hence, the method bridges the microscopic and the macroscopic description of the

The birefringence of brain tissue is measured by passing linearly polarized light through histological brain sections and by detecting local changes in the polarization state of light using a polarimeter setup (Figure 1a-b). The polarimeter is equipped with a pair of crossed polarizers, a tilting specimen stage, a quarter-wave retarder, an LED light source (with a narrow-band green wavelength spectrum), and a charge-coupled device (CCD) camera. By rotating the optical devices simultaneously around the stationary brain section and by imaging the sample at discrete rotation angles *ρ*, a sinusoidal variation of the measured light intensity (i.e., the light intensity profile) is observed for each image pixel (Figure 1c). The individual course of a light intensity profile essentially depends on the locally prevailing 3D fiber orientation described by an orientation unit vector which is defined by the in-section direction angle *φ* and the out-

demonstrated with the required precision.

356 Functional Brain Mapping and the Endeavor to Understand the Working Brain

**1.2. 3D-polarized light imaging (3D-PLI)**

anatomical connectivity of the human brain.

The amplitude of the profile quantifies the phase retardation *δ* induced to the light wave by the myelin. This phase retardation is a function of the light wavelength *λ*, the section thickness *d*, the birefringence *Δn* of the myelin, and the inclination angle *α*. The transmittance *I*<sup>0</sup> denotes the intensity of the incident light modified by local extinction effects. While *ρ* is the azimuth angle of the transmission axis of the first polarizer, *φ* represents the projection of the fiber axis onto the section plane with respect to the starting angle of the polarimeter. As a consequence, the fundamental data structure gained by 3D-PLI is a 3D vector field description of fibers and fiber tract orientations – the basis for subsequent tractography and fiber modeling.

For estimating the sinusoidal profile of the 3D-PLI signal Discrete Fourier Analysis (DFA) can be used to deduce *I* from Equation (1):

$$I = a\_0 + a\_1 \cdot \sin\left(2\rho\right) + b\_1 \cdot \cos\left(2\rho\right) \tag{2}$$

with *a*<sup>0</sup> <sup>=</sup> *<sup>I</sup>*<sup>0</sup> , *<sup>a</sup>*<sup>1</sup> <sup>=</sup> *<sup>I</sup>*<sup>0</sup> <sup>⋅</sup> sin (*δ*) <sup>⋅</sup> cos (2*φ*), *<sup>b</sup>*<sup>1</sup> <sup>=</sup> *<sup>I</sup>*<sup>0</sup> ⋅ sin (*δ*) ⋅ sin (2*φ*).

#### **1.3. Relevance and restoration of the light intensity profile**

In 3D-PLI, the light intensity profile (cf. Figure 1) represents a crucial data set in terms of fiber orientation determination, since peak position and signal amplitude of the profile are directly related to the in-section direction and out-of-section inclination of the fiber orientation, respectively. A precise and undisturbed recording of the light intensity passing through a microtome section is therefore mandatory for the reliable reconstruction of high-resolution fiber tracts. The signal quality is however influenced by several conditions. Thermal effects and electrical noise in both the light source and the CCD-electronics deteriorate the PLI signal characteristics at each image pixel. Filter inhomogeneities of the polarizers or retarder may also manipulate the intensity profile. Therefore, a standardized image calibration technique using flat fields is usually applied to all raw images prior to analysis to compensate pixel-wise for inhomogeneities across the field-of-view [22]. Depending on the section thickness the pixelby-pixel intensity is also influenced by the scatter properties of the investigated object. Another possible source of artifacts is dust on the polarizer. Although the polarimeter should be operated in a shielded construction to prevent for external light and dust particles, dust cannot completely be avoided. As a result of the rotation of the polarizer dust particles will deteriorate the measured light intensity only once within a half circle, if and only if they are located on the rotating system (cf. Figure 1). Hence, the measured intensity at the CCD array is a linear mixture of different light sources.

Recently, signal enhancement and restoration techniques for PLI images utilizing independent component analysis (ICA) have been introduced [22,23]. As a result, component maps corresponding to gray and white matter structures as well as noise and artifacts can be identified automatically using statistical analysis tools [23]. Remarkably, even in the presence of dust on the polarimeter ICA can effectively restore the original sinusoidal signal (Figure 2). After ICA filtering the noise and artifacts are removed and the sinusoidal nature of the PLI signal is restored (Figure 2).

**Figure 2.** In (a) multiple artifacts including dust contamination are present. The trajectory of a dust particle is high‐ lighted by the signal power inside the white circles along all rotation angles. After ICA filtering (b) the signal deroga‐ tion is greatly removed. At the bottom, PLI raw signals (black) are shown together with corresponding signals after artifact rejection (red). The signals were extracted from the two locations as indicated by the red dots. In comparison to the PLI raw data ICA obviously removed the signal derogation and is capable of restoring the original sinusoidal

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

359

**2. A new concept for optimal signal decomposition in polarized light**

Although blind source separation methods are successfully applied for signal separation in all kinds of neuroimaging modalities [24–29] two major problems remain: i) the method applied must be selected carefully and the separation strategy (i.e., the internal cost function) of the applied method should be optimal for the type of data. This however, is one of the reasons why many different ICA algorithms co-exist. ii) Assuming the measured signals are adequately separated into the underlying source signals (i.e., signal of interest and non-interest), identi‐ fication of the signal of interest should be performed user independently; preferably in an automatic fashion. For the latter, significant effort has to be invested in order to identify components of interest automatically from data recorded utilizing different neuroimaging techniques [29–33]. In contrast to many other ICA applications, the great advantage in PLI signal decomposition is that all profiles of the basis vectors that correspond to the signal of

signal in case of dust (signal #1) or noise (signal #2) contaminations.

**imaging**

**Figure 1.** Polarimetry at a glance. (a) Scheme of the rotating polarimeter with tilting stage (N-North, W-West, E-East, S-South). (b) Scheme of the optical fiber model. The refractive index of a negative uniaxially birefringent medium, such as a myelinated axon, is described by an elliptically shaped oblate surface, the indicatrix (gray mesh). A beam of linear‐ ly polarized light (blue trace) interacts locally with the myelin sheath of a single axon (black line) and becomes ellipti‐ cally polarized. (c) A typical 3D-PLI raw image data set consists of 18 images corresponding to equidistant rotation angles between 0° and 170°. Here, a selection of four images of a coronal section is shown, while the sketched arrow indicates one representative pixel. To obtain the fiber orientation, the measured light intensities are studied pixel-wise as a function of discrete rotation angles. The derived physical model provides a precise mathematical description of the measurement (continuous black line) and relates the sine phase to the direction angle φ and the amplitude to the inclination angle α. The highlighted data points correspond to the selected images.

the rotating system (cf. Figure 1). Hence, the measured intensity at the CCD array is a linear

Recently, signal enhancement and restoration techniques for PLI images utilizing independent component analysis (ICA) have been introduced [22,23]. As a result, component maps corresponding to gray and white matter structures as well as noise and artifacts can be identified automatically using statistical analysis tools [23]. Remarkably, even in the presence of dust on the polarimeter ICA can effectively restore the original sinusoidal signal (Figure 2). After ICA filtering the noise and artifacts are removed and the sinusoidal nature of the PLI

**Figure 1.** Polarimetry at a glance. (a) Scheme of the rotating polarimeter with tilting stage (N-North, W-West, E-East, S-South). (b) Scheme of the optical fiber model. The refractive index of a negative uniaxially birefringent medium, such as a myelinated axon, is described by an elliptically shaped oblate surface, the indicatrix (gray mesh). A beam of linear‐ ly polarized light (blue trace) interacts locally with the myelin sheath of a single axon (black line) and becomes ellipti‐ cally polarized. (c) A typical 3D-PLI raw image data set consists of 18 images corresponding to equidistant rotation angles between 0° and 170°. Here, a selection of four images of a coronal section is shown, while the sketched arrow indicates one representative pixel. To obtain the fiber orientation, the measured light intensities are studied pixel-wise as a function of discrete rotation angles. The derived physical model provides a precise mathematical description of the measurement (continuous black line) and relates the sine phase to the direction angle φ and the amplitude to the

inclination angle α. The highlighted data points correspond to the selected images.

mixture of different light sources.

358 Functional Brain Mapping and the Endeavor to Understand the Working Brain

signal is restored (Figure 2).

**Figure 2.** In (a) multiple artifacts including dust contamination are present. The trajectory of a dust particle is high‐ lighted by the signal power inside the white circles along all rotation angles. After ICA filtering (b) the signal deroga‐ tion is greatly removed. At the bottom, PLI raw signals (black) are shown together with corresponding signals after artifact rejection (red). The signals were extracted from the two locations as indicated by the red dots. In comparison to the PLI raw data ICA obviously removed the signal derogation and is capable of restoring the original sinusoidal signal in case of dust (signal #1) or noise (signal #2) contaminations.

### **2. A new concept for optimal signal decomposition in polarized light imaging**

Although blind source separation methods are successfully applied for signal separation in all kinds of neuroimaging modalities [24–29] two major problems remain: i) the method applied must be selected carefully and the separation strategy (i.e., the internal cost function) of the applied method should be optimal for the type of data. This however, is one of the reasons why many different ICA algorithms co-exist. ii) Assuming the measured signals are adequately separated into the underlying source signals (i.e., signal of interest and non-interest), identi‐ fication of the signal of interest should be performed user independently; preferably in an automatic fashion. For the latter, significant effort has to be invested in order to identify components of interest automatically from data recorded utilizing different neuroimaging techniques [29–33]. In contrast to many other ICA applications, the great advantage in PLI signal decomposition is that all profiles of the basis vectors that correspond to the signal of interest must show a sinusoidal waveform [16,22]. Since these waveforms can only vary in amplitude and phase, an automatic extraction of the signal components can be achieved utilizing this property as demonstrated in [23].

Alternatively, it has been demonstrated that by incorporating prior knowledge, i.e., by impos‐ ingtemporalorspatialconstraintsforthesourceseparationtask,decompositioncanbeeffective‐ ly improved [34–37]. Hesse and James [38], for example, showed that using different types of spatial constraints ICA can be trained to *i)* identify ocular artifacts automatically and *ii)* to detect and trace ictal activity. Liu and colleagues recently presented a reference based ICA concept to successfully target a specific genetic variation [36] in real and simulated data sets. For the investigation of signal detection in the imaging system of the retina in cats Barriga and collea‐ guesintroducedaconstrainedICAalgorithmwhichincreasesthedetectionofresponsestovisual stimuliincatsevenforlowlevelsofsignal-to-noiseratio[37].Inthenextsectionthebasicprinciple of ICA is described shortly. For a more detailed review, we refer to [39,40].

#### **2.1. A short introduction to independent component analysis (ICA)**

For 3D-PLI a linear superposition of light at the CCD camera is assumed, where each elemen‐ tary signal component refers to a distinct region in space. By applying independent component analysis (ICA) to a set of polarized light images (here a stack of 18 images at different rotation angles is used) the decomposition of the data results in spatially independent components (often called feature or basis vectors) yielding maximally (i.e., statistically) independent spatial component maps.

Let **X**=(X1, X2, …, XN)*<sup>T</sup>* be the *N*-dimensional measured PLI signal mixture, where each X*<sup>i</sup>* is one image mixture (flatted to a one-dimensional image vector with *M* pixels) detected at a corresponding rotation angle *ρ*. In spatial ICA **X** is considered to be a *N* × *M* matrix (as op‐ posed to temporal ICA, where the dimension of the matrix is transposed), with *M* = number of pixels included in the analysis and *N* = number of rotation angles reflecting *N* different instances of the signal. The contribution of each source image varies *N* times over the angles *ρ*. Similarly, **S**=(S1, S2, …, S*<sup>N</sup>* )*T* represents the *N*-dimensional true source signals. Hence, the linear relationship of mixed sources can be expressed as follows:

$$\mathbf{X} = \mathbf{A}\mathbf{S} \tag{3}$$

information and thus is maximally independent. These independent components (or spatial

independent "temporal" profiles (i.e., the contrast changing function along the rotation angle,

signal contributions from unwanted (e.g. artifact) sources. This is identical to zeroing rows in **C**, where the independent component maps of interest are stored in **C**'. The cleaning of measured

Signal separation in constrained ICA is based on incorporating prior knowledge about the underlying signals. In the study of Barriga and colleagues [37], *a priori* information is incor‐ porated into one column of the estimated mixing matrix of the spatial decomposition via an off-on constraint, meaning that the visual response signal used in the study has to follow the stimulus signal. For this, the authors filled one column of the mixing matrix with ones for the time of the stimulation, while zero values reflect off-stimulation. Similar to the concept introduced by Barriga and colleagues in 2011, it is possible to incorporate the *a priori* informa‐ tion from 3D-PLI, but with the major difference that for birefringent signals an analytical expression exists to model the expected signal of interest, which only differs in amplitude and phase. Therefore, a new algorithm motivated by the results described in [23] has been devel‐ oped, where signal separation and the identification of the underlying 3D-PLI source signals are combined into one method. The new algorithm (referred to as *constrained ICA for polarized light imaging*, *cICAP*) is based on a modified version of the Infomax principle [41] by means of incorporating prior information to the cost function of Infomax; this is similar to the work

Infomax cost function to control for the results of the decomposition. Recall that if the *i-*th

Starting from Equation (4) the natural-gradient version in Infomax is used to update the weight

(2), where the parameters *a0*, *a1* and *b1* are fitted involving all rotation angles *ρ* [16].

^ has a sinusoidal profile, the corresponding component map in **C** represents a

∆**W**= **I** + (1 - 2y)**C***\*T* **W** (6)

Signal restoration, i.e. the cleaning process, is performed by zeroing columns in **W**-1

**X**' = **A**

**2.2. Constrained ICA for polarized light imaging (***cICAP***)**

reported in [34,37]. In *cICAP* a theoretically expected function *f <sup>i</sup>*

signal of interest (cf. Equation (4-5)). The expected function *f <sup>i</sup>*

matrix **W** during the ICA learning procedure [42,43]:

data is performed by back transformation of **C**', which results in a new set of PLI data **X**'.

the associated spatially

http://dx.doi.org/10.5772/55246

(*uρ*) is implemented into the

(*uρ*) can be derived from Equation

^ <sup>⋅</sup>**C**' (5)

Optimized Signal Separation for 3D-Polarized Light Imaging

which reflect

361

maps) are stored in the rows of **C**, while in the columns of **W**-1

the basis vector) are stored.

with **W**-1 =**A**

column in **A**

with y=1 /(1 + *e*-**C***\**

).

^ .

with **A** being the unknown mixing matrix of dimension *N* × *N*. The key problem in ICA is to find an unmixing matrix **W** (similar to a pseudo-inverse **A**-1, with **W**≈**A**-1) while imposing that the sources in **S** are statistically independent, that is:

$$\mathbf{C} = \mathbf{W}\mathbf{X}.\tag{4}$$

Within ICA the *N*-dimensional data array **X** is transformed into an *N*-dimensional component space **C**, where each of the spatial component maps carries a minimum amount of mutual information and thus is maximally independent. These independent components (or spatial maps) are stored in the rows of **C**, while in the columns of **W**-1 the associated spatially independent "temporal" profiles (i.e., the contrast changing function along the rotation angle, the basis vector) are stored.

Signal restoration, i.e. the cleaning process, is performed by zeroing columns in **W**-1 which reflect signal contributions from unwanted (e.g. artifact) sources. This is identical to zeroing rows in **C**, where the independent component maps of interest are stored in **C**'. The cleaning of measured data is performed by back transformation of **C**', which results in a new set of PLI data **X**'.

$$\mathbf{X}^{\cdot} = \hat{\mathbf{A}} \cdot \mathbf{C}^{\cdot} \tag{5}$$

with **W**-1 =**A** ^ .

interest must show a sinusoidal waveform [16,22]. Since these waveforms can only vary in amplitude and phase, an automatic extraction of the signal components can be achieved

Alternatively, it has been demonstrated that by incorporating prior knowledge, i.e., by impos‐ ingtemporalorspatialconstraintsforthesourceseparationtask,decompositioncanbeeffective‐ ly improved [34–37]. Hesse and James [38], for example, showed that using different types of spatial constraints ICA can be trained to *i)* identify ocular artifacts automatically and *ii)* to detect and trace ictal activity. Liu and colleagues recently presented a reference based ICA concept to successfully target a specific genetic variation [36] in real and simulated data sets. For the investigation of signal detection in the imaging system of the retina in cats Barriga and collea‐ guesintroducedaconstrainedICAalgorithmwhichincreasesthedetectionofresponsestovisual stimuliincatsevenforlowlevelsofsignal-to-noiseratio[37].Inthenextsectionthebasicprinciple

For 3D-PLI a linear superposition of light at the CCD camera is assumed, where each elemen‐ tary signal component refers to a distinct region in space. By applying independent component analysis (ICA) to a set of polarized light images (here a stack of 18 images at different rotation angles is used) the decomposition of the data results in spatially independent components (often called feature or basis vectors) yielding maximally (i.e., statistically) independent spatial

Let **X**=(X1, X2, …, XN)*<sup>T</sup>* be the *N*-dimensional measured PLI signal mixture, where each X*<sup>i</sup>* is one image mixture (flatted to a one-dimensional image vector with *M* pixels) detected at a corresponding rotation angle *ρ*. In spatial ICA **X** is considered to be a *N* × *M* matrix (as op‐ posed to temporal ICA, where the dimension of the matrix is transposed), with *M* = number of pixels included in the analysis and *N* = number of rotation angles reflecting *N* different instances of the signal. The contribution of each source image varies *N* times over the angles *ρ*. Similarly, **S**=(S1, S2, …, S*<sup>N</sup>* )*T* represents the *N*-dimensional true source signals. Hence,

with **A** being the unknown mixing matrix of dimension *N* × *N*. The key problem in ICA is to find an unmixing matrix **W** (similar to a pseudo-inverse **A**-1, with **W**≈**A**-1) while imposing that

Within ICA the *N*-dimensional data array **X** is transformed into an *N*-dimensional component space **C**, where each of the spatial component maps carries a minimum amount of mutual

**X**=**AS** (3)

**C**=**WX**. (4)

of ICA is described shortly. For a more detailed review, we refer to [39,40].

**2.1. A short introduction to independent component analysis (ICA)**

the linear relationship of mixed sources can be expressed as follows:

the sources in **S** are statistically independent, that is:

utilizing this property as demonstrated in [23].

360 Functional Brain Mapping and the Endeavor to Understand the Working Brain

component maps.

#### **2.2. Constrained ICA for polarized light imaging (***cICAP***)**

Signal separation in constrained ICA is based on incorporating prior knowledge about the underlying signals. In the study of Barriga and colleagues [37], *a priori* information is incor‐ porated into one column of the estimated mixing matrix of the spatial decomposition via an off-on constraint, meaning that the visual response signal used in the study has to follow the stimulus signal. For this, the authors filled one column of the mixing matrix with ones for the time of the stimulation, while zero values reflect off-stimulation. Similar to the concept introduced by Barriga and colleagues in 2011, it is possible to incorporate the *a priori* informa‐ tion from 3D-PLI, but with the major difference that for birefringent signals an analytical expression exists to model the expected signal of interest, which only differs in amplitude and phase. Therefore, a new algorithm motivated by the results described in [23] has been devel‐ oped, where signal separation and the identification of the underlying 3D-PLI source signals are combined into one method. The new algorithm (referred to as *constrained ICA for polarized light imaging*, *cICAP*) is based on a modified version of the Infomax principle [41] by means of incorporating prior information to the cost function of Infomax; this is similar to the work reported in [34,37]. In *cICAP* a theoretically expected function *f <sup>i</sup>* (*uρ*) is implemented into the Infomax cost function to control for the results of the decomposition. Recall that if the *i-*th column in **A** ^ has a sinusoidal profile, the corresponding component map in **C** represents a signal of interest (cf. Equation (4-5)). The expected function *f <sup>i</sup>* (*uρ*) can be derived from Equation (2), where the parameters *a0*, *a1* and *b1* are fitted involving all rotation angles *ρ* [16].

Starting from Equation (4) the natural-gradient version in Infomax is used to update the weight matrix **W** during the ICA learning procedure [42,43]:

$$
\Delta \mathbf{W} = \left[ \mathbf{I} + (\mathbf{1} \cdot \mathbf{2} \mathbf{y}) \mathbf{C}^{\ast T} \right] \mathbf{W} \tag{6}
$$

with y=1 /(1 + *e*-**C***\** ). In Equation (6) **I** and denote the identity matrix and the learning rate, respectively, as it is used in the standard Infomax algorithm [41]. **C\*** denotes the intermediate component matrix which is identical to **C** when the learning procedure is finished. In common with the approach of [37], the mixing matrix **A** ^ is estimated by **<sup>A</sup>** ^ <sup>≔</sup>(**P**⋅**W**)-1 , with **P** being the sphering matrix of the decorrelation process within Infomax. By incorporating the *a priori* information at each iteration step, **A** ^ is updated by:

$$\stackrel{\wedge}{\mathbf{a}}\_{i}^{\wedge} = (1 \cdot c) \cdot \stackrel{\wedge}{\mathbf{a}}\_{i} + \stackrel{\wedge}{c} \cdot \stackrel{\wedge}{\cdot} f\_{i}(\mathbf{u}\_{\rho}) \tag{7}$$

*MSE* (with *ε*≫*t*) is used to prevent the algorithm of running into an endless loop. This

For the weight update described in Equation (10) the choice of the confidence value *c* and the tolerance values *t* and *ε* should be performed in a representative but independent data set,

For measuring the noise reduction and signal enhancement in a set of 3D-PLI data after ICA application we use the weighted reduced chi-squared statistic as introduced in [23]. The

> *<sup>N</sup>* (*I<sup>ρ</sup>* - *f* (*uρ*))2 *σρ*

obtained from 100 flat field images, which were independently generated for the image calibration process [22]. The weighted sum over *N* angles is further normalized by the degrees of freedom υ. The normalization ensures a good fit when the reduced chi-square value equals one, which is achieved when the squared difference between the measurement and the expected function resemble the variance of the measurement. In order to measure the signal improvement (and the reduction of noise) through the ICA process the test statistic is deter‐

weight factor *ω* that penalizes the goodness of fit, whenever a signal component is missing

The weight factor ω increases when the squared difference between the two expectation functions *f raw*(*uρ*) and *f ICA*(*uρ*) (derived from the raw and the ICA filtered PLI data sets, respectively) is large. In case of missed components of interest, the signal strength of the ICA filtered 3D-PLI signal would be largely reduced at corresponding pixel locations. Consequent‐

2 *ω* ⋅ *χICA*

*wrGOF* (*x*, *<sup>y</sup>*)= *<sup>χ</sup>raw*

involves the variance σ<sup>2</sup>

*<sup>χ</sup>* <sup>2</sup> <sup>=</sup> <sup>1</sup>

Here, *I<sup>ρ</sup>* denote the intensity measured at angle *ρ*. The variance *σ*<sup>2</sup>

*<sup>ν</sup>* ⋅ ∑ *ρ*

) <*ε* ∧ *MSE*(*a*

^ *i* , *f <sup>i</sup>*

(*uρ*)) >*t*

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

of the observation, where the statistical

used in Equation (11) was

<sup>2</sup> . (11)

<sup>2</sup> , (12)

<sup>2</sup> ) ICA. Using the ratio of these values and introducing a

(10)

363

is expressed in Equation (10):

(*uρ*) , *if kurt*(*di*

*<sup>i</sup>* , *else*

differentiation of updating *a*

*a* ^ *i \** ={

which is described later.

reduced chi-squared statistic χ<sup>2</sup>

mined before (*χraw*

with *<sup>ω</sup>* <sup>=</sup> <sup>1</sup>

*<sup>ν</sup>* ⋅ ∑ *ρ*

(1 - *c*) ⋅ *a* ^ *<sup>i</sup>* + *c* ⋅ *f <sup>i</sup>*

*2.2.1. Evaluation of the signal enhancement*

^ *i \**

*a* ^

output is weighted based on the measurement error

<sup>2</sup> ) and after (*χICA*

[23]. The weighted test statistic reads as follows:

*<sup>N</sup>* ( *f raw*(*uρ*) - *f ICA*(*uρ*))2 *σρ*

<sup>2</sup> ≥1.

where a ^ *i* refers to the *i*-th column of the estimated mixing matrix **A** ^ and *<sup>c</sup>* is a confidence parameter ranging from 0 to 1. The confidence parameter *c* in Equation (7) indicates the percentage of weighting of the prior knowledge by means of the expected function *f <sup>i</sup>* (*uρ*). In the limit of *c*→0 the original Infomax weight update is used (cf. Equation (6)). It is important to note that for each iteration step only one column in **A** ^ is updated as expressed in Equation (7) since we do not want to influence the decomposition of the remaining components during the unmixing process. The column in **A** ^ which entries are most similar to the corresponding expectation function is selected and modified according to Equation (7). As a measure for similarity the kurtosis of the deviation function *di* = *a* ^ *<sup>i</sup>* - *f <sup>i</sup>* (*uρ*) is used:

*kurt*(*di* ) <sup>=</sup> <sup>1</sup> *<sup>N</sup>* ⋅ ∑ *ρ N* ( *di*,*<sup>ρ</sup>* - *<sup>d</sup>*¯ *i σi* ) 4 - 3, (8)

where *d* ¯ *i* is the expected value of the deviation function of the *i*-th component. In the ideal case *di* will be zero (or close to zero); for any deviation at one or more corresponding angles *ρ* the measure *kurt*(*di* ) will largely be increased. For a subsequent selection of columns in **A** ^ the column profile with smallest kurtosis value as expressed by Equation (8) is used first. More‐ over, if the mean square error (*MSE*) between the updated column *a* ^ *i \** and the corresponding theoretically expectation function *f <sup>i</sup>* (*uρ*) is less than a predefined tolerance value *t* the entries of *a* ^ *i \** are fixed. Thus, the *MSE* is expressed by:

$$MSE\left(\stackrel{\scriptstyle \sim\_{\star}}{a}\_{i^{\prime}} \mid f\_{i}(\mu\_{\rho})\right) = \frac{1}{N} \quad \cdot \stackrel{\scriptstyle \sim\_{N}}{\sum} \{\stackrel{\scriptstyle \sim\_{\star}}{a}\_{i}(\rho) - \mid f\_{i}(\mu\_{\rho})\}^{2} . \tag{9}$$

The basic idea here is that a very small *MSE* (i.e., *MSE* ≤ *t*) in Equation (9) indicates that *a* ^ *i \** and *f i* (*uρ*) are very similar, which means that the *i*-th column in **A** ^ then represents a source of interest. This approach is repeated for the remaining columns in **A** ^ until no further components of interest are found. As a stopping criterion for the iteration, a predefined threshold, *ε,* for the *MSE* (with *ε*≫*t*) is used to prevent the algorithm of running into an endless loop. This differentiation of updating *a* ^ *i \** is expressed in Equation (10):

$$\begin{array}{ll} a\_{\uparrow \star} \stackrel{\scriptstyle \Lambda\_{\star}}{=} \begin{pmatrix} (1 \cdot \varepsilon) \cdot \stackrel{\scriptstyle \Lambda\_{i}}{=} \bar{a}\_{i} + \varepsilon \cdot f\_{i}(u\_{\rho}) & \text{, if } \text{kurt}(\text{d}\_{i}) < \varepsilon \wedge \text{MSE} \left(\stackrel{\scriptstyle \Lambda\_{i}}{=} \iota\_{i} \cdot f\_{i}(u\_{\rho})\right) > \text{t} \\ & \stackrel{\scriptstyle \Lambda\_{i}}{=} \iota\_{i} \end{array} \tag{10}$$

For the weight update described in Equation (10) the choice of the confidence value *c* and the tolerance values *t* and *ε* should be performed in a representative but independent data set, which is described later.

#### *2.2.1. Evaluation of the signal enhancement*

In Equation (6) **I** and denote the identity matrix and the learning rate, respectively, as it is used in the standard Infomax algorithm [41]. **C\*** denotes the intermediate component matrix which is identical to **C** when the learning procedure is finished. In common with the approach

the decorrelation process within Infomax. By incorporating the *a priori* information at each

*<sup>i</sup>* + *c* ⋅ *f <sup>i</sup>*

parameter ranging from 0 to 1. The confidence parameter *c* in Equation (7) indicates the

the limit of *c*→0 the original Infomax weight update is used (cf. Equation (6)). It is important

(7) since we do not want to influence the decomposition of the remaining components during

expectation function is selected and modified according to Equation (7). As a measure for

^ *<sup>i</sup>* - *f <sup>i</sup>*

is the expected value of the deviation function of the *i*-th component. In the ideal case

) will largely be increased. For a subsequent selection of columns in **A**

*di* will be zero (or close to zero); for any deviation at one or more corresponding angles *ρ* the

column profile with smallest kurtosis value as expressed by Equation (8) is used first. More‐

*<sup>N</sup>* ⋅ ∑ *ρ N* (*a* ^ *i \**(*ρ*) - *<sup>f</sup> <sup>i</sup>*

The basic idea here is that a very small *MSE* (i.e., *MSE* ≤ *t*) in Equation (9) indicates that *a*

of interest are found. As a stopping criterion for the iteration, a predefined threshold, *ε,* for the

^ <sup>≔</sup>(**P**⋅**W**)-1

, with **P** being the sphering matrix of

(*uρ*) (7)

^ is updated as expressed in Equation


^ *i \**

(*uρ*) is less than a predefined tolerance value *t* the entries

(*uρ*))<sup>2</sup>

^ which entries are most similar to the corresponding

(*uρ*) is used:

^ and *<sup>c</sup>* is a confidence

(*uρ*). In

^ the

^ *i \** and

and the corresponding

. (9)

^ then represents a source of

^ until no further components

^ is estimated by **<sup>A</sup>**

**\***=(1 - *c*) ⋅ **a**

^

percentage of weighting of the prior knowledge by means of the expected function *f <sup>i</sup>*

refers to the *i*-th column of the estimated mixing matrix **A**

**a** ^ *i*

362 Functional Brain Mapping and the Endeavor to Understand the Working Brain

to note that for each iteration step only one column in **A**

similarity the kurtosis of the deviation function *di* = *a*

*kurt*(*di*

) <sup>=</sup> <sup>1</sup> *<sup>N</sup>* ⋅ ∑ *ρ N* ( *di*,*<sup>ρ</sup>* - *<sup>d</sup>*¯ *i σi* ) 4

over, if the mean square error (*MSE*) between the updated column *a*

(*uρ*) are very similar, which means that the *i*-th column in **A**

interest. This approach is repeated for the remaining columns in **A**

(*uρ*)) <sup>=</sup> <sup>1</sup>

the unmixing process. The column in **A**

theoretically expectation function *f <sup>i</sup>*

are fixed. Thus, the *MSE* is expressed by:

*MSE*(*a* ^ *i \** , *f <sup>i</sup>*

of [37], the mixing matrix **A**

^ is updated by:

iteration step, **A**

where a ^ *i*

where *d* ¯ *i*

of *a* ^ *i \**

*f i*

measure *kurt*(*di*

For measuring the noise reduction and signal enhancement in a set of 3D-PLI data after ICA application we use the weighted reduced chi-squared statistic as introduced in [23]. The reduced chi-squared statistic χ<sup>2</sup> involves the variance σ<sup>2</sup> of the observation, where the statistical output is weighted based on the measurement error

$$\propto \chi^2 = \frac{1}{\nu} \quad \cdot \sum\_{\rho}^{N} \frac{(l\_{\rho} \cdot f\left(u\_{\rho}\right))^2}{\sigma\_{\rho}^{\frac{2}{\nu}}}.\tag{11}$$

Here, *I<sup>ρ</sup>* denote the intensity measured at angle *ρ*. The variance *σ*<sup>2</sup> used in Equation (11) was obtained from 100 flat field images, which were independently generated for the image calibration process [22]. The weighted sum over *N* angles is further normalized by the degrees of freedom υ. The normalization ensures a good fit when the reduced chi-square value equals one, which is achieved when the squared difference between the measurement and the expected function resemble the variance of the measurement. In order to measure the signal improvement (and the reduction of noise) through the ICA process the test statistic is deter‐ mined before (*χraw* <sup>2</sup> ) and after (*χICA* <sup>2</sup> ) ICA. Using the ratio of these values and introducing a weight factor *ω* that penalizes the goodness of fit, whenever a signal component is missing [23]. The weighted test statistic reads as follows:

$$
bar{GOF}\{\mathbf{x},\;y\} = \frac{\chi^2\_{nm}}{\omega \cdot \chi^2\_{\mathbb{C}\mathcal{A}}}\;/\tag{12}$$

with *<sup>ω</sup>* <sup>=</sup> <sup>1</sup> *<sup>ν</sup>* ⋅ ∑ *ρ <sup>N</sup>* ( *f raw*(*uρ*) - *f ICA*(*uρ*))2 *σρ* <sup>2</sup> ≥1.

The weight factor ω increases when the squared difference between the two expectation functions *f raw*(*uρ*) and *f ICA*(*uρ*) (derived from the raw and the ICA filtered PLI data sets, respectively) is large. In case of missed components of interest, the signal strength of the ICA filtered 3D-PLI signal would be largely reduced at corresponding pixel locations. Consequent‐ ly the expectation function *f ICA*(*uρ*) would differ in terms of its amplitude, while the shape of the waveform may be still sinusoidal. The assumption here is that the signal power (across the rotation angles) of the signal of interest is larger compared to noise. The restriction of ω ≥ 1 is needed in order to not artificially improve the goodness-of-fit in cases where the two expect‐ ation functions *f raw*(*uρ*) and *f ICA*(*uρ*) are very similar (e.g., at pixel locations with low noise).

#### *2.2.2. Finding the optimal parameters*

The optimal parameters *c*, *t* and *ε* for the weight update as described above (cf. Equation (10)), should be determined using a representative but independent data set. In the test data from the same post-mortem brain were used, but from a different microtome section which is not included in the subsequent analysis. The determination of the optimal parameters for *cICAP* requires the maximization of the weighted goodness-of-fit (*wrGOF*) values across all pixel locations (≈3.41⋅10<sup>5</sup> pixels). As a control measure for this requirement the minimal *wrGOF* value was checked for all trajectories (i.e., the intensity profiles at pixel location *x,y*) within the brain section. In total we performed 1225 (=352 ) *cICAP* calculations while varying both the confidence value *c* and the tolerance value *t* accordingly to the range from 0.0 to 0.7 and 10-3 to 10-<sup>7</sup> , respectively. Figure 3a shows the results of this benchmark run. The largest improvement (cf. blue dot in Figure 3a) was found for the confidence and tolerance values of *c*=0.16 and *t* =1.07⋅10-<sup>6</sup> , respectively. Once both parameters *c* and *t* are fixed the stopping criterion *ε*, which controls the number of iterations, can be determined. For this task, we investigated the smallest and largest error found in the base functions of **A** ^ by means of the *MSE* value as expressed by Equation (9). As a reference for the definition of signal components, the user independent component selection tool was used, which is based on the statistically analysis of the base functions as reported in [23]. Throughout all iterations in *cICAP* the maximal *MSE* across the columns in matrix **A** ^ that were identified as signal components was found to be well below 0.01. Figure 3b shows that after about 35 iterations the *MSE* value converges to almost a constant for both the signal (red) and noise components (blue), where both types of components can clearly be identified by means of the corresponding *MSE* value.

In Figure 4 the progress of changes in the sinusoidal profile of one exemplary basis vector is shown for all 18 iterations (blue) within *cICAP*. The algorithm converges with a very small fit error (*MSE*) of 1.07⋅10-<sup>6</sup> . The resulting sinusoidal trajectory (black) almost perfectly fits the theoretically expectation function (red).

The process of updating the weight matrix and how *a priori* knowledge is included in *cICAP* is shown in Figure 5. Once a basis vector resembles the expectation function under the criteria of *MSE*(*a* ^ *i \** ) < *t*, this basis vector is not changed in further iterations. Note, throughout the optimization the selection of the underlying signal components is automatically included in *cICAP* and no further user interaction is required.

**Figure 4.** Stepwise changes of one exemplary basis function (changing from light blue to black). After 18 iterations

**Figure 3.** Parameter benchmark run. In (a) the optimal confidence value *c* and the tolerance value *t* was determined through a test run of 1225 signal decomposition steps varying both parameters *c* and *t* 35 times each through the

ness-of-fit statistics (wrGOF) is increased at all pixel locations, the largest minimal *wrGOF* value (blue circle) is deter‐

(b) The development of the maximal mean squared error (*MSE*) is shown for 147 iterations to determine the stopping criterion ε, which controls the number of iterations. The maximal and minimal *MSE* value is plotted for components reflecting signal of interest (red) and noise components (blue), respectively. After a few iterations the largest *MSE* val‐ ue for signal components was found to be well below 1%. Therefore, setting ε to 0.01 both types of components can

, respectively. As a control measure to find the best parameters, where the good‐

with a minimal wrGOF value of about 21.

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

365

. The resulting sinusoidal trajectory

using *cICAP* the algorithm converges with a very small fit error (*MSE*) of 1.07⋅10-6

(black) almost perfectly fits the theoretically expectation function (red).

range from 0.0 to 0.7 and 10-3

clearly be identified.

to 10-7

mined. The best set of parameters was found for *c=0.16* and *t* =1.07⋅10-6

ly the expectation function *f ICA*(*uρ*) would differ in terms of its amplitude, while the shape of the waveform may be still sinusoidal. The assumption here is that the signal power (across the rotation angles) of the signal of interest is larger compared to noise. The restriction of ω ≥ 1 is needed in order to not artificially improve the goodness-of-fit in cases where the two expect‐ ation functions *f raw*(*uρ*) and *f ICA*(*uρ*) are very similar (e.g., at pixel locations with low noise).

The optimal parameters *c*, *t* and *ε* for the weight update as described above (cf. Equation (10)), should be determined using a representative but independent data set. In the test data from the same post-mortem brain were used, but from a different microtome section which is not included in the subsequent analysis. The determination of the optimal parameters for *cICAP* requires the maximization of the weighted goodness-of-fit (*wrGOF*) values across all

*wrGOF* value was checked for all trajectories (i.e., the intensity profiles at pixel location *x,y*)

both the confidence value *c* and the tolerance value *t* accordingly to the range from 0.0 to 0.7

improvement (cf. blue dot in Figure 3a) was found for the confidence and tolerance values of

criterion *ε*, which controls the number of iterations, can be determined. For this task, we

*MSE* value as expressed by Equation (9). As a reference for the definition of signal components, the user independent component selection tool was used, which is based on the statistically analysis of the base functions as reported in [23]. Throughout all iterations in *cICAP* the

found to be well below 0.01. Figure 3b shows that after about 35 iterations the *MSE* value converges to almost a constant for both the signal (red) and noise components (blue), where both types of components can clearly be identified by means of the corresponding *MSE* value.

In Figure 4 the progress of changes in the sinusoidal profile of one exemplary basis vector is shown for all 18 iterations (blue) within *cICAP*. The algorithm converges with a very small fit

The process of updating the weight matrix and how *a priori* knowledge is included in *cICAP* is shown in Figure 5. Once a basis vector resembles the expectation function under the criteria

optimization the selection of the underlying signal components is automatically included in

) < *t*, this basis vector is not changed in further iterations. Note, throughout the

investigated the smallest and largest error found in the base functions of **A**

, respectively. Figure 3a shows the results of this benchmark run. The largest

, respectively. Once both parameters *c* and *t* are fixed the stopping

. The resulting sinusoidal trajectory (black) almost perfectly fits the

pixels). As a control measure for this requirement the minimal

) *cICAP* calculations while varying

^ that were identified as signal components was

^ by means of the

*2.2.2. Finding the optimal parameters*

pixel locations (≈3.41⋅10<sup>5</sup>

to 10-<sup>7</sup>

*c*=0.16 and *t* =1.07⋅10-<sup>6</sup>

error (*MSE*) of 1.07⋅10-<sup>6</sup>

of *MSE*(*a* ^ *i \**

theoretically expectation function (red).

*cICAP* and no further user interaction is required.

and 10-3

within the brain section. In total we performed 1225 (=352

364 Functional Brain Mapping and the Endeavor to Understand the Working Brain

maximal *MSE* across the columns in matrix **A**

**Figure 3.** Parameter benchmark run. In (a) the optimal confidence value *c* and the tolerance value *t* was determined through a test run of 1225 signal decomposition steps varying both parameters *c* and *t* 35 times each through the range from 0.0 to 0.7 and 10-3 to 10-7 , respectively. As a control measure to find the best parameters, where the good‐ ness-of-fit statistics (wrGOF) is increased at all pixel locations, the largest minimal *wrGOF* value (blue circle) is deter‐ mined. The best set of parameters was found for *c=0.16* and *t* =1.07⋅10-6 with a minimal wrGOF value of about 21. (b) The development of the maximal mean squared error (*MSE*) is shown for 147 iterations to determine the stopping criterion ε, which controls the number of iterations. The maximal and minimal *MSE* value is plotted for components reflecting signal of interest (red) and noise components (blue), respectively. After a few iterations the largest *MSE* val‐ ue for signal components was found to be well below 1%. Therefore, setting ε to 0.01 both types of components can clearly be identified.

**Figure 4.** Stepwise changes of one exemplary basis function (changing from light blue to black). After 18 iterations using *cICAP* the algorithm converges with a very small fit error (*MSE*) of 1.07⋅10-6 . The resulting sinusoidal trajectory (black) almost perfectly fits the theoretically expectation function (red).

**3. Signal enhancement in 3D-PLI data sets**

**3.2. Performance on signal decomposition using** *cICAP*

old for this data set was adapted to 2.2.

across the 102 brain sections was found to be (3.32 ± 0.95) ⋅ 10<sup>5</sup>

The performance of *cICAP* was tested on 102 histological sections of one adult human brain acquired from the body donor program of the University of Rostock, Germany, in accordance with legal requirements. The brain was fixed in 4% formalin for 6 months, cryo-protected with glycerin and cut after freezing in a cryostat microtome (Polycut CM 3500, Leica, Germany). The brain was sectioned coronally with slice a thickness of 70 µm. The sections were mounted on glass slides with Aquatex© (Merck, Germany) and cover-slipped. Note that all preparation

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

367

The sections were digitized using the polarimeter setup described in section 1.2. Each brain section was imaged at 18 equidistant rotation angles of the polarimeter covering an angle range between 0° and 170° (Figure 1b). The acquired RGB-colored images have a size of 2776 × 2080 pixels with a pixel size of 64 µm × 64 µm. The intensities were sampled with a dynamic range of 14 bits per color channel. Since the light source of the polarimeter is composed of lightemitting diodes (LED) emitting a narrow-band green wavelength spectrum (central wave‐ length of 525 nm), only the green channel from the RGB color triplet was used for further

Before the images were decomposed by ICA, a standardized calibration technique using flat fields was applied, in order to compensate pixel-wise for inhomogeneities across the field of view [22]. In addition, the separation of brain tissue from the non-relevant image background was done by means of the interactive learning and segmentation toolkit (ilastik) [44]. Thus, the following calculations and statistical analyses were solely done on basis of brain tissue

To test the performance of the signal decomposition and signal restoration *cICAP* was applied to the data set described in section 3.1. The mean number of pixels used for ICA

decomposition and enhancement gained by *cICAP* were compared to the automatic component selection routine, which were recently published in [23], whereas the thresh‐

After the determination of the parameters *c*, *t* and *ε* all values were fixed within *cICAP* and both ICA algorithms were applied to the same set of PLI data for comparison. As shown in Table 1, the identified number of signal components (on average) does not vary dramatically across both algorithms. However, the smallest variability with the fewest number of iterations needed was found in the results of *cICAP*. The number of iterations needed (on average) in *cICAP* was 192, while for the Infomax algorithm using the automatic selection routine it was found to be 406. More importantly, the goodness-of-fit criterion *wrGOF* was found to be largest using *cICAP*. It is important to note that this criterion expresses the signal enhancement (or degradation) after ICA filtering and was calculated for all trajectories (i.e., at each pixel

. The results of signal

steps are conforming to preserve the integrity of the birefringent myelin sheaths.

**3.1. Signal acquisition and preprocessing**

analysis.

measurements.

**Figure 5.** Scheme of the weight matrix update and optimization procedure in *cICAP*.

### **3. Signal enhancement in 3D-PLI data sets**

### **3.1. Signal acquisition and preprocessing**

The performance of *cICAP* was tested on 102 histological sections of one adult human brain acquired from the body donor program of the University of Rostock, Germany, in accordance with legal requirements. The brain was fixed in 4% formalin for 6 months, cryo-protected with glycerin and cut after freezing in a cryostat microtome (Polycut CM 3500, Leica, Germany). The brain was sectioned coronally with slice a thickness of 70 µm. The sections were mounted on glass slides with Aquatex© (Merck, Germany) and cover-slipped. Note that all preparation steps are conforming to preserve the integrity of the birefringent myelin sheaths.

The sections were digitized using the polarimeter setup described in section 1.2. Each brain section was imaged at 18 equidistant rotation angles of the polarimeter covering an angle range between 0° and 170° (Figure 1b). The acquired RGB-colored images have a size of 2776 × 2080 pixels with a pixel size of 64 µm × 64 µm. The intensities were sampled with a dynamic range of 14 bits per color channel. Since the light source of the polarimeter is composed of lightemitting diodes (LED) emitting a narrow-band green wavelength spectrum (central wave‐ length of 525 nm), only the green channel from the RGB color triplet was used for further analysis.

Before the images were decomposed by ICA, a standardized calibration technique using flat fields was applied, in order to compensate pixel-wise for inhomogeneities across the field of view [22]. In addition, the separation of brain tissue from the non-relevant image background was done by means of the interactive learning and segmentation toolkit (ilastik) [44]. Thus, the following calculations and statistical analyses were solely done on basis of brain tissue measurements.

#### **3.2. Performance on signal decomposition using** *cICAP*

**Figure 5.** Scheme of the weight matrix update and optimization procedure in *cICAP*.

366 Functional Brain Mapping and the Endeavor to Understand the Working Brain

To test the performance of the signal decomposition and signal restoration *cICAP* was applied to the data set described in section 3.1. The mean number of pixels used for ICA across the 102 brain sections was found to be (3.32 ± 0.95) ⋅ 10<sup>5</sup> . The results of signal decomposition and enhancement gained by *cICAP* were compared to the automatic component selection routine, which were recently published in [23], whereas the thresh‐ old for this data set was adapted to 2.2.

After the determination of the parameters *c*, *t* and *ε* all values were fixed within *cICAP* and both ICA algorithms were applied to the same set of PLI data for comparison. As shown in Table 1, the identified number of signal components (on average) does not vary dramatically across both algorithms. However, the smallest variability with the fewest number of iterations needed was found in the results of *cICAP*. The number of iterations needed (on average) in *cICAP* was 192, while for the Infomax algorithm using the automatic selection routine it was found to be 406. More importantly, the goodness-of-fit criterion *wrGOF* was found to be largest using *cICAP*. It is important to note that this criterion expresses the signal enhancement (or degradation) after ICA filtering and was calculated for all trajectories (i.e., at each pixel location). If the *wrGOF* value is larger than 1, then noise and artifacts were removed (or suppressed) successfully. In cases, where the signal decomposition fails (or is not optimal) the *wrGOF* value becomes even smaller than 1 (Figure 6). To check the decomposition results for signal degradation after the application of ICA, the minimal *wrGOF* value was calculated. As illustrated in Table 1 a *wrGOF* value smaller than 1 was found in 3.57% and 0.35% of all pixels in the results from Infomax and *cICAP*, respectively.

demonstrated, where one component of interest was manually excluded before the 3D-PLI signal was reconstructed. For this section, 24% of all pixels have a *wrGOF* value below 1.

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

369

**Figure 6.** χ 2 ratios and pixel locations illustrating goodness-of-fit values. In a) χ <sup>2</sup> ratios are shown using all com‐ ponents identified by *cICAP*, where the goodness-of-fit statistics was larger than 1 at all pixel locations. In the coronal brain section pixels highlighted in red show locations with goodness-of-fit values larger than 300. b)To demonstrate the effect of being sensitive to missing signal of interest, one identified component was manually deselected. As a result, large areas of white matter structure do show a *wrGOF* value of less than 1 (highlighted

In recent years 3D-polarized light imaging (3D-PLI) has been shown to provide new insights into the organization of the human brain including mapping of the fiber anatomy [15,45,46]. Through advances in the experimental setup of the employed polarimeter, as well as in signal processing, this modality provides unique data sets to explore the 3D fiber architecture in the human brain at a submillimeter resolution [16,47,48]. Though the technique as described here is applicable solely to postmortem brain tissue, the comprehensive description of complex fiber orientations in distinct brain regions (e.g., prevalent fiber crossings) can be used to guide and evaluate fiber tractography algorithms based on diffusion MRI. By this means the fiber

in blue).

**4. Conclusion**

In general both ICA algorithms produced remarkable decomposition results, where signal enhancement (*wrGOF* ≥ 1) is evident in almost all pixels. However, we found a much better signal decomposition performance in the results of *cICAP*. This is expressed by the fact that the *wrGOF* value was found to be equal or larger than 10 in about 80% of all trajectories, while for the Infomax algorithm it was 39%. For *cICAP* we still found a goodness-of-fit value of larger or equal than 100 in 56% of all pixels. In comparison to the standard Infomax this number drops down to about 15%.


**Table 1.** Statistical analysis of ICA decomposition results. The signal components for the ICA routine *I* were selected by means of a user independent algorithm (see text).

In Figure 6a pixel locations of a representative section is shown, where the test statistic (*wrGOF*) is larger or equal 300. In this example signal enhancement was evident at all pixel locations (i.e., *wrGOF* > 1; not shown). By focusing on larger values of the test statistic (e.g., *wrGOF* ≥ 300) the figure demonstrates that the enhancement of the optical signals is ensured in particular in the region of the neo-cortex and at the white to gray matter boundary (Figure 6a). This in general is challenging, since in such regions the 3D-PLI signal has low amplitudes and therefore the measured signal has a poor signal-to-noise ratio.

With the introduction of the weighting factor (or penalty factor) as expressed by Equation (12), the test statistic *wrGOF* is sensitive to both, changes in the signal-to-noise ratio as well as in reductions of signal strength after ICA. This in particular will be the case when one or more components of interest are not selected for signal reconstruction (e.g., when the selection is performed manually). As a result, at corresponding pixel locations the *wrGOF* value will be less than 1. The idea behind this is, that the weighted test statistic *wrGOF* accounts for missing gray and white matter components and improvements in SNR, independently to the metric used for identification of components of interest. In Figure 6b such a testing scenario is demonstrated, where one component of interest was manually excluded before the 3D-PLI signal was reconstructed. For this section, 24% of all pixels have a *wrGOF* value below 1.

**Figure 6.** χ 2 ratios and pixel locations illustrating goodness-of-fit values. In a) χ <sup>2</sup> ratios are shown using all com‐ ponents identified by *cICAP*, where the goodness-of-fit statistics was larger than 1 at all pixel locations. In the coronal brain section pixels highlighted in red show locations with goodness-of-fit values larger than 300. b)To demonstrate the effect of being sensitive to missing signal of interest, one identified component was manually deselected. As a result, large areas of white matter structure do show a *wrGOF* value of less than 1 (highlighted in blue).

### **4. Conclusion**

location). If the *wrGOF* value is larger than 1, then noise and artifacts were removed (or suppressed) successfully. In cases, where the signal decomposition fails (or is not optimal) the *wrGOF* value becomes even smaller than 1 (Figure 6). To check the decomposition results for signal degradation after the application of ICA, the minimal *wrGOF* value was calculated. As illustrated in Table 1 a *wrGOF* value smaller than 1 was found in 3.57% and 0.35% of all pixels

In general both ICA algorithms produced remarkable decomposition results, where signal enhancement (*wrGOF* ≥ 1) is evident in almost all pixels. However, we found a much better signal decomposition performance in the results of *cICAP*. This is expressed by the fact that the *wrGOF* value was found to be equal or larger than 10 in about 80% of all trajectories, while for the Infomax algorithm it was 39%. For *cICAP* we still found a goodness-of-fit value of larger or equal than 100 in 56% of all pixels. In comparison to the standard Infomax this number

> *identified ICs* (mean) 7.7 ± 4.75 3.8 ± 1.19 *number of iterations (mean)* 405.5± 32.3 191.8± 23.13 *wrGOF* (median) 18.8 1910.5 wrGOF < 1 (%) 3.57 % 0.35 % wrGOF ≥ 10 (%) 38.8 % 80.0 % wrGOF ≥ 100 (%) 15.1 % 55.6 %

**Table 1.** Statistical analysis of ICA decomposition results. The signal components for the ICA routine *I* were selected by

In Figure 6a pixel locations of a representative section is shown, where the test statistic (*wrGOF*) is larger or equal 300. In this example signal enhancement was evident at all pixel locations (i.e., *wrGOF* > 1; not shown). By focusing on larger values of the test statistic (e.g., *wrGOF* ≥ 300) the figure demonstrates that the enhancement of the optical signals is ensured in particular in the region of the neo-cortex and at the white to gray matter boundary (Figure 6a). This in general is challenging, since in such regions the 3D-PLI signal has low amplitudes

With the introduction of the weighting factor (or penalty factor) as expressed by Equation (12), the test statistic *wrGOF* is sensitive to both, changes in the signal-to-noise ratio as well as in reductions of signal strength after ICA. This in particular will be the case when one or more components of interest are not selected for signal reconstruction (e.g., when the selection is performed manually). As a result, at corresponding pixel locations the *wrGOF* value will be less than 1. The idea behind this is, that the weighted test statistic *wrGOF* accounts for missing gray and white matter components and improvements in SNR, independently to the metric used for identification of components of interest. In Figure 6b such a testing scenario is

and therefore the measured signal has a poor signal-to-noise ratio.

**I) Infomax II) cICAP**

in the results from Infomax and *cICAP*, respectively.

368 Functional Brain Mapping and the Endeavor to Understand the Working Brain

drops down to about 15%.

means of a user independent algorithm (see text).

In recent years 3D-polarized light imaging (3D-PLI) has been shown to provide new insights into the organization of the human brain including mapping of the fiber anatomy [15,45,46]. Through advances in the experimental setup of the employed polarimeter, as well as in signal processing, this modality provides unique data sets to explore the 3D fiber architecture in the human brain at a submillimeter resolution [16,47,48]. Though the technique as described here is applicable solely to postmortem brain tissue, the comprehensive description of complex fiber orientations in distinct brain regions (e.g., prevalent fiber crossings) can be used to guide and evaluate fiber tractography algorithms based on diffusion MRI. By this means the fiber orientation maps provided by 3D-PLI might help to optimize the reliability of in-vivo diffusion MRI results. Precise information about the local individual fiber architecture of a patient is of particular interest in case of planning and performing a neurosurgical intervention, for instance.

**Author details**

Jürgen Dammers1

**References**

trum Jülich, Jülich, Germany

Sep;, 1(4), 245-51.

15;, 29(4), 1092-105.

human subject. Lancet. (1901).

(1967). Apr;, 4(4), 369-74.

, Lukas Breuer1

Brain, INM-1, Forschungszentrum Jülich, Jülich, Germany

the brain. Oxford University Press; (2006).

New York: Springer; (2007). page , 149-1167.

Oxford University Press; (2011). page 624.

In-vivo Neuroanatomy. Academic Press; (2009). page 490.

, Giuseppe Tabbì2

1 Institute of Neuroscience and Medicine, Medical Imaging Physics, INM-4, Forschungszen‐

2 Institute of Neuroscience and Medicine, Structural and Functional Organisation of the

[1] Mesulam, M. Foreword. In: Schmahmann JD, Pandya DN, editors. Fiber pathways of

[2] Kötter, R. Handbook of brain connectivity. In: Jirsa VK, McIntosh AR, editors. Berlin

[3] Sporns, O, Tononi, G, & Kötter, R. The human connectome: A structural description of the human brain. PLoS computational biology. Public Library of Science; (2005).

[4] Johansen-berg, H. Behrens TEJ. Diffusion MRI: From Quantitative Measurement to

[5] Jones, D. K. Diffusion MRI: Theory, Methods, and Applications. Jones DK, editor.

[6] Türe, U, Yasargil, M. G, Friedman, A. H, & Al-mefty, O. Fiber dissection technique:

[7] Klingler, J. Erleichterung der makroskopischen Präparation des Gehirns durch den

[8] Bürgel, U, Amunts, K, Hoemke, L, Mohlberg, H, Gilsbach, J. M, & Zilles, K. White matter fiber tracts of the human brain: three-dimensional mapping at microscopic resolution, topography and intersubject variability. NeuroImage. Department of Neurosurgery, RWTH Aachen University, D-52074 Aachen, Germany.; (2006). Feb

[9] Flechsig, P. Developmental (myelongenetic) localisation of the cerebral cortex in the

[10] Fink, R. P, & Heimer, L. Two methods for selective silver impregnation of degenerat‐ ing axons and their synaptic endings in the central nervous system. Brain Research.

lateral aspect of the brain. Neurosurgery. (2000). Aug;, 47(2), 417-27.

Gefrierprozess. Schweizer Archiv für Neurologie und Psychiatrie. (1935).

and Markus Axer2

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

371

However, in 3D-PLI the reconstruction of nerve fiber pathways strongly depends on the quality of the measured intensities represented by the so-called light intensity profiles or 3D-PLI signals, respectively. Hence, advanced signal processing tools are required to enable the precise determination of locally prevailing fiber orientations in form of unit vectors defined by the in-section direction angle *φ* and the out-of-section inclination angle *α*.

Independent component analysis (ICA) turned out to improve 3D-PLI signals significantly [22,23]. It was shown that ICA is capable of restoring the original birefringent signal by effectively removing noise and artifact components in the measured data. In addition, measures for the qualitative and quantitative evaluation of the 3D-PLI signals before and after the ICA filtering were introduced. In particular, the signal enhancement after ICA based denoising is large at the white to gray matter boundary, where the 3D-PLI signal is weak due to decreasing fiber density when approaching gray matter domains [22].

In *cICAP*, the "constrained ICA for polarized light imaging", the signal decomposition is optimized using a weight update, which incorporates the prior knowledge that the expected signal in 3D-PLI can be modeled utilizing the Jones calculus (cf. Equation (1)). This prior knowledge gives rise to a training of the unmixing procedure in *cICAP*, where the weight matrix update is optimal with respect to the sources of interest (cf. Equation (10)). Utilizing this approach the source separation and identification in *cICAP* is carried out automatically in one routine, where no further component selection tool is needed for subsequent extraction of the signal of interest.

With the introduction of *cICAP* an ICA-based source separation method is available which is optimal for the extraction of the sinusoidal trajectory along different rotation angles in a set of spatially mixed 3D-PLI images. The better quality of the source separation in *cICAP* is reflected by the increased goodness-of-fit statistic (*wrGOF* ), which not only takes the improvement in the signal-to-noise ratio (SNR) into account, but also accounts for signal restoration [23].

The dedicated signal processing tool *cICAP* largely contributes to the reliability of 3D-PLI signals and, hence, the accuracy of subsequent fiber tract reconstruction. As a consequence, *cICAP* plays a key role in the 3D-PLI processing chain and helps to establish the so far missing link between the microscopic and the macroscopic characterization of the anatomical connec‐ tivity of the human brain.

### **Acknowledgements**

We thank Prof. Dr. A. Wree, Institute for Anatomy, University of Rostock, for providing us with the postmortem human brain.

### **Author details**

orientation maps provided by 3D-PLI might help to optimize the reliability of in-vivo diffusion MRI results. Precise information about the local individual fiber architecture of a patient is of particular interest in case of planning and performing a neurosurgical intervention, for

However, in 3D-PLI the reconstruction of nerve fiber pathways strongly depends on the quality of the measured intensities represented by the so-called light intensity profiles or 3D-PLI signals, respectively. Hence, advanced signal processing tools are required to enable the precise determination of locally prevailing fiber orientations in form of unit vectors defined

Independent component analysis (ICA) turned out to improve 3D-PLI signals significantly [22,23]. It was shown that ICA is capable of restoring the original birefringent signal by effectively removing noise and artifact components in the measured data. In addition, measures for the qualitative and quantitative evaluation of the 3D-PLI signals before and after the ICA filtering were introduced. In particular, the signal enhancement after ICA based denoising is large at the white to gray matter boundary, where the 3D-PLI signal is weak due

In *cICAP*, the "constrained ICA for polarized light imaging", the signal decomposition is optimized using a weight update, which incorporates the prior knowledge that the expected signal in 3D-PLI can be modeled utilizing the Jones calculus (cf. Equation (1)). This prior knowledge gives rise to a training of the unmixing procedure in *cICAP*, where the weight matrix update is optimal with respect to the sources of interest (cf. Equation (10)). Utilizing this approach the source separation and identification in *cICAP* is carried out automatically in one routine, where no further component selection tool is needed for subsequent extraction of

With the introduction of *cICAP* an ICA-based source separation method is available which is optimal for the extraction of the sinusoidal trajectory along different rotation angles in a set of spatially mixed 3D-PLI images. The better quality of the source separation in *cICAP* is reflected by the increased goodness-of-fit statistic (*wrGOF* ), which not only takes the improvement in the signal-to-noise ratio (SNR) into account, but also accounts for signal restoration [23].

The dedicated signal processing tool *cICAP* largely contributes to the reliability of 3D-PLI signals and, hence, the accuracy of subsequent fiber tract reconstruction. As a consequence, *cICAP* plays a key role in the 3D-PLI processing chain and helps to establish the so far missing link between the microscopic and the macroscopic characterization of the anatomical connec‐

We thank Prof. Dr. A. Wree, Institute for Anatomy, University of Rostock, for providing us

by the in-section direction angle *φ* and the out-of-section inclination angle *α*.

370 Functional Brain Mapping and the Endeavor to Understand the Working Brain

to decreasing fiber density when approaching gray matter domains [22].

instance.

the signal of interest.

tivity of the human brain.

**Acknowledgements**

with the postmortem human brain.

Jürgen Dammers1 , Lukas Breuer1 , Giuseppe Tabbì2 and Markus Axer2

1 Institute of Neuroscience and Medicine, Medical Imaging Physics, INM-4, Forschungszen‐ trum Jülich, Jülich, Germany

2 Institute of Neuroscience and Medicine, Structural and Functional Organisation of the Brain, INM-1, Forschungszentrum Jülich, Jülich, Germany

### **References**


[11] Clarke, S, & Miklossy, J. Occipital cortex in man: organization of callosal connections, related myelo- and cytoarchitecture, and putative boundaries of functional visual areas. J. Comp. Neurol. (1990).

[25] Cong, F, Kalyakin, I, Ahuttunen-scott, T, Li, H, Lyytinen, H, & Ristaniemi, T. Singletrial based Independent Component Analysis on Mismatch Negativity in Children.

Optimized Signal Separation for 3D-Polarized Light Imaging

http://dx.doi.org/10.5772/55246

373

[26] Esposito, F, Formisano, E, Seifritz, E, Goebel, R, Morrone, R, Tedeschi, G, et al. Spa‐ tial independent component analysis of functional MRI time-series: to what extent do results depend on the algorithm used? Human brain mapping. (2002). Jul;, 16(3),

[27] Hild, K. E, & Nagarajan, S. S. Source localization of EEG/MEG data by correlating columns of ICA and lead field matrices. IEEE transactions on bio-medical engineer‐

[28] Bartlett, M. S. Information maximization in face processing. Neurocomputing. (2007).

[29] Dammers, J, Schiek, M, Boers, F, Silex, C, Zvyagintsev, M, Pietrzyk, U, et al. Integra‐ tion of amplitude and phase statistics for complete artifact removal in independent components of neuromagnetic recordings. IEEE Trans. Biomed. Eng. (2008). Oct;,

[30] Escudero, J, Hornero, R, Abásolo, D, & Fernández, A. Quantitative evaluation of arti‐ fact removal in real magnetoencephalogram signals with blind source separation. Annals of biomedical engineering. Springer Netherlands; (2011). Aug 21;, 39(8),

[31] Li, Y, Ma, Z, Lu, W, & Li, Y. Automatic removal of the eye blink artifact from EEG using an ICA-based template matching approach. Physiological measurement.

[32] Mantini, D, Franciotti, R, Romani, G. L, & Pizzella, V. Improving MEG source locali‐ zations: an automated method for complete artifact removal based on independent

[33] Nikulin, V V, Nolte, G, & Curio, G. A novel method for reliable and fast extraction of neuronal EEG/MEG oscillations on the basis of spatio-spectral decomposition. Neu‐

[34] Huang, D-S, & Mi, J-X. A New Constrained Independent Component Analysis Meth‐

[35] Hesse, C. W. On estimating the signal subspace dimension of high-density multi‐ channel magnetoencephalogram measurements. Conference proceedings : Annual International Conference of the IEEE Engineering in Medicine and Biology Society. IEEE Engineering in Medicine and Biology Society. Lyon, France. (2007). Jan;,

od. IEEE Transactions on Neural Networks. (2007). Sep;, 18(5), 1532-5.

component analysis. NeuroImage. (2008). Mar;, 40(1), 160-73.

International Journal of Neural Systems. (2010).

ing. (2009). Nov;, 56(11), 2619-26.

Aug;70(13-15):2204-17.

(2006). Apr;, 27(4), 425-36.

roImage. (2011). Apr 15;, 55(4), 1528-35.

55(10), 2353-62.

2274-86.

6228-31.

146-57.


[25] Cong, F, Kalyakin, I, Ahuttunen-scott, T, Li, H, Lyytinen, H, & Ristaniemi, T. Singletrial based Independent Component Analysis on Mismatch Negativity in Children. International Journal of Neural Systems. (2010).

[11] Clarke, S, & Miklossy, J. Occipital cortex in man: organization of callosal connections, related myelo- and cytoarchitecture, and putative boundaries of functional visual

[12] Burkhalter, A, Bernardo, K. L, & Charles, V. Development of local circuits in human visual cortex. The Journal of neuroscience- the official journal of the Society for Neu‐

[13] Lanciego, J, & Wouterlood, F. Neuroanatomical tract-tracing methods beyond 2000: what's now and next. Journal of Neuroscience Methods. (2000). Nov;, 103(1), 1-2. [14] Schmahmann, J, & Pandya, D. Fiber Pathways of the Brain. Oxford University Press;

[15] Axer, M, Grässel, D, Kleiner, M, Dammers, J, Dickscheid, T, Reckfort, J, et al. Highresolution fiber tract reconstruction in the human brain by means of three-dimen‐ sional polarized light imaging. Frontiers in neuroinformatics. (2011). Jan;, 5(34), 1-13.

[16] Axer, M, Amunts, K, Grässel, D, Palm, C, Dammers, J, Axer, H, et al. A novel ap‐ proach to the human connectome: ultra-high resolution mapping of fiber tracts in the

[17] Axer, H, Beck, S, Axer, M, Schuchardt, F, Heepe, J, Flücken, A, et al. Microstructural analysis of human white matter architecture using polarized light imaging: views

[18] Schmitt, F. O, & Bear, R. S. The optical properties of vertebrate nerve axons as related to fiber size. Journal of Cellular and Comparative Physiology. (1937). Feb;, 9(2),

[19] Schmidt, W. J. Zur Doppelbrechung des Nervenmarks. Zeitschrift für wissenschaftli‐

[20] Göthlin, G. F. Die doppelbrechenden Eigenschaften des Nervengewebes. Kungliga

[21] Jones, R. C. A New Calculus for the Treatment of Optical Systems. Journal of the Op‐

[22] Dammers, J, Axer, M, Grässel, D, Palm, C, Zilles, K, Amunts, K, et al. Signal enhance‐ ment in polarized light imaging by means of independent component analysis. Neu‐

[23] Dammers, J, Breuer, L, Axer, M, Kleiner, M, Eiben, B, Grässel, D, et al. Automatic identification of gray and white matter components in polarized light imaging. Neu‐

[24] Beckmann, C. F, & Smith, S. M. Probabilistic independent component analysis for functional magnetic resonance imaging. IEEE Transactions on Medical Imaging.

from neuroanatomy. Frontiers in neuroinformatics. (2011). Jan;, 5(28), 1-12.

areas. J. Comp. Neurol. (1990).

(2010). page 672.

261-73.

roscience. (1993). May;, 13(5), 1916-31.

372 Functional Brain Mapping and the Endeavor to Understand the Working Brain

brain. NeuroImage. (2011). Jan;, 54(2), 1091-101.

Svenska Vetenskapsakademiens. (1913). Handlingar(51).

tical Society of America. (1941). Jul;, 31(7), 500-3.

che Mikroskopie. (1923). , 41, 29-38.

roImage. (2010). Jan;, 49(2), 1241-8.

roImage. (2012). Jan 16;, 59(2), 1338-47.

(2004). Feb;, 23(2), 137-52.


[36] Liu, J, Ghassemi, M. M, Michael, A. M, Boutte, D, Wells, W, Perrone-bizzozero, N, et al. An ICA with reference approach in identification of genetic variation and associ‐ ated brain networks. Frontiers in human neuroscience. (2012). Jan;, 6(21), 1-10.

**Chapter 19**

**Causal Relationships and**

**Effective Brain Connectivity**

Guiomar Niso, Ernesto Pereda, Fernando Maestú, María Gudín, Sira Carrasco and Francisco del-Pozo

The human brain comprises approximately 1011 neurons, each of which makes about 103 synaptic connections. This huge number of connections between individual processing elements provides the fundamental scaffold for neuronal ensembles to become transiently synchronized or functionally connected [1]. A similar complex network configuration and dynamics can also be found at the macroscopic scales of systems neuroscience and brain imaging [2]. The emergence of dynamically coupled cell assemblies represents the neurophy‐ siological substrate for cognitive function such as perception, learning, thinking [3]. Under‐ standing the complex network organization of the brain on the basis of neuroimaging data

Several measures to evaluate at various scales (single cells, cortical columns, or brain areas) how the different parts of the brain communicate have been recently proposed. We can classify them, according to their symmetry, into two groups: symmetric and asymmetric measures. Symmetric measures, such as correlation, coherence, phase synchronization indexes (PLV, PLI), evaluate *functional connectivity,* i.e. statistically significant relationships between signals recorded from spatially remote neurophysiological events. On the other hand, the asymmetric ones are able to detect *effective connectivity* i.e. information flow (the influence one neuronal system exerts over another) [4], and therefore they reveal the direction of the interaction.

> © 2013 Niso et al.; licensee InTech. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use,

© 2013 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution,

distribution, and reproduction in any medium, provided the original work is properly cited.

and reproduction in any medium, provided the original work is properly cited.

represents one of the most impervious challenges for systems neuroscience.

Granger Causality (GC) belongs to this latter group.

Additional information is available at the end of the chapter

**Network Parameters in**

http://dx.doi.org/10.5772/56499

**1. Introduction**

