Multi-sensor multi-target tracking requires the solution to a challenging data association problem. The problem simplifies when a portion of the target state vector and the corresponding sensor data satisfy a particular Markovian assumption. This leads to quantifiable benefits in performance vs. complexity of the tracking solution. This paper summarizes recently-obtained technical advances in graph-based tracking and applies this to a benchmark study with respect to an advanced track-oriented multiple-hypothesis tracking solution.

Bayesian models provide powerful tools for analyzing complex time series data, but performing inference with large datasets is a challenge. Stochastic variational inference (SVI) provides a new framework for approximating model posteriors with only a small number of passes through the data, enabling such models to be fit at scale. However, its application to time series models has not been studied. In this paper we develop SVI algorithms for several common Bayesian time series models, namely the hidden Markov model (HMM), hidden semi-Markov model (HSMM), and the nonparametric HDP-HMM and HDP-HSMM. In addition, because HSMM inference can be expensive even in the minibatch setting of SVI, we develop fast approximate updates for HSMMs with durations distributions that are negative binomials or mixtures of negative binomials.

The problem of efficiently drawing samples from a Gaussian graphical model or Gaussian Markov random field is studied. We introduce the subgraph perturbation sampling algorithm, which makes use of any pre-existing tractable inference algorithm for a subgraph by perturbing this algorithm so as to yield asymptotically exact samples for the intended distribution. The subgraph can have any structure for which efficient inference algorithms exist: for example, tree-structured, low tree-width, or having a small feedback vertex set. The experimental results demonstrate that this subgraph perturbation algorithm efficiently yields accurate samples for many graph topologies.

For inference in Gaussian graphical models with cycles, loopy belief propagation (LBP) performs well for some graphs, but often diverges or has slow convergence. When LBP does converge, the variance estimates are incorrect in general. The feedback message passing (FMP) algorithm has been proposed to enhance the convergence and accuracy of inference. In FMP, standard LBP is run twice on the subgraph excluding the pseudo-FVS (a set of nodes that breaks most crucial cycles) while nodes in the pseudo-FVS use a different protocol. In this paper, we propose recursive FMP, a purely distributed extension of FMP, where all nodes use the same message-passing protocol. An inference problem on the entire graph is recursively reduced to those on smaller subgraphs in a distributed manner. One advantage of this recursive approach compared with FMP is that there is only one active feedback node at a time, so centralized communication among feedback nodes can be turned into message broadcasting from the single feedback node. We characterize this algorithm using walk-sum analysis and provide theoretical results for convergence and accuracy. We also demonstrate the performance using both simulated models on grids and large-scale sea surface height anomaly data.

Sampling inference methods are computationally difficult to scale for many models in part because global dependencies can reduce opportunities for parallel computation. Without strict conditional independence structure among variables, standard Gibbs sampling theory requires sample updates to be performed sequentially, even if dependence between most variables is not strong. Empirical work has shown that some models can be sampled effectively by going Hogwild'' and simply running Gibbs updates in parallel with only periodic global communication, but the successes and limitations of such a strategy are not well understood. As a step towards such an understanding, we study the Hogwild Gibbs sampling strategy in the context of Gaussian distributions. We develop a framework which provides convergence conditions and error bounds along with simple proofs and connections to methods in numerical linear algebra. In particular, we show that if the Gaussian precision matrix is generalized diagonally dominant, then any Hogwild Gibbs sampler, with any update schedule or allocation of variables to processors, yields a stable sampling process with the correct sample mean. "

Gaussian Graphical Models (GGMs) or Gauss Markov random fields are widely used in many applications, and the trade-off between the modeling capacity and the efficiency of learning and inference has been an important research problem. In this paper, we study the family of GGMs with small feedback vertex sets (FVSs), where an FVS is a set of nodes whose removal breaks all the cycles. Exact inference such as computing the marginal distributions and the partition function has complexity $O(k^{2}n)$ using message-passing algorithms, where k is the size of the FVS, and n is the total number of nodes. We propose efficient structure learning algorithms for two cases: 1) All nodes are observed, which is useful in modeling social or flight networks where the FVS nodes often correspond to a small number of high-degree nodes, or hubs, while the rest of the networks is modeled by a tree. Regardless of the maximum degree, without knowing the full graph structure, we can exactly compute the maximum likelihood estimate in $O(kn^2+n^2\log n)$ if the FVS is known or in polynomial time if the FVS is unknown but has bounded size. 2) The FVS nodes are latent variables, where structure learning is equivalent to decomposing a inverse covariance matrix (exactly or approximately) into the sum of a tree-structured matrix and a low-rank matrix. By incorporating efficient inference into the learning steps, we can obtain a learning algorithm using alternating low-rank correction with complexity $O(kn^{2}+n^{2}\log n)$ per iteration. We also perform experiments using both synthetic data as well as real data of flight delays to demonstrate the modeling capacity with FVSs of various sizes. We show that empirically the family of GGMs of size $O(\log n)$ strikes a good balance between the modeling capacity and the efficiency.

There has been a growing interest in exploiting contextual information in addition to local features to detect and localize multiple object categories in an image. A context model can rule out some unlikely combinations or locations of objects and guide detectors to produce a semantically coherent interpretation of a scene. However, the performance benefit of context models has been limited because most of the previous methods were tested on data sets with only a few object categories, in which most images contain one or two object categories. In this paper, we introduce a new data set with images that contain many instances of different object categories, and propose an efficient model that captures the contextual information among more than a hundred object categories using a tree structure. Our model incorporates global image features, dependencies between object categories, and outputs of local detectors into one probabilistic framework. We demonstrate that our context model improves object recognition performance and provides a coherent interpretation of a scene, which enables a reliable image querying system by multiple object categories. In addition, our model can be applied to scene understanding tasks that local detectors alone cannot solve, such as detecting objects out of context or querying for the most typical and the least typical scenes in a data set.

The problem of designing policies for in-network function computation with minimum energy consumption subject to a latency constraint is considered. The scaling behavior of the energy consumption under the latency constraint is analyzed for random networks, where the nodes are uniformly placed in growing regions and the number of nodes goes to infinity. The special case of sum function computation and its delivery to a designated root node is considered first. A policy which achieves order-optimal average energy consumption in random networks subject to the given latency constraint is proposed. The scaling behavior of the optimal energy consumption depends on the path-loss exponent of wireless transmissions and the dimension of the Euclidean region where the nodes are placed. The policy is then extended to computation of a general class of functions which decompose according to maximal cliques of a proximity graph such as the k-nearest neighbor graph or the geometric random graph. The modified policy achieves order-optimal energy consumption albeit for a limited range of latency constraints.

We consider the problem of Ising and Gaussian graphical model selection given n i.i.d. samples from the model. We propose an efficient threshold-based algorithm for structure estimation based known as conditional mutual information test. This simple local algorithm requires only low-order statistics of the data and decides whether two nodes are neighbors in the unknown graph. Under some transparent assumptions, we establish that the proposed algorithm is structurally consistent (or sparsistent) when the number of samples scales as n= Omega(J_{min}^{-4} log p), where p is the number of nodes and J_{min} is the minimum edge potential. We also prove novel non-asymptotic necessary conditions for graphical model selection.

Loopy belief propagation has been employed in a wide variety of applications with great empirical success, but it comes with few theoretical guarantees. In this paper, we analyze the performance of the max-product form of belief propagation for the weighted matching problem on general graphs. We show that the performance of max-product is exactly characterized by the natural linear programming (LP) relaxation of the problem. In particular, we first show that if the LP relaxation has no fractional optima then max-product always converges to the correct answer. This establishes the extension of the recent result by Bayati, Shah and Sharma, which considered bipartite graphs, to general graphs. Perhaps more interestingly, we also establish a tight converse, namely that the presence of any fractional LP optimum implies that max-product will fail to yield useful estimates on some of the edges. We extend our results to the weightedb-matching andr-edge-cover problems. We also demonstrate how to simplify the max-product message-update equations for weighted matching, making it easily deployable in distributed settings like wireless or sensor networks.

The problem of maximum-likelihood (ML) estimation of discrete tree-structured distributions is considered. Chow and Liu established that ML-estimation reduces to the construction of a maximum-weight spanning tree using the empirical mutual information quantities as the edge weights. Using the theory of large-deviations, we analyze the exponent associated with the error probability of the event that the ML-estimate of the Markov tree structure differs from the true tree structure, given a set of independently drawn samples. By exploiting the fact that the output of ML-estimation is a tree, we establish that the error exponent is equal to the exponential rate of decay of a single dominant crossover event. We prove that in this dominant crossover event, a non-neighbor node pair replaces a true edge of the distribution that is along the path of edges in the true tree graph connecting the nodes in the non-neighbor pair. Using ideas from Euclidean information theory, we then analyze the scenario of ML-estimation in the very noisy learning regime and show that the error exponent can be approximated as a ratio, which is interpreted as the signal-to-noise ratio (SNR) for learning tree distributions. We show via numerical experiments that in this regime, our SNR approximation is accurate.

There has been a growing interest in exploiting contextual information in addition to local features to detect and localize multiple object categories in an image. Context models can efficiently rule out some unlikely combinations or locations of objects and guide detectors to produce a semantically coherent interpretation of a scene. However, the performance benefit from using context models has been limited because most of these methods were tested on datasets with only a few object categories, in which most images contain only one or two object categories. In this paper, we introduce a new dataset with images that contain many instances of different object categories and propose an efficient model that captures the contextual information among more than a hundred of object categories. We show that our context model can be applied to scene understanding tasks that local detectors alone cannot solve.

We consider spatial graphical models on random Euclidean points, applicable for data in sensor and social networks. We establish limit laws for general functions of the graphical model such as the mean value, the entropy rate etc. as the number of nodes goes to infinity under certain conditions. These conditions require the corresponding Gibbs measure to be spatially mixing and for the random graph of the model to satisfy a certain localization property known as stabilization. Graphs such the k nearest neighbor graph and the geometric disc graph belong to the class of stabilizing graphs. Intuitively, these conditions require the data at each node not to have strong dependence on data and positions of nodes far away. Finally, it is shown that spatial mixing of the Gibbs measure on a random graph holds when a suitably defined degree-dependent (but otherwise independent) node percolation does not have a giant component.

For Gaussian graphical models with cycles, loopy belief propagation often performs reasonably well, but its convergence is not guaranteed and the computation of variances is generally incorrect. In this paper, we identify a set of special vertices called a feedback vertex set whose removal results in a cycle-free graph. We propose a feedback message passing algorithm in which non-feedback nodes send out one set of messages while the feedback nodes use a different message update scheme. Exact inference results can be obtained in O(k2n), where k is the number of feedback nodes and n is the total number of nodes. For graphs with large feedback vertex sets, we describe a tractable approximate feedback message passing algorithm. Experimental results show that this procedure converges more often, faster, and provides better results than loopy belief propagation.

The problem of composite binary hypothesis testing of Markov forest (or tree) distributions is considered. The worst-case type-II error exponent is derived under the Neyman-Pearson formulation. Under simple null hypothesis, the error exponent is derived in closed-form and is characterized in terms of the so-called bottleneck edge of the forest distribution. The least favorable distribution for detection is shown to be Markov on the second-best max-weight spanning tree with mutual information edge weights. A necessary and sufficient condition to have positive error exponent is derived.

