Resource Allocation in UAV-D2D Networks: A Scalable Heterogeneous Multi-Agent Deep Reinforcement Learning Approach

Abstract

In unmanned aerial vehicle (UAV)-assisted device-to-device (D2D) caching networks, the uncertainty from unpredictable content demands and variable user positions poses a significant challenge for traditional optimization methods, often making them impractical. Multi-agent deep reinforcement learning (MADRL) offers significant advantages in optimizing multi-agent system decisions and serves as an effective and practical alternative. However, its application in large-scale dynamic environments is severely limited by the curse of dimensionality and communication overhead. To resolve this problem, we develop a scalable heterogeneous multi-agent mean-field actor-critic (SH-MAMFAC) framework. The framework treats ground users (GUs) and UAVs as distinct agents and designs cooperative rewards to convert the resource allocation problem into a fully cooperative game, enhancing global network performance. We also implement a mixed-action mapping strategy to handle discrete and continuous action spaces. A mean-field MADRL framework is introduced to minimize individual agent training loads while enhancing total cache hit probability (CHP). The simulation results show that our algorithm improves CHP and reduces transmission delay. A comparative analysis with existing mainstream deep reinforcement learning (DRL) algorithms shows that SH-MAMFAC significantly reduces training time and maintains high CHP as GU count grows. Additionally, by comparing with SH-MAMFAC variants that do not include trajectory optimization or power control, the proposed joint design scheme significantly reduces transmission delay.

1. Introduction

Recently, unmanned aerial vehicles (UAVs) have gained popularity across various industries because of their compact size, affordability, and high flexibility [1]. These attributes enable UAVs to effectively address challenges in traditional communications, such as high implementation costs and limited flexibility in specialized scenarios [2]. In high-density hotspot areas, adding fixed ground nodes is not cost-effective. In such cases, UAV-assisted communication demonstrates high flexibility and ease of deployment [3]. Based on simulations and measurements, LTE networks can support UAV communication. Additionally, 5G massive multiple-input multiple-output (MIMO) antennas enable cellular UAVs to coexist with ground users (GUs), offering significant potential [4]. Edge caching is proposed as a key technology to alleviate traffic load in cellular networks [5]. By storing popular content closer to GUs, such as at base stations (BSs) and UAVs, it reduces content transmission delay and alleviates backhaul traffic load. Additionally, leveraging user-side cache capacity can further alleviate backhaul bottlenecks, for instance, through the implementation of device-to-device (D2D) communication [6]. In UAV-assisted D2D caching networks, collaboration between UAVs and GUs optimizes the use of limited cache space and improves system performance [7].

Due to the highly dynamic nature of UAV-assisted communication systems, achieving real-time accuracy in system modeling is challenging. Traditional optimization methods often convert the original non-convex problem into a series of convex subproblems, which are then solved iteratively [8]. While this approach enables tracking of the original problem, it often sacrifices accuracy. With an increase in the number of UAVs and GUs, computational complexity grows exponentially. Additionally, solving these problems optimally demands global information, which is not always possible in dynamic environments [8]. Recently, deep reinforcement learning (DRL) has garnered significant attention as a data-driven method. Due to its model-free nature and lack of need for real-time or complete prior information, DRL is considered an effective approach for managing task scheduling and wireless resource allocation in rapidly shifting UAV-assisted communication networks [9,10]. Compared to single-agent DRL, multi-agent systems offer several advantages: agents can learn strategies in a distributed manner through local environments, share experiences in a customized way, and maintain robustness by allowing the remaining agents to take over tasks if some agents fail [11].

Recent advancements in UAV-assisted communication have spurred numerous studies focusing on optimizing UAV trajectories, user association, and resource allocation to enhance quality-of-service (QoS) for GUs. To facilitate a better understanding of the methodologies and findings in the existing literature,

Table 1. Summary of existing works.

Work	Trajectory Design	Heterogeneous Services	Mixed Actions	Joint Optimization	Scalability
[8]	✓	✕	✓	✓	✕
[10]	✓	✕	✕	✓	✕
[16]	✓	✕	✓	✓	✕
[17]	✓	✓	✓	✓	✕
[18]	✕	✓	✓	✓	✕
[13]	✓	✕	✕	✕	✕
[14]	✓	✕	✕	✓	✕
[15]	✓	✕	✕	✓	✕
[19]	✕	✓	✕	✓	✕
Ours	✓	✓	✓	✓	✓

presents a comparative analysis of the related works. To fully leverage the high mobility of UAVs, some studies focus on trajectory design to enhance communication performance. For instance, Ref. [12] proposed a DRL solution for path planning in cellular-connected UAVs, leveraging a quantum-inspired experience replay framework to optimize flight direction while minimizing time costs and expected outage durations. Ref. [13] proposed a multi-agent deep reinforcement learning (MADRL) framework for optimizing the 3D trajectory design of cache-enabled UAVs in a dynamic wireless caching network, significantly enhancing network throughput by enabling UAVs to cooperatively learn from their experiences. On the other hand, several studies focused on the joint optimization of trajectory and resource allocation. Ref. [10] formulated a joint optimization problem involving trajectory design, user association, and power allocation, utilizing a strategic resource allocation algorithm that combined reinforcement learning and deep learning to maximize system utility for all served GUs. In [14], they jointly optimized user association, power allocation of non-orthogonal multiple access (NOMA), UAV deployment and caching placement, utilizing branch and bound (BaB) and DRL to minimize transmission delay in UAV-assisted caching networks for mixed augmented reality (AR) and multimedia applications, demonstrating lower long-term transmission delay in dynamic environments. Ref. [15] proposed a joint trajectory planning and cache management framework for UAV-assisted content delivery in intelligent transportation systems, addressing the challenges of dynamic vehicular networks through a proximal policy optimization (PPO) approach. Ref. [16] proposed a multi-UAV cooperative caching network optimization method supporting NOMA, which decoupled the problem into caching decisions and trajectory design, significantly improving cache hit probability (CHP) using the matching-DRL algorithm. However, [10,13,14] categorized continuous values into different levels of discrete values, while [15] only involved a discrete action space, limiting their ability to adapt to mixed action spaces. Furthermore, in these jobs GUs rely solely on UAVs to provide services.

There is another line of works that have explored collaboration between UAVs and ground nodes. Ref. [8] proposed a joint optimization framework for UAV trajectory design and GU access control, utilizing the MADRL algorithm to enhance throughput and fairness in air–ground coordinated communication systems. Ref. [18] employed a multi-agent deep deterministic policy gradient (MADDPG) algorithm for resource management in UAV-assisted vehicular networks. In this framework, the macro base station (MBS) and UAV collaboratively make association decisions and allocate appropriate amounts of resources to vehicles to meet heterogeneous QoS requirements while maximizing task offloading efficiency. In [17], they proposed a framework that treats the MBS and UAV as cooperative agents for content transmission. In dense multiple access cellular networks, this framework utilizes the Dual-Clip PPO algorithm to minimize transmission delay by optimizing multi-user association, cache placement, UAV trajectory, and transmission power. Ref. [19] proposed a cooperative multi-agent actor–critic (MAAC) algorithm that enables UAVs and BSs to collaborate on caching decisions, considering user mobility to maximize CHP. However, as the number of agents in the system increases, current MADRL methods face the curse of dimensionality. The cache and communication resource allocation problem is characterized by high-dimensional parameters, resulting in excessive data processing requirements and complicating the task of finding efficient configurations [20]. The algorithms proposed in [8,10,13,14,15,16,17,18,19] face scalability issues due to the exponential growth in the joint action space as the number of agents increases.

Existing studies have not simultaneously considered heterogeneous services, mixed action spaces, and system scalability, which provides a significant motivation for this research. The network dynamics and node heterogeneity complicate the analysis of communication performance for cache offloading requests. Thus, it is essential to develop a communication model for content distribution among different types of cache nodes and analyze the CHP of various communication modes. To maximize the CHP, we need to develop a multi-agent framework that includes both GUs and UAVs, allowing them to coexist and interact, which differs from [10,13,14,15,16] that considered only homogeneous agents. Typically, resource management in UAV-assisted communication systems poses a mixed-integer problem [11], requiring DRL algorithms to handle mixed action spaces. Furthermore, we need to address the scalability issue of MADRL as the number of agents increases in large-scale environments. Inspired by [21], we ensure effective learning by having each agent interact directly with a limited set of neighboring agents. The key innovations and contributions of this paper are summarized as follows:

(1) We propose a UAV-assisted heterogeneous D2D caching network, jointly optimizing caching strategies, flight trajectories, and transmission power to achieve high CHP and reduce transmission delay. Taking into account the channel model and GU mobility, we derive the successful transmission probability (STP) and CHP.

(2) We introduce a scalable heterogeneous multi-agent mean-field actor–critic (SH-MAMFAC) framework to maximize total CHP and reduce training time. This framework includes homogeneous caching users (CUs) and heterogeneous UAVs, effectively addressing the exponential increase in agent interactions in traditional MADRL and enhancing scalability.

(3) We assess the performance of different algorithms by analyzing CHP and transmission delay across various system parameters. Our comprehensive quantitative comparison includes SH-MAMFAC, MADDPG-based algorithms, and traditional centralized and independent DRL algorithms. In addition, the proposed joint design scheme significantly reduces transmission delay compared to variants without trajectory optimization or power control.

The paper is organized as follows: Section 2 introduces UAV-assisted D2D caching networks, deriving the STP and CHP. In Section 3, we define the CHP maximization problem and propose the SH-MAMFAC algorithm to address it. Section 4 covers the presentation and analysis of simulation results, and Section 5 offers conclusions. For convenience, the commonly used abbreviations in this paper are compiled in

Table A1. Summary of acronyms.

Acronym	Description
AR	Augmented reality
BaB	Branch and bound
BS	Base station
CDF	Cumulative distribution function
CHP	Cache hit probability
CU/RU	Caching user/requesting user
CTDE	Centralized training with decentralized execution
D2D	Device to device
DRL	Deep reinforcement learning
LoS/NLoS	Line of sight/non-line of sight
MADDPG	Multi-agent deep deterministic policy gradient
MADRL	Multi-agent deep reinforcement learning
MIMO	Multiple input multiple output
MBS	Macro base station
NP-SH-MAMFAC	No-power scalable heterogeneous multi-agent mean-field actor–critic
NPNT-SH-MAMFAC	No-power and no-trajectory scalable heterogeneous multi-agent mean-field actor–critic
NT-SH-MAMFAC	No-trajectory scalable heterogeneous multi-agent mean-field actor–critic
PPO	Proximal policy optimization
SH-MAMFAC	Scalable heterogeneous multi-agent mean-field actor–critic
SI-DRL	Scalable independent deep reinforcement learning
STP	Successful transmission probability
TDMA/NOMA	Time-division multiple access/non-orthogonal multiple access
UH-MADDPG	Unscalable heterogeneous multi-agent deep deterministic policy gradient
UC-DRL	Unscalable centralized deep reinforcement learning

