REPFS: Reliability-Ensured Personalized Function Scheduling in Sustainable Serverless Edge Computing

In recent years, serverless edge computing has been widely employed in the deployments of Internet-of-things (IoT) applications. Despite considerable research efforts in this field, existing works fail to jointly consider essential factors such as energy, reliability, personalized user requirements, and stochastic application executions. This oversight results in an inefficient utilization of computation and communication resources within serverless edge computing networks, subsequently diminishing the profit of service providers and degrading the quality-of-experience (QoE) of end users. In this paper, we explore the problem of reliability-ensured personalized function scheduling (REPFS) to jointly optimize the profit of service providers and the holistic QoE of end users in sustainable serverless edge computing. A personality-driven user QoE prediction method is first designed to accurately estimate the QoE of individual end users with differentiated personality types. Afterward, a deterministic function scheduling policy is developed on the problem-specific augmented non-dominated sorting genetic algorithm II (PSA-NSGA-II). Given the inherent uncertainty of application executions, a stochastic function scheduling strategy that can be easily parallelized for modern multicore scheduler platforms is also devised to accelerate solution generation for stochastic applications. Experimental results show that our deterministic function scheduling policy achieves 15% performance enhancement compared with representative multiobjective evolutionary algorithms. Furthermore, our stochastic function scheduling strategy promotes the service profit by 78% and the holistic user QoE by 118% on average compared with the developed deterministic scheduling policy.

reality, the demand for low-latency computation has been explosively growing.Edge computing offers a fantastic opportunity to meet computation requirements through deploying networking, control, computing, and storage infrastructures between end users and cloud servers [1], [2].Recently, the concept of serverless computing has attracted widespread attention from the edge computing community.It is oftentimes interchangeable with the term of function-as-a-service (FaaS), and was originally developed for the cloud computing landscape as a progressive step forward on service provisioning models [3].Motivated by the potentials of serverless computing, the investigations on the integration of serverless computing and edge computing, i.e., serverless edge computing, have become a research frontier in both academia and industry.
In serverless edge computing, a hot topic is to deal with the scheduling of latency-constrained IoT applications that contain dependent functions to maximize the profit of service providers.In practice, this is a challenging task mainly due to the personalized user requirements regarding the service payment and the application latency.However, the personalities of end users are completely neglected in all state-of-the-art works (e.g., [4], [5], [6], [7], [8]) on optimizing the profit of service providers in serverless edge computing.Besides, existing profit optimization techniques merely focus on the interests of service providers yet neglect the concerns of partial end users.From the perspective of service providers, they greedily intend to maximize the overall profit of sold services to end users.To this end, the applications with higher payments are more likely to be scheduled earlier compared with other applications with lower payments in conventional function scheduling methods.As a result, the qualityof-experience (QoE) of these end users with lower payments is oftentimes not guaranteed, thereby potentially degrading the holistic service level of a serverless edge computing network.
Furthermore, both energy and reliability issues have imposed increasingly critical challenges on serverless edge computing.On one hand, as the demand for carbon neutrality by 2050, it is imperative to take decisive action by embracing and promoting sustainable serverless edge computing [9], [10], [11], [12], [13], [14].On the other hand, it is crucial to address reliability issues in the context of serverless edge computing.While serverless architectures offer numerous advantages, they also introduce new challenges related to reliability.Ensuring reliable execution of functions in serverless edge computing environments requires robust reliability augmentation mechanisms.However, most of the existing reliability optimization techniques (e.g., [15], [16], [17], [18], [19]) are only suitable for corner cases of application executions, e.g., the best or worst cases, which potentially results in the low utilization of computation and communication resources in serverless edge computing.
To our best knowledge, there is a lack of investigation on the tradeoff between the profit of service providers and the holistic QoE of end users.Furthermore, the design of serverless edge computing networks should give careful consideration to key factors such as user personality, energy, reliability, and stochastic application executions.This paper is motivated by the need to address these gaps.We present the first study that jointly optimizes the profit of service providers and the holistic QoE of end users while considering user personality, energy, reliability, and stochastic application executions.

B. Technical Challenges
We list the pivotal technical challenges for reliability-ensured personalized function scheduling (REPFS) in sustainable serverless edge computing environments as follows.
Challenge 1. Establishing the Connection Between User QoE and Personalities: In practical scenarios, serverless edge computing should cater to the personalized user requirements, e.g., the tradeoff between the service payment and application latency.For example, certain end users may prioritize minimizing application completion time and are thus willing to pay higher service fees to achieve this objective.On the other hand, there also exist end users who prioritize cost savings and are willing to accept a longer, yet still acceptable, application completion time.Generally, the concept of personality can distinguish end users based on their unique characteristics.It is reported that personality has a significant influence on how individuals perceive and respond to various situations.Therefore, we reasonably believe that the user QoE is strongly associated with their personalities in serverless edge computing.However, forming this relationship is a multifaceted process necessitating thorough investigation into how disparate personality traits affect an individual's QoE perception.

Challenge 2. Accommodating the Stochastic Characteristic of Application Executions in the Design of Reliability-Aware Function Scheduling Policies:
The design of reliability-aware function scheduling policies in serverless edge computing confronts the challenge of dealing with the stochastic characteristic of application executions.Since real-world IoT applications often incorporate conditional commands and operations, distinct inputs to such applications would likely result in dissimilar completion time [20], [21].Traditional reliability optimization techniques commonly rely on the metric of worst-case execution time (WCET) to pessimistically estimate the actual execution time of application instances.However, this approach may result in the underutilization of computation and communication resources in serverless edge computing.
Challenge 3. Enhancing the Efficiency of Function Schedule Generation: From the perspective of service providers, it is vital to achieve a tradeoff between the service profits and holistic service levels.Given the target of this biobjective optimization, it is beneficial to leverage the power of non-dominated sorting genetic algorithm II (NSGA-II).The NSGA-II, first introduced in [22], is a prevalent evolutionary algorithm for solving common problems with multiple contradicting objectives.Nonetheless, our observations suggest that a direct utilization of the standard NSGA-II can be excessively time-consuming, largely attributed to the substantial quantity of chromosomes required for population evolution.Hence, an imperative necessity arises to design an augmented NSGA-II method specifically tailored for our problem.Furthermore, in light of the prevalence of modern multicore scheduler platforms, it is meaningful to investigate the parallelization acceleration of function scheduling policies on such platforms.

C. Contributions
In this paper, we explore the REPFS problem to balance the profit of service providers and the holistic QoE of end users with stochastic IoT applications in sustainable serverless edge computing.We make the following main contributions.r A personality-driven user QoE prediction method that can accurately estimate the QoE of individual end users under differentiated personality types.
r A deterministic function scheduling policy developed on the problem-specific augmented NSGA-II, i.e., PSA-NSGA-II, for sustainable serverless edge computing.r A stochastic function scheduling strategy that considers the uncertainty in application executions and can be implemented in parallel for multicore scheduler platforms.