We consider recovering the salient feature subset for distinguishing between two probability models from i.i.d. samples. Identifying the salient set improves discrimination performance and reduces complexity. The focus in this work is on the high-dimensional regime where the number of variables d, the number of salient variables k and the number of samples n all grow. The definition of saliency is motivated by error exponents in a binary hypothesis test and is stated in terms of relative entropies. It is shown that if n grows faster than max{ck log((d-k)/k), exp(c'k)} for constants c, c', then the error probability in selecting the salient set can be made arbitrarily small. Thus, n can be much smaller than d. The exponential rate of decay and converse theorems are also provided. An efficient and consistent algorithm is proposed when the distributions are graphical models which are Markov on trees.

We consider Gaussian multiresolution (MR) models in which coarser, hidden variables serve to capture statistical dependencies among the finest scale variables. Tree-structured MR models have limited modeling capabilities, as variables at one scale are forced to be uncorrelated with each other conditioned on other scales. We propose a new class of Gaussian MR models that capture the residual correlations within each scale using sparse covariance structure. Our goal is to learn a tree-structured graphical model connecting variables across different scales, while at the same time learning sparse structure for the conditional covariance within each scale conditioned on other scales. This model leads to an efficient, new inference algorithm that is similar to multipole methods in computational physics.

The problem of binary hypothesis testing is considered when the measurements are drawn from a Markov random field (MRF) under each hypothesis. Spatial dependence of the measurements is incorporated by explicitly modeling the influence of sensor node locations on the clique potential functions of each MRF hypothesis. The nodes are placed i.i.d. in expanding areas with increasing sample size. Asymptotic performance of hypothesis testing is analyzed through the Neyman-Pearson type-II error exponent. The error exponent is expressed as the limit of a functional over dependency edges of the MRF hypotheses for acyclic graphs. Using the law of large numbers for graph functionals, the error exponent is derived.

The problem of maximum-likelihood learning of the structure of an unknown discrete distribution from samples is considered when the distribution is Markov on a tree. Large-deviation analysis of the error in estimation of the set of edges of the tree is performed. Necessary and sufficient conditions are provided to ensure that this error probability decays exponentially. These conditions are based on the mutual information between each pair of variables being distinct from that of other pairs. The rate of error decay, or error exponent, is derived using the large-deviation principle. The error exponent is approximated using Euclidean information theory and is given by a ratio, to be interpreted as the signal-to-noise ratio (SNR) for learning. Numerical experiments show the SNR approximation is accurate.

We propose a Bayesian nonparametric approach to relating multiple time series via a set of latent, dynamical behaviors. Using a beta process prior, we allow data-driven selection of the size of this set, as well as the pattern with which behaviors are shared among time series. Via the Indian buffet process representation of the beta process predictive distributions, we develop an exact Markov chain Monte Carlo inference method. In particular, our approach uses the sum-product algorithm to efficiently compute Metropolis-Hastings acceptance probabilities, and explores new dynamical behaviors via birth/death proposals. We validate our sampling algorithm using several synthetic datasets, and also demonstrate promising unsupervised segmentation of visual motion capture data.

In this paper, we investigate the use of message-passing algorithms for the problem of finding the max-weight independent set (MWIS) in a graph. First, we study the performance of the classical loopy max-product belief propagation. We show that each fixed-point estimate of max product can be mapped in a natural way to an extreme point of the linear programming (LP) polytope associated with the MWIS problem. However, this extreme point may not be the one that maximizes the value of node weights; the particular extreme point at final convergence depends on the initialization of max product. We then show that if max product is started from the natural initialization of uninformative messages, it always solves the correct LP, if it converges. This result is obtained via a direct analysis of the iterative algorithm, and cannot be obtained by looking only at fixed points. The tightness of the LP relaxation is thus necessary for max-product optimality, but it is not sufficient. Motivated by this observation, we show that a simple modification of max product becomes gradient descent on (a smoothed version of) the dual of the LP, and converges to the dual optimum. We also develop a message-passing algorithm that recovers the primal MWIS solution from the output of the descent algorithm. We show that the MWIS estimate obtained using these two algorithms in conjunction is correct when the graph is bipartite and the MWIS is unique. Finally, we show that any problem of maximuma posteriori(MAP) estimation for probability distributions over finite domains can be reduced to an MWIS problem. We believe this reduction will yield new insights and algorithms for MAP estimation.

We consider the problem of learning multiscale graphical models. Given a collection of variables along with covariance specifications for these variables, we introduce hidden variables and learn a sparse graphical model approximation on the entire set of variables (original and hidden). Our method for learning such models is based on maximizing entropy over an exponential family of graphical models, subject to divergence constraints on small subsets of variables. We demonstrate the advantages of our approach compared to methods that do not use hidden variables (which do not capture long-range behavior) and methods that use tree-structure approximations (which result in blocky artifacts).

Compressed sensing allows perfect recovery of sparse signals (or signals sparse in some basis) using only a small number of measurements. The results in the literature have focused on the asymptotics of how many samples are required and the probability of making an error for & fixed batch of samples. We investigate an alternative scenario where observations are available in sequence and can be stopped as soon as there is reasonable certainty of correct reconstruction. This approach does not require knowing how sparse is the signal, and allows reconstruction using the smallest number of samples. Central to our sequential approach is the stopping rule. For the random Gaussian ensemble we show that a simple stopping rule gives the absolute minimum number of observations required for exact recovery, with probability one. However, for other ensembles like Bernoulli or Fourier, this is no longer true, and the rule is modified to trade off delay in stopping and probability of error. We also consider near-sparse signals and describe how to estimate the reconstruction error from the sequence of solutions. This enables stopping once the error falls below a desired tolerance. Our sequential approach to compressed sensing involves a sequence of linear programs, and we outline how such a sequence can be solved efficiently.

We present a method for sequential learning of increasingly complex graphical models for discriminating between two hypotheses. We generate forests for each hypothesis, each with no more edges than a spanning tree, which optimize an information-theoretic criteria. The method relies on a straightforward extension of the efficient max-weight spanning tree (MWST) algorithm by incorporating multivalued edge-weights. Each iteration produces nested forests with increasing number of edges; each provably optimal as compared to alternative forests. Empirical results demonstrate superior probability of error as compared to generative approaches.

The hierarchical Dirichlet process hidden Markov model (HDP-HMM) is a flexible, nonparametric model which allows state spaces of unknown size to be learned from data. We demonstrate some limitations of the original HDP-HMM formulation (Teh et al., 2006), and propose a sticky extension which allows more robust learning of smoothly varying dynamics. Using DP mixtures, this formulation also allows learning of more complex, multimodal emission distributions. We further develop a sampling algorithm that employs a truncated approximation of the DP to jointly resample the full state sequence, greatly improving mixing rates. Via extensive experiments with synthetic data and the NIST speaker diarization database, we demonstrate the advantages of our sticky extension, and the utility of the HDP-HMM in real-world applications.

Many nonlinear dynamical phenomena can be effectively modeled by a system that switches among a set of conditionally linear dynamical modes. We consider two such models: the switching linear dynamical system (SLDS) and the switching vector autoregressive (VAR) process. In this paper, we present a nonparametric approach to the learning of an unknown number of persistent, smooth dynamical modes by utilizing a hierarchical Dirichlet process prior. We develop a sampling algorithm that combines a truncated approximation to the Dirichlet process with an efficient joint sampling of the mode and state sequences. The utility and flexibility of our model are demonstrated on synthetic data, sequences of dancing honey bees, and the IBOVESPA stock index.

This paper presents recursive cavity modeling--a principled, tractable approach to approximate, near-optimal inference for large Gauss-Markov random fields. The main idea is to subdivide the random field into smaller subfields, constructing cavity models which approximate these subfields. Each cavity model is a concise, yet faithful, model for the surface of one subfield sufficient for near-optimal inference in adjacent subfields. This basic idea leads to a tree-structured algorithm which recursively builds a hierarchy of cavity models during an "upward pass" and then builds a complementary set of blanket models during a reverse "downward pass." The marginal statistics of individual variables can then be approximated using their blanket models. Model thinning plays an important role, allowing us to develop thinned cavity and blanket models thereby providing tractable approximate inference. We develop a maximum-entropy approach that exploits certain tractable representations of Fisher information on thin chordal graphs. Given the resulting set of thinned cavity models, we also develop a fast preconditioner, which provides a simple iterative method to compute optimal estimates. Thus, our overall approach combines recursive inference, variational learning and iterative estimation. We demonstrate the accuracy and scalability of this approach in several challenging, large-scale remote sensing problems.

We propose a principled framework for recursively segmenting deformable objects across a sequence of frames. We demonstrate the usefulness of this method on left ventricular segmentation across a cardiac cycle. The approach involves a technique for learning the system dynamics together with methods of particle-based smoothing as well as nonparametric belief propagation on a loopy graphical model capturing the temporal periodicity of the heart. The dynamic system state is a low-dimensional representation of the boundary, and the boundary estimation involves incorporating curve evolution into recursive state estimation. By formulating the problem as one of state estimation, the segmentation at each particular time is based not only on the data observed at that instant, but also on predictions based on past and future boundary estimates. Although this paper focuses on left ventricle segmentation, the method generalizes to temporally segmenting any deformable object.

We propose a manifold learning approach to fiber tract clustering using a novel similarity measure between fiber tracts constructed from dual-rooted graphs. In particular, to generate this similarity measure, the chamfer or Hausdorff distance is initially employed as a local distance metric to construct minimum spanning trees between pairwise fiber tracts. These minimum spanning trees are effective in capturing the intrinsic geometry of the fiber tracts. Hence, they are used to capture the neighborhood structures of the fiber tract data set. We next assume the high-dimensional input fiber tracts to lie on low-dimensional non-linear manifolds. We apply Locally Linear Embedding, a popular manifold learning technique, to define a low-dimensional embedding of the fiber tracts that preserves the neighborhood structures of the high-dimensional data structure as captured by the method of dual-rooted graphs. Clustering is then performed on this low-dimensional data structure using the k-means algorithm. We illustrate our resulting clustering technique on both synthetic data and on real fiber tract data obtained from diffusion tensor imaging.

Summary form only given. Despite significant recent progress in automatic target exploitation (ATE) and recognition (ATR), current ATE systems do not meet the requirements of modern battlefield environments. Next generation ATE systems must actively manage sensor resources, aggregate sensed information across multiple platforms and diverse signaling modalities, and adapt to increasingly agile adversaries and operating conditions. Thus, the fundamental research challenge is to develop an integrated systems theory that jointly treats information fusion, control, and adaptation using multiple, dynamic multimodal sensor platforms in resource constrained environments.

We consider the problem of state estimation for a dynamic system driven by unobserved, correlated inputs. We model these inputs via an uncertain set of temporally correlated dynamic models, where this uncertainty includes the number of modes, their associated statistics, and the rate of mode transitions. The dynamic system is formulated via two interacting graphs: a hidden Markov model (HMM) and a linear-Gaussian state space model. The HMM's state space indexes system modes, while its outputs are the unobserved inputs to the linear dynamical system. This Markovian structure accounts for temporal persistence of input regimes, but avoids rigid assumptions about their detailed dynamics. Via a hierarchical Dirichlet process (HDP) prior, the complexity of our infinite state space robustly adapts to new observations. We present a learning algorithm and computational results that demonstrate the utility of the HDP for tracking, and show that it efficiently learns typical dynamics from noisy data.

We consider the problem of computing variances in large-scale Gauss-Markov random field (GMRF) models. In our prior work we considered the short-range correlation case, and we proposed a simple low-rank method which computes approximate variances with linear complexity in the number of nodes. In addition to its low complexity, the method has good guarantees on the quality of the approximation. In this paper we extend our method and analysis using a wavelet-based multi-scale approach which is applicable to models with much longer correlation lengths.

Many estimation problems involve sensors which can be actively controlled to alter the information received and utilized in the underlying inference task. In this paper, we discuss performance guarantees for heuristic algorithms for adaptive sensor control in sequential estimation problems, where the inference criterion is mutual information. We also demonstrate the performance of our tighter online computable performance guarantees through computational simulations. The guarantees may be applied to other estimation criteria including the Cramer-Rao bound.

We consider the estimation problem in Gaussian graphical models with arbitrary structure. We analyze the Embedded Trees algorithm, which solves a sequence of problems on tractable subgraphs thereby leading to the solution of the estimation problem on an intractable graph. Our analysis is based on the recently developed walk-sum interpretation of Gaussian estimation. We show that non-stationary iterations of the Embedded Trees algorithm using any sequence of subgraphs converge in walk-summable models. Based on walk-sum calculations, we develop adaptive methods that optimize the choice of subgraphs used at each iteration with a view to achieving maximum reduction in error. These adaptive procedures provide a significant speedup in convergence over stationary iterative methods, and also appear to converge in a larger class of models.

Loopy belief propagation has been employed in a wide variety of applications with great empirical success, but it comes with few theoretical guarantees. In this paper we investigate the use of the max-product form of belief propagation for weighted matching problems on general graphs. We show that max-product converges to the correct answer if the linear programming (LP) relaxation of the weighted matching problem is tight and does not converge if the LP relaxation is loose. This provides an exact characterization of max-product performance and reveals connections to the widely used optimization technique of LP relaxation. In addition, we demonstrate that max-product is effective in solving practical weighted matching problems in a distributed fashion by applying it to the problem of self-organization in sensor networks.

We investigate the use of message-passing algorithms for the problem of finding the max-weight independent set (MWIS) in a graph. First, we study the perfor- mance of loopy max-product belief propagation. We show that, if it converges, the quality of the estimate is closely related to the tightness of an LP relaxation of the MWIS problem. We use this relationship to obtain sufficient conditions for correctness of the estimate. We then develop a modification of max-product – one that converges to an optimal solution of the dual of the MWIS problem. We also develop a simple iterative algorithm for estimating the max-weight independent set from this dual solution. We show that the MWIS estimate obtained using these two algorithms in conjunction is correct when the graph is bipartite and the MWIS is unique. Finally, we show that any problem of MAP estimation for probability distributions over finite domains can be reduced to an MWIS problem. We believe this reduction will yield new insights and algorithms for MAP estimation.

Variational methods are frequently used to approximate or bound the partition or likelihood function of a Markov random field. Methods based on mean field theory are guaranteed to provide lower bounds, whereas certain types of convex relaxations provide upper bounds. In general, loopy belief propagation (BP) provides (often accurate) approximations, but not bounds. We prove that for a class of attractive binary models, the value specified by any fixed point of loopy BP always provides a lower bound on the true likelihood. Empirically, this bound is much better than the naive mean field bound, and requires no further work than running BP. We establish these lower bounds using a loop series expansion due to Chertkov and Chernyak, which we show can be derived as a consequence of the tree reparameterization characterization of BP fixed points.

We develop an integrated, probabilistic model for the appearance and three-dimensional geometry of cluttered scenes. Object categories are modeled via distributions over the 3D location and appearance of visual features. Uncertainty in the number of object instances depicted in a particular image is then achieved via a transformed Dirichlet process. In contrast with image-based approaches to object recognition, we model scale variations as the perspective projection of objects in different 3D poses. To calibrate the underlying geometry, we incorporate binocular stereo images into the training process. A robust likelihood model accounts for outliers in matched stereo features, allowing effective learning of 3D object structure from partial 2D segmentations. Applied to a dataset of office scenes, our model detects objects at multiple scales via a coarse reconstruction of the corresponding 3D geometry.

We consider the problem of detecting and localizing a material release utilizing sparse sensor measurements. We formulate the problem as one of abrupt change detection. The problem is challenging because of the sparse sensor deployment and complex system dynamics. We restrict ourselves to propagation models consisting of diffusion plus transport according to a Gaussian puff model. We derive optimal inference algorithms, provided the model parametrization is known precisely, within a hybrid detection-localization hypothesis testing framework with linear growth in the hypothesis space. The primary assumptions are that the mean wind field is deterministically known and that the Gaussian puff model is valid. Under these assumptions, we characterize the change in performance of detection, time-to-detection and localization as a function of the number of sensors. We then examine some performance impacts when the underlying dynamical model deviates from the assumed model

We consider the problem of variance estimation in large-scale Gauss-Markov random field (GMRF) models. While approximate mean estimates can be obtained efficiently for sparse GMRFs of very large size, computing the variances is a challenging problem. We propose a simple rank-reduced method which exploits the graph structure and the correlation length in the model to compute approximate variances with linear complexity in the number of nodes. The method has a separation length parameter trading off complexity versus estimation accuracy. For models with bounded correlation length, we efficiently compute provably accurate variance estimates

We present a new framework based on walks in a graph for analysis and inference in Gaussian graphical models. The key idea is to decompose the correlation between each pair of variables as a sum over all walks between those variables in the graph. The weight of each walk is given by a product of edgewise partial correlation coefficients. This representation holds for a large class of Gaussian graphical models which we call walk-summable. We give a precise characterization of this class of models, and relate it to other classes including diagonally dominant, attractive, non-frustrated, and pairwise-normalizable. We provide a walk-sum interpretation of Gaussian belief propagation in trees and of the approximate method of loopy belief propagation in graphs with cycles. The walk-sum perspective leads to a better understanding of Gaussian belief propagation and to stronger results for its convergence in loopy graphs.

The execution of reaching movements involves the coordinated activity of multiple brain regions that relate variously to the desired target and a path of arm states to achieve that target. These arm states may represent positions, velocities, torques, or other quantities. Estimation has been previously applied to neural activity in reconstructing the target separately from the path. However, the target and path are not independent. Because arm movements are limited by finite muscle contractility, knowledge of the target constrains the path of states that leads to the target. In this letter, we derive and illustrate a state equation to capture this basic dependency between target and path. The solution is described for discrete-time linear systems and gaussian increments with known target arrival time. The resulting analysis enables the use of estimation to study how brain regions that relate variously to target and path together specify a trajectory. The corresponding reconstruction procedure may also be useful in brain-driven prosthetic devices to generate control signals for goal-directed movements.

Some applications in remote sensing require estimating a field containing a discontinuity whose exact location is a priori unknown. Such fields of interest include sea surface temperature in oceanography and soil moisture in hydrology. For the former, oceanic fronts form a temperature discontinuity, while in the latter sharp changes exist across the interface between soil types. To complicate the estimation process, remotely sensed measurements often exhibit regions of missing observations due to occlusions such as cloud cover. Similarly, water surface and ground-based sensors usually provide only an incomplete set of measurements. Traditional methods of interpolation and smoothing for estimating the fields from such potentially sparse measurements often blur across the discontinuities in the field.

We explore the application of a homotopy continuation-based method for sparse signal representation in overcomplete dictionaries. Our problem setup is based on the basis pursuit framework, which involves a convex optimization problem consisting of terms enforcing data fidelity and sparsity, balanced by a regularization parameter. Choosing a good regularization parameter in this framework is a challenging task. We describe a homotopy continuation-based algorithm to find and trace efficiently all solutions of basis pursuit as a function of the regularization parameter. In addition to providing an attractive alternative to existing optimization methods for solving the basis pursuit problem, this algorithm can also be used to provide an automatic choice for the regularization parameter, based on prior information about the desired number of non-zero components in the sparse representation. Our numerical examples demonstrate the effectiveness of this algorithm in accurately and efficiently generating entire solution paths for basis pursuit, as well as producing reasonable regularization parameter choices. Furthermore, exploring the resulting solution paths in various operating conditions reveals insights about the nature of basis pursuit solutions.

Understanding the dependency structure of a set of variables is a key component in various signal processing applications which involve data association. The simple task of detecting whether any dependency exists is particularly difficult when models of the data are unknown or difficult to characterize because of high-dimensional measurements. We review the use of nonparametric tests for characterizing dependency and how to carry out these tests with high-dimensional observations. In addition we present a method to assess the significance of the tests.

Inter-sensor communication often comprises a significant portion of energy expenditures in a sensor network as compared to sensing and computation. We discuss an integrated approach to dynamically routing measurements and models in a sensor network. Specifically, we examine the problem of tracking objects within a region wherein the responsibility for combining measurements and updating a posterior state distribution is assigned to a single sensor at any given time step. The so called leader node may change over time. Sensor nodes communicate for two reasons: firstly, measurements of target state are transmitted from sensors to the current leader node for incorporation into the state estimate model; secondly, the state model is transmitted between sensors when the leader node changes. The trade-off between these two types of communication is of primary importance to dynamic selection of the leader node. We propose an algorithm based on a dynamic programming roll-out formulation of the minimum cost problem. We obtain a cost function which can be efficiently minimized by simplifying the problem to that of an open loop feedback controller which is an upper bound to the optimal cost. We present empirical results which compare methods previously proposed in the literature to the algorithm presented here.

We describe a hierarchical probabilistic model for the detection and recognition of objects in cluttered, natural scenes. The model is based on a set of parts which describe the expected appearance and position, in an object centered coordinate frame, of features detected by a low-level interest operator. Each object category then has its own distribution over these parts, which are shared between objects. We learn the parameters of this model via a Gibbs sampler which uses the graphical model's structure to analytically average over many parameters. Applied to a database of images of isolated objects, the sharing of parts among objects improves detection accuracy when few training examples are available. We also extend this hierarchical framework to scenes containing multiple objects

This paper presents a new framework based on walks in a graph for analysis and inference in Gaussian graphical models. The key idea is to decompose correlations between variables as a sum over all walks between those variables in the graph. The weight of each walk is given by a product of edgewise partial correlations. We provide a walk-sum interpretation of Gaussian belief propagation in trees and of the approximate method of loopy belief propagation in graphs with cycles. This perspective leads to a better understanding of Gaussian belief propagation and of its convergence in loopy graphs.

Given a directed graphical model with binary-valued hidden nodes and real-valued noisy observations, consider deciding upon the maximum a-posteriori (MAP) or the maximum posterior-marginal (MPM) assignment under the restriction that each node broadcasts only to its children exactly one single-bit message. We present a variational formulation, viewing the processing rules local to all nodes as degrees-of-freedom, that minimizes the loss in expected (MAP or MPM) performance subject to such online communication constraints. The approach leads to a novel message-passing algorithm to be executed offline, or before observations are realized, which mitigates the performance loss by iteratively coupling all rules in a manner implicitly driven by global statistics. We also provide (i) illustrative examples, (ii) assumptions that guarantee convergence and efficiency and (iii) connections to active research areas.

Motivated by the problem of learning to detect and recognize objects with minimal supervision, we develop a hierarchical probabilistic model for the spatial structure of visual scenes. In contrast with most existing models, our approach explicitly captures uncertainty in the number of object instances depicted in a given image. Our scene model is based on the transformed Dirichlet process (TDP), a novel extension of the hierarchical DP in which a set of stochastically transformed mixture components are shared between multiple groups of data. For visual scenes, mixture components describe the spatial structure of visual features in an objectcentered coordinate frame, while transformations model the object positions in a particular image. Learning and inference in the TDP, which has many potential applications beyond computer vision, is based on an empirically effective Gibbs sampler. Applied to a dataset of partially labeled street scenes, we show that the TDP's inclusion of spatial structure improves detection performance, flexibly exploiting partially labeled training images.

Belief propagation (BP) is an increasingly popular method of performing approximate inference on arbitrary graphical models. At times, even further approximations are required, whether due to quantization of the messages or model parameters, from other simplified message or model representations, or from stochastic approximation methods. The introduction of such errors into the BP message computations has the potential to affect the solution obtained adversely. We analyze the effect resulting from message approximation under two particular measures of error, and show bounds on the accumulation of errors in the system. This analysis leads to convergence conditions for traditional BP message passing, and both strict bounds and estimates of the resulting error in systems of approximate BP message passing.

Automatic self-localization is a critical need for the effective use of ad hoc sensor networks in military or civilian applications. In general, self-localization involves the combination of absolute location information (e.g., from a global positioning system) with relative calibration information (e.g., distance measurements between sensors) over regions of the network. Furthermore, it is generally desirable to distribute the computational burden across the network and minimize the amount of intersensor communication. We demonstrate that the information used for sensor localization is fundamentally local with regard to the network topology and use this observation to reformulate the problem within a graphical model framework. We then present and demonstrate the utility of nonparametric belief propagation (NBP), a recent generalization of particle filtering, for both estimating sensor locations and representing location uncertainties. NBP has the advantage that it is easily implemented in a distributed fashion, admits a wide variety of statistical models, and can represent multimodal uncertainty. Using simulations of small to moderately sized sensor networks, we show that NBP may be made robust to outlier measurement errors by a simple model augmentation, and that judicious message construction can result in better estimates. Furthermore, we provide an analysis of NBP's communications requirements, showing that typically only a few messages per sensor are required, and that even low bit-rate approximations of these messages can be used with little or no performance impact.

In this paper, we present a new information-theoretic approach to image segmentation. We cast the segmentation problem as the maximization of the mutual information between the region labels and the image pixel intensities, subject to a constraint on the total length of the region boundaries. We assume that the probability densities associated with the image pixel intensities within each region are completely unknown a priori, and we formulate the problem based on nonparametric density estimates. Due to the nonparametric structure, our method does not require the image regions to have a particular type of probability distribution and does not require the extraction and use of a particular statistic. We solve the information-theoretic optimization problem by deriving the associated gradient flows and applying curve evolution techniques. We use level-set methods to implement the resulting evolution. The experimental results based on both synthetic and real images demonstrate that the proposed technique can solve a variety of challenging image segmentation problems. Futhermore, our method, which does not require any training, performs as good as methods based on training.

We introduce a new class of upper bounds on the log partition function of a Markov random field (MRF). This quantity plays an important role in various contexts, including approximating marginal distributions, parameter estimation, combinatorial enumeration, statistical decision theory, and large-deviations bounds. Our derivation is based on concepts from convex duality and information geometry: in particular, it exploits mixtures of distributions in the exponential domain, and the Legendre mapping between exponential and mean parameters. In the special case of convex combinations of tree-structured distributions, we obtain a family of variational problems, similar to the Bethe variational problem, but distinguished by the following desirable properties: i) they are convex, and have a unique global optimum; and ii) the optimum gives an upper bound on the log partition function. This optimum is defined by stationary conditions very similar to those defining fixed points of the sum-product algorithm, or more generally, any local optimum of the Bethe variational problem. As with sum-product fixed points, the elements of the optimizing argument can be used as approximations to the marginals of the original model. The analysis extends naturally to convex combinations of hypertree-structured distributions, thereby establishing links to Kikuchi approximations and variants.

