Jérémy Nadal
email: jeremy.nadal@polymtl.ca

Simon Brown
email: simon.brown@polymtl.ca

Elsa Dupraz
email: elsa.dupraz@imt-atlantique.fr

François Leduc-Primeau
email: francois.leduc-primeau@polymtl.ca

Towards an Accurate High-Level Energy Model for LDPC Decoders

Keywords: LDPC, ASIC, density evolution, energy efficiency

Estimating the energy consumption of LDPC decoders is a long and difficult task due to the large number of factors involved. Modern circuit synthesis tools can provide a satisfactory estimation of the power consumption, but this requires that the circuit be already implemented and it can take hours to provide the estimate. Currently, no accurate models are available to evaluate the decoding energy early in the design process. We propose a high-level energy model for flip-flop memory elements in LDPC architectures. The originality of the model is that it can analytically evaluate the variation of the energy due to the switching activity of the circuit gates, depending on the probability mass function (PMF) of the circuit inputs. Such PMFs are obtained through an adapted density evolution method that we propose. Therefore, the energy can be profiled for each decoding iteration and SNR value while considering several architecture choices. We illustrate the validity of the model by comparing the obtained energy estimates with measurements based on circuit simulations.

I. INTRODUCTION

As communication speeds increase and data transmission channels become more complex, the use of error-correction codes becomes increasingly critical. Low-density parity-check (LDPC) codes [START_REF] Gallager | Low-density parity-check codes[END_REF] are widely used in modern communication systems such as 5G networks [START_REF] Bae | An overview of channel coding for 5G NR cellular communications[END_REF]. For future 6G networks, preliminary studies consider that peak data rate will be up to 50 times higher that of existing 5G networks [START_REF] Akyildiz | 6g and beyond: The future of wireless communications systems[END_REF]. Therefore, there is a need for massively parallel and pipelined hardware architectures, which raises the issue of their energy consumption. This impacts the battery life of user devices and the cost of electrical consumption in base stations for the operators. Particularly, energy efficiency now becomes an important metric for studies of 6G networks [START_REF] Akyildiz | 6g and beyond: The future of wireless communications systems[END_REF], and must be taken into account when designing LDPC codes, hardware architectures, and circuits.

Developing novel energy reduction techniques requires an accurate understanding of how the code construction, decoding algorithm, and hardware architecture affect energy consumption. Previous work studied the energy consumption of coding schemes in the asymptotic regime as the block error probability goes to zero [START_REF] Ganesan | On the total power capacity of regular-ldpc codes with iterative message-passing decoders[END_REF], [START_REF] Blake | Energy, latency, and reliability tradeoffs in coding circuits[END_REF]. However, these results do not This work was supported by grant ANR-17-CE40-0020 (project EF-FECtive) of the French National Research Agency ANR and by grant 2021-NC-286323 of Fonds de recherche du Québec -Nature et technologie.

provide a way to optimize the energy consumption of finitelength codes and decoders. Reference [START_REF] Ganesan | On the total power capacity of regular-ldpc codes with iterative message-passing decoders[END_REF] also proposed a nonasymptotic energy model for regular LDPC codes. Their model however is based on post-layout simulations of specific circuit architectures and does not take the actual switching activity into account, which prevents using it to optimize the decoder architecture. As such, there is a lack of accurate analytical methods that can predict the energy consumption of an LDPC decoder before it is fully implemented. Current methods are mainly empirical, where accurate estimation of the energy consumption requires time-consuming post-synthesis or postlayout simulations to generate signal activities [START_REF] Nadal | A deeply pipelined, highly parallel and flexible LDPC decoder[END_REF].

Therefore, we propose in this paper a model to analytically estimate the energy consumption of LDPC decoders. As a first step, we focus on the energy consumption of flip-flop (FF)based components such as registers and memories, since these components represent most of the energy consumption of the decoder [START_REF] Amador | Optimum LDPC decoder: A memory architecture problem[END_REF]. The main originality of the approach is in the proposal of a clock-gated FF energy model where switching activity is taken into account. To predict the bit switching probabilities, a density evolution (DE)-based methodology is proposed. DE is a powerful analysis tool normally used for predicting the bit-error rate (BER) performance of LDPC decoders and designing good LDPC codes [START_REF] Richardson | The capacity of low-density paritycheck codes under message-passing decoding[END_REF]. The proposed method aims to adapt DE to predict the probabilistic behaviour of signals in LDPC decoder circuits.

The rest of this paper is organized as follows. Section II describes the offset min-sum (OMS) decoding algorithm and its DE analysis. Section III proposes a novel energy model for clock-gated FF based on bit switching probabilities. In Section IV, these switching probabilities are evaluated through novel DE-based methodologies adapted for most LDPC decoder architectures. The proposed energy model is then compared and validated in Section V, with post-synthesis energy measurements from a highly parallel ∆-update LDPC decoder architecture. Finally, Section VI concludes this paper.

II. LDPC DECODER

A. Decoding algorithm

We consider an Additive White Gaussian Noise (AWGN) channel described by y i = x i + w i , where y i is the channel output, x i is the i th transmitted modulated symbol, and w i is the additive Gaussian noise term of variance σ 2 . For simplicity, BPSK modulation is considered. The received and quantized channel log-likelihood ratio (LLR) L i of the i th coded bit is expressed as

L i = sat Q 2αy i σ 2 + 1 2 , ( 1
)
where y i is the received signal, σ 2 is the noise variance, α is a constant scaling factor, ⌊.⌋ is the floor operator and sat Q is a saturation operator which ensures that

L i ∈ [-Q, Q], Q ∈ N *
being the maximum LLR magnitude. Most common LDPC decoding algorithms are variants of the Sum-Product algorithm [START_REF] Kschischang | Factor graphs and the sum-product algorithm[END_REF], in which messages are iteratively sent between variable nodes (VNs) and check nodes (CNs). Let us denote by γ (ℓ) i→j the check-to-variable (C2V) message sent from CN i to VN j for decoding iteration ℓ ≥ 0, and by λ (ℓ) j→i the V2C message sent from VN j to CN i at iteration ℓ. Initially, we set γ (0) i→j = 0 ∀(i, j), and λ (0) j→i = L i ∀j. In addition, V i is the set of VN indices connected to CN i, and C j is the set of CN indices connected to VN j. In this paper, we consider the widely used OMS decoder [START_REF] Fossorier | Reduced complexity iterative decoding of low-density parity check codes based on belief propagation[END_REF], which can be formulated through the following equations:

γ (ℓ) i→j = j ′ ∈Vi\j s( λ(ℓ) j ′ →i ) × max min ∀j ′ ∈Vi\j | λ(ℓ) j ′ →i | -β, 0 , (ℓ) j→i = L j + i ′ ∈Cj \i γ (ℓ) i ′ →j , (2)
where β ∈ N is an offset to be optimized,

λ(ℓ) j→i = sat Q (λ (ℓ)
j→i ), and s(.) is the signum function. The decision on the j th bit is taken by evaluating the sign of the a-posteriori (AP) belief Λ (ℓ) j assigned to VN j and defined as Λ

(ℓ) j = λ (ℓ) j→i + γ (ℓ) i→j .

B. Density evolution

Initially developed for Sum-Product decoding in [START_REF] Richardson | The capacity of low-density paritycheck codes under message-passing decoding[END_REF], DE is an analysis tool which uses probabilistic properties of belief propagation to predict the behavior of a decoder at each iteration. That is, assuming that a given LDPC code is cyclefree, we can accurately predict the performance of a decoder at each iteration by expressing the successive probability distributions of messages exchanged in the decoder. The works presented in [START_REF] Balatsoukas-Stimming | Density evolution for minsum decoding of LDPC codes under unreliable message storage[END_REF], [START_REF] Ngassa | Density evolution and functional threshold for the noisy min-sum decoder[END_REF] extend DE to quantized OMS decoding. Both the initial DE analysis [START_REF] Richardson | The capacity of low-density paritycheck codes under message-passing decoding[END_REF] and the extensions provided in [START_REF] Balatsoukas-Stimming | Density evolution for minsum decoding of LDPC codes under unreliable message storage[END_REF], [START_REF] Ngassa | Density evolution and functional threshold for the noisy min-sum decoder[END_REF] assume that the length of the LDPC code is infinite.

It is shown in [START_REF] Richardson | The capacity of low-density paritycheck codes under message-passing decoding[END_REF] that the BER performance of the decoder is the same for any transmitted codeword. Therefore, to simplify equations, DE considers that the all-zero codeword is transmitted. In such case, the channel LLR cumulative distribution function (CDF), denoted

Φ L (k), k ∈ [-Q, Q], is Φ L (k) = 1 2 + 1 √ 4σ 2 erf k + 1 2 σ 2 α -1 , (3)
where erf(x) = 2/ √ π x 0 exp [-t 2 ]dt is the Gauss error function. The probability-mass function (PMF) of the channel LLR, denoted P L (k), can then be expressed as

P L (k) = Φ L (k) -Φ L (k -1), ∀k ∈ [-Q + 1, Q -1], P L (-Q) = Φ L (-Q) and P L (Q) = 1 -Φ L (Q).
The equation of the C2V messages PMF P (ℓ) γi→j , the V2C messages PMF λ (ℓ) j→i and the AP P (ℓ) Λ(k) can be found in [START_REF] Balatsoukas-Stimming | Density evolution for minsum decoding of LDPC codes under unreliable message storage[END_REF], [START_REF] Ngassa | Density evolution and functional threshold for the noisy min-sum decoder[END_REF] . The PMF of the saturated V2C message λ(ℓ) j→i can then be obtained by respectively summing all P (ℓ)