D. Organization
The rest of this paper is organized as follows.Section II reviews related works.Section III introduces system models.Section IV formulates our problem and outlines our approach.Section V presents our personality-driven user QoE prediction method.Section VI describes a deterministic function scheduling algorithm while Section VII exhibits our stochastic parallel function scheduling strategy.Section VIII experimentally evaluate our approach, and Section IX concludes the paper.

A. Profit-Aware Function Scheduling
Numerous research works have focused on enhancing the profitability of service providers through the design of novel function scheduling policies.For example, Bermbach et al. [4] develop an auction-based approach that enforces application functions to bid on the computation resources of serverless edge computing platforms.In such a manner, the service provider can wisely decide which functions to offload for revenue maximization.In their extended work [5], the authors develop a proof-of-concept prototype called AuctionWhisk, implemented on Apache OpenWhisk, to demonstrate the practical implementation of the auction-based approach.Tutuncuoglu et al. [6] investigate the dynamic interplay between service operators responsible for storage management and pricing, and autonomous wireless devices capable of offloading computations.In addition, a Bayesian Gaussian process bandit algorithm is devised to learn an optimal service pricing policy for serverless edge computing environments.Elgamal et al. [7] conduct an investigation on the suitable pricing of serverless applications deployed on AWS Lambda.The authors first identify key factors that influence application pricing, and then introduce an efficient algorithm that explores different solutions for function fusion and placement.Pelle et al. [8] propose a framework for latencysensitive applications over serverless resources in a hybrid edge and AWS cloud scenario.By leveraging the pricing information, the framework aims to achieve optimal resource deployment and placement for performance and cost optimization.
Despite the research efforts made in [4], [5], [6], [7], [8], an important factor that has been completely neglected is the impact of user personalities on profit optimization.In practical scenarios, users with diverse personality types often perceive distinct QoE, even when their submitted applications are configured with identical computation resources on the target serverless edge computing platform.Motivated by this observation, this paper introduces a personality-driven user QoE prediction method that accurately estimates the QoE of individual users with differentiated personality types.

B. Energy-Aware Function Scheduling
In the literature, several novel function scheduling algorithms have been investigated on improving the energy efficiency of conventional grid or renewable generations powered serverless edge computing.For example, Patros et al. [9] provide an in-depth analysis of the actual power consumption patterns in serverless computing environments, using execution traces obtained from existing literature.They also propose potential strategies that can effectively minimize the energy overhead associated with serverless execution.Aslanpour et al. [10] leverage the bin packing technique to design an energy-aware function scheduling solution.Tang et al. [11] adopt the techniques of partially observable stochastic game and Q-network to alleviate the delay and energy dissipation of IoT applications.A swarm intelligence based function scheduling algorithm is developed by Xie et al. [12] to jointly reduce the service latency, energy consumption, and monetary cost.Aiming at improving the energy efficiency of IoT societal applications, Benedict Shajulin [13] presents a dedicated three-tier architecture based on serverless blockchains.Ko et al. [14] formulate the latency-constrained and energy-efficient function offloading problem as a Markov decision process and then address it by using linear programming techniques.
Notably, the energy-aware function scheduling techniques proposed in [9], [10], [11], [12], [13], [14] predominantly prioritize energy optimization while overlooking profit maximization.While these techniques effectively address energy efficiency concerns, they neglect the crucial aspect of optimizing profits for service providers.To achieve a sustainable and economically viable serverless computing environment, it is imperative to develop function scheduling policies that strike a balance between energy efficiency and profit optimization.

C. Reliability Augmentation
While there is an abundance of research on profit optimization [4], [5], [6], [7], [8] and energy efficiency optimization [10], [11], [12], [13], [14], little attention has been given to enhancing reliability in the context of serverless edge computing.Nevertheless, existing reliability optimization techniques for traditional edge computing can serve as a guide for improving the reliability in serverless edge computing scenarios.This is because both scenarios are susceptible to common failure types, including transient faults caused by electromagnetic interference and cosmic radiation [23] and bit errors during data communication [24].In the literature, researchers have explored reliabilityaware computation offloading and assignment [15], developed frameworks to minimize the block error probability during communication [16], devised strategies to tolerate transient faults and bit errors [17].Meanwhile, researchers also have dealt with reliability-aware end-end-edge collaboration for energy minimization [18], and addressed reliability-constrained application deployment problems [19] for edge computing environments.
Unfortunately, existing reliability optimization techniques [15], [16], [17], [18], [19] may not be directly applicable to serverless edge computing due to their neglect of the stochastic application executions.In real-world IoT applications, conditional instructions and operations can lead to variations in completion time for different inputs [20], [21], [25].However, most of the existing reliability optimization techniques are designed for specific scenarios, such as the best or worst cases of application executions.This narrow focus may result in underutilization of computation and communication resources in serverless edge computing.Therefore, there is a need to develop new reliability optimization techniques that consider the stochastic nature of application executions, enabling efficient resource utilization and improved reliability in serverless edge computing environments.

A. Architecture Model
As shown in Fig. 1, we consider a common serverless edge computing network whose topological structure is depicted by an undirected connected graph G = (S, L).Specifically, S = {S 1 , S 2 , . . ., S M } is a set of total M geographically distributed edge servers.The mth (1 ≤ m ≤ M ) edge server is denoted by S m , and it is equipped with three main modules, i.e., the energy harvesting module, the energy buffering module, and the energy consumption module, as plotted in Fig. 2. The energy harvesting module automatically absorbs renewable generations from ambient environments.The energy buffering module, in a general form of lithium batteries, is designed against the uncertainty in harvested renewable energy.The energy consumption module   [28], [29], [30], [31] can obtain energy supply from either the energy harvesting module or the energy buffering module for application executions.
In the real-world construction of a serverless edge computing network, each physical edge server is generally colocated with one specific base station [26], [27].Moreover, since the number of base stations is usually much larger than that of geographically dispersed base stations, it is common to place any two edge servers at different base stations rather than the same one.The commercial edge servers, such as HPE ProLiant MicroServer Gen10 [28], Dell R230 [29], Lenovo TS250 [30] and Inspur NP3020 [31], are suitable infrastructures in practical scenarios for building serverless edge computing networks.Table II lists the main parameters of these heterogeneous edge servers.Since the study on the optimal deployment of edge servers is out the scope of this paper, we consider that the network topology and server information of our serverless edge computing network are predetermined.With this consideration, we assume that the maximal computation capacity of edge server S m is fixed as a constant Ψ m .Note that this is a common assumption and also adopted by [4], [5], [6], [7], [8], [9], [10], [11], [12], [13], [14], [15], [16], [17], [18], [19].On the other hand, during the network maintenance stage, the service provider may adjust the value of Ψ m by deploying more physical edge servers into one virtual edge server at the original location.At that moment, our scheduling algorithms in Sections VI and VII should be invoked again to redetermine the desirable function schedules for accommodating the new network.
Moreover, we consider the container-based virtualization technique to package and isolate applications running on edge servers.Traditional serverless computing platforms, such as Amazon AWS Lambda and Google Cloud Functions, stipulate the minimum volume of computation resources that should be allocated to a container [3].Considering the fixed maximal computation capacity of an edge server, we suppose that edge server S m is capable of hosting up to o m active containers simultaneously at runtime.Denote the container set as Note that the value of χ φ m is not specified in advance, but is one of the optimization variables in this paper.Since we focus on the personalized function scheduling at offline stage, the value of variable χ φ m can be set to 0 at design time.In such cases, container C φ m will not be established at runtime since it is not allocated to any computation resource.As a result, the actual number of containers that will be deployed on individual edge servers is tunable at offline stage.In addition, edge server S m can communicate with its companion S n (1 ≤ n ≤ M, m = n) through a virtual link l m,n ∈ L. The bandwidth associated with virtual link l m,n is symbolized by a constant b m,n .Following [32], the condition b m,n = b n,m holds and the data communication time incurred by two containers on an identical edge server is negligible.

