Patent Publication Number: US-10782396-B2

Title: Variational track management

Description:
RELATED APPLICATIONS 
     This application claims priority from U.S. patent application Ser. No. 14/563,824, filed 8 Dec. 2014, which is incorporated herein in its entirety. 
     BACKGROUND OF THE INVENTION 
     Field of the Invention 
     The invention relates generally to object tracking, and more specifically, to the use of variational track management in tracking moving objects. 
     Background of the Invention 
     Object tracking is the use of measurements of position, velocity, acceleration, and other kinematic parameters, to obtain an accurate determination of the positions of one or more moving objects. The field of object tracking is broad and possesses many applications, particularly in radar/sonar, robotics, and computer vision. When tracing an object, a system constructs a track that estimates the trajectory of the object over time. Where multiple objects are tracked, object tracking generally involves three steps. The first step is a data assignment process, in which various measurements are matched to existing tracks. This matching can often be supplemented by non-kinematic features of the plurality of moving objects, such as a radar cross-section in radar, or a color or shape in visual spectrum tracking, to more easily distinguish among the moving objects when the kinematic data is ambiguous. A second step is a filtering process, in which the assigned data is used to update the positions of the various objects to reduce or eliminate measurement error. Finally, a track management component utilizes various heuristics to remove or add tracks as needed due to clutter, or spurious measurements, or objects arriving in or leaving a region of interest. While there will be communication among the components performing these three functions, they generally act substantially independently. 
     SUMMARY OF THE INVENTION 
     In accordance with an aspect of the present invention, a method is provided for tracking moving objects from a set of measurements provided by an associated sensor system. An estimate of a posterior probability distribution for a plurality of track states is determined from an estimate of the posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks representing trajectories of the plurality of moving objects and the set of measurements. A new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks is determined from the set of measurements and the estimate of a posterior probability distribution for a plurality of track states. A variational lower bound, representing a lower bound for the marginal probability of the set of measurements given a model defined by the new estimate of the posterior probability distribution and the estimate of a posterior probability distribution for a plurality of track states, is determined from the new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks, the estimate of a posterior probability distribution for a plurality of track states, and the set of measurements. The steps of determining an estimate of a posterior probability distribution for a plurality of track states, determining a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks, and determining a variational lower bound are iteratively repeated until the variational lower bound is less than a threshold value. 
     In accordance with another aspect of the present invention, a system is provided for tracking a plurality of moving objects. A sensor system is configured to provide a set of measurements representing at least respective positions of the plurality of moving objects. A track state updating component is configured to determine an estimate of a posterior probability distribution for a plurality of track states from an estimate of the posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks representing trajectories of the plurality of moving objects and the set of measurements. A track assignment updating component is configured to determine a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks from the set of measurements and the estimate of a posterior probability distribution for the plurality of track states. A lower bound computation component is configured to compute a variational lower bound, representing a lower bound for the marginal probability of the set of measurements given a model defined by the new estimate of the posterior probability distribution and the estimate of a posterior probability distribution for a plurality of track states, from the new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks, the estimate of a posterior probability distribution for a plurality of track states, and the set of measurements. Each of the track state updating component, the track assignment updating component, and the lower bound computation component collectively perform an iterative determination of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks and the posterior probability distribution for a plurality of track states until the variational lower bound is less than a threshold value. 
     In accordance with yet another aspect of the present invention, a non-transitory computer readable medium stores executable instructions readable by an associated processor to perform a method for tracking at least respective positions of a plurality of moving objects from a set of measurements provided by an associated radar system. An estimate of a posterior probability distribution for a plurality of track states is determined from an estimate of the posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks representing trajectories of the plurality of moving objects and the set of measurements. The plurality of track states include a position state for each track and a metastate for each track representing an active or dormant status for the track. A new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks is determined from the set of measurements and the estimate of a posterior probability distribution for a plurality of track states. A variational lower bound, representing a lower bound for the marginal probability of the set of measurements given a model defined by the new estimate of the posterior probability distribution and the estimate of a posterior probability distribution for a plurality of track states, is determined from the new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks, the estimate of a posterior probability distribution for a plurality of track states, and the set of measurements. The steps of determining an estimate of a posterior probability distribution for a plurality of track states, determining a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks, and determining a variational lower bound are iteratively repeated until the variational lower bound is less than a threshold value. 
    
    
     
       BRIEF DESCRIPTION OF THE DRAWINGS 
       The features, objects, and advantages of the invention will become more apparent from the detailed description set forth below when taken in conjunction with the drawings, wherein: 
         FIG. 1  illustrates a system for tracking a plurality of moving objects; 
         FIG. 2  illustrates one example of a method for tracking moving objects; 
         FIG. 3  illustrates one example of a method for estimating the posterior probability distribution for a plurality of track states; 
         FIG. 4  illustrates a method for generating a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks; 
         FIG. 5  illustrates a method for calculating an estimate of the posterior probability distribution for the track assignment using a loopy belief propagation process; and 
         FIG. 6  illustrates a schematic block diagram illustrating an exemplary system of hardware components capable of implementing examples of the systems and methods disclosed in  FIGS. 1-5 . 
     
    
    
     DETAILED DESCRIPTION OF THE INVENTION 
     Most tracking algorithms treat a data association process and a filtering process as separate, alternating sequential steps. This is due to the vast number of possible data associations in situations where there are even several tracks and corresponding measurements. Data association is usually treated as a discrete optimization problem, where measurements over just a few time steps, generally no more than three, are assigned to the tracks. In the variational tracking algorithm described herein, data association and filtering are performed together to approximate the optimal state of the track given all previous data and all possible associations. To this end, a variational Bayes (VB) technique is applied to approximate the solution to a corresponding probabilistic graphical model of the tracking problem, providing an algorithm that scales linearly with the number of time steps. Current algorithms, such as multiple hypothesis trackers (MHT), instead scale exponentially with the number of time steps. With the linear scaling of the variational tracking algorithm, an approximation of the states of tens of tracks, marginalized over all possible data associations, for as many as thirty time steps, is achievable. Due to the difficulties imposed by the data association constraints, the algorithm uses a hybrid method utilizing both variational and loopy belief propagation components to make inference efficient and tractable. Further, the algorithm incorporates track management into the model itself as a series of track metastates. This allows for a significant performance gains over tradition heuristic approaches. 
       FIG. 1  illustrates a system  10  for tracking a plurality of moving objects. The system  10  utilizes a probabilistic tracking algorithm that incorporates combinatorial data association constraints and model-based track management using variational Bayes. A Bethe entropy approximation is used to incorporate data association constraints that are generally ignored in existing probabilistic tracking algorithms. The variational Bayes assignment posterior utilized in the system has an induced factorization over time with regard to assignment matrices, allowing the computational cost of the variational approach to be linear in window length as opposed to the exponential cost of existing approaches. Compared to existing methods, such as maximum a posterior inference (MAP), the current system utilizes a weaker approximation for the posterior distributions, but provides a significant increase in efficiency and accuracy. The system also utilizes new approximate inference methods on assignment matrices and a new conjugate assignment prior (CAP). 
     The system includes a sensor system  12  configured to provide a set of measurements representing at least respective positions of the plurality of moving objects. The sensor system  12  can include, for example, a camera operating in visible, infrared, or ultraviolet light, or a transceiver for receiving radio frequency or microwave radiation. In one implementation, the sensor system  12  is a radar system. To provide a starting point for the iterative analysis, an initialization component  14  is configured to initialize the posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks representing the moving objects. In the illustrated implementation, a default distribution with random perturbations is used. 
     A track state updating component  20  is configured to determine an estimate of a posterior probability distribution for a plurality of track states from an estimate of the posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks representing trajectories of the plurality of moving objects and the set of measurements. In the illustrated example, the track state updating component  20  updates, for each track, an associated position state and an associated metastate representing the activity or dormancy of a track, such that the track update component manages a status of the plurality of tracks. Many tracking algorithms unrealistically assume that the number of tracks is known apriori and fixed, with additional “wrapper logic” placed around the trackers to initiate and destroy tracks. This logic involves the use of heuristics such as M-of-N logic. By incorporating the track management into the track update component in a model-based manner, significant performance gains can be realized, as the system  10  simultaneously performs inference for track management, data association, and state estimation. 
     To this end, the track state updating component  20  includes a soft hidden Markov model (sHMM)  22  for updating the metastates of the plurality of tracks. For example, detection probabilities, representing the likelihood that the track is active at the time step, can be determined for each time step, for example, and a forward-backward algorithm, or similar analysis, can be applied to the model and the determined detection probabilities to update the metastates. For the position states, the track state updating component  20  can include a Kalman smoother  24 . As part of the position state update, the track state updating component  20  can calculate a set of pseudomeasurements from a set of measurements and the estimated posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks. The Kalman smoother  24  can be applied to the pseudomeasurements to calculate an update for the plurality of track states. 
     Once the estimate of the posterior distribution for the track states is determined, is it provided to each of a track assignment updating component  30  and a lower bound computation component  40 . The track assignment updating component  30  is configured to determine a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks from the set of measurements and the estimate of a posterior probability distribution for a plurality of track states. The illustrated implementation uses a Bethe entropy approximation to incorporate the one-to-one restrictions of the track assignment process. Accordingly, the track assignment updating component  30  utilizes a loopy belief propagation (LBP) algorithm  32  to calculate the new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks. Once the estimate of the assignment posterior distribution is updated, it is provided to each of the track state updating component  20  and the lower bound computation component  40 . The track state updating component  20  can then generate a new estimate for the posterior distribution for the track states using the improved estimate for the assignment posterior distribution. 
     The lower bound computation component  40  is configured to compute a variational lower bound, representing a lower bound for the marginal probability of the set of measurements given a model defined by the new estimate of the posterior probability distribution and the estimate of a posterior probability distribution for a plurality of track states. Each of the track state updating component  20  and the track assignment updating component  30  are responsive to the lower bound computation component  40 , such that the estimates of the posterior distributions determined at these components are accepted once the lower bound falls below a threshold value. Accordingly, the track state updating component  20 , the track assignment updating component  30 , and the lower bound computation component  40  collectively perform an iterative determination of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks and the posterior probability distribution for a plurality of track states until the variational lower bound is less than a threshold value. 
     In view of the foregoing structural and functional features described above, methods in accordance with various aspects of the present invention will be better appreciated with reference to  FIGS. 2-5 . While, for purposes of simplicity of explanation, the methods of  FIGS. 2-5  are shown and described as executing serially, it is to be understood and appreciated that the present invention is not limited by the illustrated order, as some aspects could, in accordance with the present invention, occur in different orders and/or concurrently with other aspects from that shown and described herein. Moreover, not all illustrated features may be required to implement a method in accordance with an aspect the present invention. 
       FIG. 2  illustrates one example of a method  100  for tracking moving objects from a set of measurements provided by an associated sensor system.  FIGS. 3-5  describe a specific implementation of the various steps  110 ,  130 , and  150  of  FIG. 2 . To better describe this specific implementation, it is helpful to establish a standard notation. In this implementation, the set of measurements comprises measurements of N 1  discrete time steps, referred to as frames, where N 1  is an integer greater than one, such that each object of the plurality of objects is expected to be represented by N 1  individual measurements. Specifically, for each time step, k∈N 0 , N z (k)∈N 1  measurements in a matrix Z k ={z j,k } j=1   N     z     (k)  are observed from both real objects and clutter (spurious measurements). In the illustrated example, z j,k ∈Z is a vector of position measurements in three-dimensional space. In data association, an assignment matrix A is estimated, where A ij =1 if and only if track i is associated with track j. The assignment matrix is constrained such that each track can be associated with only one measurement, and each measurement can be associated with only one track, such that: 
     
       
         
           
             
               
                 
                   
                     
                       
                         ∑ 
                         
                           i 
                           = 
                           0 
                         
                         
                           N 
                           T 
                         
                       
                       ⁢ 
                       
                         A 
                         ij 
                       
                     
                     = 
                     1 
                   
                   , 
                   
                     
                       j 
                       ∈ 
                       1 
                     
                     : 
                     
                       N 
                       Z 
                     
                   
                   , 
                   
                     
                       
                         ∑ 
                         
                           j 
                           = 
                           0 
                         
                         
                           N 
                           Z 
                         
                       
                       ⁢ 
                       
                         A 
                         ij 
                       
                     
                     = 
                     1 
                   
                   , 
                   
                     
                       i 
                       ∈ 
                       1 
                     
                     : 
                     
                       N 
                       T 
                     
                   
                   , 
                   
                     
                       A 
                       00 
                     
                     = 
                     0. 
                   
                 
               
               
                 
                   ( 
                   1 
                   ) 
                 
               
             
           
         
       
     
     The zero indices of A∈{0,1} N     i     +1xN     z     +1  are the “dummy row” and “dummy column” to represent the assignment of a measurement to clutter and the assignment of a track to a missed detection. 
     The inventors have determined that the cost functions utilized in the maximum a posteriori optimization (MAP) in multiple hypothesis tracking assume that a number of tracks, N T  is known a priori and the N Z  is random. The corresponding generative process on assignment matrices is therefore implemented by starting with a one-to-one mapping from measurement to tracks, A←I N     T     xN     T   , with each track observed with a probability P D ∈[0,1] N     T   . Only the columns of detected tracks are kept, such that A←A(⋅,d),d i ·Bernoulli (P D (i)). A Poisson number of clutter measurements is sampled: A←[A,0 N     Y     xN     c   ], N C ˜Poisson(λ). A random permutation vector, π, is used to make the order of the measurements arbitrary, A←A(⋅,π), and a dummy row and column are added to the assignment matrix to satisfy the constraints in (1). This provides a normalized prior on assignments: 
     
       
         
           
             
               
                 
                   
                     
                       P 
                       ⁡ 
                       
                         ( 
                         
                           A 
                           | 
                           
                             P 
                             D 
                           
                         
                         ) 
                       
                     
                     = 
                     
                       
                         λ 
                         
                           N 
                           c 
                         
                       
                       ⁢ 
                       
                         
                           exp 
                           ⁡ 
                           
                             ( 
                             
                               - 
                               λ 
                             
                             ) 
                           
                         
                         / 
                         
                           
                             N 
                             Z 
                           
                           ! 
                         
                       
                       ⁢ 
                       
                         
                           ∏ 
                           
                             i 
                             = 
                             1 
                           
                           
                             N 
                             T 
                           
                         
                         ⁢ 
                         
                             
                         
                         ⁢ 
                         
                           
                             
                               
                                 P 
                                 D 
                               
                               ⁡ 
                               
                                 ( 
                                 i 
                                 ) 
                               
                             
                             
                               d 
                               i 
                             
                           
                           ⁢ 
                           
                             
                               ( 
                               
                                 1 
                                 - 
                                 
                                   
                                     P 
                                     D 
                                   
                                   ⁡ 
                                   
                                     ( 
                                     i 
                                     ) 
                                   
                                 
                               
                               ) 
                             
                             
                               1 
                               - 
                               
                                 d 
                                 i 
                               
                             
                           
                         
                       
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   2 
                   ) 
                 
               
             
           
         
       
     
     It will be appreciated that the detections d, N Z , and clutter measurement count, N c , are deterministic functions of A. 
     As mentioned previously, a track model for updating a plurality of track states is implemented as over K time steps, where K is a positive integer greater than one. A set of latent states, x 1:K ∈χ K , follow a Markov process, while the measurements, z 1:K ∈Z K , are iid conditional on the track state, such that: 
                       p   ⁡     (       z     1   :   K       ,     x     1   :   K         )       =       p   ⁡     (     x   1     )       ⁢       ∏     k   =   2     K     ⁢           ⁢       p   ⁡     (       x   k     |     x     k   -   1         )       ⁢       ∏     k   =   1     K     ⁢           ⁢     p   ⁡     (       z   k     |     x   k       )                 ,           (   3   )               