λj→i (k), k > Q and k < Q, to P (ℓ) λj→i (Q) and P (ℓ) λj→i (-Q). The BER at iteration ℓ, denoted B (ℓ)
e,inf , can be calculated by evaluating the probability that Λ j is negative.

DE equations are only valid assuming a cycle-free code of infinite length, and provide asymptotic BER performance. The cycle-free assumption implies that messages are statistically independent, whereas the infinite-length assumption implies that the observed noise power within the received block averages to the channel noise power. This second assumption is what makes the DE unable to characterize the BER of finitelength codes [START_REF] Leduc-Primeau | Finite-length quasi-synchronous LDPC decoders[END_REF]. However, the observed channel error rate can be modeled as a Gaussian random variable [START_REF] Yazdani | Waterfall performance analysis of finitelength LDPC codes on symmetric channels[END_REF]. As shown in [START_REF] Leduc-Primeau | Finite-length quasi-synchronous LDPC decoders[END_REF], [START_REF] Dupraz | Noisy density evolution with asymmetric deviation models[END_REF], this model provides a good match to Monte-Carlo simulations in the waterfall region for moderately large code lengths.

III. PROPOSED ENERGY MODEL

This section proposes an energy model for clock-gated FFs, which serve to build registers and memory units. All data storage in an LDPC decoder is handled by one of these two components, and they represent most of the energy consumption of the decoder [START_REF] Amador | Optimum LDPC decoder: A memory architecture problem[END_REF]. Therefore, an accurate energy model for an LDPC decoder can be developed based on an accurate model of the FF energy consumption. Specifically, most of the energy consumed by a FF is the result of the charge and discharge of the circuit's internal capacitance as a result of inputs switching between low and high voltage values. Therefore, it is critical to consider switching activity of the decoder signals when modeling energy consumption.

In a synchronous hardware architecture, each signal S, composed of N b = ⌈log 2 (2Q+1)⌉ bits, can be associated with a random variable of PMF P

= [P (c) b0 , ..., P (c) b N b -1 ]
the length-N b probability bit vector. For sign-magnitude (SM) binary representation, the probability bit vector can be derived using

P (c) b = T P (c)
S , where T is a N b × (2Q + 1) transformation matrix. Each element T i,j of T is defined as

T i,j = ⌊ |j -Q -1| × 2 -i+1 ⌋ mod 2 if 1 ≤ i < N b , 1+s(Q+1-j) 2 if i = N b .
When the two's complement (2C) representation is used, the transformation matrix is

C, with C N b ,j = T N b ,j , ∀j, C i,j = T i,j , ∀i when j ∈ [Q + 1, 2Q + 1] and C i,j = ⌊ (2 N b + j - Q -1) × 2 -i+1 ⌋ mod 2 when j ∈ [1, Q] and i ̸ = N b .
Note that the SM and 2C representations are often used together in the same decoder architecture. Using these transformation matrices, the bit probabilities can be derived for any signal if their associated PMFs P (c) S are known. It is then possible to estimate the energy consumption based on the bits switching activity. In this section, we propose to develop such a model for clock-gated flip-flop components.

A. Flip-flop and register energy model

The energy model we propose considers FFs whose input clock ports are connected to a clock gating unit (CGU). The CGU provides the ability to disable the FF clock when no new value needs to be stored, which greatly reduces its energy consumption. To derive the energy model, the following three assumptions are considered:

(A1): The CGU enable signal state, denoted S e (c) ∈ {0, 1}, is considered fully known. This assumption is quite realistic as enable signals are scheduled by a dedicated controller unit. We also define S e (c) = 1 -S e (c).

(A2): The input bit state of the FF is modeled as a random variable. We respectively define P (c) + (r) and P (c) -(r) the probabilities that the r th input bit of a length-R flip flop array (register) is equal to 1 after ("+") and before ("-") the front edge of a given clock event c. Note that we have

P (c) -(r) = P (c-1) + (r).
This assumption is perfectly realistic since the decoder input consists of noisy channel samples, which are random by nature.

(A3): It is considered that the input bit state before and after a clock event are statistically independent. Therefore, the probability P (c) s (r) that the r th input bit switches state at the clock event number c is P

(c) s (r) = P (c) + (r)(1-P (c) -(r))+(1- P (c) + (r))P (c) -(r).
The validity of this assumption rests on the fact that under the assumption of a cycle-free Tanner graph, any two messages associated with distinct VNs and distinct CNs are independent. In practice, LDPC codes contain cycles that introduce correlations between messages. As explained in Section II-B, these correlations usually do not have a significant effect on the decoding in the waterfall region, which is of interest for evaluating energy consumption. Therefore, (A3) is reasonable as long as processing elements are associated with distinct Tanner graph nodes in two successive clock cycles.

In the absence of clock gating, we consider that the energy consumption E (c) F F (r) of the r th FF in a register for clock cycle c comes from the switching activity of the clock and input signals, regardless of the FF output state. Therefore, we have:

E (c) F F (r) = E C + E sc P (c) s (r) , ( 4
)
where E C is the energy consumed by the clock activity and E sc is the energy consumed when the FF input switches state. When adding clock gating, the energy added by the CGU is

E (c) CG = E ce S (c) e + E ce S (c) e , (5)
where E ce and E ce correspond to the CGU energy when its input-enable signal is set to 1 and 0, respectively. The parameter E ce depends on the fanout of the CGU sub-clock port. For a clock-gated register composed of W FFs, E ce depends on W since the CGU sub-clock drives W FFs. The energy consumption of such register, denoted E (c) W , is simply the sum of the energy of the W FFs and CGUs:

E (c) W = E (c) CG + S (c) W r=1 E (c) F F (r) + S (c) ē E sc W r=1 P (c) s (r) E (c) W = E (c) CG + W E C S (c) e E (c) 0 (W ) +P (c) S S (c) e E sc + S (c) ē E sc E (c) S , (6)
where

