1. Introduction

Phylogeny is the branch of biology that seeks to reconstruct evolutionary family trees. Such reconstruction can take place at various scales. For example, we could attempt to build the family tree for various present day indigenous populations in the Americas and Asia in order to glean information about the possible course of migration of humans into the Americas. At the level of species, we could seek to determine whether modern humans are more closely related to chimpanzees or to gorillas. Ultimately, we would like to be able to reconstruct the entire “tree of life” that describes the course of evolution leading to all present day species. Because the status of the “leaves” on which we wish to build a tree differs from instance to instance, biologists use the general term taxa (singular taxon) for the leaves in a general phylogenetic problem.

For example, for 4 taxa, we might seek to decide whether the tree describes the course of evolution. In such trees:

• • the arrow of time is down the page,

• • paths down through the tree represent lineages (lines of descent),

• • any point on a lineage corresponds to a point of time in the life of some ancestor of a taxon,

• • vertices other than leaves represent times at which lineages diverge,

• • the root corresponds to the most recent common ancestor of all the taxa.

In many ways, the late 1950s marked the beginning of the digital age, and with it, the beginning of a new age for the mathematics of signal processing. Highspeed analog-to-digital converters had just been invented. These devices were capable of taking analog signals like time series (think of continuous functions of time like seismograms which measure the seismic activity — the amount of bouncing — at a fixed location, or an EEG, or an EKG) and converting them to lists of numbers. These numbers were obtained by sampling the time series, that is, recording the value of the function at regular intervals, which at that time could be as fast as 300,000 times every second. (Current technology permits sampling at much higher rates where necessary.) Suddenly, reams and reams of data were being generated and new mathematics was needed for their analysis, manipulation and management.

So was born the discipline of Digital Signal Processing (DSP), and it is no exaggeration to say that the world has not been the same. In the mathematical sciences the DSP revolution has, among other things, helped drive the development of disciplines like algorithmic analysis (which was the impetus behind the creation of computer science departments), communication and information theory, linear algebra, computational statistics, combinatorics, and discrete mathematics. DSP tools have changed the face of the arts (electroacoustic music and image processing), health care (medical imaging and computed imaging), and, of course, both social and economic commerce (i.e., the internet). Suffice to say that the mathematics of DSP is one of the pillars supporting the amazing technological revolution that we are experiencing today.

1. Introduction

1.1. Background. To register two images means to align them, so that common features overlap and differences, should there be any, between the two are emphasized and- readily visible to the naked eye. We refer to the process of aligning two images as image registration.

There are a host of clinical applications requiring image registration. For example, one would like to compare two Computed Tomography (CT) scans of a patient, taken say six months ago and yesterday, and identify differences between the two, e.g., the growth of a tumor during the intervening six months (Figure 1). One could also want to align Positron Emission Tomography (PET) data to an MR image, so as to help identify the anatomic location of certain mental activation [43]. And one may want to register lung surfaces in chest Computed Tomography (CT) scans for lung cancer screening [7]. While all of these identifications can be done in the radiologist's head, the possibility always exists that small, but critical, features could be missed. Also, beyond identification itself, the extent of alignment required could provide important quantitative information, e.g., how much a tumor's volume has changed.

1. Introduction

Optical, or diffuse, tomography, refers to the use of low energy probes to obtain images of highly scattering media.

The main motivation for this line of work is, at present, the use of an infrared laser to obtain images of diagnostic value. There is a proposal to use this in neonatal clinics to measure oxygen content in the brains of premature babies as well as in the case of repeated mammography. With the discovery of highly specific markers that respond well in the optical or infrared region there are many potential applications of this emerging area; see [AI; A2].

Pure and Applied Mathematics: Two Sides of a Coin

The Bulletin of the AMS for November 1979 had a paper by L. Auslander and R. Tolimieri [3] with the delightful title “Is computing with the Finite Fourier Transform pure or applied mathematics?” This rhetorical question was answered by showing that in fact, the finite Fourier transform, and the family of efficient algorithms used to compute it, the Fast Fourier Transform (FFT), a pillar of the world of digital signal processing, were of interest to both pure and applied mathematicians.

Auslander had come of age as an applied mathematician at a time when pure and applied mathematicians still received much of the same training. The ends towards which these skills were then directed became a matter of taste.

1. Introduction

This paper is a sequel to our previous paper [3], where we considered the best sparsijying basis (BSB), and the least statistically-dependent basis (LSDB) for input data assumed to be realizations of a very simple stochastic process called the “spike process.” This process, which we will refer to as the “simple” spike process for convenience, puts a unit impulse (i.e., constant amplitude of 1) at a random location in a zero vector of length n.

1. Introduction