, while the primary symbols and variables are summarized in

Table 2. List of main symbols and variables.

Parameter	Description
I	Number of CUs.
J	Number of RUs.
F	Number of files.
$\mathbf{Q}$	UAV flight trajectory $\mathbf{Q}=\left[\mathbf{q}\left(1\right),\cdots ,\mathbf{q}\left(N\right)\right]$.
${\mathbf{C}}_{\mathrm{d}}$	Cache placement matrix of CUs.
${\mathbf{c}}_{\mathrm{u}}$	Cache placement vector of UAV.
$\mathbf{p}$	Transmission power vector of CUs.
$r_{i,j}$	Data rate from CU i to RU j.
${\mathrm{Pr}}_{i,j}$	Successful transmission probability from CU i to RU j.
$r_{\mathrm{u},j}$	Data rate from UAV to RU j.
$R_{\mathrm{d}}$	Maximum communication range for D2D transmission.
$R_{\mathrm{u}}$	Maximum communication range for UAV transmission.
${\mathrm{Pr}}_{\mathrm{LoS},\mathrm{u},j}$	LoS channel probability.
$N_0$	Noise spectral density.
$re{q}_{j,f}$	The probability of RU j requesting file f.
${\mathcal{G}}_{j,f}$	The set of CUs that meet the conditions for establishing D2D links.
H	The total cache hit probability.
$T_j$	The total transmission delay for RU j.

.

2. System Model

Figure 1 illustrates the network structure, comprising an MBS, multiple GUs, and a UAV. The GUs are categorized into CUs and RUs. The maximum communication distance between D2D devices is $R_{\mathrm{d}}$, allowing CUs to transmit files to RUs only when the distance is less than $R_{\mathrm{d}}$. There are I CUs, denoted as $\mathcal{I}=\{1,\cdots ,I\}$, each with a cache size of $S_i$ Mbits. J RUs are denoted as $\mathcal{J}=\left\{1,\cdots ,J\right\}$. The file library $\mathcal{F}$ contains F files, each of size $S_f$ Mbits.

Content popularity can be modeled as Zipf distributions, as demonstrated in various real datasets [22]. Although content popularity follows a Zipf distribution, requests exhibit significant spatio-temporal variation due to differing preferences among GUs [23]. The probability that RU j requests file f in time slot t is [23]

$$p_{j,f}\left(t\right)=p_f{\displaystyle \frac{g\left({\varphi }_{j,t},{\vartheta }_f\right)}{{\sum }_{j^{\prime }\in \mathcal{J}}g\left({\varphi }_{j^{\prime },t},{\vartheta }_f\right)}},$$

(1)

where $p_f=(1/f^{\xi })/\left({\sum }_{f\in \mathcal{F}}(1/f^{\xi })\right)$ represents the popularity of file f, and $\xi$ is the Zipf slope parameter. $g\left({\varphi }_{j,t},{\vartheta }_f\right)$ is a kernel function representing the correlation between RU j and file f, with features ${\varphi }_{j,t}$ and ${\vartheta }_f$. The request variable $re{q}_{j,f}\left(t\right)$ is defined as $re{q}_{j,f}\left(t\right)=1$ if $p_{j,f}\left(t\right)>\epsilon$ and $re{q}_{j,f}\left(t\right)=0$ otherwise, where $\epsilon$ is the file popularity threshold.

2.1. UAV Mobility Model

Assuming the UAV operates at an altitude of H with a maximum communication range of $R_{\mathrm{u}}$, it employs a time-division multiple access (TDMA) protocol to transmit content to RUs within a fixed duration of T seconds [24]. The UAV’s cache storage capacity is $S_{\mathrm{u}}$ Mbits. For clarity, we use a Cartesian coordinate system with the MBS at ${\mathbf{w}}_{\mathrm{m}}=\left[x_{\mathrm{m}},y_{\mathrm{m}},h_{\mathrm{m}}\right]$. At time slot t, the position of CU i is represented by ${\mathbf{w}}_i\left(t\right)={\left[x_i\left(t\right),y_i\left(t\right)\right]}_{i\in \mathcal{I}}$, while the position of RU j is denoted as ${\mathbf{w}}_j\left(t\right)={\left[x_j\left(t\right),y_j\left(t\right)\right]}_{j\in \mathcal{J}}$. The duration T is divided into N time slots of ${\delta }_t=T/N$, allowing the position of the UAV to remain approximately stable throughout each slot. The flight trajectory is represented as $\mathbf{Q}=\left[\mathbf{q}\left(1\right),\cdots ,\mathbf{q}\left(N\right)\right]$, where $\mathbf{q}\left(t\right)={\left[x\left(t\right),y\left(t\right)\right]}^{\mathrm{T}}$ denotes the discrete-time position at time slot t.

At the beginning of time slot t, the UAV flies in a horizontal direction, determined by the angle ${\theta }_{\mathrm{u}}\left(t\right)\in [0,2\pi ]$ with a distance of $d_{\mathrm{u}}\left(t\right)\in [0,d_{max}]$, where $d_{max}$ represents the maximum distance the UAV can travel in one time slot. The UAV’s coordinates at time slot t are represented as $x\left(t\right)=x\left(0\right)+{\sum }_{t=0}^t{d}_{\mathrm{u}}\left(t\right)cos\left({\theta }_{\mathrm{u}}\left(t\right)\right)$ and $y\left(t\right)=y\left(0\right)+{\sum }_{t=0}^t{d}_{\mathrm{u}}\left(t\right)sin\left({\theta }_{\mathrm{u}}\left(t\right)\right)$. The UAV must remain within the boundary of the target area defined by $[0,x_{max}]\times [0,y_{max}]$. Additionally, to provide reliable service to RUs, the UAV needs to return to its initial position by the end of each cycle T. This results in the following movement constraints for the UAV [17,24,25]:

$$\begin{array}{cc}\hfill & 0⩽{\theta }_{\mathrm{u}}\left(t\right)⩽2\pi ,\hfill \end{array}$$

(2a)

$$\begin{array}{cc}\hfill & 0⩽x\left(t\right)⩽x_{max},\hfill \end{array}$$

(2b)

$$\begin{array}{cc}\hfill & 0⩽y\left(t\right)⩽y_{max},\hfill \end{array}$$

(2c)

$$\begin{array}{cc}\hfill & 0⩽d_{\mathrm{u}}\left(t\right)⩽d_{max},\hfill \end{array}$$

(2d)

$$\begin{array}{cc}\hfill & \mathbf{q}\left(1\right)=\mathbf{q}\left(N\right).\hfill \end{array}$$

(2e)

2.2. Communication Model

(1) D2D communication model: In contrast to UAV communication, ground channels between D2D users typically experience additional fading with random fluctuations [26]. Assuming the channel coefficients between D2D users remain constant within a time slot, they can vary from one slot to another. The average channel power gain from CU i to RU j during time slot t is ${\beta }_{i,j}\left(t\right)={\beta }_0{d}_{i,j}^{-{\alpha }_0}\left(t\right)={\beta }_0/(∥{\mathbf{w}}_i\left(t\right)-{\mathbf{w}}_j\left(t\right)∥{)}^{{\alpha }_0}$, where $d_{i,j}\left(t\right)$ denotes the distance from CU i to RU j, and ${\beta }_0$ represents the reference channel power gain at $d_0=1$ m. ${\alpha }_0$ is the path loss exponent for the D2D link, specifically set to 2 [27]. The channel gain from CU i to RU j during time slot t is $h_{i,j}\left(t\right)=\sqrt{{\beta }_{i,j}\left(t\right)}{g}_{i,j}\left(t\right)$, where $g_{i,j}\left(t\right)$ represents the small-scale channel fading between CU i and RU j, distributed as $\mathcal{C}\mathcal{N}(0,1)$ [26]. ${\left|h_{i,j}\left(t\right)\right|}^2$ follows an exponential distribution with a mean of ${\gamma }_{i,j}$.

During time slot t, the data rate from CU i to RU j is

$$r_{i,j}\left(t\right)=B_{\mathrm{d}}{log}_2\left(1+{\displaystyle \frac{P_i\left(t\right){\left|h_{i,j}\left(t\right)\right|}^2}{I_{i,j}\left(t\right)+N_0}}\right),$$

(3)

where $P_i\left(t\right)$ denotes the transmission power of CU i at time slot t, constrained by a peak value $P_{max}$ and an average value $\overline{P}$. $B_{\mathrm{d}}$ is the bandwidth of the CUs. $I_{i,j}\left(t\right)={\sum }_{i^{\prime }=1,i^{\prime }\ne i}^I{P}_{i^{\prime }}\left(t\right){\left|h_{i^{\prime },j}\left(t\right)\right|}^2$ denotes the interference caused by all CUs except CU i and follows an exponential distribution with a mean of ${\lambda }_{i,j}$. $N_0$ represents the noise power.

The cumulative distribution function (CDF) of $r_{i,j}\left(t\right)$ is

$$\begin{array}{cc}\hfill F_{r_{i,j}\left(t\right)}\left(r\right)& =\mathrm{Pr}\left\{r_{i,j}\left(t\right)\le r\right\}\hfill \\ \hfill & =\mathrm{Pr}\left\{B_{\mathrm{d}}{log}_2\left(1+{\displaystyle \frac{P_i\left(t\right){\left|h_{i,j}\left(t\right)\right|}^2}{I_{i,j}\left(t\right)+N_0}}\right)\le r\right\}\hfill \\ \hfill & =\mathrm{Pr}\left\{{\displaystyle \frac{P_i\left(t\right){\left|h_{i,j}\left(t\right)\right|}^2}{I_{i,j}\left(t\right)+N_0}}\le 2^{{\displaystyle \frac{r}{B_{\mathrm{d}}}}}-1\right\}.\hfill \end{array}$$

(4)

Due to GU mobility, node pairs needing data transmission may not always be in contact. We assume GUs enter and exit a finite area following a Poisson random process [28]. The contact duration $T_{i,j}$ between CU i and RU j is modeled as a continuous random variable with an exponential distribution and mean ${\mu }_{i,j}$ [29]. For CU i to establish a link with RU j, the probability of successfully transmitting a file of size $S_f$ Mbits within duration $T_{i,j}$ must exceed $v_{min}$. This probability, referred to as the STP, is denoted ${\mathrm{Pr}}_{i,j}\left(t\right)$.