P (c) S = W r=1 P (c) s (r)
, and E sc is the energy consumed by the FF when the sub-clock is disabled and input bits are switching. From the above equation, it can be inferred that the total energy consumption of the clock-gated length-W register is composed of a deterministic part (E

B. Energy model of memory units

We now use the FF energy model presented previously to derive the energy consumption of memories built out of standard-cell FFs. The standard-cell-based memory (SCM) reference architecture is depicted in Figure 1. Although SCMs are not area-efficient, they can be more energy-efficient than static random-access memory (SRAM) macros when operating at low voltages [START_REF] Meinerzhagen | Benchmarking of standard-cell based memories in the sub-vt domain in 65-nm CMOS technology[END_REF]. They are thus a good choice when energy efficiency is the primary design objective.

The memory storage matrix is composed of D M rows of W M FFs, with D M being the depth of the memory and W M the data width. The i th input data bit is connected to all the FF inputs located at the i th column. All the FF outputs at a given column i are fed to a D-to-1 multiplexer. When considering clock gating, D M CGUs are instantiated, and each CGU outputs a dedicated sub-clock, generated from the main system clock, for each FF row. These sub-clocks can be enabled or disabled by the CGU depending on the write address and the write-enable signals. If the write-enable signal

E (c) M,e (D M ) = D M E ce + E sc W M r=1 P (c) s (r) . ( 7
)
It is worth mentioning that the input signal activities have a non-negligible impact on the energy consumption when the memory is not solicited (D M W M E sc in the worst case).

When the write-enable signal is activated, only one FF row is activated, and the energy of the memory unit becomes

E (c) M,e (D M ) = E ce + W M r=1 E (c) F F (r) + E (c) M,e (D M -1). (8)
A non clock-gated register can optionally be added at the output of the multiplexer. This adds E (c)

W M to the energy obtained with the above equations.

IV. DERIVATION OF SIGNAL ACTIVITIES

The previous section described a clock-gated FF energy model based on signal activities, assuming the PMFs of the decoder signals P (c) S are known. An analytical evaluation of these PMFs is however not straightforward if typical DE methods are employed. In particular, some signals in LDPC decoders are not well modeled with typical DE, such as signals associated with finding the first and second minima in the OMS decoder. Furthermore, the all-zero transmitted codeword assumption cannot be used when evaluating switching activity. This section aims to adapt the DE equations to solve these issues.

A. Signal PMFs in the processing core

Figure 2 shows the architecture of the main processing unit used in a ∆-update LDPC decoder [START_REF] Nadal | A deeply pipelined, highly parallel and flexible LDPC decoder[END_REF]. In this figure, the PMFs for the main signals of interest are labelled. Most of these PMFs can be derived using the equations presented in Section II-B. One exception concerns the signal output of the "2-min search tree" unit, whose function is to find the first and second minimum values of d c (i)

∆ = |V i | input
(ℓ) j→i ∆ = Q -| λ(ℓ) j→i |. In such case, P (ℓ) maxu,i (k) ∆ = P (ℓ) minu,i (Q -k), u ∈ {1, 2}, corresponds to the (d c (i) -r + 1) th -order statistic of the {Y (ℓ) j→i } ∀j∈Vi random variables with CDF Φ (ℓ)
Yj→i . The r th -order statistic CDF is equal to the probability that at most d c (i) -r variable values are strictly superior to k [START_REF] David | Order Statistics[END_REF]. Therefore, we have

Φ (ℓ) max1,i (k) = j∈Vi Φ (ℓ) Yj→i (k) (9)
Φ (ℓ) max2,i (k) = Φ (ℓ) max1,i (k) + j∈Vi Φ (ℓ) Yj→i (k) j ′ ∈Vi\j Φ (ℓ) Yj→i (k), (10)
where

Φ (ℓ) Yj→i (k) = 1 -Φ (ℓ) Yj→i (k) and Φ (ℓ) maxu,i (k) = k l=0 P (ℓ)
maxu,i (l) corresponds to the CDF of the first (u = 1) and second (u = 2) maximum values of the {Y (ℓ) j→i } ∀j∈Vi random variables.

Equations ( 9) and ( 10) can be used to derive the signal activities of all input registers storing the first and second minimum message value, with the exception of the internal registers in the "2-min search tree." This tree of elementary 2-min units is composed of ⌈log 2 DC MAX ⌉ levels. Each level can optionally be terminated by a pipeline register in order to increase the processing throughput of the decoder. To derive the activity of these registers, ( 9) and ( 10) must be adapted by replacing V i with V ′ i (s, u), where V ′ i (s, u) ⊆ V i is the subset of VN indices that are part of the computational tree for the u-th 2-min unit in level s of the tree.

It is worth mentioning that the proposed model can be easily adapted for other Min-Sum decoder architectures compatible with assumption (A3) introduced in Section III-A. In such case, the registers and memories are used to store the same variables, so the signal activities are similar. Therefore, adapting the model for different architectures only requires to identify where the registers are located.

B. Density evolution with random codewords

Standard DE [START_REF] Richardson | The capacity of low-density paritycheck codes under message-passing decoding[END_REF] provides the message PMFs through successive decoding iterations, under the assumption that the allzero codeword is transmitted. However, this assumption does not allow to describe the message PMFs in a real decoder, since random codewords are transmitted in practice. Therefore, standard DE cannot be straightforwardly applied to derive the signal switching activities. Instead, we consider that the bits b j affected to the VN indexes j ∈ V i are random with probability 1/2 of being 1 or 0, under the condition that the transmit bits fulfill parity check equation at the i th row of the PCM: ⊕ j ′ ∈Vi b j ′ = 0, with ⊕ being the logic exclusive-or operator.

Under this condition, it can be shown that the PMF of the C2V messages, denoted

P (ℓ) γi→j | ⊕ (k), becomes P (ℓ) γi→j | ⊕ (k) = 1 2 P γ (ℓ) i→j = k | ⊕ j ′ ∈Vi\j b j ′ = 0
+ P γ (ℓ) i→j = k | ⊕ j ′ ∈Vi\j b j ′ = 1 = 1 2 P (ℓ) γi→j (k) + P (ℓ) γi→j (-k) , ( 11
)
where

P (ℓ) γi→j (k)
is the PMF of the C2V messages when the all-zero codeword is transmitted.

Similarly, the V2C message PMF given the transmitted bits fulfill the |C j | local parity-check equations is expressed as

P (ℓ) λj→i| ⊕ (k) = 1 2 P (ℓ) λj→i (k) + P (ℓ) λj→i (-k) . ( 12
)
V. VALIDATION

The aim of this section is to estimate the parameters of the clock-gated FF model proposed in Section III, and to validate the model by comparing the predicted energy with the one estimated on post-synthesis circuits.

A. Estimation of the energy model parameters

The energy parameters of the FF model are presented in Table I. They are obtained through power estimation performed using the Cadence Genus tool on simple test circuits composed of registers and memory units. Estimations are based on switching activity files generated through postsynthesis simulations. The architectures are synthesized using TSMC CMOS 65nm GP technology. The E ce and E ce energy parameters are directly obtained through the power estimation tools. The parameter E C is then deduced by setting all FF input bits to 0 so that no activities are generated. Then, E sc and E s c are obtained by estimating the energy consumption when the input bits switch states at each clock cycle, for both cases where the clock-gating is enabled and disabled. The CGU energy when the sub-clock is enabled E ce depends on W . All W FFs being identical, E ce is expected to be linear with the fanout of the CGU. Thus it is modeled as E ce (W ) = E α W +E β , which was confirmed to closely match empirical results.

B. Validation of the proposed energy model

To validate the energy model, the check node processing element (CNPE) architecture proposed in [START_REF] Nadal | A deeply pipelined, highly parallel and flexible LDPC decoder[END_REF] and shown in Figure 2 The energy is estimated through post-synthesis simulations using the same tools presented in Section V-A. To evaluate the impact of the signal activities on the energy consumption, the energy is estimated for each decoding iteration independently. To do so, the AP samples are generated based on their expected PMFs at a given iteration and clock cycle index for 100 codewords. This method artificially removes the cycles in the code, and has negligible impact on the decoder performance evaluation as explained in Section II-B. The PMFs of the received LLRs are set according to (3) with σ 2 = 0.7079 and α = 3. Since the aims of this section is to validate the accuracy of the energy model with respect to measurements from postsynthesis simulations, the choice of the code can be arbitrary. The model must be able to predict the energy consumption independently of the type of input PMFs being considered in the CNPE. Therefore, for simplicity, we consider that the PMFs of the CNPE input signals are obtained assuming an infinite-length code based on the following protograph S P [START_REF] Dupraz | Low-latency LDPC decoding achieved by code and architecture co-design[END_REF]:

S P = 0 2 3 1 2 0 3 2 . ( 13
)
The validation and study of the energy model for the complete LDPC decoder and several codes is left for future work. Figure 3 shows the predicted energy obtained with the proposed model and the energy estimated on the CNPE circuit for registers storing the APs, identified 1 in Figure 2, the registers storing the V2C messages 2 and the registers storing the first/second minimum variables min 1 /min 2 3 . All energies are normalized per register and per iteration. For all registers, the predicted energy at each iteration follows the same trend as the circuit simulation. The energy estimated from the circuit is only between 2% (V2C register) and 8% (min2 register) higher than the predicted energy, which is acceptable for a high-level model. This slight difference can be explained by the presence of additional switching due to the signals propagating through the logic units with different delays. This additional activity, known as glitches, is not taken into account in the proposed model, since this would require information on propagation delays, which is only available after circuit synthesis. Figure 3 also shows the register energy when glitches are removed. In this case, the predicted energy is within 1% of the simulation result. This demonstrates the accuracy of the model and the validity of the proposed approach. Figure 4 shows the predicted energy obtained with the proposed model and the energy estimated on the C2V memory unit of the CPNE circuit. The predicted energy follows similar variations as the estimation, particularly during the first iterations. This shows the importance of considering the signal activity for predicting the energy consumption. After the 7 th iteration, the model slightly deviates, where the energy is up to 4% higher than the estimated ones, which is acceptable for a high-level model. The total energy of all the registers, predicted by the model, is also shown in this figure for different CNPE pipeline depths D P (CNPE). Note that increasing pipeline depth allows achieving a higher processing throughput. It can be seen that the energy consumed by the registers is of the same order of magnitude as the memory unit. Furthermore, the register energy greatly depends on the number of pipeline stages, with up to 32% energy increase between D P (CNPE) = 4 and D P (CNPE) = 8. This demonstrates the importance of considering the register energy and the number of pipeline stages when studying the energy consumption of LDPC decoders.

VI. CONCLUSION

This paper proposes a novel high-level energy model for LDPC decoders. The model predicts the energy consumption of clock-gated FFs and standard-cell-based memory units based on signal activities. The switching activity of each bit is derived through a novel DE method adapted for this energy model. The model is validated by comparing the predicted energy with the energy estimated by a power analysis tool after circuit synthesis. Future work will study how this energy model can be used to optimize the energy efficiency of LDPC decoders, for instance in terms of quantization bits, parallelism level and number of pipeline registers, as well as to guide the construction of energy-efficient codes.

SS

  (k) at a given clock event c. Equivalently, S can be represented in a bitwise format, with b i being the i th bit of S. Each bit b i is random (but not independent or identically distributed) with probability P (c) bi of being equal to 1. Let us denote P (Q)] the length-2Q+1 PMF vector of S, and P (c) b

Fig. 1 .

 1 Fig. 1. Reference clock-gated FF-based memory (example with D M = 4 and W M = 3).

  Scorresponds to an upper bound.

Fig. 2 .

 2 Fig. 2. Internal signal PMFs of a check-node processing unit

  is considered. A CNPE can process up to DC MAX = 26 messages in parallel at each clock cycle. The CNPE input corresponds to DC MAX APs of 8 bits and outputs DC MAX Λ (ℓ) j variables on N b + 1 bits, with N b = 5 being the number of bits used to represent C2V messages (including the sign).

Fig. 3 .

 3 Fig. 3. Energy per iteration of the CNPE registers.

Fig. 4 .

 4 Fig. 4. Energy per iteration of the C2V memory unit and all the registers.

  V2C message magnitudes. This output is connected to several registers and to the "C2V message memory." Therefore, the associated PMFs, denoted P

	Finding the two smallest values among the V2C message
	magnitudes {|	λ(ℓ) j→i |} ∀j∈Vi is equivalent to finding the two
	largest values among the {Y	(ℓ)
	(ℓ) min1,i and P min2,i , are critical to assess the energy (ℓ)
	consumption of the decoder.

j→i } ∀j∈Vi variables, with Y


Belaid Hamoum

Elsa Dupraz

Laura Conde-Canencia

Dominique Lavenier

Channel Model with Memory for DNA Data Storage with Nanopore Sequencing

This paper is dedicated to channel modeling and error-correction coding for DNA data storage with nanopore sequencing. We first propose a novel statistical model for DNA storage, which takes into account the memory within DNA storage error events, and follows the way nanopore sequencing works. Compared to existing channel models, the proposed model represents more accurate experimental datasets. We also propose a full error-correction scheme for DNA storage, based on a consensus algorithm and non-binary LDPC codes. Especially, we introduce a novel synchronization method which allows to eliminate remaining deletion errors after the consensus, before applying a belief-propagation LDPC decoding algorithm to correct substitution errors. This method exploits the LDPC code structure to correct deletions, and does not require adding any extra redundancy.