Data compression is a process that creates a compact data representation from a raw data source, usually with an end goal of facilitating storage or transmission. Broadly speaking, compression takes two forms, either loss less or lossy, depending on whether or not it is possible to reconstruct exactly the original datastream from its compressed version. For example, a data stream that consists of long runs of Os and Is (such as that generated by a black and white fax) would possibly benefit from simple run-length encoding, a lossless technique replacing the original datastream by a sequence of counts of the lengths of the alternating substrings of Os and Is. Lossless compression is necessary for situations in which changing a single bit can have catastrophic effects, such as in machine code of a computer program.

While it might seem as though we should always demand lossless compression, there are in fact many venues where exact reproduction is unnecessary. In particular, media compression, which we define to be the compression of image, audio, or video files, presents an excellent opportunity for lossy techniques. For example, not one among us would be able to distinguish between two images which differ in only one of the 229 bits in a typical 1024 x 1024 color image.

1. Introduction

The Fourier transform of a function on a compact Lie group computes the coefficients (Fourier coefficients) that enable its expression as a linear combination of the matrix elements from a complete set of irreducible representations of the group. In the case of abelian groups, especially the circle and its lower dimensional products (tori) this is precisely the expansion of a function on these domains in terms of complex exponentials. This representation is at the heart of classical signal and image processing (see [25; 26], for example).

The successes of abelian Fourier analysis are many, ranging from national defense to personal entertainment, from medicine to finance. The record of achievements is so impressive that it has perhaps sometimes led scientists astray, seducing them to look for ways to use these tools in situations where they are less than appropriate: for example, pretending that a sphere is a torus so as to avoid the use of spherical harmonics in favor of Fourier series - a favored mathematical hammer casting the multitudinous problems of science as a box of nails.

1. Introduction

An important activity, common to many fields of endeavor, is the act of refining high order information (detections of events, classification of objects, identification of activities, etc.) from large volumes of diverse data which is increasingly available through modern means of measurement, communication, and processing. This exploitation function winnows the available data concerning an object or situation in order to extract useful and actionable information, quite often through the application of techniques from statistical pattern recognition to the data. This may involve activities like detection, identification, and classification which are applied to the raw measured data, or possibly to partially processed information derived from it.

When new data are sought in order to obtain information about a specific situation, it is now increasingly common to have many different measurement degrees of freedom potentially available for the task.

1. Introduction and Statement of Results

The theory of matrix-valued spherical functions (see [GV; T]) gives a natural extension of the well-known theory for the scalar-valued case, see [He]. We start with a few remarks about the scalar-valued case.

The classical (scalar-valued) theory of spherical functions (put forward by Cartan and others after him) allows one to unify under one roof a number of examples that were very well known before the theory was formulated. These examples include many special functions like Jacobi polynomials, Bessel functions, Laguerre polynomials, Hermite polynomials, Legendre functions, etc.

1. Introduction

When Claude Shannon developed the mathematical theory of communication [1] he knew nothing about lasers and optical fibers. What he was mostly concerned with were communication channels using radio- and microwaves. Inherently, these channels have a narrower bandwidth than do optical fibers because of the lower carrier frequency (longer wavelength). More serious than this theoretical limitation are the practical bandwidth limitations imposed by weather and other environmental hazards. In contrast, optical fibers are a marvellously stable and predictable medium for transporting information and the influence of noise from the fiber itself can to a large degree be neglected. So, until recently there was no real need for any advanced signal processing in optical fiber communications systems. This has all changed over the last few years with the development of the internet.

Optical fiber communication became an economic reality in the early 1970s when absorption of less than 20 dB /km was achieved in optical fibers and lifetimes of more than 1 million hours for semiconductor lasers were accomplished.

1. The Impedance Matching Problem

Figure 1 shows a twin-whip HF (high-frequency) antenna mounted on a superstructure representative of a shipboard environment. If a signal generator is connected directly to this antenna, not all the power delivered to the antenna can be radiated by the antenna. If an impedance mismatch exists between the signal generator and the antenna, some of the signal power is reflected from the antenna back to the generator. To effectively use this antenna, a matching circuit must be inserted between the signal generator and antenna to minimize this wasted power.

Figure 2 shows the matching circuit connecting the generator to the antenna. Port 1 is the input from the generator. Port 2 is the output that feeds the antenna.

1. Introduction

In field after field, we are currently seeing new initiatives aimed at gathering large high-resolution three-dimensional datasets. While three-dimensional data have always been crucial to understanding the physical world we live in, this transition to ubiquitous 3-D data gathering seems novel. The driving force is undoubtedly the pervasive influence of increasing storage capacity and computer processing power, which affects our ability to create new 3-D measurement instruments, but which also makes it possible to analyze the massive volumes of data that inevitably result when 3-D data are being gathered.

Introduction

In this chapter we investigate ICA models in which the number of sources, M, may be less than the number of sensors, N: so-called non-square mixing.