where the track and measurements indices i and j have been omitted for brevity. Although the systems and methods taught herein can be used with more general models, in the illustrated example, each track independently follows a linear system, such as a Kalman filter, such that:
 
 p ( x   k   |x   k-1 )= ( x   k   |Fx   k-1   ,Q ) p ( z   k   |x   k )= ( z   k   |Hx   k   ,R )  (4)
 
     The track states can be augmented with a Markov model with an active/dormant metastate, s k , in a 1-in-N encoding to address track management, such that: 
     
       
         
           
             
               
                 
                   
                     
                       P 
                       ⁡ 
                       
                         ( 
                         
                           s 
                           
                             1 
                             : 
                             K 
                           
                         
                         ) 
                       
                     
                     = 
                     
                       
                         P 
                         ⁡ 
                         
                           ( 
                           
                             s 
                             1 
                           
                           ) 
                         
                       
                       ⁢ 
                       
                         
                           ∏ 
                           
                             k 
                             = 
                             2 
                           
                           K 
                         
                         ⁢ 
                         
                             
                         
                         ⁢ 
                         
                           P 
                           ⁡ 
                           
                             ( 
                             
                               
                                 s 
                                 k 
                               
                               | 
                               
                                 s 
                                 
                                   k 
                                   - 
                                   1 
                                 
                               
                             
                             ) 
                           
                         
                       
                     
                   
                   , 
                   
                     
                       s 
                       k 
                     
                     ∈ 
                     
                       
                         
                           { 
                           
                             0 
                             , 
                             1 
                           
                           } 
                         
                         
                           N 
                           s 
                         
                       
                       . 
                     
                   
                 
               
               
                 
                   ( 
                   5 
                   ) 
                 
               
             
           
         
       
     
     This model allows the system to handle an unknown number of tracks by making N T  arbitrarily large. The detection probability, P D , is now a function of s with a very small P D  in the dormant state and a larger P D  in the active state. Extensions with a larger number of states N S  can be implemented using the teachings provided herein. The collection of track metastates over all tracks at frame k are referred to herein as S k :={s i,k } i=1   N     T   , and likewise, X k :={x i,k } i=1   N     T   . 
     Combining the assignment process and track models to get the full model joint: 
                       p   ⁡     (       Z     1   :   K       ,     X     1   :   K       ,     A     1   :   K       ,     S     1   :   K         )       =         ∏     k   =   1     K     ⁢           ⁢     □   ⁢           ⁢     p   ⁡     (         Z   k     |     X   k       ,     A   k       )       ⁢     p   ⁡     (       X   k     |     X     k   -   1         )       ⁢     P   ⁡     (       S   k     |     S     k   -   1         )       ⁢     P   ⁡     (       A   k     |     S   k       )           =       ∏     k   =   1     K     ⁢           ⁢       P   ⁡     (       A   k     |     S   k       )       ·       ∏     i   =   1       N   T       ⁢           ⁢       p   ⁡     (       x     i   ,   k       |     x     i   ,     k   -   1           )       ⁢       P   ⁡     (       s     i   ,   k       |     s     i   ,     k   -   1           )       ·       ∏     j   =   1         N   z     ⁡     (   k   )         ⁢           ⁢           p   0     ⁡     (     z     j   ,   k       )         A   oj   k       ⁢       ∏     i   =   1       N   T       ⁢           ⁢       p   ⁡     (         z     j   ,   k       |     x     i   ,   k         ,       A   ij   k     =   1       )         A   ij   k                           ,           (   6   )               