I. INTRODUCTION

Over the last decade, the amount of generated data has been growing up exponentially, and [START_REF] Reinsel | The digitization of the world. from edge to core[END_REF] predicted that data storage needs would grow from 45 zetabytes in 2019 to 175 zetabytes by 2025. For this reason, it is necessary to find alternatives to classical storage methods (tapes, HDD, SSD, etc.). Among these, DNA data storage [START_REF] Goldman | Towards practical, high-capacity, lowmaintenance information storage in synthesized DNA[END_REF] appears as a promising solution that benefits from highly increased density and durability compared to existing storage solutions. Over the last years, improvements in DNA synthesis and sequencing techniques have made this technology more affordable, and an automated end-to-end DNA data storage device demonstrator was reported in [START_REF] Takahashi | Demonstration of End-to-End Automation of DNA Data Storage[END_REF]. This demonstrator includes the MinION, i.e., third-generation single molecule sequencing from Oxford Nanopore Technologies (ONT). This device offers a portable, real time, low-cost, and long-reads sequencing [START_REF] Leggett | A world of opportunities with nanopore sequencing[END_REF].

However, up-to-date DNA data storage systems remain very sensitive to errors introduced during synthesis and sequencing [START_REF] Heckel | A Characterization of the DNA Data Storage Channel[END_REF]. These two operations introduce not only substitution errors, like conventional memory storage systems, but also insertions and deletions at a non-negligible rate. Conventional decoding algorithms of existing error-correction solutions such as Turbo codes or Low Density Parity Check (LDPC) codes are not sufficient for DNA storage, because initially designed to efficiently protect sequences against substitution errors only.

In order to design efficient error-correction solutions for DNA storage, we first need to develop accurate statistical This work was supported by the Labex Cominlabs with funding from the French National Research Agency (ANR-10-LABX-07-01). We also acknowledge Emeline Roux for kindly providing the experimental datasets used in this work. channel models that allow for error-correction code design and simulation, without having to use very costly experimental processes. DNA storage channel models currently used for error-correction performance evaluation either assume nonrealistic independent and identically distributed (i.i.d.) errors [START_REF] Song | Codes With Run-Length and GC-Content Constraints for DNA-Based Data Storage[END_REF], or rely on Deep-Learning approaches [START_REF] Chandak | Overcoming High Nanopore Basecaller Error Rates for DNA Storage via Basecaller-Decoder Integration and Convolutional Codes[END_REF], which makes them difficult to interpret.

This paper introduces on a novel statistical DNA data storage channel model, based on a Markov chain with memory of order k, that models the MinION sequencing technology. Since this model consists of explicit probability terms, its assumptions (memory length, dependency to input sequence, etc.) are explicit as well. Further, once inferred on experimental data, this model can be easily interpreted (for instance, to identify which patterns in the input sequence are the most sensitive to specific types of errors). We show from numerical results that our model represents more accurately the experimental channel compared to the state of the art [START_REF] Song | Codes With Run-Length and GC-Content Constraints for DNA-Based Data Storage[END_REF], [START_REF] Chandak | Overcoming High Nanopore Basecaller Error Rates for DNA Storage via Basecaller-Decoder Integration and Convolutional Codes[END_REF].

The paper also addresses the design of efficient errorcorrection solutions for DNA storage. Apart from substitutions, most existing practical error-correction codes can correct only insertions [START_REF] Xue | Construction of GC-Balanced DNA With Deletion/Insertion/Mutation Error Correction for DNA Storage System[END_REF], or only deletions [START_REF] Xue | Construction of GC-Balanced DNA With Deletion/Insertion/Mutation Error Correction for DNA Storage System[END_REF], [START_REF] Abroshan | Coding for Deletion Channels with Multiple Traces[END_REF]. Existing solutions [START_REF] Lenz | Concatenated Codes for Recovery From Multiple Reads of DNA Sequences[END_REF] which can efficiently correct the three types of errors have a very high complexity, except for the i.i.d. channel model which is not realistic in the DNA storage application. To fully reconstruct the original sequence, a commonly employed solution consists in a consensus/trace-reconstruction algorithm followed by standard channel decoding to correct remaining residual errors [START_REF] Yazdi | Portable and Error-Free DNA-Based Data Storage[END_REF]- [START_REF] Lenz | Coding Over Sets for DNA Storage[END_REF]. However, because a few amount of insertions and deletion are likely to remain after the consensus, the reconstruction process may fail. In order to correct the remaining errors, it was proposed in [START_REF] Chandak | Improved read/write cost tradeoff in DNA-based data storage using LDPC codes[END_REF] to augment a conventional channel coding solution with periodical markers. These markers allow to correct insertions and deletions by resynchronizing the output sequence with the original one, while the channel code takes in charge the substitution errors.

In this work, we consider a consensus algorithm together with Non-Binary (NB) LDPC codes that fit the quaternary DNA alphabet. As an alternative to [START_REF] Chandak | Improved read/write cost tradeoff in DNA-based data storage using LDPC codes[END_REF] in order to handle the few remaining insertions and deletions, we propose a novel intermediate step between the consensus and the channel decoding. This intermediate step consists in re-synchronizing the sequence at the output of the consensus, by relying on the LDPC code structure rather than on additional markers. This represents an interesting gain in terms of coding redundancy, which may allow to reduce the expensive synthesis costs.

Section II describes the DNA data storage workflow. Section III presents existing DNA storage channel models, and Section IV introduces our channel model with memory. Section V presents the proposed scheme for error-correction. Section VI is dedicated to simulation results.

II. DNA DATA STORAGE WORKFLOW

A. Synthesis

As described in [START_REF] Hao | Current and Emerging Methods for the Synthesis of Single-Stranded DNA[END_REF], DNA synthesis consists of transforming digital data into experimental DNA strands, where each DNA strand is formed by a series of nucleotides 'A', 'C', 'G' and 'T'. Different synthesis techniques exist (Chemical, Enzymatic, Bacteria-Based, etc.), each exhibiting different constraints, in terms of e.g., sequences length, costs, etc. Data sequences used in this work were synthesized using the chemical technique, in which oligonucleotides (short DNA molecules) are synthesized and then assembled to form the ordered sequences of nucleotides. After the synthesis, we get not one but thousands of copies of the synthesized sequence.

B. Nanopore sequencing

Sequencing is the physical process where DNA strands are read. In this work, we consider the ONT MinION [START_REF] Leggett | A world of opportunities with nanopore sequencing[END_REF]. In this sequencing technique, strands go through a nanopore nucleotide by nucleotide, and an electrical signal (current intensity) is sequentially output for every group of k successive nucleotides, called k-mers. Then, a utility software called the basecaller translates the current signal at the output of the sequencer into the k nucleotides corresponding to each kmer. Multiple basecallers exist (Guppy, Albacore, Flappie,..), and they use either event segmentation or Deep-Learning approaches [START_REF] Wick | Performance of neural network basecalling tools for Oxford Nanopore sequencing[END_REF]. We considered Guppy for this work as it is the official basecaller provided by ONT. The basecaller generates fastQ files that contain the read sequences as well as related metadata.

III. EXISTING DNA DATA STORAGE CHANNEL MODELS

One major drawback of DNA data storage resides in the fact that both synthesis and sequencing introduce errors (insertions, deletions, substitutions) [START_REF] Heckel | A Characterization of the DNA Data Storage Channel[END_REF]. From an information-theoretic perspective, all the previous steps (synthesis, sequencing, basecalling) can be modeled as a channel, which takes as input a sequence composed of bases 'A', 'C', 'G', 'T', and outputs N edited replicas of the original sequence. In this section we introduce some notation and describe state-of-the-art channel models for DNA storage.

A. Notation

In what follows, we use Ins, Del, Sub, as abbreviations for Insertion, Deletion, Substitution, respectively. In addition, Match stands for "no error". In our model, X is the channel input sequence of length I, in which X i ∈ {A, C, G, T } is the base at position i. Then, kmer i = (X i-k+1 . . . X i-1 X i ) is the k-mer of length k at position i. And Y is a sequence of length I of channel events, where Y i ∈ {Ins, Del, Sub, Match} is the event at position i. We consider that a given DNA storage channel model is defined by a set of probability distributions for the successive Y i . For instance, P(Y i ) is the marginal probability of Y i , and P(Y i |kmer i , Y i-1 ) is the conditional probability of Y i with respect to kmer i and the previous event Y i-1 . Finally, we use Z to denote edited sequences at the output of the basecaller.

B. Independent and identically distributed channel model

In the literature, Y is often considered to be independent and identically distributed (i.i.d). Events Y i are further assumed to be independent from the X i , and the probabilities P(Y i ) are either inferred from sets of experimental data, or fixed arbitrarily for simulation purposes only [START_REF] Song | Codes With Run-Length and GC-Content Constraints for DNA-Based Data Storage[END_REF], [START_REF] Chandak | Improved read/write cost tradeoff in DNA-based data storage using LDPC codes[END_REF]. However, although it can simulate the correct amount of errors, this model cannot represent bursts of errors, since it assumes that Y i is independent from Y i-1 . In addition, it does not take into account the effect of kmer i onto error event Y i , and we know from several other works [START_REF] Wick | Performance of neural network basecalling tools for Oxford Nanopore sequencing[END_REF] that a statistical relation between those exists.

C. DeepSimulator

