Kenza Hamidouche
email: kenza.hamidouche@centralesupelec.fr

Walid Saad
email: walids@vt.edu

Mérouane Debbah
email: merouane.debbah@huawei.com

Breaking the Economic Barrier of Caching in Cellular Networks: Incentives and Contracts

In this paper, a novel approach for providing incentives for caching in small cell networks (SCNs) is proposed based on the economics framework of contract theory. In this model, a mobile network operator (MNO) designs contracts that will be offered to a number of content providers (CPs) to motivate them to cache their content at the MNO's small base stations (SBSs). A practical model in which information about the traffic generated by the CPs' users is not known to the MNO is considered. Under such asymmetric information, the incentive contract between the MNO and each CP is properly designed so as to determine the amount of allocated storage to the CP and the charged price by the MNO. The contracts are derived by the MNO in a way to maximize the global benefit of the CPs and prevent them from using their private information to manipulate the outcome of the caching process. For this interdependent contract model, the closed-form expressions of the price and the allocated storage space to each CP are derived. This proposed mechanism is shown to satisfy the sufficient and necessary conditions for the feasibility of a contract. Moreover, it is shown that the proposed pricing model is budget balanced, enabling the MNO to cover all the caching expenses via the prices charged to the CPs. Simulation results show that none of the CPs will have an incentive to choose a contract designed for CPs with different traffic loads.

I. INTRODUCTION

The capacity limitations of the backhaul links that connect small base stations (SBSs) to the core network in a small cell network (SCN), make it difficult for the mobile network operators (MNOs) and content providers (CPs) to ensure the required data rates for emerging, bandwidth intensive users' applications especially during peak hours [START_REF] Andrews | What will 5g be?[END_REF]. To overcome these backhaul limitations and meet users' requirements in terms of quality-of-service (QoS), distributed caching at the network edge has been recently proposed as a new promising solution [START_REF] Poularakis | Exploiting caching and multicast for 5G wireless networks[END_REF][START_REF] Blaszczyszyn | Optimal geographic caching in cellular networks[END_REF][START_REF] Blasco | Content-level selective offloading in heterogeneous networks: Multi-armed bandit optimization and regret bounds[END_REF][START_REF] Ji | Fundamental limits of caching in wireless d2d networks[END_REF].

Caching typically relies on storing the most popular files at the levels of SBSs and devices, to reduce backhaul traffic [START_REF] Poularakis | Exploiting caching and multicast for 5G wireless networks[END_REF]. To successfully deploy SCN caching solutions, MNOs require the cooperation of the CPs by sharing their content and providing this content's global [START_REF] Paschos | Wireless caching: Technical misconceptions and business barriers[END_REF]. However, although CPs can improve the QoS of their users by caching, they might be reluctant to share their content with the MNOs. This can be due to reasons such as privacy since, once the content cached at the SBSs, the MNOs can get access to all the CP's files This research was supported by ERC Starting Grant 305123 MORE (Advanced Mathematical Tools for Complex Network Engineering), the ANR project: WisePhy: Sécurité pour les communications sans fil à la couche physique, the U.S. National Science Foundation under Grants CNS-1513697 and CNS-1460316.

as well as the traffic dynamics of the users subscribed to the CP [START_REF] Paschos | Wireless caching: Technical misconceptions and business barriers[END_REF]. Thus, the MNOs must provide incentives to the CPs to share their data and cache it at the SBSs by introducing suitable economic arrangements that can be beneficial for both the MNO and the CPs [START_REF] Paschos | Wireless caching: Technical misconceptions and business barriers[END_REF][START_REF] Li | Pricing and resource allocation via game theory for a small-cell video caching system[END_REF].

Most of the existing literature on caching [START_REF] Poularakis | Exploiting caching and multicast for 5G wireless networks[END_REF][START_REF] Blaszczyszyn | Optimal geographic caching in cellular networks[END_REF][START_REF] Blasco | Content-level selective offloading in heterogeneous networks: Multi-armed bandit optimization and regret bounds[END_REF][START_REF] Ji | Fundamental limits of caching in wireless d2d networks[END_REF] has focused on determining the optimal caching policy under various network scenarios. In [START_REF] Poularakis | Exploiting caching and multicast for 5G wireless networks[END_REF], the authors proposed a caching policy that optimizes the overall energy consumption of the SBSs while accounting for the multicast opportunities. The authors in [START_REF] Blaszczyszyn | Optimal geographic caching in cellular networks[END_REF] proposed a caching policy that accounts for the geographical position of the SBSs and users. The work in [START_REF] Blasco | Content-level selective offloading in heterogeneous networks: Multi-armed bandit optimization and regret bounds[END_REF] formulated the cache placement problem as a multi armed-bandit problem in which the SBSs do not know the popularity of the files. The authors in [START_REF] Ji | Fundamental limits of caching in wireless d2d networks[END_REF] proposed a coded caching strategy for wireless device-to-device (D2D) networks.

Remarkably, none of these works have addressed the economic aspect of caching. Recently, the works in [START_REF] Li | Pricing and resource allocation via game theory for a small-cell video caching system[END_REF] and [START_REF] Chen | Caching incentive design in wireless d2d networks: A stackelberg game approach[END_REF] provides some insights on the economics of caching using game-theoretic solutions. However, these works typically assume that the players are honest and have perfect knowledge of the cellular network information, which is not a reasonable assumption in practice. In fact, the CPs can modify the outcome of the incentive mechanism in their favor by attempting to manipulate the MNO's agreements. For example, some CPs can declare inaccurate information about their traffic load so as to mislead the MNO into charging them much lower prices. Such possibilities can therefore incentivize some of the CPs to not share their truthful private information in which case the results in [START_REF] Li | Pricing and resource allocation via game theory for a small-cell video caching system[END_REF] and [START_REF] Chen | Caching incentive design in wireless d2d networks: A stackelberg game approach[END_REF] no longer hold.

The main contribution of this paper is to introduce a novel incentive mechanism for facilitating the deployment of caching in SCNs. We consider a model in which an MNO proposes agreements to the CPs to incite them to cache their content. We consider a practical scenario with asymmetric information within which the CPs can be of different, private types that are unknown to the MNO. This private information pertains to the level of generated CP traffic which corresponds to the popularity of its content. The proposed approach, based on contract theory [START_REF] Borgers | An introduction to the theory of mechanism design[END_REF], allows the MNO to define a contract for each CP in the presence of asymmetric information, by fixing a price for the allocated storage space to the CP. The goal of the MNO is to maximize the global reward of the CPs and cover the expenses of the caching process. Unlike classical contracttheoretic models [START_REF] Gao | Spectrum trading in cognitive radio networks: A contract-theoretic modeling approach[END_REF][START_REF] Duan | Cooperative spectrum sharing: a contract-based approach[END_REF][START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF], we show that the proposed model exhibits strategic interdependence between the CPs. Consequently, for this interdependent contract model, we derive the closed-form expression of the price and prove that the resulting mechanism satisfies the feasibility constraints of a contract. The designed mechanism is then shown to be budget balanced as it allows the MNO to offload its backhaul while covering the cost of caching the CPs' content. Simulation results show the performance advantage of the proposed scheme as well as its ability to incite truthfulness on the CPs.

II. SYSTEM MODEL Consider a SCN composed of a set M of macro base stations (MBSs) and a set S of SBSs deployed by an MNO. The SBSs are connected to the MBSs via capacity-limited backhaul links and serve a set of users U that are subscribed to multiple CPs from a set C. The SBSs can cache content to offload the MNO's backhaul. However, the MNOs needs the CPs' cooperation to cache their content. To motivate the CPs to participate in the caching, the MNO must offer contracts that present significantly improved QoS for the CPs' users. The contract terms between an MNO and a CP determine the price charged to the CP by the MNO and the amount of storage space offered.

We consider heterogeneous CPs with different traffic loads and content popularity. Based on this traffic load, the CPs have different incentive levels towards sharing their content with the MNOs. Naturally, there is an information asymmetry between the MNO and the CPs. The CP is aware of its users' traffic as well as its preferences while the MNO may not have that information. Consequently, the CPs may have an incentive to not reveal their correct types so as to pay lower prices to the MNO. To overcome this challenge, the MNO must specify a suitable performance-reward bundle contract (π, ρ), where π is the monetary reward that is paid by the CP, and ρ is the storage space allocated to the CP.

The goal of each CP k ∈ C is to maximize the performance of its users which depends on the amount of content that other CPs in the set C \k will be willing to cache at the SBSs. In fact, the higher the traffic load of the CPs in C\k, more storage space is allocated to cache their content. Thus, less storage space can be made available for the CP k. This will negatively impact the data rate of CP k's users as more requests need to be served via the capacity-limited backhaul. Moreover, by introducing caching, the traffic load of each SBS increases when caching highly popular files. Thus, more power is required at the SBS to serve all the requests for a highly popular file. Consequently, other CPs' users might experience large interference level from the SBSs that cache the files of CPs having a higher traffic load.

A. Transmission Data Rate

The performance r k of a CP k that results from caching its contents is measured by the transmission rate that its users experience from the serving MBS or SBS. When an SBS i serves a user j, the data rate will be:

αij(ρ(θ), θ) = Et wij log 1 + pij(ρ, θi)|hij| 2 σ 2 + I(ρ, θ) , (1)
where I(ρ, θ) = k∈S\i p kj (ρ, θ)|h kj | 2 is the interference experienced by user j from all the other SBSs. w ij is the channel bandwidth, p ij (ρ, θ) is the transmit power from SBS i to user j, |h ij | 2 is the channel gain between SBS i and user j, and σ 2 is the variance of the Gaussian noise. The vector θ = [θ 1 , ..., θ C ] represents the traffic load of the CPs. Thus, the higher the traffic load of the CPs, the higher is the interference experienced by the users served from the neighboring SBSs. Since caching is done during off-peak periods, the transmit power is averaged over the considered time period.

If the file is cached at the SBSs, then the users will experience a relatively high data rate as the content is closer to them. However, if the data is not available at the associated SBS to serving a certain user, then the SBS must fetch the user content from the MBS over the capacity-limited and congested backhaul, yielding higher delays. The data rate of a user j requesting file f from its associated SBS i can be given by:

r ij (θ) = (1 -β if (ρ(θ), θ)) min {α ij (ρ(θ), θ), α mi } + β if (ρ(θ), θ)α ij (ρ(θ), θ), (2)
where β ∈ {0, 1} S×F k is the outcome of the MNO's storage allocation ρ and F k is the cardinality of the set of files F k provided by a CP k. β depends on the caching policy ρ and the traffic load of the CPs. For instance, when a CP k has highly popular files or its willingness level to cache its content is high, it will impact the storage allocation to the other CPs C \ k. The larger the number of files that CP k wants to cache, the lower is the storage space that will be allocated to other CPs C \ k and vice versa. Each entry β if is a binary variable that equals 1 if file f is cached at SBS i and 0 otherwise. α mi is the data rate from the MBS m to SBS i and is given by:

α mi = E t w mi log 1 + p mi |h mi | 2 σ 2 + I , (3)
where I = l∈M\m p li |h li | 2 is the interference experienced by SBS i from all the other transmitting MBSs. Thus, the total rate of the users of CP k can be given by:

r k (ρ(θ), θ k ) = i∈S j∈U ki r ij (θ), (4)
where U ki ⊆ U k is the set of users that request at least one file from CP k by using SBS i, and U k is the set of users requesting files of CP k. Here, we note that, in our model, each CP will have private information that is modeled as a type of CP as discussed next.

B. Content Provider Type

We define the CP's type to be a representation of its traffic load and content popularity. In fact, when the MNO offers a contract, it must account for the generated traffic by the CPs. For example, by caching the contents of CPs with a high traffic, the MNO can serve more requests locally, thus decreasing its backhaul load considerably. Here, we consider that the number of CP types belongs to a discrete, finite space and grouped as follows:

Definition 1. There are C CPs that generate traffic over an MNO's network. The CPs' types are sorted in an ascending order and classified into K types θ 1 , ..., θ K with K ≤ C. Each type includes properties such as the willingness to cache and the global popularity of the CPs files. The types are ordered as follows: θ 1 < ... < θ k < ... < θ K , k ∈ {1, ..., K}.

Since the types are not known by the MNO, the CPs can announce wrong information about their types so that they improve the performance of their users. For example, by claiming that its content popularity is higher than it actually is, a CP k can mislead the MNO to allocate more storage space. In such a case, the CP can end up paying lower prices while also lowering the interference experienced by its users. Indeed, the truthful popularity information is necessary for the MNO to define the contracts that optimize the benefit of the CPs and cover the implementation costs of caching. Such cost includes the expenses of deploying storage devices and the required power to download the content and refresh the storage units. Here, our goal is to design contracts that incentivize the CPs to reveal the true values of their types θ to the MNO. To this end, the contracts will be designed such that no CP can profit by choosing a contract that is designed for other types.

C. Content Provider Model

The utility function of a CP k of type θ k that decides to cache a set of files F k at the operator's network is:

u k (θ) = r k (ρ k (θ), θ k ) -π k (θ), (5)
where r k (ρ) is defined in (4) and represents the valuation function regarding the rewards, which is a strictly increasing concave function of ρ k , with r(0) = 0 and r (ρ k ) > 0, r (ρ) < 0 for all ρ k . π k represents the price charged by the MNO for a storage allocation ρ k .

D. Mobile Network Operators Model

By caching the content of the CPs, the MNO will be able to reduce the traffic load on its backhaul. This benefit depends on the traffic load of the CPs as dictated by the popularity of their cached traffic. Thus, the MNO will generally prefer to cache the most popular files. By doing so, for the same storage capacity, the load can be reduced more for a CP whose files have a high popularity compared to other CPs. Thus, the cost of storage c s at the MNO can be given as a function of the traffic load of the considered CP as, c(θ) = log (1 + θ). This storage cost function increases quickly up to a certain threshold and then increases slowly. It is suitable to model the storage cost as the MNO must allocate more storage space to serve a given traffic load, and this cost becomes insignificant when the traffic load increases as some requests become redundant.

A proper utility function for the MNO can be defined as the monetary reward that is charged to the CPs minus the cost of the allocated resources by the MNO, including storage.

v k (θ) = π k (θ) -c k (θ, θ k ), (6)
where π k is the price that the operator charges CPs of type k.

The total expected utility of the operator can be given by

v = k∈C v k (θ). (7)
In the considered model, the MNO is assumed to get the CP's types directly from the CPs. Based on this information, the goal of the MNO is to determine a contract for all possible CP types that maximizes the global benefit of the CPs. At the same time, the CP ensures that its utility is nonnegative by making the prices charged to the CPs to at least cover its cost. This optimization problem can be defined as follows:

max (π k ,ρ k ) k∈C u k (ρ k (θ), θ k ) subject to v ≥ 0. (8)
In this formulation, we do not make any constraint on the participation of the CPs. Thus, when proposing the contracts resulting from solving (8), CPs may prefer not to select any of the contracts or select contracts that are not designed for their types. To analyze this economic incentive problem, next, we propose a solution based on the framework of contract theory for designing feasible contracts [START_REF] Borgers | An introduction to the theory of mechanism design[END_REF].

III. PROPOSED INCENTIVE MECHANISM FOR CACHING

The QoS achieved by the CPs' users depends on the interference from the other SBSs as shown in ( 1) and (3). Moreover, the allocated storage capacity to a given CP depends on the number of CPs that have signed contracts with the same MNO. For instance, the more storage is allocated to a CP k the lower is the available storage capacity for other CPs. Thus, in the considered model, there is interdependence between the signed contracts by the different agents.

Classical contract theory models that are used to model resource allocation problems in wireless networks such as in [START_REF] Gao | Spectrum trading in cognitive radio networks: A contract-theoretic modeling approach[END_REF][START_REF] Duan | Cooperative spectrum sharing: a contract-based approach[END_REF][START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF] cannot be applied for the analysis of caching incentive problem between CPs and an MNO defined in [START_REF] Chen | Caching incentive design in wireless d2d networks: A stackelberg game approach[END_REF]. In fact, these works assume that the contract selected by a CP does not impact the utility of other CPs or focus only on models with one MNO and one CP. Thus, none of the existing works account for the interactions between the CPs. Moreover, the revelation of misleading information by a given CP in a multiple CPs model not only impacts the MNO but also impacts other CPs, which is not considered in [START_REF] Gao | Spectrum trading in cognitive radio networks: A contract-theoretic modeling approach[END_REF][START_REF] Duan | Cooperative spectrum sharing: a contract-based approach[END_REF][START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF].

To define the most appropriate contract for the formulated problem (8), we consider the so-called truthful dominant strategy implementation. Under such contracts, the solution that maximizes the utility function of the CPs will require those CPs to reveal their private information which, in our model, pertains to the real popularity of their content. The goal of the MNO is to maximize the social welfare which effectively captures the global QoS that is experienced by the users of all CPs. Moreover, the MNO ensures that the cost of serving these users is at least covered by the price charged to the CPs. To incite the CPs to collaborate with the MNO via caching, the contract that a CP selects must be feasible in that it satisfies the following necessary and sufficient constraints: Definition 2. Ex-post Individual Rationality (IR): The contract that a CP selects should guarantee that the utility of the CP is nonnegative for any θ -k declared by the other CPs,

r k (ρ k (θ k , θ -k ), θ k ) -π k ≥ 0,
(ρ k (θ k , θ -k ), θ k )-π k ≥ r k (ρ k ( θk , θ -k ), θ k )-π k . (10)

A. Incentive Mechanism Analysis

The goal of the MNO is to determine a pricing policy that motivates the CPs to declare their real type and simultaneously participate in the caching system through a budget balanced mechanism, i.e., the MNO would not experience a negative utility and its effort is covered by the price charged to the CPs. To this end, the CPs need to declare their types for the MNO that in turn optimizes their utility while accounting for the necessary conditions for contracts feasibility. The optimization problem of the MNO can be defined as follows:

max (π k ,ρ k ) k∈C u k (ρ k (θ), θ k ) subject to (9), ( 10
), v ≥ 0. (11)
The solution of this problem consists in the determination of the components of a contract that consist in the allocated storage space and the price charged to each CP. The closedform of the contract is provided by the following theorem.

Theorem 1. The unique efficient solution of the optimization problem [START_REF] Duan | Cooperative spectrum sharing: a contract-based approach[END_REF] can be given by:

ρ * k ∈ arg max ρ k i ri(ρi( θ), θi) -ci( θ) , ∀k, (12)
π k ( θ) = max

ρ i i =k ri(ρi( θ-k ), θi) -ci( θ-k ) (a) - i =k ri(ρ * i ( θ), θi) -ci( θ) (b) , (13)
where (a) represents the maximized social welfare when CP k is not considered while in (b), CP k is considered. Moreover, θ represents the revealed type by the CPs while θ is the real type of the CPs.

Proof. The proof is provided in the Appendix.

This result shows that, in order to determine the terms of a contract with a CP k, the MNO first, allocates the storage space to CP k by solving the optimization problem [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF]. It is clear that the problem in ( 12) is NP-hard and thus it is challenging to find the optimal storage allocation. To solve [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF], we use the framework of matching theory to analyze the assignment of storage space between the MNO that acts on behalf of its SBSs and the CPs [START_REF] Gale | College admissions and the stability of marriage[END_REF]. Matching theory is a suitable framework to solve NP-hard assignment problems such as in [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF]. As stated before, the allocated storage space to a CP k depends on the allocated storage to the other CPs which is known as externalities. We solve the problem using a swap-based deferred acceptance algorithm which is guaranteed to converge to a stable outcome [START_REF] Gale | College admissions and the stability of marriage[END_REF]. Each CP i starts by requesting from the MNO, a given storage space ρ i that maximizes r i (ρ i ( θ-k ), θi )-c i ( θ-k ). After receiving all the requests from the CPs and based on the caching policy used by the MNO, it defines the accepted requests that maximize [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF] and rejects the others. This is defined while accounting for the limited storage capacity of its SBSs. If the request of a given CP is rejected, the CP decreases the amount of the storage space it requests from the MNO. The MNO accepts new requests and rejects others based on the allocation configuration that maximizes [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF]. The procedure is repeated until there does not exist a CP k that prefers to be assigned a given storage capacity ρ k and this allocation also maximizes ( 12) at the MNO.

Once the storage space is allocated to the CPs, the price paid by a CP k is found from [START_REF] Gale | College admissions and the stability of marriage[END_REF], which accounts for the impact of CP k on the utility of other CPs. This price represents the difference between the global utility achieved by all the CPs when CP k participates in the caching process, and the global utility achieved by the CPs when CP k does not participate. Note that a CP k can impact other CPs utilities in two ways. The first one is through the allocated storage space. In fact, when more storage space is allocated to CP k, less storage is available for other CPs and thus more requests of these CPs are served via the backhaul. The second is the traffic load of CP k as the transmit power of the SBSs increases by increasing the number of served requests for that CP's files. Thus, we can deduce that higher traffic load of a CP k and large amounts of allocated storage to CP k will result in an increase in the price charged by the MNO. The dependence of the price on the traffic load of the CPs, i.e., CPs type, is given next.

Corollary 1. When θ k ≥ θ l then we have π k ≥ π l .
Proof. This results follows directly from the monotonicity property of the rate function and the structure of (13).

IV. SIMULATION RESULTS

For our simulations, we consider five CPs with different traffic load levels from 1 to 5 with type 5 being the highest load. A type-1 CP is chosen with no traffic and is used as a baseline to compare the performance of our mechanism with the case in which there is no caching. We consider a set of 100 files whose popularity follows a Zipf distribution of parameter α = 0.2. The MNO has one MBS that serves all the requests that cannot be served from the SBSs' cache. The number of SBSs is 10 and the total storage capacity of the SBSs is 1 Gbits. The transmit power of the SBSs is 1 W and the bandwidth capacity to 100 MHz.