where p 0  is the clutter distribution, which is often a uniform distribution. The traditional goal in tracking is to compute p(X k |Z 1:k ), the exact computation of which is intractable due to the “combinatorial explosion” in summing out the assignments A 1:k . The methods presented herein provide an improved approximation of p(X k |Z 1:k ).
 
     To compute the posterior, P(A|Z), it is helpful to determine what conjugate priors on A are possible. Deriving approximate inference procedures is often greatly simplified if the prior on the parameters is conjugate to the complete data likelihood: p(Z, X|A). The conjugate prior for an exponential family (EF) complete likelihood can be derived as: 
                       p   ⁡     (     Z   ,     X   |   A       )       =         ∏     j   =   1       N   Z       ⁢           ⁢           p   0     ⁡     (     z   j     )         A   oj       ⁢       ∏     i   =   1       N   T       ⁢           ⁢         p   ⁡     (         z   j     |     x   i       ,       A   ij     =   1       )         A   ij       ⁢       ∏     i   =   1       N   T       ⁢           ⁢     p   ⁡     (     x   i     )                 =       ∏     i   =   1       N   T       ⁢           ⁢       p   ⁡     (     x   i     )       ⁢     exp   ⁡     (       1   T     ⁢     (     A   ⊙   L     )     ⁢   1     )               ,     
     ⁢           ⁢       L   ij     :=     log   ⁢           ⁢     p   ⁡     (         z   j     |     x   i       ,       A   ij     =   1       )           ,       L     i   ⁢           ⁢   0       :=   0     ,       L     0   ⁢           ⁢   j       :=     log   ⁢           ⁢       p   0     ⁡     (     z   j     )           ,           (   7   )               
where the matrix L∈   N     T     +1×N     z     +1  represents log likelihood contributions from various assignments. The set of log likelihood values can include contributions from at least one likelihood value associated with a non-kinematic feature, such that as a radar cross-section. This gives the following exponential family quantities: base measure
 
                 h   ⁡     (     Z   ,   X     )       =       ∏     i   =   1       N   T       ⁢           ⁢     p   ⁡     (     x   i     )           ,         
partition function g(A)=1, natural parameters η(A)=vec A, and sufficient statistics S(Z, X)=vec L. This implies the conjugate assignments prior (CAP) for P(A|χ):
 
                   CAP   (       A   ⁢        𝓍   )       :=         𝒮   ⁡     (   𝓍   )         -   1       ⁢              {     A   ∈   𝒜     }     ⁢     exp   ⁡     (       1   T     ⁢     (     𝓍   ⊙   A     )     ⁢   1     )         ,       𝒮   ⁡     (   𝓍   )       :=       ∑     A   ∈   𝒜       ⁢     exp   ⁡     (       1   T     ⁢     (     𝓍   ⊙   A     )     ⁢   1     )           ,                   (   8   )               