B. Application Model
We suppose that a group of real-time IoT applications A = {A 1 , A 2 , . . ., A I } are to execute on the target serverless edge computing network.Every application A i (1 ≤ i ≤ I) is submitted by a unique user U i , and typically modeled as a directed acyclic graph (DAG) p , w i,p , d i,p,q } is adopted to characterize the features of function f i,p .In the tuple, μ i,p ∈ (0, 1] is the activity factor exhibiting the inter-function heterogeneity in computation resource utilization.w i,p is the total number of instruction cycles, and d i,p,q is the amount of communication data from function f i,p to its direct successor f i,q .Additionally, all functions in application A i share a common deadline D i latency , a desirable response latency or completion time T * ,i latency , and a reliability goal R i goal .As indicated in [20], [21], [25], the number of instruction cycles of an application fluctuates with different inputs.As a result, both w i,p and d i,p,q are random variables considering the stochastic characteristic of application executions.Given this, we introduce a set of execution adaption variables α = {α 1 , α 2 , . . ., α I } ranging in [0,1] to capture the uncertainty in the number of instruction cycles.In this way, w i,p and d i,p,q can be represented by ( 1) and (2), respectively.
is the execution adaption variable associated with application A i .w i,p min and w i,p max are the minimal and maximal number of instruction cycles of function f i,p , respectively.Similarly, d i,p,q min and d i,p,q max are the minimal and maximal volume of communication data from function f i,p to its direct successor f i,q , respectively.Note that w i,p min , w i,p max , d i,p,q min and d i,p,q max are deemed to be constants and they have already been known at design time.Clearly, w i,p = w i,p min and d i,p,q = d i,p,q min hold if

C. Energy Model
We construct the energy model from supply and demand perspectives.Let I m harv denote the harvesting power of edge server S m at time instant t, then the renewable energy under a specific time horizon Δt is given by As aforementioned, both the energy harvesting and buffering modules can provide energy supply for application execution.
Let E m buffer represent the buffering energy sorted in the energy buffering module, the overall energy supply of edge server S m during time interval [t, t + Δt] is thus estimated by The energy dissipation of edge server S m depends not only on the computation resources allocated to individual containers, but also on the number of functions dispatched to these containers.When function f i,p is executed on container C φ m with allocated computation resource χ φ m , the energy demand incurred by function executions is expressed as [20], [21] m static is the static power, ζ m is the effective switching capacitance, and v m,φ is the supply voltage.On the other hand, let u m,n be the unit energy dissipation of communication data transfer on link l m,n , then the energy overhead of transferring data from function f i,p to all direct successors is expressed as A binary variable Λ m,i,q is introduced to indicate whether or not function f i,q is dispatched to edge server S m .If yes, Λ m,i,q is set to 1; otherwise, it is set to 0. Hence, we can readily infer the whole energy demand of edge server S m as Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
Here, ∇ i,p m,φ is also a binary variable taking the value of 1 only when function f i,p is executed on container C φ m .