In Fig. 1a, we show that the amount of content that is served via the backhaul for every CPs' users when a CP selects the contracts designed by the MNO for each type. We account for the fact that the price that can be paid by each CP is limited and the limit increases by increasing the type of the CPs. From Fig. 1a, we can observe that when high type CPs select contracts that are designed for CPs with lower traffic load, the amount of content that is served via the backhaul increases until it reaches the maximum which corresponds to the lowest type contract. Similarly, when low type CPs select the contracts designed for higher types, the amount of traffic served from the backhaul increases. This is due to the high charged price by the MNO for high type CPs and thus CPs of lower type cannot afford that which results in a larger amount of served content from the backhaul. Thus, this result validates the fact the proposed approach for cache incentive compatible as is forces each CP to choose the contract designed for its own type. In Fig. 1b, we compare the case in which all the CPs choose the contracts defined for their corresponding types by the MNO and the case in which the storage space is allocated equally by the MNO for all the CPs. Fig. 1b show that the model of equal storage allocation outperforms the proposed mechanism for low type CPs. In this case, the allocated storage space for the CPs is higher than the required storage by the CPs which appears through the utilities of type-2 CP and type-3 CP that are only 2% to 10% higher than their utilities when following the proposed mechanism. On the other hand, the proposed mechanism outperforms the equal storage space allocation model for high type CPs. In fact, the utility of high level types is higher when selecting the offered contracts by the MNO for their specific types. The utilities of type-4 and type-5 CPs are 50% to 140% higher than their achievable utilities in the model of equally allocated storage space. The proposed mechanism is more beneficial for the CPs as it allows the MNO to offer to all CPs only the amount of storage space they need. In contrast, when using the equal allocation approach, the MNO allocates to the low type CPs more than their storage requirements and insufficient space for high level type CPs.

In Fig. 1c, we show the variation of the mean utility of the CPs when increasing the total number of CPs. We consider two different values for the parameter α, 0.2 which corresponds to the case in which the files have comparable popularity, and 2 for the case in which some files are very popular while others have a very low popularity. We can observe that the mean utility for the CPs decreases by increasing the total number of CPs in the model. This is due to the data rate function that depends on the additional interference from the added CPs and the decrease of the available storage space that can be allocated to each CP. The achievable utility by a CP is up to 20% larger compared to the cases in which the CPs request a storage capacity that is larger or lower than the offered one by the MNO. Moreover, we can see that the popularity of the files impacts the mean utility of the CPs. In fact, the CPs can achieve a larger utility when the distribution of the popularity of the files is steep and thus by caching a file, a large amount of the requests can be served from the cache of the SBSs.

V. CONCLUSION In this paper, we have proposed a new incentive framework to motivate the CPs to cooperate with an MNO and cache their content at the MNO's SBSs. Based on contract theory, we have designed an incentive mechanism that allows the MNOs to offer a contract for each CP in which it sets the allocated storage for the CP and the charged price by the MNO for the caching service. This model accounts for both asymmetry of information and the interdependence between the different contracts. We have then derived the optimal pricing mechanisms and contracts that motivate the CPs to cache their content and reveal their private information. Simulations have shown the effectiveness of the proposed approach in inciting the participation of CPs in caching.

APPENDIX

We prove Theorem 1 by showing the following Lemmas.

Lemma 1. The proposed mechanism is incentive compatible.

Proof. We show this result by contradiction. Suppose that ρ is an efficient decision rule but (ρ, π) is not dominant strategy incentive compatible. Then, there exists i, θ, and θ such that:

r i (ρ i ( θi , θ -i ), θ i ) -π i (θ -i , θi ) > r i (ρ i (θ), θ i ) -π i (θ).
From [START_REF] Gale | College admissions and the stability of marriage[END_REF], this implies that

r i (ρ i ( θi , θ -i ), θ i ) -(a) > r i (ρ i (θ), θ i ) -(b), which is equivalent to ri(ρi( θi, θ-i), θi) - j =i rj(ρj( θi, θ-i), θj) -ci( θi, θ-i) > ri(ρi(θ), θi) - j =i rj(ρj(θi, θ-i), θj) -ci(θi, θ-i).
This contradicts the efficiency of ρ based on [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF] and thus, the assumption was incorrect.

Lemma 2. Truth telling is a dominant strategy under [START_REF] Gale | College admissions and the stability of marriage[END_REF].

Proof. Consider the problem of choosing the best type θi by a CP i. A best strategy for CP i solves

max θi r i (ρ i ( θ), θi ) -π i ( θ).
Substituting the payment function by the proposed mechanism [START_REF] Zhang | Contract-based incentive mechanisms for device-to-device communications in cellular networks[END_REF], we get Thus, CP i would pick a declaration θi that will lead the mechanism to pick a ρ which solves The proposed mechanism (13) will choose ρ in a way that solves the maximization problem (14) with θi = θ i . Thus, truth-telling is a dominant strategy for CP i.

Lemma 3. The proposed mechanism ( 13) is ex-post individually rational.

Proof. At the equilibrium, all the CPs are truthful and declare their real types. Thus, by replacing [START_REF] Gale | College admissions and the stability of marriage[END_REF] in the utility of a CP i, we have

ui = i ri(ρ * i (θ), θi) -ci(θ) - j =i rj(ρ * j (θ-i), θj) -cj(θ), (15)
where ρ * is the outcome that maximizes the social welfare. The CPs could have picked ρ i (θ -i ) instead of ρ i (θ) as a solution of the optimization problem, as it is one of the possible strategies. This is possible because the set of strategies, i.e., the total storage capacity of the MNO, is fixed and does not change by changing the set of participating CPs. Thus, j rj(ρ * j (θ), θj) ≥ j rj(ρ * j (θ-i), θj).

Furthermore, we know that the rate of a participating CP cannot be negative, i.e., rj(ρ * j (θ-i), θj) ≥ 0.

Therefore, i ri(ρ * i (θ), θi) ≥ j =i rj(ρ * j (θ-i), θj).

Thus, (15) is non-negative and the proposed mechanism is expost individual rational.

Lemma 4. The proposed mechanism ( 13) is weakly budgetbalanced.

Proof. Since the CPs are truth-telling at the equilibrium then we have i πi(θ) = Moreover, since the utility of a CP is a decreasing function of the number of the set of participating CPs, we have that, ∀i, j =i rj(ρj(θ-i), θj) -cj(θ-i, θj) ≥ j =i rj(ρj(θ), θj) -cj(θ, θj).

Thus, the proposed mechanism is weakly budget-balanced.

Next, we prove the provided result in Theorem 1.

Proof. Based on Lemma 1, Lemma 3 and Lemma 4, we can deduce that all the condition of the optimization problem [START_REF] Duan | Cooperative spectrum sharing: a contract-based approach[END_REF] are satisfied. Based on Lemma 2, we can deduce that the proposed mechanism [START_REF] Gale | College admissions and the stability of marriage[END_REF] is the unique efficient solution of the formulated problem [START_REF] Duan | Cooperative spectrum sharing: a contract-based approach[END_REF].

2 ρ 2 ρ 2 ρFigure 1 :

 2221 Figure 1: Numerical results for a) the amount of content served via the backhaul with respect to each CP's type, b) the utility of the CPs as a function of the CPs' type and the used storage allocation model, c) the mean utility of a CP with respect to the total number of CPs and the popularity distribution of the files.

  ( θ), θi)j =i rj(ρ( θ-i), θj) -cj(θj, θ-i) (a) + j =i rj(ρj( θ), θj) .Since (a) does not depend on θi , it is sufficient to solve max θi ri(ρi( θ), θi) + j =i rj(ρj( θ), θj) .

  ( θ), θi) + j =i rj(ρj( θ-i), θj) . (14)Under the proposed mechanism, ρ * ∈ arg max ρ ri(ρi( θ), θi) + i =j rj(ρj( θ), θj) .

  (θ-i), θj) -cj(θ-i, θj) j =i rj(ρj(θ), θj) -cj(θ, θj) .

  ∀k ∈ {1, ..., K}.[START_REF] Borgers | An introduction to the theory of mechanism design[END_REF]

	Definition 3. Incentive Compatibility (IC): A contract satisfies
	incentive compatibility constraint if each CP of type θ k prefers to reveal its real type θ k rather than another type θk , i.e.,
	r k


Mohammad Mozaffari

Walid Saad
email: walids@vt.edu

Mehdi Bennis
email: bennis@ee.oulu.fi

Mérouane Debbah
email: merouane.debbah@huawei.com

Mobile Internet of Things: Can UAVs Provide an Energy-Efficient Mobile Architecture?

In this paper, the optimal trajectory and deployment of multiple unmanned aerial vehicles (UAVs), used as aerial base stations to collect data from ground Internet of Things (IoT) devices, is investigated. In particular, to enable reliable uplink communications for IoT devices with a minimum energy consumption, a new approach for optimal mobility of the UAVs is proposed. First, given a fixed ground IoT network, the total transmit power of the devices is minimized by properly clustering the IoT devices with each cluster being served by one UAV. Next, to maintain energy-efficient communications in time-varying mobile IoT networks, the optimal trajectories of the UAVs are determined by exploiting the framework of optimal transport theory. Simulation results show that by using the proposed approach, the total transmit power of IoT devices for reliable uplink communications can be reduced by 56% compared to the fixed Voronoi deployment method. Moreover, our results yield the optimal paths that will be used by UAVs to serve the mobile IoT devices with a minimum energy consumption.

I. INTRODUCTION

The use of unmanned aerial vehicles (UAVs) as flying wireless communication platforms has received significant attention recently [START_REF] Mozaffari | Unmanned aerial vehicle with underlaid device-to-device communications: Performance and tradeoffs[END_REF]- [START_REF] Yaliniz | Efficient 3-D placement of an aerial base station in next generation cellular networks[END_REF]. On the one hand, UAVs can be used as wireless relays for improving connectivity of ground wireless devices and extending network coverage. On the other hand, UAVs can act as mobile aerial base stations to provide reliable downlink and uplink communications for ground users, and boost the capacity of wireless networks [START_REF] Yaliniz | Efficient 3-D placement of an aerial base station in next generation cellular networks[END_REF] and [START_REF] Mozaffari | Efficient deployment of multiple unmanned aerial vehicles for optimal wireless coverage[END_REF]. Compared to the terrestrial base stations, the advantage of using UAV-based aerial base stations is their ability to quickly and easily move. Furthermore, the high altitude of UAVs can enable line-of-sight (LoS) communication links to the ground users. Due to the adjustable altitude and mobility, UAVs can move towards potential ground users and establish reliable connections with a low transmit power. Hence, they can provide a cost-effective and energy-efficient solution to collect data from ground mobile users that are spread around a geographical area with limited terrestrial infrastructure.

Indeed, UAVs can play a key role in the Internet of Things (IoT) which is composed of small, battery-limited devices such as sensors, and health monitors [START_REF] Dawy | Towards massive machine type cellular communications[END_REF]- [START_REF] Hassanalieragh | Health monitoring and management using internet-of-things (IoT) sensing with cloud-based processing: Opportunities and challenges[END_REF]. These devices are typically unable to transmit over a long distance due to their energy constraints [START_REF] Lien | Toward ubiquitous massive accesses in 3gpp machine-to-machine communications[END_REF]. In such IoT scenarios, UAVs can dynamically move towards IoT devices, collect the IoT data, and transmit it to other devices which are out of the communication ranges of the transmitting IoT devices [START_REF] Lien | Toward ubiquitous massive accesses in 3gpp machine-to-machine communications[END_REF]. In this case, the UAVs play the role of moving aggregators for IoT networks. However, to effectively use UAVs for IoT communications, several challenges must be addressed such as optimal deployment and energy-efficient use of UAVs [START_REF] Mozaffari | Unmanned aerial vehicle with underlaid device-to-device communications: Performance and tradeoffs[END_REF].