where   is the set of all assignment matrices that obey the one-to-one constraints in (1). Note that χ is a function of the track metastates S. The assignment prior of (2) can be recovered in the form of the CAP distribution of (8) via the following parameter settings, with σ(⋅) denoting the logistic:
 
     
       
         
           
             
               
                 
                   
                     
                       𝓍 
                       ij 
                     
                     = 
                     
                       
                         log 
                         ⁡ 
                         
                           ( 
                           
                             
                               
                                 P 
                                 D 
                               
                               ⁡ 
                               
                                 ( 
                                 i 
                                 ) 
                               
                             
                             
                               
                                 ( 
                                 
                                   1 
                                   - 
                                   
                                     
                                       P 
                                       D 
                                     
                                     ⁡ 
                                     
                                       ( 
                                       i 
                                       ) 
                                     
                                   
                                 
                                 ) 
                               
                               ⁢ 
                               λ 
                             
                           
                           ) 
                         
                       
                       = 
                       
                         
                           
                             σ 
                             
                               - 
                               1 
                             
                           
                           ⁡ 
                           
                             ( 
                             
                               
                                 P 
                                 D 
                               
                               ⁡ 
                               
                                 ( 
                                 i 
                                 ) 
                               
                             
                             ) 
                           
                         
                         - 
                         
                           log 
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           λ 
                         
                       
                     
                   
                   , 
                   
                     
 
                   
                   ⁢ 
                   
                     i 
                     ∈ 
                     
                       1 
                       ⁢ 
                       
                         : 
                       
                       ⁢ 
                       
                         N 
                         T 
                       
                     
                   
                   , 
                   
                     j 
                     ∈ 
                     
                       1 
                       ⁢ 
                       
                         : 
                       
                       ⁢ 
                       
                         N 
                         Z 
                       
                     
                   
                   , 
                   
                     
                       𝓍 
                       
                         0 
                         ⁢ 
                         
                             
                         
                         ⁢ 
                         j 
                       
                     
                     = 
                     
                       
                         𝓍 
                         
                           i 
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           0 
                         
                       
                       = 
                       0. 
                     
                   
                 
               
               
                 
                   ( 
                   9 
                   ) 
                 
               
             
           
         
       
     
     Due of the symmetries in the prior of (9), the CAP distribution of (8) can be analytically normalized in this instance: 
     
       
         
           
             
               
                 
                   
                     
                       S 
                       ⁡ 
                       
                         ( 
                         𝓍 
                         ) 
                       
                     
                     
                       - 
                       1 
                     
                   
                   = 
                   
                     
                       P 
                       ⁡ 
                       
                         ( 
                         
                           
                             A 
                             
                               
                                 1 
                                 : 
                                 
                                   N 
                                   T 
                                 
                               
                               , 
                               
                                 1 
                                 : 
                                 
                                   N 
                                   Z 
                                 
                               
                             
                           
                           = 
                           0 
                         
                         ) 
                       
                     
                     = 
                     
                       Poisson 
                       ( 
                       
                         
                           N 
                           Z 
                         
                         ⁢ 
                         
                            
                           λ 
                           ) 
                         
                         ⁢ 
                         
                             
                         
                         ⁢ 
                         
                           
                             ∏ 
                             
                               i 
                               = 
                               1 
                             
                             
                               N 
                               T 
                             
                           
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           
                             
                               ( 
                               
                                 1 
                                 - 
                                 
                                   
                                     P 
                                     D 
                                   
                                   ⁡ 
                                   
                                     ( 
                                     i 
                                     ) 
                                   
                                 
                               
                               ) 
                             
                             . 
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   10 
                   ) 
                 
               
             
           
         
       
     
     Given that the dummy row and columns of χ are zero in (9), the normalization of equation (10) allows the CAP distribution of (8) to match (2) for the 0 assignment case. 
     Although the conjugate prior (8) allows for “computation” of the posterior, χ posterior =χ prior +L=prior+L, computing  A] or S(χ) remains difficult. As will be explained in detail below, the illustrated example utilizes a loopy belief propagation algorithm to facilitate this determination. 
     The improved variational tracker of the illustrated method enforces the factorization constraint that the posterior factorizes across assignment matrices and latent track states:
 
 p ( A   1:K   ,X   1:K   ,S   1:K   |Z   1:K )≈ q ( A   1:K   ,X   1:K   ,S   1:K )= a ( A   1:K ) q ( X   1:K   ,S   1:K ).  (11)
 
     In this formation, A can be conceptualized as the parameters and X and S can be conceptualized as the latent variables, and these two groups of variables can be factorized, giving a variational lower bound:
 
 ( q )=   q [log  p ( Z   1:K   ,X   1:K   ,A   1:K   ,S   1:K )]+ H [ q ( X   1:K   ,S   1:K )]+ H [ q ( A   1:K )],  (12)
 