D. Reliability Model
We concentrate on the transient fault at the side of edge servers during function executions and the bit error at the side of virtual links during data communication [23], [24].Let λ φ m stand for the average fault rate of container C φ m with allocated computation resource χ φ m , then it is calculated as [23] λ φ m = λ m max × 10 λ m max denotes the average fault rate corresponding to computation capacity Ψ m , and ψ φ m is a constant implying the sensitivity of fault rates to computation capacity scaling.
Generally, exponential distribution can be utilized to capture the occurrence of soft errors arising from transient faults.In this context, the execution reliability of function f i,p on container C φ m with computation capacity χ φ m yet neglecting inter-function data dependencies is modeled as [ Moreover, let ψ i,p parent = {f i, 1 , f i, 2 , . . ., f i, z } denote a set of total z direct predecessors of function f i,p , then the possibility of successfully executing all direct predecessors is given by R i, 1 nodepexe refers to the reliability of function f i, 1 with omitted container indexes for easy presentation.
suggests the conditional probability of successfully executing function f i, z when other functions f i, 1 , f i, 2 , . . ., f i, z−1 are all finished without the occurrence of soft errors.Putting together ( 9) and (10), the execution reliability of function f i,p on container C φ m considering inter-function data dependencies is hence expressed as In addition to soft errors, bit errors may happen at the side of virtual links during data communication.Let ϕ m,n denote the bit error rate at link l m,n , then the communication reliability between function f i,p deployed at edge server S m and its single direct predecessor f i, g is formulated as [ Since all communication data should be successfully transferred from direct predecessors in parent set ψ i,p parent to function f i,p , the communication reliability is inferred by In a nutshell, the overall reliability of function f i,p on container C φ m is derived by assembling execution reliability in (11) and communication reliability in (13), that is,

E. QoE Model
In our considered network, individual end users can submit their applications to edge servers for accomplishment.For each user, a shorter application completion time is desirable, but he/she will be charged more service fees by the service provider.Meanwhile, a reduced service payment is achievable at the expense of a longer application completion time yet within the deadline.However, existing works simply adopt a uniform payment model for all applications, thereby failing to capture the personalized requirements of varied end users.In this paper, we design an elastic payment model.Specifically, let Υ i payment denote the service payment of application A i , then the relationship between Υ i payment and response latency symbolled by T i latency is constructed by T Δ,i latency is a time metric to measure the violation degree of expected response latency, and can be derived by ΔT i latency = T i latency − T * ,i latency .(15) indicates that if response latency T i latency is earlier than desired finish time T * ,i latency , application A i should pay for a maximal service fee Υ max,i payment .Further, if response latency T i latency falls into the interval of (T * ,i latency , D i latency ], application A i is guaranteed to attain an acceptable execution result but the payment Υ i payment will linearly decay at the rate of coefficient δ i .However, if T i latency is later than deadline D i latency , application A i is delivered to unacceptable execution results due to timing violations, and therefore the service payment Υ i payment is set to 0. To characterize the QoE of an end user, we first define a latent variable i as where η i ∈ [0, 1] is the preference factor.At this moment, we calculate the QoE Q i of end user U i as where ϑ 1 and ϑ 2 are two parameters.Apparently, η i , ϑ 1 , and ϑ 2 play a key role in accurately estimating the QoE of an end user.Note that in our QoE model, the parameters η i , ϑ 1 , and ϑ 2 are not treated as constants, i.e., they are the variables dependent on the personalities of end users.We will discuss later the way to derive their optimal values in Section V.

A. Problem Statement
In this paper, we are dedicated to jointly optimizing the profit of service providers and the holistic QoE of end users.On one hand, since the service payment of an individual application has been defined in (15), the service profit of service providers is accordingly derived by On the other hand, the holistic user QoE can be evaluated by using the metric of mean deviation, this is, Q i max and Q i min denote the maximal and minimal values between the QoE Q i defined in (17) and the averaged QoE Meanwhile, the design constraints of reliability, inter-function dependency, energy, and server capacity should be strictly guaranteed.Formally, our problem is formulated as follows.
(21) indicates the reliability constraint, where Γ m,φ is a set of the functions scheduled on container C φ m .( 22) suggests the fulfillment of data dependency constraints for any two functions in application A i .Here, T i,q start is the start time of function f i,q , T i,p finish is the finish time of function f i,p , and T i,p,q transfer is the communication time spent on conducting data transfer between functions f i,p and f i,q .(23) imposes the energy constraint on edge server S m .( 24) restricts the allocation of computation resources to containers.To avoid wasting computation resources, (24) enforces the sum of χ 1 m , χ 2 m , . . ., χ o m m to be equal to the maximal computation capacity of edge server S m .

B. Our Approach
To solve the aforementioned problem, a novel personalitydriven user QoE prediction method is first designed in Section V to accurately estimate the QoE of individual end users under differentiated personality types.Considering the difficulty in dealing with an optimization problem with uncertainty, our approach then incorporates a deterministic function scheduling algorithm developed on the problem-specific augmented NSGA-II in Section VI.Afterward, our approach offers a stochastic function scheduling strategy in Section VII-A that takes into account the uncertainty in the number of instruction cycles and communication data volume for individual IoT applications.In order to accelerate the convergence of our strategy, a parallel function scheduling optimization is also investigated in Section VII-C for modern scheduler platforms equipped with multiprocessors or multiple graphics processing units.Through the above steps, our approach can find desirable function schedules while satisfying all design constraints.

V. PERSONALITY-DRIVEN QOE PREDICTION
In this section, we design a personality-driven user QoE prediction method.We first build the relationship between user QoE and preferences, then link user preferences with personality scores, and finally detail our QoE prediction algorithm.

A. Link User QoE With Preferences
We have modeled the QoE of an end user in (17) but introduced random variables η 1 , η 2 , . . ., η I , ϑ 1 , and ϑ 2 .To determine their optimal values, we design a questionnaire based method by linking user QoE with user preferences.Specifically, our questionnaire contains a certain number of investigation items, and each investigation item is depicted by a tuple <Participant Identity, Response Latency, Service Payment, User QoE>.The first element, "Participant Identity", records the unique identity of one selected user involved in the questionnaire.The remaining elements, "Response Latency", "Service Payment", and "User QoE", collect the user QoE under a given combination of response latency and service payment.In our scenario, the user QoE is normalized to an integer picked from the interval of [1,10], with 1 indicating the lowest QoE and 10 suggesting the highest QoE.Besides, we construct baseline and other latencypayment levels to explore the possible combinations of response latency and service payment.For the baseline latency-payment level, its response latency and service payment are both set to 100%, i.e., the combination of shortest response latency and maximal service payment.Similarly, other latency-payment levels can be easily generated.For instance, the response latency choices are increased to 110%, 120%, 130%, 140%, and 150% of the shortest response latency, while the service payment choices are scaled to 90%, 85%, 80%, 75%, and 70% of the maximal service payment.For each latency-payment level, we request individual users participating in the questionnaire to mark their corresponding QoE scores.
Let U par = {U j |1 ≤ j ≤ I par } denote a collection of total I par (0 < I par ≤ I) participants in the questionnaire.Meanwhile, let L denote the number of all possible latency-payment levels.For participant U j under the lth latency-payment level, the response latency, service payment, preference, and user QoE are represented by T j,l latency , Υ j,l payment , θ j,l , and Q j,l , respectively.At this moment, we can derive ϑ 1 , ϑ 2 , η 1 , η 2 , . . ., η I par by addressing the below mean-square-error minimization problem, that is, Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
where θ j,l = η j × T j,l latency + (1 − η j ) × Υ j,l payment .In practice, it is impractical to directly solve this mean-square-error minimization problem due to huge time overheads.Consequently, a heuristic three-step approach is adopted here.First, a group of initial preference factors η 1 , . . ., η I par is randomly generated and therefore the preference of every participant is readily derived by (16).Second, the values of ϑ 1 and ϑ 2 under fixed preference θ j,l can be determined by addressing the following linear regression problem arg min On the other hand, the preference factor η j under fixed values of ϑ 1 and ϑ 2 is calculated as In the third step, we compare the difference between randomly generated preference factors and the derived preference factors in (27).If their difference exceeds a specific threshold, another set of preference factors will be immediately produced and fed back to step one.Otherwise, the loop is hence terminated after outputting the optimal ϑ * 1 , ϑ * 2 , η * 1 , η * 2 , . . ., η * I par .

B. Link User Preferences With Personality Scores
In our context, a difficulty lies in accurately quantifying user personality types since there are too many factors influencing personality traits.In this paper, we borrow the five-factor model of personality [33] from the psychological trait theory to quantify user personality types.This model decomposes the user personalities into five distinct traits: Agreeableness (A), Conscientiousness (C), Extraversion (E), Openness (O), and Neuroticism (N).Every constituent measures a singular dimension of the user personality.For example, the trait of openness estimates the degree to which a person is intellectually curious, emotionally open, and eager for new things.The trait of extraversion can portray a person's sensitivity to the breadth of activity and the urgency from ambient environments.
Moreover, a questionnaire of ten-item personality inventory (TIPI) is conducted to precisely measure the five dimensions of personality traits.In this questionnaire, a participant is asked to answer a certain number of delicately designed questions.After the questionnaire accomplishment, we obtain a score vector associated with this participant.Let Q j = [AScore j , CScore j , EScore j , OScore j , NScore j ] represent the personality score vector of participant U j .Then, a personality score matrix Q I par ×5 = [Q j |1 ≤ j ≤ I par ] with a dimension of I par × 5 is leveraged to store the personality score vectors of all participants.Afterward, a linear mapping from personality scores to user preferences is depicted by where W I par ×1 = [η j |1 ≤ j ≤ I par ] T refers to as the preference factor matrix, and X 5×1 is a mapping matrix with the dimension of 5 × 1 that exhibits the impact of personality scores on user preferences.Evidently, once an optimal solution X * 5×1 to ( 28) is derived, the optimal preference factor matrix W * I×5 corresponding to all users under an inputted personality score matrix Q I×5 is hence constructed by This equitation indicates that we can estimate the preference factor η * i of application A i once obtaining the personality scores of user U i .Finally, the user QoE Q i,l under an arbitrary combination of T i,l latency and Υ j,l payment is therefore inferred by

C. Personality-Driven QoE Prediction Algorithm
We devise a novel personality-driven user QoE prediction approach, as outlined in Algorithm 1.The target of this algorithm is to estimate the QoE Q i of user U i under any latency-payment combination {T i latency , Υ i payment }.Specifically, lines 2-3 initialize a group of iteration flags = [ 1 , 2 , . . ., I par ] and a collection of preference factors η = [η 1 , η 2 , . . ., η I par ], respectively.Line 4 judges whether or not all the iteration flags are equal to zero.If not, line 5 estimates ϑ * 1 and ϑ * 2 by exploiting (26).Lines 6-12 iteratively derive the preference factors of all participants in the questionnaire.In this process, line 7 calculates preference factor η * j by utilizing (27).Line 8 compares the absolute difference between the randomly generated preference factor η j and the derived preference factor η * j .If the absolute difference is small enough, line 11 sets the iteration flag of participant U j to zero.Otherwise, line 12 renews matrix η = [η 1 , η 2 , . . ., η I par ] by utilizing the Latin hypercube sampling technique [34].Subsequently, a linear mapping from personality scores to user preferences is obtained by solving (28) in line 13.Based on the efforts made in lines 2-13, line 15 constructs the personality score vector S j of user U i , and line 16 derives preference factor η * i by employing (29).As a result, line 17 can predict the QoE Q i under an inputted latency-payment combination {T i latency , Υ i payment }.Finally, the entire algorithm exits before outputting QoE Q i in line 18.Refer to the online supplementary material (available at https://www.jianguoyun.com/p/DU3g1qUQ3IvkCRin4awFIAA)for the time and space complexity analysis of Algorithm 1.

VI. DETERMINISTIC FUNCTION SCHEDULING BASED ON PSA-NSGA-II
In this section, we describe our deterministic function scheduling technique.We first detail the developed PSA-NSGA-II method, and then show the deterministic function scheduling algorithm based on PSA-NSGA-II.

A. Our Tailored PSA-NSGA-II
We design an augmented NSGA-II, i.e., PSA-NSGA-II, dedicated to our studied problem, as presented below.

1) Solution Representation:
As mentioned earlier, the target of this paper is to jointly optimize the service profit and the holistic user QoE.To this end, we should optimally determine the execution order of functions, the allocation of computation resources to containers, and the dispatch of functions to containers.For facilitating chromosome representations, we introduce the function order table X order , resource-to-container allocation table X allocate , and function-tocontainer dispatch table X dispatch .Specifically, table X order has a size of I × max{ξ 1 , . . ., ξ i , . . ., ξ I }, and it used to maintain the inter-function precedence constraints.Table X allocate has a size of M × max{o 1 , . . ., o m , . . ., o M }, where each row records the resource-to-container allocation solution for an edge server.Similarly, table X dispatch has a size of I × max{ξ 1 , . . ., ξ i , . . ., ξ I }, where each element stores the index of a container preparing for running a specified function.
2) Population Initialization: In our optimization problem, we should initialize the three tables constituting an individual chromosome, i.e., the function order table X order , resource-tocontainer allocation table X allocate , and function-to-container dispatch table X dispatch .Essentially, table X order indicates the inter-function precedence constraints.Therefore, existing topological sorting techniques can be adopted to initialize X order .For every application A i , we first sort all functions in function set F i by using Kahn's topological sorting algorithm [35].Then, we orderly choose an integer number from 1 to ξ i for individual functions based on the sorting results.The above two steps can be repeatedly conducted until the functions of all applications are properly arranged.During the initialization of table X allocate , we evenly divide the whole computation resource of an edge server into multiple pieces (e.g., 10 pieces).Let Δ m denote the unit computation resource specified by edge server S m , e.g., Δ m = Ψ m /10.Then, we select an integer from 1 to 10 for element X (m,φ) dispatch is assigned to the index of one container for the execution of function f i,p stored in X (i,p) order of table X order .Note that the above operations cannot guarantee the validity of all chromosomes since the constraints of energy and reliability are completely ignored.Given this, we further infer the degree of violating either energy or reliability constraints for each chromosome by using the below two equations.
Subsequently, we can compare the constraint violating degree VioDeg(O u ) with a predefined threshold ∂.If VioDeg(O u ) ≤ ∂ holds, chromosome O u is still kept in the initial population; otherwise, it is directly eliminated from the initial population.Obviously, the initial population contains both the valid chromosomes satisfying all constraints and the chromosomes with a smaller overall constraint violation.In such a manner, the diversity of the initial population can be maintained.
3) Crossover Operator: We develop a new crossover operator based on the technique of sequential importance sampling (SIS) [36].The procedure of our crossover operator is consistent with the original SIS approach, which contains three major steps, i.e., estimation, measurement, and resampling.At the estimation step, a sample set of new crossover positions is stochastically generated based on the given probability distribution function (PDF).We assume that the crossover positions obey a Gaussian distribution N (ψ g , δ 2 g ) with mean ψ g and variance δ 2 g at gth generation of population evolution.At the measurement step, the performance of all position samples is evaluated by deriving the non-domination levels and crowding distances of resultant offspring.Suppose that there are total H position samples {w h g , x h g = {x h g,1 , . . ., x h g,K }|1 ≤ h ≤ H}, where x h g,k (1 ≤ k ≤ K) is the kth interval based on the observed data y k , and ω h g is the weight of the hth position sample x h g .Then, the weight ω h g of the hth position sample x h g is renewed up to a normalizing constant, that is, where Pro(y k |x h g ) is the posterior probability for observed data y k on position sample x h g .At the resampling step, a total of Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
are greedily reserved while the remaining samples are directly eliminated due to their poor quality.At this moment, the mean ψ g+1 and variance δ 2 g+1 utilized for the next iteration are hence updated by (35) and (36), respectively.