In [START_REF] Mozaffari | Unmanned aerial vehicle with underlaid device-to-device communications: Performance and tradeoffs[END_REF], we investigated the optimal deployment and movement of a single UAV for supporting downlink wireless communications. However, this work was restricted to a single UAV and focused on the downlink. The work in [START_REF] Han | Optimization of MANET connectivity via smart deployment/movement of unmanned air vehicles[END_REF] analyzed the optimal trajectory of UAVs to enhance connectivity of ad-hoc networks. Nevertheless, this work did not study the optimal movement of multiple UAVs acting as aerial base stations. The work in [START_REF] Mozaffari | Optimal transport theory for power-efficient deployment of unmanned aerial vehicles[END_REF] studied the optimal deployment of UAVs and UAV-users association for static ground users with the goal of meeting users' rate requirements. In [START_REF] Pang | Efficient data collection for wireless rechargeable sensor clusters in harsh terrains using UAVs[END_REF], the authors used UAVs to efficiently collect data and recharge the clusters' head in a wireless sensor network which is partitioned into multiple clusters. However, this work is limited to a static sensor network, and does not investigate the optimal deployment of the UAVs. While the energy efficiency of uplink data transmission in a machine-to-machine (M2M) communication network was investigated in [START_REF] Abuzainab | Cognitive hierarchy theory for heterogeneous uplink multiple access in the internet of things[END_REF] and [START_REF] Tu | Energy-efficient algorithms and evaluations for massive access management in cellular based machine to machine communications[END_REF], the presence of UAVs was not considered. In fact, none of the prior studies [START_REF] Mozaffari | Unmanned aerial vehicle with underlaid device-to-device communications: Performance and tradeoffs[END_REF]- [START_REF] Yaliniz | Efficient 3-D placement of an aerial base station in next generation cellular networks[END_REF], and [START_REF] Han | Optimization of MANET connectivity via smart deployment/movement of unmanned air vehicles[END_REF]- [START_REF] Tu | Energy-efficient algorithms and evaluations for massive access management in cellular based machine to machine communications[END_REF], addressed the problem of optimal deployment and mobility of UAVs for enabling reliable and energy-efficient communications for mobile IoT devices.

The main contribution of this paper is to propose a novel approach for deploying multiple, mobile UAVs for energyefficient uplink data collection from mobile, ground IoT devices. First, to minimize the total transmit power of the IoT devices, we create multiple clusters where each one is served by one of the UAVs. Next, to guarantee energyefficient communications for the IoT devices in mobile and time-varying networks, we determine the optimal paths of the UAVs by exploiting dynamic clustering and optimal transport theory [START_REF] Villani | Topics in optimal transportation[END_REF]. Using the proposed method, the total transmit power of the IoT devices required for successful transmissions is minimized by the dynamic movement of the UAVs. In addition, the proposed approach will minimize the total energy needed for the UAVs to effectively move. The results show that, using our proposed framework, the total transmit power of the devices during the uplink transmissions can be reduced by 56% compared to the fixed Voronoi deployment method. Furthermore, given the optimal trajectories for UAVs, they can serve the mobile IoT devices with a minimum energy consumption.

The rest of this paper is organized as follows. In Section II, we present the system model and problem formulation. Section III presents the optimal devices' clustering approach. In Section IV, we address the mobility of the UAVs using discrete transport theory. In Section V, we provide the simulation results, and Section VI draws some conclusions.

II. SYSTEM MODEL AND PROBLEM FORMULATION

Consider an IoT system consisting of a set L = {1, 2, ..., L} of L IoT devices deployed within a geographical area. In this system, a set K = {1, 2, ..., K} of K UAVs must be deployed to collect the data from the ground devices in the uplink. The locations of device i ∈ L and UAV j ∈ K are, respectively, given by (x i , y i ) and (x u,j , y u,j , h j ) as shown in Figure 1. We assume that devices transmit in the uplink using orthogonal frequency division multiple access (OFDMA) and UAV j can support at most M j devices simultaneously. Note that, we consider a centralized network in which the locations of devices and UAVs are known to a control center such as a central cloud server. The ground IoT devices can be mobile (e.g. smart cars) and their data availability can be intermittent (e.g. sensors). Therefore, to effectively respond to the network mobility, it is essential that the UAVs optimally move for establishing reliable and energy-efficient communications with the devices. Note that, in our analysis, without loss of generality any mobility model can be accommodated.

For ground-to-air communications, each device will typically have a LoS view towards a specific UAV with a given probability. This LoS probability depends on the environment, location of the device and the UAV, and the elevation angle between the device and the UAV [START_REF] Hourani | Modeling air-to-ground path loss for low altitude platforms in urban environments[END_REF]. One suitable expression for the LoS probability is given by [START_REF] Hourani | Modeling air-to-ground path loss for low altitude platforms in urban environments[END_REF]:

P LoS = 1 1 + ψ exp(-β [θ -ψ]) , (1)
where ψ and β are constant values which depend on the carrier frequency and type of environment such as rural, urban, or dense urban, and θ is the elevation angle. Clearly, θ = 180 π × sin -1 hj dij , where

d ij = (x i -x u,j ) 2 + (y i -y u,j ) 2 + h 2 j
is the distance between device i and UAV j.

From (1), we can observe that by increasing the elevation angle or increasing the UAV altitude, the LoS probability increases. We assume that, the necessary condition for connecting a device to a UAV is to have a LoS probability greater than a threshold ( close to 1). In other words, P LoS (θ) ≥ ε, and hence, θ ≥ P -1

LoS (ε) leading to:

d ij ≤ h j sin P -1 LoS (ε) . ( 2
)
Note that (2) shows the necessary condition for connecting the device to the UAV. Therefore, a device will not be assigned to UAVs which are located at distances greater than

hj sin(P -1 LoS (ε)) . Now, considering the LoS link, the received signal power at UAV j from device i is given by [START_REF] Hourani | Modeling air-to-ground path loss for low altitude platforms in urban environments[END_REF] (in dB):

P ij r = P t,i -10α log 4πf c c d ij -η, (3)
where P t,i is the transmit power of device i in dB, f c is the carrier frequency, α = 2 is the path loss exponent for LoS propagation, η is an excessive path loss added to the free space propagation loss, and c is the speed of light.

In our model, the transmit power of the devices must satisfy the minimum signal-to-noise-ratio (SNR) required for a successful decoding at UAVs. For quadrature phase shift keying (QPSK) modulation, the minimum transmit power of device i needed to reach a bit error rate requirement of δ is:

P ij t = Q -1 (δ) 2 R b N o 2 10 η/10 4πf c d ij c 2 , (4)
where Q -1 (.) is the inverse Q-function, N o is the noise power spectral density, and R b is the transmission bit rate. Note that, to derive (4) using ( 3), we use the bit error rate expression for QPSK modulation as

P e = Q( 2P ij r R b No ) [14]
. Our first goal is to optimally move and deploy the UAVs in a way that the total transmit power of devices to reach the minimum SNR requirement for a successful decoding at the UAVs is minimized. In fact, the objective function is:

min Cj ,µj K j=1 i∈Cj P ij t , j ∈ K, (5)
where P ij t is the transmit power of device i to UAV j, and K is the number of UAVs. Also, C j is the set of devices assigned to UAV j, and µ j is the 3D location of UAV j.

From (4), we can observe that the transmit power is directly proportional to the distance squared. Hence, minimizing the power is equivalent to minimizing the distance squared. Then, using (2), (4), and (5), and considering the constraint on the maximum number of devices that can connect to each UAV, our optimization problem can be formulated as:

C * j , µ * j = arg min Cj ,µj K j=1 i∈Cj d 2 ij , j ∈ K, (6)
s.t. C j ∩ C m = ∅, j = m, j, m ∈ K, (7)
K j=1 |C j | = L, (8)
d ij ≤ h j sin P -1 LoS (ε) , (9)
|C j | ≤ M j , (10)
where |C j | is the number of devices assigned to UAV j, L is the total number of devices, and M j is the maximum number of devices that UAV j can support. ( 7) and ( 8) guarantee that each device connects to only one UAV. Clearly, we can consider the set of devices which are assigned to a UAV as a cluster, and place the corresponding UAV at the center of the cluster. Placing a UAV at the cluster center ensures that the UAV has a minimum total distance squared to all the cluster members. Hence, to solve problem (6), we need to find L clusters, and their centers which effectively correspond to the locations of the UAVs. Note that, in a time varying network in which the location of IoT devices change, the clusters will also change. Consequently, the location of UAVs as the center of the clusters must be updated accordingly. However, moving these UAVs to the center of the clusters should be done with a minimum energy consumption. Therefore, in the mobile scenario, while finding the optimal paths of UAVs, we need to determine which UAV must go to which cluster center. Next, we present the IoT devices' clustering approach for minimizing the total transmit of the IoT devices.

III. CLUSTERING IOT DEVICES

Our first step is to optimally cluster the devices and deploy the UAVs at the center of the formed clusters so as to minimize the transmit power of the ground IoT devices.

We solve the problem in (6) by exploiting the constrained K-mean clustering approach [START_REF] Höppner | Clustering with size constraints[END_REF]. In the K-mean clustering problem, given L points in R n , the goal is to partition the points into K disjoint clusters such that the sum of distances squared of the points to their corresponding cluster center is minimized. Hence, considering ( 5) and ( 6), the total transmit power of devices is minimized by placing the UAVs in the center of the optimal clusters. This problem can be solved in two steps using an iterative process. The first step is related to the assignment, and the second step is called update.

A. Assignment Step

In the assignment step, given the location of the clusters' center (given all µ j ), we find the optimal clusters for which the total distance squared between the cluster members and their center is minimized. Therefore, in our problem, given the location of the UAVs, we will first determine the optimal assignment of the devices to UAVs which can be written as:

min Aij K j=1 L i=1 A ij ||v i -µ j || 2 , (11)
s.t. L i=1 A ij ≤ M j , j ∈ K, (12)
K j=1 A ij = 1, i ∈ L, (13)
A ij ||v i -µ j || ≤ h j sin P -1 LoS (ε) , (14)
A ij ∈ {0, 1}, (15)
where v i = (x i , y i ) is the two-dimensional location vector of device i, and A ij is equal to 1 if device i is assigned to UAV j, otherwise A ij will be equal to 0. The problem presented in ( 11) is an integer linear programming which is solved by using the cutting plane method [START_REF] Garfinkel | Integer programming[END_REF].

B. Update Step

In the update step, given the clusters obtained in the assignment step, we update the location of the UAVs which is equivalent to updating the center of the clusters. Thus, the update location of UAVs is the solution to the following optimization problem:

min (xu,j ,yu,j ,hj ) i∈Cj (x u,j -x i ) 2 + (y u,j -y i ) 2 + h j 2 , (16)
s.t. (x u,j -x i ) 2 + (y u,j -y i ) 2 + h 2 j 1 - 1 sin 2 P -1 LoS (ε) ≤ 0,
for all i ∈ C j , and j ∈ K.

Theorem 1. The solution to ( 16) is s * = (x * u,j , y * u,j , h * j ) = -P (λ) -1 Q(λ), with the vector λ that maximizes the following concave function:

max λ 1 2 Q(λ) T P (λ) -1 Q(λ) + r(λ), (18) s
.t. λ ≥ 0, (19)
where Proof: As we can see from [START_REF] Garfinkel | Integer programming[END_REF], the optimization problem is a quadratically constrained quadratic program (QCQP) whose general form is given by:

P (λ) = P o + i λ i P i , Q(λ) = Q o + i λ i Q i and r(λ) = r o + i λ i r i , with P o , Q o , r o , P i , Q i ,
min s 1 2 s T P o s + Q T o s + r o , (20) s
.t. 1 2 s T P i s + Q T i s + r i ≤ 0, i = 1, ..., |C j |. (21)
Given ( 16) and ( 17), we have:

P o =   2|C j | 0 0 0 2|C j | 0 0 0 2|C j |   , P i =   2 0 0 0 2 0 0 0 ω   , ω = 1 - 1 sin 2 (P -1 LoS (ε)) , Q o = -2 |Cj | i=1 x i -2 |Cj | i=1 y i 0 T , Q i = -2x i -2y i 0 T , r o = |Cj | i=1 x 2 i + |Cj | i=1 y 2
i , and

r i = x 2 i + y 2 i .
Note that, P o and P i are positive semidefinite matrices, and, hence, the QCQP problem in (20) is convex. Now, we write the Lagrange dual function as:

f (λ) = inf s 1 2 s T P o s + Q T o s + r o + i λ i 1 2 s T P i s + Q T i s + r i = inf s 1 2 s T P (λ)s + Q(λ) T s + r(λ) .
Clearly, by taking the gradient of the function inside the infimum with respect to s, we find s * = -P (λ) -1 Q(λ).

As a result, using s * , f (λ) = 1 2 Q(λ) T P (λ) -1 Q(λ)+r(λ). Finally, the dual of problem (20) or ( 16) will be:

max f (λ), s.t. λ ≥ 0, (22)
which proves Theorem 1.

The assignment and update steps are applied iteratively until there is no change in the location update step. Clearly, at each iteration, the total transmit power is reduced and the objective function is monotonically decreasing. Hence, the solution converges after several iterations.

In summary, for given locations of the ground IoT devices, we determined the optimal locations of the UAVs (cluster centers) for which the transmit power of the IoT devices used for reliable uplink communications is minimized. In a mobile IoT network, the UAVs must update their locations and follow the cluster centers as they evolve due to time-varying dynamics. Next, we investigate how to optimally move the UAVs to the center of the clusters.

IV. MOBILITY OF UAVS: OPTIMAL TRANSPORT THEORY

Here, we find the optimal trajectory of the UAVs to guarantee the reliable uplink transmissions of mobile IoT devices. To move along the optimal trajectories, the UAVs must spend a minimum total energy on their mobility so as to remain operational for a longer time. In the considered mobile ground IoT network, the location of the devices and their availability might change, and hence, the clusters will change. Consequently, the UAVs must frequently update their locations accordingly. Now, given the location of the cluster centers obtained in Section III, and initial locations of the UAVs, we determine which UAV should fly to which cluster center such that the total energy consumed for their mobility is minimized. In other words, given I and J , the initial and new sets of UAVs' locations, one needs to find the optimal mapping between these two sets in a way that the energy used for transportations (between two sets) is minimized.

This problem can be modeled using discrete optimal transport theory [START_REF] Villani | Topics in optimal transportation[END_REF]. In its general form, optimal transport theory deals with finding an optimal transportation plan between two sets of points that leads to a minimum transportation cost [START_REF] Villani | Topics in optimal transportation[END_REF]. These sets can be either discrete or continuous, with arbitrary distributions, and can have a general transportation cost function. The optimal transport theory was originated from the following Monge problem [START_REF] Villani | Topics in optimal transportation[END_REF]. Given piles of sands and holes with the same volume, find the best move (transport map) to completely fill up the holes with the minimum total transportation cost. In general, this problem does not necessarily have a solution as each point must be mapped to only one location. However, Kantorovich relaxed this problem by using transport plans instead of maps, in which one point can go to multiple points [START_REF] Villani | Topics in optimal transportation[END_REF].

In our model, the UAVs need to move from initial locations to the new destinations. The transportation cost for each move is the energy used by each UAV for the mobility. We model this problem based on the discrete Monge-Kantorovich problem as follows [START_REF] Xia | Motivations, ideas and applications of ramified optimal transportation[END_REF]:

min Z kl l∈J k∈I E kl Z kl , (23)
s.t. l∈J Z kl = m k , ( 24
) k∈I Z kl = m l , (25)
Z kl ∈ {0, 1}, (26)
where I and J , are the initial and new sets of UAVs' locations. Z is the |J | × |I| transportation plan matrix with each element Z kl being 1 if UAV k is assigned to location l, and 0 otherwise. E kl is the energy used for moving a UAV from its initial location with index k ∈ I to a new location with index l ∈ J . m l and m k are the number of points (UAVs) at the locations with indices l and k. The energy consumption of a UAV moving with a constant speed as a function of distance is given by [START_REF] Franco | Energy-aware coverage path planning of UAVs[END_REF]:

E(D, v) = t=D/v t=0 p(v)dt = p(v) v D, (27)
where D is the travel distance of the UAV, v is the constant speed, t is the travel time, and p(v) is the power consumption as a function of speed. As we can see from (27), energy consumption for mobility is linearly proportional to the travel distance. Using the Kantorovich-duality, the discrete optimal transport problem in (23) is equivalent to:

max ϕ,ξ k∈I ξ(k) - l∈J ϕ(l) , (28) s
.t. ϕ(l) -ξ(k) ≤ E kl , (29)
where ξ : I → R, and ϕ : J → R are unknown functions of the maximization problem. The dual problem in (28) is used to solve the primal problem in (23) by applying the complementary slackness theorem [START_REF] Villani | Topics in optimal transportation[END_REF]. In this case, the optimal solution including the optimal transport plan between I and J is achieved when ϕ(l) -ξ(k) = E kl [START_REF] Villani | Topics in optimal transportation[END_REF]. Here, to find the optimal mapping between initial set of locations and the destination set, we use the revised simplex method [START_REF] Ford | Solving the transportation problem[END_REF]. The result will be the transportation plan (Z) that optimally assigns the UAVs to the destinations. Hence, the location of the UAVs are updated according to the new destinations. Subsequently, having the destinations of each UAV at different time instances, we can find the optimal trajectory of the UAVs. As a result, given the optimal paths, the UAVs are able to serve the mobile IoT devices in an energy-efficient way.

V. SIMULATION RESULTS AND ANALYSIS

In our simulations, the IoT devices are deployed within a geographical area of size 1.2 km × 1.2 km. We consider UAVbased communications in an urban environment with ψ = 11.95, and β = 0.14 at 2 GHz carrier frequency [START_REF] Hourani | Modeling air-to-ground path loss for low altitude platforms in urban environments[END_REF]. Moreover, we use the energy consumption model for UAVs' mobility as E(D, v) = D 0.95v 2 -20.4v + 130 [START_REF] Franco | Energy-aware coverage path planning of UAVs[END_REF]. Furthermore, in a time-varying network, to capture the mobility and availability of the ground IoT devices, we generate the new devices' locations by adding zero mean Gaussian random variables to the initial devices' locations. Table I lists the simulation parameters. Note that, all statistical results are averaged over a large number of independent runs.

Figure 2 shows a snapshot of the optimal devices' clustering as well as the optimal UAVs' locations resulting from our proposed approach. In this example, 5 UAVs are used to support 100 IoT devices. We assume that each UAV has a limited number of resource blocks which can be allocated to at most 30 devices. Therefore, we have 5 clusters with a maximum size of 30 devices and 5 cluster centers corresponding to the locations of the UAVs. As we can see from Figure 2, the location of IoT devices significantly impacts the number of devices per cluster and also the optimal locations of the UAVs. In this figure, the minimum and maximum cluster sizes are 3 and 27.

Figure 3 shows the total transmit power of devices versus the number of UAVs averaged over multiple simulation runs. In this figure, the performance of the proposed approach is compared with the fixed Voronoi case which is known to be a typical deployment method for static base stations. Note that, for a fair comparison, we assume that the total number of resource blocks is fixed (L), and hence, the number of resources per UAV is L K . In other words, the maximum size of each cluster will decrease as the number of UAVs increases. In the Voronoi method, assuming a uniform distribution of devices, we fix the location of UAVs at an altitude of 500 m, and then, we assign each device to the closest UAV. However, in the proposed clustering algorithm, we find the optimal clusters and deploy the UAVs at the center of the clusters. As shown in Figure 3, the proposed method outperforms the classical Voronoi scheme as the UAVs can be placed closer to the devices. As expected, increasing the number of UAVs reduces the total transmit power of IoT devices. For instance, when the number of UAVs increases from 4 to 8, the total transmit power decreases from 77 mW to 38 mW for the proposed method, and from 115 mW to 95 mW for the Voronoi case. Figure 3 shows that our approach results in about 56% reduction in the transmit power of the IoT devices.

Figure 4 shows the trajectory of one of the UAVs in a mobile IoT scenario derived from optimal transport theory. Here, we consider 8 UAVs, and 400 devices whose locations are updated at each time by adding a Gaussian random variable with N (0, 50 m) to the previous locations. Clearly, since the locations of the devices may change over time, the optimal clusters must be updated accordingly. In Figure 4, the red dots correspond to the optimal destinations of the UAV at different times. In fact, as the clusters are changing over time, the UAV uses the proposed scheme to optimally move to one of the new cluster centers.

Figure 5 shows the energy consumed by each UAV during its mobility. In this case, we use 8 UAVs for supporting 400 ground IoT devices. We consider the network at 10 time instances during which the UAVs move at a speed of 10 m/s while updating their locations. As shown in Figure 5, for the given scenario, the total amount of energy that UAVs use for mobility is around 106 kJ. Note that, this is the minimum total energy consumption that can be achieved via the optimal transport of the UAVs. As shown, different UAVs spend different amount of energy on the mobility. Depending on the optimal clustering of devices over time, different UAVs might have different travel distances to the cluster centers. For instance, UAV 1 consumes 1.8 times more energy than UAV 3. Hence, the number of UAVs may also change over time.

Figure 6 shows the energy consumption per UAV when the number of UAVs changes. Here, we assume that, initially the UAVs are optimally deployed for a given IoT system, however, after a while some of the UAVs (q UAVs) will not be operational due to the lack of battery. Consequently, the number of UAVs decreases and the remaining UAVs must update their locations to maintain the power efficiency of the ground devices. In Figure 6, for the average case, we take the average of energy over all possible combinations of removing q UAVs among the total UAVs. However, in the worst-case scenario, we remove the q UAVs whose loss leads to the highest energy consumption for the remaining UAVs. Clearly, as more UAVs become inoperational, the energy consumption of the functioning UAVs will increase. For example, when the number of lost UAVs increases from 2 to 4, the average energy consumption per UAV increases from 1520 J to 2510 J.

VI. CONCLUSIONS

