Jean-François Duhé

Pierre Melchior

Youssef Abdelmoumen

François Roubertie

Jean-François Duhé

Stéphane Victor
email: stephane.victor@ims-bordeaux.fr

Youssef Abdelmounen

J.-F Duhé

J.-F Duhé

Modeling Thermal Systems With Fractional Models: Human Bronchus Application

Keywords: Heat equation, two-port network, thermal impedance, thermal systems, biological systems, constant phase element (CPE), fractional calculus, fractional systems, lung thermal modeling Heat equation, two-port network, thermal impedance, thermal systems, biological systems, constant phase element (CPE), fractional calculus, fractional systems, lung thermal modeling ---

Thermal modeling of systems allows heat and temperature simulations for many applications, such as refrigeration design, heat dissipation in power electronics, melting processes and bio-heat transfers. Sufficiently accurate models are especially needed in openheart surgery where lung thermal modeling will prevent pulmonary cell dying. For simplicity purposes, simple RC circuits are often used but such models are too simple and lack of precision in dynamical terms. A more complete description of conductive heat transfer can be obtained from the heat equation by means of a twoport network. The analytical expressions obtained from such circuit models are complex and nonlinear in the frequency ω. This complexity in Laplace domain is difficult to handle when it comes to control applications and more specifically during surgery, as heat transfer and temperature control of a tissue may help in reducing necrosis and preserving a greater amount of a given organ. Therefore, a frequency domain analysis of the series and shunt impedances will be presented and different techniques of approximations will be explored in order to obtain simple but sufficiently precise linear fractional transfer function models. Several approximations are proposed to model heat transfers of a human

Introduction

Thermal modeling of systems is of particular interest in applications where temperature might be critical. This includes air conditioning, industrial refrigeration, electronic device cooling, biological tissue heat losses etc. The complexity of heat transfer usually implies the use of finite element methods to solve the heat equation in a chosen region of space. However, finite element models tend to be complex to calculate and in many cases a simpler model that only takes into account temperatures around particular points of interest are precise enough for specific applications.

A typical simple model to take into account heat dissipation and its dynamics is the RC circuit. It is usually used in the domain of power electronics [START_REF] Künzi | Thermal design of power electronic circuits[END_REF], building simulation [START_REF] Parnis | Building thermal modelling using electric circuit simulation[END_REF][START_REF] Danza | A simplified thermal model to control the energy fluxes and to improve the performance of buildings[END_REF] and even to model human heat losses [START_REF] Jiang | A circuit simulating method for heat transfer mechanism in human body[END_REF][START_REF] Ismail | Electric circuit analogy of heat losses of clothed walking human body in windy environment[END_REF]. RC circuit models are also present in other similar applications, such as the measurement of bio-impedances [START_REF] Chinen | New equivalent-electrical circuit model and a practical measurement method for human body impedance[END_REF] or lithium-ion battery models [START_REF] Zhang | Comparative research on RC equivalent circuit models for lithiumion batteries of electric vehicles[END_REF]. In most of these applications, the system input is assumed a low frequency signal. In order to widen the frequency range and especially in high frequency, the thermal two-port network was introduced in [START_REF] Maillet | Thermal Quadrupoles: Solving the Heat Equation through Integral Transforms[END_REF] and such a model is directly derived from the heat equation.

The thermal two-port network models heat conduction in a single direction as a T circuit (see figure 2) with two series impedances Z 1 (s) and Z 2 (s) and a shunt impedance Z 3 (s). These impedance expressions are complex and nonlinear in ω which do not allow obtaining rational transfer functions, the latter being more suited models for control design. It can be shown in high frequency that the thermal impedance of a plane wall is given by a half-order integrator [START_REF] Battaglia | Solving an inverse heat conduction problem using a noninteger identified model[END_REF][START_REF] Malti | Thermal modeling and identification of an aluminium rod using fractional calculus[END_REF]:

lim ω→∞ Z thermal (jω) = lim ω→∞ T (jω) φ(jω) = K √ jω . (1)
Such an expression shows that the heat transfer in high frequency may be modelled by using constantphase elements (CPE) (see [START_REF] Krishna | Studies on fractional order differentiators and integrators: A survey[END_REF][START_REF] Nakagawa | Basic characteristics of a fractance device[END_REF][START_REF] Yakuphanoglu | Active and passive realization of fractance device of order 1/2[END_REF]). Constant-phase elements have particularly been used in biological and medical applications such as modeling intestine tissue [START_REF] Elwakil | Fractional-order circuits and systems: An emerging interdisciplinary research area[END_REF], porous films [START_REF] Das | Characterization of a fractional order element realized by dipping a capacitive type probe in polarizable medium[END_REF], cardiac tissue [START_REF] Magin | Modeling the cardiac tissue electrode interface using fractional calculus[END_REF] or even lung mechanics [START_REF] Copot | Fractional calculus based methods and models to characterize diffusion in the human body[END_REF], thus showing the growing importance of fractional calculus in modeling applications.

During surgery, heat transfer and temperature control of a tissue may help in reducing necrosis and preserving a greater amount of a given organ. Patients can suffer from lung injuries during cardiopulmonary bypass, such as ischemia [START_REF] Buggeskov | Lung protection strategies during cardiopulmonary bypass affect the composition of blood electrolytes and metabolites-a randomized controlled trial[END_REF]. Perfusion techniques involving some type of indirect temperature control exist [START_REF] Maltesen | Lung protection strategies during cardiopulmonary bypass affect the composition of bronchoalveolar fluid and lung tissue in cardiac surgery patients[END_REF]. Some lung injuries may be treated by applying mild hypothermia [START_REF] Wu | Effect of mild hypothermia on lung injury after cardiac arrest in swine based on lung ultrasound[END_REF]. The precision required to keep a fairly precise temperature in the tissue justifies our interest in enhancing the knowledge regarding the dynamics of thermal models for physiological scenarios.

The contributions of this paper lies on at giving different approximation models for both the series and shunt impedances. The series impedances will be approximated by means of an asymptotic model (taking into account only low and high-frequency behavior), pole-zero cells and the use of a fractional Butterworth impedance model. On the other hand, the shunt impedance will be approximated by a capacitance (its most basic model) and by filtered capacitance models: fractional slope filter and multiple pole filters. The poles for the last filter will be analyzed with and without a recursive factor imposed between the parameters. All these model approximations for the series and shunt impedances will be tested on a biological system: the heat transfer in a human bronchus.

This paper is organized as follows. The presentation of the thermal two-port network and its frequencydomain analysis will be presented in section 2. The approximation propositions will be explored in section 3 by providing simulation results. A biological application by modeling the human bronchus will be presented in section 4. Finally, section 5 presents conclusions, final remarks and research perspectives.

2 The analytical thermal two-port network

φ in φ out T in (x = 0) T out (x = L) k, ρ, c x S w Fig. 1 1D thermal system
Let us consider the heat conduction on a simple plane wall in its longitudinal direction x, as shown in figure 1where k, ρ and c respectively stand for the medium thermal conductivity, density and heat capacity. The total length is L and the cross section is S w . The heat equation is expressed as:

ρc ∂T ∂t = k∇ 2 T, (2)
and as conduction only goes in x direction, it can be written as:

ρc ∂T ∂t = k ∂ 2 T ∂x 2 . ( 3
)
The Laplace transform of equation ( 3), with null initial conditions, leads to:

sT (x, s) = k ρc ∂ 2 T (x, s) ∂x 2 , (4)
where s denotes the Laplace variable. By considering a heat flux input at x = 0 and a heat flux output at x = L, it then comes:

Qin (s) = -kS w ∂T (x,s) ∂x x=0 Qout (s) = -kS w ∂T (x,s) ∂x x=L , (5)
or even, put under a matrix form:

T in (s) Qin (s) = M T out (s) Qout (s) (6)
where

M = cosh(δL)
1 kSwδ sinh(δL) kS w sinh(δL) cosh(δL) [START_REF] Cretikos | Respiratory rate: the neglected vital sign[END_REF] with δ = s a and a = ρc k .

T in Z 1 (s) Qin Z 3 (s) Z 2 (s) Qout T out
Fig. 2 Thermal two-port network

Note that δ involves the presence of a fractional (non-integer) operator, as its time-domain representation is a half-order derivative. This matrix equation ( 6) can be represented by a T circuit model, such as represented in figure 2.

The series impedances can then be expressed by:

Z 1 (s) = Z 2 (s) = 1 kS w δ [coth(δL) -csch(δL)] (8)
and the shunt impedance by:

Z 3 (s) = 1 kS w δ csch(δL) (9)
where coth and csch respectively are the hyperbolic cotangent and cosecant functions.

Low-frequency behavior

Processes in thermal applications are usually assumed to be quasi-static, which means that the impedance behavior is only considered in low frequency, thus leading to:

lim ω→0 Z 1 (jω) = L 2kS w = R (10) lim ω→0 Z 3 (jω) = 1 ρcLS w jω = 1 jωC t ( 11
)
where

C t = ρcLS w .
As it can be seen, low-frequency heat transfers lead to a network consisting only of thermal resistances for the series impedances (which represent energy dissipation) and a capacitance for the shunt impedance (which models energy storage in the medium). This is coherent with the classic RC circuit model used to characterize heat transfer dynamics. However, it should be noted that, depending on the thermal properties of the medium, it is possible to have situations in which the transfer occurs at higher frequencies, which limits the validity of the RC model.

High-frequency behavior

Now, by considering the limits in high frequency, the impedance analysis of the thermal quadrupole, by decoupling them into gain and argument, gives:

lim ω→∞ |Z 1 (jω)| = 0, lim ω→∞ arg |Z 1 (jω)| = -45 • (12)
and

lim ω→∞ |Z 3 (jω)| = 0, lim ω→∞ arg |Z 3 (jω)| = -∞. (13)
The impedance Z 1 clearly exhibits a low-pass filter in high frequency, as shown by its zero gain in high frequency. This suggests a capacitive behavior, but the argument has to be analyzed with more care. As can be observed, the high frequency argument of this impedance is not -90 • , but its half. This allows us to deduce that the high frequency behavior of the series impedance Z 1 is a constant phase element of phase -45 • or in other words:

Z 1-HF (s) = 1 C s s 0.5 . (14)
At first, the shunt impedance Z 3 seems coherent with the capacitance model in low frequency, as it is normal for a capacitance to exhibit zero gain in high frequency. The argument, on the other hand, suggests that there is an additional filtering that appears in high frequency, which explains why this argument is not the expected -90 • . The high frequency model of Z 3 can be a combination of a capacitance and a filter, such as:

Z 3-HF (s) = 1 C t s H f ilter (s). ( 15
)
Remark : Adding a correction filter will provide a wider validity domain for the approximation of Z 3 , but it is actually impossible to design a filter that will allow an infinite frequency-range: theoretically, such a filter would have an infinite number of poles in order to obtain the argument shown in equation [START_REF] Ismail | Electric circuit analogy of heat losses of clothed walking human body in windy environment[END_REF].

Impedance approximations

It is impossible to well approximate the quadrupole by a model valid for any frequency just by using transfer function models, unless having an infinite number of parameters. However, a truncated transfer function model may be a fairly useful approximation with a wider validity compared to RC circuits. Consequently, a model going from low frequency (or even a static case) to a high but finite upper frequency will be developed. An academic example will be first used to illustrate the main ideas of section 3 and a more complex example

a 1 m 2 • s -1 k 1 W • m -1 • K -1 L 1 m S w
1 m 2 will be described in section 4. The academic model parameters are presented in table 8.

The proposed approximations will require the determination of different parameters in order to be as close as possible to the true impedance behavior. A common quadratic error criterion will be used by taking the error of the gain in dB defined as:

J(θ) = N i=1 λ 2(N -i) [|Z(jω i )| dB -|Z app (jω i , θ)| dB ] 2 ∆ i log ω N -log ω 1 , (16)
where ω = [ω 1 , . . . , ω N ] is a frequency vector where ω i ∈ [0.1, 100 rad/s] for i = 1, . . . , N and λ is a weighting coefficient (0 < λ < 1) that provides a weight for low or high frequencies. If λ = 1, all frequencies have the same weight. Unless indicated, it will be taken as unity. This criterion ( 16) measures the error in the gain curve given in dB.

The series impedance approximations will be optimized by using the Flower Pollination Algorithm (FPA) (see [START_REF] Yang | Flower pollination algorithm for global optimization[END_REF] for the complete description), which is a novel meta-heuristic optimization method. The parameters to be tuned in order to apply FPA are pretty low, which makes it faster to tune and test when compared to other algorithms. The main parameters to be tuned are the population size n (normally between 10-25) and a switching probability p which states the possibility of doing one of the two principal operations of the algorithm in an iteration: global pollination or local pollination. The most critical parameter is the switching probability, but applications found in literature tend to use a value of 0.8 [START_REF] Yang | Flower pollination algorithm for global optimization[END_REF]. A mathematical study [START_REF] Lukasik | Study of flower pollination algorithm for continuous optimization[END_REF] stated that, even though the value of p does not degrade the results when it comes to unimodal functions, a range between 0.5 and 0.8 is preferred. For this study, all optimizations were performed by using n = 10 and p = 0.8. The main steps for classic FPA are shown in Algorithm 1.

Remark : The choice of FPA as an optimization tool was mainly due to its simple structure and low parameters to be tuned. However, it should be noted that any other optimization algorithm could also be used to approximate the parameters, such as genetic algorithm or particle swarm optimization.

Algorithm 1: Flower Pollination Algorithm

Objective function f (x), x = (x 1 , x 2 , ..., x d ); Initialize a population of n flowers with random solutions Find the best solution g best in the initial population Define a switch probability p while t < M axIterations do for i = 1 : n do if rand < p then Calculate L (Levy flight distribution);

x t+1 i = x t i + L(x t i -g best ); else x t+1 i = x t i + ǫ(x t j -x t k ); end end end
3.1 Approximation of the series impedances Z 1 and Z 2

Asymptotic approximation

As it was seen in the previous section, the series impedance Z 1 and Z 2 are equal and behave like a resistance for low frequencies and like a constant phase element at high frequencies. For simplicity sake, only Z 1 is used. A simple approximation of Z 1 can be defined as the parallel combination of a resistance and a fractance, such as:

Z 1-asymp (s) = R 1 + RC s √ s , (17)
with

R = L 2kSw = 0.5 W • K -1 and C s = kSw √ a = 1 J • K -1 .
This first proposition will be called "asymptotic approximation", as it is obtained by only considering both extremes of the frequency domain. The cut-off frequency obtained at -3 dB under the static gain corresponds to (see appendix C):

ω b ≈ 0.2679 (RC s ) 2 . ( 18
)
Figures 3 and4 respectively show the gain and phase for the approximation Z 1-asymp compared with the exact expression of Z 1 over a frequency range [0.1-100 rad/s]. As expected, both curves get closer towards the limits in low and high frequencies, however, high gain and phase errors occurs in the middle range frequencies. Another important error by using the approximation Z 1-asymp , is the cut-off frequency. The true -3 dB frequency under the static gain is almost one decade greater than the one obtained with the asymptotic approximation, which means this approximation leads to a slower dynamic.