Lemma 1.

Assuming $e_1$ and $e_2$ follow independent exponential distributions, $x={\varsigma }_1{e}_1$, $y={\varsigma }_2{e}_2+1$, and ${\varsigma }_1,{\varsigma }_2>0$, the CDF of the random variable $z=x/y$ is expressed as

$$F_z\left(z\right)=1-{\displaystyle \frac{{\varsigma }_1}{{\varsigma }_1+{\varsigma }_{2}z}}exp\left(-{\displaystyle \frac{z}{{\varsigma }_1}}\right).$$

(5)

Proof. 

The proof is presented in Appendix B.    □

Let $y={\displaystyle \frac{P_i\left(t\right){\left|h_{i,j}\left(t\right)\right|}^2}{I_{i,j}\left(t\right)+N_0}}$, ${\varsigma }_1={\displaystyle \frac{P_i\left(t\right){\gamma }_{i,j}}{N_0}}$, and ${\varsigma }_2={\displaystyle \frac{{\lambda }_{i,j}}{N_0}}$. Based on Lemma 1, the CDF of $r_{i,j}\left(t\right)$ is

$$\begin{array}{c}\hfill F_{r_{i,j}\left(t\right)}\left(r\right)=F_y\left(2^{{\displaystyle \frac{r}{B_{\mathrm{d}}}}}-1\right)=1-{\displaystyle \frac{{\varsigma }_1}{{\varsigma }_1+{\varsigma }_{2}y}}exp\left(-{\displaystyle \frac{2^{\frac{r}{B_{\mathrm{d}}}}-1}{{\varsigma }_1}}\right).\end{array}$$

(6)

Through the derivation presented in Appendix C, we can obtain the STP for the transmission of file f from CU i to RU j as follows:

$$\begin{array}{cc}\hfill {\mathrm{Pr}}_{i,j}\left(t\right)& =\mathrm{Pr}\left\{T_{i,j}\ge {\displaystyle \frac{S_f}{r_{i,j}\left(t\right)}}\right\}\hfill \\ \hfill & ={\int }_0^{+\infty }{\displaystyle \frac{{\varsigma }_1}{{\mu }_{i,j}·\left({\varsigma }_1+{\varsigma }_2\left(2^{\frac{S_f}{t·B_{\mathrm{d}}}}-1\right)\right)}}\times exp\left(-\left({\displaystyle \frac{2^{\frac{S_f}{t·B_{\mathrm{d}}}}-1}{{\varsigma }_1}}+{\displaystyle \frac{t}{{\mu }_{i,j}}}\right)\right)dt.\hfill \end{array}$$

(7)

(2) UAV communication model: Due to the altitude of the UAV and the complexity of the environment, the wireless channels between the UAV and GUs, as well as between the MBS and the UAV, are subject to probabilistic line-of-sight (LoS) and non-line-of-sight (NLoS) conditions [30]. The path loss from the UAV to RU j can be described by both LoS and NLoS models as follows [31]:

$$\begin{array}{cc}\hfill & {\mathrm{PL}}_{\mathrm{LoS},\mathrm{u},j}\left(t\right)=20log\left({\displaystyle \frac{4\pi f_c{d}_{\mathrm{u},j}\left(t\right)}{c}}\right)+{\eta }_{\mathrm{LoS}},\hfill \end{array}$$

(8a)

$$\begin{array}{cc}\hfill & {\mathrm{PL}}_{\mathrm{NLoS},\mathrm{u},j}\left(t\right)=20log\left({\displaystyle \frac{4\pi f_c{d}_{\mathrm{u},j}\left(t\right)}{c}}\right)+{\eta }_{\mathrm{NLoS}},\hfill \end{array}$$

(8b)

where $d_{\mathrm{u},j}\left(t\right)$ represents the distance between the UAV and RU j at time slot t. The carrier frequency is denoted by $f_c$, and c represents the speed of light. Additionally, ${\eta }_{\mathrm{LoS}}$ and ${\eta }_{\mathrm{NLoS}}$ represent the additional losses for the LoS and NLoS scenarios, respectively.

The probability of a LoS link between the UAV and RU j is determined by their elevation angles and environmental factors [17] and can be expressed as

$${\mathrm{Pr}}_{\mathrm{LoS},\mathrm{u},j}\left(t\right)={\displaystyle \frac{1}{1+c_{1}exp\left[-c_2\left({\theta }_{\mathrm{u},j}\left(t\right)-c_1\right)\right]}},$$

(9)

where ${\theta }_{\mathrm{u},j}\left(t\right)$ represents the elevation angle between the UAV and RU j at time slot t, $c_1$ and $c_2$ are environment-dependent constants. Then, the NLoS probability is ${\mathrm{Pr}}_{\mathrm{NLoS},\mathrm{u},j}\left(t\right)=1-{\mathrm{Pr}}_{\mathrm{LoS},\mathrm{u},j}\left(t\right)$. Therefore, the average path loss from the UAV and RU j can is

$${\mathrm{PL}}_{\mathrm{u},j}\left(t\right)={\mathrm{Pr}}_{\mathrm{LoS},\mathrm{u},j}\left(t\right)\times {\mathrm{PL}}_{\mathrm{LoS},\mathrm{u},j}\left(t\right)+{\mathrm{Pr}}_{\mathrm{NLoS},\mathrm{u},j}\left(t\right)\times {\mathrm{PL}}_{\mathrm{NLoS},\mathrm{u},j}\left(t\right).$$

(10)

If the distance is sufficiently large, interference can be neglected in millimeter-wave UAV networks [13]. During time slot t, the data rate from the UAV to RU j is given by

$$r_{\mathrm{u},j}\left(t\right)=B_{\mathrm{u}}{log}_2\left(1+{\displaystyle \frac{P_{\mathrm{u}}{\left|g_{\mathrm{u},j}\right|}^2{\mathrm{PL}}_{\mathrm{u},j}\left(t\right)}{N_0}}\right),$$

(11)

where $B_{\mathrm{u}}$ denotes the bandwidth of the UAV. $P_{\mathrm{u}}$ denotes the transmission power of the UAV. ${\left|g_{\mathrm{u},j}\left(t\right)\right|}^2$ represents the Rayleigh channel gain from the UAV to RU j and follows an exponential distribution with unit mean [31].

During time slot t, the data rate from the MBS to the UAV is

$$r_{\mathrm{m},\mathrm{u}}\left(t\right)=B_{\mathrm{m}}{log}_2\left(1+{\displaystyle \frac{P_{\mathrm{m}}{\left|g_{\mathrm{m},\mathrm{u}}\right|}^2{\mathrm{PL}}_{\mathrm{m},\mathrm{u}}\left(t\right)}{N_0}}\right),$$

(12)

where $B_{\mathrm{m}}$ represents the bandwidth of the MBS. $P_{\mathrm{m}}$ represents the transmission power of the MBS. The small-scale fading between the MBS and the UAV is represented by ${\left|g_{\mathrm{m},\mathrm{u}}\left(t\right)\right|}^2$, while the path loss in this link is indicated by ${\mathrm{PL}}_{\mathrm{m},\mathrm{u}}\left(t\right)$.

2.3. Transmission Delay Model

The cache placement for CU i and the UAV at time slot t is represented by the vectors ${\mathbf{c}}_i\left(t\right)$ and ${\mathbf{c}}_{\mathrm{u}}\left(t\right)$, respectively. If a file is cached at a CU or the UAV, the corresponding entry in the cache placement vector is set to 1; otherwise, it is 0. When RU j requests file f at time slot t, it first checks if any CU can establish a D2D link. The set of CUs that meet the D2D link conditions is denoted as ${\mathcal{G}}_{j,f}\left(t\right)=\{i\mid {\mathrm{Pr}}_{i,j}\left(t\right)\ge v_{min},c_{i,f}\left(t\right)=1,d_{i,j}\left(t\right)\le R_{\mathrm{d}},i\in \mathcal{I}\}$. If no CU satisfies these conditions, the file is requested from the UAV.

The downlink wireless transmission delay from CU i to RU j for transmitting file f is

$$T_{i,j,f}\left(t\right)={\displaystyle \frac{S_f}{r_{i,j}\left(t\right)}}.$$

(13)

For UAV communication, the transmission delay consists of two components: the downlink wireless delay and the backhaul link delay [2]. The downlink wireless delay from the UAV to RU j is

$$T_{1,\mathrm{u},j,f}\left(t\right)={\displaystyle \frac{S_f}{r_{\mathrm{u},j}\left(t\right)}}.$$

(14)

While the backhaul link delay from the MBS to the UAV for sending file f is

$$T_{2,\mathrm{u},j,f}\left(t\right)={\displaystyle \frac{S_f·\left(1-c_{\mathrm{u},f}\left(t\right)\right)}{r_{\mathrm{m},\mathrm{u}}\left(t\right)}}.$$

(15)

If the content requested by RU j is already stored on the UAV, the backhaul link will not be required, that is, when $c_{\mathrm{u},f}\left(t\right)=1$, $T_{2,\mathrm{u},j,f}\left(t\right)=0$.

Thus, the overall transmission delay for RU j in time slot t is

$$\begin{array}{c}\hfill T_j\left(t\right)=\sum _{f\in \mathcal{F}}\left[{\mathbb{I}}_{j,f}\left(t\right)·\underset{i\in {\mathcal{G}}_{j,f}}{min}\left(T_{i,j,f}\left(t\right)\right)+\left(1-{\mathbb{I}}_{j,f}\left(t\right)\right)·T_{\mathrm{u},j,f}\left(t\right)\right],\end{array}$$

(16)

where ${\mathbb{I}}_{j,f}\left(t\right)=1$ indicates that file f is available from CUs via D2D, while ${\mathbb{I}}_{j,f}\left(t\right)=0$ means the file is requested from the UAV, resulting in $T_{\mathrm{u},j,f}\left(t\right)=T_{1,\mathrm{u},j,f}\left(t\right)+T_{2,\mathrm{u},j,f}\left(t\right)$.

2.4. UAV Energy Consumption Model

In this work, the energy consumption of the UAV includes two components: the energy required for propulsion during movement and the energy needed for content transmission.

The power required for a rotary-wing UAV to travel at a speed V is [32]:

$$\begin{array}{cc}\hfill P\left(V\right)=& P_0\left(1+{\displaystyle \frac{3V^3}{U^2}}\right)+P_1{\left({\left(1+{\displaystyle \frac{V^4}{4v_r^4}}\right)}^{{\displaystyle \frac{1}{2}}}-{\displaystyle \frac{V^2}{2v_r^2}}\right)}^{{\displaystyle \frac{1}{2}}}+{\displaystyle \frac{1}{2}}{d}_r\rho sA{V}^3,\hfill \end{array}$$

(17)

where $P_0$ and $P_1$ are constants associated with the blade profile power and the induced power of the UAV while in a hovering state. U signifies the rotor blade tip speed, $v_r$ is the average rotor-induced velocity in hover, $d_r$ indicates the fuselage drag ratio, and $\rho$ denotes air density. Additionally, s and A correspond to the rotor solidity and rotor disc area, respectively [25]. The power required for hovering is $P_{\mathrm{h}}=P(V=0)=P_0+P_1$.

When the UAV travels at a constant speed V, the energy used for propulsion in one time slot can be represented as

$$E\left(x\right)={\displaystyle \frac{x}{V}}P\left(V\right)+max\left\{{\delta }_t-{\displaystyle \frac{x}{V}},0\right\}·P_{\mathrm{h}},$$

(18)

where $x=∥\mathbf{q}\left(t+1\right)-\mathbf{q}\left(t\right)∥$ represents the distance covered by the UAV in one time slot.

The communication energy consumption from the UAV to RUs for transmitting file f is given by [23]

$$E_{\mathrm{u}}\left(t\right)=\sum _{j=1}^J\sum _{f=1}^F{P}_{\mathrm{u}}·min\left\{{\delta }_t,T_{\mathrm{u},j,f}\left(t\right)\right\}.$$

(19)

2.5. Cache Hit Probability Model

In time slot t, the probability that RU j obtains file f from the CUs is

$$H_{\mathrm{d},j,f}\left(t\right)=re{q}_{j,f}\left(t\right)·{\displaystyle \frac{{\sum }_{i\in {\mathcal{G}}_{j,f}}{c}_{i,f}\left(t\right)·{\mathrm{Pr}}_{i,j}\left(t\right)}{{\sum }_{i\in {\mathcal{G}}_{j,f}}{c}_{i,f}\left(t\right)}}.$$

(20)

The total probability that RUs can retrieve the file from the CUs is

$$H_{\mathrm{d}}\left(t\right)={\displaystyle \frac{{\sum }_{j=1}^J{\sum }_{f=1}^F{\mathbb{I}}_{j,f}\left(t\right)·H_{\mathrm{d},j,f}\left(t\right)}{{\sum }_{j=1}^J{\sum }_{f=1}^{F}re{q}_{j,f}\left(t\right)}}.$$

(21)

In time slot t, the probability that RU j acquires file f from the UAV is

$$H_{\mathrm{u},j,f}\left(t\right)=re{q}_{j,f}\left(t\right)·c_{\mathrm{u},f}\left(t\right).$$

(22)

The total probability that RUs can retrieve the file from the UAV is

$$H_{\mathrm{u}}\left(t\right)={\displaystyle \frac{{\sum }_{j=1}^J{\sum }_{f=1}^F(1-{\mathbb{I}}_{j,f}\left(t\right))·H_{\mathrm{u},j,f}\left(t\right)}{{\sum }_{j=1}^J{\sum }_{f=1}^{F}re{q}_{j,f}\left(t\right)}}.$$

(23)

The total CHP at time slot t is defined as the probability that the content requested by RU is present in either the CU cache or the UAV cache, expressed as

$$H\left(t\right)=1-\left(1-H_{\mathrm{d}}\left(t\right)\right)·\left(1-H_{\mathrm{u}}\left(t\right)\right).$$

(24)

3. Problem Formulation and Optimization

3.1. Problem Formulation

We propose a joint optimization framework to maximize the overall CHP by simultaneously optimizing the UAV trajectory $\mathbf{q}\left(t\right)$, the CUs’ transmission power $\mathbf{p}\left(t\right)=[P_1\left(t\right),\cdots ,P_I\left(t\right)]$, the CUs’ cache placement ${\mathbf{C}}_{\mathrm{d}}\left(t\right)=[{\mathbf{c}}_1\left(t\right),\cdots ,{\mathbf{c}}_I\left(t\right)]$, and the UAV’s cache placement ${\mathbf{c}}_{\mathrm{u}}\left(t\right)$. The optimization problems are defined as follows:

$$\begin{array}{cc}\hfill {\mathit{P}}_1:& \underset{\mathbf{q}\left(t\right),\mathbf{p}\left(t\right),{\mathbf{C}}_{\mathrm{d}}\left(t\right),{\mathbf{c}}_{\mathrm{u}}\left(t\right)}{max}{\displaystyle \frac{1}{N}}\sum _{t=1}^{N}H\left(t\right)\hfill \end{array}$$

(25a)

$$\begin{array}{cc}\hfill \mathrm{s}.\mathrm{t}.& \sum _{f=1}^F{c}_{i,f}\left(t\right)·S_f⩽S_i,\forall i,t,\hfill \end{array}$$

(25b)

$$\begin{array}{cc}\hfill & \sum _{f=1}^F{c}_{\mathrm{u},f}\left(t\right)·S_f⩽S_{\mathrm{u}},\forall t,\hfill \end{array}$$

(25c)

$$\begin{array}{cc}\hfill & 0⩽{\theta }_{\mathrm{u}}\left(t\right)⩽2\pi ,\hfill \end{array}$$

(25d)

$$\begin{array}{cc}\hfill & 0⩽x\left(t\right)⩽x_{max},\hfill \end{array}$$

(25e)

$$\begin{array}{cc}\hfill & 0⩽y\left(t\right)⩽y_{max},\hfill \end{array}$$

(25f)

$$\begin{array}{cc}\hfill & 0⩽d_{\mathrm{u}}\left(t\right)⩽d_{max},\hfill \end{array}$$

(25g)

$$\begin{array}{cc}\hfill & \mathbf{q}\left(1\right)=\mathbf{q}\left(N\right),\hfill \end{array}$$

(25h)

$$\begin{array}{cc}\hfill & \sum _{t=1}^N\left(E\left(∥\mathbf{q}\left(t\right)-\mathbf{q}(t-1)∥\right)+E_{\mathrm{u}}\left(t\right)\right)\le E_{\mathrm{u},max},\hfill \end{array}$$

(25i)

$$\begin{array}{cc}\hfill & P_i\left(t\right)\le P_{max},\forall i,t,\hfill \end{array}$$

(25j)

$$\begin{array}{cc}\hfill & {\displaystyle \frac{1}{N}}\sum _{t=1}^N{P}_i\left(t\right)\le \overline{P},\forall i,\hfill \end{array}$$

(25k)

$$\begin{array}{cc}\hfill & T_j\left(t\right)\le {\delta }_t,\forall j,t,\hfill \end{array}$$

(25l)

where $E_{\mathrm{u},max}$ denotes the total onboard energy available to the UAV. Equations (25b) and (25c) limit the maximum cached size for the CUs and UAV, respectively. Equations (25d)–(25h) impose constraints on the UAV trajectory. Equation (25i) imposes a cap on the highest amount of energy that the UAV is allowed to use. Equations (25j) and (25k) affect the power allocation for the CU. Equation (25l) states that the total transmission delay within each time slot must not exceed the duration of the time slot.

${\mathit{P}}_1$ represents a non-convex mixed-integer programming challenge, categorized as NP-hard [33]. Typically, a brute-force method is employed to identify the optimal solution. However, this approach is computationally intensive and impractical for large-scale systems [10]. Additionally, obtaining the optimal solution for the next time slot requires full knowledge of future information, meaning performance may be compromised without it. Therefore, we consider using a DRL-based approach to estimate the optimal solution when this information is unavailable.

3.2. Partially Observable Markov Game

Due to the optimization goal of maximizing the total CHP, we convert the ${\mathit{P}}_1$ into a partially observable Markov game involving $I+1$ agents, comprising I CU agents and 1 UAV agent. This can be described by a tuple $\langle \mathcal{S},\mathcal{A},\mathcal{R},\mathcal{P},\gamma \rangle$, where $\mathcal{S}$ denotes the collection of states defining the environment; $\mathcal{A}=\left\{{\mathcal{A}}_1,\cdots ,{\mathcal{A}}_{I+1}\right\}$ denotes the joint action set encompassing all agents; $\mathcal{R}=\left\{{\mathcal{R}}_1,\cdots ,{\mathcal{R}}_{I+1}\right\}$ denotes the set of reward functions that assigns rewards according to states and the actions of the agents, i.e., $\mathcal{R}:\mathcal{S}\times \mathcal{A}\mapsto {\mathbb{R}}^{I+1}$; $\mathcal{P}$ denotes the state transition function, which specifies the probability of moving to the subsequent state given the joint actions taken by all agents in the current state, i.e., $\mathcal{P}:\mathcal{S}\times {\mathcal{A}}_1\times \cdots \times {\mathcal{A}}_{I+1}\times \mathcal{S}\mapsto [0,1]$; and $\gamma$ represents the factor for discounting future rewards.

During each time slot t, the agents are able to engage with the environment, where the state of the environment is represented by $s_t$, with $s_t\in \mathcal{S}$. Consequently, agent k (we use k as the index of agents in the entire heterogeneous system, where $k\in \{1,\cdots ,I+1\}$; when $k=I+1$, the agent corresponds to the UAV; otherwise, the agent is a CU) receives local observations of the environment, denoted as $o_{k,t}$, determined by observation functions $o_{k,t}=b_k\left(s_t\right)$. Subsequently, agent k takes actions $a_{k,t}={\pi }_k\left(o_{k,t}\right)$, where ${\pi }_k(·)$ represents the strategy functions of agent k. Let ${\mathbf{o}}_t$ represent the set of observations $\{o_{1,t},\cdots ,o_{I+1,t}\}$, and ${\mathbf{a}}_t$ represent the set of actions $\{a_{1,t},\cdots ,a_{I+1,t}\}$. As illustrated in Figure 2, each agent obtains local observations and selects actions according to its policy. After taking joint action ${\mathbf{a}}_t$, the agent obtains the associated reward ${\mathbf{r}}_t=\left\{r_{1,t},\cdots ,r_{I+1,t}\right\}$ according to the reward function $\mathcal{R}$. Furthermore, based on the state transition probability function $\mathcal{P}$, the joint observation updates from ${\mathbf{o}}_t$ to ${\mathbf{o}}_{t+1}$.

3.3. Markov Decision Process