DeepSimulator is another popular tool to simulate the DNA data storage system [START_REF] Chandak | Overcoming High Nanopore Basecaller Error Rates for DNA Storage via Basecaller-Decoder Integration and Convolutional Codes[END_REF], [START_REF] Li | DeepSimulator: a deep simulator for Nanopore sequencing[END_REF]. DeepSimulator relies on a Deep-Learning approach combined with a basecaller. In a first Deep-Learning-based step, DeepSimulator takes as input a sequence of bases, and outputs electrical current levels after the sequencing. In a second step, the current levels are sent to a basecaller which transcripts the current levels into bases. Event sequences Y generated by DeepSimulator contain some memory and are statistically dependent from the input sequences X. However, after having performed a significant amount of simulations we could observe that DeepSimulator does not always reflect well the physical DNA storage channel. For instance, as shown in Figure 1, the same kind of error appears in most of the simulated sequences in the same particular position, which does not correspond to an experimental channel model. Moreover, we could also observe an inaccurate predominance of substitutions and insertions over deletions. This shows the need to develop more accurate channel models for code design and performance evaluation.

IV. CHANNEL MODEL WITH MEMORY FOR DNA DATA

STORAGE

We now introduce a novel statistical channel model for DNA data storage, that aims at solving the drawbacks just described of the existing models. The proposed model takes into account the statistical dependencies between the event sequence Y and the input sequence X. It also considers that event Y i depends on kmer i , which allows to model editions due to successive k-mer reads, according to the way the MinION sequencer works. Our model also assumes some internal memory in Y, by considering that previous event Y i-1 can affect current event Y i This allows to take burst errors into account.

A. Probability terms

In this section, we use B ∈ {A, C, G, T } to denote the obtained base after a substitution. Our model is then described by the following conditional probability distribution:

• P(Y i |kmer i , Y i-1 ) characterizes the dependency between the current event Y i , the currently read kmer i , and the previous event Y i-1 . This captures the probability of different types of events Y i ∈ {Ins, Del, Sub, Match} depending on kmer i , and allows to consider bursts of errors through the dependency to Y i-1 .

• P(L|kmer i , Y i = Ins) characterizes the insertion length L depending on the currently read kmer i . Since only insertions can be of length L > 1, the probability of L is conditioned to event Y i = Ins.

• P(B|kmer i , Y i = Sub) characterizes the probability to substitute the last base X i of kmer i by the base B, given that Y i = Sub. We assume that the probabilities P(Y i |kmer i , Y i-1 ) and P(L|kmer i , Y i = Ins) do not vary with k < i < I. On the contrary, we observed from experimental data that the probability to get an error is higher at the first position i = 1 and at the last one i = I, compared with middle positions 1 < i < I. Therefore, we allow for different probability distributions P(Y ), P(L|Y 1 = Ins) and P(Y I ), P(L|Y I = Ins) to be used when i = 1 and i = I. Finally, for 1 < i k, since no complete k-mer was observed already, we consider probabilities P (Y i |X i ), P (B|X i , Y i = Sub) and P (L|X i , Y i = Ins) that only depend on the base value of X i .

B. Training

To train our model, i.e, to estimate all the considered probability terms, we used V = 9 sets of experimental data which went through the whole DNA data storage process, as described in Section II. Each of these sets provided one reference sequence taken as channel input X, and N sequences obtained after sequencing, taken as channel outputs Z. Then, the conditional probabilities P(Y i = D|kmer i , Y i-1 = E) were estimated by counting over the N V pairs (X, Y) the number of outcomes of each event D ∈ {Ins, Del, Sub, Match}, divided by the number of outcomes of the considered (kmer i , Y i-1 = E). The other probability terms were estimated following the same approach. Obtaining experimental data is very costly, while there are 4 k+1 different combinations (kmer i , E), which grows fastly with k. In our case, when a given combination (kmer i , Y i-1 = E) was left unobserved, we estimated the corresponding probabilities by averaging over all observed combinations. At the end, note that the considered training is not only specific to a given value k, but also to a particular DNA data storage process, with fixed synthesis and sequencing techniques. In particular, considering updated or other synthesis and sequencing techniques requires to retrain the model from new sets of experimental data.

C. Channel simulator

Once we get all the probabilities from training, we can build a simulator that takes as input a given sequence and generates outputs according to the channel model. This model simulates the whole DNA data storage process and has three main advantages. First, in case of technology evolution (synthesis, sequencing, basecalling,..), the model can be retrained with new sets of experimental data (adaptability). Second, since model also takes into account the basecaller, it is faster than DeepSimulator which requires to run the basecaller during simulations (faster simulations). Third, as opposed to a blackbox approach, all our probabilities terms are explicit. The values of these probabilities can then be used, either to better understand how errors are introduced during the DNA storage process, or to incorporate them to any source and channel coding methods which would be considered in the workflow (explainability).

V. FULL RECONSTRUCTION SOLUTION

We now propose an error-correction scheme to correct insertion, deletion and substitution errors after the sequencing step. This scheme relies on three components: a consensus algorithm, a resynchronization step, and a NB-LDPC decoder. The first one aims to provide a good quality sequence based on the sequencing data redundancy, while the last two ones aim to correct residual errors. The consensus algorithm relies on the fact that known primers (short sequence of about 40 nucleotides) are added at the beginning and at the end of the sequence X. These primers consist of sequences of known bases which are mainly used to select sequences of interest through biotechnology manipulations.

A. Consensus algorithm

The consensus algorithm we consider, called CCSA, has been specifically designed for the DNA data storage system introduced in Section II. CCSA takes as input m sequences (Z 1 , • • • , Z m ) selected randomly from the set of N sequences output by the basecaller. Then, for given parameters T and L, it forms a directed graph whose nodes are given by the subsequences of length L that appear at least T times among the m input sequences. There is an edge between two nodes of the graph if the corresponding two subsequences overlap by at least d bases. Finally, a Viterbi-like algorithm is applied over the directed graph in order to select path between the start primer (p 1 ) and the end primer (p 2 ) with the highest score. Note that CCSA may output either one or several sequences of equal (highest) score of different lengths J close to I. The score is calculated from the nodes weights (number of occurrence of the subsequence over the m sequences) and edges weights (L minus overlap length between the two subsequences) over the path. A detailed description of the CCSA algorithm can be found in [START_REF] Lavenier | Constrained Consensus Sequence Algorithm for DNA Archiving[END_REF]. CCSA results reported in [START_REF] Lavenier | Constrained Consensus Sequence Algorithm for DNA Archiving[END_REF], as well as our own simulations, show that this algorithm is able to correct most of the errors introduced by the DNA storage process, although a few residual insertions, deletions, and substitutions remain. This is why an additional correction step based on NB-LDPC decoders is needed.

B. NB-LDPC Codes and decoders

LDPC codes are linear capacity-approaching blocks codes commonly used in communication systems to correct substitution errors from the channel. NB-LDPC codes are defined over Galois fields GF(q) of length q ≥ 2 [START_REF] Poulliat | Design of regular (2, dc)-LDPC codes over GF (q) using their binary images[END_REF]. In this work, we use NB-LDPC codes in GF(4) for consistency with the quaternary bases alphabet. The parity check matrix H of the code is sparse, and its non-zero elements take values in GF(q). Then, any codeword X of the code verifies HX = 0. The standard Belief-Propagation (BP) LDPC decoder [START_REF] Davey | Low density parity check codes over GF(q)[END_REF] takes as input the consensus sequence X and seeks to outputs a sequence X close to X and that verifies the condition H X = 0. However, the standard BP decoder as well as most existing LDPC decoders can only correct substitution errors. This is why we now propose a synchronization method to also correct a few amount of deletions. We only apply this synchronization step if the consensus does not output a sequence of length J = I, and by assuming that output sequences of length J < I only result from deletions.

C. NB-LDPC codes synchronization

In this part, we define the score of a given sequence X as the number of unsatisfied parity check equations (that is the number of non-zero components in the vector H X). We now describe the synchronization process when the consensus outputs a sequence X of length I -1. In this process, for a given segment length l s , we try to insert a base with arbitrary value 'A' at position 1, then at position l s + 1, then 2l s + 1, and so on. For each considered position, we compute the corresponding score. At the end, we definitively add a base at the position kl s + 1 that gives the lowest score, thus creating a new vector Xs . Even if the base value is incorrect and the base position not entirely correct either, adding this base close to the correct position should greatly reduce the score by re-synchronizing the output sequence Xs with the original codeword X. Now, if the sequence X is of lower length Iv, v bases are inserted one after each other in a greedy maneer: the first base is inserted at the position which gives the lowest score, then the process is repeated so as to insert the second base, and so on. This process allows to replace deletion errors by substitution errors which can then be corrected by the LDPC BP decoder. If after synchronization, the decoder fails (unsatisfied parity check equations remain), a new consensus will be restarted. Note that the segment length l s plays a crucial role as it addresses a tradeoff between algorithm complexity and amount of substitution errors in the resulting sequence. Finally, this synchronization technique only relies on the LDPC code structure, and does not require any additional redundancy, unlike in a solution with periodical markers [START_REF] Chandak | Improved read/write cost tradeoff in DNA-based data storage using LDPC codes[END_REF].

D. Full reconstruction solution

We now describe our full reconstruction solution. We first provide to the consensus algorithm m sequences at random among the M sequences Z that have correct primers p 1 and p 2 and length greater than I. The consensus algorithm then outputs several sequences X of length J. If the consensus outputs a sequence of length J = I, we set this sequence as Xs . Otherwise, we pick a sequence of length J < I as close as possible to I, and we apply the synchronization step described in Section V-C in order to get a sequence Xs of the correct length I. In both cases, the sequence Xs is passed through a BP decoder to correct residual substitution errors. If the sequence X is incorrect in the sense that H X = 0, we restart the full reconstruct process (consensus + synchronization + BP decoder), from another set of m random sequences provided to the consensus. This process allows to reconstruct the original sequence X with a small number of restarts, as we now evaluate through simulations.