Pole-zero approximation

In order to improve the approximation Z 1-asymp , a correction term may be added to the asymptotic approximation. This correction term will adjust the impedance for a mid-band frequency range, as the asymptotic term can properly handle both low and high-frequency limits. This means that the additional term should not contribute to the gain or to the phase in the low or the high-frequency range. A first proposition for the corrected impedance is to use pole-zero cells that will be placed inside the frequency range of interest:

Z 1-pz (s) = R 1 + RC √ s N cells i 1 + s z i 1 + s p i . ( 19
)
The parameter vector θ is given by:

θ = [p z] (20)
where p and z are vectors containing all the poles and zeros, respectively.

The number of cells N cells to be used for the correction will depend on each specific case. Each cell adds a single pole and zero to the transfer function Z 1 . A large number of cells may give pretty accurate approximations, but the thermal model may become too complex and identifiability problems may occur (large number of parameters in the transfer functions).

In order to determine an optimal number of cells to add for the pole-zero approximation, the error criteria J(θ) was optimized for N cells going from 1 to 8 in order to estimate an optimal number of poles and zeros to be added. Figure 5 shows the different values for the error criteria. As can be seen, no further improvement in the approximation is obtained for this case beyond 3 polezero cells. Even though this technique might succeed in recreating the exact frequency response, the dimension of the parameter vector of the correction term should be taken into account. In this case, dim(θ) = 6. If this thermal impedance is only a section of a global thermal model, it is possible that an approximation of this type will lead to extremely complex expressions for global transfer functions in thermal systems.

Fractional Butterworth approximation

By observing the gain plot in figure 3, another approximation can be made for the thermal impedance Z 1 . As it can be seen in the blue curve, the frequency response is considerably flat before the cut-off frequency. This suggests that this impedance behaves like a Butterworth filter [START_REF] Butterworth | Theory of filter amplifier[END_REF]. However, the high-frequency gain slope is not an integer multiple of -20 dB/dec, but -10 dB/dec. Therefore, inspired by [START_REF] Ali | Fractional order butterworth filter: Active and passive realizations[END_REF], a fractional But-terworth filter can be used as an alternative approximation model:

Z 1-BW (s) = d s α+β + as α + bs β + c . (21)
This new type of model is not a correction term multiplied by the asymptotic approximation, as it can include both the low and high frequency behaviors within its expression:

lim ω→0 Z 1-BW (jω) ≈ d c = R (22)
and

lim ω→∞ Z 1-BW (jω) ≈ d (jω) α+β . ( 23
)
After identification, it is then possible to deduce that:

α + β = 0.5 (24) c = 1 RC (25)
d = 1 C . (26)
For this type of approximation, the parameter vector will always be limited to 3 parameters:

θ = [α a b]. (27)
The fractional Butterworth filter may be a simpler solution with respect to the parameter vector dimension. However, the stability of a system with such a transfer function should be discussed a priori. The BIBO stability for fractional commensurate transfer functions was established by Matignon [START_REF] Matignon | Stability properties for generalized fractional differential systems[END_REF]. An extended criteria for incommensurate systems was later developed [START_REF] Rivero | Stability of fractional order systems[END_REF]. In [START_REF] Ali | Fractional order butterworth filter: Active and passive realizations[END_REF], it is stated that parameters a and b need to be equal or inferior to zero in order to guarantee the fractional stability of the system. The parameters obtained by using the fractional Butterworth filter are given in table 3.1.3. It is interesting to observe that parameter a is zero, which simplifies even more the transfer function of the impedance. For this particular case, the thermal impedance obtained is of the form:

Z(s) = d s 0.5 + bs β + c (28)
or even, its thermal admittance is:

Y (s) = s 0.5 d + bs β d + c d . (29)
From this last expression, impedance Z 1 (s) may actually be reinterpreted as being the parallel combination of a half-capacitance, a resistance and fractionalorder new element which can be expressed as:

Z new (s) = G 0 s α with G 0 < 0. ( 30
)
Further studies are necessary in order to analyze and give a possible physical meaning to this capacitance-like expression, particularly its negative sign.

Comparison of the approximations of the series impedance Z 1

In order to compare the previous approximations, the absolute error is defined as:

ǫ(jω) = ||Z 1-exact (jω)| dB -|Z 1-approx (jω)| dB | . (31)
Figure 6 illustrates the gains for the exact case and the proposed approximations of Z 1 whereas the error signals are plotted in figure 7. The minimum, maximum and mean values of the error obtained for the three approximations are provided in table 3.1.4. Moreover, the mean error for the asymptotic approximation is almost 3 dB and the lowest error is around 1 dB, which translates an inaccurate approximation for all mid-range frequencies. When pole-zero cells are added, the mean error is reduced by a factor of 35, and the error signal can get really close to zero, as shown by its minimum error value (which is not exactly zero, but around 10 -4 ). The Butterworth approximation offers intermediate results, as there is a rough factor of 4 between the asymptotic and the Butterworth mean error. Its peak is still acceptable, but 4 times higher than the one obtained with the polezero approximation. If one needs an extremely precise result, the complexity of the pole-zero cells might be adequate but the simpler Butterworth structure may offer accurate enough results for many applications without adding too much parameters.

3.2 Approximation of the shunt impedance Z 3

Capacitance approximation

The shunt impedance Z 3 behaves like a pure capacitance for low frequencies, but its high frequency behavior is more complex than that of the series impedance, namely:

Z 3-cap (s) = 1 C t s , (32)
with C t = akLS w = 1.

The gain frequency response of Z 3 and its capacitance approximation are shown in figure 8. The -20 dB/dec slope, seen in low frequencies, confirms the capacitance behavior of Z 3 . However, the slope of the gain curve increases as the frequencies increase. Even though a capacitance tends to behave like a shortcircuit on high frequencies, the slope of the gain is not increased and the error curve increases indefinitely in high frequency (see figure 9).

Fractional slope approximation

If the true thermal system is excited in the mid-range frequencies, it is suggested to modify the classic capacitance model to include a slope increment.

In this academic example, and by looking at the low and high frequency slopes present in figure 8, one gets:

The difference between the initial and the final slopes is of almost -60 dB/dec. In order to correct this error, a simple fractional filter is proposed as a multiplicative term to the low-frequency capacitance:

Z 3-f rac (s) = 1 C t s 1 1 + (τ s) φ . ( 35
)
This fractional correction term increases the high frequency slope magnitude by the factor φ, namely 20φ dB/dec. 1/τ indicates the breaking point from where the added correction term begins to influence. An initial estimate for φ would be 2.90 in order to be close to the final slope of -78.28 dB/dec around the higher frequency. The parameter vector for this case is:

θ f rac-slope = [τ φ]. (36)
The correction filter term in equation (35) has the following form:

G(s) = 1 1 + (τ s) φ = 1 s α + b , (37)
a filter that reminds a commensurate transfer function with order α and Matignon's stability theorem may be applied [START_REF] Matignon | Stability properties for generalized fractional differential systems[END_REF], namely:

|arg(λ i )| > α π 2 , (38)
where, for λ = s α , λ i are the roots of the characteristic polynomial in λ.

For the general structure presented in equation (37), we have:

λ i = -b ⇒ | arg(-b)| = π. (39)
Therefore, for this type of filter, the stability condition reads α < 2 for b > 0.

Multiple fractional slope approximation