where H[⋅] represents the Shannon entropy. From the variational Bayes lower bound in (12) and the model in (6), an induced factorization can be determined without further approximation:
 
     
       
         
           
             
               
                 
                   
                     
                       q 
                       ⁡ 
                       
                         ( 
                         
                           A 
                           
                             1 
                             : 
                             K 
                           
                         
                         ) 
                       
                     
                     = 
                     
                       
                         ∏ 
                         
                           k 
                           = 
                           1 
                         
                         K 
                       
                       ⁢ 
                       
                           
                       
                       ⁢ 
                       
                         q 
                         ⁡ 
                         
                           ( 
                           
                             A 
                             k 
                           
                           ) 
                         
                       
                     
                   
                   , 
                   
                     
                       q 
                       ⁡ 
                       
                         ( 
                         
                           
                             X 
                             
                               1 
                               : 
                               K 
                             
                           
                           , 
                           
                             S 
                             
                               1 
                               : 
                               K 
                             
                           
                         
                         ) 
                       
                     
                     = 
                     
                       
                         ∏ 
                         
                           i 
                           = 
                           1 
                         
                         
                           N 
                           T 
                         
                       
                       ⁢ 
                       
                           
                       
                       ⁢ 
                       
                         
                           q 
                           ⁡ 
                           
                             ( 
                             
                               
                                 x 
                                 i 
                               
                               , 
                               ′ 
                             
                             ) 
                           
                         
                         ⁢ 
                         
                           q 
                           ⁡ 
                           
                             ( 
                             
                               
                                 s 
                                 i 
                               
                               , 
                               ′ 
                             
                             ) 
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   13 
                   ) 
                 
               
             
           
         
       
     
     It can be noted that the approximate posterior on assignment matrices factorizes over time, and the approximate posterior on latent states factorizes across tracks. 
     At  110 , an estimate of a posterior probability distribution is determined for a plurality of track states from an estimate of the posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks representing trajectories of the plurality of moving objects and the set of measurements. It will be appreciated that the method is iterative, such that an estimate of the posterior probability distribution for the plurality of possible assignments will usually be available from a previous step in the method. For the first step of the method, however, it will be necessary to generate an initial distribution. To this end, an estimate of the posterior probability distribution for the plurality of possible assignments can be initialized as a default distribution with random perturbations before the method  100  begins. 
       FIG. 3  illustrates one example of a method  110  for estimating the posterior probability distribution for a plurality of track states. In the illustrated example, the updates for the track states, x i , and the metastates, s i , can be determined from the induced factorization in (13), and each update can be derived separately. Accordingly, the method begins at  112 , where a next track is selected. It will be appreciated, of course, that in a first iteration of the method, the selected “next” track will be a first track. The method  110  then proceeds to  114 , where a next time step, or frame, is selected. Again, in a first iteration for each track, the selected time step is a first time step. 
     At  116 , a pseudomeasurement is calculated for the selected time step from a set of measurements and the estimated posterior probability distribution for a plurality of possible assignments of the set of measurements to a set of tracks. The variational updates for the track states representing position, x i , can be expressed as follows, denoting equality to within an additive constant with  : 
     
       
         
           
             
               
                 
                   
                     log 
                     ⁢ 
                     
                         
                     
                     ⁢ 
                     
                       q 
                       ⁡ 
                       
                         ( 
                         
                           
                             x 
                             i 
                           
                           , 
                           ′ 
                         
                         ) 
                       
                     
                   
                   ⁢ 
                   
                     = 
                     c 
                   
                   ⁢ 
                   
                     
                       log 
                       ⁢ 
                       
                           
                       
                       ⁢ 
                       
                         p 
                         ⁡ 
                         
                           ( 
                           
                             
                               x 
                               i 
                             
                             , 
                             ′ 
                           
                           ) 
                         
                       
                     
                     + 
                     
                       
                         ∑ 
                         
                           k 
                           = 
                           1 
                         
                         K 
                       
                       ⁢ 
                       
                         
                           ∑ 
                           
                             j 
                             = 
                             1 
                           
                           
                             
                               N 
                               Z 
                             
                             ⁡ 
                             
                               ( 
                               k 
                               ) 
                             
                           
                         
                         ⁢ 
                         
                           
                             𝔼 
                             ⁡ 
                             
                               [ 
                               
                                 A 
                                 ij 
                                 k 
                               
                               ] 
                             
                           
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           log 
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           
                             𝒩 
                             ( 
                             
                               
                                 z 
                                 
                                   j 
                                   , 
                                   k 
                                 
                               
                               ⁢ 
                               
                                  
                                 
                                   
                                     Hx 
                                     
                                       i 
                                       , 
                                       k 
                                     
                                   
                                   , 
                                   R 
                                 
                                 ) 
                               
                             
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   14 
                   ) 
                 
               
             
             
               
                 
                   
                       
                   
                   ⁢ 
                   
                     
                       ⟹ 
                       
                         q 
                         ⁡ 
                         
                           ( 
                           
                             
                               x 
                               i 
                             
                             , 
                             ′ 
                           
                           ) 
                         
                       
                     
                     ∝ 
                     
                       
                         p 
                         ⁡ 
                         
                           ( 
                           
                             
                               x 
                               i 
                             
                             , 
                             ′ 
                           
                           ) 
                         
                       
                       ⁢ 
                       
                         
                           ∏ 
                           
                             k 
                             = 
                             1 
                           
                           K 
                         
                         ⁢ 
                         
                           
                             ∏ 
                             
                               j 
                               = 
                               1 
                             
                             
                               
                                 N 
                                 Z 
                               
                               ⁡ 
                               
                                 ( 
                                 k 
                                 ) 
                               
                             
                           
                           ⁢ 
                           
                             𝒩 
                             ( 
                             
                               
                                 z 
                                 
                                   j 
                                   , 
                                   k 
                                 
                               
                               ⁢ 
                               
                                 
                                    
                                   
                                     
                                       Hx 
                                       
                                         i 
                                         , 
                                         k 
                                       
                                     
                                     , 
                                     
                                       R 
                                       / 
                                       
                                         𝔼 
                                         ⁡ 
                                         
                                           [ 
                                           
                                             A 
                                             ij 
                                             k 
                                           
                                           ] 
                                         
                                       
                                     
                                   
                                   ) 
                                 
                                 . 
                               
                             
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   15 
                   ) 
                 
               
             
           
         
       
     
     Using (6), this can be shown to proportional to: 
     
       
         
           
             
               
                 
                   
                     q 
                     ⁡ 
                     
                       ( 
                       
                         
                           x 
                           i 
                         
                         , 
                         ′ 
                       
                       ) 
                     
                   
                   ∝ 
                   
                     
                       p 
                       ⁡ 
                       
                         ( 
                         
                           
                             x 
                             i 
                           
                           , 
                           ′ 
                         
                         ) 
                       
                     
                     ⁢ 
                     
                       
                         ∏ 
                         
                           k 
                           = 
                           1 
                         
                         K 
                       
                       ⁢ 
                       
                         𝒩 
                         ( 
                         
                           
                             
                               
                                 z 
                                 ~ 
                               
                               
                                 i 
                                 , 
                                 k 
                               
                             
                             ⁢ 
                             
                                
                               
                                 
                                   Hx 
                                   
                                     i 
                                     , 
                                     k 
                                   
                                 
                                 , 
                                 
                                   R 
                                   / 
                                   
                                     𝔼 
                                     ⁡ 
                                     
                                       [ 
                                       
                                         d 
                                         
                                           i 
                                           , 
                                           k 
                                         
                                       
                                       ] 
                                     
                                   
                                 
                               
                               ) 
                             
                           
                           , 
                           
                             
                               
                                 z 
                                 ~ 
                               
                               
                                 i 
                                 , 
                                 k 
                               
                             
                             := 
                             
                               
                                 1 
                                 
                                   𝔼 
                                   ⁡ 
                                   
                                     [ 
                                     
                                       d 
                                       
                                         i 
                                         , 
                                         k 
                                       
                                     
                                     ] 
                                   
                                 
                               
                               ⁢ 
                               
                                 
                                   ∑ 
                                   
                                     j 
                                     = 
                                     1 
                                   
                                   
                                     N 
                                     Z 
                                   
                                 
                                 ⁢ 
                                 
                                   
                                     𝔼 
                                     ⁡ 
                                     
                                       [ 
                                       
                                         A 
                                         ij 
                                         k 
                                       
                                       ] 
                                     
                                   
                                   ⁢ 
                                   
                                     z 
                                     
                                       j 
                                       , 
                                       k 
                                       , 
                                     
                                   
                                 
                               
                             
                           
                           , 
                           
                             
 
                           
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           
                             
                               with 
                               ⁢ 
                               
                                   
                               
                               ⁢ 
                               
                                 𝔼 
                                 ⁡ 
                                 
                                   [ 
                                   
                                     d 
                                     
                                       i 
                                       , 
                                       k 
                                     
                                   
                                   ] 
                                 
                               
                             
                             = 
                             
                               
                                 1 
                                 - 
                                 
                                   𝔼 
                                   ⁡ 
                                   
                                     [ 
                                     
                                       A 
                                       
                                         i 
                                         ⁢ 
                                         
                                             
                                         
                                         ⁢ 
                                         0 
                                       
                                       k 
                                     
                                     ] 
                                   
                                 
                               
                               = 
                               
                                 
                                   ∑ 
                                   
                                     j 
                                     = 
                                     1 
                                   
                                   
                                     N 
                                     Z 
                                   
                                 
                                 ⁢ 
                                 
                                   
                                     𝔼 
                                     ⁡ 
                                     
                                       [ 
                                       
                                         A 
                                         ij 
                                         k 
                                       
                                       ] 
                                     
                                   
                                   . 
                                 
                               
                             
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   16 
                   ) 
                 
               
             
           
         
       
     
     In addition to the pseudomeasurements, {tilde over (z)} i,k , detection probabilities, P D , representing the likelihood that the track is active at the time step, can also be determined for each time step. The system then determines at  118  if values have been calculated for all time steps. If not (N), the method returns to  114  to select a next time step. Otherwise (Y), the method advances to  120 , where a track state for the selected track is determined from the pseudomeasurements for each track. The form of the posterior, q(x i , .), is equivalent to a linear dynamic system with pseudomeasurements, {tilde over (z)} i,k , and non-stationary measurement covariance R/ [d i,k ]. Accordingly, in the illustrated implementation, the state update, q(x i , .), can be determined via a Kalman smoothing process. 
     At  122 , a track metastate for each track is determined from the calculated probabilities for each track. The posterior on the track metastates can be expressed as: 
                     log   ⁢           ⁢     q   ⁡     (       s   i     ,   ′     )         ⁢     =   c     ⁢       log   ⁢           ⁢     P   ⁡     (       s   i     ,   ′     )         +       ∑     k   =   1     K     ⁢         𝔼   q     ⁡     (     A   k     )       [         log   ⁢           ⁢     P   (       A   k     ⁢          S   k     )       ]       ⁢     =   c     ⁢       log   ⁢           ⁢     P   ⁡     (       s   i     ,   ′     )         +       ∑     k   =   1     K     ⁢       s     i   ,   k     T     ⁢     ℓ     i   ,   k               ,                   (   17   )                       ⁢           ℓ     i   ,   k       ⁡     (   s   )       :=         𝔼   ⁡     [     d     i   ,   k       ]       ⁢           ⁢     log   ⁡     (       P   D     ⁡     (   s   )       )         +       (     1   -     𝔼   ⁡     [     d     i   ,   k       ]         )     ⁢           ⁢     log   ⁡     (     1   -       P   D     ⁡     (   s   )         )             ,     
     ⁢           ⁢     s   ∈     1   ⁢     :     ⁢     N   S                   (   18   )                       ⁢         ⟹     q   ⁡     (       s   i     ,   ′     )         ∝       P   ⁡     (       s   i     ,   ′     )       ⁢       ∏     k   =   1     K     ⁢     exp   ⁡     (       s     i   ,   k     T     ⁢     ℓ     i   ,   k         )             ,             (   19   )               