VI. SIMULATION RESULTS

In this part, we first compare the numerical simulations obtained with the proposed channel model to the existing ones (i.i.d. and DeepSimulator), and with the experimental data. Our channel model was trained as described in section IV, using 9 experimental datasets called barcodes, each of them containing multiple Fastq files related to one particular reference sequence. Since the MinION sequencer reads kmers of length 6, we fixed the channel memory length to k = 6, and we performed the training over 22182 output sequences. Figure 1 shows for barcode03 four edit maps (i.i.d. channel model, DeepSimulator, our model, experimental data) containing positions of deletions, insertions, and substitutions observed for 1000 output sequences after an alignment with the reference sequence. We can observe that for the i.i.d model, errors are incorrectly uniformly distributed over all the sequence, while for DeepSimulator same errors of the same type seem to always appearing at same positions. Finally, the proposed channel model with memory seems to approach the most the experimental data, although some long-run error events are not properly taken into account.

We then evaluate the performance of our full sequence reconstruction method, by using the proposed channel model with memory. To do so, we considered barcodes 01 and 03 among the nine at our disposal, and encoded the corresponding sequences with a regular (3, 6) LDPC code in GF(4) of size (500, 1000) and code rate R = 1/2, from a PEG algorithm. We further set segment length to l s = 50 for our synchronization method. To evaluate the proposed reconstruction method, we considered three setups: (i) consensus alone, (i) consensus + NB-LDPC decoder without resynchronization, (iii) consensus + re-synchronization + NB-LDPC decoder. For each of these setups, we applied the successive reconstruction steps 900 times, and evaluated the proportion of perfectly recovered sequences. This metric is of interest in our setup, since the condition HX = 0 allows to decide that the original sequence was correctly retrieved, and therefore to stop the reconstruction loop. Figure 2 shows the proportion of correctly retrieved sequences for the two considered barcodes, with respect to the number m of input sequences to the consensus. For both barcodes, as expected, the success probability increases with m, before reaching a peak at around m = 100. The decrease in performance after this peak probably comes from the fact that the consensus has difficulties to handle too many sequences, due to the initial majority voting operation. Then, we observe two different configurations, depending on the considered barcode. For barcode01, it seems that the resynchronization and NB-LDPC decoder are not useful, in the sense that the sequences at the output of the consensus are either correct, or very far from the original data. On the opposite, for barcode03, we observe that resynchronization greatly improves the reconstruction rate, while NB-LDPC decoder alone after consensus does not help. This is probably due to the fact that for this barcode, most sequences at the output of the consensus contain a few amount of deletions introduced by an homopolymer of length 6, which are corrected by the synchronization method.

Fig. 1 .Fig. 2 .

 12 Fig. 1. Errors and types of errors observed on barcode03, for 1000 output sequences.

VII. CONCLUSION

In this work, we proposed a channel model with memory for DNA storage. Compared to existing ones, our model represents experimental data more accurately, and should allow for efficient source/channel codes design. We then introduced a full reconstruction scheme, based on a consensus algorithm and on NB-LDPC codes, and proposed a synchronization method in order to correct remaining deletions after the consensus. Future work will include improving the channel model by considering larger sets of experimental data and optimizing the proposed error-correction scheme in terms of complexity, code rate, etc.


Elsa Dupraz

Lav R Varshney

Franc ¸ois Leduc-Primeau

Power-Efficient Deep Neural Networks with Noisy Memristor Implementation

This paper considers Deep Neural Network (DNN) linear-nonlinear computations implemented on memristor crossbar substrates. To address the case where true memristor conductance values may differ from their target values, it introduces a theoretical framework that characterizes the effect of conductance value variations on the final inference computation. With only second-order moment assumptions, theoretical results on tracking the mean, variance, and covariance of the layer-bylayer noisy computations are given. By allowing the possibility of amplifying certain signals within the DNN, power consumption is characterized and then optimized via KKT conditions. Simulation results verify the accuracy of the proposed analysis and demonstrate the significant power efficiency gains that are possible via optimization for a target mean squared error.

I. INTRODUCTION

When implemented on electronic circuits, Deep Neural Networks (DNNs) demand important data transfers between memories and processors, which severely affects their power consumption and latency. As an alternative emerging paradigm, computation-in-memory consists of performing computation operations directly within the memory [START_REF] Strukov | The missing memristor found[END_REF]. In this work, we consider implementing Deep Neural Networks (DNNs) from memristor crossbars developed for analog in-memory computation [START_REF] Yakopcic | Memristor crossbar deep network implementation based on a convolutional neural network[END_REF]- [START_REF] Li | Long short-term memory networks in memristor crossbar arrays[END_REF], where the synaptic weights of neural networks are encoded in memristance values [START_REF] Kendall | Training end-to-end analog neural networks with equilibrium propagation[END_REF]. When building a computational system from a memristor crossbar, one of the main difficulties is to set up memristance values with arbitrary precision [START_REF] Liu | A memristor-based optimization framework for artificial intelligence applications[END_REF]. Hence in this work, we consider the setting where memristance values are noisy, and evaluate the effect of noise on the inference phase of the DNN.

It is shown in [START_REF] Karbasi | Noise facilitation in associative memories of exponential capacity[END_REF] that noise can sometimes help in the inference performance of neural networks by getting them unstuck from local minima. Several works in the literature also investigate the robustness of DNNs to undesired noise [START_REF] Kim | Energy-efficient neural network acceleration in the presence of bitlevel memory errors[END_REF]- [START_REF] Chatterjee | Energy-reliability limits in nanoscale feedforward neural networks and formulas[END_REF], though without considering the in-memory computation framework. As a key issue, the DNN robustness is strongly related to the power consumption of its electrical circuit implementation [START_REF] Chatterjee | Energy-reliability limits in nanoscale feedforward neural networks and formulas[END_REF]. In this work, we investigate this key issue in the context of in-memory computing, by first aiming to predict the effect of noisy memristor values onto the final DNN computation. The physical substrate of memristors and their computational and noise properties make the mathematical problem quite different from the above works [START_REF] Karbasi | Noise facilitation in associative memories of exponential capacity[END_REF]- [START_REF] Chatterjee | Energy-reliability limits in nanoscale feedforward neural networks and formulas[END_REF].

We first introduce a theoretical analysis that consists of tracking the evolution of the first and second-order moments over the successive layers of a DNNs, as functions of the noise variance of memristor values. These recursive expressions are obtained from second-order Taylor expansions of the means, variances, and covariances of the layers outputs, taking into account the noisy linear computation realized from memristor crossbars, and the non-linear activation functions at each layer. This extends our previous moment analysis of [START_REF] Dupraz | Noisy in-memory recursive computation with memristor crossbars[END_REF] which considered noisy dot-product computation from memristor crossbars. The proposed analysis is very generic as it only performs a few assumptions on the mean and variance of each memristance value, without any further assumption on their statistical distribution. The accuracy of the proposed analysis is verified using Monte Carlo simulations. Then, as a second important contribution, we propose to optimize the power consumption of a DNN memristor-based implementation. In a first step, we show how to evaluate the expected power consumption of one memristor crossbar, by applying the second-order Taylor expansion approach to electrical formula that provide the power of one specific memristor in a crossbar. In a second step, we formulate a power optimization problem under constraints on the final variance at the network output, and provide analytical solutions to this problem by relying on the Karush-Kuhn-Tucker (KKT) conditions. Our numerical results demonstrate significant power gains compared to the non-optimized case.

The remainder of this paper is organized as follows. Section II presents the noisy memristor crossbar architecture. Section III develops the moment-evolution methodology.

II. NOISY ARCHITECTURES FOR DEEP NEURAL NETWORKS A. DNN computation

In this work, we consider a feedforward NN with T layers. The network input is a vector x of length K, and the network output is a vector y of length N . Layer t ∈ 1, T of the network is composed by N t neurons and outputs a vector x (t) of length N t . We use W (t) to denote the weight matrix of layer t, where the matrix W (t) is of size N t × N t-1 . For a given activation function f , layer t performs the following linear and non-linear operations:

z (t) = W (t) x (t-1) , for all t ∈ 1, T , (1)
x (t) k = f (z (t)
k ), for all t ∈ 1, T and all k ∈ 1, N t , (2) where z (t) is a vector of length N t . In this work, we consider f to be the sigmoid function f (z) = 1 1+e -z , which is commonly considered in several NN implementations. The theoretical analysis developed in this paper can be easily extended to other activation functions f , given they are twice differentiable.

B. Memristor-based implementation of DNN

In this work, we consider the setting where the non-linear operation (2) is realized either in the CMOS layer of a CMOL architecture [START_REF] Likharev | CMOL: Devices, circuits, and architectures[END_REF], by analog computation [START_REF] Kendall | Training end-to-end analog neural networks with equilibrium propagation[END_REF], or a similarly mathematizable physical substrate. There is no fundamental difference in considering either of these implementations and this implies that in our analysis, the noise is only on the conductance values of the memristors.