According to Matignon stability theorem, the fractional order α being limited to 2, it may not be possible to go beyond this limit and to get a more suited highfrequency slope. However, a multi fractional-slope filter may also be proposed:

Z 3-mult-f rac (s) = 1 C t s N i=1 1 1 + (τ i s) ν . ( 40
)
The number of cells to be added N depends on the required additional slope. As each cell may add up to -40 dB/dec, in this case N = 2.

Reduced pole approximation

Looking at the gain curve and its slope in high frequency, an alternative correction term could be the addition of poles to take into account the slope increment.

The difference of slopes in the frequency band of interest, namely:

Ñ = |slope init -slope f inal | 20 , (41)
⌊•⌉ denoting the nearest integer function, provides the number of poles to be added. A new approximation of the shunt impedance can be expressed as:

Z 3-rec-poles (s) = 1 C t s Ñ i=1 1 1 + τ i s (42)
with parameter vector:

θ poles = [τ ]. ( 43
)
If the required number of poles is low enough, one could simply use this approximation and take each parameter τ i as an additional variable to approximate in order to get the correction term. However, there may be cases in which the number of additional poles might be too high and its determination will be too complex in terms of direct optimization (local minima convergence problem). In this case, an arbitrary relation may be established between poles as a way to reduce parameters:

τ i+1 = τ iγ 1 . ( 44
)
where all the poles are related to the first pole τ 1 .

The main advantage of this correction is that the parameter vector will always be limited to the following:

θ red = [τ 1 γ]. (45)
The parameter space will always be 2-D in this case. On the other hand, the main drawback is the loss of flexibility with respect to the general pole addition as detailed in paragraph 3.2.3.

Comparison of the approximations of the shunt impedance Z 3

The approximations obtained by using the fractional slope, the multiple fractional poles and the reduced poles are compared to the pure capacitance model.

Note that to help the convergence of the FPA algorithm for a better fit in high frequency, the weighting factor λ was fixed to 0.95 in the criterion [START_REF] Künzi | Thermal design of power electronic circuits[END_REF] for the fractional and multiple fractional slope models. . Figure 10 illustrates the gains for the exact case and the proposed approximations of Z 3 whereas the error signals are plotted in figure 11.

The minimum, maximum and mean values of the error obtained for the four approximations are provided in table 3.2.5. Even though there is an error that may go up to almost 4.16 dB by using the fractional slope, its reduction when compared to the pure capacitance model is evident. For this example, the fractional order was optimized to φ = 1.45, which is a rather surprising result. The optimization did not lead to a saturated constraint at φ = 2 even though a weight was included for higher frequencies. The obtained value of 1/τ = 8.22 rad/s is coherent with the expected the frequency band [1 -10] rad/s in which the error becomes important (see Figure 11).

Concerning the multiple fractional slope model, the optimization have led to better results. There is roughly a factor of 3 between the mean errors of this model and the simple fractional slope one. The commensurate order was found to be ν = 1.13, which means a highfrequency order of 2ν = 2.26. The additional flexibility provided by the multiple cell fractional filter permits to go beyond the limit of 2 such as imposed by Matignon's stability criterion in order to get a better fit.

The reduced pole approximation provides even more accurate results and the dimension of the research space for the optimization stays the same. There is a dramatic difference at the point in which the approximation starts to modify the curve, as the recursive pole gave a first pole located at τ 1 = 0.033 rad/s, which is two decades below the optimal value obtained by the fractional slope. This approximation may allow a more predictive approach to error reduction than the fractional slope and multiple-cell fractional filter.

Reduced pole approximation versus all pole approximation

In view of the previous results, let us consider the latter approximation given in equation ( 46). Its optimized model can be used as an initial starting point to perform a further optimization: instead of using reduced poles (which indeed has the advantage to reduce the number of parameters in the optimization procedure), this constraint defined in equation ( 44)) can be released so that all the parameters can be optimized independently. Such an approximation of the shunt impedance can be expressed as:

Z 3-all-poles (s) = 1 C t s Ñ i=1 1 1 + τ i s (46)
with the more general parameter vector:

θ poles = [τ ], (47)
where all τ i are independent.

The FPA algorithm was performed with the reduced and all poles by using a gradient descent. As expected, both approximations give really good results. Figure 12 illustrates the comparison of the associated errors for the reduced and the all pole appxoimation cases.

Z3-all-poles Z3-red-poles

Frequency (rad/s) Gain (dB) Fig. 12 Absolute errors ǫ(jω) associated with the reduced and all pole approximations of Z 3 It is indeed difficult to distinguish a further improvement for this new case by allowing the poles to be located at any position. Table 3.2.6 shows the pole locations as well as the criterion J for both cases. Even if there is a slight improvement regarding the criterion J, the non-recursive approximation increased the vector parameter size (three times more than in the recursive case) and the criterion improvement proved to be very slight. On the other hand, it can be seen that all pole approximation locations are almost the same as the ones with the reduced pole approximation. This means that, for this case, the optimal pole location for the approximation was really close to the results provided by the reduced case. It may also justify the use of a reduced pole model in more complex cases.

Application example: heat transfer in a lung bronchus

A potential application for these approximations is to model heat transfer in the human body. During surgery, heat transfer and temperature control of a tissue may help in reducing necrosis and preserving a greater amount of a given organ.

Patients can suffer from lung injuries during cardiopulmonary bypass, such as ischemia [START_REF] Buggeskov | Lung protection strategies during cardiopulmonary bypass affect the composition of blood electrolytes and metabolites-a randomized controlled trial[END_REF]. Perfusion techniques involving some type of indirect temperature control exist [START_REF] Maltesen | Lung protection strategies during cardiopulmonary bypass affect the composition of bronchoalveolar fluid and lung tissue in cardiac surgery patients[END_REF]. Some lung injuries may be treated by applying mild hypothermia [START_REF] Wu | Effect of mild hypothermia on lung injury after cardiac arrest in swine based on lung ultrasound[END_REF]. The precision required to keep a fairly precise temperature in the tissue justifies our interest in enhancing the knowledge regarding dynamic of thermal models for physiological scenarios.

For the application, the heat transfer is considered in the right middle bronchus. Under normal conditions, the nasal cavity and larynx play a significant role regarding air heating and humidification (see [START_REF] Elad | Air-conditioning in the human nasal cavity[END_REF] and [START_REF] Dias | Larynx and cervical trachea in humidification and heating of inhaled gases[END_REF]) and procedures impairing this conditioning may change thermal conditions inside the lungs. Human bronchi may be exposed to thermal stress in such conditions. For this analysis, an intermediate length of L = 0.0236 m was chosen as it is a mean value for right middle bronchus length [START_REF] Mi | Measurement and analysis of the tracheobronchial tree in chinese population using computed tomography[END_REF] and radius r ≈ 1 mm (see [START_REF] Mcfawn | Bronchial compliance and wall structure during development of the immature human and pig lung[END_REF]).

In terms of frequency, even though human respiration is not a naturally fast mechanism, it may have nonnegligible dynamics for heat transfers. Human breathing can go up to 20 breaths per minute in normal conditions, this breathing being around 12 breaths per minute for an adult in average. Slow breathing is considered in the range 0.07-0.16 Hz (see [START_REF] Russo | The physiological effects of slow breathing in the healthy human[END_REF]) and above 14 breaths per minute (or 0.23 Hz), such breathing can already be considered for some as an abnormal value (see [START_REF] Cretikos | Respiratory rate: the neglected vital sign[END_REF]). Therefore, the optimization will be carried out in the frequency range [0.01 -5] rad/s in order to go slightly over a normal condition breathing.