The SH-MAMFAC framework introduces two types of heterogeneous multi-agents. The first type includes all homogeneous CUs.

(1) Observation space: The observation space for the CU agent includes the transmission power and caching decisions of the CU from the previous time slot, as well as its own position, the positions of all RUs (RUs equipped with GPS can report their positions to the MBS; the MBS then distributes this location information to all CUs [17]), and the request distribution for each RU during the present time slot. Therefore, the observation space for the CU is

$$o_{i,t}=\left\{P_i(t-1),{\mathbf{c}}_i(t-1),{\mathbf{w}}_i\left(t\right),{\left\{{\mathbf{w}}_j\left(t\right)\right\}}_{j\in \mathcal{J}},{\left\{p_{j,f}\left(t\right)\right\}}_{j\in \mathcal{J}\cup f\in \mathcal{F}}\right\}.$$

(26)

(2) Action space: In our scenario, the action space $a_{i,t}$ is a combination of discrete and continuous components. It includes discrete sub-actions related to caching decisions ${\mathbf{c}}_i\left(t\right)$ and continuous sub-actions for transmission power $P_i\left(t\right)$. In the SH-MAMFAC framework, the actor network produces a normalized continuous action space $a_{\mathrm{con},i,t}=\left\{{\mathbf{c}}_{i,t}^{\prime },P_{i,t}^{\prime }\right\}$, which requires further processing before it can be applied in the environment. Therefore, we design a mixed-action mapping strategy to address this issue. First, we transform the output ${\mathbf{c}}_{i,t}^{\prime }$ into a probability distribution ${\widehat{\mathbf{c}}}_{i,t}^{\prime }={\widehat{\pi }}_i\left(o_{i,t};{\widehat{\mathsf{\Theta }}}_{\pi ,i}\right)$, where the fth element of ${\widehat{\mathbf{c}}}_{i,t}^{\prime }$ represents the probability of CU i caching file f. ${\widehat{\pi }}_i$ is the actor network responsible for outputting the probabilities. We then sample from this probability distribution to obtain the actual cache decision ${\mathbf{c}}_i\left(t\right)=\psi \left({\widehat{\mathbf{c}}}_{i,t}^{\prime }\right)$, with $\psi (·)$ representing the sampling process. This approach effectively maps the generated cache actions to the discrete action space [8]. Additionally, $P_{i,t}^{\prime }$ is mapped to the real value by multiplying by its maximum value $P_{max}$. Thus, it effectively addresses the issue of hybrid action space.

The second category of agents comprises the UAV, which differs from the CUs in terms of heterogeneity:

(1) Observation space: The observation space for the UAV agent includes the UAV caching decisions from the previous time slot, as well as its own position, the positions of all RUs, and the request distribution for each RU during the present time slot. Therefore, the observation space for the UAV is

$$\begin{array}{c}\hfill o_{I+1}\left(t\right)=\left\{{\mathbf{c}}_{\mathrm{u}}(t-1),\mathbf{q}\left(t\right),{\left\{{\mathbf{w}}_j\left(t\right)\right\}}_{j\in \mathcal{J}},{\left\{p_{j,f}\left(t\right)\right\}}_{j\in \mathcal{J}\cup f\in \mathcal{F}}\right\}.\end{array}$$

(27)

(2) Action space: Similar to the action space for CU agents, the action space of the UAV agent $a_{I+1,t}$ includes an integer sub-action for the caching decision ${\mathbf{c}}_{\mathrm{u}}\left(t\right)$, as well as continuous sub-actions ${\theta }_{\mathrm{u}}\left(t\right)$ and $d_{\mathrm{u}}\left(t\right)$. We transform normalized continuous output ${\mathbf{c}}_{\mathrm{u},t}^{\prime }$ into a probability distribution ${\widehat{\mathbf{c}}}_{\mathrm{u},t}^{\prime }$ and then sample from it to obtain the actual cache decision ${\mathbf{c}}_{\mathrm{u}}\left(t\right)=\psi \left({\widehat{\mathbf{c}}}_{\mathrm{u},t}^{\prime }\right)$. The continuous sub-actions ${\theta }_{\mathrm{u},t}^{\prime }$ and $d_{\mathrm{u},t}^{\prime }$ are multiplied by their maximum values $2\pi$ and $d_{max}$, respectively, to ensure effective application in the actual operation.

Rewards guide each agent toward its optimal policy, and when the reward signals at each step are aligned with the desired objectives, system performance can be enhanced [34]. Given that shaped rewards can accelerate algorithm learning compared to sparse rewards, we designed the reward element $r_{i,t+1}={\displaystyle \frac{{\sum }_{j=1}^J{\sum }_{f=1}^F{\mathbb{I}}_{j,f}\left(t\right)·c_{i,f}\left(t\right)·{\mathrm{Pr}}_{i,f}\left(t\right)·re{q}_{j,f}\left(t\right)}{{\sum }_{j=1}^J{\sum }_{f=1}^{F}re{q}_{j,f}\left(t\right)}}$ for CU i, which can be derived from (20) and (21). Similarly, based on (23), the reward element for the UAV is $r_{I+1,t+1}=H_{\mathrm{u}}\left(t\right)$. Given the shared goal of the optimization problem, all $I+1$ agents must work together to optimize the total CHP. To achieve a cooperative Markovian game, we assume that each agent receives the same immediate reward [18,34], i.e., $r_{1,t}=\cdots =r_{I+1,t}=r_t$. To reduce significant fluctuations in rewards, we define $r_t={\displaystyle \frac{1}{J}}{\sum }_{j=1}^J{r}_{j,t}$, where $r_{j,t}$ represents the reward obtained by RU j during time slot t, and is defined as