Then, the linear operation ( 1) is a dot-product computation which can be realized by a memristor crossbar. Since memristor conductance values can only be positive, the computation operation (1) is realized from two crossbars: one for the positive part, and one for the negative part of the computation [START_REF] Liu | A memristor-based optimization framework for artificial intelligence applications[END_REF]. Therefore, we let W (t) = W (t,+) -W (t,-) , where W (t,+) contains the positive components of W (t) , and W (t,-) contains the absolute values of the negative components of W (t) . The matrices W (t,+) and W (t,-) are then converted into matrices G (t,+) and G (t,-) , respectively, where G (t,+) and G (t,-) contain positive conductance values g (t,+) i,j and g (t,-) i,j . These conductance values are chosen so that the computation of the z (t) j in (1) can be equivalently performed as

z (t) j = Nt-1 i=1 g (t,+) i,j x (t) i g 0 + Nt-1 k=1 g (t,+) k,j - Nt-1 i=1 g (t,-) i,j x (t) i g 0 + Nt-1 k=1 g (t,-) k,j (3)
where g 0 is a pull-down conductance. Note that in expression (3), a part of the coefficients g (t,+) i,j and g (t,-) i,j are equal to zero, when the corresponding coefficients in W (t,+) and W (t,-) are equal to zero.

C. Noisy conductance values

When considering dot-product computation from memristor crossbars, one key issue is the difficulty in setting up conductance values g

(t,+) i,j , g (t,-) i,j
with arbitrary precision [START_REF] Liu | A memristor-based optimization framework for artificial intelligence applications[END_REF]. Therefore, we think of the conductance values as noisy, so that the deterministic quantities g [START_REF] Liu | A memristor-based optimization framework for artificial intelligence applications[END_REF]. We further assume that whenever g

(t,+) i,j , g (t,-) i,j , are replaced by random variables G (t,+) i,j , G (t,-) i,j , in
(t,+) i,j = 0, the corresponding random variable G (t,+) i,j has mean g (t,+) i,j
and variance σ 2 . In the same way, we assume that the random variables G (t,-) i,j such that g (t,-) i,j = 0 have mean g t,- i,j and the same variance value σ 2 . On the contrary, if g

(t,+) i,j = 0 (respectively g (t,-) i,j = 0), we assume that G (t,+) i,j = 0 (respectively G (t,-) i,j
= 0) as well. In order to develop an analysis which applies to a large class of problems, physical substrates, and noise environments, we do not make any further assumption beyond the second-order moments on e.g., the distribution of the random variables G (t,+) i,j and G (t,-) i,j . In addition, since the non-linear operations (2) are not realized from memristor crossbars, we assume that no additional noise is introduced during these operations. Finally, we denote the random versions of vectors x (t) and z (t) by X (t) and Z (t) .

III. PERFORMANCE ANALYSIS OF NOISY MEMRISTOR-BASED IMPLEMENTATION OF DNN

With the mathematical model of the noisy memristor-based DNN implementation in place, now we derive our theoretical performance analysis. We aim to express the MSE between the correct network output x (T ) and its noisy version X (T ) . For this, we follow the approach of [START_REF] Dupraz | Noisy in-memory recursive computation with memristor crossbars[END_REF], and provide iterative expressions of the first and second-order moments (means, variances, and covariances) of the successive random variables X (t) and Z (t) calculated at each layer t of the network. These iterative expressions are obtained from second-order Taylor expansions [START_REF] Seltman | Approximations for mean and variance of a ratio[END_REF] of the considered moments. Compared to [START_REF] Dupraz | Noisy in-memory recursive computation with memristor crossbars[END_REF] which focused on linear recursive computations, we must now take the non-linear computation (2) into account in the moment expressions. We must also consider the fact that the linear computation ( 3) is separated into one positive part and one negative part; the effect of this separation was not evaluated in [START_REF] Dupraz | Noisy in-memory recursive computation with memristor crossbars[END_REF]. Finally, in this section, we consider the viewpoint of one layer, and drop the index t in order to simplify the notation.

A. Moments after linear computation

In this part, we evaluate the second-order moments of the random vector Z after applying the linear operation (3) from memristor crossbars. We first rewrite the components Z j of Z as

Z j = Z (+) j -Z (-) j = T (+) j ∆ (+) j - T (-) j ∆ (-) j , (4)
where

T (+) j = N i=1 G (+) ij X i , ∆ (+) j = G 0 + N i=1 G (+)
ij , and likewise for T (-) j and ∆ (-) j . We then introduce the notation

E[X i ] = ν i , Var[X i ] = γ 2 i , and for all i = i, C[X i , X i ] = γ i,i
, respectively for the mean, variance, and covariance of the components X i of X. The next two theorems provide the mean, variance, and covariance of the positive components Z (+) j . Theorem 1. For all j ∈ 1, N , the second-order Taylor expansions of the mean µ

(+) j = E[Z (+) j ] and variance ρ (+) j 2 = V[Z (+) j ] of Z (+) j are given by µ (+) j = E[T (+) j ] E[∆ (+) j ] - C[T (+) j , ∆ (+) j ] E[∆ (+) j ] 2 + V[∆ (+) j ]E[T (+) j ] E[∆ (+) j ] 3 ρ (+) j 2 = V[T (+) j ] E[∆ (+) j ] 2 + 3E[T (+) j ] 2 V[∆ (+) j ] E[∆ (+) j ] 4 - 4E[T (+) j ]C[∆ (+) j , T (+) j ] E[∆ (+) j ] 3
where

E[T (+) j ] = N i=1 g (+) ij νi (5)
E[∆ (+) j ] = g0 + N i=1 g (+) ij (6) V[∆ (+) j ] = N (+) j σ 2 (7)
C[∆ (+) j , T (+) j ] = σ 2 N i=1 p (+) i,j νi (8)
V[T

(+) j ] = σ 2 N i=1 p (+) i,j (γ 2 i + ν 2 i ) + N i=1 N i =1 g (+) i,j g (+) i ,j γ i,i (9)
In the above expressions, p

(+) i,j = 1 if g (+) ij = 0, p (+)
i,j = 0 otherwise, and

N (+) j = N i=1 p (+) i,j .
Theorem 2. For all (j, j) ∈ 1, N 2 such that j = j, the covariance ρ

(+) j,j = C(Z (+) j , Z (+) j ) can be expressed as ρ (+) j,j = N i=1 N i =1 E G (+) i,j ∆ (+) j E G (+) i ,j ∆ (+) j γ i,i (10)
and the second-order Taylor expansion of E

G (+) i,j ∆ (+) j is E G (+) i,j ∆ (+) j = E[G (+) i,j ] E[∆ (+) j ] - C[G (+) i,j , ∆ (+) j ] E[∆ (+) j ] 2 + V(∆ (+) j )E[G (+) i,j ] E[∆ (+) j ] 3 where C[G (+) i,j , ∆ (+) j ] = σ 2 .
In order to obtain the mean, variance, and covariance of Z (-) j , it suffices to replace "+" by "-" in Theorems 1 and 2. The next theorem then provides the mean, variance, and covariance of the components Z j of Z.

Theorem 3. For all (j, j) ∈ 1, N 2 , and j = j, the mean

µ j = E[Z j ], variance ρ j 2 = V[Z j ],
and covariance ρ j,j = C(Z j , Z j ) are given by

µ j = µ (+) j -µ (-) j ρ 2 j = ρ (+) j 2 + ρ (-) j 2 -2C[Z (+) j , Z (-) j ] ρ j,j = ρ (+) j,j + ρ (-) j,j -C[Z (+) j , Z (-) j ] -C[Z (-) j , Z (+) j ] where C[Z (+) j , Z (-) j ] = N i=1 N i =1 E G (+) i,j ∆ (+) j E G (-) i ,j ∆ (-) j γ i,i .

B. Moments after non-linear computation

We now evaluate the means, variances, and covariances, of random variables X i = f (Z i ), after the sigmoid activation function f . Theorem 4. For all (i, i) ∈ 1, N 2 , and i = i, the secondorder Taylor expansions of the mean ν i = E[X i ] and variance γ 2 i = V[X i ] of the random variable X i , are given by

ν i = f (µ i ) + 1 2 f (µ i ) ρ 2 i γ 2 i = 1 2 g (µ i )ρ 2 i -f (µ i )f (µ i )ρ 2 i γ i,i = f (µ i )f (µ i )ρ i,i
where f (u) = 1 1+exp(-x) , and g = f 2 . This theorem can be adapted to any other twicedifferentiable activation functions (tanh, softplus, etc.), by replacing the functions f , g, and their derivatives, according to the newly considered activation function.

To summarize, the four theorems given in this section provide recursive expressions of the first and second-order moments of random vectors X (t) and Z (t) at successive layers t. In addition, the MSE at each layer can be estimated as the average of the variance terms ρ 2 j and γ 2 i . In the simulation section, we show from numerical simulations that these estimates provide good approximations of the successive MSE. At the end, the analysis we have developed allows us to investigate the effect of noisy conductance values onto the final NN output Z (T ) , depending on various parameters such as the crossbar size, the conductance noise variance, etc..

IV. POWER OPTIMIZATION OF MEMRISTOR-BASED DNNS

Having characterized the final MSE as a function of system parameters and chosen operating points, in this section we address the trade-off between the power consumption and the robustness to noise of memristor crossbars. We first evaluate the power consumption of a memristor crossbar, and then optimize this power consumption for a given target MSE at the network output. Here, we only evaluate the memory power consumption, and therefore do not take into account the non-linear activation parts, which may be realized in CMOS or similar technologies. As in Section III, we consider the viewpoint of one layer, and drop the index t in the notation.

A. Power consumption of a memristor crossbar

We assume that the variance σ 2 of conductance values