For the series impedance Z 1 , the optimal number of cells for the pole-zero correction was found to be 2, as shown in figure 13. The gain diagrams of the different approximations proposed in paragraph 3.1 (asymptotic, pole-zero, and fractional Butterworth approximations) and their associated errors are respectively shown in figures 14 and 15. Table 4 summarizes a comparison between the obtained errors.

As expected, the asymptotic approximation offers a significant error for the whole frequency range, which may indicate that a human bronchus cannot be considered in the low-frequency limit as human breathing cycles are in the mi-band frequency range. Once again, it can be seen that the pole-zero approximation offers the most accurate results. However, for this scenario, the improvement obtained by using the Butterworth is more significant than in the academic case. For this case, the pole-zero approximation has 6 parameters, which includes 3 additional parameters with respect to the Butterworth case. However, the Butterworth impedance still provides an improved estimation than the asymptotic one. It can be seen in Appendix B that for this case, the Butterworth impedance was once again reduced to the form shown in equation [START_REF] Russo | The physiological effects of slow breathing in the healthy human[END_REF], which offers further advantages in terms of parameter space (2D).

The gain diagrams of the different approximations proposed in paragraph 3.2 (capacitance, fractional slope, multiple fractional slope and pole approximations) and their associated errors are respectively shown in figures 16 and 17. Table 4 summarizes a comparison between the obtained errors. The pure capacitance model provides the worst approximation: the mean error is of 7.58 dB in the considered frequency range and it can go up to 42.05 dB as we get closer to the higher frequency. For this particular geometry, three integer poles were required to provide a good approximation. It should be noted that for a frequency of 1 rad/s (close to normal breathing), the slope of the Z 3 gain curve is already -38 dB/dec, which means that a normal breathing is already in its midband frequency range. The best results are obtained with the reduced poles approximation with a mean error of 0.14 dB.

Conclusion and final remarks

Different propositions for approximating thermal impedances were proposed and compared in terms of accuracy in the frequency ban of interest. The series impedance of a plane wall proved to have an interesting behavior in the frequency domain, as its usual simplification as a simple resistance is inaccurate beyond the lowfrequency spectrum. The presence of a half-capacitance indicates a possible additional thermal accumulation that cannot be taken into account by a simple capacitance model. Even if the resistance-half capacitance parallel model succeeds in approximating the frequency limits, there is a significative error in the mid-band frequency that can be reduced by using pole-zero cells or a fractional Butterworth expression. The pole-zero approximation proved to be better in terms of accuracy, but the impedance expression is complex and the number of parameters required for an optimal approximation may be too large. The pole-zero cells are an exclusively mathematical construct and lack physical mean-ing. The fractional Butterworth approximation may not offer such accurate results, but it still greatly improves the error when compared to the resistance-half capacitance parallel model. The parameter vector is limited to a maximum dimension of 3 and it is suggested that this type of expression may actually be reinterpreted in a circuit model as the parallel addition of additional elements to the resistance-half capacitance model. Additional physical meaning may be obtained from this model in future studies.

The shunt impedance has a more complex frequency behavior and only a low-frequency analysis was done. In low-frequency, this impedance behaves like a capacitance. However, the quick slope decreases (both in gain and phase) in high frequency, a circuit element nor a simple fractional model is sufficient. This is the reason why the approximations presented for Z 3 are limited to a specific frequency range. The fractional slope helped in reducing the errors and has little parameters to be determined. However, the approximation with recursive poles proved to be more accurate without adding more parameters to the approximation. For the example analyzed, the accuracy was not significantly improved by allowing the poles to be independent. A simple recursive pole approximation gives accurate enough results and provides a simple expression. However, it should be noted that Z 1 and Z 3 are inseparable elements of a same system. Therefore, equivalent transfer functions obtained from the T circuit will involve fractional-order derivatives as it cannot be avoided for the series impedances Z 1 and Z 2 .

All the proposed approximations for modeling the series and shunt impedance where applied in a simulation example of a human bronchus. As previously stated, the simulation results well illustrates that the pole-zero approximation well fits the series impedance Z 1 and Z 2 . For the shunt impedance Z 3 , it could be fitted by using a single integer-order pole.

Research perspectives of this study include analyzing different type of scenarios with the thermal twoport networks and its transfer function properties when applying the proposed approximation techniques. The possibility of getting transfer function models with physically meaningful parameters may improve the analysis of system identification for thermal systems. The frequency response of the asymptotic approximation for the series impedance: Thermal two-port network Figure 3 Z1 gain diagram and its asymptotic approximation

Z asym (jω) = R 1 + RC s √ jω (48)
Z1 phase diagram and its asymptotic approximation

J criterion with varying Ncells for the series impedance correction

Gain diagrams of Z1 and its approximations

please see the manuscript le for the full caption

Z3 gain diagram and its capacitance approximation please see the manuscript le for the full caption

J criterion versus pole-zero units for a human bronchus

Z1 approximation gains for a human bronchus Figure 15 please see the manuscript le for the full caption

Gain diagrams of Z3 and its approximations for a human bronchus

please see the manuscript le for the full caption

Fig. 3 Z

 3 Fig. 3 Z 1 gain diagram and its asymptotic approximation

Fig. 4 Z

 4 Fig. 4 Z 1 phase diagram and its asymptotic approximation

Fig. 5 J

 5 Fig. 5 J criterion with varying N cells for the series impedance correction

Fig. 6

 6 Fig. 6 Gain diagrams of Z 1 and its approximations

Fig. 7

 7 Fig. 7 Absolute errors ǫ(jω) associated with the different approximations of Z 1

Fig. 8 Z

 8 Fig. 8 Z 3 gain diagram and its capacitance approximation

Fig. 9

 9 Fig. 9 Absolute error ǫ(jω) associated with the capacitance approximation of Z 3

=

  -78.28 dB/dec.

Fig. 10 Fig. 11

 1011 Fig. 10 Gain diagrams of Z 3 and its approximations

Fig. 13 JFig. 14 Z

 1314 Fig.13J criterion versus pole-zero units for a human bronchus

Fig. 15

 15 Fig. 15 Absolute errors ǫ(jω) associated with the different approximations of Z 1 for a human bronchus

Fig. 16

 16 Fig.[START_REF] Künzi | Thermal design of power electronic circuits[END_REF] Gain diagrams of Z 3 and its approximations for a human bronchus

Fig. 17

 17 Fig. 17 Absolute errors ǫ(jω) associated with the different approximations of Z 3 for a human bronchus

Z 1 - 2 Z 3 τ 2 = τ γ 1 , τ 3 = τ 2γ 1 and γ = 0. 674 B

 123211674 BW 1 s 0.5 -0.859s 0.274 + i s with τ 1 = 0.033, Estimated transfer functions for human bronchus C Proof of fractional first-order system bandwidth

τ 2 = τ γ 1 , τ 3 = τ 2γ 1 and

 211 i s with τ 1 = 2.407, γ = -0.438 from which one gets the gain:|Z asym (jω)| = R (RC s ) 2 ω + RC s √ 2ω + 1(49)it can be seen that the highest gain in dB for this system is:|Z asym (jω)| dB-max = 20 log 10 (R)(50)The -3 dB frequency is given by:|Z asym (jω b )| -3dB = |Z asym (jω)| dB-max -3 dB(51)By replacing equation (49) in (51), one gets the following polynomial:(RC s ) 2 ω b + √ 2RC s √ ω b + (1 -10 0.3 ) = 0 (52)By taking the single positive and real-valued solution for ω b :ω b ≈ 0.2679 (RC s ) 2(53)which is the expression for the bandwidth.

Figures

  Figures

Figure 1 1D thermal system Figure 2

 12 Figure 1