In this paper, we have proposed a novel framework for efficiently deploying and moving UAVs to collect data from ground IoT devices. In particular, we have determined the optimal clustering of IoT devices as well as the optimal deployment and mobility of the UAVs such that the total transmit power of IoT devices is minimized while meeting a required bit error rate. To perform clustering given the limited capacity of each UAV, we have adopted the constrained size clustering approach. Furthermore, we have obtained the optimal trajectories that are used by the UAVs to serve the mobile IoT devices with a minimum energy consumption. The results have shown that by carefully clustering the devices and deploying the UAVs, the total transmit power of devices significantly decreases compared to the classical Voronoi-based deployment. Moreover, we have shown that, by intelligently moving the UAVs, they can remain operational for a longer time while serving the ground devices.

Fig. 1 :

 1 Fig. 1: System model.

  and r i given in the proof.

Fig. 2 :

 2 Fig. 2: Optimal clusters and UAVs' locations in one snapshot.

Fig. 3 :Fig. 4 :

 34 Fig. 3: Average of total transmit power vs. number of UAVs.

Fig. 5 :

 5 Fig. 5: Energy consumption of each UAV.

Fig. 6 :

 6 Fig. 6: Energy consumption vs. number of battery depleted UAVs.

Table I :

 I Simulation parameters.

		Parameter		Description			Value
		f c		Carrier frequency		2 GHz
		v			UAV speed			10 m/s
		δ		Bit error rate requirement		10 -8
				P LoS requirement		0.95
		N o		Noise power spectral density	-170 dBm/Hz
		R b		Transmission data rate		200 Kbps
		B	Transmission bandwidth per device		200 KHz
		η	Additional path loss to free space		5 dB
		600	UAV
			Cluster 1
		400	Cluster 2 Cluster 3
	Y-coordinate (m)	0 200	Cluster 4 Cluster 5
	-200
	-600 -400	-400	-200	0	200	400	600
				X-coordinate (m)

This work was supported by the U.S. National Science Foundation under Grants AST-1506297, by the Office of Naval Research (ONR) under Grant N00014-15-1-2709, and, by the ERC Starting Grant 305123 MORE (Advanced Mathematical Tools for Complex Network Engineering), and by the Academy of Finland.


Nicolas Houel
email: nicolas.houel@insep.fr.

Daniel Dinu

Antoine Faury

Didier Seyfried

Accuracy and reliability of the Myotest Pro system to evaluate a squat jump

Keywords: Accelerometer, squat jump, Myotest Pro, accuracy, reliability, human motion

come

Introduction

Human motion measurements present interest in the aim to improve the subject's performance and ability [START_REF] Bobbert | Human adjust control to initial squat depth in vertical squat jumping[END_REF] [START_REF] Runge | Determinants of musculoskeletal frailty and the risk of falls in old age[END_REF]. Height and velocity of the center of mass measured during vertical jump have been widely applied as functional approach to access muscular [START_REF] Bojsen-Moller | Muscle performance during maximal isometric and dynamic contractions is influenced by the stiffness of the tendinous structures[END_REF], coordination [START_REF] Bobbert | Human adjust control to initial squat depth in vertical squat jumping[END_REF] [2] [START_REF] Nagano | Optimal coordination of maximal effort horizontal and vertical jump motions -a computer simulation study[END_REF] and physical capacities [5 -8] in human motion and training sport. Kinematics (VICON, etc.) and kinetic systems (force plate, etc.) have been usually used to realize accurate measurements of the center of mass [START_REF] Bobbert | Human adjust control to initial squat depth in vertical squat jumping[END_REF] [START_REF] Crewther | Validating two systems for estimating force and power[END_REF]. But these methods require to reduce measurements at laboratory or very constraining conditions [START_REF] Crewther | Validating two systems for estimating force and power[END_REF] [10] [START_REF] Fong | A wireless motion sensing system using ADXL MEMs accelerometers for sports science applications[END_REF]. During the past five years, devices with small accelerometer have been increasingly used to measure human motions in sports [START_REF] Crewther | Validating two systems for estimating force and power[END_REF] [6] [START_REF] Jidovtseff | Validité et reproductibilité d'un dynamomètre inertiel basé sur l'accélérométrie[END_REF]. For example, Myotest SA proposes to measure specific human motion in sport training (squat jump, bench press, etc.) using Myotest Pro sensor [5] [9]. These studies focused on measurement accuracy of this tool during bench press motion [START_REF] Jidovtseff | Validité et reproductibilité d'un dynamomètre inertiel basé sur l'accélérométrie[END_REF] and squat jump [START_REF] Crewther | Validating two systems for estimating force and power[END_REF].

In these cases, the sensor was attached on the fitness bar, in the aim to limit the rotation motion [5] [9].

In these cases, low reliability and systematic bias have been already observed [5] [9]. Few studies estimated the accuracy and reliability of a sensor directly attached on the subject during motion like squat jumps [START_REF] Dionyssiotis | Assessment of musculoskeletal system in women with jumping mechanography[END_REF]. The aim of this study is to define the accuracy and reliability to estimate kinetics' variables of the centre of mass using the Myotest Pro during squat jump. Data of this sensor are compared with a reference force plate.

Methods

Nine male subjects (mean ± standard deviation SD: height = 179.8 ± 5.2 cm; mass = 76.02 ± 6.9 kg) voluntarily participated to this study. All subjects were asked to perform three squat jumps without arm movement (with their hands on the hips). All the squat jumps were performed on a force plate (Kisler, 500 Hz). For each subject, only the best squat jump was studied. The Myotest Pro system includes a triaxial accelerometer (500 Hz) sensor [START_REF] Crewther | Validating two systems for estimating force and power[END_REF] as well as software developed by Myotest SA. For each jump, the sensor Myotest Pro was fixed vertically at the side of hip's athlete, on the coxo-femoral joint in agreement with Myotest SA recommendations. In the training mode of the Myotest Pro software, it defines the vertical acceleration of the subject and then estimates the maximal velocity (Vmax). In order to determine velocity of the subject centre de mass on the vertical axis, numerical integration of the acceleration data were performed using trapezoidal rule. Velocity on the vertical axis of both measurement tools (force plate and Myotest Pro System) were used to calculate the vertical maximal velocity of the subject during the impulse phase of the squat jump (Vmax in m.s-1), the vertical velocity at take off (Vtoff, in m.s-1) and the flight time of the subject (t in s). The acceleration data on vertical axis of both tools were synchronized when the acceleration data is equal to zero before the subject takes off. A paired t-test was used to compare the significant difference between the kinetics' variables (Vmax, Vtoff, t) measured with the devices (Myotest sensor, force plate). Correlation coefficient (r) was used to estimate the relation between the kinetics' variables measured with both devices. A Bland and Altman test [START_REF] Bland | Comparing methods measurement: why plotting difference against standard method is misleading[END_REF] was used to define the accuracy and reliability between the kinetics' variables of both devices.

Results

Results showed no significant difference between the measurements of Vmax and Vtoff (Table 1). Significant difference was observed between the measurements of t (p < 0.05). Correlations between measurement for Vmax, Vtoff and t are respectively r > 0.92, r > 0.58 and r > 0.77 (Table 2). Bland & Altman test shows very low bias and high reliability (± 0.125 m.s-1) between tools for Vmax data (Figure 1). Bland & Altman test shows a significant under estimation of the bias for t data of the Myotest Pro system and low reliability for Vtoff data (± 0.35 m.s-1) (Figures 2 and3).

Discussion and Conclusion

Results show that the Myotest Pro system can be used to evaluate Vmax of subject's centre of mass during a squat jump with acceptable accuracy (error > 0.3 m.s-1) and reliability (bias < 0.1 m.s-1). The Myotest Pro estimate Vtoff with a small validity (95% limit of agreement > 0.8 m.s-1) and under estimate t with a significant different bias (> 0.03 s). So it cannot be used to estimate Vtoff and t of subject's centre of mass. Difference between results on Vtoff and t on both devices can be explained by the hypothesis of measurements of the Myotest Pro. If the centre of mass is currently used to evaluate the squat jump performance [1 -8], the Myotest Pro sensor can only estimate the acceleration of the point where it is fixed. In contrary to some study where soft development takes into account the position of the accelerometer to estimate the centre of mass kinetic's variables (Vmax, Vtoff, t) [START_REF] Houel | Accuracy and reliability of the Mensens system to evaluate a squat jump[END_REF]; the Myotest Pro system only estimates the hip motion. In conclusion, the Myotest Pro can be used only to estimate Vmax of the center of mass during a squat jump, or Vtoff and t of the hip where it is fixed. Force plate or other sensors [START_REF] Houel | Accuracy and reliability of the Mensens system to evaluate a squat jump[END_REF] could be preferred to estimate kinetics variables of the centre of mass.

Fig. 1 .Fig. 2 .Fig. 3 .

 123 Fig. 1. Bland and Altman plot depicting the limits of bias (green) between the two devices of measurement and the 95% limits of agreement (red) for Vmax.

Table 1 . Mean and standard deviation values of Vmax, Vtoff and t between Myotest Pro and force plate measurements.

 1

	Variable	Force plate	Myotest Pro
	Vmax	2.4±0.2	2.46±0.21
	Vtoff	2.21±0.22	2.2±0.22
	t	0.45±0.05	0.42±0.04*
	*p < 0.05 Myotest vs. force plate

Table 2 . Significant correlation between Myotest Pro and force plate measurements.

 2

	Variable	R Value
	Vmax	0.92**
	Vtoff	0.57
	t	0.74*

*p < 0.05; **p < 0.01


Stefan Mijovic

Luca Sanguinetti

Chiara Buratti

Mérouane Debbah

On the Optimum Number of Cooperating Nodes in Interfered Cluster-Based Sensor Networks

This paper presents a cooperative multiple-input multiple-output (MIMO) scheme for a wireless sensor network consisting of inexpensive nodes, organised in clusters and transmitting data towards sinks. The transmission is affected by hardware imperfections, imperfect synchronisation, data correlation among nodes of the same cluster, channel estimation errors and interference among nodes of different clusters. Within this setting, we are interested in determining the number of nodes per cluster that maximises the energy efficiency of the network. The analysis is conducted in the asymptotic regime in which the number N of sensor nodes per cluster grows large without bound. Numerical results are used to validate the asymptotic analysis in the finite system regime and to investigate different configurations. It turns out that the optimum number of sensor nodes per cluster increases with the inter-cluster interference and with the number of sinks.

I. INTRODUCTION

Wireless sensor networks (WSNs) have recently gained increasing attention as a practical technology being introduced to different applications. A considerable number of these applications require transmission of the acquired data over long distances using transmission resources available only at sensor nodes. In this situation, direct transmission from a source node to a sink over a fading channel often presents difficulties mainly due to the large amount of energy required to establish a reliable transmission, fostering an inefficient use of batteries.

Multiple-input multiple-output (MIMO) systems are well known for their capability of achieving high spectral efficiency in the presence of fading channels [START_REF] Telatar | Capacity of Multi-antenna Gaussian Channels[END_REF]. However, the need to install multiple antennas in sensor nodes can be problematic for economic and practical reasons. To extend the advantages of MIMO systems to single antenna devices, the idea of deploying a cooperative (also known as virtual or distributed) MIMO architecture appears to be very promising. Many works in the literature deal with cooperative schemes for WSNs: one of the first studies was presented in [START_REF] Dohler | A resource allocation strategy for distributed MIMO multi-hop communication systems[END_REF], where nodes cooperate to establish virtual antenna arrays; [START_REF] Del Coso | Cooperative distributed mimo channels in wireless sensor networks[END_REF] presents a multi-hop cooperative WSN, minimising end-to-end outage probability; in [START_REF] Barriac | Distributed beamforming for information transfer in sensor networks[END_REF] and [START_REF] Tessema | Distributed beamforming with close to optimal number of nodes for green wireless sensor networks[END_REF] cooperative beamforming is considered. Despite its promises, the deployment of a cooperative MIMO architecture in WSNs poses several technical challenges mainly because This research is supported by the FP7 NEWCOM# (Grant no. 318306). L. Sanguinetti was supported by the ERC StG MORE (Grant no. 305123) and by the research project 5GIOTTO funded by the University of Pisa.