The ‘extra’ sensor observations are explained as observation noise. This general approach may be called Probabilistic Independent Component Analysis (PICA) by analogy with the Probabilistic Principal Component Analysis (PPCA) model of Tipping & Bishop [I9971; ICA and PCA don't have observation noise, PICA and PPCA do.

Non-square ICA models give rise to a likelihood model for the data involving an integral which is intractable. In this chapter we build on previous work in which the integral is estimated using a Laplace approximation. By making the further assumption that the unmixing matrix lies on the decorrelating manifold we are able to make a number of simplifications. Firstly, the observation noise can be estimated using PCA methods, and, secondly, optimisation takes place in a space having a much reduced dimensionality, having order M2 parameters rather than M × N. Again, building on previous work, we derive a model order selection criterion for selecting the appropriate number of sources. This is based on the Laplace approximation as applied to the decorrelating manifold. This is then compared with PCA model order selection methods on music and EEG datasets.

Non-Gaussianity is of paramount importance in ICA estimation. Without non-Gaussianity the estimation is not possible at all (unless the independent components have time-dependences). Therefore, it is not surprising that non-Gaussianity could be used as a leading principle in ICA estimation.

In this chapter, we derive a simple principle of ICA estimation: the independent components can be found as the projections that maximize non-Gaussianity. In addition to its intuitive appeal, this approach allows us to derive a highly efficient ICA algorithm, Fast ICA. This is a fixed-point algorithm that can be used for estimating the independent components one by one. At the end of the chapter, it will be seen that it is closely connected to maximum likelihood or infomax estimation as well.

Whitening

First, let us consider preprocessing techniques that are essential if we want to develop fast ICA methods.

The rather trivial preprocessing that is used in many cases is to centre x, i.e. subtract its mean vector m = ε{x} so as to make x a zero-mean variable. This implies that s is zero-mean as well. This preprocessing is made solely to simplify the ICA algorithms: it does not mean that the mean could not be estimated. After estimating the mixing matrix A with centred data, we can complete the estimation by adding the mean vector of s back to the centred estimates of s.

An unsupervised classification algorithm is derived by modelling observed data as a mixture of several mutually exclusive classes that are each described by linear combinations of independent, non-Gaussian densities. The algorithm estimates the density of each class and is able to model class distributions with non-Gaussian structure. It can improve classification accuracy compared with standard Gaussian mixture models. When applied to images, the algorithm can learn efficient codes (basis functions) for images that capture the statistical structure of the images. We applied this method to the problem of unsupervised classification, segmentation and de-noising of images. This method was effective in classifying complex image textures such as trees and rocks in natural scenes. It was also useful for de-noising and filling in missing pixels in images with complex structures. The advantage of this model is that image codes can be learned with increasing numbers of classes thus providing greater flexibility in modelling structure and in finding more image features than in either Gaussian mixture models or standard ICA algorithms.

Introduction

Recently, Blind Source Separation by Independent Component Analysis has been applied to signal processing problems including speech enhancement, telecommunications and medical signal processing. ICA finds a linear non-orthogonal coordinate system in multivariate data determined by second- and higher-order statistics. The goal of ICA is to linearly transform the data in such a way that the transformed variables are as statistically independent from each other as possible [Jutten & Herault, 1991, Comon, 1994, Bell & Sejnowski, 1995, Cardoso & Laheld, 1996, Lee et al., 2000b]. ICA generalizes the technique of Principal Component Analysis (PCA) and, like PCA, has proven a useful tool for finding structure in data.

Introduction

Independent Component Analysis (ICA), as a signal processing tool, has shown great promise in many application domains, two of the most successful being telecommunications and biomedical engineering. With the growing awareness of ICA many other less obvious applications of the transform are starting to appear, for example financial time series prediction [Back & Weigend, 19981 and information retrieval [Isbell & Viola, 1999, Girolami, 2000a, Vinokourov & Girolami, 20001. In such applications ICA is being used as an unsupervised means of exploring and, hopefully, uncovering meaningful latent traits or structure within the data. In the case of financial time series prediction the factors which drive the evolution of the time series are hopefully uncovered, whereas within information retrieval the latent concepts or topics which generate key-words that occur in the documents are sought after. The strong assumption of independence of the hidden factors in ICA is difficult to argue for when there is limited a priori knowledge of the data, indeed it is desired that the analysis uncover injormative components which may, or may not, be independent. Nevertheless the independence assumption allows analytically tractable statistical models to be developed.

This chapter will consider how the standard ICA model can be extended and used in the unsupervised classification and visualisation of multivariate data. Prior to the formal presentation, to set the context of the remainder of the chapter, a short review of the ICA signal model and the corresponding statistical representation is given. The remaining sections propose ICA inspired techniques for the unsupervised classification and visualisation of multivariate data.