$$r_j\left(t\right)=\left\{\begin{array}{cc}\sum _{f=1}^F{H}_{\mathrm{d},j,f}\left(t\right),\hfill & \mathrm{if}{\mathbb{I}}_{j,f}\left(t\right)=1,\hfill \\ \sum _{f=1}^F{H}_{\mathrm{u},j,f}\left(t\right),\hfill & \mathrm{if}{\mathbb{I}}_{j,f}\left(t\right)=0.\hfill \end{array}\right.$$

(28)

3.4. Scalable Heterogeneous Multi-Agent Mean-Field Actor–Critic

As illustrated in Figure 3, each agent utilizes an actor network ${\pi }_k\left(o_k;{\mathsf{\Theta }}_{\pi ,k}\right)$ with weights ${\mathsf{\Theta }}_{\pi ,k}$ for decentralized execution. The actor network processes local observations to generate actions as outputs. In conventional MADRL, the action–value function, i.e., the Q-function, calculates the expected discounted rewards for given states and the actions taken by the agents. Agents optimize their behavior by learning these Q-functions. However, learning Q-functions becomes exceedingly difficult when the dimension of joint actions grows significantly [21,35]. To tackle scalability challenges in large-scale MADRL, we approximate the standard Q-function of agent k with a mean-field Q-function as [36]

$$Q_k\left(\mathbf{o},\mathbf{a}\right)={\displaystyle \frac{1}{\left|{\mathcal{K}}_{k_1}\right|}}\sum _{k_1}{Q}_k\left(\mathbf{o},a_k,a_{1,k_1}\right)+{\displaystyle \frac{1}{\left|{\mathcal{K}}_{k_2}\right|}}\sum _{k_2}{Q}_k\left(\mathbf{o},a_k,a_{2,k_2}\right),$$

(29)

where $k_1\in {\mathcal{K}}_{k_1}$ and $k_2\in {\mathcal{K}}_{k_2}$ are the index sets of neighboring CU agents and UAV agents (in this scenario, we consider a single UAV; however, this approach can be easily extended to multiple UAVs, which will be covered in our future work) associated with agent k, respectively. ${\mathcal{K}}_k={\mathcal{K}}_{k_1}\cup {\mathcal{K}}_{k_2}$ represents the set of neighboring agents for agent k determined by specific application settings. And (29) can be rewritten as

$$Q_k\left(\mathbf{o},\mathbf{a}\right)={\displaystyle \frac{1}{\left|{\mathcal{K}}_k\right|}}\sum _k{Q}_k\left(\mathbf{o},a_k,a_{1,k_1},a_{2,k_2}\right).$$

(30)

It is important to note that the pairwise approximation between an agent and its connected neighbors preserves the global interaction between any two agents. Mean-field theory simplifies the interactions of multiple agents to those of two agents, where the second agent represents the average effect of all other agents [37]. Then, Theorem 1 is founded for the Q-function:

Theorem 1.

Q-function $Q_k\left(\mathbf{o},\mathbf{a}\right)$ of agent k can be approximated by mean action ${\overline{a}}_{l,k}$ and action $a_k$, i.e.,

$$Q_k\left(\mathbf{o},\mathbf{a}\right)={\displaystyle \frac{1}{\left|{\mathcal{K}}_k\right|}}\sum _k{Q}_k\left(\mathbf{o},a_k,a_{1,k_1},a_{2,k_2}\right)\approx Q_k\left(\mathbf{o},a_k,{\overline{a}}_{l,k}\right).$$

(31)

where ${\overline{a}}_{l,k}$ ($l\in \{1,2\}$) denotes the mean action of agents belonging to class l in the neighborhood of agent k.

Proof. 

The proof is presented in Appendix D.    □

The mean-field Q-function can be updated recursively as follows:

$$Q_{k,t+1}\left(\mathbf{o},a_k,{\overline{a}}_{l,k}\right)=\left(1-{\alpha }_t\right)Q_{k,t}\left(\mathbf{o},a_k,{\overline{a}}_{l,k}\right)+{\alpha }_t\left[r_k+\gamma v_{k,t}\left({\mathbf{o}}^{\prime }\right)\right],$$

(32)

where ${\alpha }_t$ is the learning rate. The value function $v_{k,t}\left({\mathbf{o}}^{\prime }\right)$ of agent k at time t is defined as

$$v_{k,t}\left({\mathbf{o}}^{\prime }\right)=\sum _{a_k}{\pi }_{k,t}\left(a_k|{\mathbf{o}}^{\prime },{\overline{a}}_{l,k}\right){\mathbb{E}}_{a_{l,-k}\sim {\pi }_{-k,t}}\left[Q_{k,t}\left(\mathbf{o},a_k,{\overline{a}}_{l,k}\right)\right].$$

(33)

We compute the mean actions of heterogeneous agents using the following equation:

$${\overline{a}}_{l,k}={\displaystyle \frac{1}{\left|{\mathcal{K}}_{k_l}\right|}}\sum _{k_l}{a}_{l,k_l},a_{l,k_l}\sim {\pi }_{k_l,t}\left(·\mid \mathbf{o},{\overline{a}}_{l-,k_l}\right),$$

(34)

where ${\pi }_{k_l,t}$ denotes the policy of agent $k_l$ (in the neighborhood of k), and ${\overline{a}}_{l-,k_l}$ represents the previous mean action of the neighbors of class l for agent k. Then, the Boltzmann policy for each agent k is

$${\pi }_{k,t}\left(a_k|\mathbf{o},{\overline{a}}_{l,k}\right)={\displaystyle \frac{exp\left(\beta Q_{k,t}\left(\mathbf{o},a_k,{\overline{a}}_{l,k}\right)\right)}{{\sum }_{a_{k^{\prime }}\in {\mathcal{A}}_k}exp\left(\beta Q_{k,t}\left(\mathbf{o},a_{k^{\prime }},{\overline{a}}_{l,k}\right)\right)}},$$

(35)

where $\beta$ denotes the Boltzmann softmax parameter.

By iterating Equations (33)–(35), the mean actions of all agents and their corresponding policies are continuously improved. In the proof of Theorem 2, we demonstrate that this method converges to a fixed point, which is boundedly close to the Nash equilibrium. We first present three assumptions regarding mean-field updates, followed by Theorem 2.

Assumption 1.

Each action–value pair is visited infinitely often, while the rewards are bounded by a constant.

Assumption 2.

The agent’s policy is known as Greedy in the Limit with Infinite Exploration (GLIE). With the Boltzmann policy, as the temperature decreases, the policy becomes greedy with respect to the Q-function in the limit.

Assumption 3.

For each stage game $[Q_{1,t}\left(\mathbf{o}\right),\cdots ,Q_{I+1,t}\left(\mathbf{o}\right)]$ at time t and state s. For all $t,s,k\in 1,\cdots ,I+1$, the Nash equilibrium ${\pi }_{\ast }=[{\pi }_{1,\ast },\cdots ,{\pi }_{I+1,\ast }]$ is regarded as the global optimal solution or saddle point, as shown below:

$$1.{\mathbb{E}}_{{\pi }_{\ast }}\left[Q_{k,t}\left(\mathbf{o}\right)\right]⩾{\mathbb{E}}_{\pi }\left[Q_{k,t}\left(\mathbf{o}\right)\right],\forall \pi \in \mathsf{\Omega }\left({\mathsf{\Pi }}_{k_l}{\mathfrak{A}}_{k_l}\right),$$

(36)

$$2.{\mathbb{E}}_{{\pi }_{\ast }}\left[Q_{k,t}\left(\mathbf{o}\right)\right]⩾{\mathbb{E}}_{{\pi }_k}{\mathbb{E}}_{{\pi }_{-k,\ast }}\left[Q_{k,t}\left(\mathbf{o}\right)\right],\forall {\pi }_k\in \mathsf{\Omega }\left({\mathfrak{A}}_k\right),$$

(37)

$$3.{\mathbb{E}}_{{\pi }_{\ast }}\left[Q_{k,t}\left(\mathbf{o}\right)\right]⩽{\mathbb{E}}_{{\pi }_{k,\ast }}{\mathbb{E}}_{{\pi }_{-k,\ast }}\left[Q_{k,t}\left(\mathbf{o}\right)\right],\forall {\pi }_{-k}\in \mathsf{\Omega }\left({\mathsf{\Pi }}_{k^l\ne k}{\mathfrak{A}}_{k_l}\right).$$

(38)

Theorem 2.

Based on (32)–(35), under the conditions of Assumptions 1–3, updating the mean-field Q-function will converge to a fixed point that is boundedly close to the Nash equilibrium.

Proof. 

The proof is presented in Appendix E.    □

During the centralized training process, each agent stores experience $\langle \mathbf{o},\mathbf{a},\mathbf{r},{\mathbf{o}}^{\prime },{\overline{\mathbf{a}}}_l\rangle$ in a replay buffer ${\mathcal{M}}_k$ and randomly samples mini-batches from it to update their parameters. This process aims to ensure system stability and improve sample efficiency. Consider the m-th transition as an example; the critic network samples a batch of ${\mathcal{B}}_k$ transitions from ${\mathcal{M}}_k$ and minimizes the loss function to adjust the parameters of the evaluation network using the learning rate ${\rho }_Q$:

$$L\left({\mathsf{\Theta }}_{Q,k}\right)={\displaystyle \frac{1}{\left|{\mathcal{B}}_k\right|}}\sum _m\left[{\left(y_k^m-Q_k\left({\mathbf{o}}^m,a_k^m,{\overline{a}}_{l,k}^m;{\mathsf{\Theta }}_{Q,k}\right)\right)}^2\right],$$

(39)

where $y_k^m=r_k^m+\gamma Q_k^{\prime }\left({\mathbf{o}}^{m^{\prime }},a_k^{m^{\prime }},{\overline{a}}_{l,k}^{m^{\prime }};{\mathsf{\Theta }}_{Q^{\prime },k}\right)$ represents the estimated target value using the target critic network. On the other hand, the parameters of the actor network are updated using the learning rate ${\varrho }_{\pi }$ by maximizing the policy objective function:

$${\nabla }_{{\mathsf{\Theta }}_{\pi ,k}}J\left({\mathsf{\Theta }}_{\pi ,k}\right)\approx {\displaystyle \frac{1}{\left|{\mathcal{B}}_k\right|}}\sum _m\left[{\nabla }_{{\mathsf{\Theta }}_{\pi ,k}}log{\pi }_k\left({\mathbf{o}}^m;{\mathsf{\Theta }}_{\pi ,k}\right)\left(y_k^m-Q_k\left({\mathbf{o}}^m,a_k^m,{\overline{a}}_{l,k}^m;{\mathsf{\Theta }}_{Q,k}\right)\right)\right].$$

(40)

As ${\mathsf{\Theta }}_{\pi ,k}$ and ${\mathsf{\Theta }}_{Q,k}$ are updated in real time, the parameters of the target networks ${\mathsf{\Theta }}_{{\pi }^{\prime },k}$ and ${\mathsf{\Theta }}_{Q^{\prime },k}$ can also be smoothly updated as follows:

$$\begin{array}{cc}\hfill {\mathsf{\Theta }}_{{\pi }^{\prime },k}& =\kappa {\mathsf{\Theta }}_{\pi ,k}+(1-\kappa ){\mathsf{\Theta }}_{{\pi }^{\prime },k},\hfill \end{array}$$

(41a)

$$\begin{array}{cc}\hfill {\mathsf{\Theta }}_{Q^{\prime },k}& =\kappa {\mathsf{\Theta }}_{Q,k}+(1-\kappa ){\mathsf{\Theta }}_{Q^{\prime },k}.\hfill \end{array}$$

(41b)

where $\kappa$ is the update rate.

As depicted in Figure 3, the SH-MAMFAC framework is comprised of an aerial–ground collaborative environment along with $I+1$ agents. During the centralized offline training phase, agent k is provided with the observations and actions of all other agents in addition to its own local observations. The actor network for CU i processes the local observation $o_{i,t}$ and generates a normalized continuous action $a_{\mathrm{con},i,t}$. The UAV provides its observation $o_{I+1,t}$ to the actor network ${\pi }_{I+1}\left(o_{I+1};{\mathsf{\Theta }}_{\pi ,I+1}\right)$, which produces the normalized continuous action $a_{\mathrm{con},I+1,t}$. Exploration noise $\mathcal{N}$ is required to encourage exploration and prevent policies from becoming trapped in local optima. The noise follows a normal distribution with a mean of zero and a standard deviation of ${\sigma }_{\mathcal{N}}$. As the agents carry out their joint actions ${\mathbf{a}}_t$, the environment progresses to the next state $s_{t+1}$, and the agents receive rewards ${\mathbf{r}}_t$ and updated observations ${\mathbf{o}}_{t+1}$.

If the UAV breaches trajectory or energy constraints, it faces a penalty in its reward, and the affected flight may be terminated. In the considered scenario, the total CHP is jointly determined by the actions of $I+1$ agents. Learning the Q-function with observations that include information about the actions of other agents aids in training. Since each agent can perceive the actions of other agents during the offline training phase, the environment is considered stable for each agent, addressing the challenges posed by a dynamic environment. In the execution phase, the actor network relies solely on local observations, allowing each agent to independently determine its actions without needing information about other agents. The training procedure for the proposed algorithm is outlined in Algorithm 1.

Algorithm 1 SH-MAMFAC Algorithm

1: Input: Structures of actor and critic networks for agent k ($k\in \mathcal{K}=\{1,\cdots ,I+1\}$), buffer capacity $\left|\mathcal{M}\right|$, batch size $\left|\mathcal{B}\right|$, and episode duration $\widehat{T}$.   2: Output: Optimally trained parameters of actor networks for all agents.   3: Initialize actor networks ${\left\{{\mathsf{\Theta }}_{\pi ,k}\right\}}_{k\in \mathcal{K}}$, target actor networks ${\left\{{\mathsf{\Theta }}_{{\pi }^{\prime },k}\right\}}_{k\in \mathcal{K}}$, critic networks ${\left\{{\mathsf{\Theta }}_{Q,k}\right\}}_{k\in \mathcal{K}}$, target critic networks ${\left\{{\mathsf{\Theta }}_{Q^{\prime },k}\right\}}_{k\in \mathcal{K}}$.   4: Set up experience replay buffer.   5: for each episode do   6:     Establish initial positions of UAV and GUs.   7:     for each time slot t do   8:         for each agent k do   9:            Obtain observation $o_{k,t}$. 10:            Choose action $a_{k,t}={\pi }_k(o_{k,t};{\mathsf{\Theta }}_{\pi ,k})+\mathcal{N}$. 11:         end for 12:         Each agent performs its actions. 13:         if UAV violates trajectory or energy constraints then 14:            Recalibrate UAV position. 15:         end if 16:         Agents receive rewards ${\mathbf{r}}_t$ and updated observations ${\mathbf{o}}_{t+1}$. 17:         Save $({\mathbf{o}}_t,{\mathbf{a}}_t,{\mathbf{r}}_t,{\mathbf{o}}_{t+1},{\overline{\mathbf{a}}}_{l,t})$ in $\mathcal{M}$. 18:     end for 19:     if experience replay memory is full then 20:         Draw mini-batches of size $\left|\mathcal{B}\right|$ from the buffer. 21:         Adjust critic network parameters by minimizing the loss in (39). 22:         Optimize actor network by maximizing the policy objective in (40). 23:         Update target network parameters in (41). 24:     end if 25: end for

3.5. Computational Complexity Analysis

The computational complexity of the SH-MAMFAC algorithm is influenced by the designs of the actor and critic networks for both the CUs and UAV. Specifically, assuming that the actor network for CU i is comprised of $L_{\mathrm{ac},i}$ fully connected layers, and the critic network contains $L_{\mathrm{cr},i}$ fully connected layers, the number of neurons in the l-th layer of the actor network for CU i is denoted as $n_{\mathrm{ac},i,l}$, while the number of neurons in the $l^{\prime }$-th layer of the critic network for CU i is denoted as $n_{\mathrm{cr},i,l^{\prime }}$. Thus, the computational complexity for CU i at each time step is $\mathcal{O}\left({\sum }_{l=0}^{L_{\mathrm{ac},i}-1}{n}_{\mathrm{ac},i,l}{n}_{\mathrm{ac},i,l+1}+{\sum }_{l^{\prime }=0}^{L_{\mathrm{cr},i}-1}{n}_{\mathrm{cr},i,l^{\prime }}{n}_{\mathrm{cr},i,l^{\prime }+1}\right)$. Next, we analyze the computational complexity of the UAV network structure. Similarly, we assume that the UAV actor network comprises $L_{\mathrm{ac},\mathrm{UAV}}$ fully connected layers, while its critic network consists of $L_{\mathrm{cr},\mathrm{UAV}}$ fully connected layers. In the UAV actor network, the l-th layer is composed of $n_{\mathrm{ac},\mathrm{UAV},l}$ neurons, whereas the $l^{\prime }$-th layer of the UAV critic network includes $n_{\mathrm{cr},\mathrm{UAV},l^{\prime }}$ neurons. Thus, the computational complexity for the UAV at each time step is $\mathcal{O}\left({\sum }_{l=0}^{L_{\mathrm{ac},\mathrm{UAV}}-1}{n}_{\mathrm{ac},\mathrm{UAV},l}{n}_{\mathrm{ac},\mathrm{UAV},l+1}+\right.$ $\left.{\sum }_{l^{\prime }=0}^{L_{\mathrm{cr},\mathrm{UAV}}-1}{n}_{\mathrm{cr},\mathrm{UAV},l^{\prime }}{n}_{\mathrm{cr},\mathrm{UAV},l^{\prime }+1}\right)$. This paper considers a set of agents comprising I CUs and one UAV. If the process of training the neural network weights has a computational cost of W, then the total complexity of the SH-MAMFAC algorithm is given by $\mathcal{O}\left(W\right(IX+Y\left)\right)$, where $X={\sum }_{l=0}^{L_{\mathrm{ac},i}-1}{n}_{\mathrm{ac},i,l}{n}_{\mathrm{ac},i,l+1}+{\sum }_{l^{\prime }=0}^{L_{\mathrm{cr},i}-1}{n}_{\mathrm{cr},i,l^{\prime }}{n}_{\mathrm{cr},i,l^{\prime }+1}$ and $Y={\sum }_{l=0}^{L_{\mathrm{ac},\mathrm{UAV}}-1}{n}_{\mathrm{ac},\mathrm{UAV},l}{n}_{\mathrm{ac},\mathrm{UAV},l+1}+{\sum }_{l^{\prime }=0}^{L_{\mathrm{cr},\mathrm{UAV}}-1}{n}_{\mathrm{cr},\mathrm{UAV},l^{\prime }}{n}_{\mathrm{cr},\mathrm{UAV},l^{\prime }+1}$.

During the execution phase, the well-trained actor independently makes optimal decisions without relying on the actions or states of other agents, effectively reducing communication overhead. Each CU agent uses only the trained actor network and does not utilize the critic network. Thus, the execution complexity of the CU agent is defined by $\mathcal{O}\left({\sum }_{l=0}^{L_{\mathrm{ac},i}-1}{n}_{\mathrm{ac},i,l}{n}_{\mathrm{ac},i,l+1}\right)$. Similarly, the execution complexity of the UAV agent is defined as $\mathcal{O}\left({\sum }_{l=0}^{L_{\mathrm{ac},\mathrm{UAV}}-1}{n}_{\mathrm{ac},\mathrm{UAV},l}{n}_{\mathrm{ac},\mathrm{UAV},l+1}\right)$.

4. Numerical Result

4.1. Simulation Settings

In a simulation of the proposed algorithm within a $200\times 200{\mathrm{m}}^2$ square area [7], the system consists of an MBS, a UAV, and several GUs. The mobility of GUs is modeled using a reference point group mobility model in [25], with the number of CUs set to $I=50$ and the number of RUs set to $J=300$, unless otherwise specified. The height of the static MBS is $h_{\mathrm{m}}=10\mathrm{m}$, located at coordinates $[100,100,10]$. The initial flight height of the UAV is $H=100\mathrm{m}$, with a departure direction randomly generated within the range $[0,\pi /2]$, and a maximum movement distance in one time slot of $d_{\max }=10\mathrm{m}$. In each time slot t, each RU requests a file f from a set of $F=30$ files with probability $p_{j,f}\left(t\right)$, where each file has a size of $S_f=2\mathrm{Mbits}$ [7]. The file popularity follows a Zipf distribution with a slope parameter set to $\xi =0.8$ [17]. Both CUs and the UAV have limited cache capacities, capable of storing $S_i$ and $S_{\mathrm{u}}$ Mbits, respectively. The cache capacity for CUs is set to $S_i=10\mathrm{Mbits}$, while $S_{\mathrm{u}}$ varies from $2\mathrm{Mbits}$ to $20\mathrm{Mbits}$. Other system simulation parameters are summarized in

Table 3. Summary of simulation parameters.

System Parameters		Training Parameters
Parameters	Values	Parameters	Values
${\delta }_t,T$ [7,8]	0.2 s, 200 s	${\sigma }_{\mathcal{N}}$	0.1
$A,U$ [25]	0.503, 120	$\gamma$	0.95
$s,\rho$ [25]	0.05, 1.225	$\mathcal{M},{\mathcal{M}}_i$	3000
$B_{\mathrm{m}}$ [17]	20 MHz	$\mathcal{B},{\mathcal{B}}_i$	128
$d_r,v_r$ [25]	0.3, 4.03	$\kappa$	0.01
$P_{\mathrm{u}}$	3.5 W	${\varrho }_{\pi }$	${10}^{-5}$
$N_0$ [17]	−100 dBm/Hz	${\varrho }_Q$	${10}^{-4}$
$P_0,P_1$ [25]	79.86 W, 88.63 W	$L_{\mathrm{ac}},L_{\mathrm{cr}}$	5, 5
$E_{\mathrm{u},max}$ [23]	50 KJ	$n_{\mathrm{ac},1\to (L_{\mathrm{ac}}-1)}$	256, 64, 32
V [23]	15 m/s	$n_{\mathrm{cr},1\to (L_{\mathrm{cr}}-1)}$	256, 64, 32

.

We conducted experimental simulations on a laptop equipped with 16.0 GB of RAM and an Intel Core i7-11370H @ 3.30 GHz processor. The software platform used is Python 3.9 and PyTorch. ReLU is set as the activation function for each hidden layer. The Adam optimizer is utilized to update the actor and critic networks. The parameters for the proposed SH-MAMFAC algorithm’s actor and critic networks can also be found in Table 3. Additionally, the training process for our proposed algorithm consists of $\widehat{T}=6000$ episodes, with each episode containing $N=1000$ time slots.

4.2. Performance Evaluation

In this subsection, we assess the performance of the SH-MAMFAC algorithm by comparing it with the following benchmarks:

(1) Unscalable heterogeneous MADDPG (UH-MADDPG) [8]: All CUs and a UAV are considered as different types of agents, each with a decentralized actor and a centralized critic. The critic of each agent can access the actions of all other agents, but as the number of agents increases, the joint action space grows exponentially, leading to performance limitations in large-scale environments due to high computational and communication overhead.

(2) Unscalable centralized DRL (UC-DRL) [11]: A single CU is selected as the central agent, responsible for receiving collective observation information from all CUs and making decisions for the entire group. As the number of agents increases, the communication overhead for information gathering grows significantly. Additionally, communication delays in large-scale environments may result in outdated state information, affecting the timeliness and accuracy of decision-making.

(3) Scalable independent DRL (SI-DRL) [38]: Each agent makes decisions and updates its policy based solely on its own local observations, without considering the actions or observations of other agents. This method can scale to large agent environments but lacks mechanisms for collaboration and information sharing between agents. Additionally, it may face non-stationarity issues, and the convergence of the learning policy cannot be guaranteed.

This section evaluates and compares the training efficiency and computational performance of four DRL methods. Figure 4 illustrates the cumulative reward trends over 6000 episodes, where the solid lines represent the moving average over 100 episodes, and the shaded areas indicate the corresponding variance. The markers on the lines show the convergence within 6000 episodes. Additionally, we recorded and compared the computation time per episode, the number of episodes required to converge, and the total computation time. The results are summarized in

Table 4. Computational performance of different algorithms.

Algorithm	Episodic Time (s)	Number of Episodes	Total Time (h)
SH-MAMFAC	1.38	1600	0.61
UH-MADDPG	1.79	3000	1.49
UC-DRL	2.62	5500 *	4.00 *
SI-DRL	0.96	1100	0.29

* Fail to complete training within 4 h.

.

As depicted in Figure 4, all DRL algorithms demonstrate similar low episode rewards during the initial training phase, attributed to the agents exploring the environment without optimized control strategies. As learning progresses, the agents accumulate valuable experience, which is utilized for network updates, leading to refined strategies and an increase in rewards across all DRL methods. Notably, the simplest method, SI-DRL, exhibits the highest variance and unstable learning behavior, resulting in the lowest cumulative rewards. This instability stems from agents training independently, relying solely on individual observations and actions, while disregarding the actions of others, which in turn leads to a non-stationary environment. In contrast, both UC-DRL and UH-MADDPG address the issue of non-stationarity by incorporating information from all agents, thereby enhancing training performance. However, these methods raise concerns regarding privacy, as they necessitate the sharing of local observations and actions from all agents. Furthermore, UC-DRL becomes infeasible in large-scale environments (e.g., when $I=50$), as it suffers from the curse of dimensionality, resulting in an inability to complete training within a 4-h period, with the process halting at approximately 5500 episodes. SH-MAMFAC, leveraging a mean-field approximation, efficiently captures the average interaction effects among agents, thus optimizing multi-agent collaboration and significantly improving scalability, ultimately achieving the highest reward levels.

The computational performance of the various methods is compared as presented in Table 4, revealing that SI-DRL exhibits the shortest per-episode computation time. This is primarily due to its requirement to process only local observations and action inputs. In contrast, UC-DRL demonstrates the longest computation time, largely attributed to its higher computational overhead and latency issues. As previously mentioned, UC-DRL failed to complete the training within 4 h, and the total computation time shown in the table corresponds only to 5500 episodes. In comparison, UH-MADDPG reduces training time effectively by employing a centralized training with decentralized execution (CTDE) framework. During the execution phase, agents can make decisions independently, resulting in more efficient utilization of computational resources and further reducing training time. In simpler environments without coordination, SI-DRL achieves convergence the fastest [37]. The proposed SH-MAMFAC method, through efficient cooperation among agents, better adapts to the dynamic nature of communication environments and achieves convergence within 1600 episodes, taking only 0.61 h. Its performance surpasses both UH-MADDPG and UC-DRL.

Figure 5 compares the CHP of the proposed SH-MAMFAC algorithm with all baseline methods under different RU densities. At the initial stage, with lower RU density, the CHP of all schemes increases due to sufficient system cache space and bandwidth. However, when the number of RUs exceeds 300, cache and transmission resources become strained, leading to a sharp drop in CHP for most strategies, indicating that traditional methods can effectively serve up to around 300 RUs. Notably, SH-MAMFAC consistently achieves the highest CHP across all RU counts. For instance, at $J=300$, SH-MAMFAC’s CHP reached 91.36%, exceeding UH-MADDPG by 4.42%, UC-DRL by 15.84%, and SI-DRL by 38.99%. When $J=500$, SH-MAMFAC’s CHP decreased to 81.32%, with a decline of 10.04%, significantly lower than the declines of UH-MADDPG (24.84%) and UC-DRL (26.86%), and slightly lower than SI-DRL (11.62%). This can be attributed to MAMFAC’s multi-agent cooperative learning, which enables efficient resource allocation and cache utilization. In contrast, baseline methods suffer from more severe performance degradation at high RU densities. UH-MADDPG adopts a CTDE framework, which results in an exponentially increasing joint action space as the number of agents grows, leading to increased computational and communication complexity, thereby hindering effective learning and decision making. UC-DRL employs centralized training, with a single CU acting as the central agent, but as the network scales up, it struggles with managing large volumes of global information, causing increased system delay and reduced performance. In SI-DRL, agents train and make decisions independently, lacking resource coordination, leading to lower CHP. Furthermore, LRU (the least recently used (LRU) algorithm is a widely adopted cache eviction method that targets the item in the cache which has been accessed the least recently for removal) lacks an intelligent learning mechanism and relies solely on historical requests, making it unable to adapt to dynamic content and complex network environments, resulting in consistently low CHP. In summary, SH-MAMFAC effectively addresses the scalability issues inherent in traditional CTDE methods in large-scale multi-agent environments, making it suitable for real-world large-scale networks that require servicing a substantial number of GUs and handling high data requests, thus improving system CHP and QoS.

Figure 6 presents a comparison of transmission delay for various algorithms under different bandwidth settings. As the CU bandwidth increases, the transmission delay for all strategies shows a decreasing trend. This indicates that a higher bandwidth improves the data transmission rate, thereby reducing transmission delay. The increase in bandwidth allows more data to be transmitted per unit time, enhancing the STP of D2D communication. This enables RUs to more frequently obtain the required files directly from CUs via D2D communication without frequently accessing the UAV or MBS, thereby reducing transmission delay. SH-MAMFAC demonstrates the lowest transmission delay across all bandwidth settings, due to its efficient cooperation among numerous CU agents, allowing RUs to acquire content more effectively from CUs and ensuring stable and superior performance under varying bandwidth conditions. UH-MADDPG, due to its reliance on joint action space, fails to fully utilize the potential of increased bandwidth, resulting in a less significant reduction in transmission delay compared to SH-MAMFAC. UC-DRL depends on a central agent for decision making, and the high communication overhead and delay required for information collection limit its performance improvement with increasing bandwidth. Thus, UC-DRL faces challenges in achieving efficient real-time decision making in large-scale networks. SI-DRL, although having good scalability with independent decision making for each agent, lacks a cooperation mechanism, leading to poor delay performance. Additionally, LRU, due to its inability to adapt to network dynamics and bandwidth changes, consistently has the highest transmission delay.

In Figure 7, the optimized trajectory of the UAV is presented. To clearly illustrate the network layout and communication mechanism, only 5 CUs and 20 RUs are shown. The trajectory of a UAV is obtained using the proposed SH-MAMFAC algorithm. CUs are represented by blue circles, RUs by green crosses, and the UAV’s flight path is indicated by a red line, illustrating its trajectory over a service duration of $T=200$ s. The D2D communication coverage area is depicted by cyan dashed circles. Some RUs are within the coverage of CUs and are primarily served by them. However, RUs outside the D2D coverage are served by either the UAV or the MBS, with those served by the UAV marked as green stars in the figure. This highlights the important role of UAVs in enhancing communication networks, particularly for RUs facing connectivity challenges. Additionally, due to limited battery capacity, the UAV must return after serving some RUs, resulting in only a subset of RUs being able to connect with it.

Figure 8 illustrates the transmission delay comparison of different algorithms under different UAV flight altitudes. The flight altitude of the UAV significantly affects the channel conditions and path loss, impacting communication performance. We observe that when the UAV is at a lower altitude, transmission delay increases with altitude. This is attributable to the reduced likelihood of LoS transmission and the greater average path loss associated with higher altitudes. As H continues to increase, the probability of LoS transmission improves and becomes more prevalent. During this phase, the average path loss decreases with increasing altitude, leading to a reduction in delay. When the probability of LoS transmission reaches a high level, for example, at $H=160$ m, the effect of altitude on the likelihood of LoS and NLoS transmissions becomes minimal. Consequently, the average path loss rises with H, leading to an increase in delay. It can be observed that SH-MAMFAC dynamically adjusts its strategy to effectively adapt to changes in the channel environment caused by altitude variations, maintaining the lowest transmission delay across various flight altitudes. This makes it particularly suitable for complex network environments requiring dynamic altitude adjustments.

Figure 9 compares the CHPs of different algorithms under varying UAV cache sizes. The results show that as the cache size increases, the CHPs of all algorithms improve significantly. This improvement is mainly due to larger cache capacity allowing more popular files to be stored, increasing the probability of requested files being cached and thus boosting the CHP. Specifically, the increase in CHP is especially notable as the cache size grows from smaller values. For instance, SH-MAMFAC shows a 42.48% increase in CHP when the cache size expands from 2 Mbits to 6 Mbits. This is because each additional unit of cache allows the system to store more frequently requested files, significantly enhancing the CHP. However, as the cache size further expands to larger values, the additional space is often used for less frequently accessed files, slowing down the CHP improvement. SH-MAMFAC consistently achieves the highest CHP across all cache sizes, demonstrating its ability to optimize caching strategies more effectively and maintain good scalability and robustness in large-scale network environments.

To demonstrate the advantages of air–ground cooperative networks and joint optimization, we implement SH-MAMFAC variants without trajectory optimization or power control. These strategies are named no-power SH-MAMFAC (NP-SH-MAMFAC), no-trajectory SH-MAMFAC (NT-SH-MAMFAC), and no-power and no-trajectory SH-MAMFAC (NPNT-SH-MAMFAC). When trajectory optimization is not performed, UAVs follow the circular trajectory defined in [39]. Figure 10 shows the transmission delay of different schemes at maximum transmission power $P_{max}$. Clearly, as $P_{max}$ increases, the delay of all four schemes decreases. SH-MAMFAC achieves the lowest transmission delay by jointly optimizing UAV trajectory and power control. NP-SH-MAMFAC optimizes only the trajectory without power control, while NT-SH-MAMFAC optimizes only power control without trajectory optimization. At low transmission power ($P_{max}<1.7$ W), NP-SH-MAMFAC significantly outperforms NT-SH-MAMFAC in terms of delay. However, when $P_{max}>1.7$ W, NT-SH-MAMFAC performs better due to its power control, which allows further delay reduction at higher power levels. NPNT-SH-MAMFAC, lacking both power and trajectory optimization, shows the highest delay across all power settings, confirming the importance of joint optimization for performance improvement. These results indicate that trajectory optimization is more critical at low transmission power, as optimizing the UAV path effectively covers user areas, reduces transmission distance and path loss, and lowers delay. At higher transmission power, power control becomes more crucial, as dynamically adjusting transmission power improves efficiency, reduces interference, and significantly decreases transmission delay. SH-MAMFAC, by optimizing both trajectory and power control, leverages their synergy to ensure optimal communication performance under varying power conditions.

5. Conclusions

This paper addressed the joint power allocation, trajectory design, and caching placement problem in UAV-assisted D2D caching networks. We proposed a mixed-integer programming problem to maximize the total CHP, incorporating constraints on user demands and limited resources. To tackle scalability issues in large-scale UAV-assisted caching systems, we introduced the SH-MAMFAC algorithm, which employed mean-field MADRL to model interaction dynamics in heterogeneous multi-agent systems, including both homogeneous CUs and a UAV. In this framework, each agent makes autonomous decisions based on its local observations after centralized training. The simulation results demonstrated that our algorithm significantly reduced training time for large-scale multi-agent systems and offered strong deployability. Compared to traditional DRL algorithms, our approach effectively improved the CHP and reduced the transmission delay. Furthermore, the SH-MAMFAC algorithm outperformed the NP-SH-MAMFAC, NT-SH-MAMFAC, and NPNT-SH-MAMFAC algorithms in terms of delay, highlighting the necessity of joint optimization for enhancing network performance.

While SH-MAMFAC demonstrates excellent performance in addressing environmental dynamics, enhancing service satisfaction for GUs and improving deployment capabilities in real heterogeneous networks, there are still areas for further improvement. Specifically, the UAV-D2D network must continuously adapt and dynamically reconfigure for optimal performance, which presents significant challenges for centralized training in environments requiring real-time learning and adaptation. To address these challenges, we need to explore mechanisms for distributed decision making, system scalability, and overall efficiency. By enabling each CU or UAV agent to make autonomous decisions based on local observations and interactions rather than relying on centralized control, we can significantly enhance the system’s flexibility and adaptability. This distributed approach reduces communication overhead and the computational burden of centralized control. It also ensures effective scaling of the network as the number of agents increases. We will continue to explore these potential improvement directions in future work.