Figure 9 please

 9 Figure 9

Table 1

 1 Academic example simulation parameters

	Parameter	Value

Table 2 Z 1

 21

		(s) approximation parameters for the fractional
	Butterworth filter
	a	b	α	J(θ)(dB 2 )
	0	-0.859	0.226	0.790

Table 3

 3 Extreme and mean values of the absolute errors obtained for the Z 1 approximations

	Approximation Max (dB) Mean (dB) Min (dB)
	Asymptotic	4.99	2.84	0.96
	Pole-zero	0.39	0.08	0.00
	Butterworth	1.37	0.78	0.04

Table 4

 4 Extreme and mean values of the absolute errors obtained for the Z 3 approximations

	Approximation	Max (dB) Mean (dB) Min (dB)
	Capacitance	35.40	5.69	0.0004
	Frac slope	4.16	1.32	0.0100
	Multi frac slope	1.76	0.39	0.0137
	Red poles	1.48	0.18	0.0027

Table 5

 5 Pole location and J criterion for Z 3 approximations

		Recursive	Non-recursive
		9.91	9.94
	Poles	30.16	30.34
		98.31	92.16
	J(dB 2 )	0.033	0.0286

Table 6

 6 Extreme and mean values of the absolute errors obtained for the Z 1 approximations of a human bronchus

	Approximation Max (dB) Mean (dB) Min (dB)
	Asymptotic	4.99	3.05	1.0829
	Pole-zero	0.26	0.05	0.0001
	Butterworth	1.34	0.80	0.0151

Table 7

 7 Error signal peak and mean values for Z 3 (s) approximations in a human bronchus

	Approximation	Peak (dB) Mean (dB) Min (dB)
	Capacitance	42.05	7.58	0.0031
	Frac slope	8.66	1.57	0.0494
	Multi frac slope	3.43	0.54	0.0177
	Red poles	1.39	0.14	0.0000

Table 8

 8 Competing interests: The authors have no relevant financial or non-financial interests to disclose. Code availability (software application or custom code): Not applicable Authors' Contributions -in relation with your submission: J.F. Duhé has developed the results with S. Victor. He also made the simulations. S. Victor, P. Melchior, Y. Abdelmoumen and F. Roubertie have contributed to the biological and physiological results. Ethics approval (include appropriate approvals or waivers): Estimated models for the academic example of section 3

	6 Declarations
	Funding: Not applicable
	Conflicts of interest/

Table 9

 9 Estimated models for the bronchus applications of section 4 = 0.010, z 2 = 0.112, z 3 = 0.981 p 1 = 0.014, p 2 = 0.166, p 3 = 0.528

	Z 1-asymp		150240 1 + 2.5158	√	s
	Z 1-zp		150240 1 + 2.5158	√	s	3 1	1 + s z i 1 + s p i
	z 1 Z 1-BW Z 3-cap	59720.2 s 0.5 -0.406s 0.25 + 0.3975 11869 s
	Z 3-f rac		11869 s	1 1 + (3.3873s) 1.36
	Z 3-mul-f rac	11869 s	1 1 + (0.734s) 1.12	1 1 + (2.954s) 1.12


domain: info spi

text:
Melan Vijayaratnam

Marco Cagnazzo

Giuseppe Valenzise

Anthony Trioux

Michel Kieffer

Towards zero-latency video transmission through frame extrapolation

Keywords: Extrapolation, low-latency video delivery, video deep learning

INTRODUCTION

Ultra low-latency video delivery is an essential feature in many applications involving interactions among humans (e.g., video conferencing, virtual and augmented reality) or between humans and machines (e.g., teleoperation of unmanned vehicles or robots, remote surgery, etc.). In these scenarios, the Glass-to-Glass (G2G) latency, intended as the delay between the acquisition of a video frame by an agent and its display by a second (remote) agent [START_REF] Bachhuber | On the Minimization of Glass-to-Glass and Glass-to-Algorithm Delay in Video Communication[END_REF], plays a major role in the overall quality of experience perceived by users [START_REF] Brunnström | Latency impact on quality of experience in a virtual reality simulator for remote control of machines[END_REF].

G2G latency consists of acquisition, coding, buffering at transmitter, transmission, buffering at receiver, decoding, and display delays. In the past decades, significant efforts have been devoted to optimizing each of these individual sources of delay. Nevertheless, the minimum achievable latency is still constrained by technological and physical limits (the most evident being the speed of light), which represent a hard lower bound beyond which latency cannot be further reduced. With these clear limitations in mind, in this paper, we investigate how we can further reduce the G2G latency to be arbitrary low. Given the uncompressible physical delays in video delivery, a possible way to achieve this goal is to predict video frames that have not been received yet using those already available at the receiver.

The use of prediction/extrapolation to compensate latency is not new. It has been already employed in a number of applications, in-cluding virtual reality [START_REF] Garcia-Agundez | An evaluation of extrapolation and filtering techniques in head tracking for virtual environments to reduce cybersickness[END_REF], tactile interfaces [START_REF] Henze | Software-reduced touchscreen latency[END_REF], cloud gaming [START_REF] Lee | Outatime: Using speculation to enable low-latency continuous interaction for mobile cloud gaming[END_REF], etc. Nevertheless, to the best of our knowledge, this simple idea of using extrapolation has not been explored and analyzed yet for the case of latency compensation in video transmission. In this paper, we consider for the first time the application of extrapolation to ultra-low or even zero latency video communication. To this end, we analyze two possible schemes to integrate extrapolation together with a conventional video codec, i.e., extrapolating from the original frames at the encoder side before coding and transmission, or at the decoder side based on the available decoded (quantized) frames. In both cases, if the prediction horizon is large enough, all sources of delay (including the time needed for extrapolation) can be compensated and a predicted frame can be displayed at the receiver while it is acquired at the transmitter.

Since the extrapolated frames are in general different from the true ones, latency compensation produces a loss of fidelity. In this paper, using different recent learning-based extrapolation algorithms, we characterize the trade-off between latency and fidelity. Our aim is to study the feasibility of latency compensation based on frame extrapolation. Our preliminary results on a driving car dataset indicate that the effectiveness of this approach depends significantly on the content and that there is still large space for improvement in video extrapolation techniques.

RELATED WORK

Extrapolation as a means to compensate latency has been applied in many fields. In virtual reality headsets, predicting future head movements provides a low-cost solution to reduce latency, which is one of the main causes of motion sickness [START_REF] Garcia-Agundez | An evaluation of extrapolation and filtering techniques in head tracking for virtual environments to reduce cybersickness[END_REF]. To provide a seamless experience with respect to touchscreen devices, the authors of [START_REF] Henze | Software-reduced touchscreen latency[END_REF] extrapolate the finger's movement to compensate software latency. In cloud gaming, it has been proposed to render speculative frames of future possible outcomes, delivering them to the client ahead of time [START_REF] Lee | Outatime: Using speculation to enable low-latency continuous interaction for mobile cloud gaming[END_REF]. Our work relies on similar ideas, but applies them to a generic video coding architecture.