4) Mutation Operator:
We design a new mutation operator, namely logistic mutation operator, borrowed from the wellknown logistic population model (LPM).In our mutation operation, the native LPM cannot be directly utilized to guide the chromosome mutation since it has a strict upper boundary (i.e., natural resource capacity).The literature [37] develops an improved version of LPM for bandwidth allocation in data centers but it is quite different from our optimization objectives and unique constraints.Inspired by the literature [37], we put forward a logistic growth model for mutation operations, as outlined in ( 37)- (39).
In (37), β 1 is a constant, r k g,u is the mutation factor for the kth interval in O u g , and τ u g denotes the mutation coefficient of chromosome O u g .dr k g,u /dg represents the derivative of mutation control factor r k g,u with respect to iteration counter g.In (38), β 2 is also a constant, and dτ u g /dg represents the derivative of mutation coefficient τ u g with respect to iteration counter g.In (39), ω u g is the non-domination level of chromosome O u g , while τ max bgt and τ unit bgt are orderly the whole and unit mutation capacity specified as constants.Accordingly, the mutation probability Pro(x k g,u ) of position x k g,u can be derived by (40).
g th and g max are sequentially the threshold and maximum of iteration counter g, and β 3 a constant.Normalize(r k g,u ) is a normalization function used to scale r k g,u into [0,1].Fig. 3 shows the basic idea of our LPM based mutation control method.At the early stage of population evolution, mutation factor has a slow start speed with the generation counter.This is mainly due to the consistency considerations for the operation in the population initialization phase.Recall that we retain a group of chromosomes with smaller overall constraint violations in the population initialization phase.If the mutation probabilities of these chromosomes take larger values, the constraint violation degrees of mutated chromosomes are extremely likely to deteriorate.On the contrary, a slow start of mutation speed can potentially avoid this shortcoming and therefore facilitate the maintenance of initial population diversity.At the middle stage of population evolution, mutation factor presents an accelerated growing trend and later keeps smooth without fluctuation.As a result, the diversity from one generation to the next is significantly boosted, which prevents the evolutionary procedure from falling into local optimal solutions.At the late stage of population evolution, mutation factor is set to a small constant in order to hasten the convergence of the entire evolutionary procedure.
5) Selection: In our PSA-NSGA-II, the current population is composed of parents preserved by the last generation and offspring generated via genetic operators.It is clear that not all chromosomes are guaranteed to satisfy all design constraints due to the retention of infeasible solutions in the population initialization phase and the stochastic nature of the entire evolutionary process.On one hand, the chromosomes satisfying all design constraints are certainly desirable since they represent feasible solutions to our studied problem.On the other hand, the chromosomes with a small constraint violation may also carry useful information that guides the direction of population evolution.Considering the two aspects, we rank all chromosomes in descending order according to their non-dominated levels and crowding distances.Meanwhile, the constraint violating degree of every infeasible chromosome is derived by (33).As a result, a certain number of chromosomes satisfying all constraints and slightly violating the design constraints are selected to constitute a new generation.The variable ratio method is adopted here to make the proportion of retained infeasible individuals decrease linearly with the population iterations.