G (+) ij and G (-) ij
is fixed by internal properties of the memristors [START_REF] Starzyk | Memristor crossbar architecture for synchronous neural networks[END_REF], [START_REF] James | Analog neural computing with super-resolution memristor crossbars[END_REF]. Therefore, in order to make the computation more robust to noise, we consider coefficients c j such that target conductance values are now given by g(+) i,j = c j g (+)

i,j , g(-)

i,j = c j g (-)

i,j , and g0 = c j g 0 . This does not change the computation operation (3), since the coefficient c j appears as a factor in both parts of the ratio for Z j . The corresponding random variables G(+) i,j and G(-) i,j then have means c j g

(+)
i,j and c j g (-) i,j , respectively, and both have variance σ 2 . Note that it is physically possible for each component Z j to have its own parameter c j , which in turn leads to many degrees of freedom in the optimization.

Then, for the positive part of the computation, the power P (+) i,j consumed by the memristor at position (i, j) is a random variable which can be expressed as

P (+) i,j = G(+) i,j U 2 i,j
, where U i,j is the voltage across memristor at position (i, j). We can note that

U i,j = (X i -Z (+) j ) = X i - N i =1 G(+) i ,j X i G0 + N i =1 G(+) i ,j
. [START_REF] Likharev | CMOL: Devices, circuits, and architectures[END_REF] Since P i,j is a random variable, we are interested in its expectation E[P (+) i,j ], which can be rewriten as

E[P (+) i,j ] = E G(+) i,j Ṽ (+) i,j ∆(+) j ( 12
)
where

Ṽ (+) i,j = N i =1 G(+) i ,j (X i -X i ), and ∆(+) j = G0 + N i =1 G(+) i ,j .
The following theorem gives the second-order Taylor expansion of E[P (+) i,j ]. Theorem 5. For all (i, j) ∈ 1, N 2 , the second-order Taylor expansion of E[P (+) i,j ] is given by

E[P (+) i,j ] = E[ G(+) ij ]E[V (+) i,j 2 ] E[∆ (+) j ] 2 - 4E[V (+) i,j ]E[ G(+) i,j ]C[V (+) i,j , ∆ (+) j ] E[∆ (+) j ] 3 - 2E[V (+) i,j ] 2 C[ G(+) i,j , ∆ (+) j ] E[∆ (+) j ] 3 + 3E[ G(+) ij ]E[V (+) i,j ] 2 V[∆ (+) j ] E[∆ (+) j ] 4 where E[ G(+) ij ] = c j g (+) i,j , C( G(+) i,j , ∆ (+) j ) = σ 2 , V(∆ (+) j ) = N (+) j σ 2
, and

E[V (+) i,j ] = cj N i =1 g (+) i ,j (νi -ν i ) E[∆ (+) j ] = cj g0 + N i =1 g (+) i ,j C[V (+) i,j , ∆ (+) j ] = σ 2 N i =1 p (+) i,j (νi -ν i ) V[V (+) i,j ] = (i ,k) i =i,k =i,k = i c 2 j g (+) i ,j g (+) k,j (γ 2 i -γ i,i -γ i,k + γ i ,k ) + N i =1 p (+) i,j (σ 2 + c 2 j g (+) i ,j 2
)(γ 2 i + γ 2 i -2γ i,i ) + 2σ 2 (νi -ν i ) 2 .
Finally, for one crossbar, we express the expected total power as

E[P ] = i,j E[P (+) i,j ] + E[P (-) i,j ]
, where E[P (-) i,j ] is also obtained from Theorem 5 by replacing "+" by "-" in the various expressions. The average total power can then be evaluated for the full DNN architecture, by summing over all T crossbars.

B. Power optimization

In this subsection, we use ρj to denote the variance of component Z j , calculated given that the conductance values were multiplied by c j . In addition, we can show that ρj = ρ 2 j /c 2 j , where ρ 2 j was given in Theorem 3. Then, at each layer of the network, both the average power E[P ] and the variances ρ2 j depend on the parameters c j . In this section, we propose to optimize these parameters so as to minimize E[P ] under a certain performance constraint expressed onto the targeted variance of each component Z j of Z (t) . We treat layers one after the other, assuming that previous layers from 1 to (t -1) were already optimized when considering layer t. Therefore, we define the following optimization problem for layer t:

min cj i,j E[P (+) i,j ] + E[P (-) i,j ] s.t. ∀j, ρ2 j ≤ ρ, c j ≥ 0 ( 13
)
where ρ is a constraint on the maximum variance on each component Z j at the output of the linear part. We further constrain that c j ≥ 0 in order to get positive conductance values. Then, the Lagrangian of the constrained optimization problem ( 13) is given by

L = i,j E[P (+) i,j ] + E[P (-) i,j ] + j α j (ρ 2 j -ρ) -β j c j (14)
and can be rewritten as

L = N j=1 F 1,j c j + F 2,j c j + α j ρ 2 j c 2 j -ρ -β j c j (15)
where the terms F 1,j and F 2,j can be identified from the expressions of E[P i,j ] given in Theorem 5, and we show that F 1,j ≥ 0 and F 2,j ≤ 0. The Lagrangian is separable with respect to j, and its derivative with respect to c j is

∂L ∂c j = F 1,j - F 2,j c 2 j -2α j ρ 2 j c 3 j -β j . ( 16
)
We now apply the KKT conditions in order to solve the optimization problem. The constraint c j ≥ 0 is necessarily inactive (c j > 0, β j = 0) since c j = 0 does not allow us to satisfy ρ2 j ≤ ρ. In addition, considering an inactive constraint ρj < ρ with α j = 0 does not lead to a solution for c j . Therefore, the constraint ρj ≤ ρ is active, which leads to

c j = ρ 2 j /ρ. (17)
Finally, by calculating α j from the condition ∂L ∂cj = 0, the second-order derivative of the Lagrangian is given by

∂ 2 L ∂c 2 j = 1 c 2 j 3F 1,j c j - F 2,j c j > 0, (18)
which allows us to conclude that (17) provides a strict global minimum. Therefore, the values of c j in (17) allow us to minimize the power consumption of a crossbar, while satisfying the variance conditions ρ2 j ≤ ρ.

V. SIMULATION RESULTS

In this part, we first evaluate the accuracy of the proposed theoretical analysis. We consider a Neural Network structure with 7 layers, with the following number of neurons per layer: (100, 100, 200, 150, 120, 80, 10), and with sigmoid activation functions. We generate the weight matrices W (t) at random and uniformly between 0 and 10, and then convert them into conductance values g (t,+) i,j and g (t,-) i,j , with g 0 = 10. We also generate input vectors x at random, with components distributed uniformly between -5 and 5. In Figure 1, we compare the variance measured from Monte-Carlo simulations and evaluated with the theoretical analysis of Section III, with respect to the noise variance σ 2 , after one layer and after the seven layers. After one layer, we observe that the theoretical variance value accurately predicts the MSE measured from Monte Carlo simulations: the two curves are visually indistinguishable. In addition, as expected, the MSE increases with σ 2 . Then, after seven layers, we also observe that the theoretical variance values are close to the MSE measured from simulations, although there is a small gap which is probably due to error accumulation over the successive Finally, we note that after one layer, the MSE is higher or very close to the corresponding noise variance σ 2 . On the opposite, after seven layers, the MSE is actually smaller than the noise variance value σ 2 on memristors. This shows that memristor crossbar computation can be inherently robust to noise on conductance values. We then evaluate the power optimization method proposed in Section IV. To do so, we consider two standard datasets that are IRIS [START_REF] Fisher | The use of multiple measurements in taxonomic problems[END_REF] and MNIST [START_REF] Lecun | The MNIST database of handwritten digits[END_REF]. For IRIS (input length 4), we consider a NN with 2 layers and respectively 50 and 10 neurons on each layer, and for MNIST (input length 784), we consider a NN with 3 layers, and respectively 200, 50, 10, neurons on each layer. Both networks integrate sigmoid activation functions, and were trained over 2/3 of the dataset, with 100 epochs and learning rate η = 0.1. For both networks, we set σ = 10 -2 , and we optimize the parameters c j with the method of Section IV, for various target variance values ρ. The results are shown in Figure 2, where we also display the total power in the non-optimized case (this non-optimized total power is for comparison and does not necessarily allow to achieve the corresponding value of ρ). We observe a clear power gain after optimization, which is due to the fact that the optimized values c j allow to exactly reach the target variance ρ. In addition, as expected, the total power decreases as ρ increases, and the NN for MNIST requires more power than the NN for IRIS.

VI. CONCLUSION

In this paper, we considered a DNN implementation from noisy memristor crossbars, and we developed a theoretical analysis to predict the effect of noise onto the final inference computation performance. We further proposed a method to optimize the memristor crossbars power consumption, and demonstrated significant power gains compared to the nonoptimized case. Future works will include the study of other neural network architectures such as convolutional neural networks and recurrent neural networks.

  Section IV addresses power optimization. Section V provides simulation results. In what follows, E[•] refers to the expectation, V[•] to the variance, and C[•, •] to the covariance. The notation 1, Ψ indicates the set {1, 2, . . . , Ψ}.

Fig. 1 .

 1 Fig. 1. Comparison between MSE measured from Monte Carlo simulations and evaluated with the theoretical analysis developed therein. The red and blue curve (one layer) are superimposed.

Fig. 2 .

 2 Fig.2. Optimized and non-optimized total power calculated over the T layers, for two NNs. The non-optimized case does not necessarily achieve the final target variance.

This work was supported in part by the "Make our Planet Great Again" Initiative of the Thomas Jefferson Fund, by grant ANR-17-CE40-0020 of the French National Research Agency ANR (project EF-FECtive), and by IVADO grant PRF-2019-4784991664.