From a video coding point of view, all recent video coding standards from ITU-T ISO/IEC (e.g., HEVC, SHVC, or VVC) have lowlatency profiles that optimize the predictive coding structure. The MPEG-I VVC [START_REF] Bross | General video coding technology in responses to the joint call for proposals on video compression with capability beyond HEVC[END_REF] has added some features enabling low latency at a signaling level such as the Gradual Decoding Refresh (GDR). GDR redistributes the bitrate of the intra-frames over many frames, smoothing the bandwidth consumption, and thus the transmission latency, over time. To reduce the transmission delay, the 5G ultrareliable low latency communications (URLLC) mode targets radio access delays of less than one millisecond [START_REF] Ali | URLLC for 5G and Beyond: Requirements, Enabling Incumbent Technologies and Network Intelligence[END_REF]. Despite these de- velopments, as pointed out in [START_REF] Sharma | Toward Tactile Internet in Beyond 5G Era: Recent Advances, Current Issues, and Future Directions[END_REF], current video systems are still subject to G2G delays from 50 to 400 ms. In order to satisfy the requirements of latency-sensitive applications, Bachhuber et al. [START_REF] Bachhuber | On the Minimization of Glass-to-Glass and Glass-to-Algorithm Delay in Video Communication[END_REF] analyzed all delay components in a video communication system. Frame skipping and frame preemption mechanisms are proposed to further reduce the G2G latency. Neglecting the buffering and network transmission delays, direct communication with a latency of about 20 ms is achieved.

These works optimize specific sources of delay in the processing pipeline, which intrinsically limits the maximum latency reduction that can be achieved. Instead, in this paper we consider the global G2G latency, with the goal to arbitrarily reduce it.

LATENCY COMPENSATION VIA EXTRAPOLATION

We propose to compensate latency using frame extrapolation. This extrapolation can be applied at either the encoder or at the decoder side, as illustrated in Fig. 1. For the sake of simplicity, we assume here a low-delay encoder configuration, where the encoding order of frames is identical to the display order. Nevertheless, the proposed scheme generalizes to any GOP structure. In the following, we describe the two compensation schemes in detail.

Extrapolation at the decoder side

In this scheme, no modification of the transmitter is required with respect to a standard transmission pipeline: the acquired frames are compressed and transmitted to the receiver, where they arrive with some latency. Figure 1(a) illustrates the extrapolation at the decoder side with an example. In this figure, In represents the n-th frame of the video. We assume in this example that the G2G latency δ is equal to two frames (corresponding to a time of 2/f seconds, where f is the frame rate). Therefore the decoder receives Ĩn-2 (the compressed version of In-2) when the encoder is already acquiring In. In order to compensate this latency, the decoder runs a frame extrapolation algorithm F which takes as input a given number k of available decoded frames, and produces a prediction În of frame n as:

În = F({ Ĩn-h , Ĩn-h-1 , . . . , Ĩn-h-k+1 }; h). (1)
The k input frames are also referred to as context frames. For example, in Figure 1, k = 2 context frames are considered. The parameter h is the temporal horizon of the extrapolation mechanism. It determines the amount of latency we want to compensate for. The larger is h, the more difficult it is to get a reliable prediction, as we show in the experimental results. The parameter h determines the trade-off between latency compensation and quality degradation. In the example we have h = 2, meaning that we want to fully compensate for the latency δ. While the context frames k are dependent of the extrapolators used, the temporal horizon h is chosen depending on the application, there is no relationships between these two variables. With these settings, the extrapolation algorithm produces:

În = F({ Ĩn-2, Ĩn-3}; 2). (2)
Therefore, the estimation În of frame In, extrapolated from compressed frames, can be displayed at the decoder while frame In is acquired by the encoder.

Extrapolation at the encoder side

The extrapolation at the encoder sides is illustrated in Figure 1(b), following up the example used in the previous section. The extrapolation method F(•; h), is used this time with the acquired frames: In, In-1, . . ., I n-k+1 used as context frames. Likewise, the prediction depends on the temporal horizon h, i.e., we compute

Î(n + h)=F({In, In-1, . . . , I n-k }; h). (3)
As for the previous example, we have k=2 and h=2 here. Then, the extrapolated frames are compressed, transmitted, decoded, and displayed. In this case, only the transmitter has to be adapted and a plain, standard receiver is employed. The extrapolation performed at the encoder provides a preemptive transmission latency compensation. In the considered example, the compressed, extrapolated frame Ĩn is displayed while the frame In is acquired by the transmitter.

context k horizon h Copylast 1 1 FlowNet2 + warping 2 1 MCNet 5 5 SDCNet 2 1
Table 1. Frame extrapolation methods and their number of context and default prediction horizon.

With either the decoder-side or encoder-side extrapolation schemes, it is possible to compensate the same amount of latency. This potentially includes also the time to perform extrapolation, which might be not negligible, depending on the adopted method. In a multicast scenario, the extrapolation at the encoder has an advantage compared to the extrapolation at the decoder, as extrapolation is run only once, and more complex techniques with possibly more powerful hardware may be employed. Nevertheless, this also means that the latency compensation is the same for all the receivers. This may be suboptimal when the G2G latency differs among receivers. In such case, the extrapolation at the decoder gives more flexibility. Furthermore, one could also consider a variant of the proposed approaches where extrapolation is performed both at the encoder and the decoder side. In this case, extrapolation at encoder is used to compensate for some minimum, common latency; then, each receiver can further compensate latency according to its needs. Finally, latency compensation comes with a cost in terms of fidelity loss. In the following experiments, we study this latency-distortion trade-off using some off-the-shelf video extrapolators.

EXPERIMENTS AND DISCUSSION

This section evaluates the feasibility of frame extrapolation as an effective latency compensation tool. Specifically, we characterize the loss in quality for different amounts of compensated latency, using different off-the-shelf extrapolators to implement F. All the experiments are conducted using the HEVC HM codec implementation (16.24) [START_REF] Rosewarne | High Efficiency Video Coding (HEVC) Test Model 16 (HM 16) Encoder Description Update 13 | MPEG[END_REF] using the encoder intra main.cfg configuration.

Choice of the extrapolators

As discussed in [START_REF] Gao | Disentangling propagation and generation for video prediction[END_REF], video prediction methods can be motion-based, pixel-based, or fusion-based. Motion-based methods seek to compute the motion of image pixel intensities which corresponds to the motion of objects in a scene. Deep learning methods in this category compute optical flow with neural networks. Combined with a warping operation, next frames are inferred from previous images. Pixel-based methods seek to generate each pixel from scratch, i.e., they do not use explicit motion representations such as an optical flow or motion vectors. Fusion-based methods combine both motion and pixel-based methods. In this work, we select one method from each category. FlowNet-2 [START_REF] Ilg | Flownet 2.0: Evolution of optical flow estimation with deep networks[END_REF] is a motion-based technique combined with a warping operation, which can effectively predict the next frame. MCNet [START_REF] Villegas | Decomposing motion and content for natural video sequence prediction[END_REF] is a pixel-based technique built upon the Encoder-Decoder Convolutional Neural Network and Convolutional LSTM for pixel-level prediction, which independently capture the spatial layout of an image and the corresponding temporal dynamics. Finally, SDCNet [START_REF] Reda | SDCNet: Video prediction using spatially-displaced convolution[END_REF] is a fusion-based approach that models moving appearances with both convolutional kernels and vectors as optical flow. Some characteristics of the considered extrapolation methods are reported in Table 1. Each extrapolator is designed to work with a certain number k of context frames, and a certain prediction horizon h. The latter might be different (smaller) than the desired latency to compensate. In this case, we simply apply the extrapolation recursively, i.e., we progressively add extrapolated frames to the context.

We also include a simple frame-copy extrapolation, dubbed Copylast. This method just copies the last available frame at the encoder/decoder, and is reported as a very low complexity baseline.

Datasets