when (18) follows from (2), above. If P(s i , .) follows a Markov chain, then the form for q(s i , .) is the same as a hidden Markov model (HMM) with emission log likelihoods l i,k ∈[   − ] N     s   . Therefore, the metastate posterior q(s i , .) update can be implemented via a forward-backward inference algorithm on the calculated probabilities P D . The algorithm works in an online fashion using a sliding window. The system then determines at  124  if track states have been updated for all tracks. If not (N), the method returns to  112  to select a next track. Otherwise (Y), the method terminates, returning the posterior probability distribution for a plurality of track states.
 
     Returning to  FIG. 2 , at  130 , a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks is determined from the set of measurements and the estimate of a posterior probability distribution for a plurality of track states.  FIG. 4  illustrates a method  130  for generating a new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks. At  132 , a next time step is selected. It will be appreciated, of course, that in a first iteration of the method, the selected next time step will be a first time step. The method  130  then proceeds to  134 , where a next track is selected. Again, in a first iteration for each time step, the selected track is a first track. 
     At  136 , expected sufficient statistics, E q(Ak)  [Lk], are calculated for the selected track and time step from a set of likelihood values. At  138 , an expected value for a track state associated with the selected track at the selected time step, E q(S     k     ) [χ k ], is calculated from the expected sufficient statistics for that time step. For example, χ k  can be calculated using (9) above. It is then determined at  140  if all of the tracks have been evaluated for the selected time step. If not (N), the method returns to  134  to select a next track. 
     Otherwise, the method advances to  142 , where a new value for the distribution parameters, χ′ k , are calculated from the expected values for the track states for that time step and the expected sufficient statistics, such that χ′ k ←E[χ k ]+E[L k ]. At  144 , the new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks from distribution parameters. Exact updates for the assignment matrix under the lower bound yields a product of CAP distributions: 
     
       
         
           
             
               
                 
                   
                     q 
                     ⁡ 
                     
                       ( 
                       
                         A 
                         
                           1 
                           : 
                           K 
                         
                       
                       ) 
                     
                   
                   = 
                   
                     
                       ∏ 
                       
                         k 
                         = 
                         1 
                       
                       K 
                     
                     ⁢ 
                     
                       CAP 
                       ( 
                       
                         
                           A 
                           k 
                         
                         ⁢ 
                         
                           
                              
                             
                               
                                 
                                   𝔼 
                                   
                                     q 
                                     ⁡ 
                                     
                                       ( 
                                       
                                         X 
                                         k 
                                       
                                       ) 
                                     
                                   
                                 
                                 ⁡ 
                                 
                                   [ 
                                   
                                     L 
                                     k 
                                   
                                   ] 
                                 
                               
                               + 
                               
                                 
                                   𝔼 
                                   
                                     q 
                                     ⁡ 
                                     
                                       ( 
                                       
                                         S 
                                         k 
                                       
                                       ) 
                                     
                                   
                                 
                                 ⁡ 
                                 
                                   [ 
                                   
                                     𝓍 
                                     k 
                                   
                                   ] 
                                 
                               
                             
                             ) 
                           
                           . 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   20 
                   ) 
                 
               
             
           
         
       
     
     To obtain a tractable algorithm, a loopy belief propagation algorithm can be used to compute    q(A     k     ) [A k ], for use in (16) and (19), above. The CAP distribution from (8) can be represented as a factor graph: 
                   CAP   (         A   ⁢        𝓍   )       ∝       ∏     i   =   1       N   T       ⁢         f   i   R     ⁡     (       A   i     ⁢   ′     )       ⁢       ∏     j   =   1       N   Z       ⁢         f   j   C     ⁡     (     A   ⁢           ⁢     ′   j       )       ⁢       ∑     i   =   0       N   T       ⁢       ∑     j   =   0       N   Z       ⁢       f   ij   S     ⁡     (     A   ij     )                     ,             (   21   )               
with f i   R (v):= {Σ i=0   N     x    v j =1}(R for row factors), f i   C (v):= {Σ i=0   N     T    v i =1} (C for column vectors), and f ij   S (v):=exp(χ ij v). Reparameterization methods are applied to convert this factor graph to a pairwise factor graph, where derivation of a Bethe free energy is facilitated. The Bethe free entropy is:
 
                       H   β     ⁡     [     q   ⁡     (   A   )       ]       :=         ∑     i   =   1       N   T       ⁢       ∑     j   =   0       N   Z       ⁢     H   ⁡     [     q   ⁡     (       r   i     ,     A   ij       )       ]           +       ∑     j   =   1       N   Z       ⁢       ∑     i   =   0       N   T       ⁢     H   ⁡     [     q   ⁡     (       c   j     ,     A   ij       )       ]           -       ∑     i   =   1       N   T       ⁢       N   Z     ⁢     H   ⁡     [     q   ⁡     (     r   i     )       ]           -       ∑     j   =   1       N   Z       ⁢       N   T     ⁢     H   ⁡     [     q   ⁡     (     c   j     )       ]           -       ∑     i   =   1       N   T       ⁢       ∑     j   =   1       N   Z       ⁢     H   ⁡     [     q   ⁡     (     A   ij     )       ]                     (   22   )                       ⁢       =         ∑     i   =   1       N   T       ⁢     H   ⁡     [     q   ⁡     (       A   i     ⁢   ′     )       ]         +       ∑     j   =   1       N   Z       ⁢     H   ⁡     [     q   ⁡     (     A   ⁢           ⁢     ′   j       )       ]         -       ∑     i   =   1       N   T       ⁢       ∑     j   =   1       N   Z       ⁢     H   ⁡     [     q   ⁡     (     A   ij     )       ]               ,             (   23   )               