B. Deterministic Function Scheduling Using PSA-NSGA-II
We develop a deterministic function scheduling approach based on the PSA-NSGA-II, as detailed in Algorithm 2. Line 1 sets the initial iteration counter g to the value of 0. Line 2 then selects the chromosomes with smaller constraint violating degrees to constitute initial population P g = {P 1  g , P 2 g , . . ., P U g }.Line 3 initializes the operation flag for three tables of function order, function-to-container dispatch, and resource-to-container allocation.After the above initialization operations, the algorithm enters into the population evolution procedure in lines 4-36.Specifically, line 6 calls the modulo function Modulo(flag++, 3) to derive the remainder component of operation flag divided by 3. If the result is 1, the evolutionary operation is applied to the function order table; otherwise, the evolutionary operation is imposed on function-to-container dispatch and resource-tocontainer allocation tables.Accordingly, line 7 generates a collection of cut-off samples {η g,x flag |1 ≤ x ≤ U } based on gaussian distribution N (ψ g flag , δ g flag ).Line 8 sets two chromosome indices u and z to the values of 1 and U , respectively.If flag = 1 holds, the cut-off sample indices x 1 and x 2 are orderly set to u and z; otherwise, the two cut-off sample indices take an identical cut-off sample (lines 9-13).
Subsequently, line 14 exploits the splitting function Split(P u g , X g,u flag , η g,x 1 flag ) to divide table X g,u flag in chromosome P u g into two parts of X g,u 1 flag and X g,u 2 flag based on the cut-off sample η g,x 1 flag .Likewise, splitting function Split(P z g , X g,z flag , η g,x 2 flag ) is involved in line 15 to disunite table X g,z flag in chromosome P z g into two fractions of X g,z 1  flag and X g,z 2 flag based on the cut-off sample η g,x 2 flag .After division operations in lines 14-15, lines 16-17 employ the connection functions Connect(P u g , P z g , X g,u 1 flag , X g,z 2 flag ) and Connect(P u g , P z g , X g,z 1 flag , X g,u 2 flag ) to construct offsprings of chromosomes P u g and P z g .Line 18 updates chromosome indices u and z for the crossover operations of other chromosomes.Line 19 derives gaussian distribution N (ψ g+1 flag , δ g+1 flag ) of crossover samples for the next generation of population by using ( 35)- (36).If u = z holds, the crossover operations are finished (lines 19-20) while the mutation operations immediately start to run in lines 22-33.For every chromosome P u g , line 23 calculates mutation probability vectors for three tables X g,u order , X g,u allocate and X g,u dispatch by using (40).Line 24 randomly generates mutation probability thresholds σ 1 , σ 2 and σ 3 falling into [0,1].For every position x k 1 g,u in table X g,u order , if mutation probability Pro(x k 2 g,u ) is larger than threshold σ 1 , mutation order function MutateOrder(X g,u flag , x k 1 g,u ) is applied to position x k 1 g,u (lines 25-27).Similarly, for every position ) is applied to position x k 2 g,u (lines 28-30).For every position x k 3 g,u in table X g,u allocate , if mutation probability Pro(x k 3 g,u ) is larger than threshold σ 3 , mutation allocation function MutateAllocate(X g,u allocate , x k 3 g,u ) is applied to position x k 3 g,u (lines 31-33).
Once the current population evolution is accomplished, the algorithm then enters into the competitive selection stage.For evaluating the quality of individual chromosomes, line 34 uses Algorithm 1 to derive the objective values defined in (20).Afterward, line 35 adopts the non-dominated ranking method [22] to sort all chromosomes in descending order.Line 36 infers the constraint violating degrees of each chromosome by using (33).Through the efforts made in lines 34-36, line 37 invokes the population generation function MakeNewPop(P g , τ g ) to construct the next generation of population P g .As mentioned earlier, we leverage the variable ratio method to make the proportion (i.e., variable τ g ) of retained infeasible individuals in the population decreases linearly with the number of population iterations.After several rounds of population iterations, the algorithm is terminated after outputting non-dominated solution set P in line 39.Refer to the online supplementary material (available at https://www.jianguoyun.com/p/DU3g1qUQ3IvkCRin4awFIAA)for the time and space complexity analysis of Algorithm 2.

VII. STOCHASTIC PARALLEL FUNCTION SCHEDULING
In this section, we describe our stochastic parallel function scheduling strategy.We first introduce the design flow of stochastic function scheduling, then present the Monte Carlo Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

A. Design Flow of Stochastic Function Scheduling
As mentioned earlier, we introduce a set of execution adaptation variables α = {α 1 , α 2 , . . ., α I } to capture the stochastic property of application executions in real-world edge networks.Since 0 ≤ α i ≤ 1 holds for ∀i ∈ [1, I], the total number of execution cycles of function f i,p ranges from w i,p min to w i,p max , and the corresponding communication data from function f i,p to its direct successor f i,q fluctuates between d i,p,q min and d i,p,q max , as shown in ( 1)- (2).Unlike existing approaches only suitable for corner cases, our technique enables the design of stochastic scheduling scheme according to the status quo of application executions.
Fig. 4 exhibits the design flow of our stochastic function scheduling scheme.As illustrated in this figure, the execution adaptation variables α = {α 1 , α 2 , . . ., α I } are iteratively generated to derive a deterministic function schedule that confines the deadline missing ratio of all applications to predefined design requirements.In every round of search iteration, Algorithm 2 developed on the problem-specific augmented NSGA-II is utilized to produce a deterministic function schedule under designated values of execution adaption variables α = {α 1 , α 2 , . . ., α I }.Afterward, the Monte Carlo simulation technique is exploited to appraise the performance of current execution adaption variables in terms of application deadline missing ratio.If the deadline missing ratio falls into an acceptable range, the search procedure is hence immediately terminated after outputting an uncertainty-aware function schedule.Otherwise, the search procedure will update the values of execution adaption variables α = {α 1 , α 2 , . . ., α I } and move to the next round of iteration.

B. Tuning Execution Adaption Variables
During every iteration, a certain number of Monte Carlo samples of varied inputs of individual applications should be generated based on their probability distribution.As pointed out in [20], [21], [38], Gaussian distribution or exponential distribution can approximately capture the stochastic characteristic of application executions.Thus, we adopt the Gaussian distribution and exponential distribution in subsequent experimental evaluations to generate the actual number of instruction cycles of an application to be executed.In addition, as indicated in [39], a total of 10000 Monte Carlo samples are sufficient enough to cover a wide range of the entire sample space.Therefore, we generate 10000 Monte Carlo samples (i.e., application instances) for each target application.Once Monte Carlo samples under a given probability distribution are generated, we can evaluate the performance of the current function schedules by using the time metric of deadline missing ratio.Specifically, the deadline missing ratio of a function schedule is inferred by the rate of the number of desirable samples satisfying the deadlines of all applications to the total number of generated samples.If the inferred deadline missing ratio is less than a threshold, the current function schedules and associated execution adaptation variables are hence the preferent scheduling solutions and the whole iteration process exits.Otherwise, the proposed policy updates the execution adaptation variables α = {α 1 , α 2 , . . ., α I } for the next round of search iteration.
We should point out that our stochastic function scheduling approach is not restricted to the Gaussian and exponential distributions, i.e., it fits any designated types of probability distributions.This observation stems from the fact that, once the required Monte Carlo samples are generated, we can evaluate and then tune the current function schedules, irrespective of the underlying probability distribution from which these samples are drawn.

C. Parallel Function Scheduling Optimization
Essentially, the execution adaptation variables capture the stochastic property of application executions, thus, can be utilized to estimate the actual execution cycles of individual functions in each application.Given this fact, we can derive the optimal values of execution adaption variables α = {α 1 , α 2 , . . ., α I } by utilizing a step-by-step search technique.Initially, the values of α 1 , α 2 , . . ., α I are all set to 0, and the step length is assumed to be a fixed value falling into [0,1].For a better understanding, we take variable α 1 ∈ α as an example to illustrate the manner to decide its optimal value.Suppose that the step length is 0.1, then variable α 1 can pick one value from [0, 0.1, 0.2, . . ., 1.0].For each possible value of variable α 1 , one round of Monte Carlo simulation is then performed.The first value of variable α 1 at which the function schedule satisfies the deadline missing ratio constraint is considered to Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.be the one that captures the stochastic attribute of application executions.Likewise, the values of remaining variables α 2 , α 3 , . . ., α I can also be obtained by conducting a step-bystep search.Practically, this step-by-step search program could be conducted in parallel since it allows searching for the optimal values of these random variables from the ranges of [0, 0.5] and [0.5, 1] simultaneously or even recursively.Therefore, our search method could significantly accelerate the entire search process for modern scheduler platforms equipped with multiprocessors or multiple graphics processing units, as illustrated in Fig. 5.

VIII. EVALUATION
A. Experimental Settings 1) Hardware Platform: We randomly choose 600 out of 3233 base stations Shanghai Telecom dataset [26], [40] to construct the target network.For each of the 600 base stations, it is supposed that an edge server randomly selected from Table II has already been placed at the same location.The average fault rate corresponding to the maximum computation capacity of an edge server falls into [10 −8 , 10 −6 ], and the sensitivity of fault rates to computation capacity scaling of a container is selected from [1,10] [41].The harvesting power of an edge server is generated by [42] where N (t) obeys a standardized Gaussian distribution and γ ∈ [1,30] is a random variable.The bandwidth between any two base stations falls into the range of [500, 5000]KB/s, and it is equally shared by all virtual links over the two base stations.The average fault rate of a virtual link falls into [10 −7 , 10 −5 ].
2) Application Information: We sample a total of 2119 DAG applications with varied depth and breadth sizes from Alibaba cluster trace [32], [43].For each function in a DAG, its minimum and maximal number of execution cycles are scaled to the intervals of [4 × 10 3 , 5 × 10 5 ] and [6 × 10 6 , 7 × 10 8 ], respectively.Accordingly, the minimum and maximal amount of communication data from one function to its direct successor could be scaled to the interval of [500, 1500] megabytes.For the Gaussian distribution of an application, the mean and variance are simulated as follows.