of the large amount of signalling packets required to enable cooperation among different nodes.

Many works in the literature deal with the optimisation of cluster size in non-cooperative WSNs (see, e.g., [START_REF] Abusaimeh | Energy-aware optimization of the number of clusters and cluster-heads in wsn[END_REF]). Some works address this issue considering cooperative schemes: [START_REF] Tessema | Distributed beamforming with close to optimal number of nodes for green wireless sensor networks[END_REF] derives a close-to-optimal number of nodes and a selection method for distributed beamforming; in [START_REF] Shin | Cooperative communication with joint optimization of cluster size and resource allocation in wireless sensor networks[END_REF], the optimal cluster size minimizing the outage probability under a Rayleigh fading channel is derived for a cooperative WSN. The aforementioned works, as many others in the literature, do not account for many realistic aspects of WSNs, such as hardware limitation and imperfect synchronisation.

To overcome the issue of overhead in cooperative WSN, this work presents a simple mechanism, which does not require any signalling and is suitable for cheap sensor nodes with limited hardware capabilities. More precisely, we consider a WSN organised in clusters, wherein nodes of each cluster cooperate to transmit data to one or more sinks located in the same cluster [START_REF] Agrakhed | Cluster based energy efficient qos routing in multi-sink wireless multimedia sensor networks[END_REF]. This multi-cluster scenario is analysed by employing the Wyner model [START_REF] Wyner | Shannon-theoretic approach to a Gaussian cellular multipleaccess channel[END_REF], which simplifies the analysis and allows for intuitive interpretation of results. Furthermore, in order to facilitate the mathematical computation, we perform the analysis of the asymptotic regime where the number of sensor nodes grows without bound. Results, confirmed by simulations, show that the asymptotic analysis and subsequent optimisation are valid even for relatively small number of sensor nodes.

The main contribution of this work is the analysis and the optimisation of a cooperative WSN affected by inter-cluster interference, imperfect synchronisation and channel estimation, hardware impairments and data correlation within clusters. We provide insights into the most energy efficient cluster configuration, i.e., number of nodes and sinks per cluster, and the system performance under previously mentioned realistic phenomena.

The remainder of this work is organised as follows. 1 Section II describes system and signal models whereas the problem formulation is illustrated in Section III. The asymptotic analysis is presented in Section IV. Numerical results are given in Section V while conclusions are drawn in Section VI.

II. SYSTEM AND SIGNAL MODELS

Consider a WSN composed of L clusters, each consisting of N sensor nodes and K sinks (see Fig. 1). All sensor nodes and sinks are equipped with a single antenna and operate over a bandwidth B. A double index notation is used to refer to each node or sink in a given cluster. Under this convention, we refer to node ni as "node n in cluster i" and to sink kl as "sink k in cluster l".

A. Channel model

We consider a block flat-fading channel with coherence time T and assume that the transmission takes place according to the time-division duplex (TDD) protocol shown in Fig. 2, with T s being the time required to transmit a symbol and M being the number of symbols transmitted within a frame. As seen, the transmission phase is preceded by a training phase during which pilots, of length τ symbols, are sent by sinks with 0 ≤ τ ≤ M . Pilots enable sensors to estimate the channels. The TDD protocol is assumed to be matched to the coherence time (i.e., M T s ≤ T ). Therefore, the channel can be considered as reciprocal and the sensors can make use of pilot-based estimates for data transmission. Let us denote as h nikl the channel coefficient between node ni and sink kl and assume that

h nikl = d nikl w nikl (1)
where w nikl ∼ CN (0, 1) is the small-scale fading channel and d nikl accounts for the pathloss. For the sake of compactness, we call

h ikl = [h 1ikl . . . h N ikl ] T , H ikl = diag {h 1ikl . . . h N ikl }, w ikl = [w 1ikl . . . w N ikl ]
T and W ikl = diag {w 1ikl . . . w N ikl } to denote the channel/fading vector/matrix between all nodes in cluster i and sink kl.

To facilitate the analysis, the Wyner model is employed [START_REF] Wyner | Shannon-theoretic approach to a Gaussian cellular multipleaccess channel[END_REF]: • Interference from all the other clusters is neglected; where 0 < α < 1 is the parameter defined by the proximity between the clusters (see Fig. 1).

• All

B. Signal model

We assume that sensors within a cluster l measure the same parameter such that the transmit data vector

s l = [s 1l , . . . , s N l ] T has correlation matrix C = E[s l s H l ], with [C] n,n = E[|s nl | 2 ] = 1 ∀n and [C] n,m = E[s * nl s ml ] = c ∀n = m,
where c defines the level of correlation between data of different nodes within the same cluster. On the other hand, data from different clusters are assumed to be uncorrelated.

We denote by v nl the precoding coefficient used by node nl and assume that it is computed as:

v nl = 1 √ λ nl K k=1 h nlkl (2)
where λ nl is chosen such that the following average constraint

E{|v nl | 2 } = 1 is satisfied. Therefore, from (1) it follows that λ nl = K k=1 d nlkl .
In a more compact form, we may write

V l = diag{v 1l , . . . , v N l } as V l = Λ -1/2 l K k=1 H lkl (3)
with

Λ -1/2 l = diag{λ -1/2 1l , . . . , λ -1/2 N l }.
Let Vl be an estimate of V l and assume that hardware impairments (such as nonlinearities in amplifiers, clock drifts, I/Q imbalance in mixers, finite-precision ADCs and so forth) affect transmission. Similar to [START_REF] Zhang | A General Framework for Transmission with Transceiver Distortion and Some Applications[END_REF]- [START_REF] Björnson | Massive MIMO Systems With Non-Ideal Hardware: Energy Efficiency, Estimation, and Capacity Limits[END_REF], we model the hardware impairments as a reduction of the original signal by a factor (1ǫ 2 ) (where ǫ is related to error vector magnitude) and replacing such a loss with Gaussian distortion noise with same power. We denote by Φ l the matrix that describes imperfect synchronisation among nodes of cluster l. Then, the signal received at sink kl takes the form:

y kl = h H lkl Φ l p (1 -ǫ 2 ) Vl s l + η l + L i=1,i =l y ikl + n kl (4)
where p is the transmit power, n kl ∼ CN (0, σ 2 ) is the thermal noise while η l = √ pǫ Vl ξ l accounts for non-ideal hardware.

We assume ξ l ∼ CN (0, I N ) such that the distortion noise at sensor nl is distributed as η nl ∼ CN (0, pǫ 2 |v nl | 2 ). The term y ikl accounts for the interference generated by cluster i at sink kl given by

y ikl = h H ikl Φ i p (1 -ǫ 2 ) Vi s i + η i . (5)

III. PROBLEM STATEMENT AND ENERGY-EFFICIENCY OPTIMISATION

One of the most common way to define the EE is as a benefit-cost ratio, where the service quality is compared with the associated energy costs. The objective of this work is to solve the following problem:

max N ∈Z+ EE l = 1 -τ M K k=1 r kl P T l . ( 6
)
where r kl denotes the achievable rate (in bit/s) at sink kl and P T l accounts for the total consumed power (in Joule/s) in cluster l. The factor 1 -τ M accounts for pilot overhead. The total consumed power is computed as

P T l = 1 - τ M P TX l + τ M P CE + N P SEN + KP SINK (7)
where P TX l accounts for the power consumption of the sensor nodes in transmission state, P CE of the pilot transmission phase whereas P SEN and P SINK are constant quantities accounting for the fixed power consumption required by each sensor and sink, respectively, for running the circuitry.

In the sequel, we show how to model and compute all the terms in EE l when the transmission is affected by intercluster interference, imperfect synchronisation and channel state information (CSI) as well as hardware impairments.

A. Pilot-based Estimation of Precoding Coefficients

Observe that d nlkl corresponds to the long-term average channel attenuation, which changes in time some orders of magnitudes slower that the fast fading component w nikl . In practice, this means that d nlkl is constant for a sufficiently large number of reception phases to be accurately estimated at the sensor. For this reason, in all subsequent discussions we assume that the quantities {d nlkl ; ∀k, l} are known at sensor nl. Therefore, we are only left with the estimation of K k=1 h nlkl . In the training phase, we assume sinks of the same cluster to be perfectly synchronised. This can be justified by the fact that sinks can be complex devices employing sophisticated synchronisation mechanisms. Moreover, relatively low number of sinks per cluster guarantees the feasibility of the procedure and limits the overhead. The pilot signal transmitted by sinks in cluster l can be represented by a deterministic vector u l ∈ C τ ×1 with elements of power p τ . Therefore, we have that P CE = Kp τ /µ SINK , where 0 < µ SINK ≤ 1 accounts for the transceiver efficiency of sinks. We assume that pilot sequences used in different clusters are mutually orthogonal and that the pilot reuse factor is such that the so-called pilot contamination effect is negligible. The collective received signal x nl ∈ C τ ×1 at sensor nl is given by

x nl = K k=1 h nlkl u l + n nl = ν nl u l + n nl (8)
where n nl ∼ CN (0, ς 2 I N ) represents the additive noise at node nl during the pilot signalling. To keep the complexity of nodes at a tolerable level, we employ the least-squares estimator of ν nl defined as

νnl = 1 τ p τ u H l x nl . (9)
The variance of the estimation error is given by E{|ν nlνnl | 2 } = ς 2 τ pτ . Plugging (8) into (9) and using (2) yields

Vl = V l + E l ( 10
)
where

E l ∼ CN 0, 1 λ nl ς 2
τ pτ I N is the diagonal estimation error matrix.

Observe that a single pilot signal (i.e., τ = 1) from all sinks would be sufficient to estimate the precoding coefficients at all sensors. This is a consequence of the adoption of the precoding scheme in [START_REF] Dohler | A resource allocation strategy for distributed MIMO multi-hop communication systems[END_REF], which requires sensor nl to have only knowledge of the composite channel K k=1 h nlkl . Different precoding schemes based on knowledge of {h nlkl ; k = 1, . . . , K} would require τ ≥ K. This might not be possible when K is relatively large.

B. Synchronisation Errors

In any distributed system (such as the one considered in this work), nodes within a cluster cannot be assumed to be perfectly synchronised in time. The cause for imperfect synchronisation ranges from hardware to communication protocol limitations. The transmissions of nodes will be dispersed in time. Considering as a reference time the target time instant of transmission, t 0 , the actual transmission instants of each node can be modelled as a random variable t nl , uniformly distributed in range [-t max /2, t max /2], where t max is the maximum synchronisation error represented as a fraction of symbol time T s . Assuming an OFDM-based system, i.e., WiFibased WSN, the dispersion in time domain can be represented by phase shift in frequency domain [START_REF] Balan | AirSync: Enabling Distributed Multiuser MIMO With Full Spatial Multiplexing[END_REF] φ nl = e j2πt nl . The overall effect of synchronisation error is expressed as a complex diagonal matrix Φ l = diag{φ 1l , . . . , φ N l } which multiplies the channel matrix H ll . Synchronisation error only introduces phase shift implying that the elements of matrix Φ l are complex numbers with unit modulus. Matrix Φ l is unknown to the transmitters, like the channel estimation error, therefore these effects cannot be compensated.