where the pairwise conversion uses constrained auxiliary variables r i :=A i  and c j :=A j  and the implied relations H[q(r i , A ij )]=H[q(r i )]+H[q(A ij |r i )]=H[q(r i )]=H[q(A i .)].
 
     The method  130  uses an altered variational lower bound,    β (q), which replaces the entropy H[q(A k )] with H β [q(A k )] 2 . It will be appreciated that    β (q)   (q) with respect to q(X i:K , S 1:K ) which implies that the state posterior updates under the old lower bound in (16) and (19) remain unchanged under the new lower bound. To get the update equations for q(A k ), the new lower bound is expressed in terms of q(A 1:K ), such that: 
     
       
         
           
             
               
                 
                   
                     
                       ℒ 
                       β 
                     
                     ⁡ 
                     
                       ( 
                       q 
                       ) 
                     
                   
                   ⁢ 
                   
                     = 
                     c 
                   
                   ⁢ 
                     
                   ⁢ 
                   
                     
                       
                         𝔼 
                         q 
                       
                       [ 
                       
                         
                           log 
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           
                             p 
                             ( 
                             
                               
                                 Z 
                                 
                                   1 
                                   : 
                                   K 
                                 
                               
                               ⁢ 
                               
                                  
                                 
                                   
                                     X 
                                     
                                       1 
                                       : 
                                       K 
                                     
                                   
                                   , 
                                   
                                     A 
                                     
                                       1 
                                       : 
                                       K 
                                     
                                   
                                 
                                 ) 
                               
                             
                             ] 
                           
                         
                         + 
                         
                           
                             𝔼 
                             q 
                           
                           [ 
                           
                             log 
                             ⁢ 
                             
                                 
                             
                             ⁢ 
                             
                               P 
                               ( 
                               
                                 A 
                                 
                                   1 
                                   : 
                                   K 
                                 
                               
                                
                             
                             ⁢ 
                             
                               S 
                               
                                 1 
                                 : 
                                 K 
                               
                             
                           
                           ) 
                         
                       
                       ] 
                     
                     + 
                     
                         
                     
                     ⁢ 
                     
                       ( 
                       24 
                       ) 
                     
                   
                 
               
             
             
               
                 
                     
                   ⁢ 
                   
                     
                       ∑ 
                       
                         k 
                         = 
                         1 
                       
                       K 
                     
                     ⁢ 
                     
                       
                         H 
                         β 
                       
                       ⁡ 
                       
                         [ 
                         
                           q 
                           ⁢ 
                           
                             ( 
                             
                               A 
                               k 
                             
                             ) 
                           
                         
                         ] 
                       
                     
                   
                 
               
             
             
               
                 
                   
                     
                       
                         = 
                         c 
                       
                       ⁢ 
                         
                       ⁢ 
                       
                         
                           ∑ 
                           
                             k 
                             = 
                             1 
                           
                           K 
                         
                         ⁢ 
                         
                           
                             
                               𝔼 
                               
                                 q 
                                 ⁡ 
                                 
                                   ( 
                                   
                                     A 
                                     k 
                                   
                                   ) 
                                 
                               
                             
                             [ 
                             
                               
                                 1 
                                 T 
                               
                               ⁢ 
                               
                                 ( 
                                 
                                   
                                     
                                       A 
                                       k 
                                     
                                     ⊙ 
                                     
                                       
                                         𝔼 
                                         
                                           q 
                                           ⁡ 
                                           
                                             ( 
                                             
                                               X 
                                               k 
                                             
                                             ) 
                                           
                                         
                                       
                                       ⁡ 
                                       
                                         [ 
                                         
                                           L 
                                           k 
                                         
                                         ] 
                                       
                                     
                                   
                                   + 
                                   
                                     
                                       𝔼 
                                       
                                         q 
                                         ⁡ 
                                         
                                           ( 
                                           
                                             S 
                                             k 
                                           
                                           ) 
                                         
                                       
                                     
                                     ⁡ 
                                     
                                       [ 
                                       
                                         𝓍 
                                         k 
                                       
                                       ] 
                                     
                                   
                                 
                                 ) 
                               
                             
                             ) 
                           
                           ⁢ 
                           1 
                         
                       
                     
                     ] 
                   
                   + 
                   
                       
                   
                   ⁢ 
                   
                     ( 
                     25 
                     ) 
                   
                 
               
             
             
               
                 
                     
                   ⁢ 
                   
                     
                       ∑ 
                       
                         k 
                         = 
                         1 
                       
                       K 
                     
                     ⁢ 
                     
                       
                         H 
                         β 
                       
                       ⁡ 
                       
                         [ 
                         
                           q 
                           ⁢ 
                           
                             ( 
                             
                               A 
                               k 
                             
                             ) 
                           
                         
                         ] 
                       
                     
                   
                 
               
             
             
               
                 
                   
                     = 
                     c 
                   
                   ⁢ 
                     
                   ⁢ 
                   
                     
                       ∑ 
                       
                         k 
                         = 
                         1 
                       
                       K 
                     
                     ⁢ 
                     
                       
                         𝔼 
                         
                           q 
                           ⁡ 
                           
                             ( 
                             
                               A 
                               k 
                             
                             ) 
                           
                         
                       
                       [ 
                       
                         
                           log 
                           ⁢ 
                           
                               
                           
                           ⁢ 
                           
                             CAP 
                             ( 
                             
                               
                                 A 
                                 k 
                               
                               ⁢ 
                               
                                  
                                 
                                   
                                     
                                       𝔼 
                                       
                                         q 
                                         ⁡ 
                                         
                                           ( 
                                           
                                             X 
                                             k 
                                           
                                           ) 
                                         
                                       
                                     
                                     ⁡ 
                                     
                                       [ 
                                       
                                         L 
                                         k 
                                       
                                       ] 
                                     
                                   
                                   + 
                                   
                                     
                                       𝔼 
                                       
                                         q 
                                         ⁡ 
                                         
                                           ( 
                                           
                                             S 
                                             k 
                                           
                                           ) 
                                         
                                       
                                     
                                     ⁡ 
                                     
                                       [ 
                                       
                                         𝓍 
                                         k 
                                       
                                       ] 
                                     
                                   
                                 
                                 ) 
                               
                             
                             ] 
                           
                         
                         + 
                         
                           
                             
                               H 
                               β 
                             
                             ⁡ 
                             
                               [ 
                               
                                 q 
                                 ⁡ 
                                 
                                   ( 
                                   
                                     A 
                                     k 
                                   
                                   ) 
                                 
                               
                               ] 
                             
                           
                           . 
                           
                               
                           
                           ⁢ 
                           
                             ( 
                             26 
                             ) 
                           
                         
                       
                     
                   
                 
               
             
           
         
       
     
     This corresponds to the Bethe free energy of the factor graph described in (21), with  [L k ]+ [χ k ] as the CAP parameter. Therefore, the estimate of the posterior probability distribution for the track assignment can be determined using loopy belief propagation. 
       FIG. 5  illustrates a method  150  for calculating an estimate of the posterior probability distribution for the track assignment using a loopy belief propagation process. At  152 , the distribution parameters, χ k , at for each time step, are conditioned, with the columns normalized in log scale and values below a threshold value clipped at the threshold value (e.g., −25). At  154 , initial row and column quantities, μ and v, are initialized to random values to break symmetry. Specifically, the row and column quantities are defined as; 
                       μ   ij   R     :=     msg       f   i   R     →     A   ij           ,       μ   ij   C     :=     msg       f   i   C     →     A   ij           ,       v   ij   R     :=     msg       A   ij     →     f   i   R           ,     
     ⁢       v   ij   C     :=         msg       A   ij     →     f   i   C         ⁢           ⁢       v   ij   R     ⁡     (   x   )         =         μ   ij   C     ⁡     (   x   )       ⁢       f   ij   S     ⁡     (   x   )             ,         μ   ij   R     ⁡     (   1   )       =       ∏     k   ≠   j       ⁢       v   ik   R     ⁡     (   0   )           ,     
     ⁢         μ   ij   R     ⁡     (   0   )       =       ∑     l   ≠   j       ⁢         v   il   R     ⁡     (   1   )       ⁢       ∏       k   ≠   j     ,   l       ⁢       v   ik   R     ⁡     (   0   )               ,           (   27   )               
where all of the messages form functions in {0,1}→   + . This formulation exploits the fact that there is only one non-zero value in the row A i .
 
     At  156 , the row and column values are updated using the distribution parameters, χ k . It can be noted from (27) that: 
     
       
         
           
             
               
                 
                   
                     
                       
                         μ 
                         ij 
                         R 
                       
                       ⁡ 
                       
                         ( 
                         1 
                         ) 
                       
                     
                     = 
                     
                       
                         
                           
                             ∏ 
                             
                               k 
                               = 
                               0 
                             
                           
                           
                             N 
                             Z 
                           
                         
                         ⁢ 
                         
                           
                             
                               v 
                               ik 
                               R 
                             
                             ⁡ 
                             
                               ( 
                               0 
                               ) 
                             
                           
                           / 
                           
                             
                               
                                 v 
                                 ij 
                                 R 
                               
                               ⁡ 
                               
                                 ( 
                                 0 
                                 ) 
                               
                             
                             ⟹ 
                             
                               
                                 μ 
                                 ~ 
                               
                               ij 
                               R 
                             
                           
                         
                       
                       := 
                       
                         
                           
                             
                               μ 
                               ij 
                               R 
                             
                             ⁡ 
                             
                               ( 
                               0 
                               ) 
                             
                           
                           
                             
                               μ 
                               ij 
                               R 
                             
                             ⁡ 
                             
                               ( 
                               1 
                               ) 
                             
                           
                         
                         = 
                         
                           
                             
                               
                                 ∑ 
                                 
                                   l 
                                   = 
                                   0 
                                 
                                 
                                   N 
                                   Z 
                                 
                               
                               ⁢ 
                               
                                 
                                   
                                     v 
                                     il 
                                     R 
                                   
                                   ⁡ 
                                   
                                     ( 
                                     1 
                                     ) 
                                   
                                 
                                 
                                   
                                     v 
                                     il 
                                     R 
                                   
                                   ⁡ 
                                   
                                     ( 
                                     0 
                                     ) 
                                   
                                 
                               
                             
                             - 
                             
                               
                                 
                                   v 
                                   ij 
                                   R 
                                 
                                 ⁡ 
                                 
                                   ( 
                                   1 
                                   ) 
                                 
                               
                               
                                 
                                   v 
                                   ij 
                                   R 
                                 
                                 ⁡ 
                                 
                                   ( 
                                   0 
                                   ) 
                                 
                               
                             
                           
                           ∈ 
                           
                             ℝ 
                             + 
                           
                         
                       
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   28 
                   ) 
                 
               
             
           
         
       
     
     The ratio of the messages to the row factors can then be expressed as:
 