Mean(
For the exponential distribution of an application, its rate parameter is also generated by Mean(A i ).The reliability thresholds of individual applications are randomly picked from the interval of [0.7, 0.9999].The desired finish time and the deadline of an application are set as follows.
T * ,i latency = In ( 44)-( 45), 1 , 2 , 3 , and 4 are coefficients related to the number of functions for execution.We set 1 = 2 = 0.3 × I i=1 ξ i and 3 = 4 = 0.7 × I i=1 ξ i .The 2000 GHz and 2000 KB/s are the estimated computation and commutation capacities on average of the target network, respectively.The upper boundary on deadline missing ratio of applications is set to 10%.We recruit 120 participants to build the personality-driven user QoE prediction model.

B. Benchmarking Algorithms
We compare our PSA-NSGA-II method with representative multiobjective evolutionary algorithms, as described below.
r B-NSGA-II [44] is a basic NSGA-II version with conven- tional crossover and mutation operations.In B-NSGA-II, the crossover and mutation probabilities are specified as constants, i.e., they are not allowed to be altered during the entire evolutionary procedure.r c-DPEA [45] is a newfangled evolutionary algorithm that adopts dual populations to balance the solution diversity and algorithm convergence.In c-DPEA, two complementary populations are generated and simultaneously participated in the evolutionary procedure.
r DLS-MOEA [46] is a multiobjective evolutionary algo- rithm with dual local search.As its name indicates, one significant improvement is to perform a dual local search in both objective and decision spaces.Upon the novel dual local search mechanism, an archive based non-dominated solution generator is developed to enforce the population evolution toward the Pareto frontier.r MOALO [47] is a multiobjective ant lion optimizer that incorporates the popular roulette wheel approach to guide promising ants (i.e., solutions) towards desirable searching regions of multiobjective optimization problems.
r MOEA/D [48] is a decomposition based multiobjective evolutionary algorithm.In MOEA/D, an initial problem is wisely partitioned into several scalar optimization subproblems that could be simultaneously addressed for alleviating computational complexity.
r MaOEA/IGD [49] is a competitive indicator based mul- tiobjective evolutionary algorithm.It adopts the inverted generational distance (IGD) as an indicator to pick out propitious solutions to maintain better algorithm convergence and population diversity.Besides, we compare our stochastic scheduling policy with the following profit-aware function scheduling algorithms.
r AUCT [5] is an auction-based approach allows service providers to make informed decisions regarding function offloading in order to maximize revenue.
r BAYE [6] is a pricing optimization algorithm based on the Bayesian Gaussian process bandit techniques to learn an optimal service pricing policy.
r CostLess [7] leverages the constrained shortest path tech- nique to find a profit maximization solution under specific application latency.
r PAWS [8] adopts the real-world pricing scheme of AWS for latency-sensitive applications over serverless resources in a hybrid edge and cloud scenario.The executable programs or pseudocodes, and appropriate parameters of the above benchmarking algorithms are taken from their original studies.We implement them on a machine configured with 64 GB memory and Intel Xeon W-10855 M 6-core processor.

C. Results and Analysis
We below present the results and analysis of our function scheduling algorithms under the Gaussian distribution of stochastic application executions.Due to page limit, the evaluations under the exponential distribution of stochastic application executions are put into the supplementary material (available at https://www.jianguoyun.com/p/DU3g1qUQ3IvkCRin4awFIAA).
1) Evaluation for Our User QoE Prediction Method: Fig. 6 plots six confusion matrixes of our user QoE prediction method under a varied number of training participants.The diagonal elements in each confusion matrix imply the QoE prediction accuracy while other elements indicate the prediction errors.As demonstrated, the prediction accuracy markedly increases with the number of participants involved in the training process.When 90% of the total participants join in the training process, our QoE prediction method has 89.8% possibility on average to generate a correct prediction result for remaining participants.
2) Evaluation for Deterministic Function Scheduling: In this set of experiments, we compare our PSA-NSGA-II based deterministic algorithm with six benchmarking multiobjective algorithms.For thorough investigations, we first randomly choose 600 out of 3233 base stations as one test network and then repeat this operation to construct 20 varied network structures.Table IV lists the hypervolume of different algorithms, where each data point is the averaged result of 100 runs of individual algorithms.A larger hypervolume is preferred here because it indicates the corresponding algorithm yields non-dominated solutions closer to the Pareto frontier.As observed, the hypervolume of our augmented NSGA-II method is 8.3%, 8.2%, 15.0%, 14.8%, 7.3%, and 10.4% on average larger than that of benchmarking algorithms B-NSGA-II, c-DPEA, DLS-MOEA, MOALO, MOEA/D, and MaOEA/IGD, respectively.III gives a comparison between our deterministic and stochastic scheduling policies.In this group of comparison experiments, we adopt three fixed settings of execution adaption variables α = {α 1 , α 2 , . . ., α I } for our deterministic function scheduling algorithm, i.e., α = 0, α = 0.5, and α = 1.In contrast, our stochastic function scheduling scheme permits every element α i ∈ α to pick any value between 0 and 1.We select the median data point of the two sorted non-dominated solution sets for comparison.As observed, our stochastic scheme outperforms the deterministic algorithm in terms of both service profit and holistic user QoE regardless of network structure variations.Fig. 7 illustrates the convergence of our stochastic parallel function scheduling scheme across a variety of network structures.We adopt the consolidation ratio, a common metric introduced in [50], to quantify the convergence towards the Pareto optimal front.This ratio is defined as the percentage of old solutions that remain non-dominated in the current generation.It becomes evident that the convergence of our scheduling scheme is signified when the consolidation ratio increases to 1, indicative of the absence of new non-dominated solutions.We use colored circles to mark the instances where the consolidation ratios for specific network structures reach 1 for the first time.As can be seen from the figure, our scheduling scheme consistently identifies these colored circles across various network structures, although the precise positions of these markers vary with network structures.
We also evaluate the profits of service providers achieved by our stochastic function scheduling scheme as well as profitaware scheduling algorithms AUCT, BAYE, CostLess and PAWS.Note that the latter four algorithms are originally designed for single-objective optimization problems and therefore may not be directly comparable to our biobjective optimization technique.Therefore, we compare the maximum gain of service providers that can be achieved by the five algorithms under a given holistic QoE.As seen from Table V, our stochastic function scheduling scheme outperforms the profit-aware scheduling algorithms AUCT, BAYE, CostLess and PAWS.For example, when the threshold on normalized holistic QoE is set to 0.7 (0.8), Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