The considered learning-based extrapolation methods require data to be trained on. In our experiments, we have trained the unsupervised methods SDCNet and MCNet on the Caltech Pedestrian dataset [START_REF] Dollár | Pedestrian detection: A benchmark[END_REF], which contains 128419 images from 65 different video sequences at 30 fps with the resolution 640×480. All the frames are center cropped into 640 × 448. This dataset provides us with both quantity and variety and reflects a real-case scenario. For the supervised method FlowNet2 [START_REF] Ilg | Flownet 2.0: Evolution of optical flow estimation with deep networks[END_REF], we need ground-truth optical flow, which is not available on the previous dataset. Therefore, we use the weights pre-trained on the MPI Sintel dataset [START_REF] Butler | A Naturalistic Open Source Movie for Optical Flow Evaluation[END_REF]. Sintel is an action movie and as such contains many fast movements that are difficult for traditional non learning-based methods. The test experiments (along with rate-distortion analysis shown in next section) are done on the DriveSeg [START_REF] Ding | MIT driveseg (manual) dataset for dynamic driving scene segmentation[END_REF] Manual scene. Similar to the training set, the dataset has been captured from a moving vehicle during continuous daylight driving through a crowded city street. It comprises 5000 frames taken at 30 fps. We center crop the frames of the sequence into 640 × 448.

Study of the distortion-fidelity trade-off

To evaluate the scheme where extrapolation is performed at the decoder, the DriveSeg test set at 30 fps is compressed using HEVC with different quantization parameters QP ∈ {22, 27, 32, 37} for all the selected methods. Rate-Distortion (RD) points are computed considering an extrapolation horizon h ∈ {1, . . . , 5} to compensate a latency δ∈{33 ms, 67 ms, 100 ms, 133 ms, 167 ms}. For the case of extrapolation at the encoder, the same test set is considered, with the same an horizon h ∈ {1, . . . , 5} and the same QPs for compressing extrapolated frames with HEVC. Table 2. Quantitative results on DriveSeg scene: results for the extrapolation at the decoder side and gain/loss (in parenthesis) obtained with the extrapolation scheme at encoder side.

Fig. 3. Qualitative results for a latency compensation of 100 ms. The extrapolation at decoder side scheme is used. The HEVC codec uses a quantization parameter QP = 32. Images taken from the Kitti [START_REF] Geiger | Vision meets robotics: The kitti dataset[END_REF] dataset.

The distortion between displayed and original frames is evaluated using three objective metrics: PSNR in the YCbCr color space [START_REF] Sullivan | Objective quality metric and alternative methods for measuring coding efficiency[END_REF], SSIM [START_REF] Wang | Image quality assessment: from error visibility to structural similarity[END_REF], and VMAF [START_REF] Bampis | Spatiotemporal Feature Integration and Model Fusion for Full Reference Video Quality Assessment[END_REF]. This yields a rate-distortion curve for each possible extrapolation horizon (latency compensation) h. We compare all curves to the case h = 0, i.e., without any extrapolation. The average quality losses over different rates, expressed by the negative Bjøntegaard delta metric (BDPSNR, BDSSIM and BD-VMAF), are reported in Table 2, where the first numbers are the BD metric values for the extrapolation at the decoder side, and the numbers in parentheses are the differences with respect to this scheme when using extrapolation at the encoder. We can observe that extrapolation at the encoder/decoder produces essentially the same quality losses for a given amount of compensated latency.

Extrapolation at encoder is performed on original frames, and may therefore be more efficient. Nevertheless, for the considered values of QP , extrapolation at the decoder is still efficient, even using decoded context frames. In both cases, the fidelity loss introduced by the extrapolation dominates that induced by compression.

Figure 2 illustrates visually the latency-distortion trade-off in terms of BDPSNR. As expected, the larger the amount of latency we compensate, the higher the quality loss we incur. For typical G2G latency values of 100 ms observed in remote control applications [START_REF] Marai | Smooth and Low Latency Video Streaming for Autonomous Cars During Handover[END_REF][START_REF] Sharma | On the Latency in Vehicular Control using Video Streaming over Wi-Fi[END_REF], on the test dataset used in this work, the PSNR loss ranges between 12.3 and 8.3 dB with the considered extrapolation techniques.

Discussion

We have verified quantitatively the expected trade-off between fidelity and latency compensation for a specific dataset. This aspect is key for the feasibility of latency compensation, as distortion is the unavoidable consequence of arbitrarily reducing latency. The effectiveness of latency compensation depends significantly on the interpolation operator F. The proposed scheme employing recent frame prediction techniques entails PSNR drops in excess of 8 dB for 100 ms reduction. This motivates the study of better future frame prediction algorithms in the future.

The fidelity metrics alone might not be sufficient to explain the feasibility of latency compensation through extrapolation. Figure 3 shows the results of compensating 100 ms on a frame of the Kitti dataset [START_REF] Geiger | Vision meets robotics: The kitti dataset[END_REF], which similar to the training set. The interpolation of SDCNet has clear visual deformations compared to the original frame. Nevertheless, we can observe that the position of the bike is approximately aligned with the original frame. On the other hand, Copylast (which only uses the last available frame and corresponds to not compensating any latency) produces a lag between the actual position of the bike and the displayed one. We can imagine that, depending on the application (e.g., teleoperation), the SDCNet prediction brings valuable information (the bike position), even though such quality metrics as PSNR are not able to catch this aspect. This example highlights the essential objective of latency reduction, i.e., to gain more knowledge about the future representation of the dynamic scene by anticipating next images. This question is not only essential at the human perceptual level, but also at the machine level. We can imagine this future understanding as key components for the machine to anticipate behaviors.

CONCLUSION

This paper introduces a tool that allows us to compensate latency in a video transmission scheme at the price of additional complexity but also of degradation of the frame fidelity. Depending on the application, many configurations are possible extrapolation at the encoder, extrapolation at the decoder, or both at the encoder and decoder. The degradation is essentially caused by the extrapolation approaches. The goal of this paper is not to propose a better extrapolation approach to reduce this loss, but rather demonstrating the applicability of the latency compensation relying on such tool. Future works may concern the improvement of existing extrapolating methods for such task.

Fig. 1 .

 1 Fig. 1. An example of latency compensation with (a) extrapolation at the decoder and (b) extrapolation at the encoder. Ĩ and Î indicate decoded (quantized) and predicted (extrapolated) frames, respectively. F is the extrapolation function.

Fig. 2 .

 2 Fig. 2. Evolution of the latency compensation-distortion trade-off in terms of PSNR.

  +0.02) 14.21 (+0.02) 15.61 (+0.01) 0.19 (+0.00) 0.35 (+0.00) 0.43 (+0.00) 46.52 (-0.14) 65.44 (-0.09) 71.65 (-0.07) MCNet 8.34 (-0.15) 12.26 (-0.28) 14.79 (-0.29) 0.06 (+0.00) 0.17 (-0.01) 0.28 (-0.02) 30.28 (-0.94) 50.86 (-2.49) 64.54 (-3.08) FlowNet-2 + warp 7.35 (+0.03) 11.23 (-0.03) 13.17 (-0.04) 0.05 (+0.00) 0.15 (+0.00) 0.24 (+0.00) 25.42 (-0.20) 51.15 (-0.03) 65.29 (+0.23)

			-BDPSNR ↓			-BDSSIM ↓			-BDVMAF ↓
	latency compensa-	33	100	167	33	100	167	33	100	167
	tion (ms)
	Copylast 10.78 (SDCNet 5.20 (+0.14)	8.28 (+0.10)	9.74 (+0.09)	0.05 (-0.02)	0.15 (-0.07)	0.24 (-0.13) 17.27 (-0.13) 36.00 (-0.13) 46.63 (-0.09)

ACKNOWLEDGMENTS

This work was funded by the ANR AAPG2020 national fund (ANR-20-CE25-0014).