Joint Multi-User Tracking and Signal Detection in Reconfigurable Intelligent Surface-Assisted Cell-Free ISAC Systems

1 Join t Multi-User T rac king and Signal Detection in Recongurable In telligen t Surface-Assisted Cell-F ree ISA C Systems W eifeng Zh u, Jun yuan Gao, Sh uo w en Zhang, Meixia T ao, and Liang Liu Abstract—This pap er inv estigates the cell-free multi-user in tegrated sensing and communication (ISAC) system, where m ultiple base stations collaboratively track the users and de- tect their signals. Moreov er, recongurable intelligen t surfaces (RISs) are deplo yed to serv e as additional reference nodes to o vercome the line-of-sigh t blockage issue of mobile users for accomplishing seamless sensing. Due to the high-speed user mobilit y , the m ulti-user tracking and signal detection p erfor- mance can b e signicantly deteriorated without elab orated online user kinematic state updating principles. T o tackle this c hallenge, we rst manage to establish a probabilistic signal mo del to comprehensively c haracterize the in terdependencies among user states, transmit signals, and received signals during the tracking pro cedure. Based on the Ba y esian problem formu- lation, we further propose a nov el hybrid v ariational message passing (HVMP) algorithm to realize computationally ecient join t estimation of user states and transmit signals in an online manner, which integrates VMP and standard MP to derive the p osterior probabilities of estimated v ariables. F urthermore, the Ba yesian Cramér-Rao b ound is pro vided to c haracterize the p erformance limit of the multi-user tracking problem, which is also utilized to optimize RIS phase proles for tracking p erfor- mance enhancemen t. Numerical results demonstrate that the prop osed algorithm can signican tly impro v e b oth trac king and signal detection performance o ver the representativ e Bay esian estimation counterparts. Index T erms—Integrated sensing and comm unication (ISA C), recongurable in telligent surface (RIS), multi-user trac king, signal detection, h ybrid v ariational message passing. I. In tro duction The adven t of the six-generation (6G) wireless netw orks marks a signicant revolution to encompass adv anced en- vironmen tal sensing with information transmission under the common cellular infrastructure through integrated sensing and communication (ISAC) [1]–[3]. Under the dual-functional integration, a v ariet y of emerging intel- ligen t applications including autonomous driving, low- altitude economy , and smart man ufacturing [2]–[4] can b e supported by 6G cellular systems without em b edding extra infrastructures. As the basic issue in ISAC, the acquisition of accurate p osition information is usually realized by measuring the line-of-sight (LOS) c hannel W eifeng Zhu, Jun yuan Gao, Shuo wen Zhang, and Liang Liu are with Departmen t of Electrical and Electronic Engineering, The Hong Kong Polytec hnic Universit y , Hong Kong SAR, China (email: {eee- wf.zhu, jun yuan.gao, shuo wen.zhang, liang-eie.liu}@p olyu.edu.hk). Meixia T ao is with School of Information Science and Electronic Engineering, Shanghai Jiao T ong Universit y , Shanghai 200240, China (email: mxtao@sjtu.edu.cn). parameters tow ards the BS [5]–[7], which includes the dela y , Doppler frequency , angle of arriv al (AO A), etc. Ho w ev er, the LOS path to wards the BS can often b e blo c k ed by v arious obstacles in practical scenarios, leading to missing p osition measurements. T o ov ercome the reliance on the LOS c hannel exis- tence, recongurable intelligen t surface (RIS) emerges as a promising technique in ISAC systems by reshaping electromagnetic propagation en vironmen ts [8]–[11]. In par- ticular, RISs can b e widely deploy ed to provide seamless surv eillance cov erage by creating virtual LOS channels along the user-RIS-BS paths, which bypass environmen t obstacles via con trollable reections. By regarding RISs as additional anchors, our preliminary w ork [12] rigourously sho w the feasibilit y to acquire the position information from these virtual LOS channels, when the LOS blo ck- age o ccurs. On the other hand, RISs are also capable of mitigating the wireless interference and fading issue under adv anced phase prole design, which signican tly enhances the communication quality [8], [9]. Thanks to these attractive b enets, the RIS b ecomes one of the indisp ensable basic hardw are facilities in 6G ISAC systems [13]. In the context of RIS-assisted ISAC, early researches mainly fo cus on the p ersp ectives of theoretical per- formance analysis and optimization, where the sensing p erformance is c haracterized b y sensing signal-to-noise ratio (SNR) [14] and Cramér-Rao b ound (CRB) [15]– [17]. F rom the p ersp ective of algorithm design, the work [18] proposes a ngerprin t-based metho d to lo calize indoor users based on their unique received signal strength (RSS) patterns. How ev er, the high-resolution RSS information is hard to be obtained in the outdoor scenario, which restricts the applications. Alternativ ely , the w ork [19] prop oses a trilateration-based metho d b y utilizing the range information from b oth direct and reecting links, while the angular measurement capabilit y of RISs is not exca v ated with limited lo calization resolution. Then, the w orks [20], [21] leverage the maximum likelihoo d methods to estimate b oth of the distance and A OA to w ards the RIS for lo calization, while their complexities are prohibitively in tensiv e in the multi-user scenario. In our work [12], an ecien t subspace-based algorithm and Discrete F ourier T ransform (DFT)-based RIS phase prole are designed to facilitate the multi-user lo calization problem. How ev er, these works perform lo calization in eac h time slot in- 2 dep enden tly , which neglects the temp oral correlation of user kinematic states in adjacent frames. As such, in the mobilit y case, these metho ds may only pro vide limited estimation p erformance on p osition information. T o promote high-accuracy p osition tracking of mobile users, it is b enecial to leverage b oth of the historical and fresh observ ations to acquire the p osition information in eac h time slot. F ollowing this guideline, recent w orks [22], [23] prop ose to apply the famous Kalman lter (KF) tec hnique [24] to improv e the trac king p erformance. In [22], an extended KF method is prop osed to track the mobile user based on the range measurements tow ards the BS and RIS. Then, the work [23] prop oses a manifold- based unscented KF metho d for user tracking based on dela y-Doppler-angle measurements. Ho wev er, these tw o w orks are limited to the single-user case and assume that the LOS path to the BS alwa ys exists. Alternatively , the w ork [25] prop oses a Bay eisan tracking algorithm for the m ulti-user tracking problem in the pure LOS blockage scenario based on only the p osition-related measurements to w ards the RIS. In practice, the blo ckage conditions of direct and reecting links are usually random and unkno wn, implying that the received signals contain a com bination of measuremen ts from these tw o t ypes of links. In contrast to the aforemen tioned works, this w ork aims to provide a comprehensive in vestigation on the general RIS-assisted cell-free ISAC system with unknown link blo c kages, which has the following three dierences. First, it is necessary to iden tify the blo ckage conditions of b oth direct and reecting links in each time slot, which then enables us to extract and combine the kinematic mea- suremen ts in receiv ed signals to accomplish the tracking task. Second, compared with single-BS systems [22], [23], [25], there exist strong correlations among the reecting links to all BSs through the identical RIS in cell-free ISA C systems, which share the same position-related parameters to w ards the RIS. Thus, this work prop oses to directly estimate the position-related parameters tow ards eac h RIS from all the related reecting links instead of treating eac h reecting link individually . Third, note that the ab ov e w orks [12], [18]–[23], [ 25] all apply the pilot signals to the sensing task, leading to signicantly increased sensing o v erhead in the multi-user scenario. T o address this issue, w e are motiv ated to leverage the data signals to facilitate the sensing task for b o osting b oth of the sensing and comm unication p erformance. In this work, we studies the joint multi-user tracking and signal detection problem in the RIS-assisted cell-free ISA C system. The main distinctions and contributions of this work are summarized as follows: 1) F or the considered system, w e rst prop ose an ISA C uplink transmission proto col for realizing join t m ulti-user trac king and communication. Then, we establish a nov el probabilistic signal mo del that includes geometric constraints and statistical dep en- dencies in the received signals for the considered system. Subsequen tly , we formulate the joint m ulti- user trac king and signal detection problem under the Ba y esian inference framework, enabling to realize either maxim um a p osterior (MAP) or minimum mean square error (MMSE) estimation. 2) Under the graphical representation of the prob- abilistic model, we devise a no v el h ybrid v aria- tional message passing (HVMP) algorithm for the join t estimation of user states and transmit signals. Sp ecically , the HVMP algorithm integrates VMP and MP to realize ecien t calculations of the desired p osterior probabilities and side information for fu- ture tracking, b oth of which are given in explicit expressions. Compared with existing works, the blo c kage condition of direct and reecting links can also b e automatically estimated, and the p osition- related measurements of all av ailable links can b e deriv ed and exibly combined to facilitate the joint estimation. 3) W e provide the Bay esian Cramér-Rao b ound (BCRB) to characterize the performance limit of the considered problem. Meanwhile, we formulate a w eigh ted BCRB minimization problem to optimize the RIS phase prole, whic h enhances the track- ing and detection p erformance ov er the traditional heuristic RIS phase prole design. 4) Extensive n umerical results are provided to v ali- date the p erformance sup eriority of the proposed approac h to represen tativ e b enchmarks. It is veried that the HVMP can pro vide consisten tly centimeter- lev el tracking accuracy even if there are a large n um b er of tracking users. Moreov er, the millimeter- lev el tracking accuracy can p oten tially b e accom- plished b y eectiv ely combining measurements from b oth direct and reecting links when the blo ckage condition b ecomes less serious. The rest of the pap er is organized as follo ws. Sec- tion I I in tro duces the system model of the considered RIS-assisted ISAC system and Section I I I pro vides the Ba y esian problem formulation of the considered prob- lem. In Section IV, the prop osed HVMP algorithm is dev elop ed, whose message passing equations are sp ecied. In Section V, the BCRB is obtained for p erformance analysis and a RIS phase prole optimization problem is studied. Finally , Section VI ev aluates the p erformance of the proposed m ulti-user tracking and signal detection algorithm and Section VI I concludes this work. Notations: In this pap er, vectors and matrices are denoted b y b oldface lo w er-case letters and boldface upp er- case letters, resp ectively . F or a complex matrix X , [ X ] N , M denote the sub-matrix with ro ws and columns in N and M , resp ectively , while X T and X H denote its transp ose, and conjugate transp ose, resp ectively . F or a complex n um b er x ∈ C , its phase is denoted as ∠ x . F urther, E a [ · ] denotes the expectation operation o v er random v ector a , ⊙ denotes the Hadamard pro duct, and ⊗ denotes the Kroneck er pro duct. Gaussian distributions o v er the real random v ector x ∈ R N × 1 and complex 3 RIS  Us er         RIS Control le r Obs ta cl es        Us er                 RIS                  BS                               BS                                Fig. 1. Illustration of the system model. Multiple RISs are deploy ed in the system to assist the sensing and communication service of users. random vector y ∈ C N × 1 are denoted by N ( x ; ¯ x , Σ Σ Σ x ) and C N ( y ; ¯ y , Σ Σ Σ y ) , resp ectively . The von Misés (VM) distri- bution ov er the v ariable θ is denoted by V M ( θ ; µ, κ ) = (2 π I 0 ( κ )) − 1 exp( κ cos( θ − µ )) , with I 0 ( κ ) denoting the zeroth-order Bessel function. I I. System Mo del W e consider an uplink cell-free ISAC system with G BSs eac h with M B an tennas, K single-an tenna mobile users, and R RISs each with M I reecting elements, as sho wn in Fig. 1. In the tw o-dimension (2D) system, the co ordinate of each BS g ∈ G = { 1 , 2 , . . . , G } is denoted as p B g = [ p B x,g , p B y,g ] T ∈ R 2 × 1 and the co ordinate of each RIS r ∈ R = { 1 , 2 , . . . , R } is denoted as p I r = [ p I x,r , p I y,r ] T ∈ R 2 × 1 . In each time slot t , the p osition and velocity of eac h user k are denoted as p U k,t = [ p U x,k,t , p U y,k ,t ] T ∈ R 2 × 1 and ˙ p U k,t = [ ˙ p U x,k,t , ˙ p U y,k ,t ] T ∈ R 2 × 1 , resp ectively . All the users keep transmitting orthogonal frequency division m ultiplexing (OFDM) signals in the uplink for V time slots, while the BSs need to dynamically track the users and detect their transmit signals simultaneously . A. Mobilit y Mo del W e assume that the motion of each user follows the discrete-time state model [26]. Dene the state vector for all users at each time slot t as ψ ψ ψ t = [ ψ ψ ψ T 1 ,t , . . . , ψ ψ ψ T K,t ] T ∈ R 4 K × 1 with ψ ψ ψ k,t = [( p U k,t ) T , ( ˙ p U k,t ) T ] T ∈ R 4 × 1 . As such, the state transition equation [26]–[28] can b e written by 1 ψ ψ ψ t = F ψ ψ ψ t − 1 + ω ω ω t − 1 , ∀ t, (1) where the transition matrix F is dened as F = I K ⊗ F 0 with F 0 =  I 2 ∆ T I 2 0 I 2  ∈ R 4 K × 4 K . (2) In the mobilit y model (1), ω ω ω t − 1 ∈ R 4 K × 1 denotes the mo- tion noise follo wing the Gaussian distribution N ( 0 , Q t − 1 ) with Q t − 1 = blkdiag ( Q 1 ,t − 1 , . . . , Q K,t − 1 ) , and ∆ T denote the time in terv al betw een t w o adjacen t trac king time slots. 1 In practice, the velocity of each user k usually cannot keep perfectly constant in the whole time slot t . T o address this issue, the velocity uncertain ty is modeled as the motion noise term in the discrete-time state model. It is assumed that the motion noise vectors w t and w t ′ are mutually indep endent for t  = t ′ . F rom the ab o v e mobilit y mo del, w e can nd that the historical user state measuremen ts can b e leveraged to help rene the m ulti- user trac king performance, which then enhances the signal detection p erformance as well. B. Signal Mo del In the ab ov e system, each user transmits Q OFDM sym b ols at each time slot, while the num b er of sub- carriers is N and the sub-carrier spacing is ∆ f Hz, suc h that the bandwidth is B = N ∆ f Hz. Thus, there are totally QN resource elements in the uplink trans- mission phase, contributing to the resource element set E = { ( q , n ) | q = 1 , . . . , Q, n = 1 , . . . , N } . T o av oid the in ter-sym b ol interference, a cyclic prex (CP) comprised b y J OFDM samples is inserted at the b eginning of each transmit OFDM sym b ol. As suc h, the p erio d of eac h OFDM symbol is ∆ t = N + J N ∆ f . A t eac h BS g , the received signal at the ( q , n ) -th resource element in each time slot t can b e given by y ( q ) g,n,t = K X k =1 s ( q ) k,n,t h ( q ) k,g ,n,t + z ( q ) g,n,t , ∀ g , n, q, t, (3) where s ( q ) k,n,t is the data signal at the ( q , n ) -th resource elemen t transmitted from the user k , h ( q ) k,g ,n,t ∈ C M B × 1 is the channel b etw een the user k and the BS g at the ( q , n ) - th resource elemen t, and z ( q ) g,n,t ∈ C M B × 1 is the bac kground noise v ector follo wing complex Gaussian distribution of C N ( z ( q ) g,n,t ; 0 , σ 2 z I ) . C. Channel Mo del In the signal mo del (3), each of the eective channel v ectors, i.e., { h ( q ) k,g ,n,t } , is comprised by the direct user-BS c hannel and cascaded user-RIS-BS channels. Sp ecically , for the ( q , n ) -th resource elemen t at time slot t , let h UB , ( q ) k,g ,n,t , h UI , ( q ) k,r,n,t , and H IB g,r,n denote the c hannels from user k to the BS g , from user k to the RIS r , and from the RIS r to the BS g , respectively . Then, the corresp onding c hannel v ector h ( q ) k,g ,n,t is denoted as h ( q ) k,g ,n,t = h UB , ( q ) k,g ,n,t + R X r =1 H IB g,r,n diag ( ϕ ϕ ϕ ( q ) r,t ) h UI , ( q ) k,r,n,t . (4) where ϕ ϕ ϕ ( q ) r,t = [ ϕ ( q ) 1 ,r,t , · · · , ϕ ( q ) M I ,r,t ] T ∈ C M I × 1 is the phase prole vector of the RIS r during the q -th OFDM symbol duration in the time slot t with | ϕ ( q ) m I ,r,t | = 1 , ∀ m I , r , t, q . In this w ork, w e assume that the user-BS channels, user-RIS c hannels, and RIS-BS channels all follow the binary LOS c hannel mo del 2 , under which each c hannel is either a LOS c hannel in the case without blo ckage or zero in the case with blo ckage. Sp ecically , the RIS r -BS g , user k -BS g , and user k -RIS r channels can b e given b y H IB g,r,n = α IB g,r β IB g,r e − ȷζ F τ IB g,r a B ( θ B g,r ) a T I ( θ I g,r ) , (5) 2 The small-scale fading component generally has muc h weak er power than the LOS path (smaller than -20 dB) [29], [30] and thus can b e eliminiated, due to severe path loss and reection loss. 4 h UB , ( q ) k,g ,n,t = α UB k,g ,t β UB k,g ,t e − ȷ ( ζ F ( n − 1) τ UB k,g ,t + ζ T ( q − 1) ν UB k,g ,t ) × a B ( θ UB k,g ,t ) , (6) h UI , ( q ) k,r,n,t = α UI k,r,t β UI k,r,t e − ȷ ( ζ F ( n − 1) τ UI k,r,t + ζ T ( q − 1) ν UI k,r,t ) × a I ( θ UI k,r,t ) , (7) resp ectiv ely , where ζ F = 2 π ∆ f and ζ T = 2 π ∆ t . In the ab o v e, α ∈ { 0 , 1 } , β , τ , ν , and θ are the blo c kage indicator, complex channel gain, dela y prole, A OA (angle of departure (AOD)), resp ectively , whose sup erscripts and subscripts are omitted for clarity . The geometric relationship b etw een the p osition-related parameters ( τ , ν , and θ ) and user states is detailed provided in Section I I I-B. Consider that the uniform linear array (ULA) with d spacing is emplo yed for b oth BSs and RISs, the steering v ectors a i ( θ ) , ∀ i ∈ { B , I } with the AO A (A OD) θ can b e expressed by a i ( θ ) = [1 , e − ȷζ S cos θ , . . . , e − ȷζ S ( M i − 1) cos θ ] T ∈ C M i × 1 , (8) where ζ S = 2 πd λ . Based on the mobility mo del (1), signal mo del (3), and the ab ov e channel mo del, the BSs aim to achiev e t w o goals: tracking the user p ositions, i.e., estimating { p U k,t } K k =1 , ∀ t , and decoding user messages, i.e., estimating s ( q ) k,n,t , ∀ k , n, q , t . Remark 1: Note that there are quite a lot of works considering lo calization and tracking problem in wireless net w orks [5], [6], [31]. Compared to these works, we hav e t w o new challenges. First, we consider the case when the LOS c hannels for user-BS and user-RIS links can b e blo ck ed, due to v arious obstacles. Thus, it is also essen tial to estimate the blo ckage conditions of all LOS c hannels during the tracking pro cedure, which enables to extract the user state information. Second, in [5], [6], [25], [31], eac h user transmits known pilot signals to the BS, whic h only needs to fo cus on the tracking job. On the con trary , in our considered system, the BS can utilize the data signal to p erform trac king. Therefore, the BSs need to estimate both the position-related information (angle, delay , and Doppler frequency) and the unknown user messages simultaneously from their received signals. The ab ov e tw o new challenges call for totally new signal pro cessing solutions in our setup. I I I. Problem F orm ulation In this section, w e rst introduce an ISA C protocol such that BSs can simultaneously trac k the users and detect their signals based on their receiv ed signals. Then, we form ulate the join t tracking and signal detection problem under this proto col. A. ISA C Proto col Under the prop osed proto col, all the resource elements are divided into t wo sets, i.e., Set I, whic h is termed as the ISA C set E I with E I = Q I N I resource elements, and Set I I, which is termed as the communication set E C = E \ E I with E C = ( QN − Q I N I ) resource elements. Sp ecically , all the BSs lev erage received signals in the ISAC Set to T im e d o m ain Frequ en cy d o m ain ISA C reso u rce el em en t Com m u n i cati o n reso u rce el em en t                                                 res o u rce el em en t T im e d o m ain     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜                            (a) T ransmission signal structure for uplink users. T im e d o m ain Frequ en cy d o m ain ISA C reso u rce el em en t Com m u n i cati o n reso u rce el em en t                                                 res o u rce el em en t T im e d o m ain     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜     󰇛  󰇜                            (b) Reection pattern of the phase prole for each RIS. Fig. 2. An example illustration of the ISAC proto col with Q 1 = 4 , I = 3 , ∆ Q = 9 , and ∆ N = 3 . Then, we hav e N I = { 1 , 4 , 7 } and Q I = { 1 , 2 , 3 , 4 } ∪ { 10 , 11 , 12 , 13 } ∪ { 19 , 20 , 21 , 22 } . join tly trac k the users and detect their signals at eac h time slot t . Then, since the user channels are kno wn with their estimated AO A/A OD, delay , and Doppler information, BSs only need to detect the signals in the communication set. In the rest of this pap er, we merely fo cus on the ISA C set, b ecause signal detection with known channels in the comm unication set is standard. In this w ork, we consider the ISAC set given by E I = { ( q , n ) | q ∈ Q 0 , n ∈ N 0 } , whic h has a com b-lik e pattern sho wn in Fig. 2(a) to facilitate precise estimation on both dela ys and Doppler frequencies. Sp ecically , w e hav e N 0 = { ς ˘ n | ς ˘ n = 1 + ( ˘ n − 1)∆ N , ˘ n = 1 , . . . , N I } and Q 0 = S I i =1 Q i , resp ectiv ely , where ∆ N ≥ 1 and Q i = { ι i, ˘ q | ι i, ˘ q = ˘ q + ( i − 1)∆ Q, ˘ q = 1 , . . . , Q 1 } , ∀ i with Q I = I Q 1 , Q 1 > 1 , and ∆ Q ≫ 1 . 3 Next, w e introduce the user transmission strategy in the ISA C set. Similar to [32], [33], we adopt repetition co ding at the user side. Sp ecically , we ha v e s ( q ) k,n,t = ˜ s k,t , ∀ k , t, ( q , n ) ∈ E I . (9) The reason to adopt repetition co ding in the ISA C set is as follo ws. F rom the model (3) – (7), we require to estimate E I indep enden t signals and one delay-Doppler pair from E I signal measurements with indep endent co ding, which is generally an underdetermined problem and intractable. Th us, rep etition co ding is motiv ated to address this issue b y reducing the num b er of estimated signals, contributing to an ov erdetermined problem. Note that rep etition coding ma y reduce the spectral eciency (SE). How ever, this SE loss is inevitable if BSs ha v e to jointly track the users and detect their signals. Moreov er, repetition co ding is merely used in the ISAC set, while we do not need rep etition coding an ymore for the communication set since the channels are known for signal detection. Moreo v er, for extracting the AO A information of user- RIS links, w e follow our preliminary w ork [12] to design the 3 The spacing ∆ Q is usually set to b e large enough to facilitate the Doppler frequency estimation. 5  3rd o rde r tensor                     4th orde r tensor                           T im e S pa ti a l Facto ri ze t h e t i m e d i m en s i o n i n t o t w o new d i m en s i o n s D i rect l i n k Refl ect i n g link Fig. 3. The tensor representation of the received signal. RIS phase prole with ϕ ϕ ϕ ( ι i, ˘ q ) r,t = ˜ ϕ ϕ ϕ ( ˘ q ) r,t , ∀ r, t, i, ˘ q = 1 , . . . , Q 1 and ˜ ϕ ϕ ϕ ( ˘ q 1 ) r,t  = ˜ ϕ ϕ ϕ ( ˘ q 2 ) r,t , ∀ ˘ q 1  = ˘ q 2 . As shown in Fig. 2(b), this design implies that the RIS has the common RIS phase pattern in each subset Q i , ∀ i . This also explains why we consider Q 1 > 1 for eac h Q i , since the AO A of eac h user- RIS link can never b e extracted with Q 1 = 1 [12]. Remark 2: Note that the proposed ISAC protocol diers from the traditional OFDM transmission proto col only in terms of transmit signals in the ISAC sub-blo ck. Sp ecically , w e replace the pilot signals in the traditional OFDM proto col with information signals to enhance the sp ectral eciency . This also means that the ISA C protocol can b e easily applied to existing communication systems. In Section VI, we will show that the prop osed ISA C proto col has the approaching tracking performance to the pilot-based transmission paradigm. F or the ISAC set, the received signals can b e expressed as a 3rd order tensor Y Y Y g,t ∈ C M B × Q I × N I . Under the ab o v e design, w e can further factorize the time dimension in to tw o new dimensions, res ulting in a 4th order tensor represen tation, as sho wn in Fig. 3. F or conv enience, we dene the vector form of the received signal as y g,t = v ec ( Y Y Y g,t ) = h ( y g,ς 1 ,t ) T , . . . , ( y g,ς N I ,t ) T i T , ∀ g , t, (10) where y g,ς ˘ n ,t =   y ( ι 1 , 1 ) g,ς ˘ n ,t  T , . . . ,  y ( ι 1 ,Q 1 ) g,ς ˘ n ,t  T , . . . ,  y ( ι I , 1 ) g,ς ˘ n ,t  T , . . . ,  y ( ι I ,Q 1 ) g,ς ˘ n ,t  T  T , ∀ g , ˜ q , t. (11) Sp ecically , the signal vector (10) can b e re-expressed by y g,t = K X k =1  w UB k,g ,t  a F ( τ UB k,g ,t ) ⊗ ˘ a T ( ν UB k,g ,t ) ⊗ a T ( ν UB k,g ,t ) ⊗ a B ( θ UB k,g ,t )  + R X r =1 w UI k,g ,r,t  a F ( τ UIB k,g ,r,t ) ⊗ ˘ a T ( ν UI k,r,t ) ⊗ ˘ a I ( θ I r , θ UI k,r,t , ν UI k,r,t ) ⊗ a B ( θ I g,r )   + z g,t = K X k =1  w UB k,g ,t a UB k,g ,t + R X r =1 w UI k,g ,r,t a UI k,g ,r,t  + z g,t = A g,t w g,t + z g,t , (12) where w UB k,g ,t = √ P ˜ s k,t α UB k,g ,t β UB k,g ,t and w UI k,g ,r,t = √ P ˜ s k,t α IB g,r α UI k,r,t β IB g,r β UI k,r,t are eective signals along the user k -BS g and user k -RIS r -BS g links, resp ectively; τ UIB k,g ,r,t = τ UI k,r,t + τ IB g,r is the total delay of the user k -RIS r -BS g link; a UB k,g ,t = a F 0 ( τ UB k,g ,t ) ⊗ a T 2 ( ν UB k,g ,t ) ⊗ a T 1 ( ν UB k,g ,t ) ⊗ a S ( θ UB k,g ,t ) is the eectiv e steering vector for the user k - BS g c hannel and a UI k,g ,r,t is similarly dened; A g,t = [ A 1 ,g,t , . . . , A K,g ,t ] and w g,t = [( w 1 ,g,t ) T , . . . , ( w K,g ,t ) T ] T are the eective steering matrix and eective signal vector, resp ectiv ely , with A k,g ,t = [ a UB k,g ,t , a UI k,g , 1 ,t . . . , a UI k,g ,R,t ] and w k,g ,t = [ w UB k,g ,t , w UI k,g , 1 ,t . . . , w UI k,g ,R,t ] T . Then, the aggre- gated signal matrix from all BSs can b e expressed b y Y t = [ y 1 ,t , . . . , y G,t ] . In (12), ˘ a I ( · ) , a F ( · ) , and a T ( · ) / ˘ a T ( · ) are the RIS-related, delay-related, and Doppler-related steering vectors dened by 4 ˘ a I ( φ, θ, ν ) = ( Ψ Ψ Ψ T r,t ( a I ( φ ) ⊙ a I ( θ ))) ⊙ a T ( ν ) , = ( Ψ Ψ Ψ T r,t ˜ a I ( φ, θ ))) ⊙ a T ( ν ) , (13) a F ( τ ) =  1 , e − ȷ ˘ ζ F τ , . . . , e − ȷ ˘ ζ F ( N 0 − 1) τ  T , (14) a T ( ν ) =  1 , e − ȷζ T ν , . . . , e − ȷζ T ( Q 1 − 1) ν  T , (15) ˘ a T ( ν ) =  1 , e − ȷ ˘ ζ T ν , . . . , e − ȷ ˘ ζ T ( I − 1) ν  T , (16) resp ectiv ely , where Ψ Ψ Ψ r = [ ˜ ϕ ϕ ϕ (1) r,t , . . . , ˜ ϕ ϕ ϕ ( Q 1 ) r,t ] , ∀ r , ˜ a I ( φ, θ ) = a I ( φ r ) ⊙ a I ( θ ) , ˘ ζ F = ζ F ∆ N , and ˘ ζ T = ζ T ∆ Q . B. Probabilistic Signal Representation With the ab ov e, w e establish the probabilistic signal mo del in the considered system, which is then applied to the Bay esian problem formulation. W e rst apply the geometric constraints to build the statistical relationship b et w een p osition-related parameters and user states. Us- ing the user r -BS g channel, ∀ r, g , as an example, we hav e p ( θ UB k,g ,t | ψ ψ ψ k,t ) = δ  θ UB k,g ,t − arccos(( e B g ) T e UB k,g ,t )  , (17) p ( τ UB k,g ,t | ψ ψ ψ k,t ) = δ  τ UB k,g ,t − c − 1 0 d UB k,g ,t  , (18) p ( ν UB k,g ,t | ψ ψ ψ k,t ) = δ  ν UB k,g ,t − λ − 1 ( ˙ p U k,t ) T e UB k,g ,t ,  , (19) where d UB k,g ,t = || p U k,t − p B g || 2 , ∀ k , t , is the distance betw een user k and the BS g at the time slot t , e UB k,g ,t = ( d UB k,g ,t ) − 1 ( p U k,t − p B g ) is the directional vector from the BS g to the user k at the time slot t , and e B g ∈ R 2 × 1 is the array directional vector of the BS g . The statistical relationships for user k -RIS r and RIS r -BS g channels are similar to (17) – (19). Note that the Doppler frequency for eac h RIS-BS channel is equal to zero as sho wn in (5), due to the static nature of RISs and BSs. Based on the mobilit y mo del (1), the state of eac h user k is mo deled to ev olv e according to a rst-order Marko v pro cess as p ( ψ ψ ψ k,t | ψ ψ ψ k,t − 1 ) = N ( ψ ψ ψ k,t ; F 0 ψ ψ ψ k,t − 1 , Q k,t − 1 ) , ∀ k , t. (20) W e denote ϑ ϑ ϑ t = [( ϑ ϑ ϑ 1 ,t ) T , . . . , ( ϑ ϑ ϑ K,t )] T as the p osition- related parameter v ector for all users in the time slot t . Therein, the p osition-related parameter vector of eac h user k is dened by ϑ ϑ ϑ k,t = [( ϑ ϑ ϑ UB k, 1 ,t ) T , . . . , ( ϑ ϑ ϑ UB k,G,t ) T , 4 When Q 1 is not very large, we observe that a T 1 ( ν ) ≈ 1 Q 1 , which can b e leveraged to simplify the signal mo del. F rom our sim ulation trials, we nd that there has little inuence on the joint estimation performance b y adopting this approximation. 6 ( ϑ ϑ ϑ UI k, 1 ,t ) T , . . . , ( ϑ ϑ ϑ UI k,R,t ) T ] T , ∀ k , t, (21) where ϑ ϑ ϑ c k,g ,t = [ θ c k,g ,t , τ c k,g ,t , ν c k,g ,t ] T , ∀ c ∈ { UB , UI } are the p osition-related parameter vector for each link. W e also denote the link indicator vector as α α α t = [ α α α T 1 ,t , . . . , α α α T G,t ] T with α α α g,t = [ α α α T 1 ,g,t , . . . , α α α T K,g ,t ] T and α α α k,g ,t = [ α UB k,g ,t , α UI k,g , 1 ,t , . . . , α UI k,g ,R,t ] T . By following (12), the conditional probability for the received signal at the BS g in the ISAC sub-blo ck of each time slot t ov er the p osition-related parameters of all users can b e given b y p ( y g,t | ϑ ϑ ϑ t , ˜ s t , α α α g,t ) = C N ( y g,t ; A g,t w g,t , σ 2 z I ) . (22) The probability density function (PDF) of the transmit signal is also denoted as P 0 ( ˜ s k,t ) , ∀ k , t . T o p erform online estimation, we can derive the joint probability of the receiv ed signals, transmit signals, states, and p osition- related parameters from the rst time slot to the t -th time slot as p ( { Y v , ϑ ϑ ϑ v , ψ ψ ψ v , ˜ s v , α α α v } t v =1 ) = t Y v =1 p ( y g,v | ϑ ϑ ϑ v , ˜ s v , α α α g,v ) ×  K Y k =1 p ( ϑ ϑ ϑ k,v | ψ ψ ψ k,v ) p ( ψ ψ ψ k,v | ψ ψ ψ k,v − 1 ) P 0 ( ˜ s k,v )  , (23) where p ( ψ ψ ψ k, 1 | ψ ψ ψ k, 0 ) = p ( ψ ψ ψ k, 1 ) denotes the prior informa- tion of the user states at the rst time slot derived from the initial estimation of user states. C. Ba y esian Problem F ormulation Under the probabilistic mo del, w e can follo w the Ba y es’ rule to derive the p osterior probabilities of the state and transmit signal of each user k in the time slot t as p ( ψ ψ ψ k,t |{ Y v } t v =1 ) =  \ ψ ψ ψ k,t p ( { Y v , ϑ ϑ ϑ v , ψ ψ ψ v , ˜ s v , α α α v } t v =1 ) p ( { Y v } t v =1 ) , (24) p ( ˜ s k,t |{ Y v } t v =1 ) =  \ ˜ s k,t p ( { Y v , ϑ ϑ ϑ v , ψ ψ ψ v , ˜ s v , α α α v } t v =1 ) p ( { Y v } t v =1 ) . (25) After that, w e can apply the criterion of MAP or MMSE to perform the online estimation at the time slot t . Mean while, the p osition-related parameters and link indi- cators can also b e derived simultaneously . How ever, it is intractable to derive the explicit expression of (24) – (25) due to high-dimensional in tegrals. Though the p opular KF-based metho ds [24] can approximately derive a solution, their p erformance suers severe deterioration in the case with limited tracking resources. As such, this w ork resorts to leveraging MP and VMP techniques under the Bay esian approximate inference framework to design an eective and ecient estimation algorithm. IV. Join t Multiuser T racking and Signal Detection Algorithm In this section, we prop ose a computationally ecient HVMP algorithm for the joint multi-user tracking and signal detection problem based on the built probabilistic mo del. A. Graphical Representation Before the introduction of the prop osed algorithm, we rst establish the graphical representation of factor graph                                                               󰇗     󰇗                          BS  BS            󰇗     󰇗                                     󰇗       󰇗                   󰇗       󰇗                   󰇗     󰇗       󰇗     󰇗                               RIS  Fig. 4. F actor graph of the probabilistic model. based on the join t probability (23) as shown in Fig. 4. Therein, squares and white circles denote the factor no des and v ariable no des, resp ectively . Moreov er, we prop ose to perform bidirectional message passing o v er the (black) edges for v ariable up date in each time slot, and p erform only forward message passing o ver the (red) edges for state up date, due to the online estimation. Note that v ariational message passing is p erformed along edges connected with the blue-v eined square and standard message passing is realized along the other edges. In the following, we in tro duce the details of the prop osed algorithm. B. Hybrid V ariational Message P assing In general, the prop osed HVMP algorithm integrates b oth the v ariational message passing and standard mes- sage passing to iteratively up date the extrinsic message on eac h parameter, which are utilized to derive the p osterior probability of user states and transmit signals after conv ergence. Sp ecically , the VMP is applied to estimate p osition-related parameters from the received signals, while the MP is utilized to rene the estimation of states and transmit signals. By emplo ying suc h a h ybrid of VMP and MP rules, the prop osed algorithm can up date all the extrinsic messages with tractable expressions to realize computationally ecient online estimation for user states and transmit signals. W e rst introduce the bidirectional message passing equations sequentially and then giv e the forw ard message passing equations for online estimation. By leveraging the v ariational techniques in the message passing procedure, w e consider to mo del distributions with b oth of the p osition and v elo cit y b y Gaussian distributions, while the v on Misés (VM) distribution [34], [35] is used for modeling p osition-related parameters. 1) Messages b etw een ˙ p k,t and ˙ f k,t : The bidirectional messages in the j -th iteration can b e given b y ε ˙ p k,t → ˙ f k,t ( j ) = N ( ˙ p k,t ; ˙ − → p k,t , ˙ − → Q k,t ) , (26) ε ˙ p k,t ← ˙ f k,t ( j ) ∝  p ( { ν UB k,g ,t } , { ν UI k,r,t }| ˙ p k,t ) × G  g =1 ε ˙ f k,t ← ν UB k,g ,t ( j ) R  r =1 ε ˙ f k,t ← ν UI k,r,t ( j ) 7 = N ( ˙ p k,t ; ˙ ← − p k,t ( j ) , ˙ ← − Q k,t ( j )) , (27) where ˙ − → p k,t and ˙ − → Q k,t are the mean and co v ariance of the v elo cit y from the online estimation in the last time slot t − 1 , ˙ ← − p k,t ( j ) and ˙ ← − Q k,t ( j ) are the up dated ones from the Doppler frequency estimation in the j -th iteration of the HVMP algorithm. The explicit expression of ˙ ← − p k,t ( j ) and ˙ ← − Q k,t ( j ) is hard to b e deriv ed due to the complicate relationship b et w een the Doppler frequencies and the v elo cit y , and we resort to the Gaussian-Newton metho d and T aylor series expansion to obtain ˙ ← − p k,t ( j ) and ˙ ← − Q k,t ( j ) . The detailed deriv ation is provided in App endix A. 2) Messages b etw een p k,t and f k,t : The bidirectional messages in the j -th iteration are derived by ε p k,t → f k,t ( j ) = N ( p k,t ; − → p k,t , − → Q k,t ) , ∀ j, (28) ε p k,t ← f k,t ( j ) ∝  p ( { p UB k,g ,t } , { p UI k,r,t }| p k,t ) G  g =1 ε f k,t ← p UB k,g ,t ( j ) × R  r =1 ε f k,t ← p UI k,r,t ( j ) = N ( p k,t ; ← − p k,t ( j ) , ← − Q k,t ( j )) , (29) where ← − p k,t ( j ) =  ← − Q k,t ( j )  − 1  G  g =1  ← − Q UB k,g ,t ( j )  − 1 ← − p UB k,g ,t ( j ) (30) + R  r =1  ← − Q UI k,r,t ( j )  − 1 ← − p UI k,r,t ( j )  , ← − Q k,t ( j ) =  G  g =1  ← − Q UB k,g ,t ( j )  − 1 + R  r =1  ← − Q UI k,r,t ( j )  − 1  − 1 . (31) In the ab o v e messages, − → p k,t ( − → Q k,t ) is the mean (cov ari- ance) of the p osition from the online estimation in the last time slot t − 1 and the explicit expression of the message ε f k,t ← p UB k,g ,t ( j ) ( ε f k,t ← p UI k,r,t ( j ) ) is given in (38). 3) Messages b etw een ˜ s k,t and l k,t : The message from ˜ s k,t to l k,t is obtained based on the prior distribution of the transmit signal, i.e., ε ˜ s k,t → l k,t ( j ) = P 0 ( ˜ s k,t ) , and the message from l k,t to ˜ s k,t is derived by ε ˜ s k,t ← l k,t ( j ) ∝  p ( { w UB k,g ,t } , { w UI k,r,t }| ˜ s k,t ) G  g =1 ε l k,t ← w UB k,g ,t ( j ) × R  r =1 ε l k,t ← w UI k,r,t ( j ) = C N ( ˜ s k,t ; ← − s k,t ( j ) , ← − ρ k,t ( j )) , ∀ j, (32) where ← − s k,t ( j ) = ( ← − ρ k,t ( j )) − 1  G  g =1 ( ← − ϖ UB k,g ,t ( j )) − 1 ← − w UB k,g ,t ( j ) (33) + R  r =1 ( ← − ϖ UI k,r,t ( j )) − 1 ← − w UI k,r,t ( j )  , ← − ρ k,t =  ( ← − ϖ UB k,g ,t ( j )) − 1 + R  r =1 ( ← − ϖ UI k,r,t ( j )) − 1  − 1 , (34) with ε l k,t ← w UB k,g ,t ( j ) = C N ( w UB k,g ,t ; ← − w UB k,g ,t ( j ) , ← − ϖ UB k,g ,t ( j )) and ε l k,t ← w UI k,r,t ( j ) = C N ( w UI k,r,t ; ← − w UI k,r,t ( j ) , ← − ϖ UI k,r,t ( j )) . 4) Messages b etw een ˙ f k,t and ν UB k,g ,t ( ν UI k,r,t ): The bidi- rectional message passing equations can b e written as ε ˙ f UB k,t → ν UB k,t ( j ) ∝ Z p ( { ν UB k,g ,t } , { ν UI k,r,t }| ˙ p k,t ) · Y g ′  = g ε g k,t ← ν UB k,g ′ ,t ( j ) R Y r =1 ε g k,t ← ν UI k,r,t ( j ) = V M ( − ζ T ν UB k,t ; − → µ UB k,g ,t ( ν, j ) , − → κ UB k,g ,t ( ν, j )) , (35) ε ˙ f UB k,t ← ν UB k,g ,t ( j ) = ε ν UB k,g ,t ← py t ( j ) = V M ( − ζ T ν UB k,t ; ← − µ UB k,g ,t ( ν, j ) , ← − κ UB k,g ,t ( ν, j )) , (36) where − → µ UB k,g ,t ( ν, j ) and − → κ UB k,g ,t ( ν, j )) are also derived by the Gauss-Newton metho d as detailed in App endix A; ε ˙ f UB k,t ← ν UB k,g ,t ( j ) ( ε ˙ f UI k,t ← ν UI k,r,t ( j ) ) is derived from the v aria- tional message passing. The messages b etw een ˙ f k,t and ν UI k,r,t can b e similarly obtained to (35) and (36). 5) Messages b etw een f k,t and p UB k,g ,t ( p UI k,r,t ): Due to the Gaussian message passing along the edges connected with no de f k,t , we can simply derive the message ε f k,t → p UB k,g ,t ( j ) ( ε f k,t → p UI k,r,t ( j ) ) similar to (29) as ε f k,t → p UB k,g ,t ( j ) ∝  p ( { p UB k,g ,t } , { p UI k,r,t }| p k,t ) ε p k,t → f k,t ( j ) ×  g ′  = g ε f k,t ← p UB k,g ′ ,t ( j ) R  r =1 ε f k,t ← p UI k,r,t ( j ) = N ( p UB k,g ,t ; − → p UB k,g ,t ( j ) , − → Q UB k,g ,t ( j )) , (37) and the message ε f k,t ← p UB k,g ,t ( j ) ( ε f k,t ← p UI k,r,t ( j ) ) can b e obtained as ε f k,t ← p UB k,g ,t ( j ) = Z p ( θ UB k,g ,t , τ UB k,g ,t | p UB k,g ,t ) × ε f UB k,g ,t ← θ UB k,g ,t ( j ) ε f UB k,g ,t ← τ UB k,g ,t ( j ) = N ( p UB k,g ,t ; ← − p UB k,g ,t ( j ) , ← − Q UB k,g ,t ( j )) , (38) where the mean ← − p UB k,g ,t ( j ) and cov ariance ← − Q UB k,g ,t ( j ) are deriv ed based on second-order T a ylor series expansion and Gaussian approximation based on the up dated messages on θ UB k,g ,t and τ UB k,g ,t in the part of v ariational message passing. The detailed deriv ation of (38) is omitted due to page limits. Similarly , the messages betw een f k,t and p UI k,r,t can b e obtained by following the ab ov e calculations. 6) Messages b etw een l k,t and w UB k,g ,t ( w UI k,g ,r,t ): The bidirectional messages can b e calculated by ε l k,t → w UB k,g ,t ( j ) ∝  p ( { w UB k,g ,t } , { w UI k,g ,r,t }| ˜ s k,t ) ×  g ′  = g ε l k,t ← w UB k,g ′ ,t ( j ) R  r =1 ε l k,t ← w UI k,g ,r,t ( j ) , (39) ε l k,t ← w UB k,g ,t ( j ) = C N ( w UB k,g ,t ; ← − w UB k,g ,t ( j ) , ← − ϖ UB k,g ,t ( j )) , (40) where ← − w UB ,j k,t and ← − ϖ UB ,j k,t are derived from the v ariational message passing. Messages of the v ariable w UI k,r,t is similarly deriv ed by following (39) and (40). 7) Messages betw een f UB k,g ,t and θ UB k,g ,t /τ UB k,g ,t (b et ween f UI k,r,t and θ UI k,r,t /τ UI k,r,t ): By follo wing (35), w e can also obtain the messages ε f UB k,g ,t → θ UB k,g ,t ( j ) and ε f UB k,g ,t → τ UB k,g ,t ( j ) in the form of VM distributions. Alternatively , the mes- sages ε f UB k,t ← θ UB k,t ( j ) = ε θ UB k,g ,t ← py t ( j ) and ε f UB k,t ← τ UB k,t ( j ) = 8 ε τ UB k,t ← py t ( j ) are obtained from the conv erged results for the v ariational message passing in the i -th outer iteration of HVMP . The message passing for user-RIS paths can b e similarly obtained. 8) Messages b etw een py t and ϑ ϑ ϑ UB k,g ,t ( ϑ ϑ ϑ UI k,r,t ) via v ari- ational message passing: F rom the signal mo del (12), this part aims to solve the multidimensional line sp ectra estimation problem [36] to derive the extrinsic messages o v er the v ariable vector ϑ ϑ ϑ t . Under the Ba y esian infer- ence framework, the standard sum-pro duct algorithm is in tractable to b e employ ed in this part, since all the v ariables ϑ ϑ ϑ t are coupled with complicated and non-linear dep endencies in the received signals { y g,t } . Motiv ated b y the v ariational line sp ectra inference metho ds in [25], [34], [36], the VMP is employ ed in this part to realize ecien t deriv ation of the desired extrinsic messages. How ever, the metho ds in [ 34], [36] cannot b e directly applied due to the b eamspace signal structure in (12), where the receiv ed signal at the RIS cannot b e fully observed and is pro jected in to a lo w-dimensional subspace { Ψ Ψ Ψ T r } , as sho wn in (13). On the other hand, the v ariational metho d in [25] concentrates on the single-BS case and cannot fully exploit the dep endencies of each elemen t in ϑ ϑ ϑ UB k,g ,t ( ϑ ϑ ϑ UI k,r,t ) b y optimizing each element individually via co ordinate descen t. T o address the ab ov e c hallenges, we prop ose to generalize the v ariational inference-based metho d to deal with the b eamspace multidimensional line spectra estimation problem. The objective of VMP is to appro ximate the PDF of the p osterior probabilit y ov er these v ariables by a surrogate PDF with the minimal Kullbac k-Leibler (KL) divergence 5 . Due to the fact that ln p Y Y Y ( Y t ) = D KL ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y || p ϑ ϑ ϑ, W W W , α α α | Y Y Y ) + L ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ) , L ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ) = E χ ϑ ϑ ϑ, W W W ,α α α | Y Y Y  ln p ϑ ϑ ϑ, W W W , α α α, Y Y Y ( ϑ ϑ ϑ t , W t , α α α t , Y t ) χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ( ϑ ϑ ϑ t , W t , α α α t | Y t )  , w e can turn to maximize the evidence low er b ound (ELBO) L ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ) given that p Y Y Y ( Y t ) is constant for the received signal Y t . Note that the v ariables { w UB k,g ,t } and { ϑ ϑ ϑ UB k,g ,t } ( { w UI k,g ,r,t } , { ϑ ϑ ϑ UI k,r,t } ) are all considered to b e indep enden t in the VMP part. Therefore, the surrogate PDF can b e factorized in to the pro duct of a set of indep enden t factors as χ ϑ ϑ ϑ, W W W , α α α | Y Y Y = G Y g =1 χ W W W | Y Y Y ( w g,t | Y t ) χ α α α | Y Y Y ( α α α g,t | Y t ) × K Y k =1 G Y g =1 ( χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ UB k,g ,t | Y t ) R Y r =1 χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ UI k,r,t | Y t )) ! , (41) where χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ k,t | Y t ) = Q 3 i =1 V M ([ ϑ ϑ ϑ k,t ] i ; [ ˆ ϑ ϑ ϑ k,t ] i , [ ˆ κ κ κ k,t ] i ) , χ W W W | Y Y Y ( w g,t | Y t ) is modeled b y the Gaussian distribution, and χ α α α | Y Y Y ( α α α g,t | Y t ) = Q i δ ([ α α α g,t ] i − [ ˆ α α α g,t ] i ) . T o maximize the ELBO, w e resort to the ecien t block coordinate de- scen t (BCD) metho d to optimize eac h factor in (41), which can s ignican tly outperform the CD metho d adopted in [25] when the user num b er is large. In eac h outer iteration 5 The KL divergence of p 1 ( x ) to p 2 ( x ) (within the same feasible set X ) is dened as D KL ( p 1 || p 2 ) = ∫ X p 1 ( x ) ln p 1 ( x ) p 2 ( x ) dx . of the HVMP algorithm, we p erform multiple inner iter- ations of VMP to up date these factors until conv ergence. In the following, we shall in tro duce the updating rules for these factors sequentially . F or conv enience, the indices of b oth outer and inner iterations are omitted. First, maximizing L ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ) w.r.t. the factor χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ UB k,g ,t | Y t ) giv es ln χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ UB k,g ,t | Y t ) = ln ε ˙ f k,t → ν UB k,g ,t + ln ε f UB k,g ,t → θ UB k,g ,t + ln ε f UB k,g ,t → τ UB k,g ,t + ℜ{ ( η η η UB k,g ,t ) H a UB k,g ,t } , (42) where η η η UB k,g ,t = 2 σ 2 z  ( ˆ w UB k,g ,t ) ∗ y g,t −  k ′  = k M 2 ( w UB k,g ,t , w UB k ′ ,g,t ) a UB k ′ ,g,t − K  k ′′ =1 R  r =1 M 2 ( w UB k,g ,t , w UI k ′′ ,g,r,t ) a UI k ′′ ,g,r,t  . (43) In (43), M 2 ( X, Y ) = E X Y [ X H Y ] denotes the exp ectation of X H Y ov er the up dated PDF in the last inner itera- tion. Similarly , maximizing L ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ) w.r.t. the factor χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ UI k,r,t | Y t ) also gives ln χ ϑ ϑ ϑ | Y Y Y ( ϑ ϑ ϑ UI k,r,t | Y t ) = ln ε ˙ f k,t → ν UI k,g ,t + ln ε f UI k,g ,t → θ UI k,g ,t + ln ε f UI k,g ,t → τ UI k,g ,t + G X g =1 ℜ{ ( η η η UI k,g ,r,t ) H a UI k,g ,r,t } − ℜ{ ( ˜ η η η UI k,r,t ) H a I ( θ UI k,r,t ) } , (44) where η η η UI k,g ,r,t = 2 σ 2 z  − K  k ′ =1 M 2 ( w UI k,g ,r,t , w UB k ′ ,g,t ) a UB k ′ ,g,t −  k ′′  = k R  r =1 M 2 ( w UI k,g ,r,t , w UI k ′′ ,g,r,t ) a UI k ′′ ,g,r,t + ( ˆ w UI k,g ,r,t ) ∗ y g,t  , and eac h m I -th elemen t in ˜ η η η UI k,r,t is obtained b y (45) at the top of the next page. Note that the last term exists in (44) due to the beamspace signal structure by the cascade c hannel of the user-RIS-BS links. Moreov er, the v ariable v ector ϑ ϑ ϑ UB k,g ,t is only statistically dependent on the receiv ed signal at the BS g , while ϑ ϑ ϑ UI k,r,t is statistically dep endent on the received signal at all BSs, which contributes to the hierarchical signal pro cessing op eration. The v ariable v ector ϑ ϑ ϑ UB k,g ,t in (42) is then rened via the Newton metho d [36] as ˆ ϑ ϑ ϑ UB k,g ,t = ˆ ϑ ϑ ϑ UB , old k,g ,t − ( ▽ ▽ ▽ 2 G B ( ˆ ϑ ϑ ϑ UB , old k,g ,t )) − 1 ▽ ▽ ▽ G B ( ˆ ϑ ϑ ϑ UB , old k,g ,t ) , (46) where G B ( ϑ ϑ ϑ UB k,g ,t ) is dened as ( 42), and ▽ ▽ ▽ G B ( ϑ ϑ ϑ UB , old k,g ,t ) ( ▽ ▽ ▽ 2 G B ( ϑ ϑ ϑ UB , old k,g ,t )) denotes the gradien t (Hessen matrix) at ϑ ϑ ϑ UB , old k,g ,t . Similarly , the Newton metho d is also utilized to rene the estimates of ϑ ϑ ϑ UI k,r,t with (44). Then, the concentration parameter v ectors ˆ κ κ κ UB k,g ,t = [ ˆ κ UB k,g ,t ( θ ) , ˆ κ UB k,g ,t ( τ ) , ˆ κ UB k,g ,t ( ν )] T and ˆ κ κ κ UI k,r,t = [ ˆ κ UI k,r,t ( θ ) , ˆ κ UI k,r,t ( τ ) , ˆ κ UI k,r,t ( ν )] T can also be up dated b y following [36]. Next, maximizing L ( χ ϑ ϑ ϑ, W W W , α α α | Y Y Y ) w.r.t. the factor χ W W W | Y Y Y ( w g,t | Y t ) gives χ W W W | Y Y Y ( w g,t | Y t ) = C N ([ w g,t ] ˆ A g,t ; [ ˆ w g,t ] ˆ A g,t , [ ˆ C g,t ] ˆ A g,t , ˆ A g,t ) × Y j / ∈ ˆ A δ ([ w g,t ] j ) , (47) 9 [ ˜ η η η UI k,r,t ] m I =  M B D I QM I σ 2 z  G g =1 M 2 ( w UI k,g ,r,t , w UI k,g ,r,t ) , m I = 1 , 2 M B D I Q 2 σ 2 z  G g =1  M I i = m I M 2 ( w UI k,g ,r,t , w UI k,g ,r,t )[ ¯ Ψ Ψ Ψ ∗ r,t ¯ Ψ Ψ Ψ T r,t ] i,i − m I +1 , m I  = 1 , (45) where ˆ A g,t is the supp ort set of the vector ˆ w g,t and [ ˆ w g,t ] ˆ A g,t is the vector containing the elemen ts in the subset ˆ A g,t . In (47), the mean and co v ariance of [ w g,t ] ˆ A g,t is given b y [ ˆ w g,t ] ˆ A g,t = σ − 2 z [ ˆ C g,t ] ˆ A g,t , ˆ A g,t [ ˆ A g,t ] H : , ˆ A g,t y g,t , [ ˆ C g,t ] ˆ A g,t , ˆ A g,t = σ 2 z ([ J g,t ] ˆ A g,t , ˆ A g,t + σ 2 z diag ([ ˘ β β β g,t ] ˆ A g,t )) − 1 , where [ ˘ β β β g,t ] i = [ β β β g,t ] − 1 i , ∀ i and J g,t = E [ A H g,t A g,t ] . Since χ α α α | Y Y Y ([ α α α g,t ] i | Y t ) = δ ([ α α α g,t ] i − [ ˆ α α α g,t ] i ) , we prop ose to determine [ ˆ α α α g,t ] via hypothesis test given b y log  p ([ ˆ w g,t ] i | [ ˆ α α α g,t ] i = 1) p ([ ˆ w g,t ] i | [ ˆ α α α g,t ] i = 0)  H 1 ≷ H 0 γ i , (48) where γ i is the pre-determined threshold obtained by empirical results. Compared with the greedy iterative searc h metho d in [34], the hypothesis test-based metho d can determine all of these indicator v ariables in parallel without an iterative searching pro cess to nd the lo cal optimal solution and is thus computationally ecien t in large-scale systems. Based on (48), the supp ort set ˆ A g,t can also b e updated. By iteratively executing the ab ov e calculation pro cedure, the estimation of { ϑ ϑ ϑ t , W t , α α α t } can b e rened. After the VMP algorithm con v erges, w e adopt the v ariable ϵ θ UB k,g ,t as an example to show the deriv ation of the extrinsic message by ε θ UB k,g , t ← py t = V M ( θ UB k,g ,t ; ← − µ UB k,g ,t ( θ ) , ← − κ UB k,g ,t ( θ )) , (49) where ← − κ UB k,g ,t ( θ ) e ȷ ← − µ UB k,g ,t ( θ ) = ˆ κ UB k,g ,t ( θ ) e − ȷζ S ˆ θ UB k,g ,t − − → κ UB k,g ,t ( θ ) e ȷ − → µ UB k,g ,t ( θ ) , (50) since the family of VM distributions is closed under the m ultiplication and division op erations. Similarly , extrinsic messages on the other p osition-related parameters can also b e derived. 9) F orward messages: After the state estimation in the curren t time slot t via the ab o v e bidirectional message passing is conv erged, w e can up date the forward messages to the time slot t + 1 for online user tracking and signal detection. These forward messages can b e expressed b y ε ˙ u k,t +1 → ˙ p k,t +1 = ε ˙ p k,t → ˙ u k,t +1 = ε ˙ p → u k,t +1 ∝ ε ˙ u k,t → ˙ p k,t ( ∞ ) ε ˙ p k,t ← ˙ f k,t ( ∞ ) = N ( ˙ p k,t +1 ; ˙ − → p k,t +1 , ˙ − → Q k,t +1 ) , (51) ε p k,t → u k,t +1 ∝ ε u k,t → p k,t ε p k,t ← f k,t ( ∞ ) = N ( p k,t +1 ; ˆ p k,t +1 , ˆ Q k,t +1 ) , (52) ε u k,t +1 → p k,t +1 ∝ ε ˙ p → u k,t +1 ε p k,t → u k,t +1 = N ( p k,t +1 ; − → p k,t +1 , − → Q k,t +1 ) . (53) Since all the ab o ve messages in the deriv ation of (51) – (53) are in the family of Gaussian distributions, explicit expressions of these messages for velocities and p ositions are easy to b e obtained and thus omitted here. In summary , we iteratively p erform bidirectional mes- sage passing to up date the states and position-related pa- Algorithm 1 The HVMP Algorithm. Input: p ( ψ ψ ψ 0 ) and the received signals { Y t } ; Output: States { ˆ ψ ψ ψ t } and signals { ˆ s t } ; 1: for t = 1 to V do 2: for j = 1 to J out do 3: Up date ε ˙ p k,t → ˙ f k,t ( j ) and ε p k,t → f k,t ( j ) in (26) and (28), resp ectively; 4: Up date ε f k,t → p UB k,g ,t ( j ) / ε f k,t → p UI k,r,t ( j ) in (37); 5: Up date ε p UB k,g ,t → f UB k,g ,t ( j ) = ε f k,t → p UB k,g ,t ( j ) and ε p UI k,r,t → f UI k,r,t ( j ) = ε f k,t → p UI k,r,t ( j ) ; 6: Up date ε ˙ f k,t → ν UB k,g ,t ( j ) / ε ˙ f k,t → ν UI k,r,t ( j ) in (35); 7: Up date ε f UB k,g ,t → θ UB k,g ,t ( j ) / ε f UI k,r,t → θ UI k,r,t ( j ) similar to (35); 8: Up date ε f UB k,g ,t → τ UB k,g ,t ( j ) / ε f UI k,r,t → τ UI k,r,t ( j ) similar to (35); 9: Up date ε l k,t → w UB k,g ,t ( j ) / ε l k,t → w UI k,r,t ( j ) in (39); 10: P erform VMP in Sec. IV-B8) (with maximal J inn inner iterations) to up date ε θ UB k,g ,t ← py t / ε θ UI k,r,t ← py t ( j ) , ε τ UB k,g ,t ← py t ( j ) / ε τ UI k,r,t ← py t ( j ) , ε ν UB k,g ,t ← py t ( j ) / ε ν UI k,r,t ← py t ( j ) , and ε w UB k,g ,t ← py t ( j ) / ε w UI k,r,t ← py t ( j ) ; 11: Up date ε ˙ f k,t ← ν UB k,g ,t ( j ) / ε ˙ f k,t ← ν UI k,r,t ( j ) in (36); 12: Up date ε f UB k,g ,t ← θ UB k,g ,t ( j ) / ε f UI k,r,t ← θ UI k,r,t ( j ) similar to (36); 13: Up date ε f UB k,g ,t ← τ UB k,g ,t ( j ) / ε f UI k,r,t ← τ UI k,r,t ( j ) similar to (36); 14: Up date ε l k,t ← w UB k,g , t ( j ) / ε l k,t ← w UI k,r,t ( j ) in (40); 15: Up date ε f k,t ← p UB k,g ,t ( j ) / ε f k,t ← p UI k,r,t ( j ) in (38); 16: Up date ε ˙ p k,t ← ˙ f k,t ( j ) , ε p k,t ← f k,t ( j ) , and ε ˜ s k,t ← l k,t ( j ) in (27), (29), and (32), resp ectively; 17: end for 18: Up date ε ˙ u k,t +1 → ˙ p k,t +1 in (51); 19: Up date ε p k,t → u k,t +1 in (52); 20: Up date ε u k,t +1 → p k,t +1 in (53); 21: end for rameters based on the prior information from the time slot t − 1 until HVMP conv erges, and then utilize the conv erged estimation results to rene the prior information of user states in the next time slot t + 1 via forward message passing, which can realize the online tracking and signal detection for all users with div erse and high mobilit y . The whole algorithm is outlined in Algorithm 1. Remark 3: It is highlighted that the prop osed HVMP algorithm lev erages MP techniques to eectively combine measuremen ts from LOS and reecting links, where the blo c kage conditions of dieren t links, i.e., { α α α g,t } , are auto- matically determined during eac h iteration. Therefore, the HVMP algorithm has sup erior exibilit y and accuracy by iterativ ely estimating the blo ckage condition of each link during the parameter estimation pro cedure. Remark 4: By considering BS co op eration, the fronth ual links for aggregating receiv ed signals at each BS may ha v e nite capacity due to the hardware limit [37]. Th us, the 10 quan tize-and-forw ard (QF) scheme [37] can b e adopted and the prop osed algorithm can b e applied to solv e the considered problem by combining the dequantization tec hnique in [37]. How ev er, the utilized quantization sc heme also has a large inuence on the joint estimation p erformance, which is out of the scop e for this work and can b e left as future works. V. P erformance Analysis and Optimization A. P erformance Analysis via Bay esian Cramér-Rao Bound In this w ork, we adopt the Ba yesian Cramér-Rao Bound to characterize the tracking p erformance of users, which is widely used in the estimation problem with random parameters. W e dene the estimated parameter vector ˜ ψ ψ ψ t = [ ψ ψ ψ T t , ˜ s T t ] T ∈ C 5 K × 1 . Then, the Bay esian information matrix (BIM) for the discrete-time tracking nonlinear ltering metho d [26] can b e given by the follo wing recursive deriv ation of B t = B mea t + Ξ Ξ Ξ 22 t − Ξ Ξ Ξ 21 t ( B t − 1 + Ξ Ξ Ξ 11 t ) − 1 Ξ Ξ Ξ 12 t , (54) where B mea t = E y , ˜ ψ ψ ψ t  − ∂ 2 log p ( y t | ˜ ψ ψ ψ t ) ∂ ˜ ψ ψ ψ t ∂ ˜ ψ ψ ψ H t  , (55) Ξ Ξ Ξ ij t = E ˜ ψ ψ ψ t , ˜ ψ ψ ψ t − 1   − ∂ 2 log p ( ˜ ψ ψ ψ t | ˜ ψ ψ ψ t − 1 ) ∂ ˜ ψ ψ ψ t + i − 2 ∂ ˜ ψ ψ ψ H t + j − 2   , ∀ i, j ∈ { 1 , 2 } . (56) The detailed deriv ation of the measurement BIM B mea t are omitted here due to page limits. F or the initialization of the BCRB calculation, w e ha v e B 0 = E ˜ ψ ψ ψ 0  − ∂ 2 log p ( ˜ ψ ψ ψ 0 ) ∂ ˜ ψ ψ ψ 0 ∂ ˜ ψ ψ ψ H 0  . Then, for eac h tracking time slot t , we can obtain the estimation result ¯ ψ ψ ψ t satises E h ( b ψ ψ ψ t − ˜ ψ ψ ψ t )( b ψ ψ ψ t − ˜ ψ ψ ψ t ) H i ⪰ B − 1 t , ∀ t. (57) Therefore, the BCRB serves as the low er b ound for the parameter estimation p erformance, whic h is p oten tial to help optimize the RIS phase prole during the tracking pro cedure. B. RIS Phase Prole Optimization During the tracking pro cedure, the RIS phase prole can b e adaptively optimized to realize high-accuracy trac king and signal detection p erformance. As introduced in the last subsection, the BCRB can be leveraged as the p erformance metric for the state estimation p erformance and the signal detection p erformance. As suc h, we can form ulate an optimization problem with the ob jective function of the BCRB ov er the v ariable of the RIS phase prole { ς ( q ) m I ,r = ∠ ϕ ( q ) m I ,r } . In con trast to directly minimizing the BCRB, we consider to minimize the general weigh ted BCRB and the optimization problem can b e expressed b y ( P1 ) min { ς ς ς r } T race { diag ( v ) B − 1 t } , (58) where v = [ v 1 , . . . , v 5 K ] T is the v ector to adjust the w eight on eac h estimated parameter. 6 T o solv e Problem (P1), w e 6 F or example, w e may only care about the estimation accuracy of the user p ositions and transmit signals, meaning that w e can set the weigh ts related to velocities to be zeros. lev erage the gradien t-based algorithm that can obtain a lo cal optimal solution due to the non-conv ex objective function. In detail, in each iteration e , the up dating rule for the RIS phase prole is given b y ς ς ς ( e ) r = ς ς ς ( e − 1) r − ℵ ( e ) ▽ ▽ ▽ F ( { ς ς ς r } ) , ∀ e, r, (59) where the step size ℵ ( e ) is determined by the Armijo bac ktrac king line searc h. The gradien t-based method will iterativ ely up date the RIS phase prole until the stopping criterion P R r =1 ||ℵ ( e ) ▽ ▽ ▽ r F ( { ς ς ς r } ) || < ϵ is satised. VI. Sim ulation Results This section provides simulation results to verify the eectiv eness of the prop osed algorithm for the RIS- assisted ISAC system. In the system setup, we consider that there are G = 2 BSs deploy ed at p B 1 = [0 , 0] T m and p B 2 = [90 , 0] T m, R = 2 RISs deploy ed at p I 1 = [20 , 40] T m and p I 2 = [60 , 40] T m, and K = 3 users, whic h mo ve in the square region with the center [45 , − 5] T m and the side length 50 m. W e also consider that av erage v elo cities of these three users are ¯ υ 1 = 40 m/s, ¯ υ 2 = 30 m/s, ¯ υ 3 = 15 m/s, where ¯ υ k = T − 1 P T t =1 || ˙ p U k,t || 2 , ∀ k . Eac h BS is equipp ed with M B = 6 antennas and each RIS is equipp ed with M I = 48 elements. W e consider that the bandwidth is N ∆ f = 10 MHz with N = 12 sub carriers and N 0 = { 1 , 2 , . . . , 12 } with ∆ N = 1 . The trac king interv al is set as ∆ T = 0 . 02 s and the ISAC sub-blo c k is set with Q 1 = 12 , I = 10 , and ∆ Q = 200 . The path loss factor is mo deled by β = λ 4 πd [12], [25] for each LOS link. W e consider that each user transmits Gaussian signals if not sp ecied. All the simulation results are a v eraged o ver 100 channel realizations. The track- ing p erformance is ev aluated by the metric of the ro ot MSE (RMSE): RMSE = T − 1 P T t =1 || ˆ ψ ψ ψ t − ψ ψ ψ t || 2 , and the signal detection p erformance is ev aluated by the MSE: MSE = T − 1 P T t =1 || ˆ s t − ˜ s t || 2 2 [38]. F or comparison, w e adopt three representativ e metho ds as b enchmarks: 1) Benchmark 1: MUSIC with spatial smo othing [39] is utilized to estimate the user states from the received signal Y t and then KF is applied to rene the state estimation; 2) Benchmark 2: the v ariational inference algorithm only trac ks the p osition of eac h user without Doppler frequency (v elo cit y) estimation, which is similar to that in [25]; 3) Benchmark 3: we consider to p erform m ulti-user trac king under the pilot-based transmission protocol with- out p erforming signal detection in the pilot phase, whic h can serve as one low er b ound on multi-user trac king. First, we ev aluate the conv ergence p erformance in Fig. 5 under dierent transmit p ow ers. W e can observe that the prop osed HVMP algorithm can quickly conv erge to the stationary solutions within several outer iterations. Compared with Benchmark 2, the HVMP algorithm can consisten tly provide signican tly smaller estimation errors for b oth states and transmit signals. In particular, it is observed that the initial estimation error of HVMP , kno wn as the lo cation prediction error, is also small and approac hing the con verged estimation error, indicating the eectiv eness of exploiting the full state of each user for 11 1 5 10 15 20 10 -2 10 -1 10 0 1 5 10 15 20 10 -4 10 -3 10 -2 10 -1 10 0 10 1 1 5 10 15 20 10 -1 10 0 Fig. 5. The conv ergence p erformance of the prop osed algorithm: (left) p osition estimation, (middle) signal detection, (right) velocity estimation. 0 10 20 30 10 -2 10 -1 10 0 10 1 10 2 0 10 20 30 10 -4 10 -3 10 -2 10 -1 10 0 10 1 Fig. 6. The p erformance comparison under dierent transmit p ow ers: (left) p osition estimation, (righ t) signal detection. trac king p erformance enhancement. It is also v alidated that Benc hmark 3 performs as a low er b ound for the state estimation of the prop osed scheme, which utilize pilot signals to enhance the parameter estimation p erformance. In the following, we only ev aluate the performance of the lo cation and transmit signal estimation, which act as the k ey p erformance metrics for our considered problem. W e then ev aluate the impact of the RIS phase prole design on the p erformance of the prop osed algorithm and b enc hmarks in Fig. 6. Three kinds of RIS phase prole design are considered: random phase prole, DFT phase prole, and optimized phase prole. F or the DFT phase prole, w e consider that the ⌊ Q 1 K ⌋ DFT co dew ords with the closest distances to ˜ a I ( ˜ θ UI k,t ) are selected for each user k at each RIS r , where ˜ θ UI k,r,t is the predicted AO A from user k to the RIS r . W e can nd that all of the prop osed algo- rithm and b enchmarks can b enet from the optimized RIS phase prole to enhance the tracking and signal detection p erformance, v alidating the eectiveness of the prop osed optimization metho d. Due to the elab orated up dating rules and full exploitation of user states, the prop osed HVMP algorithm can signicantly outp erform the tw o represen tativ e b enchmarks and approac h the p erformance 2 4 8 12 16 20 24 10 -3 10 -2 10 -1 10 0 10 1 10 2 2 4 8 12 16 20 24 10 -4 10 -2 10 0 10 2 Fig. 7. The performance comparison under dieren t num b ers of OFDM sym bols Q 1 in eac h group of the ISAC sub-blo ck with P = 30 dBm: (left) p osition estimation, (right) signal detection. 1 4 7 10 10 -3 10 -2 10 -1 10 0 10 1 1 4 7 10 10 -5 10 -4 10 -3 10 -2 10 -1 10 0 10 1 Fig. 8. The performance comparison under dieren t num b ers of users K with P = 30 dBm: (left) p osition estimation, (right) signal detection. limit of BCRB. Moreo ver, w e can observe that Benc hmark 2 can provide sup erior p erformance to Benc hmark 1 even without Doppler frequency estimation, since Benc hmark 1 via MUSIC suers from p o or p erformance on the state estimation in the case with inadequate time samples (snapshots). The p erformance of the prop osed algorithm with dif- feren t num b ers of OFDM symbols Q 1 in each group of the ISAC sub-blo ck is given in Fig. 7. Note the condition Q 1 ≥ 2 needs to b e satised since we cannot extract the A O A information of user-RIS links with xed RIS phase prole during each tracking time slot, i.e., Q 1 = 1 . F rom Fig. 7, the tracking and signal detection errors are rst reduced rapidly by enlarging Q 1 , since we can get muc h more accurate AO A information for each user-RIS link to enhance the lo cation estimation p erformance. When Q 1 exceeds a certain v alue, the prop osed algorithm turns to ha v e saturated tracking p erformance, since the distance estimation errors dep enden t on the bandwidth b ecome the limiting factor for the tracking p erformance at this momen t. This result also rev eals that we only require a mo derate num b er of OFDM symbols to accomplish the high-precision trac king task, whic h can sa v e the trac king o v erhead to promote the data transmission p erformance. W e also ev aluate the p erformance of the prop osed algorithm under dierent num b ers of users in Fig. 8, where the a v erage v elo cities of dieren t users are in 12 0 0.2 0.4 0.6 0.8 1 10 -3 10 -2 10 -1 10 0 0 0.2 0.4 0.6 0.8 1 10 -6 10 -5 10 -4 10 -3 10 -2 10 -1 10 0 Fig. 9. The performance comparison under dierent probabilities of LOS blo ckage for user-BS links with P = 30 dBm and M B = 12 : (left) p osition estimation, (righ t) signal detection. the region [10 , 40] m/s. With K increasing, b oth of the t w o b enchmarks are shown to provide signicantly dete- riorated tracking and signal detection p erformance and gradually fail to accurately track the lo cations of users. In contrast, the prop osed HVMP algorithm suers from m uc h milder tracking p erformance degradation and can pro vide consistent cen timeter-lev el tracking p erformance when K ≤ 10 . These results verify the robust tracking p erformance of the prop osed algorithm in the scenario with a wide range of user num b ers, implying that the HVMP algorithm is suitable to the practical scenario with a v arying num ber of users. Compared with the tracking p erformance, the prop osed algorithm suers from larger estimation error on the transmit signals due to more sev ere interference b etw een users, indicating the necessity of increasing the trac king o verhead to enhance b oth of the c hannel estimation and signal detection p erformance. The impact of the LOS blo ckage for the user-BS links on the p erformance of the prop osed algorithm is shown in Fig. 9. When the LOS blockage b ecomes more severe with a larger probability of LOS blo c kage, the m ulti- user tracking and signal detection p erformance of the prop osed algorithm and Benchmark 2 is signicantly degraded due to a smaller av erage num b er of p osition- related measurements in eac h tracking time slot. Ho w ev er, Benc hmark 1 has slightly degraded p erformance when the probabilit y of LOS blo ckage increases, since we observe that p osition-related measurements of user-BS links ob- tained b y MUSIC are still unreliable and pro vide only mild renemen t of the localization performance in the case with inadequate snapshots. Therefore, the prop osed algorithm can also mak e b etter use of the multi-path information o v er existing schemes enhance the tracking p erformance with even millimeter-level lo cation precision. Finally , we give an example of the tracking tra jectory deriv ed b y the prop osed algorithm and benchmarks in Fig. 10. It is observ ed that track ed trajectories of the prop osed algorithm keep approaching the true tra jectories, while the tw o b enchmarks usually uctuate during the tracking pro cedure. Moreov er, track ed tra jectories of b enchmarks ha v e muc h larger deviations at the intersection region of the tra jectories of dierent users, due to their severe 20 25 30 35 40 45 50 55 60 65 70 -30 -25 -20 -15 -10 -5 0 5 10 15 20 33 35 -15 -13 Fig. 10. The track ed tra jectory under the transmit p ow er P = 10 dBm. 0 20 40 60 80 100 120 10 -2 10 -1 10 0 10 1 0 20 40 60 80 100 120 10 -5 10 0 Fig. 11. The tracking RMSE in each time slot under the transmit power P = 10 dBm. m utual interferences. The tracking and signal detection p erformance across dieren t time slots during the trac king pro cedure is also ev aluated in Fig. 11. Here, we consider the setup where LOS links are a v ailable in the tracking time slot t ∈ [21 , 60] . It is observed that the tracking p erformance is sligh tly impro v ed with extra LOS measure- men ts, while the signal detection p erformance is dramati- cally enhanced b y exploiting the multi-path diversit y . The ab o v e results also v alidate that the prop osed algorithm can eectively combine the measurements from LOS and reecting links to bo ost the joint estimation p erformance, making it app ealing to the practical system. The future w ork can extend this system to the co ded system by exploiting channel co ding to b o ost the p erformance. VI I. Conclusions This w ork studied the joint multi-user tracking and signal detection problem in the RIS-assisted ISA C sys- tems. W e rst established a comprehensive probabilistic mo del for the received signals during the trac king pro cess and then prop osed a nov el HVMP algorithm with com- putational eciency to simultaneously track user states (p osition and velocity) and detect the transmit signals 13 in an online manner. W e further analyzed the tracking p erformance via the p erformance limit of BCRB and prop osed to optimize the RIS phase prole b y minimizing the weigh ted BCRB. Simulation results were pro vided to v alidate that the prop osed algorithm can realize high- accuracy tracking and signal detection p erformance for users with diverse and high mobility . App endix A Deriv ation of the Message (27) The message in (27) can b e detailed expressed by (60), whose explicit expression is hard to obtain. By further utilizing the Gaussian appro ximation, the message can b e obtained as (27), where the mean ˙ ← − p k,t ( j ) and co v ariance ˙ ← − Q k,t ( j ) are obtained by the Gauss-Newton metho d. By denoting the exp onential term in (60) as F 1 ( ˙ p k,t ) , the gradien t and the Hessian matrix calculated in eac h iteration as (61) and (62), where ← − e UB k,g ,t ( j ) = ( ← − d UB k,g ,t ( j )) − 1 ( ← − p UB k,g ,t ( j ) − p B g ) is the up dated directional v ector for the user k -BS g link in the j th iteration of HVMP and ← − d UB k,g ,t ( j ) is the up dated distance for the user k -BS g link with ← − p UB k,g ,t ( j ) up dated in (38). Then, the up dated v alue of ˙ ← − p k,t ( j, ι vel ) and ˙ ← − Q k,t ( j, ι vel ) in eac h iteration can b e given by ˙ ← − p k,t ( j, new ) = ˙ ← − p k,t ( j, old ) −  ▽ ▽ ▽ 2 F 1 ( ˙ ← − p k,t ( j, old ))  − 1 ▽ ▽ ▽ F 1 ( ˙ ← − p k,t ( j, old )) , (63) ˙ ← − Q k,t ( j, new ) = −  ▽ ▽ ▽ 2 F 1 ( ˙ ← − p k,t ( j, new )  − 1 . (64) After the Gauss-Newton metho d con v erges, we can nally obtain the mean and cov ariance in the Gaussian message (27). References [1] A. Liu et al., “A survey on fundamental limits of integrated sensing and communication,” IEEE Commun. Surveys T uts., vol. 24, no. 2, pp. 994–1034, 2nd Quat. 2022. [2] F. Liu, C. Masouros, A. P . Petropulu, H. Griths, and L. Hanzo, “Joint radar and comm unication design: Applications, state-of- the-art, and the road ahead,” IEEE T rans. Commun., vol. 68, no. 6, pp. 3834–3862, Jun. 2020. [3] L. Zheng, M. Lops, Y. C. Eldar, and X. W ang, “Radar and communication co existence: An ov erview: A review of recent methods,” IEEE Signal Pro cess. Mag., vol. 36, no. 5, pp. 85–99, Sep. 2019. [4] L. Liu, S. Zhang, and S. Cui, “Lev eraging a v ariet y of anc hors in cellular netw ork for ubiquitous sensing,” IEEE Commun. Mag., vol. 62, no. 9, pp. 98–104, Sep. 2024. [5] L. Dong, “Coop erative localization and tracking of mobile ad hoc networks,” IEEE T rans. Signal Pro cess., vol. 60, no. 7, pp. 3907–3913, Jul. 2012. [6] W. Y uan, Z. W ei, S. Li, J. Y uan, and D. W. K. Ng, “Integrated sensing and comm unication-assisted orthogonal time frequency space transmission for vehicular netw orks,” IEEE J. Sel. T opics Signal Pro cess., vol. 15, no. 6, pp. 1515–1528, Nov. 2021. [7] W. Zh u, S. Zhang, and L. Liu, “Joint transmission and compres- sion optimization for netw orked sensing with limited-capacity fronthaul links,” IEEE T rans. Wireless Commun., v ol. 24, no. 8, pp. 6643–6657, A ug. 2025. [8] W. T ang et al., “Wireless comm unications with recongurable intelligen t surface: Path loss mo deling and exp erimental mea- surement,” IEEE T rans. Wireless Commun., vol. 20, no. 1, pp. 421–439, Jan. 2021. [9] Q. W u, S. Zhang, B. Zheng, C. Y ou, and R. Zhang, “In telligent reecting surface-aided wireless communications: A tutorial,” IEEE T rans. Commun., vol. 69, no. 5, pp. 3313–3351, Ma y 2021. [10] Y. Liu et al., “Recongurable in telligent surfaces: Principles and opportunities,” IEEE Commun. Surv eys T uts., vol. 23, no. 3, pp. 1546–1577, 3rd Quart. 2021. [11] W. Zhu, J. Gao, S. Zhang, and L. Liu, “Recongurable intelli- gent surface-assisted multiuser tracking and signal detection in ISAC,” in Pro c. Wireless Commun. Signal Pro cess., Oct. 2025, pp. 1–6. [12] W. Zhu, Q. W ang, S. Zhang, B. Di, L. Liu, and Y. C. Elada, “Intelligen t reecting surface based localization of mixed near- eld and far-eld targets,” IEEE T rans. Wireless Commun., pp. 1–1, 2026. [13] S. P . Chepuri, N. Shlezinger, F. Liu, G. C. Alexandrop oulos, S. Buzzi, and Y. C. Eldar, “In tegrated sensing and commu- nications with recongurable intelligen t surfaces: F rom signal modeling to pro cessing,” IEEE Signal Pro cess. Mag., vol. 40, no. 6, pp. 41–62, Sept. 2023. [14] R. Liu, M. Li, Y. Liu, Q. W u, and Q. Liu, “Joint transmit wa veform and passiv e b eamforming design for RIS-aided DFRC systems,” IEEE J. Sel. T opics Signal Pro cess., vol. 16, no. 5, pp. 995–1010, Aug. 2022. [15] A. Elzanaty , A. Guerra, F. Guidi, and M.-S. Alouini, “Re- congurable intelligen t surfaces for lo calization: Position and orientation error b ounds,” IEEE T rans. Signal Pro cess., vol. 69, pp. 5386–5402, 2021. [16] X. Song, J. Xu, F. Liu, T. X. Han, and Y. C. Eldar, “Intelligen t reecting surface enabled sensing: Cramér-rao b ound optimiza- tion,” IEEE T rans. Signal Pro cess., v ol. 71, pp. 2011–2026, 2023. [17] X. Zhang, L. Liu, S. Zhang, W. Zhu, and H. Zhang, “Optimiz- ing rate-CRB p erformance for b eyond diagonal recongurable intelligen t surface enabled ISAC,” IEEE Commun. Lett., pp. 1–1, 2025. [18] H. Zhang, H. Zhang, B. Di, K. Bian, Z. Han, and L. Song, “Met- alocalization: Recongurable in telligent surface aided multi-user wireless indoor lo calization,” IEEE T rans. Wireless Commun., vol. 20, no. 12, pp. 7743–7757, Dec. 2021. [19] Q. W ang, L. Liu, S. Zhang, B. Di, and F. C. M. Lau, “A heterogeneous 6G netw orked sensing architecture with active and passiv e anchors,” IEEE T rans. Wireless Comm un., vol. 23, no. 8, pp. 9502–9517, Aug. 2024. [20] Y. Han, S. Jin, C.-K. W en, and T. Q. S. Quek, “Lo calization and channel reconstruction for extra large RIS-assisted massive MIMO systems,” IEEE J. Sel. T opics Signal Process., v ol. 16, no. 5, pp. 1011–1025, Aug. 2022. [21] M. Rahal, B. Denis, M. F. Keskin, B. Uguen, and H. W ymeersch, “Ris-enabled nlos near-eld join t p osition and v elocity estima- tion under user mobilit y ,” IEEE J. Sel. T opics Signal Process., vol. 18, no. 4, pp. 633–645, May 2024. [22] M. Ammous and S. V alaee, “Positioning and trac king using recongurable intelligen t surfaces and extended kalman lter,” in IEEE V eh. T ech. Conf., Jun. 2022, pp. 1–6. [23] P . Zheng, X. Liu, and T. Y. Al-Naouri, “LEO- and RIS- empowered user tracking: A riemannian manifold approach,” IEEE J. Sel. Areas Commun., vol. 42, no. 12, pp. 3445–3461, Dec. 2024. [24] G. Einic k e and L. White, “Robust extended kalman ltering,” IEEE T rans. Signal Process., vol. 47, no. 9, pp. 2596–2599, Sept. 1999. [25] B. T eng, X. Y uan, and R. W ang, “V ariational ba yesian multiuser tracking for recongurable intelligen t surface-aided MIMO- OFDM systems,” IEEE J. Sel. Areas Commun., vol. 41, no. 12, pp. 3752–3767, Dec. 2023. [26] P . Tic havsky , C. Mura v chik, and A. Nehorai, “P osterior Cramér- Rao b ounds for discrete-time nonlinear ltering,” IEEE T rans. Signal Pro cess., vol. 46, no. 5, pp. 138 6–1396, May 1998. [27] P . Djuric and others., “P article ltering,” IEEE Signal Pro cess. Mag., vol. 20, no. 5, pp. 19–38, Sep. 2003. [28] Y. Bar-Shalom, F. Daum, and J. Huang, “The probabilistic data association lter,” IEEE Control Syst. Mag., v ol. 29, no. 6, pp. 82–100, Dec. 2009. [29] M. R. Akdeniz et al., “Millimeter wa ve channel modeling and cellular capacity ev aluation,” IEEE J. Sel. Areas Commun., vol. 32, no. 6, pp. 1164–1179, Jun. 2014. 14 ε ˙ p k,t ← ˙ f k,t ( j ) ∝  G  g =1 δ  ν UB k,g ,t − λ − 1 ζ T ( ˙ p U k,t ) T e UB k,g ,t  V M ( − ζ T ν UB k,t ; ← − µ UB k,g ,t ( ν, j ) , ← − κ UB k,g ,t ( ν, j )) × R  r =1 δ  ν UI k,r,t − λ − 1 ζ T ( ˙ p U k,t ) T e UI k,r,t  V M ( − ζ T ν UI k,r,t ; ← − µ UI k,r,t ( ν, j ) , ← − κ UI k,r,t ( ν, j )) , ∝ G  g =1 V M  − λ − 1 ζ T ( ˙ p U k,t ) T e UB k,g ,t ; ← − µ UB k,g ,t ( ν, j ) , ← − κ UB k,g ,t ( ν, j )  R  r =1 V M  − λ − 1 ζ T ( ˙ p U k,t ) T e UI k,r,t ; ← − µ UI k,r,t ( ν, j ) , ← − κ UI k,t,t ( ν, j )  , (60) ▽ ▽ ▽ F 1 ( ˙ p k,t ) = − G  g =1 ˘ ζ T ← − κ UB k,g ,t ( ν, j ) λ ← − e UB k,g ,t ( j ) sin  ˘ ζ T λ ( ˙ p U k,t ) T ← − e UB k,g ,t ( j ) + ← − µ UB k,g ,t ( ν, j )  − R  r =1 ˘ ζ T ← − κ UI k,r,t ( ν, j ) λ ← − e UI k,r,t ( j ) sin  ˘ ζ T λ ( ˙ p U k,t ) T ← − e UI k,r,t ( j ) + ← − µ UI k,t ( ν, j )  , (61) ▽ ▽ ▽ 2 F 1 ( ˙ p k,t ) = − G  g =1 ˘ ζ 2 T ← − κ UB k,g ,t ( ν, j ) λ 2 ← − e UB k,g ,t ( j )( ← − e UB k,g ,t ( j )) T cos  ˘ ζ T λ ( ˙ p U k,t ) T ← − e UB k,g ,t ( j ) + ← − µ UB k,g ,t ( ν, j )  − G  g =1 ˘ ζ 2 T ← − κ UI k,r,t ( ν, j ) λ 2 ← − e UI k,r,t ( j )( ← − e UI k,r,t ( j )) T cos  ˘ ζ T λ ( ˙ p U k,t ) T ← − e UI k,r,t ( j ) + ← − µ UI k,r,t ( ν, j )  , (62) [30] C. Han, A. O. Bicen, and I. F. Akyildiz, “Multi-ray channel modeling and wideband c haracterization for wireless communi- cations in the terahertz band,” IEEE T rans. Wireless Commun., vol. 14, no. 5, pp. 2402–2412, May 2015. [31] H. Han, T. Jiang, and W. Y u, “Activ e sensing for multiuser beam tracking with recongurable intelligen t surface,” IEEE T rans. Wireless Comm un., v ol. 24, no. 1, pp. 540–554, Jan. 2025. [32] E. Nassa ji, M. Bashir, and D. T ruhachev, “Unsourced random access ov er fading channels via data rep etition, p ermutation, and scrambling,” IEEE T rans. Commun., vol. 70, no. 2, pp. 1029–1042, F eb. 2022. [33] A. V em, K. R. Naray anan, J.-F. Chamberland, and J. Cheng, “A user-independent successive in terference cancellation based coding scheme for the unsourced random access gaussian c han- nel,” IEEE T rans. Commun., vol. 67, no. 12, pp. 8258–8272, Dec. 2019. [34] M.-A. Badiu, T. L. Hansen, and B. H. Fleury , “V ariational bay esian inference of line spectra,” IEEE T rans. Signal Process., vol. 65, no. 9, pp. 2247–2261, May 2017. [35] J. Gao et al., “Integrated massive communication and target localization in 6g cell-free netw orks,” IEEE T rans. Wireless Commun., pp. 1–1, 2025. [36] Q. Zhang, J. Zhu, N. Zhang, and Z. Xu, “Multidimensional v ariational line sp ectra estimation,” IEEE Signal Pro cess. Lett., vol. 27, pp. 945–949, Ma y 2020. [37] W. Zhu, M. T ao, X. Y uan, F. Xu, and Y. Guan, “Co op erative multi-cell massive access with temporally correlated activity ,” IEEE T rans. Wireless Comm un., v ol. 23, no. 5, pp. 4347–4363, May 2024. [38] D. Guo, S. Shamai, and S. V erdu, “Mutual information and minimum mean-square error in gaussian channels,” IEEE T rans. Infor. Theory , vol. 51, no. 4, pp. 1261–1282, Apr. 2005. [39] T.-J. Shan, M. W ax, and T. Kailath, “On spatial smo othing for direction-of-arriv al estimation of coheren t signals,” IEEE T rans. Acoustics, Speech, Signal Process., v ol. 33, no. 4, pp. 806–811, Aug. 1985.