We develop and analyze methods for computing provably optimal maximum a posteriori probability (MAP) configurations for a subclass of Markov random fields defined on graphs with cycles. By decomposing the original distribution into a convex combination of tree-structured distributions, we obtain an upper bound on the optimal value of the original problem (i.e., the log probability of the MAP assignment) in terms of the combined optimal values of the tree problems. We prove that this upper bound is tight if and only if all the tree distributions share an optimal configuration in common. An important implication is that any such shared configuration must also be a MAP configuration for the original distribution. Next we develop two approaches to attempting to obtain tight upper bounds: a) a tree-relaxed linear program (LP), which is derived from the Lagrangian dual of the upper bounds; and b) a tree-reweighted max-product message-passing algorithm that is related to but distinct from the max-product algorithm. In this way, we establish a connection between a certain LP relaxation of the mode-finding problem and a reweighted form of the max-product (min-sum) message-passing algorithm.

Automatic self-calibration of ad-hoc sensor networks is a critical need for their use in military or civilian applications. In general, self-calibration involves the combination of absolute location information (e.g. GPS) with relative calibration information (e.g. estimated distance between sensors) over regions of the network. We formulate the self-calibration problem as a graphical model, enabling the application of nonparametric belief propagation (NBP), a recent generalization of particle filtering, for both estimating sensor locations and representing location uncertainties. NBP has the advantage that it is easily implemented in a distributed fashion, can represent multi-modal uncertainty, and admits a wide variety of statistical models. This last point is particularly appealing in that it can be used to provide robustness against occasional high-variance (outlier) noise. We illustrate the performance of NBP using Monte Carlo analysis on an example network.

