%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%  File         : 	chapter.tex
%%  Description  :      latex file for fg98 paper
%%  Author	 :	Baback Moghaddam
%%  Date         :    	1/12/96
%%  Note         :      adapted from ieee.tex format file
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\documentstyle[12pt,psfig,doublespace,smallcaptionfont,fullpage]{article}
\renewcommand{\baselinestretch}{1.394}


%\textheight      9.5in
%\textwidth       7in
%\oddsidemargin    -0.25 in
%\evensidemargin  -0.25 in
%\topmargin       -.85 in
%\headheight       0.0 in
%\headsep          .5 cm
%\footheight       0 in
%\flushbottom

\newcommand{\axis}[1]{\mbox{$ \hat{#1}$}}
\input{psfig}
\special{header=/v/lib/ps/lprep71.procs}
\pssilent
\newcommand{\dg}{\mbox{$^{\circ}$}}
\def\vect#1{{\bf #1}}


\title{{Probabilistic, View-based and Modular \\  Models for Face Recognition}}
\author{Baback Moghaddam and Alex Pentland \\
Perceptual Computing Section, The Media Laboratory, \\
Massachusetts Institute of Technology \\ 20 Ames Street, Cambridge MA
02139, U.S.A.\\ E-mail: {\tt baback@merl.com, sandy@media.mit.edu}}

\date{}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{document}

\begin{singlespace}
\maketitle
\end{singlespace}

\begin{singlespace}


\begin{abstract}

In this chapter we present several computational models for face
processing: a Bayesian or probabilistic approach to quantifying a
perceptually-inspired measure of facial similarity (for purposes of
face recognition), view-based or view-tuned models for recognition at
off-frontal views and lastly a modular facial description in terms of
key salient features such as the eyes, nose and the mouth. All three
of these methodologies rely, to varying degrees, on low-dimensional
descriptions of appearance as obtained by eigenspace decompositions
and in each its own way each builds on and extends the ``eigenface''
technique introduced by Turk \& Pentland \cite{Turk91}.


\end{abstract}

\end{singlespace}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\section{Introduction}
\label{sec:intro}


%The face is our primary focus of attention in social intercourse,
%playing a major role in conveying identity and emotion.
%The human ability to recognize faces is remarkable; we can
%recognize thousands of faces learned throughout our lifetime and
%identify familiar faces at a glance even after years of separation.
%This skill is quite robust, despite large changes in the visual
%stimulus due to viewing conditions, expression, aging, and
%distractions such as glasses or changes in hairstyle or facial hair.
%As a consequence the visual processing of human faces has fascinated
%philosophers and scientists for centuries, including such figures as
%Aristotle and Darwin.

%Aside from providing insight on high-level visual processing,
%computational models of face recognition will
%also to practical machine vision applications.  Computers which
%recognize faces will be applied to a wide variety of problems,
%including criminal identification, security systems, image and film
%processing, and human-computer interaction.  For example, the ability
%to model a particular face and distinguish it from a large number of
%stored face models will make it possible to vastly improve criminal
%identification.  Model-based coding of faces will benefit from the
%ability to recognize the subjects or to simply detect the presence of
%faces. 

Developing a computational model for recognition of natural objects
like faces is quite difficult, because they are complex,
multidimensional, and meaningful visual stimuli.  They are a natural
class of objects, and stand in stark contrast to sine wave gratings,
the ``blocks world'', and other artificial stimuli used in human and
computer vision research.  Thus unlike most early visual functions,
for which we may construct detailed models of retinal or striate
activity, recognition of natural objects such as faces is a high level
task for which computational approaches can currently only suggest
broad constraints on the corresponding neural activity.

Beginning with Pentland 1986 \cite{Pentland_ai86}, therefore, our
laboratory began to focus its research towards developing early,
preattentive pattern recognition capabilities that can deal with the
complexities of natural objects.  Although recognition of natural
objects is a high level visual problem, the laws of physics, biology,
and the geometry of the human body imposes quite a bit of structure
imposed on the task \cite{Pentland_ai86}.  We can capitalize
on this structure by finding ways to directly employ 2-D image
measurements to detect and recognize salient objects and classes.

The general approach we have developed is to attempt to describe the
{\sl range} of 2-D appearances of objects to be recognized.  To obtain
such an ``appearance-based'' representation, one must first transform
the image into a low-dimensional coordinate system that preserves the
general perceptual quality of the target object's image.  The
necessity for such a transformation is to address the `curse of
dimensionality': the raw image data has so many degrees of freedom
that it would require millions of examples to learn the range of
appearances directly.

Once a low-dimensional representation of the target class (face, eye,
hand, car, etc.) has been obtained, then standard methods can be used
to {\sl learn} the {\sl range of appearance} that the target exhibits
in the new, low-dimensional coordinate system.  Over the last decade
this approach has been used to recognize fractal textures and motions
(Pentland 1984, Heeger and Pentland 1985
\cite{Pentland_pami84,Heeger_on85}), decompose objects into their
component parts (Pentland 1989 \cite{Pentland_nc89}), read human lips
(Mase and Pentland 1990 \cite{Mase_90lips}), recognize hand
gestures (Darrell and Pentland 1993 \cite{Darrell93}),
recognize facial expressions (Essa and Pentland 1995
\cite{Essa_iccv95}), and even read American Sign Language
(Starner and Pentland 1995 \cite{Starner_fg95}).

And, of course, it has been applied to the problem of recognizing
human faces (Turk and Pentland 1991 \cite{Turk91}, Moghaddam and
Pentland 1997 \cite{Moghaddam_pami97}, Moghaddam, Wahid and Pentland
1998 \cite{Moghaddam_fg98} and Moghaddam, Jebara, and Pentland 1998
\cite{Moghaddam_icpr98}).  Face detection, recognition, and interpretation
is the most fully developed example of this general approach.  In part
this is because of the overwhelming importance of faces, but it is
also because faces are prime candidates for an appearance-based
approach to visual recognition: they present themselves quite
consistently in expected positions and orientations, their general
configuration of facial features is generally stable across faces, and
even changes in the appearance of a face are to some degree
predictable, such as the effect of aging or wearing sunglasses.

In this chapter we will confine ourselves almost completely to face
recognition, with only passing reference to recognition of other
objects.  The terminology, and most of the examples, will be drawn
from our Lab's most recent face recognition work, described in
Moghaddam and Pentland 1997 \cite{Moghaddam_pami97} and Moghaddam,
Jebara, and Pentland 1998 \cite{Moghaddam_icpr98}. However the reader
is requested to keep in mind that the approach described here can
equally well be applied to detection and recognition of other objects.



\subsection{Appearance-Based Detection and Recognition}

What do we mean by `the range of appearances of the human face'?
The `range of appearances' is precisely the probability density
function (PDF) of the image data for the target class.
For instance, given several examples of a target class $\Omega$ in
such a low-dimensional representation, it is straightforward to model
the probability distribution function $P({\bf x}|\Omega)$ of its
image-level features ${\bf x}$ as a mixture of Gaussians, thus
obtaining a low-dimensional, parametric {\sl appearance model} for the
target class \cite{moghaddam95}.  Once the target classes' probability
distribution function (PDF) has been learned, we can use Bayes' rule
to perform maximum {\sl a postieriori} (MAP) detection and
recognition.  

The result is a very simple, neural-net-like representation of the
target class's appearance, which can be used to detect occurances of
the class, to compactly describe its appearance, and to efficiently
compare different examples from the same class.  We have shown that
this method is very powerful for detection and recognition of human
faces, hands, and facial expressions \cite{Moghaddam_pami97}.  Other
researchers have used extensions of this basic method to recognize
industrial objects and household items \cite{Murase95}.

The use of parametric appearance models to characterize the PDF of an
object's appearance in the image is related to the simpler idea of a
view-based representation, as advocated by Ullman \cite{Ullman_pami91} and
Poggio \cite{Poggio_nature90}.  As originally developed, their idea of
view-based recognition was to accurately describe the spatial
structure of the target object by interpolating between
previously-seen views.  However in order to describe natural objects
such as faces or hands, which display a wide range of structural and
non-rigid variation, it is necessary to extend the notion of ``view''
to include characterizing the range of geometric and feature
variation, as well as the likelihood associated with such variation.
That is, one must use an appearance-based approach such as the one
described here, instead of the simpler idea of a view-based approach.

We typically use such the Karhunen-Lo\`{e}ve transform (KLT) (also
called Principal Components Analysis (PCA)) as the
dimensionality-reducing preprocessing transform.  This is because the
Karhunen-Lo\`{e}ve transform is known to provide an optimally-compact
linear basis (with respect to RMS error) for a given class of signal.
This transform has also been used by pioneers in face recognition
research, such as Kohonen \cite{kohonen} or Abdi \cite{Abdi_1988}.
However our use of the same preprocessing step leads to a false
impression of similarity.  Previous researchers have used the KLT as a
feature extraction step, which is followed by a simple classification
algorithm.  In contrast, we are using the KLT to facilitate learning
of the range of appearances (the PDF), which we then use to make
maximum {\sl a posteriori} (MAP) estimates for target detection and
recognition.

In this chapter we will first address the problem of {\sl detecting}
faces and facial features, and describe how our method of learning the
range of facial/feature appearances is used to accomplish this first,
all-important step.  We will then describe how to the same models of
facial appearance can be used for facial recognition, and report on
the robustness and accuracy of the combined detection/recognition
method.  Finally, we will discuss how recognition is extended to
included variation in head orientation, and how facial features can be
usefully included in the face recognition proccess.



\section{Visual Attention and Object Detection}

Visual attention is the process of restricting higher-level processing
to a subset of the visual field, referred to as the {\em
focus-of-attention} (FOA). The critical component of visual attention
is the {\em selection} of the FOA.  In humans this process is not
based purely on bottom-up processing and is in fact goal-driven.  The
measure of interest or {\em saliency} is modulated by the behavioral
state and the demands of the particular visual task that is currently
active.

Palmer \cite{Palmer83} has suggested that visual attention is the
process of locating the object of interest and placing it in a {\em
canonical} (or object-centered) reference frame suitable for
recognition (or template matching). We have developed a computational
technique for automatic object recognition, which is in accordance
with Palmer's model of visual attention.  The system uses a
probabilistic formulation for the estimation of the position and scale
of the object in the visual field and remaps the FOA to an
object-centered reference frame, which is subsequently used for
recognition and verification.

At a simple level the underlying mechanism of attention during a
visual search task can be based on a spatiotopic {\em saliency} map
$S(i,j)$ which is a function of the image information in a local
region $R$
\begin{equation}
S(i,j) = f\left[\{I(i+r,j+c): (r,c) \in R\}\right]
\end{equation}
For example saliency maps have been constructed which employ
spatio-temporal changes as cues for foveation \cite{Anderson85} or
other low-level image features such as local symmetry for detection of
interest points \cite{Reisfeld90}. However bottom-up techniques based
on low-level features lack {\em context} with respect to high-level
visual tasks such as object recognition. In a recognition task, the
selection of the FOA is driven by higher-level goals and therefore
requires internal representations of an object's appearance and a
means of comparing candidate objects in the FOA to the stored object
models.



\begin{figure}[tb]
\begin{center}
\begin{tabular}{cccc}
\psfig{figure=Images/white2_orig.ps,height=1.25in,clip=t} &
\psfig{figure=Images/white2_ms_eigsearch.ps,height=1.25in,clip=t} &
\psfig{figure=Images/white_scalewarp.ps,height=1.25in,clip=t} &
\psfig{figure=Images/white_eigsearch.ps,height=1.25in,clip=t} \\
%\psfig{figure=Images/hand2_orig.ps,height=1.0in,clip=t} &
%\psfig{figure=Images/hand2_detect.ps,height=1.0in,clip=t}\\
(a) & (b)
% & (c) & (d)
\end{tabular}
\end{center}
\vskip -4mm \caption[]{(a) input image, (b) face detection, 
 (c) face centering, (d) facial feature detection
}
\label{fig:face_hand_detections}
\end{figure}

Specifically, in an object-based visual search the saliency map is a
function of the degree of match between a candidate object in a local
image region and an internal model of the target object.  In
view-based recognition (as opposed to 3D geometric or invariant-based
recognition), the saliency can be formulated in terms of visual
similarity using a variety of metrics ranging from simple template
matching scores to more sophisticated measures using, for example,
robust statistics for image correlation \cite{Brunelli94}.  In this
chapter, however, we are primarily interested in saliency maps which
have a {\em probabilistic} interpretation as object-class membership
functions or {\em likelihoods}. These likelihood functions are learned
by applying density estimation techniques in complementary subspaces
obtained by an eigenvector decomposition. Our approach to this
learning problem is {\em view-based} --- {\it i.e.,} the learning and
modeling of the visual appearance of the object from a (suitably
normalized and preprocessed) set of training imagery.
Figure~\ref{fig:face_hand_detections} shows examples of the automatic
selection of FOA for detection of human faces. In each case, the
target object's probability distribution was {\em learned} from
training views and then subsequently used in computing likelihoods for
detection.  The face representation is based on the visual appearance
(normalized grayscale image). The maximum likelihood (ML) estimates of
position and scale are shown in the figure by the cross-hairs and
bounding box, respectively.




\subsection{Object Detection}

The standard detection paradigm in image processing is that of
normalized correlation or template matching.  However this approach is
only optimal in the simplistic case of a {\em deterministic} signal
embedded in additive white Gaussian noise. When we begin to consider a
target {\em class} detection problem --- {\it e.g,} finding a generic
human face or a human hand in a scene --- we must incorporate the
underlying probability distribution of the object. Subspace methods
and eigenspace decompositions are particularly well-suited to such a
task since they provide a compact and {\em parametric} description of
the object's appearance and also automatically identify the {\em
degrees-of-freedom} of the underlying statistical variability.

In particular, the eigenspace formulation leads to a powerful
alternative to standard detection techniques such as template matching
or normalized correlation.  The reconstruction error (or residual) of
the eigenspace decomposition (referred to as the
``distance-from-face-space'' in the context of the work with
``eigenfaces'' \cite{Turk91}) is an effective indicator of similarity.
The residual error is easily computed using the projection
coefficients and the original signal energy. This detection strategy
is equivalent to matching with a linear combination of {\em
eigentemplates} and allows for a greater range of distortions in the
input signal (including lighting, and moderate rotation and scale). In
a statistical signal detection framework, the use of eigentemplates
has been shown to yield superior performance in comparison with
standard matched filtering \cite{Kumar82}\cite{Pentland94}.

In \cite{Pentland94} we used this formulation for a modular eigenspace
representation of facial features where the corresponding residual ---
referred to as ``distance-from-{\em feature}-space'' or DFFS --- was
used for localization and detection.  Given an input image, a saliency
map was constructed by computing the DFFS at each pixel.  When using
$M$ eigenvectors, this requires $M$ convolutions (which can be
efficiently computed using an FFT) plus an additional local energy
computation.  The global minimum of this distance map was then
selected as the best estimate of the location of the target.

In this chapter we will show that the DFFS can be interpreted as an
estimate of a marginal component of the probability density of the
object and that a complete estimate must also incorporate a second
marginal density based on a complementary ``distance-{\em
in}-feature-space'' (DIFS).  Using our estimates of the object
densities, we formulate the problem of target detection from the point
of view of a maximum likelihood (ML) estimation problem.
Specifically, given the visual field, we estimate the position (and
scale) of the image region which is most representative of the target
of interest.  Computationally this is achieved by sliding an
$m$-by-$n$ observation window throughout the image and at each
location computing the {\em likelihood} that the local subimage
$\vect{x}$ is an instance of the target class $\Omega$ --- {\it i.e.,}
$P(\vect{x}|\Omega)$.  After this probability map is computed, we
select the location corresponding to the highest likelihood as our ML
estimate of the target location.  Note that the likelihood map can be
evaluated over the entire parameter space affecting the object's
appearance which can include transformations such as scale and
rotation.



\section{Eigenspace Methods for Visual Modeling}


In recent years, computer vision research has witnessed a growing
interest in eigenvector analysis and subspace decomposition methods.
In particular, eigenvector decomposition has been shown to be an
effective tool for solving problems which use high-dimensional
representations of phenomena which are intrinsically low-dimensional.
This general analysis framework lends itself to several closely
related formulations in object modeling and recognition which employ
the {\em principal modes} or characteristic {\em degrees-of-freedom}
for description.  The identification and parametric representation of
a system in terms of these principal modes is at the core of recent
advances in physically-based modeling \cite{Pentland91},
correspondence and matching \cite{Sclaroff95}, and parametric
descriptions of shape \cite{Cootes92}.

Eigenvector-based methods also form the basis for data analysis
techniques in pattern recognition and statistics where they are used
to extract low-dimensional subspaces comprised of statistically
uncorrelated variables which tend to simplify tasks such as
classification. The Karhunen-Loeve Transform (KLT) \cite{Loeve55} and
Principal Components Analysis (PCA) \cite{Jolliffe86} are examples of
eigenvector-based techniques which are commonly used for
dimensionality reduction and feature extraction in pattern
recognition.

In computer vision, eigenvector analysis of {\em imagery} has been
used for characterization of human faces \cite{Kirby90} and automatic
face recognition using ``eigenfaces'' \cite{Turk91}\cite{Pentland94}.
More recently, principal component analysis of imagery has also been
applied for robust target detection \cite{Pentland94}\cite{Burl94},
nonlinear image interpolation \cite{Bregler94}, visual learning for
object recognition \cite{Murase95}\cite{Weng94}, as well as visual
servoing for robotics
\cite{Nayar94}.

Specifically, Murase \& Nayar \cite{Murase95} used a low-dimensional
{\em parametric} eigenspace for recovering object identity and pose by
matching views to a spline-based hypersurface.  Nayar {\it et al.}
\cite{Nayar94} have extended this technique to visual feedback control
and servoing for a robotic arm in ``peg-in-the-hole'' insertion tasks.
Pentland {\it et al.} \cite{Pentland94} proposed a view-based
multiple-eigenspace technique for face recognition under varying pose
as well as for the detection and description of facial features.
Similarly, Burl {\it et al.} \cite{Burl94} used Bayesian
classification for object detection using a feature vector derived
from principal component images. Weng \cite{Weng94} has proposed a
visual learning framework based on the KLT in conjunction with an
optimal linear discriminant transform for learning and recognition of
objects from 2D views.


\subsection{Probabilistic Eigenspaces}

However, these authors (with the exception of \cite{Pentland94}) have
used eigenvector analysis primarily as a dimensionality reduction
technique for subsequent modeling, interpolation, or classification.
In contrast, our methods use an eigenspace decomposition as an
integral part of an efficient technique for probability density
estimation of high-dimensional data.

%In this section we present our recent work using eigenspace
%decompositions for object representation and modeling.
Our learning
method estimates the {\em complete} probability distribution of the
object's appearance using an eigenvector decomposition of the image
space.
%of the sample covariance
%matrix of a set of training views.
The desired target density is decomposed into two components: the
density in the principal subspace (containing the
traditionally-defined principal components) and its orthogonal
complement (which is usually discarded in standard PCA).  We have
derived the form for an optimal density estimate for the case of
Gaussian data and a near-optimal estimator for arbitrarily complex
distributions in terms of a Mixture-of-Gaussians density model
\cite{Moghaddam_pami97}.

We note that this learning method differs from {\em supervised} visual
learning with function approximation networks \cite{Poggio90} in which
a hypersurface representation of an input/output map is automatically
learned from a set of training examples.  Instead, we use a
probabilistic formulation which combines the two standard paradigms of
{\em unsupervised} learning --- PCA and density estimation --- to
arrive at a computationally feasible estimate of the class conditional
density function.
% $P(\vect{x}|\Omega)$ for an object based on its
%(high-dimensional) imagery $\vect{x}$. These density estimates can
then used for maximum likelihood detection of faces and facial
features, as well as Bayesian modeling for recognition.


The key to our approach to automatic visual learning is density
estimation. However, instead of applying estimation techniques
directly to the original high-dimensional space of the imagery, we use
an eigenspace decomposition to yield a computationally feasible
estimate. Specifically, the eigenspace analysis is applied to a set of
training views of the object in order to identify a principal subspace
which captures the {\em intrinsic} dimensionality of the data. The
component of the complete density in this lower-dimensional subspace
is then estimated using a suitable parametric form. In addition, we
implicitly model the component of the distribution in the {\em
orthogonal} subspace.  The complete density estimate can be
efficiently computed from the lower-dimensional principal components.
Our density estimate is shown to be {\em optimal} in the case of
Gaussian-distributed training data. 
% We also formulate an approximate
%(near-optimal) density estimate for more realistic data with arbitrary
%and multimodal distributions.



\subsubsection{Principal Component Imagery}

Given a training set of $m$-by-$n$ images $\{I^t\}_{t=1}^{N_T}$, we
can form a training set of vectors $\{\vect{x}^t\}$, where $\vect{x}
\in {\cal R}^{N=m n}$, by lexicographic ordering of the pixel elements
of each image $I^t$. The basis functions for the KLT \cite{Loeve55}
are obtained by solving the eigenvalue problem
\begin{equation}
\Lambda ~=~ \Phi^T \Sigma \Phi 
\end{equation}
where $\Sigma$ is the covariance matrix, $\Phi$ is the eigenvector
matrix of $\Sigma$ and $\Lambda$ is the corresponding diagonal matrix
of eigenvalues. The unitary matrix $\Phi$ defines a coordinate
transform (rotation) which {\em decorrelates} the data and makes
explicit the {\em invariant subspaces} of the matrix operator
$\Sigma$.  In PCA, a partial KLT is performed to identify the
largest-eigenvalue eigenvectors and obtain a principal component
feature vector $\vect{y} = \Phi_{\small M}^T ~\vect{\tilde{x}}$, where
$\vect{\tilde{x}} = \vect{x} - \vect{\bar{x}}$ is the mean-normalized
image vector and $\Phi_M$ is a submatrix of $\Phi$ containing the
principal eigenvectors. PCA can be seen as a linear transformation
$\vect{y} = {\cal T(\vect{x})}: {\cal R}^N \rightarrow {\cal R}^M$
which extracts a lower-dimensional subspace of the KL basis
corresponding to the maximal eigenvalues.  These principal components
preserve the major linear correlations in the data and discard the
minor ones.\footnote{In practice the number of training images $N_T$
is far less than the dimensionality of the imagery $N$, consequently
the covariance matrix $\Sigma$ is singular. However, the first $M <
N_T$ eigenvectors can always be computed (estimated) from $N_t$
samples using, for example, a Singular Value Decomposition
\cite{Golub89}.}

\begin{figure}[tb]
\begin{center}
\begin{tabular}{cc}
\raisebox{.1in}{\psfig{figure=Images/featurespace1.ps,height=2in,clip=t}} &
\psfig{figure=Images/eigspectrum_divide.ps,height=1.8in,clip=t} \\
(a) & (b) 
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) Decomposition into the principal subspace $F$ and its
orthogonal complement $\bar{F}$ for a Gaussian density, (b) a typical
eigenvalue spectrum and its division into the two orthogonal
subspaces.}
\label{fig:unimodal}
\end{figure}

By ranking the eigenvectors of the KL expansion with respect to their
eigenvalues and selecting the first $M$ principal components we form
an orthogonal decomposition of the vector space ${\cal R}^N$ into two
mutually exclusive and complementary subspaces: the principal subspace
(or feature space) $F=\{\Phi_i\}_{i=1}^M$ containing the principal
components and its orthogonal complement
$\bar{F}=\{\Phi_i\}_{i=M+1}^N$. This orthogonal decomposition is
illustrated in Figure~\ref{fig:unimodal}(a) where we have a
prototypical example of a distribution which is embedded entirely in
$F$. In practice there is always a signal component in $\bar{F}$ due
to the minor statistical variabilities in the data or simply due to
the observation noise which affects every element of $\vect{x}$.


In a partial KL expansion, the residual reconstruction error is
defined as
\begin{equation}
\epsilon^2(\vect{x}) ~=~ \sum_{i=M+1}^{N} y_i^2  ~=~ || \vect{\tilde{x}} ||^2
~-~ \sum_{i=1}^{M} y_i^2
\label{eq:dffs}
\end{equation}
and can be easily computed from the first $M$ principal components and
the $L_2$ norm of the mean-normalized image $\vect{\tilde{x}}$.
Consequently the $L_2$ norm of every element $\vect{x} \in {\cal R}^N$
can be decomposed in terms of its projections in these two subspaces.
We refer to the component in the orthogonal subspace $\bar{F}$ as the
``distance-from-feature-space'' (DFFS) which is a simple Euclidean
distance and is equivalent to the residual error
$\epsilon^2(\vect{x})$ in Eq.(\ref{eq:dffs}). The component of
$\vect{x}$ which lies {\em in} the feature space $F$ is referred to as
the ``distance-in-feature-space'' (DIFS) but is generally not a
distance-based norm, but can be interpreted in terms of the
probability distribution of $y$ in $F$.



\subsubsection{Density Estimation in  Eigenspace}
\label{sec:density}
\label{sec:gaussian}

One difficulty with probabilistic visual modeling is that the
intensity or intensity difference vectors are very high-dimensional,
with ${\Delta} \in {\cal R}^N$ and $N=O(10^4)$. Therefore we typically
lack sufficient independent training observations to compute reliable
2nd-order statistics for the likelihood densities ({\it i.e.,}
singular covariance matrices will result). Even if we were able to
estimate these statistics, the computational cost of evaluating the
likelihoods is formidable. Furthermore, this computation would be
highly inefficient since the {\em intrinsic} dimensionality or major
degrees-of-freedom of the data vectors each class is likely to be
significantly smaller than $N$.

Recently, an efficient density estimation method was proposed by
Moghaddam \& Pentland \cite{Moghaddam_pami97} which divides the vector
space ${\cal R}^N$ into two complementary subspaces using an
eigenspace decomposition. This method relies on a Principal Components
Analysis (PCA) \cite{Jolliffe86} to form a low-dimensional estimate of
the complete likelihood which can be evaluated using only the first
$M$ principal components, where $M<<N$.  This decomposition is
illustrated in Figure~\ref{fig:unimodal} which shows an orthogonal
decomposition of the vector space ${\cal R}^N$ into two mutually
exclusive subspaces: the principal subspace $F$ containing the first
$M$ principal components and its orthogonal complement $\bar{F}$,
which contains the residual of the expansion.  The component in the
orthogonal subspace $\bar{F}$ is the so-called
``distance-from-feature-space'' (DFFS), a Euclidean distance
equivalent to the PCA residual error.  The component of $\Delta$ which
lies {\em in} the feature space $F$ is referred to as the
``distance-in-feature-space'' (DIFS) and is a {\em Mahalanobis}
distance for Gaussian densities.


%\begin{figure}[tb]
%\begin{center}
%\begin{tabular}{c}
%\psfig{figure=/v/users/baback/Text/Chapter/Ship/Images/featurespace1.ps,height=2in,clip=t}
%\\
%(a)
%\\
%\psfig{figure=/v/users/baback/Text/Chapter/Ship/Images/eigspectrum_divide.ps,height=2in,clip=t} \\
% (b) 
%\end{tabular}
%\end{center}
%\vskip -4mm 
%\caption{(a) Decomposition of ${\cal R}^N$ into the principal subspace
%           $F$ and its orthogonal complement $\bar{F}$ for a Gaussian
%           density, (b) a typical eigenvalue spectrum and its division
%           into the two orthogonal subspaces.}
%\label{fig:unimodal}
%\end{figure}


As shown in \cite{Moghaddam_pami97}, the complete likelihood estimate can
be written as the product of two independent marginal Gaussian
densities
\begin{equation}
\begin{array}{ccl}
\hat{P}({\vec{x}}|
\Omega) & = &  
\left[ \frac{{\rm exp}\left(-\frac{1}{2} {\displaystyle \sum_{i=1}^M} \frac{y_i^2}{\lambda_i}
 \right)}{(2\pi)^{M/2}~{\displaystyle \prod_{i=1}^M} \lambda_i^{1/2}} \right]
\cdot \left[
\frac{{\rm exp}\left(-\frac{\epsilon^2({\vec{x}})}{2\rho} \right)}
{(2\pi\rho)^{(N-M)/2}} \right] \\
& & \\
& = &  P_F({\vec{x}}|\Omega) ~\hat{P}_{\bar{F}}({\vec{x}}|\Omega) 
\end{array}
\label{eq:P_hat}
\end{equation}
where $P_F({\vec{x}}|\Omega)$ is the true marginal density in
$F$, $\hat{P}_{\bar{F}}({\vec{x}}|\Omega)$ is the estimated
marginal density in the orthogonal complement $\bar{F}$, $y_i$ are the
principal components and $\epsilon^2({\vec{x}})$ is the residual
(or DFFS).  The optimal value for the weighting parameter $\rho$ is 
then found to be simply the average of the $\bar{F}$ eigenvalues
\begin{equation}
\rho ~=~ \frac{1}{N-M} \sum_{i=M+1}^{N} \lambda_i
\label{eq:rho}
\end{equation}
We note that in actual practice, the majority of the $\bar{F}$
eigenvalues are unknown but {\em can} be estimated, for example, by
fitting a nonlinear function to the available portion of the
eigenvalue spectrum and estimating the average of the eigenvalues
beyond the principal subspace.




\subsection{Maximum Likelihood Detection}

The density estimate $\hat{P}(\vect{x}|\Omega)$ can be used to compute
a local measure of target saliency at each spatial position $(i,j)$ in
an input image based on the vector $\vect{x}$ obtained by the
lexicographic ordering of the pixel values in a local neighborhood $R$
\begin{equation}
S(i,j;\Omega) = \hat{P}({\vect{x}} | \Omega) ~, ~~~ \vect{x} ~=~
\downarrow \left[ \{I(i+r,j+c) : (r,c) \in R\} \right]
\end{equation}
where $\downarrow[\bullet]$ is the operator which converts a subimage
into a vector. The ML estimate of position of the target
$\Omega$ is then given by
\begin{equation}
(i,j)^{\rm ML} ~=~ {\rm argmax} ~S(i,j;\Omega)
\end{equation}  
This ML formulation can be extended to estimate object scale with {\em
multiscale} saliency maps. The likelihood computation is performed (in
parallel) on linearly scaled versions of the input image
$I^{(\sigma)}$ corresponding to a pre-determined set of scales
$\{\sigma_1,\sigma_2,\cdots \sigma_n\}$
\begin{equation}
S(i,j,k;~\Omega) = \hat{P} \left( \downarrow \{I^{(\sigma_k)}(\sigma_k
i + r, \sigma_k j + c)
: (r,c) \in R\}  ~|~ \Omega \right)
\end{equation}
where the ML estimate of the spatial and scale indices is defined by
\begin{equation}
(i,j,k)^{\rm ML} ~=~ {\rm argmax} ~S(i,j,k ~;~\Omega)
\end{equation}  


\begin{figure}[t]
\begin{center}
\begin{tabular}{c}
\psfig{figure=/v/users/baback/Text/ICCV95/New_Images/saliency.ps,height=2.5in,clip=t}
\end{tabular}
\end{center}
\caption{Target saliency map S(i,j) showing the probability of a left eye
pattern over the input image.}
\label{fig:saliency}
\end{figure}


One important factor of variability in the appearance of the object in
grayscale imagery is that of lighting and contrast. However, one can
normalize for global illumination changes (as well as the linear
response characteristics of the CCD camera) by normalizing each
subimage $\vect{x}$ by its mean and standard deviation. This lighting
normalization is performed both during training (density estimation)
and also in the operational mode ({\em e.g.,} in detection).

This maximum likelihood detection framework can be viewed as a
Bayesian formulation of some neural network approaches to target
detection.  Perhaps the most closely related is the neural network
face detector of Sung \& Poggio \cite{Sung94} which is essentially a
trainable nonlinear binary pattern classifier.  They too learn the
distribution of the object class with a Mixture-of-Gaussians model
(using an elliptical k-means algorithm instead of EM). Instead of
likelihoods, however, input patterns are represented by a set of
distances to each mixture component (similar to a combination of the
DIFS and DFFS), thus forming a feature vector indicative of the
overall class membership. In addition, Sung \& Poggio explicitly model
the ``not-class'' by learning the distribution of nearby {\em
non-face} patterns. The set of distances to both classes are then used
to train a neural network to discriminate between face and non-face
patterns (similar to computing a likelihood ratio in MAP).

\section{Bayesian Model of Facial Similarity}

Current approaches to image matching for visual object recognition and
image database retrieval often make use of simple image similarity
metrics such as Euclidean distance or normalized correlation, which
correspond to a standard template-matching approach to recognition
\cite{Brunelli_Poggio}. For example, in its simplest form, the
similarity measure $S(I_1,I_2)$ between two images $I_1$ and $I_2$ can
be set to be inversely proportional to the norm $||I_1 - I_2||$. Such
a simple formulation suffers from a major drawback: it does not
exploit knowledge of which type of variations are critical (as opposed
to incidental) in expressing similarity.  In this chapter, we
formulate a {\em probabilistic} similarity measure which is based on
the probability that the image intensity differences, denoted by
$\Delta = I_1 - I_2$, are characteristic of typical variations in
appearance of the {\em same} object.  For example, for purposes of
face recognition, we can define two classes of facial image
variations: {\em intrapersonal} variations $\Omega_I$ (corresponding,
for example, to different facial expressions of the {\em same}
individual) and {\em extrapersonal} variations $\Omega_E$
(corresponding to variations between {\em different} individuals). Our
similarity measure is then expressed in terms of the probability
\begin{equation}
S(I_1,I_2) ~=~ P(\Delta \in \Omega_I) ~=~ P(\Omega_I | \Delta)
\end{equation}
where $P(\Omega_I | \Delta)$ is the {\em a posteriori} probability
given by Bayes rule, using estimates of the likelihoods $P(\Delta |
\Omega_I)$ and $P(\Delta | \Omega_E)$ which are derived from training
data using an efficient subspace method for density estimation of
high-dimensional data \cite{Moghaddam_pami97}. This Bayesian (MAP)
approach can also be viewed as a generalized nonlinear extension of
Linear Discriminant Analysis (LDA) \cite{Weng_pami96,Etemad_icassp96}
or ``FisherFace'' techniques \cite{Fisherfaces_pami97} for face
recognition. Moreover, our nonlinear generalization has distinct
computational/storage advantages over these linear methods for large
databases.


\subsection{Analysis of Intensity Differences}

We now consider the problem of characterizing the type of differences
which occur when matching two images in a face recognition task. We
define two distinct and mutually exclusive classes: $\Omega_I$
representing {\em intrapersonal} variations between multiple images of
the same individual ({\it e.g.,} with different expressions and
lighting conditions), and $\Omega_E$ representing {\em extrapersonal}
variations which result when matching two different individuals. We
will assume that both classes are Gaussian-distributed and seek to
obtain estimates of the likelihood functions $P({\Delta} | \Omega_I)$
and $P({\Delta} | \Omega_E)$ for a given intensity difference
$\Delta = I_1 - I_2$.

Given these likelihoods we can define the similarity score
$S(I_1,I_2)$ between a pair of images directly in terms of the
intrapersonal {\em a posteriori} probability as given by Bayes rule:
\begin{equation}
\begin{array}{ccl}
S & = & P(\Omega_I|{\Delta}) \\

  & = & {\displaystyle \frac{P({\Delta}|\Omega_I)
P(\Omega_I)}{P({\Delta}|\Omega_I) P(\Omega_I) + P({
\Delta}|\Omega_E) P(\Omega_E)}}

\end{array}
\label{eq:similarity}
\end{equation}
where the priors $P(\Omega)$ can be set to reflect specific operating
conditions ({\it e.g.,} number of test images {\it vs.} the size of
the database) or other sources of {\em a priori} knowledge regarding
the two images being matched. Additionally, this particular Bayesian
formulation casts the standard face recognition task (essentially an
$M$-ary classification problem for $M$ individuals) into a {\em
binary} pattern classification problem with $\Omega_I$ and
$\Omega_E$. This much simpler problem is then solved using the maximum
{\em a posteriori} (MAP) rule --- {\it i.e.,} two images are
determined to belong to the same individual if $P(\Omega_I|{\Delta}) >
P(\Omega_E|{\Delta})$, or equivalently, if $S(I_1,I_2) > \frac{1}{2}$.


We note that the Bayesian classification of identity outlined above is
perhaps closely related to human ``categorical perception'' and the
{\it a posteriori} probability itself is perhaps a more meaningful and
accurate computational model of {\em perceptual similarity} as judged
by humans.




%\section{Density Estimation with Eigenspace}

%In this section we present our recent work using eigenspace
%decompositions for object representation and modeling. Our learning
%method estimates the {\em complete} probability distribution of the
%object's appearance using an eigenvector decomposition of the image
%space.
%%of the sample covariance
%%matrix of a set of training views.
%The desired target density is decomposed into two components: the
%density in the principal subspace (containing the
%traditionally-defined principal components) and its orthogonal
%complement (which is usually discarded in standard PCA).  We derive
%the form for an optimal density estimate for the case of Gaussian data
%and a near-optimal estimator for arbitrarily complex distributions in
%terms of a Mixture-of-Gaussians density model.

%We note that this learning method differs from {\em supervised} visual
%learning with function approximation networks \cite{Poggio90} in which
%a hypersurface representation of an input/output map is automatically
%learned from a set of training examples.  Instead, we use a
%probabilistic formulation which combines the two standard paradigms of
%{\em unsupervised} learning --- PCA and density estimation --- to
%arrive at a computationally feasible estimate of the class conditional
%density function.
%% $P(\vect{x}|\Omega)$ for an object based on its
%%(high-dimensional) imagery $\vect{x}$.

%%These density estimates are
%%then used for maximum likelihood detection of faces and articulated
%%hands in natural images.


%%The key to our approach to automatic visual learning is density
%%estimation. However, instead of applying estimation techniques
%%directly to the original high-dimensional space of the imagery, we use
%%an eigenspace decomposition to yield a computationally feasible
%%estimate. Specifically, the eigenspace analysis is applied to a set of
%%training views of the object in order to identify a principal subspace
%%which captures the {\em intrinsic} dimensionality of the data. The
%%component of the complete density in this lower-dimensional subspace
%%is then estimated using a suitable parametric form. In addition, we
%%implicitly model the component of the distribution in the {\em
%%orthogonal} subspace.  The complete density estimate can be
%%efficiently computed from the lower-dimensional principal components.
%%Our density estimate is shown to be {\em optimal} in the case of
%%Gaussian-distributed training data.  We also formulate an approximate
%%(near-optimal) density estimate for more realistic data with arbitrary
%%and multimodal distributions.


%Specifically, given a set of training images
%$\{\vect{x}^t\}_{t=1}^{N_T}$, from an object class $\Omega$, we wish
%to estimate the class membership or {\em likelihood} function for this
%data --- {\it i.e.,} $P(\vect{x}|\Omega)$. In this section, we examine
%two density estimation techniques for visual learning of
%high-dimensional data. The first method is based on the assumption of
%a Gaussian distribution while the second method generalizes to
%arbitrarily complex distributions using a Mixture-of-Gaussians density
%model. Before introducing these estimators we briefly review
%eigenvector decomposition as commonly used in PCA.



%\subsection{Principal Component Imagery}

%Given a training set of $m$-by-$n$ images $\{I^t\}_{t=1}^{N_T}$, we
%can form a training set of vectors $\{\vect{x}^t\}$, where $\vect{x}
%\in {\cal R}^{N=m n}$, by lexicographic ordering of the pixel elements
%of each image $I^t$. The basis functions for the KLT \cite{Loeve55}
%are obtained by solving the eigenvalue problem
%\begin{equation}
%\Lambda ~=~ \Phi^T \Sigma \Phi 
%\end{equation}
%where $\Sigma$ is the covariance matrix, $\Phi$ is the eigenvector
%matrix of $\Sigma$ and $\Lambda$ is the corresponding diagonal matrix
%of eigenvalues. The unitary matrix $\Phi$ defines a coordinate
%transform (rotation) which {\em decorrelates} the data and makes
%explicit the {\em invariant subspaces} of the matrix operator
%$\Sigma$.  In PCA, a partial KLT is performed to identify the
%largest-eigenvalue eigenvectors and obtain a principal component
%feature vector $\vect{y} = \Phi_{\small M}^T ~\vect{\tilde{x}}$, where
%$\vect{\tilde{x}} = \vect{x} - \vect{\bar{x}}$ is the mean-normalized
%image vector and $\Phi_M$ is a submatrix of $\Phi$ containing the
%principal eigenvectors. PCA can be seen as a linear transformation
%$\vect{y} = {\cal T(\vect{x})}: {\cal R}^N \rightarrow {\cal R}^M$
%which extracts a lower-dimensional subspace of the KL basis
%corresponding to the maximal eigenvalues.  These principal components
%preserve the major linear correlations in the data and discard the
%minor ones.\footnote{In practice the number of training images $N_T$
%is far less than the dimensionality of the imagery $N$, consequently
%the covariance matrix $\Sigma$ is singular. However, the first $M <
%N_T$ eigenvectors can always be computed (estimated) from $N_t$
%samples using, for example, a Singular Value Decomposition
%\cite{Golub89}.}

%\begin{figure}[tb]
%\begin{center}
%\begin{tabular}{cc}
%\raisebox{.1in}{\psfig{figure=Images/featurespace1.ps,height=2in,clip=t}} &
%\psfig{figure=Images/eigspectrum_divide.ps,height=1.8in,clip=t} \\
%(a) & (b) 
%\end{tabular}
%\end{center}
%\vskip -4mm \caption{(a) Decomposition into the principal subspace $F$ and its
%orthogonal complement $\bar{F}$ for a Gaussian density, (b) a typical
%eigenvalue spectrum and its division into the two orthogonal
%subspaces.}
%\label{fig:unimodal}
%\end{figure}

%By ranking the eigenvectors of the KL expansion with respect to their
%eigenvalues and selecting the first $M$ principal components we form
%an orthogonal decomposition of the vector space ${\cal R}^N$ into two
%mutually exclusive and complementary subspaces: the principal subspace
%(or feature space) $F=\{\Phi_i\}_{i=1}^M$ containing the principal
%components and its orthogonal complement
%$\bar{F}=\{\Phi_i\}_{i=M+1}^N$. This orthogonal decomposition is
%illustrated in Figure~\ref{fig:unimodal}(a) where we have a
%prototypical example of a distribution which is embedded entirely in
%$F$. In practice there is always a signal component in $\bar{F}$ due
%to the minor statistical variabilities in the data or simply due to
%the observation noise which affects every element of $\vect{x}$.


%In a partial KL expansion, the residual reconstruction error is
%defined as
%\begin{equation}
%\epsilon^2(\vect{x}) ~=~ \sum_{i=M+1}^{N} y_i^2  ~=~ || \vect{\tilde{x}} ||^2
%~-~ \sum_{i=1}^{M} y_i^2
%\label{eq:dffs}
%\end{equation}
%and can be easily computed from the first $M$ principal components and
%the $L_2$ norm of the mean-normalized image $\vect{\tilde{x}}$.
%Consequently the $L_2$ norm of every element $\vect{x} \in {\cal R}^N$
%can be decomposed in terms of its projections in these two subspaces.
%We refer to the component in the orthogonal subspace $\bar{F}$ as the
%``distance-from-feature-space'' (DFFS) which is a simple Euclidean
%distance and is equivalent to the residual error
%$\epsilon^2(\vect{x})$ in Eq.(\ref{eq:dffs}). The component of
%$\vect{x}$ which lies {\em in} the feature space $F$ is referred to as
%the ``distance-in-feature-space'' (DIFS) but is generally not a
%distance-based norm, but can be interpreted in terms of the
%probability distribution of $y$ in $F$.


%\subsection{Gaussian Densities}
%\label{sec:gaussian} 

%We begin by considering an optimal approach for estimating
%high-dimensional Gaussian densities.  We assume that we have
%(robustly) estimated the mean $\vect{\bar{x}}$ and covariance $\Sigma$
%of the distribution from the given training set
%$\{\vect{x}^t\}$.\footnote{In practice, a full rank $N$-dimensional
%covariance $\Sigma$ can not be estimated from $N_T$ independent
%observations when $N_T<N$. But as we shall see our estimator does not
%require the full covariance, but only its first $M$ principal
%eigenvectors where $M<N_T$.} Under this assumption, the likelihood of
%an input pattern $\vect{x}$ is given by
%\begin{equation}
%P(\vect{x}|\Omega) = \frac{\mbox{exp} \left[ -\frac{1}{2}
%(\vect{x}-\vect{\bar{x}})^T \Sigma^{-1} (\vect{x}-\vect{\bar{x}})
%\right] } {(2\pi)^{N/2} ~|\Sigma|^{1/2}}  
%\end{equation}		
%The sufficient statistic for characterizing this likelihood is the
%{\em Mahalanobis} distance
%\begin{equation}
%d(\vect{x}) ~=~ \vect{\tilde{x}}^T \Sigma^{-1} \vect{\tilde{x}} 
%\label{eq:mahalanobis}
%\end{equation}
%where $\vect{\tilde{x}} = \vect{x} - \vect{\bar{x}}$.  However,
%instead of evaluating this quadratic product explicitly, a much more
%efficient and robust computation can be performed, especially with
%regard to the matrix inverse $\Sigma^{-1}$.  Using the eigenvectors
%and eigenvalues of $\Sigma$ we can rewrite $\Sigma^{-1}$ in the
%diagonalized form
%\begin{equation}
%\begin{array}{ccl}
%d(\vect{x}) & = & \vect{\tilde{x}}^T \Sigma^{-1} \vect{\tilde{x}} \\
%  & = & \vect{\tilde{x}}^T \left[ \Phi \Lambda^{-1} \Phi^T \right] \vect{\tilde{x}} \\
%  & = & \vect{y}^T \Lambda^{-1} \vect{y}   
%\end{array}
%\end{equation}
%where $\vect{y} = \Phi^T \vect{\tilde{x}}$ are the new variables
%obtained by the change of coordinates in a KLT. Because of the
%diagonalized form, the {\em Mahalanobis} distance can also be
%expressed in terms of the sum
%\begin{equation}
%d(\vect{x}) ~=~ \sum_{i=1}^{N} \frac{y_i^2}{\lambda_i}
%\label{eq:mahalanobis_decoupled}
%\end{equation}
%In the KLT basis, the {\em Mahalanobis} distance in
%Eq.(\ref{eq:mahalanobis}) is conveniently {\em decoupled} into a
%weighted sum of uncorrelated component energies. Furthermore, the
%likelihood becomes a {\em product} of independent separable Gaussian
%densities. Despite its simpler form, evaluation of
%Eq.(\ref{eq:mahalanobis_decoupled}) is still computationally infeasible due to
%the high-dimensionality. We therefore seek to {\em estimate}
%$d(\vect{x})$ using only $M$ projections.  Intuitively, an obvious
%choice for a lower-dimensional representation is the principal
%subspace indicated by PCA which captures the major degrees of
%statistical variability in the data.\footnote{We will see shortly that
%given the typical eigenvalue spectra observed in practice ({\em e.g.,}
%Figure~\ref{fig:unimodal}(b)), this choice is optimal for a
%different reason: it minimizes the information-theoretic {\em
%divergence} between the true density and our estimate of it.}
%Therefore, we divide the summation into two independent parts
%corresponding to the principal subspace $F =
%\{\Phi_i\}_{i=1}^{M}$ and its orthogonal complement
%$\bar{F}=\{\Phi_i\}_{i=M+1}^{N}$
%\begin{equation}
%d(\vect{x}) ~=~ \sum_{i=1}^{M} \frac{y_i^2}{\lambda_i} ~+~
%\sum_{i=M+1}^{N} \frac{y_i^2}{\lambda_i}
%\label{eq:d}
%\end{equation}
%We note that the terms in the first summation can be computed by
%projecting $\vect{x}$ onto the $M$-dimensional principal subspace $F$.
%The remaining terms in the second sum $\{y_i\}_{i=M+1}^{N}$, however,
%can not be computed explicitly in practice because of the
%high-dimensionality. However, the {\em sum} of these terms is
%available and is in fact the DFFS quantity $\epsilon^2(\vect{x})$
%which can be computed from Eq.(\ref{eq:dffs}).  Therefore, based on
%the available terms, we can formulate an estimator for $d(\vect{x})$
%as follows
%\begin{equation}
%\begin{array}{ccl}
%	\hat{d}(\vect{x}) & = & {\displaystyle \sum_{i=1}^{M}
%\frac{y_i^2}{\lambda_i}} ~~+~~ 	\frac{1}{\rho} \left[ {\displaystyle
%\sum_{i=M+1}^{N} y_i^2} \right]
%\\
%	   & = & {\displaystyle \sum_{i=1}^{M} \frac{y_i^2}{\lambda_i} ~~+~~
%	\frac{\epsilon^2(\vect{x})}{\rho}} 
%\label{eq:d_hat}
%\end{array}
%\end{equation}
%where the term in the brackets is $\epsilon^2(\vect{x})$, which as we
%have seen can be computed using the first $M$ principal components.
%We can therefore write the form of the likelihood estimate based on
%$\hat{d}(\vect{x})$ as the product of two marginal and independent
%Gaussian densities
%\begin{equation}
%\begin{array}{ccl}
%\hat{P}(\vect{x}|\Omega) & = &  
%\left[ \frac{{\rm exp}\left(-\frac{1}{2} {\displaystyle \sum_{i=1}^M} \frac{y_i^2}{\lambda_i}
% \right)}{(2\pi)^{M/2}~{\displaystyle \prod_{i=1}^M} \lambda_i^{1/2}} \right]
%\cdot \left[
%\frac{{\rm exp}\left(-\frac{\epsilon^2(\vect{x})}{2\rho} \right)}
%{(2\pi\rho)^{(N-M)/2}} \right]\\
% & & \\
%& = & P_F(\vect{x}|\Omega) ~\hat{P}_{\bar{F}}(\vect{x}|\Omega) 
%\end{array}
%\label{eq:P_hat}
%\end{equation}
%where $P_F(\vect{x}|\Omega)$ is the true marginal density in $F$-space
%and $\hat{P}_{\bar{F}}(\vect{x}|\Omega)$ is the estimated marginal
%density in the orthogonal complement $\bar{F}$-space.  The optimal
%value of $\rho$ can now be determined by minimizing a suitable cost
%function $J(\rho)$. From an information-theoretic point of view, this
%cost function should be the Kullback-Leibler divergence or {\em
%relative entropy} \cite{Cover94} between the true density
%$P(\vect{x}|\Omega)$ and its estimate $\hat{P}(\vect{x}|\Omega)$
%\begin{equation}
%J(\rho) ~=~ \int P(\vect{x}|\Omega)
%~\log\frac{P(\vect{x}|\Omega)}{\hat{P}(\vect{x}|\Omega)} ~d\vect{x}
%~=~ {\rm E}\left[ \log \frac{P(\vect{x}|\Omega)}{\hat{P}(\vect{x}|\Omega)}\right]
%\end{equation}
%Using the diagonalized forms of the {\em Mahalanobis} distance
%$d(\vect{x})$ and its estimate $\hat{d}(\vect{x})$ and the fact that
%${\rm E}[y_i^2] =\lambda_i$ , it can be easily shown that
%\begin{equation}
%J(\rho) ~=~ \frac{1}{2} \sum_{i=M+1}^{N} \left[ \frac{\lambda_i}{\rho}
%- 1 + \log \frac{\rho}{\lambda_i} \right]
%\end{equation}
%The optimal weight $\rho^*$ can be then found by minimizing this cost
%function with respect to $\rho$. Solving the equation $\frac{\partial
%J}{\partial \rho} = 0$ yields
%\begin{equation}
%\rho^* ~=~ \frac{1}{N-M} \sum_{i=M+1}^{N} \lambda_i
%\label{eq:rho_star}
%\end{equation}
%which is simply the arithmetic average of the eigenvalues in the
%orthogonal subspace $\bar{F}$.
%\footnote{Cootes {\sl et al.} \cite{Cootes94} have used a similar
%decomposition of the Mahalanobis distance but instead use an ad-hoc
%parameter value of $\rho = \frac{1}{2}\lambda_{M+1}$ as an
%approximation.}  In addition to its optimality, $\rho^*$ also results
%in an {\em unbiased} estimate of the {\em Mahalanobis} distance ---
%{\it i.e.,} ${\rm E}[\hat{d}(\vect{x};\rho^*)] = {\rm
%E}[d(\vect{x})]$.  This derivation shows that once we select the
%$M$-dimensional principal subspace $F$ (as indicated, for example, by
%PCA), the optimal estimate of the sufficient statistic
%$\hat{d}(\vect{x})$ will have the form of Eq.(\ref{eq:d_hat}) with
%$\rho$ given by Eq.(\ref{eq:rho_star}).

%It is interesting to consider the minimal cost $J(\rho^*)$ 
%\begin{equation}
%J(\rho^*) ~=~ \frac{1}{2} ~\sum_{i=M+1}^{N} \log
%\frac{\rho^*}{\lambda_i}
%\label{eq:J_star}
%\end{equation}
%from the point of view of the $\bar{F}$-space eigenvalues $\{\lambda_i
%: i=M+1,\cdots,N\}$. It is easy to show that $J(\rho^*)$ is minimized
%when the the $\bar{F}$-space eigenvalues have the {\em least} spread
%about their mean $\rho^*$. This suggests a strategy for selecting the
%principal subspace: choose ${F}$ such that the eigenvalues associated
%with its orthogonal complement $\bar{F}$ have the least absolute
%deviation about their mean. In practice, the higher-order eigenvalues
%typically decay and stabilize near the observation noise variance.
%Therefore this strategy is usually consistent with the standard PCA
%practice of discarding the higher-order components since these tend to
%correspond to the ``flattest'' portion of the eigenvalue spectrum (see
%Figure~\ref{fig:unimodal}(b)). In the limit, as the
%$\bar{F}$-space eigenvalues become exactly equal, the divergence
%$J(\rho^*)$ will be zero and our density estimate
%$\hat{P}(\vect{x}|\Omega)$ approaches the true density
%$P(\vect{x}|\Omega)$.

%We note that in most applications it is customary to simply discard
%the $\bar{F}$-space component and simply work with
%$P_{F}(\vect{x}|\Omega)$. However, the use of the DFFS metric or
%equivalently the marginal density $P_{\bar{F}}(\vect{x}|\Omega)$ is
%critically important in formulating the likelihood of an observation
%$\vect{x}$ --- especially in an object detection task --- since there
%are an infinity of vectors which are {\em not} members of $\Omega$
%which can have likely ${F}$-space projections.  Without
%$P_{\bar{F}}(\vect{x}|\Omega)$ a detection system can result in a
%significant number of false alarms.



%\subsection{Multimodal Densities}
%\label{sec:mixture}

%In the previous section we assumed that probability density of the
%training images was Gaussian. This lead to a likelihood estimate in
%the form of a product of two independent multivariate Gaussian
%distributions (or equivalently the sum of two {\em Mahalanobis}
%distances: DIFS $+$ DFFS). In our experience, the distribution of
%samples in the feature space is often accurately modeled by a single
%Gaussian distribution.  This is especially true in cases where the
%training images are accurately aligned views of similar objects seen
%from a standard view ({\it e.g.,} aligned frontal views of human faces
%at the same scale and orientation). However, when the training set
%represents multiple views or multiple objects under varying
%illumination conditions, the distribution of training views in
%$F$-space is no longer unimodal. In fact the training data tends to
%lie on complex and non-separable low-dimensional manifolds in image
%space. One way to tackle this multimodality is to build a view-based
%(or object-based) formulation where separate eigenspaces are used for
%each view
%\cite{Pentland94}.  Another approach is to capture the complexity of
%these manifolds in a universal or {\em parametric} eigenspace using
%splines
%\cite{Murase95}, or local basis functions \cite{Bregler94}.

%\begin{figure}[tb]
%\begin{center}
%\begin{tabular}{c}
%\psfig{figure=Images/featurespace2.ps,height=2in,clip=t}
%\\
%\end{tabular}
%\end{center}
%\vskip -4mm \caption{Decomposition into the principal
%subspace $F$ and its orthogonal complement $\bar{F}$ for an arbitrary
%density.}
%\label{fig:multimodal}
%\end{figure}

%If we assume that the $\bar{F}$-space components are Gaussian and
%independent of the principal features in $F$ (this would be true in
%the case of pure observation noise in $\bar{F}$) we can still use the
%separable form of the density estimate $\hat{P}(\vect{x}|\Omega)$ in
%Eq.(\ref{eq:P_hat}) where $P_{F}(\vect{x}|\Omega)$ is now an {\em
%arbitrary} density $P(\vect{y})$ in the principal component vector
%$\vect{y}$. Figure~\ref{fig:multimodal} illustrates the decomposition,
%where the DFFS is the residual $\epsilon^2(\vect{x})$ as before.  The
%DIFS, however, is no longer a simple {\em Mahalanobis} distance but
%can nevertheless be interpreted as a ``distance'' by relating it to
%$P(\vect{y})$ --- {\em e.g.,} as ${\rm DIFS} = -\log P(\vect{y})$.


%The density $P(\vect{y})$ can be estimated using a parametric mixture
%model.  Specifically, we can model arbitrarily complex densities using
%a Mixture-of-Gaussians
%\begin{equation}
%P(\vect{y} | \Theta) ~=~ \sum_{i=1}^{N_c} ~\pi_i ~g(\vect{y};
%\vect{\mu}_i, \Sigma_i)
%\end{equation}
%where $g(\vect{y}; \vect{\mu}, \Sigma)$ is an $M$-dimensional Gaussian
%density with mean vector $\vect{\mu}$ and covariance $\Sigma$, and the
%$\pi_i$ are the mixing parameters of the components, satisfying $\sum
%\pi_i = 1$.  The mixture is completely specified by the parameter
%$\Theta =
%\{\pi_i, \vect{\mu}_i,\Sigma_i \}_{i=1}^{Nc}$. Given a training set
%$\{\vect{y}^t\}_{t=1}^{N_T}$ the mixture parameters can be estimated
%using the ML principle
%\begin{equation}
%\Theta^* ~=~ {\rm argmax} ~\left[ \prod_{t=1}^{N_T} ~P(\vect{y}^t |
%\Theta) \right]
%\end{equation}
%This estimation problem is best solved using the
%Expectation-Maximization (EM) algorithm \cite{Dempster77} which
%consists of the following two-step iterative procedure:
%\begin{itemize}
%\item	E-step:
%	\begin{equation}
%	h_i^{k}(t) ~=~ \frac{\pi_i^k g(\vect{y}^t; \mu_i^k,
%\Sigma_i^k)}{\displaystyle \sum_{j=1}^{N_c} \pi_j^k g(\vect{y}^t; \mu_j^k , \Sigma_i^k)}
%	\end{equation}
%\item	M-step:
%	\begin{equation}
%	\pi_i^{k+1} ~=~ \frac{\displaystyle \sum_{t=1}^{N_T}
%h_i^{k}(t)}{\displaystyle \sum_{i=1}^{N_c} \sum_{t=1}^{N_T} h_i^{k}(t)}
%	\end{equation}	
%	\begin{equation}
%	\mu_i^{k+1} ~=~ \frac{\displaystyle \sum_{t=1}^{N_T}
%h_i^{k}(t) \vect{y}^t}{\displaystyle \sum_{t=1}^{N_T} h_i^{k}(t)}
%	\end{equation}	
%	\begin{equation}
%	\Sigma_i^{k+1} ~=~ \frac{\displaystyle \sum_{t=1}^{N_T}
%h_i^{k}(t) (\vect{y}^t - \mu_i^{k+1})(\vect{y}^t -
%\mu_i^{k+1})^T}{\displaystyle \sum_{t=1}^{N_T} h_i^{k}(t)}
%	\end{equation}	
%\end{itemize}
%The E-step computes the {\em a posteriori} probabilities $h_i(t)$
%which are the {\em expectations} of ``missing'' component labels
%$z_i(t)=\{0,1\}$ which denote the membership of $\vect{y}^t$ in the
%$i$-th component. Once these expectations have been computed, the
%M-step maximizes the joint likelihood of the data {\em and} the
%``missing'' variables $z_i(t)$.  The EM algorithm is monotonically
%convergent in {\em likelihood} and is thus guaranteed to find a local
%maximum in the total likelihood of the training set. Further details
%of the EM algorithm for estimation of mixture densities can be found
%in \cite{Redner84}.

%Given our operating assumptions --- that the training data is
%$M$-dimensional (at most) and resides solely in the principal subspace
%$F$ with the exception of perturbations due to white Gaussian
%measurement noise, or equivalently that the $\bar{F}$-space component
%of the data is itself a separable Gaussian density --- the estimate of
%the complete likelihood function $P(\vect{x}|\Omega)$ is given by
%\begin{equation}
%\hat{P}(\vect{x}|\Omega) ~=~ P(\vect{y}|\Theta^*) ~
%\hat{P}_{\bar{F}}(\vect{x}|\Omega)
%\label{P_hat_em}
%\end{equation}
%where $\hat{P}_{\bar{F}}(\vect{x}|\Omega)$ is a Gaussian component density
%based on the DFFS, as before.





\section{Face Detection and Recognition}

In this section we will present several examples of our face detection
and recognition systems including: ML detection of faces and facial
features ({\it e.g.,} eyes) used for facial alignment, recognition
using ``eigenfaces'' on large databases and recognition using the
Bayesian similarity measure.


\begin{figure}[tb]
\begin{center}
\begin{tabular}{cc}
\psfig{figure=Images/db_trains.ps,height=2.5in,clip=t}
&
\psfig{figure=Images/db_detects.ps,height=2.5in,clip=t} \\
(a) & (b)
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) Examples of facial feature training templates and
(b) the resulting typical detections.}
\label{fig:db_features}
\end{figure}

%%\subsection{Faces}
%\label{sec:faces}

Over the years, various strategies for facial feature detection have
been proposed, ranging from edge map projections \cite{Kanade73}, to
more recent techniques using generalized symmetry operators
\cite{Reisfeld90} and multilayer perceptrons \cite{Vincent92}. In any
robust face processing system this task is critically important since
a face must be first geometrically normalized by aligning its features
with those of a stored model before recognition can be attempted.

The eigentemplate approach to the detection of facial features in
``mugshots'' was proposed in \cite{Pentland94}, where the DFFS metric
was shown to be superior to standard template matching for target
detection. The detection task was the estimation of the position of
facial features (the left and right eyes, the tip of the nose and the
center of the mouth) in frontal view photographs of faces at fixed
scale. Figure~\ref{fig:db_features} shows examples of facial feature
training templates and the resulting detections on the MIT Media
Laboratory's database of 7,562 ``mugshots''. 



\begin{figure}[t]
\begin{center}
\begin{tabular}{c}
\raisebox{1.5in}{(a)} \psfig{figure=Images/db_roc_dffs_ml.ps,height=3in,clip=t} \\
\\
\raisebox{0.55in}{(b)} \raisebox{0in}{         } \psfig{figure=Images/thresholds.ps,height=1.5in,clip=t} \\
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) Detection performance of an SSD, DFFS and a ML detector,
(b) geometric interpretation of the detectors.}
\label{fig:db_roc}
\end{figure}


We have compared the detection performance of three different
detectors on approximately 7,000 test images from this database: a
sum-of-square-differences (SSD) detector based on the average facial
feature (in this case the left eye), an eigentemplate or DFFS detector
and a ML detector based on $S(i,j;\Omega)$ as defined in section
\ref{sec:gaussian}.  Figure~\ref{fig:db_roc}(a) shows the {\em
receiver operating characteristic} (ROC) curves for these detectors,
obtained by varying the detection threshold independently for each
detector.  The DFFS and ML detectors were computed based on a
5-dimensional principal subspace.  Since the projection coefficients
were unimodal a Gaussian distribution was used in modeling the true
distribution for the ML detector as in section
\ref{sec:gaussian}.  Note that the ML detector exhibits the best
detection vs. false-alarm tradeoff and yields the highest detection
rate (95\%). Indeed, at the {\em same} detection rate the ML
detector has a false-alarm rate which is nearly 2 orders of magnitude
lower than the SSD.  

Figure~\ref{fig:db_roc}(b) provides the geometric intuition regarding
the operation of these detectors.  The SSD detector's threshold is
based on the {\em radial} distance between the average template (the
origin of this space) and the input pattern.  This leads to
hyperspherical detection regions about the origin.  In contrast, the
DFFS detector measures the orthogonal distance to $F$, thus forming
planar acceptance regions about $F$.  Consequently to accept valid
object patterns in $\Omega$ which are very different from the mean,
the SSD detector must operate with high thresholds which result in
many false alarms.  However, the DFFS detector can not discriminate
between the object class $\Omega$ and non-$\Omega$ patterns in
$F$. The solution is provided by the ML detector which incorporates
both the $\bar{F}$-space component (DFFS) and the $F$-space likelihood
(DIFS). The probabilistic interpretation of Figure~\ref{fig:db_roc}(b)
is as follows: SSD assumes a {\em single} prototype (the mean) in
additive white Gaussian noise whereas the DFFS assumes a {\em uniform}
density in $F$. The ML detector, on the other hand, uses the complete
probability density for detection.


\begin{figure}[tb]
\begin{center}
\begin{tabular}{c}
\psfig{figure=Images/ms_head_search.ps,width=5in}
\end{tabular}
\end{center}
\vskip -4mm \caption{Examples of multiscale face detection.}
\label{fig:ms_head_search}
\end{figure}


We have incorporated and tested the multiscale version of the ML
detection technique in a face detection task.  This multiscale head
finder was tested on the ARPA FERET database where in 97\% of 2,000
images the face, eyes, nose, and mouth were correctly detected and
localized to within one pixel error.  Figure~\ref{fig:ms_head_search} shows
examples of the ML estimate of the position and scale on these images.
The multiscale saliency maps $S(i,j,k ; \Omega)$ were computed based
on the likelihood estimate $\hat{P}(\vect{x}|\Omega)$ in a
10-dimensional principal subspace using a Gaussian model (section
\ref{sec:gaussian}). Note that this detector is able to localize the
position and scale of the head despite variations in hair style and
hair color, as well as presence of sunglasses.  Illumination
invariance was obtained by normalizing the input subimage $\vect{x}$
to a zero-mean unit-norm vector.

\subsection{Using ML detection for attention and alignment}


\begin{figure}[tb]
\begin{center}\hfil
\psfig{figure=Images/face_processor2.ps,width=6.5in,clip=t}\hfil
\end{center}
\vskip -4mm \caption{The face processing system.}
\label{fig:face_processor}
\end{figure}

We have also used the multiscale version of the ML detector as the
{\em attentional} component of an automatic system for recognition and
model-based coding of faces.  The block diagram of this system is
shown in Figure~\ref{fig:face_processor} which consists of a two-stage
object detection and alignment stage, a contrast normalization stage,
and a feature extraction stage whose output is used for both
recognition and coding.  Figure~\ref{fig:white_detect} illustrates the
operation of the detection and alignment stage on a natural test image
containing a human face. The function of the face finder is to locate
regions in the image which have a high likelihood of containing a
face.


\begin{figure}[tb]
\begin{center}
\begin{tabular}{cccc}
\psfig{figure=Images/white.ps,height=1.1in,clip=t} &
\psfig{figure=Images/white_ms_eigsearch.ps,height=1.1in,clip=t} &
\psfig{figure=Images/white_scalewarp.ps,height=1.1in,clip=t} &
\psfig{figure=Images/white_eigsearch.ps,height=1.1in,clip=t} \\
(a) & (b) & (c) & (d)
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) original image, (b) position and scale estimate, (c) normalized
head image, (d) position of facial features.}
\label{fig:white_detect}
\end{figure}

The first step in this process is illustrated in
Figure~\ref{fig:white_detect}(b) where the ML estimate of the position
and scale of the face are indicated by the cross-hairs and bounding
box. Once these regions have been identified, the estimated scale and
position are used to normalize for translation and scale, yielding a
standard ``head-in-the-box'' format image
(Figure~\ref{fig:white_detect}(c)).  A second feature detection stage
operates at this fixed scale to estimate the position of 4 facial
features: the left and right eyes, the tip of the nose and the center
of the mouth (Figure~\ref{fig:white_detect}(d)).  Once the facial
features have been detected, the face image is warped to align the
geometry and shape of the face with that of a canonical model.  Then
the facial region is extracted (by applying a fixed mask) and
subsequently normalized for contrast. The geometrically aligned and
normalized image (shown in Figure~\ref{fig:white_normed}(a)) is then
projected onto a custom set of eigenfaces to obtain a feature vector
which is then used for recognition purposes as well as facial image
coding.  


\begin{figure}[tb]
\begin{center}
\begin{tabular}{ccc}
\psfig{figure=Images/white_normed.ps,height=1.35in,clip=t} &
\psfig{figure=Images/white_recons.ps,height=1.35in,clip=t} &
\psfig{figure=Images/white_jpeg.ps,height=1.35in,clip=t} \\
(a) & (b) & (c)  
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) aligned face, (b) eigenspace reconstruction (85 bytes) (c)
JPEG reconstruction (530 bytes).}
\label{fig:white_normed}
\end{figure}

Figure~\ref{fig:white_normed} shows the normalized facial image
extracted from Figure~\ref{fig:white_detect}(d), its reconstruction
using a 100-dimensional eigenspace representation (requiring only 85
bytes to encode) and a comparable non-parametric reconstruction
obtained using a standard transform-coding approach for image
compression (requiring 530 bytes to encode). This example illustrates
that the eigenface representation used for recognition is also an
effective {\em model-based} representation for data compression. The
first 8 eigenfaces used for this representation are shown in
Figure~\ref{fig:feret_8eigfaces}.


\begin{figure}
\begin{center}
\begin{tabular}{c}
\psfig{figure=Images/feret_8eigfaces.ps,width=4.5in,clip=t}
\end{tabular}
\end{center}
\vskip -4mm \caption{The first 8 eigenfaces.}
\label{fig:feret_8eigfaces}
\end{figure}


\subsection{Eigenface Recognition}

Figure~\ref{fig:feret_photobook} shows the results of a similarity
search in an image database tool called Photobook
\cite{Pentland94_photobook}. Each face in the database was
automatically detected and aligned by the face processing system in
Figure~\ref{fig:face_processor}. The normalized faces were then
projected onto a 100-dimensional eigenspace. The image in the upper
left is the one searched on and the remainder are the ranked nearest
neighbors in the FERET database.  The top three matches in this case
are images of the same person taken a month apart and at different
scales.  The recognition accuracy (defined as the percent correct
rank-one matches) on a database of 155 individuals is 99\%
\cite{Moghaddam94}.


\begin{figure}[tb]
\begin{center}
\begin{tabular}{c}
\psfig{figure=Images/feret_photobook2.ps,width=4.5in,clip=t} 
\end{tabular}
\end{center}
\vskip -4mm \caption{Photobook: FERET face database.}
\label{fig:feret_photobook}
\end{figure}


\subsection{Eigenface Recognition on Large Databases}

In order to have an estimate of the recognition performance on much
larger databases, we have conducted tests on a database of $7,562$
images of approximately $3,000$ people.  The images were collected in
a small booth at a Boston photography show, and include men, women,
and children of all ages and races.  Head position was controlled by
asking people to take their own picture when they were lined up with
the camera. Two LEDs placed at the bottom of holes adjacent to the
camera allowed them to judge their alignment; when they could see both
LEDs then they were correctly aligned.

The eigenfaces for this database were approximated using a principal
components analysis on a representative sample of 128 faces.
Recognition and matching was subsequently performed using the first 20
eigenvectors. 

To assess the average recognition rate, 200 faces were selected at
random, and a nearest-neighbor rule was used to find the most-similar
face from the entire database.  If the most-similar face was of the
same person then a correct recognition was scored. In this experiment
the eigenvector-based recognition system produced a recognition
accuracy of 95\%.  This performance is somewhat surprising because the
database contains wide variations in expression, and has relatively
weak control of head position and illumination.  In a {\em
verification} task, our system yielded a false rejection rate of 1.5\%
at a false acceptance rate of 0.01\%.


\subsection{Recognition using Edge-based Features}

We have also extended the normalized eigenface representation into an
edge-based domain for facial description. We simply run the normalized
facial image through a Canny edge detector to yield an edge map
Unfortunately binary edge maps, are highly uncorrelated with each
other due to their sparse nature, and therefore result in a very
high-dimensional principal subspace.  Therefore, to reduce the
intrinsic dimensionality, we {\em induced} spatial correlation via a
{\em diffusion} process on the binary edge map, which effectively
broadens and ``smears'' the edges, yielding a continuous-valued edge
map
as shown in Figure~\ref{fig:combined_eigfaces}(a). Such an edge
map is simply an alternative representation which imparts mostly {\em
shape} (as opposed to texture) information and has the advantage of
being less susceptible to illumination changes. The recognition rate
of a pure edge-based normalized eigenface representation (on a FERET
database of 155 individuals) was found to be 95\% which is surprising
considering that it utilizes what appears to be (to humans at least) a
rather impoverished representation.  The slight drop in recognition
rate is most likely due to the increased dimensionality of this
representation space and its greater sensitivity to expression
changes, {\it etc.}


\begin{figure}[tb]
\begin{center}
\begin{tabular}{c}
\raisebox{.75in}{(a)} \psfig{figure=Images/dual.ps,height=1.5in,clip=t}\\
\raisebox{.75in}{(b)} \psfig{figure=Images/dual_eig.ps,height=1.5in,clip=t}\\
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) Examples of combined texture/edge-based face
representations and (b) few of the resulting eigenvectors.}
\label{fig:combined_eigfaces}
\end{figure}

Interestingly, we can combine both texture and edge-based
representations of the object by simply performing a KL expansion on
the augmented images shown in Figure~\ref{fig:combined_eigfaces}.  The
resulting eigenvectors conveniently decorrelate the joint
representation and provide a basis set which optimally spans both
domains simultaneously.  With this bimodal representation, the
recognition rate was found to be 97\%. Though still less than a
normalized grayscale representation, we believe a bimodal
representation can have distinct advantages for tasks other than
recognition, such as detection and image interpolation.

\subsection{Bayesian Face Recognition}

To test our Bayesian recognition strategy we used a collection of
images from the FERET face database.  This collection of images
consists of hard recognition cases that have proven difficult for all
face recognition algorithms previously tested on the FERET database.
The difficulty posed by this dataset appears to stem from the fact
that the images were taken at different times, at different locations,
and under different imaging conditions. The set of images consists of
pairs of frontal-views (FA/FB) and are divided into two subsets: the
``gallery'' (training set) and the ``probes'' (testing set). The
gallery images consisted of 74 pairs of images (2 per individual) and
the probe set consisted of 38 pairs of images, corresponding to a
subset of the gallery members. The probe and gallery datasets were
captured a week apart and exhibit differences in clothing, hair and
lighting (see Figure~\ref{fig:dupes}).


\begin{figure}[tb]
\begin{center}
\begin{tabular}{cc}
\psfig{figure=/v/users/baback/Text/CVPR96/Images/dupes_gallery.ps,width=2.5in,clip=t}
&
\psfig{figure=/v/users/baback/Text/CVPR96/Images/dupes_probes.ps,width=2.5in,clip=t}\\
(a) & (b)
\end{tabular}
\end{center}
\vskip -4mm
\caption{Examples of FERET frontal-view image pairs used
for (a) the Gallery set (training) and (b) the Probe set (testing).}
\label{fig:dupes}
\end{figure}

%\Section{Initial Image Alignment}

%\begin{figure}[tb]
%\begin{center}
%\begin{tabular}{c}
%\psfig{figure=/v/users/baback/Text/ECCV96/aligner.ps,width=5in,clip=t}
%\end{tabular}
%\end{center}
%\vskip -4mm 
%\caption{The face alignment system}
%\label{fig:face_processor}
%\end{figure}

%\begin{figure}[tb]
%\begin{center}
%\begin{tabular}{c}
%\psfig{figure=/u/baback/Text/ICCV95/New_Images/feret_8eigfaces.ps,width=5in,clip=t}
%\end{tabular}
%\end{center}
%\vskip -4mm 
%\caption{Standard Eigenfaces.}
%\label{fig:feret_8eigfaces}
%\end{figure}

Before we can apply our matching technique, we need to perform an
affine alignment of these facial images. For this purpose we have used
an automatic face-processing system which extracts faces from the
input image and normalizes for translation, scale as well as slight
rotations (both in-plane and out-of-plane). This is achieved using the
maximum-likelihood detection and alignment method described earlier
which is summarized in Figure~\ref{fig:face_processor}. All the faces
in our experiments were geometrically aligned and normalized in this
manner prior to further analysis.




\begin{figure*}[tbh]
\begin{center}
\begin{tabular}{cc}
\psfig{figure=/v/users/baback/Text/CVPR96/Images/intra_extra_dist.ps,width=3in,clip=t}
&
\raisebox{.35in}{\psfig{figure=/v/users/baback/Text/FG98/Images/new_distribs.ps,width=3in,clip=t}}
\\ (a) & (b)
\end{tabular}
\end{center}
\vskip -4mm \caption{(a) Distribution of the two classes in the first
3 principal components (circles for $\Omega_I$, dots for $\Omega_E$)
and (b) schematic representation of the two distributions showing
orientation difference between the corresponding principal
eigenvectors.}
\label{fig:distribs}
\end{figure*}


\subsubsection{Comparison with Eigenface Matching}

As a baseline comparison, we first used an eigenface matching
technique for recognition \cite{Turk91}. The normalized images from
the gallery and the probe sets were projected onto a 100-dimensional
eigenspace and a nearest-neighbor rule based on a Euclidean distance
measure was used to match each probe image to a gallery image. We note
that this method corresponds to a generalized template-matching method
which uses a Euclidean norm type of similarity $S(I_1,I_2)$, which is
restricted to the principal component subspace of the data. A few of
the lower-order eigenfaces used for this projection are shown in
Figure~\ref{fig:feret_8eigfaces}. We note that these eigenfaces represent
the principal components of an entirely different set of images ---
{\it i.e.,} none of the individuals in the gallery or probe sets were
used in obtaining these eigenvectors. In other words, neither the
gallery nor the probe sets were part of the ``training set.''  The
rank-1 recognition rate obtained with this method was found to be 84\%
(64 correct matches out of 76), and the correct match was always in
the top 10 nearest neighbors.  Note that this performance is better
than or similar to recognition rates obtained by any algorithm tested
on this database, and that it is lower (by about 10\%) than the
typical rates that we have obtained with the FERET database
\cite{Moghaddam94}.  We attribute this lower performance to the fact
that these images were selected to be particularly challenging. In
fact, using an eigenface method to match the first views of the 76
individuals in the gallery to their second views, we obtain a higher
recognition rate of 89\% (68 out of 76), suggesting that the gallery
images represent a less challenging data set since these images were
taken at the same time and under identical lighting conditions.



\subsubsection{Intra/Extrapersonal-Based Matching}

For our probabilistic algorithm, we first gathered training data by
computing the intensity differences for a training subset of 74
intrapersonal differences (by matching the two views of every
individual in the gallery) and a random subset of 296 extrapersonal
differences (by matching images of {\em different} individuals in the
gallery), corresponding to the classes $\Omega_I$ and $\Omega_E$,
respectively.
%An example of each of these two types of differences is shown
%in Figure~\ref{fig:intra_extra_warps}.

It is interesting to consider how these two classes are distributed,
for example, are they linearly separable or embedded distributions?
One simple method of visualizing this is to plot their mutual
principal components --- {\it i.e.,} perform PCA on the {\em combined}
dataset and project each vector onto the principal eigenvectors. Such
a visualization is shown in Figure~\ref{fig:distribs}(a) which is a 3D
scatter plot of the first 3 principal components. This plot shows what
appears to be two completely enmeshed distributions, both having
near-zero means and differing primarily in the amount of scatter, with
$\Omega_I$ displaying smaller intensity differences as expected.  It
therefore appears that one can not reliably distinguish low-amplitude
extrapersonal differences (of which there are many) from intrapersonal
ones.

However, direct visual interpretation of Figure~\ref{fig:distribs}(a)
is very misleading since we are essentially dealing with
low-dimensional (or ``flattened'') hyper-ellipsoids which are
intersecting near the origin of a very high-dimensional space. The key
distinguishing factor between the two distributions is their relative
orientation. Fortunately, we can easily determine this relative
orientation by performing a separate PCA on each class and computing
the dot product of their respective first eigenvectors. This analysis
yields the cosine of the angle between the major axes of the two
hyper-ellipsoids, which was found to be $124\dg$, implying that the
orientation of the two hyper-ellipsoids is quite
different. Figure~\ref{fig:distribs}(b) is a schematic illustration of
the geometry of this configuration, where the hyper-ellipsoids have
been drawn to approximate scale using the corresponding eigenvalues.

\begin{figure}[tb]
\begin{center}
\begin{tabular}{c}
\psfig{figure=/u/baback/Text/FG98/Images/idiff_intra8.ps,width=5in,clip=t}
\\
(a)\\
\psfig{figure=/u/baback/Text/FG98/Images/idiff_extra8.ps,width=5in,clip=t}\\
(b)
\end{tabular}
\end{center}
\vskip -4mm
\caption{``Dual'' Eigenfaces: (a) Intrapersonal, (b) Extrapersonal}
\label{fig:idiff_intra_extra}
\end{figure}


\subsubsection{Dual Eigenfaces}

We note that the two mutually exclusive classes $\Omega_I$ and
$\Omega_E$ correspond to a ``dual'' set of eigenfaces as shown in
Figure~\ref{fig:idiff_intra_extra}. Note that the intrapersonal
variations shown in Figure~\ref{fig:idiff_intra_extra}-(a) represent
subtle variations due mostly to expression changes (and lighting)
whereas the extrapersonal variations in
Figure~\ref{fig:idiff_intra_extra}-(b) are more representative of
general eigenfaces which code variations such as hair color, facial
hair and glasses. Also note the overall qualitative similarity of the
extrapersonal eigenfaces to the standard eigenfaces in
Figure~\ref{fig:feret_8eigfaces}. This suggests the basic intuition that
intensity differences of the extrapersonal type span a larger vector
space similar to the volume of facespace spanned by standard
eigenfaces, whereas the {\em intrapersonal} eigenspace corresponds to
a more tightly constrained subspace. It is the representation of this
intrapersonal subspace that is the critical part of formulating a
probabilistic measure of facial similarity. In fact our experiments
with a larger set of FERET images have shown that this intrapersonal
eigenspace alone is sufficient for a simplified {\em maximum
likelihood} measure of similarity (see Section~\ref{sec:feret96}).


Finally, we note that since these classes are not linearly separable,
simple linear discriminant techniques ({\it e.g.,} using hyperplanes)
can not be used with any degree of reliability. The proper decision
surface is inherently nonlinear (quadratic, in fact, under the
Gaussian assumption) and is best defined in terms of the {\em a
posteriori} probabilities --- {\it i.e.,} by the equality $P(\Omega_I
| \Delta) = P(\Omega_E | {\Delta})$. Fortunately, the optimal
discriminant surface is automatically implemented when invoking a MAP
classification rule.


Having analyzed the geometry of the two distributions, we then
computed the likelihood estimates $P({\Delta}|\Omega_I)$ and
$P({\Delta}|\Omega_E)$ using the PCA-based method outlined in
Section~\ref{sec:density}. We selected principal subspace dimensions
of $M_I = 10$ and $M_E = 30$ for $\Omega_I$ and $\Omega_E$,
respectively.  These density estimates were then used with a default
setting of equal priors, $P(\Omega_I) = P(\Omega_E)$, to evaluate the
{\em a posteriori} intrapersonal probability $P(\Omega_I | \Delta)$
for matching probe images to those in the gallery. 

Therefore, for each probe image we computed probe-to-gallery
differences and sorted the matching order, this time using the {\em a
posteriori} probability $P(\Omega_I | {\Delta})$ as the
similarity measure. This probabilistic ranking yielded an improved
rank-1 recognition rate of 89.5\%. Furthermore, out of the 608
extrapersonal warps performed in this recognition experiment, only 2\%
(11) were misclassified as being intrapersonal --- {\it i.e.,} with
$P(\Omega_I | {\Delta}) > P(\Omega_E | {\Delta})$.


\subsubsection{The 1996 FERET  Competition Results}
\label{sec:feret96}


\begin{figure}[t]
\vspace{-10mm}
\begin{center}
\begin{tabular}{c}
\hspace{-9mm}\psfig{figure=/v/users/baback/Text/FG98/Images/fafb96.eps,width=5in,clip=t}
\\
\end{tabular}
\end{center}
\vskip -6mm 
\caption{Cumulative recognition rates for frontal FA/FB views for the
competing algorithms in the FERET 1996 test. The top curve (labeled
``MIT Sep 96'') corresponds to our Bayesian matching technique. Note
that second placed is standard eigenface matching (labeled ``MIT Mar
95'').}
\label{fig:feret96}
\end{figure}



The inter/extrapersonal approach to recognition has produced a
significant improvement over the accuracy we obtained using a standard
eigenface nearest-neighbor matching rule.  The probabilistic
similarity measure was used in the September 1996 FERET competition
(with subspace dimensionalities of $M_I = M_E = 125$) and was found to
be the top-performing system by a typical margin of 10-20\% over the
other competing algorithms \cite{Phillips_cvpr97} (see
Figure~\ref{fig:feret96}). Figure~\ref{fig:feret96_fafb} shows the
performance comparison between standard eigenfaces and the Bayesian
method from this test. Note the 10\% gain in performance afforded by
the new Bayesian similarity measure.
%  Similarly,
%Figure~\ref{fig:feret96_dups} shows the recognition results for
%``duplicate'' images which were separated in time by up to 6 months (a
%much more challenging recognition problem) which shows a 30\%
%improvement in recognition rate with Bayesian matching. 
Thus we note
that 
%in both cases (FA/FB and duplicates)
the new probabilistic
similarity measure has effectively {\em halved} the error rate of
eigenface matching.


We have recently experimented with a more simplified probabilistic
similarity measure which uses only the {\em intrapersonal} eigenfaces
with the intensity difference $\Delta$ to formulate a {\em maximum
likelihood} (ML) matching technique using
\begin{equation}
S' ~ = ~ P(\Delta  |  \Omega_I)
\label{eq:simple_similarity}
\end{equation}
instead of the {\em maximum a posteriori} (MAP) approach defined by
Equation~\ref{eq:similarity}. Although this simplified measure has not
yet been officially FERET tested, our own experiments with a database
of size 2000 have shown that using $S'$ instead of $S$ results in only
a minor (2\%) deficit in the recognition rate while cutting the
computational cost by a factor of 1/2 (requiring a single eigenspace
projection as opposed to two).

\begin{figure}[t]
\vspace{-10mm}
\begin{center}
\begin{tabular}{c}
\hspace{-12mm}\psfig{figure=/v/users/baback/Text/FG98/Images/MIT9596fafb.eps,width=5in,clip=t}
\end{tabular}
\end{center}
\vskip -8mm 
\caption{Cumulative recognition rates for frontal FA/FB views with
standard eigenface matching and the newer Bayesian similarity metric.}
\label{fig:feret96_fafb}
\end{figure}

%\begin{figure}[t]
%\vspace{-10mm}
%\begin{center}
%\begin{tabular}{c}
%\hspace{-12mm}\psfig{figure=/v/users/baback/Text/FG98/Images/MIT9596dups.eps,width=5in,clip=t}
%\end{tabular}
%\end{center}
%\vskip -8mm 
%\caption{Cumulative recognition rates for frontal duplicate views with
%standard eigenface matching and the newer Bayesian similarity metric.}
%\label{fig:feret96_dups}
%\end{figure}






\section{View-Based Face Recognition}


The problem of face recognition under general viewing conditions
(change in pose) can also be approached using an eigenspace
formulation. There are essentially two ways of approaching this
problem using an eigenspace framework. Given $N$ individuals under $M$
different views, one can do recognition and pose estimation in a
universal eigenspace computed from the combination of $NM$ images. In
this way a single ``parametric eigenspace'' will encode both identity
as well as pose.  Such an approach, for example, has recently been
used by Murase and Nayar \cite{Murase95} for general 3D object
recognition.

Alternatively, given $N$ individuals under $M$ different views, we can
build a ``view-based'' set of $M$ distinct eigenspaces, each capturing
the variation of the $N$ individuals in a common view.  The view-based
eigenspace is essentially an extension of the eigenface technique to
multiple sets of eigenvectors, one for each combination of scale and
orientation.  One can view this architecture as a set of parallel
``observers'' each trying to explain the image data with their set of
eigenvectors (see also Darrell and Pentland \cite{Darrell93}).
In this view-based, multiple-observer approach, the first step is to
determine the location and orientation of the target object by
selecting the eigenspace which best describes the input image.  This
can be accomplished by calculating the likelihood estimate using each
viewspace's eigenvectors and then selecting the maximum.



\begin{figure}[t]
\centerline{\hfill
\psfig{figure=Images/wh_4.ps,width=5in,clip=t}\hfill}
\vskip -4mm \caption[]{Some of the images used to test accuracy at
face recognition despite wide variations in head orientation.  Average
recognition accuracy was 92\%, the orientation error had a standard deviation of $15^{\circ}$.}
\label{fig:orient}
\end{figure}


The main advantage of the parametric eigenspace method is its
simplicity. The encoding of an input image using $n$ eigenvectors
requires only $n$ projections.  In the view-based method, $M$
different sets of $n$ projections are required, one for each view.
However, this does not imply that a factor of $M$ times more
computation is necessarily required.  By progressively calculating the
eigenvector coefficients while pruning alternative viewspaces, the
cost of using $M$ eigenspaces can be greatly reduced.

The key difference between the view-based and parametric
representations can be understood by considering the geometry of
facespace. In the high-dimensional vector space of an input image,
multiple-orientation training images are represented by a set of $M$
distinct regions, each defined by the scatter of $N$ individuals.
Multiple views of a face form non-convex (yet connected) regions in
image space \cite{Bichsel94}. Therefore the resulting ensemble is a
highly complex and non-separable manifold.


\begin{figure}[t]
\begin{center}
\begin{tabular}{cc}
\begin{tabular}{c}
Training View\\
\psfig{figure=Images/orient_train_recon.ps,width=2.5in,clip=t}
\\
Testing View \\
\psfig{figure=Images/orient_novel_recon.ps,width=2.5in,clip=t}
\end{tabular} &
\raisebox{-1.1in}{\psfig{figure=Images/viewbased.ps,width=2.5in,clip=t}} \\
(a) & (b)
\end{tabular}
\end{center}
\vskip -4mm \caption[]{(a) parametric vs. view-based eigenspace reconstructions for a
training view and a novel testing view. The input image is shown in
the left column. The middle and right columns correspond to the
parametric and view-based reconstructions, respectively. All
reconstructions were computed using the first 10 eigenvectors. (b) a
schematic representation of the two approaches.}
\label{fig:orient_recons}
\end{figure}

The parametric eigenspace attempts to describe this ensemble with a
projection onto a single low-dimensional linear subspace
(corresponding to the first $n$ eigenvectors of the $NM$ training
images). In contrast, the view-based approach corresponds to $M$
independent subspaces, each describing a particular region of the
facespace (corresponding to a particular view of a face). The relevant
analogy here is that of modeling a complex distribution by a single
cluster model or by the union of several component clusters.
Naturally, the latter (view-based) representation can yield a more
accurate representation of the underlying geometry.


This difference in representation becomes evident when considering the
quality of reconstructed images using the two different methods.
Figure~\ref{fig:orient_recons} compares reconstructions obtained with
the two methods when trained on images of faces at multiple
orientations. In Figure~\ref{fig:orient_recons}(a) we see first an
image in the training set, followed by reconstructions of this image
using first the parametric eigenspace and then the view-based
eigenspace.  Note that in the parametric reconstruction neither the
pose nor the identity of the individual is adequately captured. The
view-based reconstruction, on the other hand, provides a much better
characterization of the object. Similarly, in
Figure~\ref{fig:orient_recons}(b) we see a novel view ($+68\dg$) with
respect to the training set ($-90\dg$ to $+45\dg$).  Here, both
reconstructions correspond to the nearest view in the training set
($+45\dg$) but the view-based reconstruction is seen to be more
representative of the individual's identity.  Although the quality of
the reconstruction is not a direct indicator of the recognition power,
from an information-theoretic point-of-view the multiple eigenspace
representation is a more accurate representation of the signal
content.  
%In the following sections we will experimentally show that
%this does lead to an advantage in terms of recognition performance.

We have evaluated the view-based approach with data similar to that
shown in Figure \ref{fig:orient}.  This data consists of 189 images
consisting of nine views of 21 people.  The nine views of each person
were evenly spaced from $-90^{\circ}$ to $+90^{\circ}$ along the
horizontal plane.
In the first series of experiments the {\em interpolation} performance
was tested by training on a subset of the available views \{$\pm 90\dg,
\pm 45\dg, 0\dg$\} and testing on the intermediate views
\{$\pm 68\dg, \pm23\dg$\}.  A 90\% average recognition rate was
obtained.
A second series of experiments tested the {\em extrapolation}
performance by training on a range of views ({\it e.g.,} $-90\dg$ to
$+45\dg$) and testing on novel views outside the training range ({\it
e.g.,} $+68\dg$ and $+90\dg$). For testing views separated by $\pm
23\dg$ from the training range, the average recognition rates were
83\%.  For $\pm 45\dg$ testing views, the average recognition rates
were 50\% (see \cite{Pentland94} for further details).



\begin{figure}[t]
\begin{center}
\begin{tabular}{cc}
\raisebox{.2in}{\psfig{figure=Images/db_features_single.ps,width=1in,clip=t}} &
\raisebox{.0in}{\psfig{figure=Images/wf_recogrates.ps,height=3.4in,clip=t}}
\\
(a) & (b)
\end{tabular}
\end{center}
\vskip -4mm \caption[]{(a) facial eigenfeature regions, (b) recognition rates for
eigenfaces, eigenfeatures and the combined modular representation.}
\label{fig:wf_recogrates}
\end{figure}

\section{Modular Descriptions for Recognition}

The eigenface recognition method is easily extended to facial features
as shown in Figure~\ref{fig:wf_recogrates}(a).  Eye-movement studies
indicate that these particular facial features represent important
landmarks for fixation, especially in an attentive discrimination
task~\cite{Yarbus67}.
%  Therefore 
This leads to an improvement in recognition performance by
incorporating an additional layer of description in terms of facial
features. This can be viewed as either a modular or layered
representation of a face, where a coarse (low-resolution) description
of the whole head is augmented by additional (higher-resolution)
details in terms of salient facial features.

%Pentland {\em et al.} \cite{Pentland94} demonstrated a robust
%technique for automatically detecting important facial landmarks (a
%probabilistic interpretation of this technique along with an improved
%version is presented later in this chapter). With such an ability, we
%can automatically generate a modular representation of a face. 

The utility of this layered representation (eigenface plus
eigenfeatures) was tested on a small subset of our large face
database.  We selected a representative sample of 45 individuals with
two views per person, corresponding to different facial expressions
(neutral vs.  smiling).  These set of images was partitioned into a
training set (neutral) and a testing set (smiling).  Since the
difference between these particular facial expressions is primarily
articulated in the mouth, this feature was discarded for recognition
purposes.

Figure~\ref{fig:wf_recogrates}(b) shows the recognition rates as a
function of the number of eigenvectors for eigenface-only,
eigenfeature-only and the combined representation.  What is surprising
is that (for this small dataset at least) the eigenfeatures alone were
sufficient in achieving an (asymptotic) recognition rate of 95\%
(equal to that of the eigenfaces). More surprising, perhaps, is the
observation that in the lower dimensions of eigenspace, eigenfeatures
outperformed the eigenface recognition.  Finally, by using the
combined representation, we gain a slight improvement in the
asymptotic recognition rate (98\%). A similar effect was reported by
Brunelli and Poggio \cite{Brunelli93} where the cumulative normalized
correlation scores of templates for the face, eyes, nose and mouth
showed improved performance over the face-only templates.



A potential advantage of the eigenfeature layer is the ability to
overcome the shortcomings of the standard eigenface method. A pure
eigenface recognition system can be fooled by gross variations in the
input image (hats, beards, etc.).  Figure~\ref{fig:wf_retrieval}(a)
shows additional testing views of 3 individuals in the above dataset
of 45.  These test images are indicative of the type of variations
which can lead to false matches: a hand near the face, a painted face,
and a beard.  Figure~\ref{fig:wf_retrieval}(b) shows the nearest
matches found based on standard eigenface matching. Neither of the 3
matches correspond to the correct individual.  On the other hand,
Figure~\ref{fig:wf_retrieval}(c) shows the nearest matches based on
the eyes and nose, and results in correct identification in each
case. This simple example illustrates the potential advantage of a
modular representation in disambiguating low-confidence eigenface
matches.



\begin{figure}[t]
\begin{center}
\begin{tabular}{c}
\raisebox{.5in}{(a)} \psfig{figure=Images/wf_retrieval_test.ps,height=1in,clip=t} \\
\raisebox{.5in}{(b)} \psfig{figure=Images/wf_retrieval_w.ps,height=1in,clip=t} \\
\raisebox{.5in}{(c)} \psfig{figure=Images/wf_retrieval_f.ps,height=1in,clip=t} 
\end{tabular}
\end{center}
\vskip -4mm \caption[]{(a) Test views, (b) Eigenface matches, (c) Eigenfeature matches.}
\label{fig:wf_retrieval}
\end{figure}






%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{Discussion}

In this chapter we have described an eigenspace density estimation
technique for unsupervised visual learning which exploits the {\em
intrinsic} low-dimensionality of the training imagery to form a
computationally simple estimator for the complete likelihood function
of the object. Our estimator is based on a subspace decomposition and
can be evaluated using only the $M$-dimensional principal component
vector.  We derived the form for an optimal estimator and its
associated expected cost for the case of a Gaussian density.  In
contrast to previous work on learning and characterization --- which
uses PCA primarily for dimensionality reduction and/or feature
extraction --- our method uses the eigenspace decomposition as an
integral part of estimating {\em complete} density functions in
high-dimensional image spaces. 

These density estimates were then used in a maximum likelihood
formulation for target detection.  The multiscale version of this
detection strategy was demonstrated in applications in which it
functioned as an attentional subsystem for object recognition. The
performance was found to be superior to existing detection techniques
in experiments with large numbers of test data.  We have also shown
that the same representation can be extended to multiple head poses,
to incorporate edge or texture features, and to utilize facial
features such as eye or nose shape.  Each of these extensions has
provided additional robustness and generality to the core idea of
detection and recognition using probabilistic appearance models.

%We note that from a probabilistic perspective, the class conditional
%density $P(\vect{x}|\Omega)$ is the most important object
%representation to be learned. This density is the critical component
%detection, recognition, prediction, interpolation and general
%pinference.  For example, having learned these densities for several
%object classes $\{\Omega_1,\Omega_2,\cdots, \Omega_n\}$, one can
%invoke a Bayesian framework for classification and recognition:
%\begin{equation}
%P(\Omega_i | \vect{x}) ~=~ \frac{P(\vect{x}|\Omega_i) P(\Omega_i)}
%{\displaystyle \sum_{j=1}^{n} ~P(\vect{x}|\Omega_j) P(\Omega_j)}
%\end{equation} 
%where now a maximum {\em a posteriori} (MAP) classification rule can
%be used for object/pose identification.

%Such a framework is also important in detection. In fact, the ML
%detection framework can be extended using the notion of a
%``not-class'' $\bar{\Omega}$, resulting in {\em a posteriori} saliency
%maps of the form
%\begin{equation}
%S(i,j,k; \Omega) ~=~ P(\Omega | \vect{x}) ~=~ \frac{P(\vect{x}|\Omega) P(\Omega)}
%{P(\vect{x}|\bar{\Omega}) P(\bar{\Omega}) ~+~ P(\vect{x}|\Omega)
%P(\Omega)}
%\end{equation} 
%where now a maximum {\em a posteriori} (MAP) rule can be used to
%estimate the position and scale of the object. One difficulty with
%such a formulation is that the ``not-class'' $\bar{\Omega}$ is, in
%practice, too broad a category and is therefore multimodal and very
%high-dimensional. One possible approach to this problem is to use ML
%detection to identify the particular subclass of $\bar{\Omega}$ which
%has high likelihoods ({\it e.g.,} false alarms) and then to estimate
%this distribution and use it in the MAP framework. This can be viewed
%as a probabilistic approach to learning using positive as well as {\em
%negative} examples. 

%In fact, such a MAP framework can be viewed as a Bayesian formulation
%of some neural network approaches to target detection.  Perhaps the
%most closely related is the neural network face detector of Sung \&
%Poggio \cite{Sung94} which is essentially a trainable nonlinear binary
%pattern classifier.  They too learn the distribution of the object
%class with a Mixture-of-Gaussians model (using an elliptical k-means
%algorithm instead of EM). Instead of likelihoods, however, input
%patterns are represented by a set of distances to each mixture
%component (similar to a combination of the DIFS and DFFS), thus
%forming a feature vector indicative of the overall class
%membership. In addition, Sung \& Poggio explicitly model the
%``not-class'' by learning the distribution of nearby {\em non-face}
%patterns. The set of distances to both classes are then used to train
%a neural network to discriminate between face and non-face patterns
%(similar to computing a likelihood ratio in MAP).  Another recent
%example of a neural network technique for object detection which also
%utilizes negative examples is the face-finder system of Rowley {\it et
%al.} \cite{Rowley95}. The experimental results obtained with these
%methods clearly demonstrate the need for incorporating negative
%examples in building robust detection systems.






%\footnotesize{

\begin{singlespace}

\bibliography{Biblio}
\bibliographystyle{plain}
%\bibliographystyle{latex8}

\end{singlespace}

% }






\end{document}