We are interested in observing separately real and imaginary part, R{φ nl } = cos (2πt nl ) and I{φ nl } = sin (2πt nl ). First and second order moments of the two random variables are given by:

E[R{φ nl }] = sinc (t max ) (11) E[|R{φ nl }| 2 ] = 1 2 (1 + sinc (2t max )) (12)
and

E[I{φ nl }] = 0 (13) E[|I{φ nl }| 2 ] = 1 2 (1 -sinc (2t max )) . ( 14
)
where sinc(x) = sin(πx) πx . The effect of imperfect synchronisation is twofold: it reduces the useful signal and introduces an interference-like term due to unmatched phases.

γ kl = p(1-ǫ 2 )dll K A p(1-ǫ 2 )d ll K B + p(1-ǫ 2 )d ll K C(α) + p K ς 2 τ pτ D(α) + pǫ 2 d ll K E(α) + σ 2 (15) EE l = 1 -τ M B K k=1 log 2 (1 + γ kl ) 1 -τ M p µSEN 1 K tr{ K k=1 W H lkl K k=1 W lkl } + ς 2 τ pτ N d ll + τ M K pτ µSINK + N P SEN + KP SINK (16) γ l (α) = 1 -ǫ 2 d ll c sinc 2 (t max ) (1 -ǫ 2 ) d ll (1 -sinc (2t max )) + ς 2 τ pτ (1 + α + α 2 ) + ǫ 2 d ll (K + 1) + Kd ll (α + α 2 ) (17)
C. Energy Efficiency Plugging ( 3) and ( 10) into ( 4) and exploiting properties of Wyner model (described in Section II-A) one gets 2

y kl = p (1 -ǫ 2 ) d ll K w H lkl R{Φ l }W lkl s l + p (1 -ǫ 2 ) d ll K w H lkl I{Φ l }W lkl s l + p (1 -ǫ 2 ) d ll K w H lkl Φ l K m=1,m =k W lml s l + p (1 -ǫ 2 ) d ll w H lkl Φ l E l s l + ǫ pd ll w H lkl Φ l 1 K K m=1 W lml + E l ξ l + y ikl + y jkl + n kl . ( 18
)
The achievable rate at sink kl is thus given by r kl = B log 2 (1 + γ kl ) where γ kl is Signal to Interference and Noise Ratio (SINR) computed as in [START_REF]Texas Instruments CC2530 -Second Generation System-on-Chip Solution for 2[END_REF] with

A = w H lkl R{Φ l }W lkl CW H lkl R{Φ l }w lkl (19) B = w H lkl I{Φ l }W lkl CW H lkl I{Φ l }w lkl (20) C(α) = w H lkl Φ l m =k W lml C m =k W H lml Φ H l w lkl + αw H ikl Φ i K k=1 W iki C K k=1 W H iki Φ H i w ikl + α 2 w H jkl Φ j K k=1 W jkj C K k=1 W H jkj Φ H j w jkl (21) D(α) = w H lkl w lkl + αw H ikl w ikl + α 2 w H jkl w jkl (22) E(α) = w H lkl K k=1 W lkl K k=1 W H lkl w lkl + αw H ikl K k=1 W iki K k=1 W H iki w ikl + α 2 w H jkl K k=1 W jkj K k=1 W H jkj w jkl ( 23
)
2 Terms y ikl and y jkl can be expanded in a similar way.

being coefficients depending on fading and synchronisation error only. Regarding the power consumption, starting from (4) and using assumptions from II-A, we obtain

P TX l = p µ SEN 1 K tr{ K k=1 W H lkl K k=1 W lkl } + ς 2 τ p τ N d ll . (24)
where 0 < µ SEN ≤ 1 accounts for the transceiver efficiency of sensor nodes. Putting all the above results together, the EE takes the form in [START_REF]IEEE Standard for Local and metropolitan area networks-Part 15.4: Low-Rate Wireless Personal Area Networks (LR-WPANs)[END_REF]. At this point one can find EE-optimising N by performing an exhaustive search over the integer set, i.e., performing Monte Carlo simulations of ( 16) for each possible value of N . However, to gain more intuitive insights in the solution of ( 6) we perform the asymptotic analysis, described in the next section.

IV. ASYMPTOTIC ANALYSIS

The analysis is conducted in the regime in which the number of sensors N is infinitely large. Lemma 1. If N grows without bound, then 1 N γ lγ l → 0 almost surely with γ l given by (17). Moreover, we have that 1 N P TX l -P TX l → 0 almost surely with

P TX l = p µ SEN 1 + ς 2 τ p τ 1 Kd ll . (25)
Proof. The results easily follow using simple statistical arguments and asymptotic results.

Lemma 1 shows that the SINR and the transmit power increase linearly with N . Although valid for N growing without bounds, next we use this result for a system with a large but finite number of sensors. This yields

EE l = 1 -τ M BK log 2 (1 + γ l (α) N ) κ + N θ l (26)
where we have defined (for notational compactness)

θ l = 1 - τ M p µ SEN 1 + ς 2 τ p τ 1 Kd ll + P SEN (27)
and We now look at the EE-optimal value of N when α is given.

κ = KP SINK + τ M K p τ µ SIN K . (28)
Lemma 2. For α given, the value of N maximising (26) is given by

N ⋆ = e (z ⋆ +1) -1 γ l ( 29
)
where

z ⋆ = W γ l (α) κ θ l e - 1 e (30)
and W (x) is the Lambert function defined by the equation t = W (t)e W (t) for any t ∈ C.

Proof. The proof relies on using the same augments of Theorem 2 in [START_REF] Bjornson | Optimal design of energy-efficient multi-user MIMO systems: Is massive MIMO the answer?[END_REF]. Let EE = g log(1+bN ) c+dN denote the objective function in (26). Note that ∂EE/∂N = 0 if and only if

1 ln(2) b(c + dN ) 1 + bN -d log(1 + bN ) = 0 (31)
or, equivalently,

bc -d 1 + bN = d ln(1 + bN ) -1 . (32)
Plugging z = ln(1 + bN ) -1 into (32) yields bc de -1 e = ze z whose solution is eventually found to be z ⋆ = W ( bc de -1 e ) where W (•) is the Lambert function. Since z ⋆ = ln(1 + bN )-1, the result in (29) follows.

A close inspection of (29) reveals that N ⋆ increases with α. This is because higher interference level brings to lower SINR, with ensuing reduction of the achievable rate. This can only be compensated by increasing N ⋆ . Another interesting observation is that N ⋆ must increase with K. This is due to the fact that when K increases, the energy consumption of sinks becomes the dominant component of the overall energy consumption, due to the high transmit power of pilots and complex circuitry of sinks (P SINK ). In this setting, increasing N , while being beneficial for the achievable rate, does not affect the total energy consumption too much, implying that EE-optimal N should be higher.

V. NUMERICAL RESULTS

Monte Carlo simulations have been used to validate the analysis above. Results are obtained averaging over 1000 realisations of all the random variables, such as fading, noise, etc. The network parameters are given in Table I. Those related to the energy consumption are taken from the datasheet of a popular WSN device, TI CC2530 [START_REF]Texas Instruments CC2530 -Second Generation System-on-Chip Solution for 2[END_REF], while PHY layer parameters are inspired by the IEEE 802.15.4 standard [START_REF]IEEE Standard for Local and metropolitan area networks-Part 15.4: Low-Rate Wireless Personal Area Networks (LR-WPANs)[END_REF]. Unless otherwise specified, the level of correlation of data within a cluster is c = 0.8 and the maximum synchronisation error is t max = 0.5 (symbol time). Fig. 3 plots the energy efficiency as a function of N for K = 3 and different values of α. The curves obtained from the asymptotic results of Lemma 1 closely follow the ones obtained through Monte Carlo simulations. This proves that the asymptotic analysis is accurate even for a relatively low number of sensor nodes. As expected, the energy efficiency decreases with α due to the higher level of interference.

Fig. 4 illustrates the EE-optimal value of sensor nodes as a function of α for K = 1, 3 and 5. Firstly, notice that the closed form solution obtained through (29) is very close to the values obtained through simulations. The difference between the two is higher for higher values of K. However, for high values of K, i.e., K ≥ 5, a wider range of values of N achieves similar energy efficiency (see Fig. 5) such that N ⋆ computed by (29) achieves energy efficiency very close to the optimal one. Secondly, as predicted at the end of Section IV, N ⋆ increases with α and K.

Fig. 5 reports energy efficiency as a function of N for α = 0.5 and different values of K and t max . We can observe a large performance gap between the case where synchronisation is perfect, t max = 0, and the case where the maximum synchronisation error is t max = 0.5. The global maximum of the energy efficiency is obtained for K = 1, while the value of N ⋆ depends on the level of desynchronisation among nodes. Fig. 5 also gives an insight about the solution of the reverse problem: if N is given, find the EE-maximising value of K. These remarks would not change for other values of α and c.

VI. CONCLUSIONS

In this paper, we focused on the maximisation of the energy efficiency in a cooperative MIMO scheme for WSN in which low-cost sensor nodes are organised in clusters and transmit data to sinks. Wyner model was considered for analytic tractability. The analysis took into account hardware limitations that typically affect low-cost devices, imperfect synchronisation and data correlation among nodes of the same cluster, imperfect CSI as well as inter-cluster interference. The asymptotic analysis was used to compute a closed form expression for the EE-optimal number of nodes per cluster, N ⋆ . This allowed to get some insights on how N ⋆ is affected by the network parameters. In particular, it turned out that N ⋆ must increase with the inter-cluster interference level as well as with the number of sinks in each cluster. We also provided an insight on how imperfect synchronisation affects the performance, however this aspect is yet to be analysed in details. Monte Carlo simulations were eventually used to validate the asymptotic analysis, which is proved to be accurate even for a relatively low number of sensor nodes.

Fig. 1 .

 1 Fig. 1. Reference scenario and Wyner model.

Fig. 2 .

 2 Fig. 2. Time division duplex protocol.

5 Average α = 1 Asymptotic α = 0 5 Asymptotic α = 1 Fig. 3 .

 510513 Fig. 3. Energy efficiency vs. N for K = 3 and α = 0, 0.5 and 1.

Average K = 1 Average K = 3 Average K = 5 Asymptotic K = 1 Asymptotic K = 3 Asymptotic K = 5 Fig. 4 .K = 1 , 5 K = 5 , 5 K = 9 , 5 K = 1 ,K = 9 , tmax = 0 Fig. 5 .

 13513541555951905 Fig.4. EE-optimal value of sensor nodes N ⋆ vs. α for K = 1, 3 and 5.

Matrices and vectors are denoted by uppercase boldface and lowercase boldface letters respectively, {•} H denotes the Hermitian operator, tr{•} denotes the trace of a matrix, diag{•} denotes a diagonal matrix while Z + represents the set of all strictly positive integer numbers.