{tilde over ( v )} ij   R   :=v   ij   R (1)/ v   ij   R (0)=(μ ij   C (1)/μ ij   C (0))exp(χ ij )∈   +   (29)
 
     The relations in (27)-(29) can be applied to the column messages to acquire a message passing update scheme in terms of message ratios: 
     
       
         
           
             
               
                 
                   
                     
                       
                         
                           μ 
                           ~ 
                         
                         ij 
                         R 
                       
                       ⁡ 
                       
                         ( 
                         1 
                         ) 
                       
                     
                     = 
                     
                       
                         
                           ∑ 
                           
                             l 
                             = 
                             0 
                           
                           
                             N 
                             Z 
                           
                         
                         ⁢ 
                         
                           
                             v 
                             ~ 
                           
                           il 
                           R 
                         
                       
                       - 
                       
                         
                           v 
                           ~ 
                         
                         ij 
                         R 
                       
                     
                   
                   , 
                   
                     
                       
                         v 
                         ~ 
                       
                       ij 
                       R 
                     
                     = 
                     
                       
                         exp 
                         ⁡ 
                         
                           ( 
                           
                             𝓍 
                             ij 
                           
                           ) 
                         
                       
                       
                         
                           μ 
                           ~ 
                         
                         ij 
                         C 
                       
                     
                   
                   , 
                   
                     
                       
                         μ 
                         ~ 
                       
                       ij 
                       C 
                     
                     = 
                     
                       
                         
                           ∑ 
                           
                             l 
                             = 
                             0 
                           
                           
                             N 
                             T 
                           
                         
                         ⁢ 
                         
                           
                             v 
                             ~ 
                           
                           il 
                           C 
                         
                       
                       - 
                       
                         
                           v 
                           ~ 
                         
                         ij 
                         C 
                       
                     
                   
                   , 
                   
                     
 
                   
                   ⁢ 
                   
                     
                       
                         v 
                         ~ 
                       
                       ij 
                       C 
                     
                     = 
                     
                       
                         
                           exp 
                           ⁡ 
                           
                             ( 
                             
                               𝓍 
                               ij 
                             
                             ) 
                           
                         
                         
                           
                             μ 
                             ~ 
                           
                           ij 
                           R 
                         
                       
                       . 
                     
                   
                 
               
               
                 
                   ( 
                   30 
                   ) 
                 
               
             
           
         
       
     
     Once the row and column values from the loopy belief propagation have been updated, the lower Bethe free energy is calculated at  158 , for example, via the relation shown in (26) above. At  160 , it is determined if the calculated free energy is below a threshold value. If not (N), the method  150  returns to  156  to update the row and column values. If so (Y), at  162  the marginal distributions for the track assignment are calculated by normalizing the product of the incoming messages to each variable, such that:  [A ij ]=P(A ij =1)=σ(χ ij −log {tilde over (μ)} ij   R −log {tilde over (μ)} ij   C ). 
     Returning to  FIG. 4 , once the loopy belief propagation algorithm has finished for a given time step, it is determined at  146  if marginal distributions have been updated for all time steps. If not (N), the method returns to  132  to select a next time step. Otherwise (Y), the method terminates, returning the posterior probability distribution for the track assignment. 
     Returning to  FIG. 2 , once the posterior distributions for the track states and the track assignment have been determined, the method advances to  160 , where a variational lower bound is calculated from the new estimate of the posterior probability distribution for the plurality of possible assignments of the set of measurements to the set of tracks, the estimate of a posterior probability distribution for a plurality of track states, and the set of measurements. The variational lower bound represents a lower bound for the marginal probability of the set of measurements given a model defined by the new estimate of the posterior probability distribution and the estimate of a posterior probability distribution for a plurality of track states. For example, the variational lower bound can be calculated using the relation in (12). 
     At  170 , it is determined if the calculated variational lower bound is less than a threshold value. If not (N), the method returns to  110  to update the estimate of the posterior distribution of the track states using the new estimate for the posterior distribution of the track assignment. If so (Y), the method terminates and returns the current estimates for the track state and track assignment distributions. 
       FIG. 6  is a schematic block diagram illustrating an exemplary system  200  of hardware components capable of implementing examples of the systems and methods disclosed in  FIGS. 1-5 . The system  200  can include various systems and subsystems. The system  200  can be a personal computer, a laptop computer, a workstation, a computer system, an appliance, an application-specific integrated circuit (ASIC), a server, a server blade center, a server farm, etc. 
     The system  200  can includes a system bus  202 , a processing unit  204 , a system memory  206 , memory devices  208  and  210 , a communication interface  212  (e.g., a network interface), a communication link  214 , a display  216  (e.g., a video screen), and an input device  218  (e.g., a keyboard and/or a mouse). The system bus  202  can be in communication with the processing unit  204  and the system memory  206 . The additional memory devices  208  and  210 , such as a hard disk drive, server, stand-alone database, or other non-volatile memory, can also be in communication with the system bus  202 . The system bus  202  interconnects the processing unit  204 , the memory devices  206 - 210 , the communication interface  212 , the display  216 , and the input device  218 . In some examples, the system bus  202  also interconnects an additional port (not shown), such as a universal serial bus (USB) port. 
     The processing unit  204  can be a computing device and can include an application-specific integrated circuit (ASIC). The processing unit  204  executes a set of instructions to implement the operations of examples disclosed herein. The processing unit can include a processing core. 
     The additional memory devices  206 ,  208  and  210  can store data, programs, instructions, database queries in text or compiled form, and any other information that can be needed to operate a computer. The memories  206 ,  208  and  210  can be implemented as computer-readable media (integrated or removable) such as a memory card, disk drive, compact disk (CD), or server accessible over a network. In certain examples, the memories  206 ,  208  and  210  can comprise text, images, video, and/or audio, portions of which can be available in formats comprehensible to human beings. 
     Additionally or alternatively, the system  200  can access an external data source or query source through the communication interface  212 , which can communicate with the system bus  202  and the communication link  214 . 
     In operation, the system  200  can be used to implement one or more parts of n object tracking system using variational tracking. Computer executable logic for implementing the system resides on one or more of the system memory  206 , and the memory devices  208 ,  210  in accordance with certain examples. The processing unit  204  executes one or more computer executable instructions originating from the system memory  206  and the memory devices  208  and  210 . The term “computer readable medium” as used herein refers to a medium that participates in providing instructions to the processing unit  204  for execution, and can include either a single medium or multiple non-transitory media operatively connected to the processing unit  204 . 
     The invention has been disclosed illustratively. Accordingly, the terminology employed throughout the disclosure should be read in an exemplary rather than a limiting manner. Although minor modifications of the invention will occur to those well versed in the art, it shall be understood that what is intended to be circumscribed within the scope of the patent warranted hereon are all such embodiments that reasonably fall within the scope of the advancement to the art hereby contributed, and that that scope shall not be restricted, except in light of the appended claims and their equivalents.