We consider the problem of enforcing a sparsity prior in underdetermined linear problems, which is also known as sparse signal representation in overcomplete bases. The problem is combinatorial in nature, and a direct approach is computationally intractable, even for moderate data sizes. A number of approximations have been considered in the literature, including stepwise regression, matching pursuit and its variants, and, recently, basis pursuit (/spl lscr//sub 1/) and also /spl lscr//sub p/-norm relaxations with p<1. Although the exact notion of sparsity (expressed by an /spl lscr//sub 0/-norm) is replaced by /spl lscr//sub 1/ and /spl lscr//sub p/ norms in the latter two, it can be shown that under some conditions these relaxations solve the original problem exactly. The seminal paper of D.L. Donoho and X. Huo (see Stanford Univ. Tech. report: http://www-sccm.stanford.edu/pub/sccm/sccm02-17.pdf) establishes this fact for /spl lscr//sub 1/ (basis pursuit) for a special case where the linear operator is composed of an orthogonal pair. We extend their results to a general underdetermined linear operator. Furthermore, we derive conditions for the equivalence of /spl lscr//sub 0/ and /spl lscr//sub p/ problems, and extend the results to the problem of enforcing sparsity with respect to a transformation (which includes total variation priors as a special case). Finally, we describe an interesting result relating the sign patterns of solutions to the question of /spl lscr//sub 1/-/spl lscr//sub 0/ equivalence.

Automatic self-calibration of ad-hoc sensor networks is a critical need for their use in military or civilian applications. In general, self-calibration involves the combination of absolute location information (e.g. GPS) with relative calibration information (e.g. time delay or received signal strength between sensors) over regions of the network. Furthermore, it is generally desirable to distribute the computational burden across the network and minimize the amount of inter-sensor communication. We demonstrate that the information used for sensor calibration is fundamentally local with regard to the network topology and use this observation to reformulate the problem within a graphical model framework. We then demonstrate the utility of nonparametric belief propagation (NBP), a recent generalization of particle filtering, for both estimating sensor locations and representing location uncertainties. NBP has the advantage that it is easily implemented in a distributed fashion, admits a wide variety of statistical models, and can represent multi-modal uncertainty. We illustrate the performance of NBP on several example networks while comparing to a previously published nonlinear least squares method.

Belief propagation (BP) is an increasingly popular method of perform- ing approximate inference on arbitrary graphical models. At times, even further approximations are required, whether from quantization or other simplified message representations or from stochastic approxima- tion methods. Introducing such errors into the BP message computations has the potential to adversely affect the solution obtained. We analyze this effect with respect to a particular measure of message error, and show bounds on the accumulation of errors in the system. This leads both to convergence conditions and error bounds in traditional and approximate BP message passing.

We describe a threedimensional geometric hand model suitable for vi- sual tracking applications. The kinematic constraints implied by the model's joints have a probabilistic structure which is well described by a graphical model. Inference in this model is complicated by the hand's many degrees of freedom, as well as multimodal likelihoods caused by ambiguous image measurements. We use nonparametric belief propaga- tion (NBP) to develop a tracking algorithm which exploits the graph's structure to control complexity, while avoiding costly discretization. While kinematic constraints naturally have a local structure, self occlusions created by the imaging process lead to complex interpenden- cies in color and edgebased likelihood functions. However, we show that local structure may be recovered by introducing binary hidden vari- ables describing the occlusion state of each pixel. We augment the NBP algorithm to infer these occlusion variables in a distributed fashion, and then analytically marginalize over them to produce hand position esti- mates which properly account for occlusion events. We provide simula- tions showing that NBP may be used to refine inaccurate model initializa- tions, as well as track hand motion through extended image sequences. 1 Introduction Accurate visual detection and tracking of threedimensional articulated objects is a chal- lenging problem with applications in humancomputer interfaces, motion capture, and scene understanding [1]. In this paper, we develop a probabilistic method for tracking a geometric hand model from monocular image sequences. Because articulated hand mod- els have many (roughly 26) degrees of freedom, exact representation of the posterior dis- tribution over model configurations is intractable. Trackers based on extended and un- scented Kalman filters [2, 3] have difficulties with the multimodal uncertainties produced by ambiguous image evidence. This has motived many researchers to consider nonparamet- ric representations, including particle filters [4, 5] and deterministic multiscale discretiza- tions [6]. However, the hand's high dimensionality can cause these trackers to suffer catas- trophic failures, requiring the use of models which limit the hand's motion [4] or sophisti- cated prior models of hand configurations and dynamics [5, 6]. An alternative way to address the high dimensionality of articulated tracking problems is to describe the posterior distribution's statistical structure using a graphical model. Graph- Figure 1: Projected edges (left block) and silhouettes (right block) for a configuration of the 3D structural hand model matching the given image. To aid visualization, the model is also projected following rotations by 35 (center) and 70 (right) about the vertical axis. ical models have been used to track viewbased human body representations [7], con- tour models of restricted hand configurations [8], viewbased 2.5D "cardboard" models of hands and people [9], and a full 3D kinematic human body model [10]. Because the variables in these graphical models are continuous, and discretization is intractable for threedimensional models, most traditional graphical inference algorithms are inapplica- ble. Instead, these trackers are based on recently proposed extensions of particle filters to general graphs: mean field Monte Carlo in [9], and nonparametric belief propagation (NBP) [11, 12] in [10]. In this paper, we show that NBP may be used to track a threedimensional geometric model of the hand. To derive a graphical model for the tracking problem, we consider a redun- dant local representation in which each hand component is described by its own three dimensional position and orientation. We show that the model's kinematic constraints, including selfintersection constraints not captured by joint angle representations, take a simple form in this local representation. We also provide a local decomposition of the likelihood function which properly handles occlusion in a distributed fashion, a significant improvement over our earlier tracking results [13]. We conclude with simulations demon- strating our algorithm's robustness to occlusions. 2 Geometric Hand Modeling Structurally, the hand is composed of sixteen approximately rigid components: three pha- langes or links for each finger and thumb, as well as the palm [1]. As proposed by [2, 3], we model each rigid body by one or more truncated quadrics (ellipsoids, cones, and cylin- ders) of fixed size. These geometric primitives are well matched to the true geometry of the hand, allow tracking from arbitrary orientations (in contrast to 2.5D "cardboard" mod- els [5, 9]), and permit efficient computation of projected boundaries and silhouettes [3]. Figure 1 shows the edges and silhouettes corresponding to a sample hand model configu- ration. Note that only a coarse model of the hand's geometry is necessary for tracking. 2.1 Kinematic Representation and Constraints The kinematic constraints between different hand model components are well described by revolute joints [1]. Figure 2(a) shows a graph describing this kinematic structure, in which nodes correspond to rigid bodies and edges to joints. The two joints connecting the phalanges of each finger and thumb have a single rotational degree of freedom, while the joints connecting the base of each finger to the palm have two degrees of freedom (cor- responding to grasping and spreading motions). These twenty angles, combined with the palm's global position and orientation, provide 26 degrees of freedom. Forward kinematic transformations may be used to determine the finger positions corresponding to a given set of joint angles. While most modelbased hand trackers use this joint angle parameteriza- tion, we instead explore a redundant representation in which the ith rigid body is described by its position qi and orientation ri (a unit quaternion). Let xi = (qi, ri) denote this local description of each component, and x = {x1, . . . , x16} the overall hand configuration. Clearly, there are dependencies among the elements of x implied by the kinematic con- (a) (b) (c) (d) Figure 2: Graphs describing the hand model's constraints. (a) Kinematic constraints (EK ) de- rived from revolute joints. (b) Structural constraints (ES) preventing 3D component intersections. (c) Dynamics relating two consecutive time steps. (d) Occlusion consistency constraints (EO). straints. Let EK be the set of all pairs of rigid bodies which are connected by joints, or equivalently the edges in the kinematic graph of Fig. 2(a). For each joint (i, j) EK , define an indicator function K (x i,j i, xj ) which is equal to one if the pair (xi, xj ) are valid rigid body configurations associated with some setting of the angles of joint (i, j), and zero otherwise. Viewing the component configurations xi as random variables, the following prior explicitly enforces all constraints implied by the original joint angle representation: pK(x) K (x i,j i, xj ) (1) (i,j)EK Equation (1) shows that pK (x) is an undirected graphical model, whose Markov structure is described by the graph representing the hand's kinematic structure (Fig. 2(a)). 2.2 Structural and Temporal Constraints In reality, the hand's joint angles are coupled because different fingers can never occupy the same physical volume. This constraint is complex in a joint angle parameterization, but simple in our local representation: the position and orientation of every pair of rigid bodies must be such that their component quadric surfaces do not intersect. We approximate this ideal constraint in two ways. First, we only explicitly constrain those pairs of rigid bodies which are most likely to intersect, corresponding to the edges ES of the graph in Fig. 2(b). Furthermore, because the relative orientations of each finger's quadrics are implicitly constrained by the kinematic prior pK (x), we may detect most intersections based on the distance between object centroids. The structural prior is then given by 1 ||q p i - qj || > i,j S (x) S (x (x i,j i, xj ) S i,j i, xj ) = (2) 0 otherwise (i,j)ES where i,j is determined from the quadrics composing rigid bodies i and j. Empirically, we find that this constraint helps prevent different fingers from tracking the same image data. In order to track hand motion, we must model the hand's dynamics. Let xt denote the i position and orientation of the ith hand component at time t, and xt = {xt1, . . . , xt16}. For each component at time t, our dynamical model adds a Gaussian potential connecting it to the corresponding component at the previous time step (see Fig. 2(c)): 16 pT xt | xt-1 = N xt - xt-1; 0, i i i (3) i=1 Although this temporal model is factorized, the kinematic constraints at the following time step implicitly couple the corresponding random walks. These dynamics can be justified as the maximum entropy model given observations of the nodes' marginal variances i. 3 Observation Model Skin colored pixels have predictable statistics, which we model using a histogram distribu- tion pskin estimated from training patches [14]. Images without people were used to create a histogram model pbkgd of nonskin pixels. Let (x) denote the silhouette of projected hand configuration x. Then, assuming pixels are independent, an image y has likelihood p p skin(u) C (y | x) = pskin(u) pbkgd(v) (4) pbkgd(u) u(x) v(x) u(x) The final expression neglects the proportionality constant p v bkgd(v), which is inde- pendent of x, and thereby limits computation to the silhouette region [8]. 3.1 Distributed Occlusion Reasoning In configurations where there is no selfocclusion, pC (y | x) decomposes as a product of local likelihood terms involving the projections (xi) of individual hand components [13]. To allow a similar decomposition (and hence distributed inference) when there is occlu- sion, we augment the configuration xi of each node with a set of binary hidden variables zi = {zi } = 0 if pixel u in the projection of rigid body i is occluded (u) u. Letting zi(u) by any other body, and 1 otherwise, the color likelihood (eq. (4)) may be rewritten as 16 p z 16 i(u) p skin(u) C (y | x, z) = = p p C (y | xi, zi) (5) bkgd(u) i=1 u(xi) i=1 Assuming they are set consistently with the hand configuration x, the hidden occlusion variables z ensure that the likelihood of each pixel in (x) is counted exactly once. We may enforce consistency of the occlusion variables using the following function: 0 if x = 1 (x j occludes xi, u (xj ), and zi(u) j , zi ; x (6) (u) i) = 1 otherwise Note that because our rigid bodies are convex and nonintersecting, they can never take mutually occluding configurations. The constraint (xj, zi ; x (u) i) is zero precisely when pixel u in the projection of xi should be occluded by xj, but zi is in the unoccluded state. (u) The following potential encodes all of the occlusion relationships between nodes i and j: O (x (x ; x ; x i,j i, zi, xj , zj ) = j , zi(u) i) (xi, zj(u) j ) (7) u These occlusion constraints exist between all pairs of nodes. As with the structural prior, we enforce only those pairs EO (see Fig. 2(d)) most prone to xj occlusion: pO(x, z) O (x i,j i, zi, xj , zj ) (8) (i,j)EO z y i(u) xi Figure 3 shows a factor graph for the occlusion relationships between xi and its neighbors, as well as the observation potential pC (y | xi, zi). x u k The occlusion potential (xj, zi ; x (u) i) has a very Figure 3: Factor graph showing weak dependence on xi, depending only on p(y | xi, zi), and the occlusion con- whether xi is behind xj relative to the camera. straints placed on xi by xj , xk. Dashed lines denote weak dependencies. The 3.2 Modeling Edge Filter Responses plate is replicated once per pixel. Edges provide another important hand tracking cue. Using boundaries labeled in training images, we estimated a histogram pon of the response of a derivative of Gaussian filter steered to the edge's orientation [8, 10]. A similar histogram poff was estimated for filter outputs at randomly chosen locations. Let (x) denote the oriented edges in the projection of model configuration x. Then, again assuming pixel independence, image y has edge likelihood p 16 p z 16 i(u) p on(u) on(u) E (y | x, z) = = p p E (y | xi, zi) (9) off (u) poff(u) u(x) i=1 u(xi) i=1 where we have used the same occlusion variables z to allow a local decomposition. 4 Nonparametric Belief Propagation Over the previous sections, we have shown that a redundant, local representation of the geometric hand model's configuration xt allows p (xt | yt), the posterior distribution of the hand model at time t given image observations yt, to be written as 16 p xt | yt pK(xt)pS(xt)pO(xt, zt) pC(yt | xt, zt)p , zt) i i E (yt | xti i (10) zt i=1 The summation marginalizes over the hidden occlusion variables zt, which were needed to locally decompose the edge and color likelihoods. When video frames are observed, the overall posterior distribution is given by p (x | y) p xt | yt pT (xt | xt-1) (11) t=1 Excluding the potentials involving occlusion variables, which we discuss in detail in Sec. 4.2, eq. (11) is an example of a pairwise Markov random field: p (x | y) i,j (xi, xj) i (xi, y) (12) (i,j)E iV Hand tracking can thus be posed as inference in a graphical model, a problem we propose to solve using belief propagation (BP) [15]. At each BP iteration, some node i V calculates a message m (x ij j ) to be sent to a neighbor j (i) {j | (i, j) E}: mn (x mn-1 (x ij j ) j,i (xj , xi) i (xi, y) ki i) dxi (13) xi k(i)\j At any iteration, each node can produce an approximation ^ p(xi | y) to the marginal distri- bution p (xi | y) by combining the incoming messages with the local observation: ^ pn(xi | y) i (xi, yi) mn (x ji i) (14) j(i) For treestructured graphs, the beliefs ^ pn(xi | y) will converge to the true marginals p (xi | y). On graphs with cycles, BP is approximate but often highly accurate [15]. 4.1 Nonparametric Representations For the hand tracking problem, the rigid body configurations xi are sixdimensional con- tinuous variables, making accurate discretization intractable. Instead, we employ nonpara- metric, particlebased approximations to these messages using the nonparametric belief propagation (NBP) algorithm [11, 12]. In NBP, each message is represented using either a samplebased density estimate (a mixture of Gaussians) or an analytic function. Both types of messages are needed for hand tracking, as we discuss below. Each NBP message update involves two stages: sampling from the estimated marginal, followed by Monte Carlo ap- proximation of the outgoing message. For the general form of these updates, see [11]; the following sections focus on the details of the hand tracking implementation. The hand tracking application is complicated by the fact that the orientation component ri of xi = (qi, ri) is an element of the rotation group SO(3). Following [10], we represent orientations as unit quaternions, and use a linearized approximation when constructing den- sity estimates, projecting samples back to the unit sphere as necessary. This approximation is most appropriate for densities with tightly concentrated rotational components. 4.2 Marginal Computation BP's estimate of the belief ^ p(xi | y) is equal to the product of the incoming messages from neighboring nodes with the local observation potential (see eq. (14)). NBP approximates this product using importance sampling, as detailed in [13] for cases where there is no selfocclusion. First, M samples are drawn from the product of the incoming kinematic and temporal messages, which are Gaussian mixtures. We use a recently proposed multi- scale Gibbs sampler [16] to efficiently draw accurate (albeit approximate) samples, while avoiding the exponential cost associated with direct sampling (a product of d M Gaussian mixtures contains M d Gaussians). Following normalization of the rotational component, each sample is assigned a weight equal to the product of the color and edge likelihoods with any structural messages. Finally, the computationally efficient "rule of thumb" heuris- tic [17] is used to set the bandwidth of Gaussian kernels placed around each sample. To derive BP updates for the occlusion masks zi, we first cluster (xi, zi) for each hand component so that p (xt, zt | yt) has a pairwise form (as in eq. (12)). In principle, NBP could manage occlusion constraints by sampling candidate occlusion masks zi along with rigid body configurations xi. However, due to the exponentially large number of possible occlusion masks, we employ a more efficient analytic approximation. Consider the BP message sent from xj to (zi, xi), calculated by applying eq. (13) to the occlusion potential (x ; x u j , zi(u) i). We assume that ^ p(xj | y) is well separated from any candidate xi, a situation typically ensured by the kinematic and structural constraints. The occlusion constraint's weak dependence on xi (see Fig. 3) then separates the message computation into two cases. If xi lies in front of typical xj configurations, the BP message j,i(u)(zi ) is uninformative. If x (u) i is occluded, the message approximately equals j,i(u)(zi = 0) = 1 = 1) = 1 - Pr [u (x (u) j,i(u)(zi(u) j )] (15) where we have neglected correlations among pixel occlusion states, and where the prob- ability is computed with respect to ^ p(xj | y). By taking the product of these messages k,i(u)(zi ) from all potential occluders x (u) k and normalizing, we may determine an ap- proximation to the marginal occlusion probability i Pr[z = 0]. (u) i(u) Because the color likelihood pC (y | xi, zi) factorizes across pixels u, the BP approximation to pC (y | xi) may be written in terms of these marginal occlusion probabilites: p p skin(u) C (y | xi) i + (1 - ) (16) (u) i(u) pbkgd(u) u(xi) Intuitively, this equation downweights the color evidence at pixel u as the probability of that pixel's occlusion increases. The edge likelihood pE(y | xi) averages over zi similarly. The NBP estimate of ^ p(xi | y) is determined by sampling configurations of xi as before, and reweighting them using these occlusionsensitive likelihood functions. 4.3 Message Propagation To derive the propagation rule for nonocclusion edges, as suggested by [18] we rewrite the message update equation (13) in terms of the marginal distribution ^ p(xi | y): ^ pn-1(x mn (x i | y) dx ij j ) = j,i (xj , xi) i (17) x mn-1 (x i ji i) Our explicit use of the current marginal estimate ^ pn-1(xi | y) helps focus the Monte Carlo approximation on the most important regions of the state space. Note that messages sent 1 2 1 2 Figure 4: Refinement of a coarse initialization following one and two NBP iterations, both without (left) and with (right) occlusion reasoning. Each plot shows the projection of the five most significant modes of the estimated marginal distributions. Note the difference in middle finger estimates. along kinematic, structural, and temporal edges depend only on the belief ^ p(xi | y) follow- ing marginalization over occlusion variables zi. Details and pseudocode for the message propagation step are provided in [13]. For kine- matic constraints, we sample uniformly among permissable joint angles, and then use forward kinematics to propagate samples from ^ pn-1(xi | y) /mn-1 (x ji i) to hypothesized configurations of xj. Following [12], temporal messages are determined by adjusting the bandwidths of the current marginal estimate ^ p(xi | y) to match the temporal covariance i. Because structural potentials (eq. (2)) equal one for all state configurations outside some ball, the ideal structural messages are not finitely integrable. We therefore approximate the structural message m (x ij j ) as an analytic function equal to the weights of all kernels in ^ p(xi | y) outside a ball centered at qj, the position of xj.

In many applications of graphical models arising in computer vision, the hidden variables of interest are most naturally specified by continuous, non-Gaussian distributions. There exist inference algorithms for discrete approximations to these continuous distributions, but for the high-dimensional variables typically of interest, discrete inference becomes infeasible. Stochastic methods such as particle filters provide an appealing alternative. However, existing techniques fail to exploit the rich structure of the graphical models describing many vision problems. Drawing on ideas from regularized particle filters and belief propagation (BP), this paper develops a nonparametric belief propagation (NBP) algorithm applicable to general graphs. Each NBP iteration uses an efficient sampling procedure to update kernel-based approximations to the true, continuous likelihoods. The algorithm can accommodate an extremely broad class of potential functions, including nonparametric representations. Thus, NBP extends particle filtering methods to the more general vision problems that graphical models can describe. We apply the NBP algorithm to infer component interrelationships in a parts-based face model, allowing location and reconstruction of occluded features.

An information-theoretic method for multiphase image segmentation, in an active contour-based framework is proposed. Our approach is based on nonparametric density estimates, and is able to solve problems involving arbitrary probability densities for the region intensities. This is achieved by maximizing the mutual information between the region labels and the image pixel intensities, in order to segment up to 2/sup m/ regions using m curves. The method does not require any prior training regarding the regions of interest, but rather learns the probability densities during the evolution process. We present some illustrative experimental results, demonstrating the power of the proposed segmentation approach.

The problem of approximating the product of several Gaussian mixture distributions arises in a number of contexts, including the nonparametric belief propagation (NBP) inference algorithm and the training of prod- uct of experts models. This paper develops two multiscale algorithms for sampling from a product of Gaussian mixtures, and compares their performance to existing methods. The first is a multiscale variant of pre- viously proposed Monte Carlo techniques, with comparable theoretical guarantees but improved empirical convergence rates. The second makes use of approximate kernel density evaluation methods to construct a fast approximate sampler, which is guaranteed to sample points to within a tunable parameter (cid:15) of their true probability. We compare both multi- scale samplers on a set of computational examples motivated by NBP, demonstrating significant improvements over existing methods.

Efficient computation of extensions of banded, partially known covariance matrices is provided by the classical Levinson algorithm. One contribution of this paper is the introduction of a generalization of this algorithm that is applicable to a substantially broader class of extension problems. This generalized algorithm can compute unknown covariance elements in any order that satisfies certain graph-theoretic properties, which we describe. This flexibility, which is not provided by the classical Levinson algorithm, is then harnessed in a second contribution of this paper, the identification of a multiscale autoregressive (MAR) model for the maximum-entropy (ME) extension of a banded, partially known covariance matrix. The computational complexity of MAR model identification is an order of magnitude below that of explicitly computing a full covariance extension and is comparable to that required to build a standard autoregressive (AR) model using the classical Levinson algorithm.

We present a tree-based reparameterization (TRP) framework that provides a new conceptual view of a large class of algorithms for computing approximate marginals in graphs with cycles. This class includes the belief propagation (BP) or sum-product algorithm as well as variations and extensions of BP. Algorithms in this class can be formulated as a sequence of reparameterization updates, each of which entails refactorizing a portion of the distribution corresponding to an acyclic subgraph (i.e., a tree, or more generally, a hypertree). The ultimate goal is to obtain an alternative but equivalent factorization using functions that represent (exact or approximate) marginal distributions on cliques of the graph. Our framework highlights an important property of the sum-product algorithm and the larger class of reparameterization algorithms: the original distribution on the graph with cycles is not changed. The perspective of tree-based updates gives rise to a simple and intuitive characterization of the fixed points in terms of tree consistency. We develop interpretations of these results in terms of information geometry. The invariance of the distribution, in conjunction with the fixed-point characterization, enables us to derive an exact expression for the difference between the true marginals on an arbitrary graph with cycles, and the approximations provided by belief propagation. More broadly, our analysis applies to any algorithm that minimizes the Bethe free energy. We also develop bounds on the approximation error, which illuminate the conditions that govern their accuracy. Finally, we show how the reparameterization perspective extends naturally to generalizations of BP (e.g., Kikuchi (1951) approximations and variants) via the notion of hypertree reparameterization.

We propose a shape-based approach to curve evolution for the segmentation of medical images containing known object types. In particular, motivated by the work of Leventon, Grimson, and Faugeras, we derive a parametric model for an implicit representation of the segmenting curve by applying principal component analysis to a collection of signed distance representations of the training data. The parameters of this representation are then manipulated to minimize an objective function for segmentation. The resulting algorithm is able to handle multidimensional data, can deal with topological changes of the curve, is robust to noise and initial contour placements, and is computationally efficient. At the same time, it avoids the need for point correspondences during the training phase of the algorithm. We demonstrate this technique by applying it to two medical applications; two-dimensional segmentation of cardiac magnetic resonance imaging (MRI) and three-dimensional segmentation of prostate MRI.

We propose a novel non-parametric technique for source localization with passive sensor arrays. Our approach involves formulation of the problem in a variational framework where regularizing sparsity constraints are incorporated to achieve super-resolution and noise suppression. Compared to various source localization schemes, our approach offers increased resolution, significantly reduced sidelobes, and improved robustness to limitations in data quality and quantity. We demonstrate the effectiveness of the method on simulated data.

We propose a method for edge-preserving regularized reconstruction in coherent imaging systems. In our framework, image formation from measured data is achieved through the minimization of a cost function, which includes nonquadratic regularizing constraints for suppressing noise artifacts, while preserving the object boundaries in the reconstruction. The cost function we use effectively deals with the complex-valued and random-phase nature of the scattered field, which is inherent in many coherent systems. We solve the challenging optimization problems posed in our framework by a novel extension of half-quadratic regularization methods. We present experimental results from three coherent imaging applications: digital holography, synthetic aperture radar, and medical ultrasound. The proposed technique produces images where coherent speckle artifacts are effectively suppressed, and boundaries between different regions in the scene are preserved.

We present a novel information theoretic approach to image segmentation. We cast the segmentation problem as the maximization of the mutual information between the region labels and the image pixel intensities, subject to a constraint on the total length of the region boundaries. We assume that the probability densities associated with the image pixel intensities within each region are completely unknown a priori, and we formulate the problem based on nonparametric density estimates. Due to the nonparametric structure, our method does not require the image regions to have a particular type of probability distribution, and does not require the extraction and use of a particular statistic. We solve the information-theoretic optimization problem by deriving the associated gradient flows and applying curve evolution techniques. We use fast level set methods to implement the resulting evolution The evolution equations are based on nonparametric statistics, and have an intuitive appeal. The experimental results based on both synthetic and real images demonstrate that the proposed technique can solve a variety of challenging image segmentation problems.

In this paper, we introduce a novel approach for simultaneous restoration and segmentation of blurred noisy images by approaching a variant of the Mumford-Shah functional from a curve evolution perspective. In particular, by viewing the active contour as the set of discontinuities in the image, we derive a gradient flow to minimize an extended Mumford-Shah functional where the known blurring function is incorporated as part of the data fidelity term. Each gradient step involves solving a discrete approximation of the corresponding partial differential equation to obtain a smooth and deblurred estimate of the observed image without blurring across the curve. The experimental results based on both synthetic and real images demonstrate that the proposed method segments and restores the blurred images effectively. We conclude that our work is an edge-preserving image restoration technique that couples segmentation, denoising, and deblurring within a single framework. In addition, this framework provides an intellectual connection between regularization theory (used to solve the deblurring inverse problem) and the theory of curve evolution.

We describe a method for computing provably exact maximum a poste- riori (MAP) estimates for a subclass of problems on graphs with cycles. The basic idea is to represent the original problem on the graph with cy- cles as a convex combination of tree-structured problems. A convexity argument then guarantees that the optimal value of the original problem (i.e., the log probability of the MAP assignment) is upper bounded by the combined optimal values of the tree problems. We prove that this upper bound is met with equality if and only if the tree problems share an opti- mal configuration in common. An important implication is that any such shared configuration must also be the MAP configuration for the original problem. Next we develop a tree-reweighted max-product algorithm for attempting to find convex combinations of tree-structured problems that share a common optimum. We give necessary and sufficient conditions for a fixed point to yield the exact MAP estimate. An attractive feature of our analysis is that it generalizes naturally to convex combinations of hypertree-structured distributions.

Reviews a significant component of the rich field of statistical multiresolution (MR) modeling and processing. These MR methods have found application and permeated the literature of a widely scattered set of disciplines, and one of our principal objectives is to present a single, coherent picture of this framework. A second goal is to describe how this topic fits into the even larger field of MR methods and concepts-in particular, making ties to topics such as wavelets and multigrid methods. A third goal is to provide several alternate viewpoints for this body of work, as the methods and concepts we describe intersect with a number of other fields. The principle focus of our presentation is the class of MR Markov processes defined on pyramidally organized trees. The attractiveness of these models stems from both the very efficient algorithms they admit and their expressive power and broad applicability. We show how a variety of methods and models relate to this framework including models for self-similar and 1/f processes. We also illustrate how these methods have been used in practice.

While computer vision systems often assume simple illumination models, real-world illumination is highly complex, consisting of reflected light from every direction as well as distributed and localized primary light sources. One can capture the illumination incident at a point in the real world from every direction photographically using a spherical illumination map. This paper illustrates, through analysis of photographically-acquired, high dynamic range illumination maps, that real-world illumination shares many of the statistical properties of natural images. In particular, the marginal and joint wavelet coefficient distributions, directional derivative distributions, and harmonic spectra of illumination maps resemble those documented in the natural image statistics literature. However, illumination maps differ from standard photographs in that illumination maps are statistically non-stationary and may contain localized light sources that dominate their power spectra. Our work provides a foundation for statistical models of real-world illumination that may facilitate robust estimation of shape, reflectance, and illumination from images.

We propose a model-based curve evolution technique for segmentation of images containing known object types. In particular, motivated by the work of Leventon et al. (2000), we derive a parametric model for an implicit representation of the segmenting curve by applying principal component analysis to a collection of signed distance representations of the training data, The parameters of this representation are then calculated to minimize an objective function for segmentation. We found the resulting algorithm to be computationally efficient, able to handle multidimensional data, robust to noise and initial contour placements, while at the same time, avoiding the need for point correspondences during the training phase of the algorithm. We demonstrate this technique by applying it to two medical applications.

We present extensions to our previous work in modelling dynamical processes. The approach uses an information theoretic criterion for searching over subspaces of the past observations, combined with a nonparametric density characterizing its relation to one-step-ahead prediction and uncertainty. We use this methodology to model handwriting stroke data, specifically signatures, as a dynamical system and show that it is possible to learn a model capturing their dynamics for use either in synthesizing realistic signatures and in discriminating between signatures and forgeries even though no forgeries have been used in constructing the model. This novel approach yields promising results even for small training sets.

We develop a tree-based reparameterization framework that pro(cid:173) vides a new conceptual view of a large class of iterative algorithms for computing approximate marginals in graphs with cycles. It includes belief propagation (BP), which can be reformulated as a very local form of reparameterization. More generally, we consider algorithms that perform exact computations over spanning trees of the full graph. On the practical side, we find that such tree reparameterization (TRP) algorithms have convergence properties superior to BP. The reparameterization perspective also provides a number of theoretical insights into approximate inference, in(cid:173) cluding a new characterization of fixed points; and an invariance intrinsic to TRP /BP. These two properties enable us to analyze and bound the error between the TRP /BP approximations and the actual marginals. While our results arise naturally from the TRP perspective, most of them apply in an algorithm-independent manner to any local minimum of the Bethe free energy. Our re(cid:173) sults also have natural extensions to more structured approxima(cid:173) tions [e.g. , 1, 2].

In this work, we first address the problem of simultaneous image segmentation and smoothing by approaching the Mumford-Shah paradigm from a curve evolution perspective. In particular, we let a set of deformable contours define the boundaries between regions in an image where we model the data via piecewise smooth functions and employ a gradient flow to evolve these contours. Each gradient step involves solving an optimal estimation problem for the data within each region, connecting curve evolution and the Mumford-Shah functional with the theory of boundary-value stochastic processes. The resulting active contour model offers a tractable implementation of the original Mumford-Shah model (i.e., without resorting to elliptic approximations which have traditionally been favored for greater ease in implementation) to simultaneously segment and smoothly reconstruct the data within a given image in a coupled manner. Various implementations of this algorithm are introduced to increase its speed of convergence. We also outline a hierarchical implementation of this algorithm to handle important image features such as triple points and other multiple junctions. Next, by generalizing the data fidelity term of the original Mumford-Shah functional to incorporate a spatially varying penalty, we extend our method to problems in which data quality varies across the image and to images in which sets of pixel measurements are missing. This more general model leads us to a novel PDE-based approach for simultaneous image magnification, segmentation, and smoothing, thereby extending the traditional applications of the Mumford-Shah functional which only considers simultaneous segmentation and smoothing.

In this work, we approach the classic Mumford-Shah problem from a curve evolution perspective. In particular we let a given family of curves define the boundaries between regions in an image within which the data are modeled by piecewise smooth functions plus noise as in the standard Mumford-Shah functional. The gradient descent equation of this functional is then used to evolve the curve. Each gradient descent step involves solving a corresponding optimal estimation problem which connects the Mumford-Shah functional and our curve evolution implementation with the theory of boundary-value stochastic processes. The resulting active contour model, therefore, inherits the attractive ability of the Mumford-Shah technique to generate, in a coupled Mumford-Shah a smooth reconstruction of the image and a segmentation as well. We demonstrate applications of our method to problems in which data quality is spatially varying and to problems in which sets of pixel measurements are missing. Finally, we demonstrate a hierarchical implementation of our model which leads to a fast and efficient algorithm capable of dealing with important image features such as triple points.

We present an estimation-theoretic approach to curve evolution for the Mumford-Shah problem. By viewing an active contour as the set of discontinuities in the Mumford-Shah problem, we may use the corresponding functional to determine gradient descent evolution equations to deform the active contour. In each gradient descent step, we solve a corresponding optimal estimation problem, connecting the Mumford-Shah functional and curve evolution with the theory of boundary-value stochastic processes. In employing the Mumford-Shah functional, our active contour model inherits its attractive ability to generate, in a coupled manner, both a smooth reconstruction and a segmentation of the image. Next, by generalizing the data fidelity term of the original Mumford-Shah functional to incorporate a spatially varying penalty, we extend our method to problems in which data quality varies across the image and to images in which sets of pixel measurements are missing. This more general model leads us to a novel PDE-based approach for simultaneous image magnification, segmentation, and smoothing, thereby extending the traditional applications of the Mumford-Shah functional which only considers simultaneous segmentation and smoothing.

Multiresolution representations play an important role in image processing and computer vision, as well as in modeling stochastic processes. We have developed a semi-parametric class of non-Gaussian multiscale statistical processes defined by random cascades on wavelet trees. This model class is rich enough to accurately capture the remarkably regular non-Gaussian features of natural images, but sufficiently structured to permit estimation of the underlying state variables. We showed that our models accurately fit both the marginal and joint histograms of wavelet coefficients from natural images. We developed a Newton-like method for exact MAP state estimation that exploits fast algorithms for tree estimation, and hence is very efficient. Applications of this algorithm to denoising of both 1D signals and natural images were presented. The GSM-tree model class is related to a number of previous approaches to image coding and denoising.

We present the embedded trees algorithm, an iterative technique for estimation of Gaussian processes defined on arbitrary graphs. By exactly solving a series of modified problems on embedded span(cid:173) ning trees, it computes the conditional means with an efficiency comparable to or better than other techniques. Unlike other meth(cid:173) ods, the embedded trees algorithm also computes exact error co(cid:173) variances. The error covariance computation is most efficient for graphs in which removing a small number of edges reveals an em(cid:173) bedded tree. In this context, we demonstrate that sparse loopy graphs can provide a significant increase in modeling power rela(cid:173) tive to trees, with only a minor increase in estimation complexity.

We introduce a family of first-order multidimensional ordinary differential equations (ODEs) with discontinuous right-hand sides and demonstrate their applicability in image processing. An equation belonging to this family is an inverse diffusion everywhere except at local extrema, where some stabilization is introduced. For this reason, we call these equations "stabilized inverse diffusion equations" (SIDEs). Existence and uniqueness of solutions, as well as stability, are proven for SIDEs. A SIDE in one spatial dimension may be interpreted as a limiting case of a semi-discretized Perona-Malik equation. In an experiment, SIDE's are shown to suppress noise while sharpening edges present in the input signal. Their application to image segmentation is also demonstrated.

This paper addresses the problem of both segmenting and reconstructing a noisy signal or image. The work is motivated by large problems arising in certain scientific applications, such as medical imaging. Two objectives for a segmentation and denoising algorithm are laid out: it should be computationally efficient and capable of generating statistics for the errors in the reconstruction and estimates of the boundary locations. The starting point for the development of a suitable algorithm is a variational approach to segmentation (Shah 1992). This paper then develops a precise statistical interpretation of a one dimensional (1-D) version of this variational approach to segmentation. The 1-D algorithm that arises as a result of this analysis is computationally efficient and capable of generating error statistics. A straightforward extension of this algorithm to two dimensions would incorporate recursive procedures for computing estimates of inhomogeneous Gaussian Markov random fields. Such procedures require an unacceptably large number of operations. To meet the objective of developing a computationally efficient algorithm, the use of previously developed multiscale statistical methods is investigated. This results in the development of an algorithm for segmenting and denoising which is not only computationally efficient but also capable of generating error statistics, as desired.

A nonlinear diffusion process known to be effective for image segmentation is analyzed in 1-D. It is shown that it optimally solves certain edge detection problems. A fast implementation of the algorithm is introduced.

Computing the linear least-squares estimate of a high-dimensional random quantity given noisy data requires solving a large system of linear equations. In many situations, one can solve this system efficiently using the conjugate gradient (CG) algorithm. Computing the estimation error variances is a more intricate task. It is difficult because the error variances are the diagonal elements of a complicated matrix. This paper presents a method for using the conjugate search directions generated by the CG algorithm to obtain a converging approximation to the estimation error variances. The algorithm for computing the error variances falls out naturally from a novel estimation-theoretic interpretation of the CG algorithm. The paper discusses this interpretation and convergence issues and presents numerical examples.

We describe a new region based approach to active contours for segmenting images composed of two or three types of regions characterizable by a given statistic. The essential idea is to derive curve evolutions which separate two or more valves of a pre-determined set of statistics computed over geometrically determined subsets of the image. Both global and local image information is used to evolve the active contour. Image derivatives, however, are avoided, thereby giving rise to a further degree of noise robustness compared to most edge based snake algorithms.

We recently developed a multiscale recursive estimation procedure for the estimation of large-scale dynamic systems. The procedure propagates multiscale models for the estimation errors more efficiently than the Kalman filter's propagation of the error covariances, with a resulting computational complexity of 𝒪(N) and 𝒪(N3/2), where N is the number of variables estimated, for 1-D and 2-D dynamic systems, respectively. To further reduce the computational cost, we introduce in this paper a new class of reduced-order spatially-interpolated multiscale models, and demonstrate their use in remote.

A novel region-based approach to snakes is introduced in this paper for the segmentation of images composed of two or three types of regions where each region may be distinguished by a given statistic. The basic idea behind this technique is to formulate curve evolutions which separate two or more values of a predetermined set of statistics computed over geometrically determined subsets of the image data. Our methodology provides a natural framework for incorporating both global and local image information in the active contour motion while avoiding the use of image derivatives. As such, this technique possesses a robustness to noise which is noncharacteristic of most edge-based snake algorithms.

Statistically self-similar (SSS) processes can be used to describe a variety of physical phenomena, yet modeling these phenomena has proved challenging. Most of the proposed models for SSS and approximately SSS processes have power spectra that behave as 1/f/sup /spl gamma//, such as fractional Brownian motion (fBm), fractionally differenced noise, and wavelet-based syntheses. The most flexible framework is perhaps that based on wavelets, which provides a powerful tool for the synthesis and estimation of 1/f processes, but assumes a particular distribution of the measurements. An alternative framework is the class of multiresolution processes proposed by Chou et al. (1994), which has already been shown to be useful for the identification of the parameters of fBm. These multiresolution processes are defined by an autoregression in scale that makes them naturally suited to the representation of SSS (and approximately SSS) phenomena, both stationary and nonstationary. Also, this multiresolution framework is accompanied by an efficient estimator, likelihood calculator, and conditional simulator that make no assumptions about the distribution of the measurements. We show how to use the multiscale framework to represent SSS (or approximately SSS) processes such as fBm and fractionally differenced Gaussian noise. The multiscale models are realized by using canonical correlations (CC) and by exploiting the self-similarity and possible stationarity or stationary increments of the desired process. A number of examples are provided to demonstrate the utility of the multiscale framework in simulating and estimating SSS processes.

The multiscale autoregressive (MAR) framework was introduced to support the development of optimal multiscale statistical signal processing. Its power resides in the fast and flexible algorithms to which it leads. While the MAR framework was originally motivated by wavelets, the link between these two worlds has been previously established only in the simple case of the Haar wavelet. The first contribution of this paper is to provide a unification of the MAR framework and all compactly supported wavelets as well as a new view of the multiscale stochastic realization problem. The second contribution of this paper is to develop wavelet-based approximate internal MAR models for stochastic processes. This will be done by incorporating a powerful synthesis algorithm for the detail coefficients which complements the usual wavelet reconstruction algorithm for the scaling coefficients. Taking advantage of the statistical machinery provided by the MAR framework, we will illustrate the application of our models to sample-path generation and estimation from noisy, irregular, and sparse measurements.

Few fast statistical signal processing algorithms exist for large problems involving non-stationary processes and irregular measurements. A previously introduced class of multiscale autoregressive models indexed by trees admits signal processing algorithms which can efficiently deal with problems of this type. In this paper we provide a novel and efficient algorithm for translating any second-order prior model to a multiscale autoregressive prior model so that these efficient signal processing algorithms may be applied.

We present an estimation-theoretic approach for reducing motion effects in SAR imaging. Our approach may be viewed as a multi-dimensional matched-filter, whose parameters are determined by the relative motion between the SAR antenna and the target. In contrast to similar multi-dimensional matched filter methods, we propose a fast and easily implementable solution.

In order to be efficient and robust, distributed computing applications must accommodate new data sources and alter their computational structure automatically. We demonstrate the utility of mobile agents for addressing such challenges in a distributed, real-time application. We describe our implementation of a distributed information fusion system based on Bayesian networks using the D'Agent mobile agent system. The Bayesian networks infer identity for clusters of vehicles using information from distributed sensors and databases.

Variational methods have been employed with considerable success in computer vision, particularly for surface reconstruction problems. Formulations of this type require the solution of computationally complex Euler–Lagrange partial differential equations (PDEs) to obtain the desired reconstructions. Further, the calculation of reconstruction error covariances for such approaches are usually neglected. In this paper we describe a computationally efficient multiscale approach to surface reconstruction which differs fundamentally from other multiresolution methods that are used to solve the Euler–Lagrange PDEs. Instead, we interpret the variational problem as a statistical estimation problem in order to define a nearby, but slightly different , multiscale estimation problem that admits efficient solutions for both surface reconstruction and the calculation of error statistics. In particular, the membrane and thin-plate variational models for surfaces are interpreted as 1/ f 2 prior statistical models for the surface and its gradients, respectively. Such 1/ f 2 behavior is then achieved using a recently introduced class of multiresolution models that admits algorithms with constant per-pixel computational complexity.

In this paper, we investigate the problems of anomaly detection and localization from noisy tomographic data. These are characteristic of a class of problems that cannot be optimally solved because they involve hypothesis testing over hypothesis spaces with extremely large cardinality. Our multiscale hypothesis testing approach addresses the key issues associated with this class of problems. A multiscale hypothesis test is a hierarchical sequence of composite hypothesis tests that discards large portions of the hypothesis space with minimal computational burden and zooms in on the likely true hypothesis. For the anomaly detection and localization problems, hypothesis zooming corresponds to spatial zooming - anomalies are successively localized to finer and finer spatial scales. The key challenges we address include how to hierarchically divide a large hypothesis space and how to process the data at each stage of the hierarchy to decide which parts of the hypothesis space deserve more attention. For the latter, we pose and solve a nonlinear optimization problem for a decision statistic that maximally disambiguates composite hypotheses. With no more computational complexity, our optimized statistic shows substantial improvement over conventional approaches. We provide examples that demonstrate this and quantify how much performance is sacrificed by the use of a suboptimal method as compared to that achievable if the optimal approach were computationally feasible.

Magnetic resonance imaging (MRI) has become a widely used research and clinical tool in the study of the human brain. The ability to accurately segment the MRI data set into homogeneous regions such as gray matter, white matter, and cerebro spinal fluid aids in morphological quantification of brain features. The large amount of data associated with typical MRI brain scans makes completely manual segmentation prohibitive on a large scale. We develop an estimation-theoretic interpretation of the segmentation problem which leads to a computationally efficient, statistically based recursive technique for its solution. Being statistically based, the method also provides associated measures of uncertainty of the resulting estimates, which are useful both for evaluation of the estimates as well as their combination with other sources of information.

To abate the enormous costs incurred in the transmission and storage of SAR data, we present a segmentation driven compression technique using hierarchical stochastic modeling within a multiscale framework. Our approach to SAR image compression is unique in that we exploit the multiscale stochastic structure inherent in SAR imagery. This structure is well captured by a set of scale auto-regressive models that accurately characterize the evolution in scale. We thus use the local evolution in scale of SAR imagery to generate a segmentation map which is then used in tandem with the corresponding models to provide a robust, hierarchical compression technique.

Magnetic resonance imaging (MRI) has become a widely used research and clinical tool in the study of the human brain. The ability to robustly and accurately quantify repeatable morphological measures from such data is aided by the ability to accurately segment the MRI data set into homogeneous regions such as gray matter, white matter, and cerebro spinal fluid. The large amount of data associated with typical MRI scans makes completely manual segmentation prohibitive on a large scale. In this paper an efficient approach to the segmentation of such MR imagery is presented. The approach uses an estimation-theoretic interpretation of the segmentation problem to develop a computationally efficient, statistically-based recursive technique for its solution. Being statistically based, the method also provides associated measures of uncertainty of the resulting estimates, which are extremely important both for evaluation of the estimates as well as their combination with other sources of information.

This paper covers the design of multiscale stochastic models that can be used to fuse measurements of a random field or random process provided at multiple resolutions. Such sensor fusion problems arise in a variety of contexts, including many problems in remote sensing and geophysics. An example, which is used in this paper as a vehicle to illustrate our methodology, is the estimation of variations in hydraulic conductivity as required for the characterization of groundwater flow. Such a problem is typical in that the phenomenon to be estimated cannot be measured at fine scales throughout the region of interest, but instead must be inferred from a combination of measurements of very different types, including point measurements of hydraulic conductivity at irregular collections of points and indirect measurements that provide only coarse and nonlocal information about the conductivity field. Fusion of such disparate and irregular measurement sets is a challenging problem, especially when one includes the objective of producing, in addition to estimates, statistics characterizing the errors in those estimates. In this paper, we show how modeling a random field at multiple resolutions allows for the natural fusion (or assimilation) of measurements that provide information of different types and at different resolutions. The key to our approach is to take advantage of the fast multiscale estimation algorithms that efficiently produce both estimates and error variances even for very large problems. The major innovation required in our case, however, is to extend the modeling of random fields within this framework to accommodate multiresolution measurements. In particular to take advantage of the fast algorithms that the models admit, we must be able to model each nonlocal measurement as the measurement of a single variable of the multiresolution model at some appropriate resolution and scale. We describe how this can be done and illustrate its effectiveness for an ill-posed inverse problem in groundwater hydrology.

We use a natural pixel-type representation of an object, originally developed for incomplete data tomography problems, to construct nearly orthonormal multiscale basis functions. The nearly orthonormal behavior of the multiscale basis functions results in a system matrix, relating the input (the object coefficients) and the output (the projection data), which is extremely sparse. In addition, the coarsest scale elements of this matrix capture any ill conditioning in the system matrix arising from the geometry of the imaging system. We exploit this feature to partition the system matrix by scales and obtain a reconstruction procedure that requires inversion of only a well-conditioned and sparse matrix. This enables us to formulate a tomographic reconstruction technique from incomplete data wherein the object is reconstructed at multiple scales or resolutions. In case of noisy projection data we extend our multiscale reconstruction technique to explicitly account for noise by calculating maximum a posteriori probability (MAP) multiscale reconstruction estimates based on a certain self-similar prior on the multiscale object coefficients. The framework for multiscale reconstruction presented can find application in regularization of imaging problems where the projection data are incomplete, irregular, and noisy, and in object feature recognition directly from projection data.

We present efficient multiscale approaches to the segmentation of natural clutter, specifically grass and forest, and to the enhancement of anomalies in synthetic aperture radar (SAR) imagery. The methods we propose exploit the coherent nature of SAR sensors. In particular, they take advantage of the characteristic statistical differences in imagery of different terrain types, as a function of scale, due to radar speckle. We employ a class of multiscale stochastic processes that provide a powerful framework for describing random processes and fields that evolve in scale. We build models representative of each category of terrain of interest (i.e., grass and forest) and employ them in directing decisions on pixel classification, segmentation, and anomalous behaviour. The scale-autoregressive nature of our models allows extremely efficient calculation of likelihoods for different terrain classifications over windows of SAR imagery. We subsequently use these likelihoods as the basis for both image pixel classification and grass-forest boundary estimation. In addition, anomaly enhancement is possible with minimal additional computation. Specifically, the residuals produced by our models in predicting SAR imagery from coarser scale images are theoretically uncorrelated. As a result, potentially anomalous pixels and regions are enhanced and pinpointed by noting regions whose residuals display a high level of correlation throughout scale. We evaluate the performance of our techniques through testing on 0.3-m resolution SAR data gathered with Lincoln Laboratory's millimeter-wave SAR.

Recently, a class of multiscale stochastic models has been introduced in which random processes and fields are described by scale-recursive dynamic trees. A major advantage of this framework is that it leads to an extremely efficient, statistically optimal algorithm for least-squares estimation. In certain applications, however, estimates based on the types of multiscale models previously proposed may not be adequate, as they have tended to exhibit a visually distracting blockiness. We eliminate this blockiness by discarding the standard assumption that distinct nodes on a given level of the multiscale process correspond to disjoint portions of the image domain; instead, we allow a correspondence to overlapping portions of the image domain. We use these so-called overlapping-tree models for both modeling and estimation. In particular, we develop an efficient multiscale algorithm for generating sample paths of a random field whose second-order statistics match a prespecified covariance structure, to any desired degree of fidelity. Furthermore, we demonstrate that under easily satisfied conditions, we can "lift" a random field estimation problem to one defined on an overlapped tree, resulting in an estimation algorithm that is computationally efficient, directly produces estimation error covariances, and eliminates blockiness in the reconstructed imagery without any sacrifice in the resolution of fine-scale detail.

The successful application of a previously developed tree-based multiscale estimation framework to large-scale static estimation problems has proven the statistical flexibility and the computational efficiency of the framework. However, the estimation of processes evolving in time remains a considerable challenge. We address this challenge by investigating multiscale models for the steady-state estimation of dynamic systems. In particular, we build such models for 1-D and 2-D heat diffusion processes, making use of a canonical correlations realization technique. The problem of the assimilation of satellite altimetric measurements of the dynamic ocean surface over time is also addressed.

The general problem of joint detection of linearly dependent users in an uncoded multiple access (MA) system is N-P hard. We look to exploit the existing structure in our problem so that low complexity algorithms may be devised to yield the optimal solution. Advantage is taken over the design of user signatures in a typical MA communication system. By imposing a hierarchical cross-correlation structure on the user signature waveforms, the receiver design problem is reduced so that it is no longer N-P complete. A tree joint detection algorithm which takes advantage of such a cross-correlation structure is presented. The tree detector gives the optimal estimate with an extremely low computational complexity that is typically low-order-polynomial in the number of users. This is an enormous savings in computations over the O(2/sup K/) computations needed if the signatures did not exhibit any structure.

This work presents a method for segmenting images based on gradients in the intensity function. Past approaches have centered on formulating the problem in the context of variational calculus as the minimization of a functional involving the image intensity and edge functions. Computational methods for finding the minima of such variational problems are prone to two shortfalls: they are often computationally intensive and almost always incapable of computing error statistics associated with the segmentation. Using a particular variational formulation as a starting point, this paper presents a derivation of an associated statistical formulation using multiscale models. The result is an algorithm which is fast and capable of computing error statistics.

A fast multiscale optimal interpolation algorithm has been adapted to the mapping of hydrographic and other oceanographic data. This multiscale algorithm produces solution and error estimates which are consistent with those obtained from exact least-squares methods, but at a small fraction of the computational cost. Problems whose solution would be completely impractical using exact least-squares, that is problems with tens or hundreds of thousands of measurements and estimation grid points, can easily be solved on a small workstation using the multiscale algorithm. Contrary to methods previously proposed for solving large least-squares problems, the multiscale approach provides error statistics while permitting long-range correlations, using all measurements, and permitting arbitrary measurement locations.

We present an efficient multiscale approach to the segmentation of natural clutter, specifically grass and forest, in synthetic aperture imagery (SAR) and to the enhancement of anomalous image regions therein. The methods we propose exploit the coherent nature of SAR sensors. In particular, they characterize the scale-to-scale statistical differences in imagery of various terrain categories due to radar speckle. To achieve this, we employ a recently introduced class of multiscale stochastic processes that provide a powerful framework for describing random processes and fields that evolve in scale. We build models representative of each relevant category of terrain and use them to direct subsequent decisions on pixel classification, segmentation, and anomaly presence.

Anomaly characterization from tomographic measurements is of interest in a wide range of fields. In this paper we address the problems of single anomaly detection and localization which we formulate as hypothesis testing problems. While the optimal hypothesis test is easy to formulate, it is computationally infeasible due to the overwhelming number of hypotheses which must be considered. We propose the multiscale hypothesis test (MSHT) as an efficient suboptimal alternative. We show how to find decision statistics to achieve maximal composite hypothesis distinguishability for the composite hypothesis tests which comprise the MSHT.

In this work, we develop statistically based algorithms to reconstruct binary polygonal objects from sparse and noisy tomographic-based observation data. Traditional approaches to the reconstruction of geometric objects from projection data often lead to highly nonlinear estimation problems. To avoid the difficulties associated with such nonlinear problems, we first examine the problem of reconstruction of an object based on knot location measurements, i.e., measurements of the locations of abrupt change in the projections. The ties between this problem and that of multitarget radar tracking enable us to develop a sequential hypothesis-testing algorithm requiring only the solution of a series of linear estimation problems. In particular, data association hypotheses are generated, under each of which the inversion is linear. The complexity of the association possibilities are kept in check through the use of constraints on the reconstruction imposed by the tomography problem. The solution of this first problem is then used as an initialization to a more complete reconstruction which, while utilizing all the projection data, is nonlinear. We demonstrate that the estimates provided by the first, efficient algorithm are of good quality on their own, and, when combined with a fully nonlinear inversion, produce excellent object estimates.

Includes bibliographical references (p. 34-36).

We describe a variational framework for the tomographic reconstruction of an image from the maximum likelihood (ML) estimates of its orthogonal moments. We show how these estimated moments and their (correlated) error statistics can be computed directly, and in a linear fashion from given noisy and possibly sparse projection data. Moreover, thanks to the consistency properties of the Radon transform, this two-step approach (moment estimation followed by image reconstruction) can be viewed as a statistically optimal procedure. Furthermore, by focusing on the important role played by the moments of projection data, we immediately see the close connection between tomographic reconstruction of nonnegative valued images and the problem of nonparametric estimation of probability densities given estimates of their moments. Taking advantage of this connection, our proposed variational algorithm is based on the minimization of a cost functional composed of a term measuring the divergence between a given prior estimate of the image and the current estimate of the image and a second quadratic term based on the error incurred in the estimation of the moments of the underlying image from the noisy projection data. We show that an iterative refinement of this algorithm leads to a practical algorithm for the solution of the highly complex equality constrained divergence minimization problem. We show that this iterative refinement results in superior reconstructions of images from very noisy data as compared with the classical filtered back-projection (FBP) algorithm.

The authors represent the standard ramp filter operator of the filtered-back-projection (FBP) reconstruction in different bases composed of Haar and Daubechies compactly supported wavelets. The resulting multiscale representation of the ramp-filter matrix operator is approximately diagonal. The accuracy of this diagonal approximation becomes better as wavelets with larger numbers of vanishing moments are used. This wavelet-based representation enables the authors to formulate a multiscale tomographic reconstruction technique in which the object is reconstructed at multiple scales or resolutions. A complete reconstruction is obtained by combining the reconstructions at different scales. The authors' multiscale reconstruction technique has the same computational complexity as the FBP reconstruction method. It differs from other multiscale reconstruction techniques in that (1) the object is defined through a one-dimensional multiscale transformation of the projection domain, and (2) the authors explicitly account for noise in the projection data by calculating maximum a posteriori probability (MAP) multiscale reconstruction estimates based on a chosen fractal prior on the multiscale object coefficients. The computational complexity of this maximum a posteriori probability (MAP) solution is also the same as that of the FBP reconstruction. This result is in contrast to commonly used methods of statistical regularization, which result in computationally intensive optimization algorithms.

A significant problem in oceanographic remote sensing is the dense gridding or smoothing of sparsely sampled altimetric data. The smoothing of altimetric measurements has an application much broader than just the regular production of elevation maps for oceanographers, however. In particular, the ability to estimate ocean circulation patterns from altimetric data can serve as an important measure for the verification of ocean acoustic tomographic results. The authors present a multiscale technique capable of extremely efficient interpolation of altimetric data: about 256000 estimates and estimation error variances are computed in one minute on a Sun Sparc-10. They also demonstrate how similar techniques may be used to directly estimate the surface gradient and biases in the geoid-model error.

The joint detection of all users in a multiple access (MA) communication system in which user transmissions are correlated has been shown in recent literature to enhance the system performance relative to that achieved without joint detection. Over the past several years the area of low complexity joint detectors has received much attention. This paper explains the problem of multiple access joint detection in geometrical terms. Geometric interpretation leads to the proposal of an alternating projection joint detection algorithm (APJD). Due to some similarities between our APJD and the multistage joint detector (MJD) of Varansi and Aazhang (1990), the MJD is also discussed. The APJD is guaranteed to converge and a proof is given. The geometric interpretation of the MA joint detection problem allows for the exploration of determining, a priori, the error probability of a joint detector and user waveform set in the absence of noise. Simulations offer empirical characterization of the error behavior of both detectors.

We propose a best basis algorithm for signal enhancement in white Gaussian noise. We base our search of best basis on a criterion of minimal reconstruction error of the underlying signal. We subsequently compare our simple error criterion to the Stein (1981) unbiased risk estimator, and provide a substantiating example to demonstrate its performance. A review is also given of noise removal by thresholding and of wavepacket orthonormal bases.

Heterogeneities of hydraulic conductivity across multiple spatial scales can significantly affect the flow of groundwater, yet the sparsity of subsurface observations does not allow accurate reconstruction at all such scales. Instead, a common approach is to estimate an effective parameter which predicts large scale variations in flow. We build a new framework upon maximum a posteriori inversions to directly address the problems of scale. This framework consists of (1) a fractal prior model for conductivity which is an autoregressive (AR) process evolving from coarse to fine scale, and (2) a measurement model in which each observation is well-approximated by a linear functional of the AR process at some scale. This framework efficiently incorporates multiple measurement sources, which consist of point measurements of hydraulic conductivity and a piezometric head. The multi-resolution framework produces a conductivity estimate with spatially varying resolution tailored to the measurement sampling geometry.

A class of multiscale stochastic models has been introduced in which Gaussian random processes are described by scale-recursive dynamics that are indexed by the nodes of a tree. One of the primary reasons the framework is useful is that it leads to an extremely fast, statistically optimal algorithm for least-squares estimation in the context of 2-D images. We refine this approach to estimation by eliminating the visually distracting blockiness that has been observed in the previous work. We eliminate the blockiness by discarding the standard assumption that distinct nodes at a given level of our tree correspond to disjoint portions of the image domain; as a consequence of this simple idea, a given image pixel may now correspond to several tree nodes. We develop tools for systematically building overlapping-tree multiscale representations of prespecified statistics, and we develop a corresponding estimation algorithm for this processes. In this way, we achieve nearly optimal estimation results, we generate corresponding error covariance information, and we eliminate blockiness without sacrificing the resolution of fine-scale detail.

In this work, a novel method of multiscale geometric feature extraction and object recognition is developed. In particular, the new representation should have the following characteristics. First, the coarse scale features should have a geometric interpretation so that the overall geometry of the object is discernible from just these features. Second, the presence of fine scale detail should not change the coarse scale representation. These two goals are not achieved by current techniques which are based on error as measured by the L/sup 2/ norm. Two methods to accomplish these goals are presented. In the first, morphological filtering and wavelet networks are used. In the second, the correlation criteria of the matching pursuit algorithm of Mallat and Zhang (1993) is modified to obtain a variable, high resolution matching pursuit.

A recently developed multiresolution estimation framework offers the possibility of highly efficient statistical analysis, interpolation, and smoothing of extremely large data sets in a multiscale fashion. This framework enjoys a number of advantages not shared by other statistically-based methods. In particular, the algorithms resulting from this framework have complexity that scales only linearly with problem size, yielding constant complexity load per grid point independent of problem size. Furthermore these algorithms directly provide interpolated estimates at multiple resolutions, accompanying error variance statistics of use in assessing resolutionlaccuracy tradeoffs and in detecting statistically significant anomalies, and maximum likelihood estimates of parameters such as spectral power law coefficients. Moreover, the efficiency of these algorithms is completely insensitive to irregularities in the sampling or spatial distribution of measurements and to heterogeneities in measurement errors or model parameters. For these reasons this approach has the potential of being an effective tool in a variety of remote sensing problems. In this paper, we demonstrate a realization of this potential by applying the multiresolution framework to a problem of considerable current interest-the interpolation and statistical analysis of ocean surface data from the TOPEXPOSEIDON altimeter

A class of multiscale stochastic models based on scale-recursive dynamics on trees has previously been introduced. Theoretical and experimental results have shown that these models provide an extremely rich framework for representing both processes which are intrinsically multiscale, e.g., 1/f processes, as well as 1D Markov processes and 2D Markov random fields. Moreover, efficient optimal estimation algorithms have been developed for these models by exploiting their scale-recursive structure. The authors exploit this structure in order to develop a computationally efficient and parallelizable algorithm for likelihood calculation. They illustrate one possible application to texture discrimination and demonstrate that likelihood-based methods using the algorithm achieve performance comparable to that of Gaussian Markov random field based techniques, which in general are prohibitively complex computationally.

The estimation of dynamically evolving ellipsoids from noisy lower-dimensional projections is examined. In particular, this work describes a model-based approach using geometric reconstruction and recursive estimation techniques to obtain a dynamic estimate of left-ventricular ejection fraction from a gated set of planar myocardial perfusion images. The proposed approach differs from current ejection fraction estimation techniques both in the imaging modality used and in the subsequent processing which yields a dynamic ejection fraction estimate. For this work, the left ventricle is modeled as a dynamically evolving three-dimensional (3-D) ellipsoid. The left-ventricular outline observed in the myocardial perfusion images is then modeled as a dynamic, two-dimensional (2-D) ellipsoid, obtained as the projection of the former 3-D ellipsoid. This data is processed in two ways: first, as a 3-D dynamic ellipsoid reconstruction problem; second, each view is considered as a 2-D dynamic ellipse estimation problem and then the 3-D ejection fraction is obtained by combining the effective 2-D ejection fractions of each view. The approximating ellipsoids are reconstructed using a Rauch-Tung-Striebel smoothing filter, which produces an ejection fraction estimate that is more robust to noise since it is based on the entire data set; in contrast, traditional ejection fraction estimates are based only on true frames of data. Further, numerical studies of the sensitivity of this approach to unknown dynamics and projection geometry are presented, providing a rational basis for specifying system parameters. This investigation includes estimation of ejection fraction from both simulated and real data.

Proposes a technique for discrimination of fractal fields with different fractal dimensions, directly from the noisy and sparse tomographic projection data. This application is motivated from the medical field, where a change in fractal dimension is used to differentiate between normal and abnormal conditions in many different contexts, including diagnosis of liver abnormalities. The conventional method for discrimination of fractal fields from tomographic data is based on the calculation of the slope of the power spectra of the corresponding projections. This method, derived from the Radon transform results, breaks down in case the projection data are sparse and/or noisy. In order to avoid any restrictions on the duality and quantity of the projection data, we formulate our discrimination problem in a discrete hypothesis testing framework, the solution to which is given by the maximum-log-likelihood discrimination rule. The problem of discriminating fractal fields through likelihood calculations is, however, complicated by the fact that inverses and determinants of large, full, and generally ill conditioned fractal-field data covariance matrices are required. We show that these complications in the likelihood calculations can be removed by a transformation to the multiscale framework. The multiscale data covariance matrices are sparse and in addition, are naturally partitioned into ill conditioned coarsest scale approximation blocks and relatively well conditioned multiscale detail blocks. We simplify our likelihood calculations by using the class of multiscale stochastic models defined on trees to realize accurate approximations of the detail block of the data covariance matrices.>

Large data assimilation problems present a number of challenges, many of which the authors' multiscale estimation framework is capable of addressing. They demonstrate the application of their estimation framework to an oceanographic data assimilation problem of considerable interest: the smoothing of Topex/Poseidon ocean altimetric data. they produce surface estimates and error statistics on a 512/spl times/512 grid in one minute on a Sparc-10. Using a generic surface reconstruction problem they demonstrate alternative uses of the multiscale framework: higher order state models, combining estimates from an ensemble of trees, and overlapping tree models.>

In this paper, we discuss two interesting instantiations of the moment problem in image processing. The first involves the estimation of moments of an image indirectly from projections, and the reconstruction of the image from these moments. The second relates the reconstruction of binary polygons from moments to well-known algorithm in array signal processing. Through these examples, we place the moment problem into a geometric perspective and illustrate how this perspective leads to a number of interesting practical applications in image processing and other fields.>

Discusses a class of multiresolution models of random fields based on a generalization of the midpoint deflection construction of the 1D Brownian motion. The authors then present least squares (LS) algorithms for the estimation of parameters which define these models and hence provide a framework for synthesizing and analyzing images with fractal-like properties such as those found in statistical representation of natural terrain and other geophysical phenomena. The authors also briefly discuss possible applications of this modeling framework to target detection in images.>

In the computation of dense optical flow fields, spatial coherence constraints are commonly used to regularize otherwise ill-posed problem formulations, providing spatial integration of data. We present a temporal, multiframe extension of the dense optical flow estimation formulation proposed by Horn and Schunck (1981) in which we use a temporal coherence constraint to yield the optimal fusing of data from multiple frames of measurements. Conceptually, standard Kalman filtering algorithms are applicable to the resulting multiframe optical flow estimation problem, providing a solution that is sequential and recursive in time. Experiments are presented to demonstrate that the resulting multiframe estimates are more robust to noise than those provided by the original, single-frame formulation. In addition, we demonstrate cases where the aperture problem of motion vision cannot be resolved satisfactorily without the temporal integration of data enabled by the proposed formulation. Practically, the large matrix dimensions involved in the problem prohibit exact implementation of the optimal Kalman filter. To overcome this limitation, we present a computationally efficient, yet near-optimal approximation of the exact filtering algorithm. This approximation has a precise interpretation as the sequential estimation of a reduced-order spatial model for the optical flow estimation error process at each time step and arises from an estimation-theoretic treatment of the filtering problem. Experiments also demonstrate the efficacy of this near-optimal filter.

A new approach to regularization methods for image processing is introduced and developed using as a vehicle the problem of computing dense optical flow fields in an image sequence. The solution of the new problem formulation is computed with an efficient multiscale algorithm. Experiments on several image sequences demonstrate the substantial computational savings that can be achieved due to the fact that the algorithm is noniterative and in fact has a per pixel computational complexity that is independent of image size. The new approach also has a number of other important advantages. Specifically, multiresolution flow field estimates are available, allowing great flexibility in dealing with the tradeoff between resolution and accuracy. Multiscale error covariance information is also available, which is of considerable use in assessing the accuracy of the estimates. In particular, these error statistics can be used as the basis for a rational procedure for determining the spatially-varying optimal reconstruction resolution. Furthermore, if there are compelling reasons to insist upon a standard smoothness constraint, the new algorithm provides an excellent initialization for the iterative algorithms associated with the smoothness constraint problem formulation. Finally, the usefulness of the approach should extend to a wide variety of ill-posed inverse problems in which variational techniques seeking a "smooth" solution are generally used.

The authors describe and demonstrate a hierarchical reconstruction algorithm for use in noisy and limited-angle or sparse-angle tomography. The algorithm estimates an object's mass, center of mass, and convex hull from the available projections, and uses this information, along with fundamental mathematical constraints, to estimate a full set of smoothed projections. The mass and center of mass estimates are made using a least squares estimator derived from the principles of consistency of the Radon transform. The convex hull estimate is produced by first estimating the positions of support lines of the object from each available projection and then estimating the overall convex hull using prior shape information. Estimating the position of two support lines from a single projection is accomplished using a generalized likelihood ratio technique for estimating jumps in linear systems. Results for simulated objects in a variety of measurement situations are shown, and several possible extensions to this work are discussed.

An overview is provided of the several components of a research effort aimed at the development of a theory of multiresolution stochastic modeling and associated techniques for optimal multiscale statistical signal and image processing. A natural framework for developing such a theory is the study of stochastic processes indexed by nodes on lattices or trees in which different depths in the tree or lattice correspond to different spatial scales in representing a signal or image. In particular, it is shown how the wavelet transform directly suggests such a modeling paradigm. This perspective then leads directly to the investigation of several classes of dynamic models and related notions of multiscale stationarity in which scale plays the role of a time-like variable. The investigation of models on homogeneous trees is emphasized. The framework examined here allows for consideration, in a very natural way, of the fusion of data from sensors with differing resolutions. Also, thanks to the fact that wavelet transforms do an excellent job of 'compressing' large classes of covariance kernels, it is seen that these modeling paradigms appear to have promise in a far broader context than one might expect.>

The authors introduce a class of multiscale stochastic processes which are Markov in scale and which are characterized by dynamic state models evolving in scale. The models for these processes are motivated by the theory of multiscale representations and the wavelet transform. The authors formulate an optimal estimation problem based on these models, which has potential applications to sensor fusion problems where there exist data from sensors of differing resolution, and provide an efficient algorithm based on the wavelet transform. They give examples applying these models to first-order Gauss-Markov processes.>

In this paper we present several algorithms for reconstructing 2D convex sets given support line measurements for which the angles are known precisely but the lateral displacements are noisy. We extend the algorithms given in a previous paper by explicitly incorporating prior information about the shape of the objects to be reconstructed. We develop the Scale-Invariant algorithms, which incorporate prior shape information by defining prior probabilities on support vectors, where a support vector is a vector formed from the lateral displacements of a particular set of support lines of an object. We also develop the Ellipse-Based algorithms, which either assume or jointly estimate the parameters of an ellipse, given prior distributions that favor ellipses. In order to relate the support vector prior probability to the expected shape of an object we develop a vector decomposition called the Size/Shape/Shift decomposition, which helps to provide insight into the detailed geometric relationship between support vectors and 2D convex objects. We then use the maximum a posteriori criterion to determine the specific form of the support vector estimator. The computations involve a quadratic programming optimization stage, which is used to determine one component of the decomposition, and either a line search or a conjugate gradient stage, which is used to determine the remaining components. The performance of the algorithms is demonstrated using simulated support line measurements of an ellipse.

article Stability and stabilizability of discrete event dynamic systems Share on Authors: Cüneyt M. Özveren Massachusetts Institute of Technology, Cambridge Massachusetts Institute of Technology, CambridgeView Profile , Alan S. Willsky Massachusetts Institute of Technology, Cambridge Massachusetts Institute of Technology, CambridgeView Profile , Panos J. Antsaklis Univ. of Notre Dame, Notre Dame, IN Univ. of Notre Dame, Notre Dame, INView Profile Authors Info & Claims Journal of the ACMVolume 38Issue 3July 1991 pp 729–751https://doi.org/10.1145/116825.116855Online:01 July 1991Publication History 92citation799DownloadsMetricsTotal Citations92Total Downloads799Last 12 Months6Last 6 weeks0 Get Citation AlertsNew Citation Alert added!This alert has been successfully added and will be sent to:You will be notified whenever a record that you have chosen has been cited.To manage your alert preferences, click on the button below.Manage my AlertsNew Citation Alert!Please log in to your account Save to BinderSave to BinderCreate a New BinderNameCancelCreateExport CitationPublisher SiteGet Access

Efficient recursive smoothing algorithms are developed for isotropic random fields that can be obtained by passing white noise through rational filters. The estimation problem is shown to be equivalent to a countably infinite set of 1-D separable two-point boundary value smoothing problems. The 1-D smoothing problems are solved using a Markovianization approach followed by a standard 1-D smoothing algorithm. The desired field estimate is then obtained as properly weighted sum of the 1-D smoothed estimates. The 1-D two-point boundary value problems are also shown to have the same asymptotic properties and yield a stable spectral factorization of the power spectrum of the isotropic random fields.>

Algorithms are proposed for reconstructing convex sets given noisy support line measurements. It is observed that a set of measured support lines may not be consistent with any set in the plane. A theory of consistent support lines which serves as a basis for reconstruction algorithms that take the form of constrained optimization algorithms is developed. The formal statement of the problem and constraints reveals a rich geometry that makes it possible to include prior information about object position and boundary smoothness. The algorithms, which use explicit noise models and prior knowledge, are based on maximum-likelihood and maximum a posteriori estimation principles and are implemented using efficient linear and quadratic programming codes. Experimental results are presented. This research sets the stage for a more general approach to the incorporation of prior information concerning the estimation of object shape.>

The authors describe and demonstrate a hierarchical reconstruction algorithm for use in noisy and limited-angle or sparse-angle tomography. The algorithm estimates the object's mass, center of mass, and convex hull from the available projections, and uses this information, along with fundamental mathematical constraints, to estimate a full set of smoothed projections. The mass and center of mass are estimated using a maximum-likelihood (ML) estimator derived from the principles of consistency of the Radon transform. The convex hull estimate is produced by first estimating the positions of support lines of the object from each available projection and then estimating the overall convex hull using ML or maximum a posteriori (MAP) techniques. The position of two support lines from a single projection is estimated using either a generalized likelihood ratio technique for estimating jumps in linear systems or a support-width penalty method that uses Akaike's model-order estimation technique.>

Presents a method to reconstruct images from finite sets of noisy projections which are available only over limited or sparse angles. The method solves a constrained optimization problem to find a maximum a posteriori (MAP) estimate of the full 2-D Radon transform of the object, using prior knowledge of object mass, center of mass, and convex support, and information about fundamental constraints and smoothness of the Radon transform. This efficient primal-dual algorithm consists of an iterative local relaxation stage which solves a partial differential equation in Radon-space, followed by a simple Lagrange multiplier update stage. The object is reconstructed using convolution backprojection applied to the Radon transform estimate.>

A novel linear maximum-entropy method (MEM) spectral-estimation algorithm for 2-D isotropic random fields is presented. This procedure differs from pervious 2-D MEM algorithms by the fact that maximum advantage is taken of the symmetries implied by isotropy. It is shown that the isotropic MEM problem has a linear solution and that it is equivalent to the problem of constructing the optimal linear filter for estimating the underlying isotropic field at a point on the boundary of a disk of radius R, given noisy measurements of the field inside the disk. A fast algorithm for computing the estimation filter is then used to obtain the MEM spectral estimate.>

A new algorithm for the hierarchical aggregation of singularly perturbed finite-state Markov processes is derived. The approach taken bridges the gap between conceptually simple results for a relatively restricted class of processes and the significantly more complex results for the general case. The critical role played by (almost) transient states is exposed, resulting in a straightforward algorithm for the construction of a sequence of aggregate generators associated with various time scales. These generators together provide a uniform asymptotic approximation of the original probability transition function.

An efficient algorithm for computing the maximum-likelihood estimates of multiple signals observed by an array of sensors is presented. The algorithm provides estimates of parameters related to the directional patterns of the sources as well as estimates of the location parameters of the sources. Furthermore, the algorithm is equally applicable to wideband sources and narrowband sources and does not require a knowledge of the statistical properties of the signals.>

Sampling theorems are developed for isotropic random fields and their associated Fourier coefficient processes. A wave-number-limited isotropic random field z(r) is considered whose spectral density function is zero outside a disk of radius B centered at the origin of the wavenumber plane. z(r) can be reconstructed in the mean-square sense from its observation on the countable number of circles with radius r/sub i/=i pi /B, i in N, or of radius r/sub i/=a/sub i,n//B, i in N, where a/sub 1,n/ denotes the ith zero of the nth-order Bessel function J/sub n/(x), and n is arbitrary.>

Martingale decomposition techniques are used to derive Markovian models for the error in smoothed estimates of processes described by linear models driven by white noise. These models, together with some simple Hilbert space decomposition ideas, provide a simple unified framework for examining a variety of problems involving the efficient assimilation of spatial data, which we refer to as mapping problems. Algorithms for several different mapping problems are derived. A specific example of map updating for a two-dimensional random field is included.

The problem of detecting, locating and characterizing objects in a 2D cross-section from noisy projection data has been considered recently [1-3], in which objects are characterized by a finite number of parameters, which are estimated directly from noisy projection measurements. In this paper, the problem of maximum likelihood (ML) estimation of those parameters characterizing the geometry of an object (e.g. size and orientation) is considered, and estimation performance is investigated.

This paper considers the problem of observing a 2D function via its 1D projections (Radon transform); it presents a framework for detecting, locating and describing objects contained within a 2D cross-section by using noisy measurements of the Radon transform directly, rather than post-processing a reconstructed image. This framework offers the potential for significant improvements in applications where (1) attempts to perform an initial inversion with insufficient measurement data result in severely degraded reconstructions, and (2) the ultimate goal of the process is to obtain several specific pieces of information about the cross-section. To illustrate this perspective, we focus our attention on the problem of obtaining maximum-likelihood (ML) estimates of the parameters characterizing a single random object situated within a deterministic background medium, and we investigate the performance, robustness, and computational structure of the ML estimation procedure.

We consider a class of matrix signal processes that are received in the presence of multiplicative observation noise. By examining the differential version of the observation, we are able to derive finite-dimensional optimal detection-estimation equations that involve a linear filter with gain computed on-line using the incoming observations. An example involving the detection of an actuator failure on a rotating rigid body is considered.

A wide variety of continuous- and discrete-time estimation problems on the circleS^1are considered with the aid of Fourier series analysis. Measurement and diffusion update equations are derived for the conditional expectation of certain functions of the parameter to be estimated, and we investigate the use of Fourier series to obtain easily implemented optimal estimation equations. A variety of important examples--phase tracking, frequency demodulation, and phase demodulation in the presence of oscillator instabilities, additive noise, Rayleigh fading, or any combination of these--are considered.

The practical implementation of the infinite-dimensional optimal estimation results presented in Part I of this series is considered. Several techniques are described in detail. Included among these is the so-called "assumed density" approximation technique. Finite-dimensional suboptimal filtering equations based on this method are derived for several of the phase-tracking/demodulation problems studied in Part I. Finally, these techniques are applied to a phase tracking problem of importance in navigation systems such as Omega, and simulation results are reported that favorably compare a system designed using these techniques to an optimal phase-lock loop and an optimal linear system.