Paper deep dive
Adaptive Distributed Physical-Layer Authentication and Attack Detection in 6G Non-Terrestrial Networks via Causal Meta-Learning
Parsa Rajabi, Mohammad Reza Abedi, Nader Mokari, Paeiz Azmi, Halim Yanikomeroglu
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 92%
Last extracted: 9/10/2026, 2:53:19 AM
Summary
This paper proposes SAFA-MZ, a secure adaptive framework for distributed physical-layer authentication (DPLA) in 6G non-terrestrial networks (NTNs). It addresses challenges like Doppler shifts and channel variations by using a multi-feature fingerprint (spatial, angular, Doppler-delay, etc.) and a structural causal model (SCM) to capture environmental relations. The framework employs Model-Agnostic Meta-Learning (MAML) with Invariant Risk Minimization (IRM) for fast adaptation to unseen environments. A two-stage authentication scheme uses local recognition first, activating Time-Difference-of-Arrival (TDOA) localization via a Graph Attention Network (GAT) only when necessary, reducing backhaul overhead. Simulations show 92% accuracy and 96% AUC.
Entities (13)
Relation Signals (12)
SAFA-MZ → appliesto → Non-Terrestrial Networks
confidence 95% · This paper proposes a secure adaptive framework for authentication in multi-zone networks (SAFA-MZ), a causal meta-learning framework for distributed PLA (DPLA) in NTNs.
SAFA-MZ → uses → Structural Causal Model
confidence 95% · Second, we formulate a structural causal model (SCM) to capture the relations among design choices, environmental factors, extracted features, and authentication outcomes.
SAFA-MZ → uses → MAML
confidence 95% · Third, we develop a model-agnostic meta-learning (MAML) strategy with invariant risk minimization (IRM) and causal consistency regularization...
SAFA-MZ → uses → IRM
confidence 95% · Third, we develop a model-agnostic meta-learning (MAML) strategy with invariant risk minimization (IRM)...
SAFA-MZ → uses → Multi-Feature Fingerprint
confidence 95% · First, we design a multi-feature fingerprint that combines spatial, angular, combiner, subspace, and Doppler-delay features. The fingerprint is adaptive and distributed...
Multi-Feature Fingerprint → combines → Angular Features
confidence 90% · First, we design a multi-feature fingerprint that combines spatial, angular, combiner, subspace, and Doppler-delay features.
Multi-Feature Fingerprint → combines → Spatial Features
confidence 90% · First, we design a multi-feature fingerprint that combines spatial, angular, combiner, subspace, and Doppler-delay features.
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Physical-layer authentication (PLA) in non-terrestrial networks (NTNs) is challenged by severe Doppler shifts, long delays, and fast channel variations, which cause distribution shifts and degrade conventional learning methods. Existing PLA schemes often rely on single features or generalize poorly to unseen environments. This paper proposes a secure adaptive framework for authentication in multi-zone networks (SAFA-MZ), a causal meta-learning framework for distributed PLA (DPLA) in NTNs. First, we design a multi-feature fingerprint that combines spatial, angular, combiner, subspace, and Doppler-delay features. The fingerprint is adaptive and distributed, as it fuses heterogeneous physical-layer features and measurements from multiple aerial nodes. Second, we formulate a structural causal model (SCM) to capture the relations among design choices, environmental factors, extracted features, and authentication outcomes. Third, we develop a model-agnostic meta-learning (MAML) strategy with invariant risk minimization (IRM) and causal consistency regularization for fast adaptation to unseen NTN environments with few labeled samples. Fourth, we propose a two-stage authentication scheme that performs local recognition and activates time-difference-of-arrival (TDOA) localization with a graph attention (GAT) network only when needed, which reduces backhaul overhead. Simulations show that SAFA-MZ achieves 92% accuracy and 96% AUC, outperforming centralized deep learning and single-feature baselines across diverse environments.
Tags
Links
- Source: https://arxiv.org/abs/2609.09511v1
- Canonical: https://arxiv.org/abs/2609.09511v1
Trouble viewing inline? Open PDF directly →
Full Text
75,548 characters extracted from source content.
Expand or collapse full text
1 Adaptive Distributed Physical-Layer Authentication and Attack Detection in 6G Non-Terrestrial Networks via Causal Meta-Learning Parsa Rajabi , Mohammad Reza AbediStudent Member, IEEE, Nader MokariSenior Member, IEEE, Paeiz AzmiSenior Member, IEEE, Halim YanikomerogluFellow, IEEE Abstract—Physical-layerauthentication(PLA)innon- terrestrial networks (NTNs) is challenged by severe Doppler shifts, long delays, and fast channel variations, which cause distribution shifts and degrade conventional learning methods. Existing PLA schemes often rely on single features or generalize poorly to unseen environments. This paper proposes a secure adaptive framework for authentication in multi-zone networks (SAFA-MZ), a causal meta-learning framework for distributed PLA (DPLA) in NTNs. First, we design a multi-feature fingerprint that combines spatial, angular, combiner, subspace, and Doppler-delay features. The fingerprint is adaptive and distributed, as it fuses heterogeneous physical-layer features and measurements from multiple aerial nodes. Second, we formulate a structural causal model (SCM) to capture the relations among design choices, environmental factors, extracted features, and authentication outcomes. Third, we develop a model-agnostic meta-learning (MAML) strategy with invariant risk minimization (IRM) and causal consistency regularization for fast adaptation to unseen NTN environments with few labeled samples. Fourth, we propose a two-stage authentication scheme that performs local recognition and activates time-difference- of-arrival (TDOA) localization with a graph attention (GAT) network only when needed, which reduces backhaul overhead. Simulations show that SAFA-MZ achieves 92% accuracy and 96% AUC, outperforming centralized deep learning and single-feature baselines across diverse environments. IndexTerms—Distributedphysical-layerauthentication (DPLA), non-terrestrial networks (NTN), causal inference, meta-learning, MIMO, graph attention network. I. INTRODUCTION P HYSICAL-layer authentication (PLA) verifies transmit- ter identity using channel characteristics and hard- ware impairments. Machine learning methods have been widely used for PLA [1], including channel prediction-based, classification-based, and similarity-based schemes [2]. Com- pared with threshold-based hypothesis testing [3]–[6], artificial intelligence (AI)-based authentication avoids fixed thresholds and improves robustness in dynamic environments [7], [8]. This work has been submitted to the IEEE for possible publication. Copyright may be transferred without notice, after which this version may no longer be accessible. P. Rajabi, MR. Abedi, N. Mokari, and P. Azmi are with the Department of Electrical and Computer Engineering, Tarbiat Modares University, Tehran, 14115-111, Iran (e-mail: parsa_rajabi@modares.ac.ir; Mohammadreza_abedi@modares.ac.ir;nader.mokari@modares.ac.ir; pazmi@modares.ac.ir). H. Yanikomeroglu is with the Department of Systems and Com- puter Engineering, Carleton University, Ottawa, K1S 5B6, Canada (e-mail: halim@sce.carleton.ca). Existing PLA methods rely on single features, tag- embedding [5], or hybrid schemes [9]. In contrast, this paper uses multiple features for dynamic non-terrestrial networks (NTNs). Security in NTN is challenged by time-varying chan- nels and hardware impairments [10], which degrade authen- tication. To address spoofing in such settings, we propose a multi-feature distributed PLA (DPLA) framework. DPLA architectures, including semi-distributed, decision-level fusion, and fully distributed designs, offer different tradeoffs among accuracy, overhead, and robustness [11]. Prior works have studied DPLA for millimeter-wave multiple-input multiple- output (mmWave MIMO) systems using beam-pattern devi- ations [12], dynamic weighted voting in smart grids [13], and missed-detection bounds for ultra-reliable low-latency communications (URLLC) [14]. Deep learning-based fingerprint identification can be vulner- able to adversarial perturbations [15] and cross-receiver distri- bution shifts caused by heterogeneous receiver hardware [16]. To improve robustness under distribution shifts and adversarial variations, we combine causal learning and meta-learning. Causal learning extracts invariant relations and reduces spuri- ous correlations, while meta-learning enables fast adaptation to unseen environments with few labels. Related studies have used structural causal models (SCMs) for radio-frequency (RF) fingerprint purification [17], invariant graph learning [18], counterfactual fairness [19], reinforcement learning (RL) [20]–[22], 6G explainability [23], wireless sustainability anal- ysis [24], and semantic communication [25]. Meta-learning has also shown effectiveness in internet of things (IoT)-based multi-user PLA [26] and low-Earth orbit (LEO) satellites [27], and sparse fingerprint-based localization using graph-based spatial aggregation [28]. Pilot contamination in multi-user uplinks degrades chan- nel estimates and physical features when users share pilots. Similar problems appear in cell-free massive MIMO [29]. Hence, we use a normalized minimum mean-square error (MMSE) combiner, interference-aware features, and causal meta-learning. Related attacks on channel state information (CSI)-based PLA are studied in [30], while pilot contamination in mmWave grant-free massive machine-type communications is addressed in [31]. This paper proposes a secure adaptive framework for au- thentication in multi-zone networks (SAFA-MZ). The frame- work combines SCM, backdoor adjustment, invariant risk min- imization (IRM), and model-agnostic meta-learning (MAML) arXiv:2609.09511v1 [eess.SP] 8 Sep 2026 2 to learn environment-invariant authentication features and adapt rapidly to unseen NTN conditions. It also estimates the conditional average treatment effects (CATEs) of key design variables. For efficiency, local multi-feature fingerprint recognition is performed first, and time-difference-of-arrival (TDOA)-based distributed verification is activated only when local authentication is unreliable. In the simulations, SAFA- MZ is compared with representative baselines, including a node-level deep learning method without causal meta-learning or distributed verification, and single-feature authentication schemes that rely on individual physical-layer features and do not employ distributed verification. Finally, physically grounded causal interpretations are provided to explain the influence of environmental factors and design variables on the authentication performance. The contributions of this paper are summarized as follows: • We propose a dual-mode DPLA framework for NTN that uses local multi-feature authentication and activates TDOA-based attack localization only when needed. This reduces overhead and improves robustness and security. • We develop a causal-invariant fingerprinting method that uses multiple features under pilot contamination, includ- ing directional similarity, beam concentration, combiner consistency, subspace deviation, Doppler-delay character- istics, SINR, and interference-to-noise ratio (INR). This improves robustness against environmental changes and spurious correlations in dynamic NTN channels. • We integrate causal representation learning and meta- learning for fast adaptation to new environments with minimal samples. Additionally, we propose a graph-based scheme for coordinated attack detection among aerial nodes and provide causal interpretations for authentica- tion performance. NOTATION We use bold lowercase letters (e.g., y) for vectors. The Euclidean norm is ∥ · ∥, and the complex modulus is | · |. Expectation and probability are denoted by E[·] and P(·), respectively. 1· is the indicator function. The symbol ◦ denotes the Hadamard product. The M × M unitary discrete Fourier transform (DFT) matrix is F M . It transforms the signal from the antenna domain to the angular (beam) domain. I. SYSTEM MODEL As illustrated in Fig. 1, we consider a set K≜ k 1 ,k 2 ,...,k |K| of aerial nodes (e.g., high-altitude platform stations (HAPS), LEO, or geostationary orbit (GEO) satel- lites). Each node k is equipped with M k antennas. The nodes are spatially distributed and connected via a backhaul link to a central fusion center (FC). A set U≜ u 1 ,u 2 ,...,u |U| of single-antenna users, including legitimate users and ille- gitimate transmitters (attackers), transmits uplink signals. For a given user u, the subset K u ⊆ K denotes the nodes that can receive its signal (i.e., the serving nodes). Due to the relative motion between the nodes and the users, the channel experiences severe Doppler shifts and long propagation delays, Fig. 1: System model with multiple HAPSs/satellites and a fusion center. making DPLA challenging. Therefore, we consider a two- phase authentication framework. In the local authentication phase, legitimate users are registered at the FC using their physical-layer signatures. During uplink transmission, each observation (either legitimate user or attacker) is locally compared with the enrolled profiles. Matched transmissions are accepted, while suspicious ones are forwarded to the distributed verification phase, where multiple nodes cooperate and the FC makes the final decision. A. Uplink Signal Model and Feature Extraction The received uplink snapshot at node k over block n is denoted by y k (n) ∈ C M k and modeled using the following explicit multi-user representation for any target user u: y k (n) = X u∈U k (n) p p u (n) h k,u (n)x u (n) + n k (n) = p p u (n) h k,u (n)x u (n) + X u ′ ∈U k (n) u ′ ̸=u p p u ′ (n) h k,u ′ (n)x u ′ (n) + n k (n), (1) where U k (n) denotes the set of users simultaneously active toward node k in block n, p u (n) is the uplink transmit power of user u, h k,u (n)∈ C M k is the corresponding uplink channel vector, x u (n) ∈ C is the transmitted pilot/data symbol with E[|x u (n)| 2 ] = 1, and n k (n) ∼ CN (0,σ 2 k I M k ) is additive noise (n k (n) ∈ C M k ), where σ 2 k denotes the noise variance per antenna at node k. For each active target user u∈U k (n), node k applies a user-specific linear combiner w k,u (n)∈ C M k to the received vector y k (n). The resulting scalar output is y (out) k,u (n) = w H k,u (n)y k (n).(2) When multiple users are simultaneously active over the same time-frequency resource, the users u̸= u ′ contribute inter-user interference (IUI); if only one user is active, the interference term vanishes. Assuming pilot resources are available, let ˆ h k,u (n)∈ C M k denote the estimated uplink channel of active user u at node k in block n. Let P u (n) denote the set of users reusing the same pilot as user u in block n. Then the pilot-aided channel estimate at node k can be expressed as ˆ h k,u (n) = h k,u (n) + X u ′ ∈P u (n) u ′ ̸=u α u,u ′ h k,u ′ (n) + e k,u (n), (3) 3 where the second term represents pilot contamination caused by other users sharing the same pilot resource, α u,u ′ ∈ C denotes the contamination coefficient, and e k,u (n) collects the residual estimation noise and modeling mismatch. To suppress both noise and IUI, we adopt a normalized linear MMSE-type combiner w k,u (n) = R (k) y (n) −1 ˆ h k,u (n) R (k) y (n) −1 ˆ h k,u (n) .(4) Here, the received auto-correlation matrix is defined as R (k) y (n)≜ E y k (n)y H k (n)|h k,u (n) u∈U k (n) ,(5) where the expectation is taken over user symbols and noise, conditioned on the instantaneous channel vectors and under mutually uncorrelated user symbols, can be expanded as R (k) y (n) =p u (n)h k,u (n)h H k,u (n)(6) + X u ′ ∈U k (n) u ′ ̸=u p u ′ (n)h k,u ′ (n)h H k,u ′ (n) + σ 2 k I M k . Therefore, the auto-correlation matrix explicitly captures the desired signal component, the aggregate IUI covariance, and the thermal noise term. In practice, R (k) y (n) is estimated from N c received samples as ˆ R (k) y (n) = 1 N c N c X i=1 y (i) k (n)y (i)H k (n),(7) where y (i) k denotes the i-th received snapshot at node k. The normalization in (4) ensures ∥w k,u (n)∥ = 1. Compared with maximum-ratio combining (MRC), the MMSE combiner is more robust in dense multi-user scenarios because it accounts for the spatial covariance of both noise and co-channel interference. Moreover, its computation relies on the instantaneous received statistics and the channel estimate of the target stream, without requiring any higher-layer user identity information. The uplink channel is modeled as a time-varying Rician fading process. Since the user velocity is assumed negligible relative to the aerial node motion, the dominant Doppler shift is induced by node k and is approximated by f (k,u) d = v k c f c cos θ (k,u) elev ,(8) where v k is the velocity of node k, c = 3× 10 8 m/s is the speed of light, f c is the carrier frequency, and θ (k,u) elev ∈ [0,π/2] is the elevation angle between node k and user u. Under a narrowband common-Doppler approximation, the uplink channel evolves as h k,u (n)≈ h LoS k,u (n)e j2πf (k,u) d nT s + h NLoS k,u (n),(9) where h LoS k,u (n) ∈ C M k denotes the slowly varying line- of-sight component and h NLoS k,u (n) ∈ C M k represents the residual diffuse component. The dominant deterministic phase rotation is captured through the LoS component, while the residual diffuse component is absorbed into h NLoS k,u (n). This formulation is especially useful in the present setting because multi-user interference is governed not only by user activity and transmit powers, but also by the temporal evolution and spatial overlap of the corresponding channel vectors. During the enrollment phase, node k computes a set of reference statistics for each legitimate user u based on N enr observations. The first quantity is the mean channel vector μ (k,u) h ≜ 1 N enr N enr X n=1 ˆ h k,u (n),(10) whereμ (k,u) h ∈ C M k represents the average spatial signature of legitimate user u at node k, capturing its mean angle-of- arrival structure and average channel gain. The second quantity is the enrollment channel correlation matrix R (k,u) h ≜ 1 N enr N enr X n=1 ˆ h k,u (n) ˆ h H k,u (n),(11) and the corresponding reference subspace projection matrix is P (k,u) ref ≜ U (k,u) ref U (k,u) ref H ,(12) where the columns of U (k,u) ref ∈ C M k ×r k are the r k dominant eigenvectors of R (k,u) h , with r k ≤ M k . The matrix P (k,u) ref cap- tures the dominant spatial subspace occupied by the legitimate user’s channel realizations. In NTN scenarios, this subspace may drift slowly due to node motion and geometry variation, and hence periodic re-enrollment may be required. The reference combiner associated with user u at node k is constructed from the enrollment statistics as w (ref) k,u = R (k) y,ref −1 μ (k,u) h R (k) y,ref −1 μ (k,u) h ,(13) where R (k) y,ref denotes the received-signal covariance matrix estimated from the enrollment observations of user u at node k. This reference combiner is consistent with a normalized linear MMSE-type combiner in (4) and provides a baseline spatial filtering profile for subsequent authentication. For angular-domain analysis, we employ the unitary DFT matrix F M k ∈ C M k ×M k to transform antenna-domain channel vectors into the beam domain. Let N (k,u) main ⊆ 1, 2,...,M k denote the set of dominant DFT bins corresponding to the main lobe of legitimate user u as observed at node k. These bins are identified during enrollment, for example, as those containing a prescribed fraction of the total beam-domain energy. Since the dominant arrival direction depends on the NTN geometry, N (k,u) main may evolve with the elevation angle and should be updated when needed. Moreover, we define the aggregate post- combining IUI power as I k,u (n)≜ X u ′ ∈U k (n) u ′ ̸=u p u ′ (n) w H k,u (n) ˆ h k,u ′ (n) 2 .(14) For each active target user u∈U k (n) and each observation block n, node k extracts the feature vector f k,u (n) = [f k,u,1 (n),...,f k,u,d f (n)] T ,(15) 4 where d f = 8 is the feature dimension considered in this paper. We define the elements of f k,u (n)∈ R d f as follows: f k,u,1 (n) = ˆ h H k,u (n)μ (k,u) h ∥ ˆ h k,u (n)∥μ (k,u) h ∥ .(16) f k,u,1 (n) is the directional similarity between the current estimated channel and the enrolled mean channel. It measures the alignment of the instantaneous spatial signature with the legitimate reference signature. f k,u,2 (n) = P i∈N (k,u) main [F M k ˆ h k,u (n)] i 2 P M k i=1 [F M k ˆ h k,u (n)] i 2 .(17) f k,u,2 (n) is the main-lobe concentration ratio in the beam domain. It quantifies how much of the channel energy remains concentrated inside the enrolled dominant angular support. f k,u,3 (n) = w H k,u (n)w (ref) k,u ∥w k,u (n)∥w (ref) k,u ∥ .(18) f k,u,3 (n) is the combiner similarity between the current nor- malized MMSE combiner and the enrollment-time reference combiner. It reflects the consistency of the spatial filtering direction across time. f k,u,4 (n) = w H k,u (n)y k (n) ∥w k,u (n)∥y k (n)∥ .(19) f k,u,4 (n) is a normalized projection score obtained from the combiner output magnitude relative to the received-snapshot norm. It captures the instantaneous compatibility between the combiner and the received signal. f k,u,5 (n) = 1− ∥P (k,u) ref ˆ h k,u (n)∥ 2 ∥ ˆ h k,u (n)∥ 2 .(20) f k,u,5 (n) is the subspace deviation, namely the fraction of channel energy lying outside the enrolled legitimate-user subspace. Larger values indicate stronger deviation from the enrolled spatial structure. f k,u,6 (n) = exp −α (k,u) d ˆτ k,u (n)− τ (ref) k,u (n) 1 + ˆ f (k,u) d (n)−f (k,u) d,ref (n) σ (k,u) f 2 .(21) f k,u,6 (n) is the Doppler-delay consistency score, which jointly measures the agreement of the observed transmission with the enrolled geometric and kinematic signature of user u. Here, ˆ f (k,u) d (n) and ˆτ k,u (n) denote the estimated Doppler shift and propagation delay, respectively, while f (k,u) d,ref (n) and τ (ref) k,u (n) are their expected reference values. Moreover, σ (k,u) f > 0 is a scaling factor for the Doppler mismatch, and α (k,u) d > 0 is a decay constant that controls how quickly the score drops with increasing delay error. f k,u,7 (n) = p u (n) w H k,u (n) ˆ h k,u (n) 2 I k,u (n) + σ 2 k ∥w k,u (n)∥ 2 .(22) f k,u,7 (n) is the estimated post-combining SINR, measuring the target user’s separability from co-channel interferers. f k,u,8 (n) = I k,u (n) σ 2 k ∥w k,u (n)∥ 2 .(23) f k,u,8 (n) denotes the post-combining INR, which measures the aggregate IUI relative to the noise floor and captures inter- ference variations, particularly in dense NTN access scenarios. Since the normalized MMSE combiner in (4) satisfies ∥w k,u (n)∥ = 1, the terms ∥w k,u (n)∥ 2 in (22) and (23) can be omitted. The first six features measure spatial, subspace, and geometry consistency with the enrolled legitimate profile, while the last two capture interference. Thus, the feature vector is interference-aware and better suited to dense multi-user NTN uplinks than channel-similarity-only features. B. Authentication 1) Local phase: At coherence block n, each node k independently tests observed feature vector f k,u (n) against enrolled user identity u ′ for local authentication. The local soft score, interpreted as a legitimacy score with respect to user u ′ , is computed from the Mahalanobis distance as d M f k,u (n);u ′ = q f k,u (n)−μ (k,u ′ ) f T Σ (k,u ′ ) f −1 f k,u (n)−μ (k,u ′ ) f , (24) whereμ (k,u ′ ) f ∈ R d f and Σ (k,u ′ ) f ∈ R d f ×d f denote the mean vector and covariance matrix of the legitimate feature distribution for enrolled user u ′ at node k, respectively. During enrollment, N enr legitimate feature snapshots are collected for each enrolled user u ′ at node k. Based on these enrollment samples, the sample mean and covariance are estimated as ˆ μ (k,u ′ ) f = 1 N enr N enr X i=1 f k,u ′ (i),(25) ˆ Σ (k,u ′ ) f = 1 N enr − 1 N enr X i=1 f k,u ′ (i)− ˆ μ (k,u ′ ) f f k,u ′ (i)− ˆ μ (k,u ′ ) f T . (26) In practice, d M (f k,u (n);u ′ ) is evaluated by replacingμ (k,u ′ ) f and Σ (k,u ′ ) f with their sample estimates in (25)–(26). The Mahalanobis distance normalizes the feature space using the spread and correlation structure of legitimate sam- ples. Hence, it is scale-invariant and accounts for inter-feature correlations. A small value of d M f k,u (n);u ′ indicates strong consistency with the legitimate profile of user u ′ , whereas a large value suggests a potential attacker or identity mismatch. This distance is then mapped to a local soft score defined as p k,u,u ′ (n)≜ 1 1 + e d M (f k,u (n);u ′ )−τ k ,(27) 5 where τ k > 0 is the local decision threshold at node k. Moreover, p k,u,u ′ (n) ≈ 1 when d M f k,u (n);u ′ ≪ τ k , corresponding to a close match with the legitimate profile of user u ′ , whereas p k,u,u ′ (n)≈ 0 when d M f k,u (n);u ′ ≫ τ k , indicating an attacker or an identity mismatch. 2) Distributed Verification and Attack Detection: The pro- posed framework operates in two modes. In the local authenti- cation mode, each node independently compares the observed feature vector against the enrolled legitimate-user profiles. We define U enr as the set of enrolled users. For an observation f k,u (n) at node k, let u ⋆ k (n) = arg min u ′ ∈U enr d M f k,u (n);u ′ ,(28) denote the best-matching enrolled identity, and define the corresponding minimum distance as d k,min (n) = min u ′ ∈U enr d M f k,u (n);u ′ .(29) Let τ ac > 0 denote the local acceptance threshold. If d k,min (n) ≤ τ ac , the observation at node k is considered locally consistent with the enrolled profile of user u ⋆ k (n). However, if there exists at least one participating node k ∈K u such that d k,min (n) > τ ac , the transmission is flagged as suspicious and the system switches to the adaptive distributed verification mode. In the distributed verification (attack detection) mode, multi- ple nodes collaboratively contribute to the final authentication decision. The FC collects from each participating node k ∈K u the local soft score p k,u,u ⋆ k (n), the associated best-matching identity u ⋆ k (n), and the time-of-arrival (TOA) estimate ˆτ k,u (n). The TOA measurement is modeled as ˆτ k,u (n) = ∥r u (n)− r k (n)∥ c + ε k,u (n),(30) where r u (n) ∈ R 3 is the unknown transmitter position, r k (n)∈ R 3 is the known position of node k, c is the speed of light, and ε k,u (n) denotes the TOA estimation error. A model for ε k,u (n) is zero-mean Gaussian noise with variance σ 2 toa . We assume that the aerial nodes are time-synchronized through the backhaul network. Using the TOA measurements collected from the nodes in K u , the FC forms TDOA observations and computes a position estimate ˆ r u (n) for the transmitter. Based on this estimate, two location-consistency scores are constructed. The first score quantifies the stability of the estimated position over time: f tri,u (n) = exp − Tr d Cov( ˆ r u (n)) 2σ 2 pos ! ,(31) where d Cov( ˆ r u (n)) denotes the empirical covariance of recent position estimates and σ 2 pos is a scaling parameter. A small value of σ pos makes the score more sensitive to temporal fluctuations in the estimated position. The second score measures the consistency between the estimated transmitter location and the reference location as- sociated with the claimed or best-matching identity: f dist,u,u ⋆ (n) = exp − ∥ ˆ r u (n)− r (ref) u ⋆ (n)∥ 2 2σ 2 dist ! ,(32) where r (ref) u ⋆ (n) is the reference or claimed location associated with identity u ⋆ , and σ 2 dist controls the tolerance to spatial mismatch. A smaller σ dist imposes a stricter consistency requirement between the estimated and claimed positions. At the FC, the final decision is obtained through a scal- able graph-based fusion mechanism rather than direct feature concatenation. Specifically, among all participating nodes in K u , the FC first identifies the L most reliable nodes based on their local verification evidence. It then constructs a graph over the selected nodes using delay, Doppler, and geometric consistency information, performs graph message passing, and generates a graph-level representation through attention- based readout. Finally, the global triangulation-consistency and identity-location consistency scores are appended to the graph- level representation before the final output layer. The detailed formulation is provided in Section I-D. C. Environment and Design Variables for NTN The environmental and link-level context for user u ob- served at node k is represented by Z k,u ≜ γ k,u , M k , σ θ,k , K (k,u) r , H (k) alt , θ (k,u) elev , v k . (33) Here, γ k,u denotes the received SNR of user u at node k, expressed in dB; σ θ,k is the angular spread observed at node k, which is typically small in NTN links dominated by strong LoS components; K (k,u) r is the Rician factor of the link between user u and node k; H (k) alt ∈ R + is the altitude of node k, measured in km. The variables in Z k,u characterize the operating environment and the link geometry. In particular, low elevation angles generally lead to longer propagation distances and may increase path loss and Doppler sensitivity, while higher node velocities increase the Doppler shift according to (8). Moreover, the Rician factor, angular spread, pilot reuse factor, and analog-to-digital converter (ADC) resolution affect the quality of the extracted physical-layer features and, consequently, the authentication performance. We distinguish these environmental variables from control- lable design variables. Let T ∈ T denote the treatment, i.e., the controllable system configuration selected by the designer, where T is the set of admissible treatment values. Depending on the experiment, T may represent a scalar design choice, such as the ADC resolution, T = N (k) bit , the pilot reuse factor, T = K (k) p , the top-L selection budget, T = L, or the node orbit/platform type, e.g., HAPS, LEO, or GEO. More generally, T may be a vector of design parameters T = N (k) bit , K (k) p , L, orbit (k) T .(34) Each realization of T corresponds to a particular system configuration whose effect on the authentication performance metrics, such as TPR, FPR, AUC, or detection accuracy, can be evaluated. Unlike the environmental variables Z k,u , which describe the propagation, mobility, and hardware con- text experienced by the link, the treatment variables represent design choices that can be actively controlled by the system designer. Depending on the considered NTN deployment, such treatments can be configured either globally at the network level or locally on a per-node basis. 6 I. PROPOSED METHOD Our method has three main goals in NTN. First, it estimates how design choices in T causally affect authentication under severe Doppler, long delays, and fast channel variations, using causal learning to separate true effects from spurious correlations and account for environmental confounders. This also reveals whether the effect is direct or mediated through physical-layer features. Second, it trains a meta-learner that quickly adapts to new NTN environments, while incorporating the learned causal knowledge to improve robustness and consistency. Third, if the meta-learner cannot authenticate a user, a TDOA-based attack detection module is activated. A graph attention network then uses only the top-L most reliable nodes. This improves scalability and reduces complexity while preserving informative spatial observations. A. SCM for NTN The feature vector is obtained via (16)–(23), yet we need to discover the causal relationships between the environment, treatments, and the feature vector. For each node k and user u, the data-generating process is given by the SCM: f k,u (n) = g f,k (T,Z k,u ,U f,k ),(35) where g f,k (·) is an unknown deterministic function rep- resenting the underlying structural mechanisms, which are approximated in practice using neural network models in (36). The variable U f,k denotes exogenous random variable that is assumed to be mutually independent of T and Z k,u . It captures all unobserved sources of randomness affecting the feature extraction and authentication outcome, including thermal noise, channel estimation errors, and hardware im- pairments. The unknown function g f,k in (35) is implemented indirectly through the neural network components Φ φ for node k. The structural mechanisms are learned indirectly through the neural components. Hence, we use a neural network as a feature extractor Φ φ which maps the raw input X k,u (n) to an approximated feature vector b f k,u (n). b f k,u (n) = Φ φ (X k,u (n)),(36) The raw input X k,u (n) includes y k (n), ˆ h k,u (n), ˆ f (k,u) d (n), and ˆτ k,u (n). In summary, in reality the features are obtained via the relation in (16)–(23). These can be considered as an unknown function of T, Z k,u , and U f,k , as defined in (35). Thus, we estimate this function using (36). Generally, if the approximated feature ˆ f k,u (n) is similar to extracted feature f k,u (n), the transmission is more likely to be authenticated as legitimate. Furthermore, we find the identity by (28). We use Y k,u (n) ∈ 0, 1 as the binary local authentica- tion outcome at node k, where Y k,u (n) = 0 indicates that the observed transmission is authenticated as legitimate and Y k,u (n) = 1 otherwise. In practice, we produce a local authentication outcome score b Y k,u (n) at node k with the method in I-B2 by using estimated feature vector ˆ f k,u (n) instead of f k,u (n): b Y k,u (n) = ( 0, if d k,min (n)≤ τ ac , 1, otherwise. (37) Fig. 2: Causal graph of the data-generating process. As illustrated in Fig. 2, the corresponding causal graph includes the directed edges: Z k,u → f k,u , Z k,u → Y k,u (n), T → f k,u (n), T → Y k,u (n), and f k,u (n) → Y k,u (n). In the assumed causal graph, both Z k,u and T act as confounders for the causal relationship f k,u (n)→ Y k,u (n), since each directly affects both the cause (f k,u (n)) and the outcome (Y k,u (n)). We further allow statistical dependence between Z k,u and T, reflecting that system design choices may depend on the operating environment. Given the causal graph (Fig. 2), two causal relationships must be distinguished: 1) For the effect of f k,u (n) on Y k,u (n): Both Z k,u and T are direct confounders. Hence, to identify the causal effect f k,u (n) → Y k,u (n), we must condition on both Z k,u and T (i.e., include them as control variables). Failure to condition on either one opens a backdoor path and biases the estimate. 2) For the effect of T on Y k,u (n): There is no direct edge from Z k,u to T; the statistical dependence between them arises from an unobserved confounder U f,k . This creates the backdoor path T ← U f,k → Z k,u → Y k,u (n). Conditioning on Z k,u alone blocks this path because Z k,u is a descendant of U f,k and a parent of Y . Therefore, for the causal effect T→ Y k,u (n), it is sufficient to condition only on Z k,u , not on T itself (as that would condition on the treatment and introduce bias). The unobserved confounder U f,k does not need to be measured; its backdoor path is blocked by the observable Z k,u . In NTN scenarios, there exist unobserved confounders (e.g., network congestion, operator policies) that affect both the treatment T and the environmental variables Z k,u . These unobserved confounders create spurious backdoor paths. Al- though Z k,u themselves are not direct confounders (they have no edge to T), they act as sufficient adjustment variables be- cause conditioning on them blocks all backdoor paths between T and Y k,u (n). Therefore, we assume that the set Z k,u satisfies the backdoor criterion relative to (T,Y k,u (n)). Furthermore, we assume positivity (overlap) and the stable unit treatment value assumption (SUTVA). Under SUTVA and the backdoor criterion, the causal effect of T on Y k,u (n) is identifiable from observational NTN data. Let t 0 ,t 1 ∈T be two treatment values of interest. For each node k and user u, the interventional mean is defined as μ k,u (t)≜ E b Y k,u (n)| do(T = t) ,(38) and the corresponding average causal effect (ACE) is ACE k,u (t 1 ,t 0 )≜ μ k,u (t 1 )− μ k,u (t 0 ).(39) 7 Under the backdoor assumption with respect to Z k,u , the interventional mean can be identified via μ k,u (t) = E Z k,u h E b Y k,u (n)| T = t,Z k,u i .(40) Under the standard consistency, conditional exchangeability, and positivity assumptions, the interventional mean μ k,u (t) is identifiable from observational data via the backdoor formula; the proof is provided in Appendix A-A of the Supplemental Material. The conditional average treatment effect (CATE) at a given environment realization z k,u is defined as CATE k,u (t 1 ,t 0 | z k,u )≜ E b Y k,u (n)| do(T = t 1 ),Z k,u = z k,u − E b Y k,u (n)| do(T = t 0 ),Z k,u = z k,u .(41) Under the backdoor condition, this can be written as CATE k,u (t 1 ,t 0 | z k,u ) = E b Y k,u (n)| T = t 1 ,Z k,u = z k,u − E b Y k,u (n)| T = t 0 ,Z k,u = z k,u .(42) In NTN scenarios, CATE analysis is particularly useful for characterizing how the effectiveness of a given system design varies with environmental factors such as the elevation angle θ (k,u) elev and node velocity v k . Let f k,u (t ′ ) denote the potential feature vector at node k for user u under treatment T = t ′ at time index n. For notational simplicity, the time index n is omitted throughout this section. This vector serves as the mediator according to the causal graph. Let b Y k,u (t, f k,u (t ′ )) denote the model-predicted potential outcome at time index n when the treatment is set to t and the mediator is set to f k,u (t ′ ). We also omit n from b Y k,u (t, f k,u (t ′ )) for simplicity. Thus, the natural direct effect (NDE) is defined as NDE k,u (t 1 ,t 0 ) = E h b Y k,u t 1 , f k,u (t 0 ) i − E h b Y k,u t 0 , f k,u (t 0 ) i ,(43) and the natural indirect effect (NIE) is given by NIE k,u (t 1 ,t 0 ) = E h b Y k,u t 0 , f k,u (t 1 ) i − E h b Y k,u t 0 , f k,u (t 0 ) i ,(44) The total effect is decomposed as ACE k,u (t 1 ,t 0 ) = NDE k,u (t 1 ,t 0 ) + NIE k,u (t 1 ,t 0 ),(45) which is formally established in Appendix A-C of the Sup- plemental Material. Given N s samplescollectedfrommultipleNTN environmentsandnodes,wedenotethedatasetas (Z (i) k,u , T (i) , f (i) k,u , b Y (i) k,u (n)) N s i=1 ,whereeachsample corresponds to an observation at node k, user u, sample i, and time index n. To simplify the formulation we omit time index n. We define the nuisance models as e t (z)≜ P(T (i) = t| Z (i) k,u = z),(46) m t (z)≜ E[ b Y (i) k,u (n)| T (i) = t,Z (i) k,u = z].(47) The doubly robust (DR) estimator (also known as aug- mented inverse probability weighting, AIPW) of the interven- tional mean μ k,u (t) is given by bμ DR k,u (t) = 1 N N X i=1 " m t (Z (i) k,u ) + 1T (i) = t e t (Z (i) k,u ) b Y (i) k,u (n)− m t (Z (i) k,u ) # . (48) This estimator is consistent under standard ignorability, overlap, and correct specification of either the propensity model e t (z) or the outcome model m t (z); the proof is given in Appendix A-B of the Supplemental Material. For compar- ison, the inverse probability of treatment weighting (IPTW) estimator is defined as bμ IPTW k,u (t) = 1 N s N s X i=1 1T (i) = t b Y (i) k,u (n) e t (Z (i) k,u ) .(49) B. Meta-Learning with Causal Invariance for NTN Let E = e 1 ,e 2 ,...,e |E| denote the set of meta-training NTN environments. Each environment e ∈ E corresponds to a distinct operating condition characterized by context variables Z e k,u . For each environment e and observing node k, we define a meta-learning task indexed by the pair (e,k), with training dataset D e,k train and validation dataset D e,k val . Each dataset contains samples from multiple users. We use a shared authentication model parameterized by φ, as defined in (35). Φ φ is trained to capture causal and environment-invariant structure across tasks. For each task (e,k), one inner-loop adaptation step is performed as φ ′ e,k = φ− α∇ φ L e,k task (φ;D e,k train ),(50) where α is the inner-loop learning rate and L e,k task denotes the task-specific authentication loss. The meta-parameters are then updated using the validation losses across tasks: φ← φ− β X e,k ∇ φ L e,k task (φ ′ e,k ;D e,k val ),(51) where β is the outer-loop learning rate. To encourage Φ φ to extract environment-invariant causal features, we add an invariant risk minimization (IRM) penalty: R IRM (φ) = X e,k ∇ w R e,k (w,φ) w=1 2 ,(52) where R e,k (w,φ) denotes the empirical risk for task (e,k) under a scalar classifier parameter w used only to construct the IRM penalty. If the learned representation is invariant across environments, the IRM penalty vanishes; this statement is formally proved in Appendix A-D of the Supplemental Material. To further promote causal stability, we regularize the task- specific DR estimates of the interventional mean. Let bμ DR e,k (t) 8 denote the DR estimate of the interventional mean under treatment t computed from task (e,k), and let ̄μ DR −e,−k (t) = 1 |I|− 1 X (e ′ ,k ′ )∈I\(e,k) bμ DR e ′ ,k ′ (t),(53) where I = (e,k) : e ∈ E,k ∈ K is the set of all task indices. We define the DR consistency regularizer as R DR = X t∈T X (e,k)∈I bμ DR e,k (t)− ̄μ DR −e,−k (t) 2 .(54) This term encourages, but does not strictly enforce, stability of interventional responses across environments and nodes. The final meta-training objective is min φ X e,k L e,k task φ ′ e,k ;D e,k val + λ IRM R IRM (φ) + λ DR R DR , (55) where λ IRM and λ DR control the strength of the invariance and causal-consistency regularization terms, respectively. C. Counterfactual Design Policy for NTN We learn a policy π : Z k,u 7→ T that selects the design configuration for each environment as observed by node k. Thus, for a given Z k,u , the applied treatment is π(Z k,u )∈T . The goal is to maximize the interventional performance under the induced policy, defined as J (π) = E Z k,u μ k,u (π(Z k,u )) , where μ k,u (t) is the interventional mean defined in (38). We estimate J (π) from logged data collected from all nodes using the DR off-policy estimator: b J DR (π) = 1 N s N s X i=1 m π(Z (i) k,u ) (Z (i) k,u ) + 1 N s N s X i=1 1T (i) = π(Z (i) k,u ) e T (i) (Z (i) k,u ) b Y (i) k,u (n)− m T (i) (Z (i) k,u ) , (56) where m t (·) is the outcome model in (47), e t (·) is the propensity model in (46). This estimator is DR, i.e., it remains consistent if either the outcome model or the propensity model is correctly specified. Its off-policy nature is especially valuable in NTN, as it enables policy evaluation and optimiza- tion without live deployment, thereby reducing the cost and operational risk of HAPS and satellite platforms. D. Scalable Top-L Graph-Based Fusion at the FC In attack detection mode, the FC collects local authentica- tion evidence from nodes in K u and applies a top-L graph- based fusion strategy. For each node k ∈ K u , it first forms a node-level evidence vector as x k,u (n)≜ h p k,u,u ⋆ k (n), f T k,u (n), ˆτ k,u (n), ˆ f (k,u) d (n), r T k (n) i T , (57) where x k,u (n) ∈ R d x and d x ≜ d f + 6. The FC evaluates the reliability of the local evidence from node k through a trainable function ρ k,u (n)≜ g θ ρ x k,u (n) , where g θ ρ is a learnable node-scoring network that maps R d x to R and is parameterized by θ ρ . Hence, the FC selects the L nodes: K (L) u (n)≜ TopL k∈K u ρ k,u (n) ,(58) where L ∈ N is a design parameter satisfying 1 ≤ L ≤ |K u |. To avoid ambiguity between node labels and matrix indices, we enumerate the selected nodes as K (L) u (n) = κ 1 ,κ 2 ,...,κ L . A graph is then constructed over the selected nodes G (L) u (n)≜ V (L) u (n), E (L) u (n) , where V (L) u (n) and E (L) u (n) are the set of nodes and the set of edges of the graph, respectively. In this paper, we consider V (L) u (n) = K (L) u (n). The graph is complete on the selected nodes, i.e., (κ i ,κ j ) ∈ E (L) u (n), ∀i,j ∈1,...,L, i̸= j. For the selected nodes, the corresponding node features are simply obtained by evaluating (57) at k = κ i , i.e., x κ i ,u (n)≜ x k,u (n) k=κ i , i = 1,...,L.(59) For each pair of selected nodes (κ i ,κ j ), the FC further defines the edge feature vector e i,j,u (n)∈ R 3 as following: e i,j,u (n)≜ ∆τ i,j,u (n), ∆f i,j,u (n), ∆r i,j (n) T ,(60) where ∆τ i,j,u (n) = ˆτ κ i ,u (n) − ˆτ κ j ,u (n), ∆f i,j,u (n) = ˆ f (κ i ,u) d (n)− ˆ f (κ j ,u) d (n), and ∆r i,j (n) =∥r κ i (n)− r κ j (n)∥. Using the pairwise delay and Doppler consistency, the FC constructs a weighted adjacency matrix A (L) u (n) ∈ R L×L , whose (i,j)-th element is given by [A (L) u (n)] i,j ≜ exp − |∆τ i,j,u (n)| σ τ − |∆f i,j,u (n)| σ d , i̸= j, 0,i = j, (61) where σ τ > 0 and σ d > 0 are design parameters controlling the sensitivity to delay and Doppler mismatch, respectively. Next, the FC stacks the selected node features into the matrix X (L) u (n)∈ R L×d x , defined as X (L) u (n)≜ x T κ 1 ,u (n) x T κ 2 ,u (n) . . . x T κ L ,u (n) .(62) These features are first mapped into a latent space through a trainable input embedding function H (0) u (n)≜ φ in X (L) u (n) , where H (0) u (n) ∈ R L×d h , φ in : R d x → R d h is applied row- wise, and d h ∈ N is the hidden embedding dimension. Graph message passing is then performed over L G layers with H (ℓ+1) u (n) = φ (ℓ) e D − 1 2 u (n) e A (L) u (n) e D − 1 2 u (n)H (ℓ) u (n)W (ℓ) , (63) for ℓ = 0, 1,...,L G − 1, where e A (L) u (n)≜ A (L) u (n) + I L , and I L is the L× L identity matrix, e D u (n) ∈ R L×L is the diagonal degree matrix associated with e A (L) u (n), W (ℓ) is the trainable weight matrix of the ℓ-th layer, and φ (ℓ) (·) is the 9 corresponding nonlinear activation function. The final graph- layer output is denoted by H (L G ) u (n) = s T 1,u (n) s T 2,u (n) . . . s T L,u (n) ,(64) where s i,u (n) ∈ R d h is the final embedding associated with the selected node κ i and H (L G ) u (n)∈ R L×d h . The FC aggregates the selected node embeddings via attention-based readout, where the attention logit of the i-th node is a i,u (n)≜ v T a tanh W a s i,u (n) = exp a i,u (n) P L j=1 exp a j,u (n) , (65) where W a and v a are trainable attention parameters. The graph-level fused representation is then obtained as F (u) g (n)≜ L X i=1 α i,u (n)s i,u (n),(66) where F (u) g (n)∈ R d h . To preserve the geometric consistency information already introduced in (31) and (32), the FC augments the graph-level representation with the triangulation consistency score and the identity-location consistency score, which is represented by ̄ F (u) g (n)∈ R d h +2 and defined as ̄ F (u) g (n)≜ h F (u) g (n) T , f tri,u (n), f dist,u,u ⋆ (n) i T .(67) The final authentication score generated by the FC is S (u) FC (n)≜ s θ FC ̄ F (u) g (n) , where s θ FC is decision neural network. Here, lower values of S (u) FC (n) indicate stronger evidence of legitimacy, whereas larger values indicate stronger evidence for rejection. The overall FC-side parameter set is θ FC ≜ n θ ρ , φ in , W (ℓ) L G −1 ℓ=0 , W a , v a , w o , b o o ,(68) where w o ∈ R d h +2 and b o ∈ R. Finally, the FC makes the authentication decision by D (u) final (n) = 1 n S (u) FC (n) > τ FC o ,(69) where τ FC ∈ (0, 1) is the decision threshold. Moreover, D (u) final (n) = 0 denotes acceptance as a legitimate transmission and D (u) final (n) = 1 denotes rejection. The proposed fusion strategy balances scalability and ex- pressiveness by limiting graph construction and message pass- ing to the L most reliable nodes, thereby reducing computa- tional cost when |K u | is large. Meanwhile, the graph captures pairwise consistency in delay, Doppler, and geometry, enabling reliable authentication in dynamic NTN environments. E. Verification Score and ROC Metrics The receiver operating characteristic (ROC) curve shows the tradeoff between correctly accepting legitimate transmissions and wrongly accepting attacker transmissions as the threshold τ FC changes. Here, the positive class corresponds to legitimate Fig. 3: Summary of the proposed method and algorithms. users. For a given τ FC , the true positive rate (TPR) and false positive rate (FPR) are defined as TPR(τ FC ) = TP (τ FC ) TP (τ FC ) + FN (τ FC ) ,(70) FPR(τ FC ) = FP (τ FC ) FP (τ FC ) + TN (τ FC ) .(71) Here, TP (τ FC ) is the number of legitimate transmissions correctly accepted (D (u) final (n) = 0), FN (τ FC ) is the number of legitimate transmissions wrongly rejected (D (u) final (n) = 1), FP (τ FC ) is the number of attacker transmissions wrongly accepted (D (u) final (n) = 0), and TN (τ FC ) is the number of attacker transmissions correctly rejected (D (u) final (n) = 1). The area under the ROC curve (AUC) is defined as AUC = Z 1 0 TPR(FPR)dFPR = Z TPR(τ FC ) d FPR(τ FC ) dτ FC dτ FC ,(72) where the absolute value ensures a nonnegative area under the ROC curve. Under the convention in (69), increasing τ FC generally makes acceptance less strict, since a transmission is accepted when S (u) FC (n)≤ τ FC . F. Summary of the Algorithms As shown in Fig. 3, the proposed framework consists of four stages: data collection, causal estimation, meta-learning, and authentication. First, NTN observations are collected across environments and nodes, including variables, treatments, mea- surements, features, and labels. Second, the causal module estimates interventional quantities such as DR means, ACE, and CATE. Third, the meta-learning module optimizes φ using MAML with causal-invariance and DR-consistency regular- ization to improve cross-environment generalization. Finally, the trained model produces the FC authentication score and accept/reject decision, while the graph-based attack detection module is activated if the claimed identity is not legitimate. Algorithms 1 and 2 summarize the causal estimation and meta-learning procedures, respectively. The causal estimates guide training through invariance and DR consistency regular- izers, encouraging the model to exploit stable causal relation- ships rather than environment-specific correlations. The final model is evaluated on unseen NTN environments. IV. COMPUTATIONAL COMPLEXITY ANALYSIS We analyze the computational complexity of the proposed framework for both offline training and online inference. 10 Algorithm 1 Causal Effect Estimation Input: T , D causal = n Z (i) k,u ,T (i) ,f (i) k,u ,Y (i) k,u o N s i=1 Output: bμ DR (t), [ ACE(t 1 ,t 0 ), and \ CATE(t 1 ,t 0 | z k,u ) 1: Estimate propensity scores e t (z) via (46), ∀t∈T . 2: Estimate outcome models m t (z) via (47), ∀t∈T . 3: Compute DR means bμ DR (t) via (48), ∀t∈T . 4: Compute ACEs [ ACE(t 1 ,t 0 ) via (39), ∀(t 1 ,t 0 )∈T 2 . 5: Compute \ CATE(t 1 ,t 0 | z k,u ) via (41), ∀(t 1 ,t 0 )∈T 2 and z k,u . 6: return bμ DR (t) t∈T , [ ACE(t 1 ,t 0 ) t 0 ,t 1 ∈T , 7:and \ CATE(t 1 ,t 0 | z k,u ) t 0 ,t 1 ∈T ,z k,u ∈Z k,u . Algorithm 2 Meta-Learning with Causal Invariance Input: D e,k train ,D e,k val , α,β, λ IRM ,λ DR , and outputs from Algorithm 1 Output: Learned initialization φ ⋆ 1: Initialize φ randomly. 2: while not converged do 3:Compute task-adapted parameters φ ′ e,k using (50), ∀e∈E,k ∈K. 4:Evaluate validation losses L e,k task (φ ′ e,k ;D e,k val ), ∀e∈E,k ∈K. 5:Compute R IRM (φ) and R DR (φ) using (52) and (54). 6:Update φ using the outer update in (51). 7: end while 8: return φ ⋆ ← φ. For each node k, feature extraction in (16)–(23) is dom- inated by FFT operations and subspace projection, resulting in per-sample complexity O(M k logM k + M k r k ). The local Mahalanobis scoring has negligible complexity compared with feature extraction. Hence, the total local-processing complex- ity over N s samples is O N s P |K| k=1 (M k logM k + M k r k ) . Let M≜ max k M k denote the maximum antenna number. Under practical LoS-dominant NTN conditions with r k = O(1), this reduces to O(N s |K|M logM ), while in rich scat- tering it becomes O(N s |K|M 2 ). The causal estimation module fits the propensity and out- come models over a dataset of size N s . Since |T| is fixed and small, the overall causal-estimation complexity scales as O(N s |T|). For meta-learning, let C ML denote the complexity of one forward-backward pass, where C ML =O(d x d h + d 2 h ). The total training complexity over n iter outer iterations be- comes O(n iter N s (C ML +|T|)). If the feature extractor is trained end-to-end, the local feature-extraction cost is ad- ditionally included. Since d x , d h , and |T| remain bounded in practical implementations, the meta-learning and causal- estimation terms scale approximately linearly with n iter N s . For one authentication attempt, each node performs channel estimation, feature extraction, and local scoring with com- plexity O(M k logM k + M k r k ). If distributed verification is activated, Top-L node selection requires O(|K u | log|K u |), while graph construction and message passing require O(L 2 ) and O(L G L 2 d h ), respectively. Optional TDOA localization adds O(L 3 ). However, since L, L G , and d h are fixed small values in our implementation, and TDOA is only triggered when local authentication fails, the dominant online complex- ity is O P k∈K u (M k logM k + M k r k ) . For homogeneous node dimensions M k ≈ M and bounded r k , this reduces to O(|K u |M logM ). The local-only baseline employs a conventional deep learning-based authentication model using the same local TABLE I: Final computational complexity comparison. SchemeComplexity ProposedOffline: O(N s |K|M logM + n iter N s ) Online: O(|K u |M logM ) Local-onlyOffline: O(N s |K|M logM ) Online: O(|K u |M logM ) Typical parameter ranges N s ≈ 10 3 −10 4 , |K|≈ 2−20 |K u |≈ 2−10, M ≈ 32−256 n iter ≈ 10−50 feature extraction pipeline, but without causal estimation, meta-learning adaptation, graph-based fusion, or TDOA- based attack detection. It still performs local feature extrac- tion and local scoring. Thus, its offline complexity is ap- proximately O(N s |K|M logM ), while its online complexity is O(|K u |M logM ). Under LoS-dominant NTN conditions, the proposed framework has simplified offline complexity O(N s |K|M logM + n iter N s ), and dominant online complex- ity O(|K u |M logM ). Thus, the proposed framework intro- duces additional offline computation for causal estimation and meta-learning, while the online cost is dominated by local processing and remains scalable. Table I summarizes the dominant computational complexity comparison. The online expressions reported in the table correspond to the dominant terms after treating L, L G , d h , and conditional TDOA localization as lower-order contributions. V. SIMULATION RESULTS We first generate a multi-environment dataset and train the model. Then, we test it under unseen conditions and analyze the features of the proposed multi-feature fingerprint. A. Simulation Setup 1) Data Generation: We generate a dataset covering HAPS and satellite scenarios. Six scenarios are considered: five for training and one for testing. Each scenario includes ten users with IDs 1 to 10, and each user performs ten movements. The objective is user-ID recognition from the extracted fingerprint features. The uplink channel follows a Rician fading model with carrier frequency f c = 2 GHz. For each task, the envi- ronmental variables are sampled from the ranges in Table I. We generate 10 samples per user per scenario. Thus, the first- phase training set contains 500 samples, i.e., 5 scenarios× 10 users × 10 samples, while the test set contains 100 samples, i.e., 1 scenario × 10 users × 10 samples. The physical feature vector f k,u (n) is computed using (16)–(23). Channel estimation includes AWGN with power σ 2 k = 10 −12 W, TOA noise with standard deviation σ toa = 5 ns, and Doppler noise with standard deviation 0.5 Hz. Legitimate users and attackers are located within a 50 × 50 km ground area centered at the origin. For each environment, N enr = 100 enrollment snapshots per legitimate user are collected to estimateμ (k,u) h , R (k,u) h , and P (k,u) ref . The reference subspace uses the two dominant eigenvectors, i.e., r k = 2. The reference Doppler and delay values are denoted by f (k,u) d,ref and τ (ref) k,u , respectively. 11 TABLE I: Ranges of environment variables for task generation. ParameterSymbolRange / Values Number of antennas M k 16, 64, 128, 256 SNR at node kγ k,u −10, −5,..., 30 dB Angular spreadσ θ,k 2 ◦ , 5 ◦ , 10 ◦ Rician factorK (k,u) r 5, 10, 20 ADC resolutionN (k) bit 4, 6, 8 bits Pilot reuse groupsK (k) p 2, 3, 4 Elevation angleθ (k,u) elev 30 ◦ , 60 ◦ , 90 ◦ HAPS-specific HAPS altitudeH (k) alt 20, 25 km HAPS velocityv k 15, 20 m/s LEO/GEO-specific Satellite altitudeH (k) alt 500, 1200, 2000, 35786 km Satellite velocityv k 7200, 7500, 7800, 3070 m/s 2) Causal Meta-Learning for Local Authentication: The causal meta-learning framework in Algorithms 1 and 2 is trained on a set of |E| = 5 distinct environments, corre- sponding to the scenarios in the data generation code. Each environment contains 10 users (either legitimate or attacker). The feature dimension is d f = 8. For each task, we randomly sample 5 support samples and 5 query samples per user in a few-shot setting. The inner learning rate is α = 0.01 with 5 inner gradient steps. The outer learning rate is β = 0.001. The meta-learner is trained for n iter = 160 outer epochs, with 20 episodes per epoch. The hyperparameters λ IRM and λ DR are set to 0.1 and 0.05, respectively. The neural network Φ φ , used as the feature extractor is a 4-layer MLP with hidden size 64 and output dimension d f , one output per feature. 3) Distributed Verification and Attack Detection: For the second-phase experiments, i.e., distributed verification using TDOA and location-consistency scores, we consider |K| = 4 nodes located at (±50,±50) km horizontally. Each node is placed at altitude 20 km, has M k = 16 antennas, and moves with velocity v k = 15 m/s, as used in the verification code. The legitimate user and the attacker are placed on the ground within the (±50,±50) km region. The training set consists of 160 samples for legitimate users and 160 samples for attackers. The test set contains 100 samples for each class. The TDOA estimation uses 10 consecutive time steps. The position-stability score uses σ pos = 5 m, and the distance- consistency score uses σ dist = 500 m. The FC classifier is a two-layer MLP with hidden units (64, 32) and a sigmoid output. It is trained with the Adam optimizer for 1000 epochs, without early stopping. The final decision threshold is τ FC = 0.5. All results are averaged over n MC = 10 independent Monte Carlo runs. The SNR sweep varies the transmit power from 10 −4 to 10 2 W, corresponding to an SNR range from −10 to 30 dB. 4) Test Scenario: In the final test, the system is evaluated on unseen environment combinations. The meta-learner is adapted to the new environment using 10 support samples per user via the inner update and tested on 100 query samples. In distributed verification, the FC collects local soft scores p k,u,u ⋆ k (n) and TDOA measurements ˆτ k,u (n) from all nodes. The reported metrics are AUC and authentication accuracy. Fig. 4: Average ROC curves for single-feature classifiers, baseline fusion, and the proposed all-feature fusion method across the considered NTN scenarios. B. Multi-Feature Comparison and Feature Importance Fig. 4 compares the AUC of different methods. The pro- posed method achieves the highest AUC by using all eight features, learning the causal relationships, and adapting to new environments via meta-learning. In contrast, the deep-learning baseline may overfit environment-specific patterns and perform worse in unseen NTN conditions. The relevance of PLA features depends on the NTN sce- nario, link geometry, interference level, and channel estimation quality. Overall, geometry- and spatial-domain features, i.e., f k,u,1 , f k,u,5 , and f k,u,6 , are generally reliable across most NTN deployments, whereas interference-aware features, such as f k,u,7 and f k,u,8 , become increasingly important in dense multi-user settings with IUI or pilot contamination. Feature importance varies with propagation conditions, platform type, system parameters, and attacker capability. In sparse LoS-dominant scenarios, angular, subspace, and geometry-dependent features, i.e., f k,u,1 , f k,u,2 , f k,u,5 , and f k,u,6 , are typically the most discriminative, whereas remov- ing f k,u,7 or f k,u,8 causes only limited degradation due to weak interference. In dense NTN access scenarios, post- combining and interference-aware features become more im- portant; hence, f k,u,7 , f k,u,8 , f k,u,3 , and f k,u,5 are more influential, while the snapshot-dependent projection feature f k,u,4 is usually less robust under noise and interference. The dominant features also depend on the NTN platform. For LEO links, high mobility increases the importance of the delay–Doppler feature f k,u,6 , together with spatial features f k,u,1 , f k,u,3 , and f k,u,5 . For GEO links, the weaker Doppler effect makes f k,u,1 , f k,u,2 , and f k,u,5 dominant, although the delay component of f k,u,6 can still be useful. For HAPS links, stable geometry and shorter propagation distance make angular and subspace features, particularly f k,u,1 , f k,u,2 , and f k,u,5 , highly informative. System conditions further affect the ranking. High SNR and large antenna arrays improve channel estimation and angular resolution, strengthening f k,u,1 , f k,u,2 , and f k,u,5 . At low SNR, with few antennas, strong NLoS propagation, large angular spread, or pilot contamination, instantaneous and beam-domain features become less stable, while f k,u,3 , f k,u,5 , f k,u,6 , f k,u,7 , and f k,u,8 provide more robust evidence. 12 TABLE I: Qualitative importance of the considered physical-layer authentication features under different NTN environments. Feature Name (No- tation) Typical NTN Dense Ac- cess Sparse LoSLEOGEO/HAPSMain Interpretation Channel Similarity (f k,u,1 ) HighMedium– High Very HighHighVery HighCaptures spatial consistency with the enrolled user. Main-Lobe Energy (f k,u,2 ) MediumLow– Medium HighMediumHighUseful when dominant angular support is stable; weaker under scattering, beam spread, or contamination. Combiner Similar- ity (f k,u,3 ) Medium– High HighMediumHighMediumMore important when spatial filtering consistency mat- ters, especially in dense interference-limited settings. ProjectionScore (f k,u,4 ) LowLowLow– Medium LowLowUsually the least robust feature due to sensitivity to noise, interference, and snapshot variability. Subspace Deviation (f k,u,5 ) Very HighHighVery HighHighVery HighOne of the most informative features, reflecting deviation from the enrolled spatial structure. Delay–Doppler (f k,u,6 ) Very HighHighHighVery HighMediumStrong when geometry and motion are discriminative; particularly important for LEO links. SINR (f k,u,7 )HighVery HighLow– Medium HighMediumCritical in dense access and strong-IUI conditions be- cause it reflects user separability after combining. INR (f k,u,8 )MediumVery HighLowMediumLow– Medium Mainly useful in interference-heavy scenarios where il- legitimate activity alters the interference level. TABLE IV: Scenario-dependent qualitative feature importance for physical-layer authentication in NTN environments. ScenarioPrimarySecondaryLessLeastNotes Typical NTN Environmentf k,u,6 , f k,u,5 , f k,u,1 f k,u,3 , f k,u,7 f k,u,2 , f k,u,8 f k,u,4 Geometry consistency and spatial structure are usually the strongest authentication cues. Dense Multi-User Accessf k,u,7 , f k,u,8 , f k,u,3 f k,u,5 , f k,u,6 , f k,u,1 f k,u,2 f k,u,4 Interference-aware features dominate under strong IUI and pilot contamination. Sparse LoS-Dominantf k,u,1 , f k,u,5 , f k,u,6 , f k,u,2 f k,u,3 f k,u,7 f k,u,4 , f k,u,8 Spatial, angular, and geometric consistency are highly discriminative under weak interference. LEO Linksf k,u,6 , f k,u,1 f k,u,3 , f k,u,5 f k,u,7 , f k,u,2 f k,u,4 High mobility enhances the value of delay– Doppler cues. GEO/HAPS Linksf k,u,1 , f k,u,2 , f k,u,5 f k,u,6 , f k,u,3 f k,u,7 , f k,u,8 f k,u,4 Stable geometry reduces Doppler relevance and emphasizes spatial consistency. Low-SNRf k,u,6 , f k,u,7 f k,u,5 , f k,u,3 f k,u,1 , f k,u,2 f k,u,4 Spatial features degrade, while geometric and interference-aware cues become more useful. Attacker location and capability also influence feature dis- crimination. When the attacker is far from the legitimate user, f k,u,6 is highly effective due to distinct geometric signatures. However, for a nearby or more capable attacker that partially mimics delay–Doppler characteristics, authentication relies more on spatial consistency, combiner similarity, subspace deviation, and post-combining SINR, represented by f k,u,1 , f k,u,3 , f k,u,5 , and f k,u,7 , respectively. Overall, a qualitative ranking for typical NTN scenarios is f k,u,6 ≈ f k,u,5 ≈ f k,u,1 > f k,u,3 ≈ f k,u,7 > f k,u,2 > f k,u,8 > f k,u,4 .(73) In dense and interference-limited scenarios, the ranking shifts toward interference-aware and combiner-based features as f k,u,7 ≈ f k,u,8 ≈ f k,u,3 > f k,u,5 ≈ f k,u,6 ≈ f k,u,1 > f k,u,2 > f k,u,4 .(74) In sparse LoS-dominant environments, the ranking becomes f k,u,1 ≈ f k,u,5 ≈ f k,u,6 ≈ f k,u,2 > f k,u,3 > f k,u,7 > f k,u,8 ≈ f k,u,4 .(75) These rankings complement empirical feature-importance analysis and show feature roles under different NTN condi- tions. Tables I and IV summarize them from feature- and scenario-centric perspectives. Table V groups the features into spatial, combiner-based, subspace-based, geometry-based, and interference-aware categories for easier interpretation. TABLE V: Categorization of the considered physical-layer features. CategoryFeature Notation(s) Spatial / Angular Consistencyf k,u,1 , f k,u,2 Combiner / Receiver-Side Consistency f k,u,3 , f k,u,4 Subspace Structure Deviationf k,u,5 Geometry / Kinematicsf k,u,6 Interference-Awaref k,u,7 , f k,u,8 Some redundancy exists among the proposed features. Specifically, f k,u,1 , f k,u,2 , and f k,u,5 all characterize spatial consistency and may therefore be correlated, particularly in strong LoS conditions with small angular spread. Similarly, f k,u,7 and f k,u,8 are both interference-related; however, f k,u,7 is usually more informative because it reflects the effective SINR, whereas f k,u,8 only measures interference relative to noise and is mainly useful in dense access scenarios. To quantitatively validate the above observations, we per- form a leave-one-feature-out ablation study. Let AUC all de- note the AUC obtained using the complete features and let AUC −f i denote the AUC obtained after removing the i-th feature. The contribution of feature f i is derived as ∆AUC i = AUC all − AUC −f i . A large value of ∆AUC i indicates that the removed feature is important for authentication, whereas a value close to zero indicates that the feature has limited marginal contribution under the considered environment. Figs. 5a and 5b show environment-dependent feature impor- tance via leave-one-feature-out AUC and accuracy degradation in NTN scenarios. In sparse LoS settings, removing f k,u,1 , 13 (a)(b) Fig. 5: Impact of each feature across scenarios: (a) AUC and (b) accuracy. (a)(b) (c)(d) Fig. 6: AUC and accuracy analyses: (a),(b) feature-removal impacts on AUC and accuracy; (c),(d) feature-importance heatmaps for AUC and accuracy. f k,u,5 , or f k,u,6 leads to the largest degradation, highlighting the role of stable geometry and spatial structure. In dense or pilot-contaminated scenarios, f k,u,3 , f k,u,7 , and f k,u,8 become more influential due to post-combining separability and inter- ference effects. The results also show that f k,u,6 dominates in sparse LEO links because of strong Doppler signatures, while GEO scenarios depend more on spatial-domain features. The ablation results show no universally dominant feature. This motivates causal meta-learning over the full feature set. Figs. 6a and 6b report the average ∆AUC i and accuracy variations across NTN conditions. f k,u,4 is generally the least influential, while f k,u,6 is more relevant in LEO, low-SNR, and far-attacker scenarios. The interference-aware features f k,u,7 and f k,u,8 are most useful in dense, wide-angular- spread, and pilot-contaminated environments, whereas f k,u,1 and f k,u,5 remain strong in sparse LoS, GEO, and large-array settings. Figs. 6c and 6d further confirm spatial dominance in sparse LoS and interference dominance in dense scenarios. VI. CONCLUSION This paper proposed SAFA-MZ, a causal meta-learning framework for DPLA in NTNs under pilot contamination. SAFA-MZ builds an adaptive and distributed physical-layer fingerprint using multi-node features and adjusts the fusion process to changing NTN environments. It combines causal analysis with meta-learning to improve robustness under dis- tribution shifts and enable fast adaptation with few labeled samples. A two-stage authentication mechanism performs lo- cal authentication first and activates TDOA-based distributed verification only when needed to reduce overhead. Simulation results showed that SAFA-MZ achieves 92% accuracy and 96% AUC, outperforming baseline and single-feature methods. APPENDIX A PROOFS A. Backdoor identification of μ k,u (t) Proposition 1. Under consistency, conditional exchangeability Y k,u (t) ⊥ T | Z k,u , and positivity P(T = t | Z k,u = z) > 0, the interventional mean is identifiable as μ k,u (t) = E[Y k,u | do(T = t)] = E Z k,u [E[Y k,u | T = t,Z k,u ]]. Proof. By iterated expectation and using conditional ex- changeability, we have μ k,u (t) =E[Y k,u (t)] = E Z k,u [E[Y k,u (t)| Z k,u ]] =E Z k,u [E[Y k,u (t)| T = t,Z k,u ]]. By consistency, whenever T = t we have Y k,u (t) = Y k,u , hence E[Y k,u (t) | T = t,Z k,u ] = E[Y k,u | T = t,Z k,u ]. Combining the above equalities yields μ k,u (t) = E Z k,u [E[Y k,u | T = t,Z k,u ]]. Positivity guarantees that the conditional expectation on the right-hand side is well-defined on the support of Z k,u . B. Double robustness of the estimator in (48) Theorem 1. Define the population estimating function: ψ t ≜ m t (Z k,u ) + 1T = t e t (Z k,u ) Y k,u − m t (Z k,u ) . Assume positivity and bounded moments of Y k,u . Then E[ψ t ] = μ k,u (t). Consequently, the sample estimator in (48) is consis- tent under standard regularity conditions (e.g., i.i.d. sampling and fixed or cross-fitted nuisance estimators). Proof. We prove the result in the two standard cases. a) Case 1: the outcome model is correct.: E[ψ t | Z k,u = z] = m t (z) + 1 e t (z) E 1T = t Y k,u − m t (z) | Z k,u = z = m t (z) + 1 e t (z) P(T = t| Z k,u = z) E [Y k,u − m t (z)| T = t,Z k,u = z]. The last conditional expectation is zero, and therefore E[ψ t | Z k,u = z] = m t (z) = E[Y k,u | T = t,Z k,u = z]. Taking expectation over Z k,u gives E[ψ t ] = E Z k,u [E[Y k,u | T = t,Z k,u ]] = μ k,u (t), where the last equality follows from the proof above. 14 b) Case 2: the propensity model is correct.: E[ψ t | Z k,u = z] = m t (z) + 1 e t (z) E[1T = tY k,u | Z k,u = z] − m t (z) e t (z) E[1T = t| Z k,u = z] = m t (z) + P(T = t| Z k,u = z) E[Y k,u | T = t,Z k,u = z] e t (z) − m t (z) e t (z) P(T = t| Z k,u = z) = E[Y k,u | T = t,Z k,u = z]. Taking expectation over Z k,u again yields E[ψ t ] = μ k,u (t). Thus, the estimator is DR, remaining consistent if either the propensity or outcome model is correctly specified. C. Correct mediation decomposition Proof. Add and subtract the same intermediate term E[Y k,u (t 1 , f k,u (t 0 ))] to prove the decomposition: ACE k,u (t 1 ,t 0 ) = E h b Y k,u (t 1 , f k,u,t 1 )− b Y k,u (t 0 , f k,u (t 0 )) i = E h b Y k,u (t 1 , f k,u,t 1 )− b Y k,u (t 1 , f k,u (t 0 )) i + E h b Y k,u (t 1 , f k,u (t 0 ))− b Y k,u (t 0 , f k,u (t 0 )) i = NIE k,u (t 1 ,t 0 ) + NDE k,u (t 1 ,t 0 ). D. Zero IRM penalty under an invariant representation Proposition 2. Assume there exists a representation Φ φ ⋆ such that for every triple (e,k,u), w = 1 is a stationary point of the risk R e,k,u (w◦ Φ φ ⋆ ). Then R IRM (φ ⋆ ) = 0. Proof. Substituting φ ⋆ into (52), if w = 1 is stationary for every (e,k,u), then each gradient term is zero. Hence, all summands vanish and R IRM (φ ⋆ ) = 0. REFERENCES [1] T. M. Hoang, A. Vahid, H. D. Tuan, and L. Hanzo, “Physical layer authentication and security design in the machine learning era,” IEEE Communications Surveys & Tutorials, vol. 26, no. 3, p. 1830–1860, 2024. [2] J. Zhang, Y. Guo, Y.-W. P. Hong, and S. Tomasin, “Model-driven learning-based physical layer authentication for mobile wi-fi devices,” IEEE Transactions on Information Forensics and Security, vol. 21, p. 1497–1511, 2026. [3] A. V. Guglielmi, L. Crosara, and S. Tomasin, “Design of physical layer challenge response authentication with RIS and multiple-antenna devices,” IEEE Transactions on Information Forensics and Security, vol. 21, p. 3677–3688, 2026. [4] Y. Zhang, Y. Dang, M. Kaveh, Z. Yan, R. Jantti, and Z. Han, “Physical layer challenge-response authentication between ambient backscatter devices,” IEEE Transactions on Information Forensics and Security, vol. 21, p. 1512–1527, 2026. [5] Z. Chang, W. Xiong, N. Xie, and D. Niyato, “Improving wireless se- curity with phase-tag physical-layer authentication,” IEEE Transactions on Mobile Computing, vol. 24, no. 6, p. 4716–4732, 2025. [6] Y. Zhou, Y. Wang, Z. Ma, P. Fan, and M. Xiao, “Physical layer authen- tication for UAV communications under rayleigh and rician channels,” IEEE Transactions on Wireless Communications, vol. 24, no. 4, p. 2722–2733, 2025. [7] X. Cheng, R. Meng, X. Xu, H. Gao, P. Zhang, and D. Niyato, “APEG: Adaptive physical layer authentication with channel extrapolation and generative AI,” IEEE Transactions on Information Forensics and Secu- rity, vol. 21, p. 1257–1272, 2026. [8] G. Han, Y. Hu, Y. He, T. Zhang, and F. Li, “Multi-source trust evaluation using physical layer authentication and reinforcement learning for dis- tributed AUV swarms in underwater data collection,” IEEE Transactions on Cognitive Communications and Networking, 2026. [9] N. Xie, J. Zhang, Q. Zhang, H. Tan, A. X. Liu, and D. Niyato, “Hybrid physical-layer authentication,” IEEE Transactions on Mobile Computing, vol. 23, no. 2, p. 1295–1311, 2024. [10] O. T. Demir and E. Bjornson, “Impact of phase noise and power amplifier non-linearities on downlink cell-free massive mimo-ofdm systems,” in 2025 International Conference on Smart Applications, Communications and Networking (SmartNets), 2025, p. 1–6. [11] T. Zhang, Y. Huo, L. Ma, and E. Guo, “Distributed physical layer authentication: Overview and opportunities,” IEEE Communications Magazine, vol. 62, no. 1, p. 134–140, 2024. [12] P. Zhang, K. Han, Y. Zhang, Y. Shen, F. Xiao, and X. Jiang, “Distributed physical layer authentication framework exploiting array pattern feature for mmwave mimo systems,” IEEE Transactions on Mobile Computing, vol. 24, no. 7, p. 6430–6445, 2025. [13] Y. Li, Y. Huo, T. Zhang, Z. Zhou, Q. Gao, T. Yan, Y. Yang, and T. Jing, “Distributed physical layer authentication with dynamic soft voting for smart distribution grids,” IEEE Transactions on Information Forensics and Security, vol. 20, p. 1807–1821, 2025. [14] H. Forssell and R. Thobaben, “Worst-case detection performance for distributed simo physical layer authentication,” IEEE Transactions on Communications, vol. 70, no. 1, p. 485–499, 2022. [15] J. Ma, J. Zhang, G. Shen, A. Marshall, and C.-H. Chang, “Adversarial attacks against deep learning-based radio frequency fingerprint identi- fication,” IEEE Transactions on Mobile Computing, vol. 25, no. 6, p. 7831–7844, 2026. [16] L. Yang, Q. Li, L. Wen, and J. Yang, “Exploiting radio frequency fin- gerprints for device identification: Tackling cross-receiver challenges in the source-data-free scenario,” IEEE Transactions on Mobile Computing, vol. 25, no. 6, p. 7800–7813, 2026. [17] P. Tang, G. Ding, Y. Xu, Y. Jiao, Y. Song, and G. Wei, “Causal learning for robust specific emitter identification over unknown channel statistics,” IEEE Transactions on Information Forensics and Security, vol. 19, p. 5316–5329, 2024. [18] S. Job, X. Tao, T. Cai, H. Xie, and J. Yong, “Causal neighbour- hood learning for invariant graph representations,” arXiv preprint arXiv:2602.17934, 2026. [19] Y. Wang and Z. Luo, “Causal discovery and causal inference based coun- terfactual fairness in machine learning,” in ICASSP 2023 - 2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2023, p. 1–5. [20] A. Roy, S. Banerjee, J. Sadasivan, A. Sarkar, and S. Dey, “Causally- aware reinforcement learning for joint communication and sensing,” IEEE Transactions on Machine Learning in Communications and Net- working, vol. 3, p. 552–567, 2025. [21] M. Chen, X. Zhang, D. Zhang, and J. Zhang, “A multi-agent policy op- timization algorithm integrating causal inference,” in 2025 International Conference on Computer, Internet of Things and Smart City (CIoTSC), 2025, p. 1–5. [22] J. Serra, J. Moysen, J. J. Hernandez-Carlon, and P. Dini, “Causal reinforcement learning for generalizable radio resource management in 6g,” IEEE Communications Magazine, 2026. [23] M. Arana-Catania, A. Sonee, A.-M. Khan, K. Fatehi, Y. Tang, B. Jin, A. Soligo, D. Boyle, R. Calinescu, P. Yadav, H. Ahmadi, A. Tsourdos, W. Guo, and A. Russo, “Explainable reinforcement and causal learning for improving trust to 6g stakeholders,” IEEE Open Journal of the Communications Society, vol. 6, p. 4101–4125, 2025. [24] L. Mata, M. Sousa, P. Vieira, M. Paula Queluz, and A. Rodrigues, “Integrating machine learning and observational causal inference for enhanced spectral and energy efficiency in wireless networks,” IEEE Access, vol. 13, p. 142 603–142 622, 2025. [25] C. Kurisummoottil Thomas, W. Saad, and Y. Xiao, “Causal semantic communication for digital twins: A generalizable imitation learning approach,” IEEE Journal on Selected Areas in Information Theory, vol. 4, p. 698–717, 2023. [26] W. Zhao, Y. Huang, Y. Guo, Y. Sun, L. Chen, and Y. Chen, “Tempo- ral meta-metric learning for multiuser physical-layer authentication in mobile iiot,” IEEE Internet of Things Journal, vol. 13, no. 4, p. 6755– 6765, 2026. 15 [27] Y. Zhou, L. Lei, X. Zhao, L. You, Y. Sun, and S. Chatzinotas, “Decompo- sition and meta-drl based multi-objective optimization for asynchronous federated learning in 6g-satellite systems,” IEEE Journal on Selected Areas in Communications, vol. 42, no. 5, p. 1115–1129, 2024. [28] W. Yan, F. Yin, J. Gao, A. Wang, Y. Tian, and R. Chen, “Attentional graph meta-learning for indoor localization using extremely sparse fingerprints,” IEEE Transactions on Mobile Computing, vol. 25, no. 6, p. 8718–8734, 2026. [29] A. K. Kocharlakota, S. A. Vorobyov, and R. W. Heath, “Pilot contam- ination aware transformer for downlink power control in cell-free mas- sive mimo networks,” IEEE Transactions on Wireless Communications, vol. 25, p. 9656–9671, 2026. [30] X. Zeng, C. Wang, Z. Li, L. Jin, D. W. K. Ng, D. Niyato, K. J. Kim, and N. Al-Dhahir, “Intelligent physical layer authentication based on complex-valued neural networks: Defending against pilot contamination and clone attacks,” IEEE Transactions on Wireless Communications, vol. 25, p. 14 486–14 502, 2026. [31] Y. Wang, Y. Wang, T. Wang, and J. Cheng, “Joint channel estimation, user activity identification, and pilot contamination attack detection for mmwave grant-free massive mtc networks: A three-dimensional compressive sensing-based approach,” IEEE Transactions on Commu- nications, vol. 73, no. 4, p. 2584–2601, 2025.