TABLE V MAXIMUM PROFIT (DOLLAR) OF SERVICE PROVIDERS UNDER THE CONSTRAINT OF A GIVEN HOLISTIC QOE
our stochastic function scheduling scheme can boost the profit by up to 61.4% (75.6%).
To evaluate the parallelization capacity of our stochastic function scheduling scheme, we set different parallelization levels by allocating a varied number of cores from 1 to 6 for algorithm execution on the Intel Xeon W-10855 M 6-core processor platform.The level #1 is allocated to a single core, which indicates no parallelization is applied to the acceleration Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply. of our stochastic function scheduling scheme.Oppositely, level #6 is allocated to all six cores, which suggests the highest parallelization level is applied to the acceleration of our stochastic function scheduling scheme.Other parallelization levels from #2 to #4 are allocated to 2, 3, and 4 cores, respectively.Accordingly, Fig. 8 plots the running time of assorted parallelization levels.We see that our stochastic parallel function scheduling scheme exhibits an appealing advantage in reducing time overheads.As expected, as the parallelization level increases, the resultant time overhead drops significantly.For example, the time overhead of parallelization level #6 is up to six times shorter than that of the no parallelization policy.

IX. CONCLUSION
In this paper, a reliability-ensured personalized function scheduling technique is designed for stochastic IoT applications in sustainable serverless edge computing.Our technique contains a personality-driven user QoE prediction method, an augmen[?tedNSGA-II based deterministic function scheduling policy, and a stochastic parallel function scheduling strategy.Experimental results show that our technique can well balance the service profit of the target network and the holistic QoE of end users while satisfying all design constraints.
We acknowledge that providing theoretical guarantees on optimization performance and convergence speed is undeniably a pivotal facet of algorithmic research.However, achieving such guarantees could be significantly challenging for complex and multifaceted optimization problems, such as the stochastic parallel function scheduling in this study.The pursuit of theoretical guarantees will be a prime objective in our subsequent investigations.

Fig. 6 .
Fig. 6.Confusion matrixes of our QoE prediction method under a varied number of training participants.

TABLE II PARAMETERS
OF CONSIDERED EDGE SERVERS (31))measures the degree of only violating energy constraints, where E m demand (O u ) is the energy demand of edge server S m derived by chromosome O u .(32) measures the degree of merely violating reliability constraints, where R i,p softbit (O u ) denotes the reliability of function f i,p obtained by chromosome O u .Putting the two equations together, the overall degree of violating both energy and reliability constraints for chromosome O u is hence calculated as

TABLE III COMPARISON
BETWEEN OUR DETERMINISTIC AND STOCHASTIC FUNCTION SCHEDULING METHODSTABLE IV HYPERVOLUME OF DETERMINISTIC SCHEDULING ALGORITHMS AND THE IMPROVEMENT OF OUR PSA-NSGA-II METHOD3) Evaluation for Stochastic Function Scheduling: Table