Adrien Morel Symme

David Gibus
email: david.gibus@univ-smb.fr

Adrien Badel
email: adrien.morel@univ-smb.fr

Exploring the power-limit of multi-modes multielectrodes vibration energy harvesters

Keywords: Vibrations, Energy harvesting, Piezoelectricity, Coupled dynamics, Multi-modes, Wideband, Resonance

In this article, we propose a model for multimodes multi-electrodes piezoelectric energy harvester, and study its dynamics when associated with an electrical interface. We show that a multi-modes harvester allows to overcome the traditional power limit defined for single-mode harvesters, thanks to the effect of the second mode on the first one. We investigate this phenomenon with a resistive load, both theoretically and experimentally, for a piezoelectric energy harvester made with an inverted mass and four pairs of piezoelectric electrodes. Multi-modes resonances allow to fully exploit the potential of multi-modes energy harvesters, and open ways to design large bandwidth large power density energy harvesters.

I. INTRODUCTION

Piezoelectric energy harvesters (PEH) have been widely used in order to harness energy from ambient vibrations and supply low-power electronic systems. PEH are usually made of a mechanical resonator that efficiently converts the ambient vibration in a strain applied on piezoelectric electrode(s), as long as the vibration frequency remains relatively close to the resonance frequency of the first mode [START_REF] Toprak | Piezoelectric energy harvesting: State-of-the-art and challenges[END_REF]. The mechanical [START_REF] Gibus | Strongly coupled piezoelectric cantilevers for broadband vibration energy harvesting[END_REF] and electrical [START_REF] Lallart | Synchronous inversion and charge extraction (SICE): a hybrid switching interface for efficient vibrational energy harvesting[END_REF][START_REF] Morel | A unified N-SECE strategy for highly coupled piezoelectric energy scavengers[END_REF] optimization of the harvested power with such PEH has been widely investigated in the literature. An expression of the power limit that can be scavenged with such PEH has been derived, and illustrates that the maximum power that can be harvested with a PEH is proportional to its quality factor, to its effective mass, and to the squared vibration amplitude [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF][START_REF] Halvorsen | Architecture-independent power bound for vibration energy harvesters[END_REF][START_REF] Liao | Optimal power, power limit and damping of vibration based piezoelectric power harvesters[END_REF]. In this paper, we prove both theoretically and experimentally that the traditional power limit can be overcome, thanks to multi-modes multielectrodes PEH.

Multi-modes PEH have been investigated in the last decades to obtain relatively large power bandwidth PEH and to harvest energy from broadband vibrations [START_REF] Erturk | Modeling of Piezoelectric Energy Harvesting from an L-shaped Beam-mass Structure with an Application to UAVs[END_REF][START_REF] Ferrari | Piezoelectric multifrequency energy converter for power harvesting in autonomous microsystems[END_REF][START_REF] Tang | A multiple-degree-of-freedom piezoelectric energy harvesting model[END_REF]. For instance, [START_REF] Erturk | Modeling of Piezoelectric Energy Harvesting from an L-shaped Beam-mass Structure with an Application to UAVs[END_REF] proposed an L-shaped energy harvester whose second natural frequency is relatively close to the first one ( 𝑓 ! /𝑓 " ≈ 2 ). [START_REF] Ferrari | Piezoelectric multifrequency energy converter for power harvesting in autonomous microsystems[END_REF] designed a PEH made of an array of cantilevers with different frequency responses that are all fixed on the same rigid base. Such approach allows to easily bring closer the modes frequencies and to widen the harvesting bandwidth at the cost of a larger size and lower power density. [START_REF] Tang | A multiple-degree-of-freedom piezoelectric energy harvesting model[END_REF] proposed the analysis of a N degree-offreedom (DOF) PEH made with multiple small masses fixed on a single large mass. It is notably shown that it is possible to dimension such structure in order to obtain multiple close peaks in the power response. To date, however, most multimodes PEH analysis in the literature do not take into account the influence of the electrical interface on the PEH dynamics, and the coupling between the resonance modes that might appear due to this electrical interface.

In this paper, we propose to analyze multi-modes multiports PEH with coupled electrical interfaces. We first show that the traditional power limit, defined for single-mode single-port PEH, can be overcome thanks to multi-modes resonances. Thereafter, we prove, both theoretically then experimentally, that the effect of closer modes on single-mode resonance can appear with a well-chosen resistive interface. This proves that multi-modes PEH can be used not only to widen the harvesting bandwidth, but also to increase the power density of PEH.

II. DYNAMICS OF MULTI-MODES MULTI-ELECTRODES PEH

A. Electromechanical model

In this paper, we consider a resonant structure with 𝑁 resonance modes and covered by 𝑀 piezoelectric electrodes driven by ambient acceleration (Fig. 1). Assuming that the structure behaves linearly, the 𝑛 #$ mode of the structure (∀𝑛 ∈ [1, 𝑁])can be associated with an equivalent stiffness 𝐾 % , an equivalent linear damper 𝐷 % , an equivalent mass 𝑀 % , and an equivalent force factor 𝐵 & %

. The 𝑚 #$ piezoelectric electrode (∀𝑚 ∈ [1, 𝑀]) impacts the 𝑛 #$ mode dynamics with an electrically-induced force 𝐹 %' = 𝛼 %' 𝑣 ' , with 𝑣 ' being the voltage across the 𝑚 #$ piezoelectric electrode and 𝛼 %' the piezoelectric coefficient (force factor) between the 𝑚 #$ piezoelectric electrode and the 𝑛 #$ resonance mode. Note that 𝐵 & % ≠ 𝑀 % as detailed in [START_REF] Erturk | Issues in mathematical modeling of piezoelectric energy harvesters[END_REF]. Thus, the dynamics of the 𝑛 #$ mode is governed by [START_REF] Toprak | Piezoelectric energy harvesting: State-of-the-art and challenges[END_REF],

𝑀 % 𝑥̈% + 𝐷 % 𝑥̇% + 𝐾 % 𝑥 % + > 𝛼 %' 𝑣 ' = -𝐵 & % 𝛾 ( ')" (1)
with 𝑥 % the displacement of the inertial mass of the resonator due to the 𝑛 #$ mode, and 𝛾 the base acceleration of the structure. The 𝑚 #$ piezoelectric electrode is under strain, which generates electrical charges and contributes to the evolution of the 𝑚 #$ voltage, 𝑣 ' . Since all modes contribute to the strain of the 𝑚 #$ piezoelectric electrode, the equation governing the evolution of the voltage is given by (2),

𝑖 ' = 𝐶 * ' 𝑣̇' -> 𝛼 %' 𝑥̇% + %)" (2)
with 𝐶 * ' and 𝑖 ' the capacitance and the electrical current flowing out from the 𝑚 #$ piezoelectric electrode, respectively. The dynamics of a 𝑁 -modes 𝑀 -electrodes

N-modes structure Electrode 1 Electrode 2 Electrode 3 Electrode 4 !! ! " !# ! $
structures can be expressed thanks to 𝑁 mechanical equations [START_REF] Toprak | Piezoelectric energy harvesting: State-of-the-art and challenges[END_REF] and 𝑀 electrical equations [START_REF] Gibus | Strongly coupled piezoelectric cantilevers for broadband vibration energy harvesting[END_REF], and the equivalent circuit model of such harvester is shown in Fig. 2.

B. Boundaries of the harvested power

The extracted power from the harvester 𝑃 ,-# can be calculated from the total power extracted from the 𝑀 electrodes (3).

𝑃 ,-# = > 𝑃 ,-#,' ( ')" = > 1 𝑇 E 𝑖 ' 𝑣 ' / 0 𝑑𝑡 ( ')" (3)
𝑇 being the period of the periodic excitation. Replacing 𝑖 ' by its expression (2) yields [START_REF] Morel | A unified N-SECE strategy for highly coupled piezoelectric energy scavengers[END_REF].

𝑃 ,-# = 1 𝑇 > E H𝐶 * ' 𝑣̇' -> 𝛼 %' 𝑥̇% + %)" I 𝑣 ' / 0 𝑑𝑡 ( ')" (4)
Since 𝑣 ' is a 𝑇-periodic function, ∫ 𝑣̇'𝑣 ' / 0 𝑑𝑡 = 0. Thus, the expression of the harvested power becomes [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF].

𝑃 ,-# = - 1 𝑇 > > E 𝛼 %' 𝑣 ' / 0 𝑥̇%𝑑𝑡 ( ')" + %)" (5)
Combining (1) with [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF] and considering that 𝑥 % is a 𝑇periodic function (∫ 𝑥̇%𝑥 % / 0 𝑑𝑡 = ∫ 𝑥̈%𝑥̇% / 0 𝑑𝑡 = 0) yields [START_REF] Halvorsen | Architecture-independent power bound for vibration energy harvesters[END_REF].

𝑃 ,-# = 1 𝑇 > E 𝑥̇% L𝐵 & % 𝛾 + 𝐷 % 𝑥̇%M / 0 𝑑𝑡 + %)" (6)
In the case of a harmonic excitation, 𝛾(𝑡) = 𝛾 ' cos (𝜔𝑡).

Since 𝑥 % is a 𝑇-periodic function 𝑥 % = ∑ [𝑎 1 % cos(𝑘𝜔𝑡) + 23 1)" 𝑏 1 % sin(𝑘𝜔𝑡)],
with 𝑎 1 % and 𝑏 1 % the coefficients giving the inphase and out-of-phase components of the 𝑘 #$ harmonic of 𝑥 % , respectively. The derivative of 𝑥 % , 𝑥̇%, can be expressed by

𝑥̇% = 𝜔𝑘 ∑ [-𝑎 1 % sin(𝑘𝜔𝑡) + 𝑏 1 % cos(𝑘𝜔𝑡)] 23 1)"
. Combining these expressions with (6) and using Parseval equality lead to the power expression [START_REF] Liao | Optimal power, power limit and damping of vibration based piezoelectric power harvesters[END_REF].

𝑃 !"# = 𝜔 2 & '𝛾 $ 𝑎 % & 𝐵 ' ! -𝜔𝐷 & & 𝑘 ( .𝑎 ) & ( + 𝑏 ) & ( 1 *+ ),% 2 - &,% (7)
𝑃 ,-# is maximized if the optimality conditions given by ( 8) are respected,

⎩ ⎨ ⎧∀𝑛, 𝑎 " % = 𝛾 ' 𝐵 & ! 2𝜔𝐷 % = 𝛾 ' 𝐵 & ! 𝑄 % 2Ω % 𝑀 % 𝜔 % ! ∀𝑘 > 1, ∀𝑛, 𝑎 1 % = 0 ∀𝑘 > 0, ∀𝑛, 𝑏 1 % = 0 (8)
with 𝑄 % = 𝐷 % 4" ^𝑀% 𝐾 % the quality factor of the 𝑛 #$ mode, 𝜔 % = ^𝐾% /𝑀 % the resonance pulsation of the 𝑛 #$ mode, and Ω % = 𝜔/𝜔 % the normalized vibration frequency. If the conditions given by ( 8) are respected, the harvested power is equal to the power limit of the harvester 𝑃 56' , whose expression is given by ( 9),

𝑃 56' = 𝛾 ' ! 4 > 𝐵 & ! ! 𝑄 % 𝑀 % 𝜔 % + %)" = > 𝑃 56' % + %)" (9)
with 𝑃 56' % the single-mode traditional power limit valid for the 𝑛 #$ resonance mode [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF]. Some conclusions can be drawn from ( 7), ( 8) and ( 9):

• The optimal displacement waveform is a sine wave (in the case of harmonic excitation) whose phaseshift with the excitation is equal to 90°, and whose amplitude depends on the considered mode and on the excitation of amplitude. • The higher harmonics of the mechanical displacement always decrease the harvested power. This means that electrical non-linearities (with circuits such as SECE and SSHI [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF]) or mechanical non-linearities (with multi-stable harvesters, or stiffening/hardening behaviors [START_REF] Jia | Review of nonlinear vibration energy harvesting: Duffing, bistability, parametric, stochastic and others[END_REF]) might reduce the harvested power and unable to reach the power limit. However, in most case, the impact of the higher harmonics of the mechanical displacement remains relatively negligible (first-harmonic assumption [START_REF] Morel | A unified N-SECE strategy for highly coupled piezoelectric energy scavengers[END_REF][START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF]). • The power limit is equal to the sum of all the power limits associated with each mode. This means that the simultaneous exploitation of the multiple modes of the energy harvester allows to harvest more power than the traditional power limit that is associated with a single mode [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF]. This phenomenon where the harvested power can go beyond the single-mode power limit is referred as multi-modes resonance.

C. Optimal voltages across the electodes

Because they are both periodic functions, the mechanical displacement associated with the 𝑛 #$ mode and the piezoelectric voltage across the 𝑚 #$ electrode can be expressed as the sums of in-phase and out-of-phase components:

𝑥 % = ∑ [𝑎 1 % cos(𝑘𝜔𝑡) + 𝑏 1 % sin(𝑘𝜔𝑡)] 23 1)" and 𝑣 ' = ∑ [𝑐 1 ' cos(𝑘𝜔𝑡) + 𝑑 1 ' sin(𝑘𝜔𝑡)] 23 1)"
. Combining these expressions with (1) and replacing 𝑎 1 % and 𝑏 1 % by their optimal values given by ( 8) lead to optimality conditions on the electrodes voltages components, 𝑐 1 ' and 𝑑 1 ' .

⎩ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎧ ∀𝑛, > 𝛼 %' 𝑐 " ' ( ')" = 𝛾 ' 𝐵 & ! 𝑄 % (Ω % ! -1) 2Ω % ∀𝑛, > 𝛼 %' 𝑑 " ' ( ')" = - 𝐵 & % 𝛾 ' 2 ∀𝑘 > 1, ∀𝑛, > 𝛼 %' 𝑐 1 ' = > 𝛼 %' 𝑑 1 ' = 0 ( ')" ( ')" (10)
From [START_REF] Tang | A multiple-degree-of-freedom piezoelectric energy harvesting model[END_REF], one can draw the following conclusions:

• In order to attain the power limit 𝑃 56' with 𝑁 resonance modes, 2𝑁 conditions should be met by the piezoelectric voltages across all electrodes. In order to respect all these conditions, the number of electrodes should be equal or greater than the number of modes considered (𝑀 > 𝑁).

! ! " ! 1/% ! & " ' # ! ' # $ ' # % ! $ " $ 1/% $ ! % " % 1/% % … … … … ( %! ) ! ( %$ ) $ ( %" ) " ( $! ) ! ( $$ ) $ ( $" ) " ( !! ) ! ( !$ ) $ ( !" ) " ∑ ( %" +̇% & %'! -" . ( " ) " ∑ ( %! +̇% & %'! -! . ( ! ) ! ∑ ( %$ +̇% & %'! -$ . ( $ ) $ … +̇% !! +̇!
Input vibration Mechanical structure (! modes) M piezoelectric electrodes

• The third equation of [START_REF] Tang | A multiple-degree-of-freedom piezoelectric energy harvesting model[END_REF] shows that the harmonics linked with a single electrode can be compensated by the harmonics linked with another electrode, in order to reach the power limit even if the electrical interface is non-linear. This mean that some electrodes can be used to suppress the harmonic influences of other electrodes, avoiding harmonic generation on the mechanical displacement.

III. ANALYSIS OF THE HARVESTED POWER ON A RESISTOR

It is not always possible to respect the optimality conditions given by ( 8) and [START_REF] Tang | A multiple-degree-of-freedom piezoelectric energy harvesting model[END_REF] when the vibration frequency is far from one of the mode resonance frequency. To investigate this, a harvester of two modes (𝑁 = 2) and two electrodes ( 𝑀 = 2 ) is considered. The mechanical and electrical dynamics of such system are given by ( 11) and [START_REF] Jia | Review of nonlinear vibration energy harvesting: Duffing, bistability, parametric, stochastic and others[END_REF],

⎩ ⎨ ⎧𝑀 " 𝑑 ! 𝑥 " 𝑑𝑡 ! + 𝐷 " 𝑑𝑥 " 𝑑𝑡 + 𝐾 " 𝑥 " + 𝛼 "" 𝑣 " + 𝛼 "! 𝑣 ! = 𝐵 & " 𝛾 𝑀 ! 𝑑 ! 𝑥 ! 𝑑𝑡 ! + 𝐷 ! 𝑑𝑥 ! 𝑑𝑡 + 𝐾 ! 𝑥 ! + 𝛼 !" 𝑣 " + 𝛼 !! 𝑣 ! = 𝐵 & ! 𝛾 (11)
b

𝑖 " = 𝑣 " 𝑅 " = 𝐶 * " 𝑣̇" -𝛼 "" 𝑥̇" -𝛼 !" 𝑥̇! 𝑖 ! = 𝑣 ! 𝑅 ! = 𝐶 * ! 𝑣̇! -𝛼 "! 𝑥̇" -𝛼 !! 𝑥̇! ( 12
)
where 𝑅 " and 𝑅 ! correspond to the resistive loads connected to the electrodes 1 and 2, respectively. Note that in the particular case where the two electrodes are in parallel, the two currents 𝑖 " and 𝑖 ! can be summed and the total current is: 𝑣/𝑅 = d𝐶 * " + 𝐶 * ! e 𝑣̇-𝛼 "" 𝑥̇" -𝛼 !" 𝑥̇! -𝛼 "! 𝑥̇" -𝛼 !! 𝑥̇! From [START_REF] Jia | Review of nonlinear vibration energy harvesting: Duffing, bistability, parametric, stochastic and others[END_REF], the relation in frequency domain between the voltage across each electrode and the displacement of each mode can be obtained,

𝑣 ' = d𝛼 "' 𝑥 " + 𝛼 !' 𝑥 ! e 𝐶 * ' f𝜀 7 " + 𝑗𝜀 8 ' i (13)
with 𝜀 8 " and 𝜀 7 ' being the dimensionless electricalinduced damping and stiffness associated with the 𝑚 #$ electrode, respectively [START_REF] Morel | Resistive and reactive loads' influences on highly coupled piezoelectric generators for wideband vibrations energy harvesting[END_REF]. Note that the expressions of 𝜀 8 " and 𝜀 7 ' remain identical to the single-mode single-electrode case [START_REF] Morel | Resistive and reactive loads' influences on highly coupled piezoelectric generators for wideband vibrations energy harvesting[END_REF].

⎩ ⎪ ⎪ ⎨ ⎪ ⎪ ⎧ 𝜀 8 " = d𝑅 ' 𝐶 * ' 𝜔e 1 + d𝑅 ' 𝐶 * ' 𝜔e ! 𝜀 7 " = d𝑅 ' 𝐶 * ' 𝜔e ! 1 + d𝑅 ' 𝐶 * ' 𝜔e ! ( 14
)
From ( 11) and ( 13), the expressions of the modal displacements in frequency domain can be found by solving the following system:

j 𝐵 & " 𝛾 = [𝐴 " + 𝑗𝐵 " ]𝑥 " + [𝐶 + 𝑗𝐷]𝑥 ! 𝐵 & ! 𝛾 = [𝐴 ! + 𝑗𝐵 ! ]𝑥 ! + [𝐶 + 𝑗𝐷]𝑥 " (15)
with 𝐴 1 , 𝐵 1 , 𝐶 and 𝐷 (𝑘 ∈ {1,2}) being given by [START_REF] Pereira | Operating frequency prediction of piezoelectric dc-dc converters[END_REF].

⎩ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎧ 𝐴 1 = 𝐾 1 -𝜔 ! 𝑀 1 + 𝛼 1" ! 𝐶 * " 𝜀 7 " + 𝛼 1! ! 𝐶 * ! 𝜀 7 ! 𝐵 1 = 𝜔𝐷 1 + 𝛼 1" ! 𝐶 * " 𝜀 8 " + 𝛼 1! ! 𝐶 * ! 𝜀 8 ! 𝐶 = 𝛼 "" 𝛼 !" 𝐶 * " 𝜀 7 " + 𝛼 "! 𝛼 !! 𝐶 * ! 𝜀 7 ! 𝐷 = 𝛼 "" 𝛼 !" 𝐶 * " 𝜀 8 " + 𝛼 "! 𝛼 !! 𝐶 * ! 𝜀 8 ! (16)
The expression of the apparent power of the 𝑚 #$ electrode can be expressed from [START_REF] Morel | Resistive and reactive loads' influences on highly coupled piezoelectric generators for wideband vibrations energy harvesting[END_REF]:

𝑆 ' = 𝑣 ' 𝑣 ' * 2𝑅 ' (17)
The harvested power with each electrode 𝑃 ' is equal to the real part of 𝑆 ' . Its analytical expression is too long to be given here, but can be easily obtained from ( 15) and ( 17). In all the following simulations, the values of the dimensional parameters ( 𝑀 " , 𝑀 ! , 𝐷 " , 𝐷 ! , … ) are chosen to obtain 𝑄 " = 𝑄 ! = 50, 𝜔 ! = 2𝜔 " = 200 𝑟𝑎𝑑. 𝑠 4" , and 𝑃 :;< " = 2𝑃 56' ! . The values of the four normalized coupling 𝑘 ' ! | %' = 𝛼 %' ! / (𝐾 % 𝐶 * ' ) corresponding to the electromechanical coupling between the 𝑛 #$ resonance mode with the 𝑚 #$ piezoelectric electrode, are given for each simulation.

Figure 3 shows the maximum total harvested power (with the two electrodes) with optimized resistances 𝑅 " and 𝑅 ! as a function of the vibration frequency, with 𝑘 ' ! | "" = 𝑘 ' ! | !! = 0.2 , and 𝑘 ' ! | "! = 0 . The blue curve shows the maximum power when the coupling between the first electrode and the second mode is null (𝑘 ' ! | !" = 0). In this case, the power does not go beyond the power limit 𝑃 :;< " . In the other hand, when a non-null coupling is introduced in the simulation (𝑘 ' ! | !" = 0.1), the maximum power goes beyond 𝑃 56' " from about 18%. Such simple example illustrates how a bi-modal resonance can appear thanks to multi-modes coupling. The increase of power around this bi-modal resonance is due to the charges generated from the mechanical displacement of the second mode that is in-phase with the mechanical displacement of the first mode. Even though the second mode displacement remains small (because the bi-modal resonance frequency is far from the second mode resonance frequency), the small voltage generated by this displacement adds up with the already-large voltage generated by the displacement of the first mode. Therefore, the harvested power is greatly increased, because it is quadratically related to the total voltage. Note that the closest the two modes frequencies are, the greater the power around the bi-modal resonance.

! !"# $ ! !"# % Bi-modal resonance due to ! ! " | "# # # # "
Figure 4 shows that such bi-modal resonance can be obtained even if 𝑘 ' ! | !" = 𝑘 ' ! | "! = 0 , as long as the two electrodes are connected in parallel. Indeed, when the electrodes are connected in parallel, the electromechanical coupling 𝑘 ' ! | !! leads to a non-null coupling between the two electrodes with the second mode, which enables the existence of a bi-modal resonance. Interestingly, the maximum achievable harvested power with the two parallel-connected electrodes can be greater than the maximum achievable harvested power with the two electrodes when they are not connected in parallel (blue curve in Fig. 3).

IV. EXPERIMENTATIONS ON A MULTI-MODES MULTI-ELECTRODES HARVESTER

A. Custom-made energy harvester

In order to validate the proposed concepts, a multi-modes multi-electrodes piezoelectric energy harvester has been designed. As shown in Fig. 5, the harvester is made with a cantilever beam of steel covered by four pairs of electrodes of soft piezoelectric ceramics PZT (its fabrication being detailed in [START_REF] Gibus | Design of a vibration energy harvester with two close resonant frequencies[END_REF]). As shown in Fig. 5, the inertial mass of the harvester is inverted in order to bring closer the resonance frequencies of the first and second modes.

B. Experimental set-up

The designed prototype has been fixed on an electromagnetic shaker. As shown in Fig. 6, the acceleration of the electromagnetic shaker is sensed thanks to an accelerometer, whose output is sent to a dSpace control board. The input of the power amplifier (that drives the electromagnetic shaker) is also sent from the same dSpace board and the amplitude of this signal is adjusted internally in dSpace thanks to a PI controller, in order to always maintain the acceleration amplitude to 0.1m.s -2 . The displacement and speed of the tip of the harvester are sensed thanks to a vibrometer laser, and the voltage across the electrodes is sensed thanks to a unity-gain amplifier whose output is sent to dSpace. The electrodes are also connected to a programmable resistance decade box whose value is directly sent from the dSpace board.

C. Experimental characterization of the prototype

In order to characterize each mode from each piezoelectric port, each pair of piezoelectric electrodes has been successively connected to 20 resistive loads while all the other electrodes were short-circuited, and the PEH has been excited with 60 vibration frequencies. The measured displacement as well as the measured power have been used to find the characteristics of each mode seen from each piezoelectric port thanks to a fitting algorithm on Matlab. From these characterizations, the following parameters of the multi-modes multi-electrodes PEH have been identified: first mode, while the fourth pair of electrodes exhibits the largest coupling with the second mode, which is coherent with the literature on the subject [START_REF] Gibus | Strongly coupled piezoelectric cantilevers for broadband vibration energy harvesting[END_REF]. Because of the inverted mass of the harvester, the resonance frequencies of the first and second modes are relatively close ( 𝜔 " = 217.4 𝑟𝑎𝑑. 𝑠 4" , 𝜔 ! = 332 𝑟𝑎𝑑. 𝑠 4" ), and the power limits of both modes are of the same order of magnitude.

D. Experimental results

Figure 7 shows the maximum power with each mode (with the optimal resistive load maximizing the harvested power) as a function of the vibration frequency. As shown in Fig. 7, the power limit of the first mode 𝑃 56'

" is only obtained with the first pair of electrodes, because its coupling is large enough and respects the condition 𝑘 ' ! | "" > 2/𝑄 " [START_REF] Morel | Resistive and reactive loads' influences on highly coupled piezoelectric generators for wideband vibrations energy harvesting[END_REF]. Similarly, the only pair of electrodes that allows to reach the power limit of the second mode (𝑃 56' ! ) is the fourth one, because 𝑘 ' ! | != > 2/𝑄 ! [START_REF] Morel | Resistive and reactive loads' influences on highly coupled piezoelectric generators for wideband vibrations energy harvesting[END_REF].

Figure 7: a. Maximum harvested power (on a resistive load of optimal value) and b. value of the optimal resistive load maximizing the power for each pair of electrodes, with 𝛾 # = 0.1𝑚. 𝑠 %" .

Figure 7 shows that the proposed multi-modes multielectrodes harvester allows to harvest a relatively large power, up to 250 µW (resp. 150µW) for the first (resp. second) mode resonance. Another advantage of the proposed harvester can be seen from Fig. 7.b: for a given vibration frequency, the optimal resistive load maximizing the harvested power depends on which pair of electrodes is connected to the electrical circuit. Therefore, it is also possible to optimize the harvested power by adjusting the configurations of the connected electrodes to the electrical circuit, as demonstrated in [START_REF] Du | Piezoelectric vibration energy harvesting: A connection configuration schema to increase operational range and output power[END_REF].

When the four pair of electrodes are connected in parallel, the equivalent electromechanical coupling with the first mode 𝑘 ' ! | " and second mode 𝑘 ' ! | ! can be predicted using [START_REF] Pereira | Operating frequency prediction of piezoelectric dc-dc converters[END_REF].

𝑘 ' ! | % = w∑ x 𝑘 ' ! | %' 𝐶 * ' ( ')" y ! ∑ 𝐶 * ' ( ')" (16)
The coupling of the first mode 𝑘 ' ! | " is equal to 0.078, while the coupling of the second mode 𝑘 ' ! | ! is equal to 0.064. This means that both modes respect the condition 𝑘 ' ! | 6 𝑄 6 > 2 and that the harvested power with an optimal resistive load should reach both 𝑃 56' " and 𝑃 56' ! [START_REF] Morel | Resistive and reactive loads' influences on highly coupled piezoelectric generators for wideband vibrations energy harvesting[END_REF]. Figure 8 shows the harvested power in this configuration (four pairs of electrodes connected in parallel) as a function of the vibration frequency. Figure 8 confirms that both power limit 𝑃 56' " (around 𝜔 " ) and 𝑃 56' ! (around 𝜔 ! ) can (almost) be reached. Figure 9 shows the maximum harvested power (on an optimal resistance) around the first mode resonance frequency, with more resistance values (20 in Fig. 8, 60 in Fig. 9) and a smaller frequency step (0.5 Hz in Fig. 8, 0.2 Hz in Fig. 9) 9). Around 35.5 Hz, the harvested power goes beyond the traditional power limit 𝑃 56' " and grows to 294µW (20% more than 𝑃 56' " ). This power peak is due to the existence of a bi-modal resonance, and is wellpredicted by the model proposed in section III of this paper ( 14)-(15) (Fig. 3 and Fig. 4). This phenomenon shows that the influence of higher resonance modes can be positive, even though these influences are often seen as negative in the piezoelectric-based conversion [START_REF] Pereira | Operating frequency prediction of piezoelectric dc-dc converters[END_REF] and piezoelectric energy harvesting literatures [START_REF] Morel | A unified N-SECE strategy for highly coupled piezoelectric energy scavengers[END_REF] [START_REF] Brenes | Maximum power point of piezoelectric energy harvesters: a review of optimality condition for electrical tuning[END_REF]. Such phenomenon allows to expand the power limit of multi-modes structures with littleto-no increase of the harvester's volume, which may lead to the design of harvesters exhibiting higher power densities.

V. CONCLUSION

In this paper, we propose a model for multi-modes multielectrodes piezoelectric energy harvesters. Based on this model, we prove that the harvested power can be greater than the traditional power limit thanks to multi-modes resonance. Our model as well as the multi-modes resonance phenomenon are validated experimentally, thanks to a 4electrodes 2-modes piezoelectric energy harvester based on an inverted mass. Such model opens new ways to design strong power-density and large bandwidth piezoelectric energy harvesters thanks to multi-modes resonances and to clever electrode placement/choice. Experimental results

Numerical results

! !"# $

+20%

Bi-modal resonance

Figure 1 :

 1 Figure 1: Generic illustration of a multi-modes structure with four piezoelectric electrodes.

Figure 2 :

 2 Figure 2: Equivalent circuit model of a multi-modes multi-electrodes energy harvester

Figure 3 :

 3 Figure 3: Maximum harvested power (on optimized resistances 𝑅 ! and 𝑅 " ) with 𝑘 # " | !! = 𝑘 # " | "" = 0.2 and 𝑘 # " | !" = 0 . The blue curve shows the maximum power with no coupling between the second mode and first electrode (𝑘 # " | "! = 0) while the yellow curve shows the maximum power with 𝑘 # " | "! = 0.1.

Figure 4 :

 4 Figure 4: Maximum harvested power (on optimized resistance 𝑅 ) with 𝑘 # " | !! = 𝑘 # " | "" = 0.2 and 𝑘 # " | !" = 𝑘 # " | "! = 0 , when the piezoelectric electrodes are connected in parallel.

Figure 5 :

 5 Figure 5: Experimental realization of a multi-modes multi-electrodes energy harvester. Dimensions are indicated in millimeters.

Figure 6 :

 6 Figure 6: Experimental setup for the test of the multi-modes multi-electrodes energy harvester.

Figure 8 :

 8 Figure 8: Harvested power with the four pairs of electrodes connected in parallel as a function of the vibration frequency.

Figure 9 :

 9 Figure9: Maximum harvested power (on a resistive load of optimal value) with the four pairs of electrodes in parallel, around the first resonance mode, illustrating the bi-modal resonance phenomenon.

Figure 9

 9 Figure9reveals a particularly interesting phenomenon that was predicted by equation[START_REF] Ferrari | Piezoelectric multifrequency energy converter for power harvesting in autonomous microsystems[END_REF]. Around 35.5 Hz, the harvested power goes beyond the traditional power limit 𝑃 56'

  𝐶 * ' ) corresponds to the expedient electromechanical coupling between the 𝑛 #$ resonance mode and the 𝑚 #$ pair of piezoelectric electrodes. Because of the low-coupling of some combinations of electrode/mode (e.g., second electrode with second mode), the low identified values of 𝛼 and 𝑘 ' ! might be approximative. We can observe that the first pair of electrodes exhibits the largest coupling with the

						PEH
					Accelerometer
					Electromagnetic	Vibration "	_	+
					shaker	Power	Voltage
						amplifier	follower	Programmable
					Acceleration	! !"# $	% #	#	resistance decade box
					Link
	Normalized power P/P 1 lim [-]	0 0.2 0.4 0.6 0.8 1 1.2	Bi-modal resonance due to ! ! " | "" ! !"# % ! !"# $	!	Matlab Simulink	Control board (dSPACE)
	TABLE I. Parameters ELECTRODES PIEZOELECTRIC ENERGY HARVESTER IDENTIFIED PARAMETERS OF THE MULTI-MODES MULTI-Values Units 𝒎 = 𝟏 𝒎 = 𝟐 𝒎 = 𝟑 𝒎 = 𝟒 𝛼 !# 2.3 1.3 0.54 0.070 mN.m -1 𝛼 "# -2.9 -2.4 -7 -12 mN.m -1 𝐶 $ # 39 41 34 31 nF 𝑘 # " | !# 0.087 0.026 0.0055 1.0 × 10 %& -𝑘 # " | "# 0.0034 0.0022 0.023 0.073 -𝑀 ! 33 g 𝑀 " 580 g 𝐵 ' ! 178 g 𝐵 ' " 735 g 𝑄 ! 45 -𝑄 " 40 -𝜔 ! 217.4 rad.s -1 𝜔 " 332 rad.s -1 𝑃 ()# ! 245 µW 𝑃 ()# " 144 µW 𝑘 ' ! | %' = 𝛼 %' ! /(𝐾 % 15 20 25 30 35 Frequency [Hz] Inverted mass in steel Cantilever beam in steel Electrodes 1Electrodes 2 Electrodes 3 Electrodes 4 60 5 5 0.4 65

ACKNOWLEDGMENT A.M. thanks M.-Y. Durand for the fruitful discussions on dynamics and resonances.

This project has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No 862289.


Alfredo Gonzalez-Perez
email: alfredo.gonzalez@campusroslagen.se

UV-LED water disinfection for small-scale decentralized systems

Water disinfection can be achieved with different techniques using chlorine, ozone, ultraviolet light or filtration. Among the different processes ultraviolet has become a technology of choice due to their proven efficiency inactivating protozoa, bacteria, and viruses. Ultraviolet light involves high process safety, ensuring that no harmful byproduct will be formed during the application. The main germicidal UV sources are low-pressure mercury lamps, Excimer lamps, UV-LEDs, and other like medium-pressure mercury lamps and Xenon lamps that allow for pulsed operation (Bergman 2021). Those are used for water, air and surface disinfection among other uses.

In recent years, a global treaty to protect human health and the environment from the adverse effects of mercury has emerged. With the implementation of the regulation by the Minamata Convention on Mercury that enforce the elimination of mercury-based technologies and their substitution for alternative ones. Currently water purification plants all over the world depend on mercury-based UV lamps for their disinfection units (allowed by the Minamata convention on Mercury for the water business), which maybe required for replacement when alternative technologies like those on UV-LED become efficient enough for the treatment of large water volumes and regulation incorporate guidance adapted to their operating requirements.

Technology replacement Germicidal Activity Wavelenghts

Motivation Acknowledgement

Centralized water treatment systems can treat large volumes of water at high rates to accommodate all residential, business, and industrial uses and are mostly used in urban areas. A centralized approach has multiple advantages; effectively remove practically any range of raw water turbidity along with harmful pathogens, including bacteria, virus, and protozoa. Some disadvantages like the capital cost as well as operating and maintenance costs for a centralized system can be significant. In addition, the fast growth of many cities makes a centralized approach unsuitable leading to some suggestions for a shift towards decentralization have been proposed. This approach opens for the possibility to treat and reuse water at the local level as well as testing different technologies without affecting the entire city providing a more sustainable approach with increasing resilience.

In rural areas as well as remote areas away from traditional municipality services a small to medium size decentralized water supply and treatment is preferred. A decentralized system that could be point-of-entry (POE) or point-of-use (POU) a treatment unit can be installed at the individual home, small community, or business to provide drinking water and wastewater treatment. The advantage of a low initial capital investments and reduced operating and maintenance costs is a plus. In general, the shift towards decentralization have advantages for both small-medium size water treatment as well as large scale water treatment at the economic, operational and technology level adding the possibility to use emerging technologies like the UV-LED based one, and potentially more efficient, that are not ready for large scale utilities.

Decentralized water treatment Key References

(

UV-LED Challenges and opportunities

•

  UV LEDs of higher output are needed and more efficiency. • Studies indicate that combined wavelengths of 280/365nm and 280/405nm are the best for disinfection. Use of dual-wavelenghts for efficent disinfection and flashing modes. • Reactor designs with new geometries. Some excellent designs away from traditional mercury-lamps are already commertially avaliable. • New regulations need to be develop taking into account the characteristics of UV-LED technology to facilitate the use on municipality owned water treatment plants. Low-Pressure (low-output): Small/Compact systems low doses Low-Pressure (high-output): Small/Compact systems greater doses Medium-Pressure (high-output): High doses and high flow volumesThe author wants to acknowledge the Water Centre for Innovation (UCV -Utvecklingscentrum för vatten) at Campus Roslagen AB (Norrtälje Municipality) for support and discussions on water decentralization and regulation.

  1) Lui, GY., et al. 2014 Point-of-use water disinfection using ultraviolet and visible light-emitting diodes. Science of the total environment 493, 185-196. (2) Bergman, RS. 2021 Germicidal UV sources and systems. Photochem and Photobiol. 97, 466-470. (3) Song, K., et al. 2016 Application of ultraviolet light-emitting diodes (UV-LEDs) for water disinfection: A review. Water Research 94, 341-349. (4) Chevremont, A.-C.C., et al. 2012. Multivariate optimization of fecal bioindicator inactivation by coupling UV-A and UV-C LEDs. Desalination 285, 219-225. (5) Chevremont, A.C., et al. 2012. Effect of coupled UV-A and UV-C LEDs on both microbiological and chemical pollution of urban wastewaters. Sci. Total Environ. 426, 304-310. (6) Song, K., et al. 2019. Microorganisms inactivation by wavelength combinations of ultraviolet light-emitting diodes (UV-LEDs). Sci. Total Environ. 665, 1103-1110. (7) Woo, H., et al. 2019. Efficacy of Inactivation of Human Enteroviruses by Dual-Wavelength Germicidal Ultraviolet (UV-C) Light Emitting Diodes (LEDs). Water 11, 1131. (8) Xiao, Y., et al. 2018. Impact of UVA pre-radiation on UVC disinfection performance: Inactivation, repair and mechanism study. Water Res. 141, 279-288. (9) Matafonova, G., and Batoev, V. Dual-wavelength light radiation for synergistic water disinfection. Sci Total Environ. 806 (2022) 151233. Small and medium-size water treatment units can already benefit form the use of UV-LEDs lamps for their water disifection units. This alternative is optimal for decentralized water treatment units (Lui et al 2014). UV-LEDs provide with new advantages not present in traditional mercury-based UV lamps. The use of several individual waveleghts, flashing modes and new reactor designs are some of the advantages that will drive the transition to UV-LED standarization on water treatment. (Song et al 2016). The table show a selection of multiple wavelengths combinations (UV-C, UV-B and UV-A) targeting different microorganisms on different mediums. A recent review by Matafonova and Batoev (Matafonova et al. 2022) summarize the state of the art on disinfection using dual-wavelength light radiation. Below we show some key parameters for uv disinfection on using UV-mercury lamps.


Alfredo Gonzalez-Perez
email: alfredo.gonzalez@campusroslagen.se

Angel Macia

Water disinfection with micelle-clay composites

Water disinfection can be achieved with different techniques using chlorine, ozone, ultraviolet light or filtration. Among the different processes filtration is a good option to avoid the formation disinfection by-products. The use of micelle-clay composites has been proven to be efficient for the removal of bacteria, viruses and protozoan.

A small layer of modified clay was used for a preliminary filtration test, based on the micelle-clay composite, with a known solution of E. Coli. The bacteria were counted after filtration to assess the entrapment on the test filter. Preliminary results on filtration in in the laboratory with the Micelle-clay composite (sphere-like micelles) indicate the removal of bacteria in good agreement with the results previously published using other quaternary ammonium surfactants in their composition. Note that hexadecylpydinium cation is a close structural relative to octadecyltrimethylammonium and hexadecyldimethylbenzylammonium that had been used by other authors [START_REF] Benitez | Modified compositions of micelle-clay and liposome-clay composites for optimal removal from water of bacteria and hydrophobic neutral chemicals[END_REF] for the removal of bacteria from water. More studies will follow these preliminary results to investigate the optimal composition for the maximization of the bacterial removal and different flow rates as well as test on the regeneration of the used filter material to increase their life-spam.

Preliminary results

Motivation

Acknowledgement

Summary: We prepared new organoclay and micelle-clay composites using bentonite and hexadecylpyridinium bromide (C16PyBr) for bacterial removal from water. Using different surfactant concentrations, we ensure that monomeric, spherical micelles and rod-like micelles are present for the clay modification. The modification process was monitored for several days to ensure the saturation and the modified clay was extracted via centrifugation and posterior drying. Preliminary tests for bacterial removal agree with [START_REF] Benitez | Modified compositions of micelle-clay and liposome-clay composites for optimal removal from water of bacteria and hydrophobic neutral chemicals[END_REF] and indicate that bacteria can be remove from a stock solution opening the door for further systematic experiments on bacterial removal and filter recovery. The method can be used for small-scale drinking water production units and have potential uses in general for small to medium scale decentralized systems. Methods:

The organoclay and micelle-clay composites where prepared using the protocol described by [START_REF] Nir | Water purification from perchlorate by micelle-clay complex: laboratory and pilot experiments[END_REF]. using the nanoclay (hydrophilic bentonite), montmorillonite K30 and Bentonite MED powder separately and C16PyBr as a cation in monomeric and micellar form (spherical and rod-like). The concentration range for C16PyBr was chosen following the data from Gonzalez-Perez et al 2006.

  This work has been carried out within the strategic innovation program BioInnovation, a joint venture by VINNOVA, FORMAS and the Swedish Energy Agency (Vinnova project 2022-00677 Biobased composite filters for small scale water purification units). A G-P acknowledges financial support from VINNOVA through the project reference 2021-03838 Removal of Macrolide Antibiotics from Water using Clays.

Figure 1 :

 1 Figure 1: Logarithm of the critical micellization concentration in molar fraction units against the number of carbon atoms in the amphiphile chain of n-alkylpyridinium halides. Solid squares and circles correspond to the measured first and second cmc's, respectively.

Figure 1 .

 1 Temporal following of the conductivity is shown on Figure 2.

Figure 2 :

 2 Figure 2: Conductivity vs time at different C16PyBr-Motmorillonite concentrations in monomeric and micellar form.


Fabrice Kwiatkowski

Laurent Serlet

Andrzej Stos

Oncogenetic pedigrees: relation between design and ability to predict mutation

Keywords: oncogenetics, risk modelisation, family pedigree, simulation MSC Classification 62P10, 92D10, 92C50, 92-10, 92-08

We consider the risk of a disease caused by the presence within the genome of one deleterious mutation or two interacting mutations. For an individual, the probability of being mutated/doubly mutated can be estimated knowing the phenotype of the other members in the family pedigree. We study the performance of this process as a function of the size, the shape of the family tree and the parameters of the model. We carry out simulations using the parameters pertaining to breast/ovarian cancer in BRCA-mutated families.

1 Introduction :

Diseases caused by genetic mutations are numerous. One of the most representative cases is likely breast/ovarian cancer favored by BRCA mutations. But other mutations are suspected for many cancers as well as in other types of diseases like cardiovascular ones. Neurodegenerative diseases may soon join the list of pathologies associated to mutations. Even if sequencing the genome -or, as far as we are concerned, part of it -has become faster and much less expensive, it cannot yet be systematically performed for various reasons and its analysis requires to have fully identified the mutations involved. Therefore we still need reliable computer tools first to evaluate a familial predisposition considering its medical data, i.e. the pedigree, second to limit the genetic inquiry to a reduced number of susceptibility genes and third, to estimate the individual risk to inherit a deleterious genetic mutation.

Scoring systems (BRACAPRO [START_REF] Berry | Validation, sensitivity of genetic testing of BRCA1/BRCA2, and prevalence of other breast cancer susceptibility genes[END_REF], Manchester [START_REF] Evans | Addition of pathology and biomarker information significantly improves the performance of the Manchester scoring system for BRCA1 and BRCA2 testing[END_REF], Eisinger [START_REF] Eisinger | Identification et prise en charge des prédispositions héréditaires aux cancers du sein et de l'ovaire (mise à jour 2004)[END_REF]…) have been built in the context of HBOC (Hereditary Breast/Ovarian Cancer) and various refinements have been proposed [START_REF] Bonaiti | A new scoring system for the diagnosis of BRCA1/2 associated breast-ovarian cancer predisposition[END_REF]. These scoring methods only require a few additions, so their implementation is easy. However, the amount of family information they use is limited and their performance is not optimal. Nowadays, family history is becoming better and better reported, going back several generations and stored in databases. Relatively complex calculations are therefore easy to perform. One can in particular calculate the probability of mutation knowing all the family information which is the best predictor of the mutation risk. This method has long been described by [START_REF] Elston | A general model for the genetic analysis of pedigree data[END_REF]. The prognosis can of course be confirmed or reversed after sequencing, in the case of mutations that are clearly identified. And in all cases, this risk calculation is a valuable tool in the personalization of screening and medical monitoring.

In this context of an ever-increasing availability of data and of calculation possibilities, the aim of this article is to optimize the design of this method and to examine its performance using simulations.

A few words about the underlying model first. Since the prognosis method we are discussing is based only on probability calculations, it does not require a detailed knowledge of the nature of the mutation concerned (or the mutations concerned). We will assume, as in most of the known cases, that the deleterious effect of the mutation consists of a strong increase in penetrance i.e. the probability of declaring the disease at some point in one's life and an earlier age of onset. These effects are evidently seen in the case of cancer, but we aim in this work to establish a general mathematical framework that is not restricted to a specific disease. In what follows, the disease caused by the mutations will be called and noted K.

Whether the antecedents are collected by the oncogeneticist especially for this prognosis purpose or whether they are pre-existing in a database, it is necessary to question the desirable period of this family exploration, that is the size and shape of the family tree to take into account, provided of course that such a tree can be completed in relation to the individual concerned. Several issues naturally arise:

-Is the prognosis always better when a larger number of parents are included? -Is it desirable to go back the generations as much as possible? -Does taking siblings into account in the ancestral line improve the result? Are cousins useful? -Is the prognosis easier in families with several cases of illness? -Do descendants provide information as valuable as ascendants to evaluate the mutation risk? -What sensitivity and specificity can be typically achieved by this kind of method and what is their variability according to the parameters of the model?

We will answer these questions, with answers that are sometimes surprising, in Section 3. Before that, we present the modelling of the data and the principle of the computational algorithm in the following Section 2. We end this paper in Section 4 with a discussion about the limitations and possible extensions of our work.

Model and data description

The basis for constructing a random model is that individuals from the same family have genotypes stochastically linked to each other by Mendel's laws (50% chance of inheriting a mutation from one of its parents). But these genotypes are hidden and all statistically usable data (phenotypes) are random expressions of these genotypes. Mathematically, this is a hidden Markov process indexed by a tree. Two difficulties then arise:

• The combinatorial explosion in the number of possible genotypes makes unrealistic -except for small families -any approach based on the exhaustive examination of all possible genotypes, i.e. a brute force method. • Because the considered mutations are relatively rare in frequency in the general population, the effect to be analyzed remains narrow and it must be distinguished from the random noise which here takes the form of sporadic onset of the disease, generally at an older age.

In the following numerical calculations, the probability of the disease occurrence (as a function of the age), depending on whether a mutation is involved or not, will have the same distribution as what is observed for breast/ovarian cancers with or without a BRCA mutation: In addition, we only consider the case of heterozygous mutations, knowing that homozygous mutations are generally lethal during embryogenesis. We therefore consider that the genome of each individual is a random quantity having two possible states {wild-type, mutated} with the following associated probabilities: We assume that all the parameters of the model are known:

• the mutation frequency in the general population, 𝑝 "

• the penetrance of K for non-mutated subjects, 𝑝 # • the penetrance of K for the mutated subjects, 𝑝 $ However, these parameters can be varied in order to consider penetrances linked to less deleterious mutations or even interactions between mutations and polymorphisms that are very poorly penetrating alone. The methods for evaluating these different contexts will be discussed in Section 4.

Here, we assume that the data available per individual is the year of birth, whether or not the disease occurred, and if so, the age at reporting. Such an approach, although it appears succinct, has the advantage of being easily generalizable. To quantify the performance of the algorithms used, we are going to perform data simulations according to the following model: for each family pedigree, we will first generate the (hidden) genotypes: for any individual, his genotype is obtained randomly but according to that of his parents if they are present; if they do not exist in the family pedigree, we will use the frequency 𝑝 " . We then randomly assign birth years to the family members so that ages are consistent with the generations. Finally, the phenotypes are generated by taking into account the genotypes, the value of the parameters p0 and p1 and the laws described at Figure 1.

Principle of the mutation risk computation

The calculation of the risk of mutation considering the disease cases in the family by the conditional probability calculation, as reported below, has been described long ago [START_REF] Elston | A general model for the genetic analysis of pedigree data[END_REF]. However, its implementation for each specific case is complex. For this reason, practitioners rather use scoring methods which only require making a few additions or subtractions but without being able to claim the same performance [Bonaïti, 2014].

Let us denote by F the set of phenotypes of all family members -known family data -and denote by G the compilation of the (unknown) genotypes of individuals, in particular G(i) will denote the genotype of the i th individual. For the single mutation model which we are working with until stated otherwise, we will note 𝐺(𝑖) = 1 when individual i is mutated and 𝐺(𝑖) = 0 otherwise.

An elementary probability calculation from Bayes' theorem shows that:

In each of these two expressions, the sum is extended to all possible genotypes satisfying the mentioned constraint on the genotype 𝑔(𝑖) of individual i. Explicit expressions for the conditional probability of the phenotype F knowing that the genotype is g and for the probability that 𝐺 = 𝑔 are given in Section 2.4.3.

Note that the calculation of the numerator and the denominator require either scanning all possible genotypes or proceeding by a Monte-Carlo method. The first method, called "by brute-force", is possible for small families and requires the calculation of the likelihood of each genotype g. Note, moreover, that a fraction of genotypes are impossible, therefore of zero probability; the summation is de facto reduced to eligible genotypes. The Monte-Carlo method consists in simulating a series of genotypes, either by direct process or by a Metropolis algorithm [START_REF] Metropolis | Equation of state calculations by fast computing machines[END_REF]. With this method, formula (1) is rewritten:

where G1, G2, ... forms a series of random draws among the set of possible genomes, according to Mendel's laws and the contextual hypotheses mentioned above. The 1A function refers to the binary indicator function which is equal to 1 if A is carried out and 0 otherwise. Numerically, it is necessary to ensure that there is a sufficient number N of simulations so that both the numerator and especially the denominator are estimated correctly.

The (conditional) probability of mutation obtained by formula (1) is a quantification of the risk of mutation of an individual but also indirectly of the family risk. One may prefer a "mutated/non-mutated" binary predictor, which will be obtained by thresholding: for a threshold level s, the individual i is prognosticated as mutated if 𝑃( 𝐺(𝑖) = 1 | 𝐹) > 𝑠. From a certain point of view, one may regret this conversion in binary prognosis of a quantitative risk which is a priori more informative. However, this conversion then makes it possible to analyze the predictor with the usual performance indicators, i.e. in terms of sensitivity (percentage of well predicted among the mutated) and specificity (percentage of well predicted among the non-mutated). We will present the traditional Receiver Operating Characteristics (ROC) curve and in the rest of this chapter, the different predictors will be compared using the relative position of their ROC curves and the corresponding areas under curve (AUC). We can also consider evaluating the performance on the whole family tree, which amounts to averaging the predictors on all the individuals of the pedigree. The results are generally better, but less pertinent so we will leave them aside to focus on the individual case.

Standardization of family trees

One of our main purposes is to address the influence of the size and shape of the family tree considered. In a real situation, the family tree is constrained in its form by the number of children born of the different couples and by the ability to obtain reliable information about the phenotype of each individual. This last point applies particularly to older generations, because they have had time to present the disease. As a result, the size and especially the shape can vary a lot.

To facilitate comparisons, we are going to use a generic form of regular pedigree, with two parameters which are the size (number of children per couple) and the height (number of generations), and to answer specific questions later, we will define other typical shapes. We call regular pedigree with ng ³ 2 generations and ne ³ 1 children per couple and we denote by Reg(ng, ne) a pedigree whose skeleton is the ancestral line of an individual numbered 0 which goes back to ng generations including the ancestors of two genders and to this skeleton, we add ne children for each couple. Here is a diagram of a Reg(4, 3) type pedigree: For comparisons, we can fix ng and vary ne ³ 1, or fix ne ³ 1 and vary ng ³ 2. An intuitive idea --maybe naïve --is that the more individuals are taken into account (richness of phenotype), the more the performance of the predictor increases, probably up to a certain limit. A related but more subtle question is to assess the importance of fixed-size tree shapes: is it better to have more generations but fewer individuals per generation than the opposite? In particular, we will compare the performance of the predictors for two pedigrees both comprising 15 individuals. The first is Reg(4, 1) which we will call "Grt Grd Par" because it involves all the great-grandparents. Here it is: The second that we call "Wide 3 gen." also includes 15 members but only 3 generations because of the addition of uncles / aunts and cousins. In the case of the regular tree Reg(ng, ne) with ng ³ 3 and ne ³ 2, individual 0 has uncles (or aunts) but his cousins do not appear. In general, it might be interesting to isolate the influence of individuals a little further from the direct ascending line of individual 0. We will compare the family "Uncles", a family with 11 individuals, to the family "Cousins" with 19 individuals. This latter strictly contains the family "Uncles" but 8 cousins are added, drawn in green in the diagram below: Studying the influence of the shape, we also study the case of an asymmetric family pedigree in the sense that only one of the parents of individual 0 is entered. For example, we compare the predictors obtained for the following two families of almost identical sizes: In the previous diagrams, individual 0, for which the mutation risk is evaluated, is positioned on the right and his risk is calculated according to individuals of the same or previous generations represented at his left.

One might question to what extent descendants can also provide valid information about an individual's genotype. Of course, this can be informative only if an individual is old enough to have children and these children should not be too young themselves. This may limit its practical value, but this problem remains an interesting question. In particular, for a fixed pedigree size, do ascendants give more information than descendants?

Recall that to obtain our numerical results we simulate families with inevitably a choice on the law of the birth year by generation. Here, the years of birth are chosen so that the (potential) age of the individuals is uniform in [30,50], [50,70] and [70,90] respectively for the first 3 generations. We use the word "potential" because the individual may in fact have died earlier, either from the disease K studied or from another cause. It is essential to take age into account because the likelihood of having contracted the disease K before a given date depends on how long the individual has lived. Consequently, "young" individuals, for example of the first generation, a priori provide little information, which suggests that prediction with ascendants is more efficient than that with descendants. To carry out this comparison, we will compare the predictors associated with pedigrees "Grd Par." = Reg(3, 1) including 7 individuals and "Grt Grd Par." = Reg(4, 1) containing 15 individuals respectively to those calculated with the following pedigrees:

The importance of conditioning

Nowadays, the individual mutation diagnosis is primarily offered to families considered at risk and a family is considered as such if one or even several cases of sickness K have occurred. This induces a strong bias on data in present databases such as the ones collected by anti-cancer centers. In particular, estimating some parameters in these databases will lead to values radically different from the ones for the general population. In the future, larger and less biased databases will probably emerge and risk calculation will be computed systematically as a prevention approach. Consequently, in the following numerical experiments we will consider both the case of unconstrained family and the case of families conditioned to have (at least) a certain number of cases of illness. The performances and ranking of the predictors can change significantly between the two cases.

Notation and formulas used in software

To ensure reproducibility of our work, our demonstration software (in Python) can be downloaded at the address: https://drive.uca.fr/f/a07b89c2df52460686bf/?dl=1

With minor adaptations this software can execute most of the experiments that follow.

To complete the mathematical description of the risk computation algorithm and make the software clear, we give in the present subsection the necessary notation and the explicit formulas for the laws (likelihood) of the random objects involved.

Coding the family tree structure

Mathematically a family tree is a set of 𝑛 (related) individuals, numbered from 1 to 𝑛, for whom data have been collected. The genealogical relations between these individuals are modeled -both in the formulas below and in the companion software-by two vectors (𝑓(𝑖), 1 ≤ 𝑖 ≤ 𝑛) and (𝑚(𝑖), 1 ≤ 𝑖 ≤ 𝑛) where, for 1 ≤ 𝑖 ≤ 𝑛, the value 𝑓(𝑖) [resp. 𝑚(𝑖)] indicates the index of the father [resp. mother] of the individual of index 𝑖, when these parents belong to the family tree i.e. the set of individuals for whom data are available. Otherwise, we set 𝑓(𝑖) = -1 if the father of the individual 𝑖 is out of the tree and we write 𝑖 ∈ 𝑈 < . Similarly, 𝑖 ∈ 𝑈 " if 𝑚(𝑖) = -1, meaning the mother of 𝑖 is out of the tree.

Genotype: law and generation

For a given family tree, the genotype 𝐺 = (𝐺(𝑖), 1 ≤ 𝑖 ≤ 𝑛) is the collection of the genotypes of all individuals in the tree. For 1 ≤ 𝑖 ≤ 𝑛, we have 𝐺(𝑖) = 1 or 0 according to the individual being mutated or wild-type, respectively. The only parameter useful in the generation and the expression of the law of the genotype is 𝑝 " , the probability for an individual drawn at random in the population to be mutated. This is used for the genotype of individuals whose parents are out of the tree. For the other individuals, the probability of being mutated is expressed conditionally to the state of their parents, as prescribed by Mendel laws. To this end, we introduce the 3 × 2 matrix which gives the probability for an individual to be wild-type (column 0) or mutated (column 1) according to the number of mutated parents: 0 for line 0, 1 for line 1, 2 for line 2. Recall Table 1 for justifications. Taking into account the previous remarks, the law of the family genotype 𝐺 is given by: Note that this probability can be zero since the matrix 𝑀 contains a null term corresponding to the probability for an individual to be mutated when neither of his parents are.

Phenotype: law and generation

The family phenotype is the collection of the phenotypes of the 𝑛 individuals under consideration and, in the present work, these phenotypes are reduced to the occurrence or not of a cancer and, in the positive case, the age of onset (which in reality is more often the age of diagnosis). The phenotype of one individual is thus the age of diagnosis between 0 and 99 or the value 200 meaning by convention that no cancer has occurred.

We use the two parameters introduced earlier: 𝑝 $ [resp. 𝑝 # ] is the probability that a mutated [resp. wild-type] individual develops a cancer throughout his life. Given that an individual does develop a cancer, the age of diagnosis is given by one of two probability laws (𝑑 B (𝑡), 0 ≤ 𝑡 ≤ 99) depending on whether the individual is mutated 𝑔 = 1 or not 𝑔 = 0. The curves we use in practice for our numerical studies are drawn in Figure 1. The ages (𝑎(𝑖), 1 ≤ 𝑖 ≤ 𝑛) of the family members are randomly generated throughout the tree in a coherent manner as specified before. For old generations, age has to be understood as potential age and is set to 99.

With the previous notations, the conditional law of the phenotype 𝐹 knowing the genotype 𝐺 = 𝑔 can be expressed as M = 0 @ 1 0 0.5 0.5 0.25 0.75

1 A = (M [i, j], 0  i  2, 0  j  1) < l a t e x i t s h a 1 _ b a s e 6 4 = " + Y J B 5 O h Q P B Z c B r P E + B c u E 2 g E M R 4 = " > A A A D S H i c j V F d S x w x F L 0 z 2 l b t h 6 u + K L 5 c u r S s I E t m W 7 E g g t Q X X w Q F V 4 W d R T L Z 7 B q d n Z l m M o V l 2 a f + u / 4 B 8 R / U t 7 6 V v v U m G 6 F V S p t h k n P P v e c k N 0 m K V J W G s d s g n J l 9 8 v T Z 3 P z C 8 x c v X y 3 W l p Z P y 7 z S Q r Z F n u b 6 P O G l T F U m 2 0 a Z V J 4 X W v J h k s q z 5 H r f 5 s 8 + S 1 2 q P D s x o 0 J 2 h 3 y Q q b 4 S 3 B B 1 U f t y u B u n s m 8 a G C d y o L I x 1 5 q P J m M h J h j h W 2 Q Y x 8 i a W x b S 7 I L W N N r e Q o x l 1 v M S j L U a X J q N X W w c d t T m V R c 3 M d 5 B R v a f U K F b W p Z i U 3 w 1 X a I N v K j V W Z O 5 g Y 9 B 5 E E d / D j K a z c Q Q w 9 y E F D B E C R k Y A i n w K G k r w M R M C i I 6 8 K Y O E 1 I u b y E C S y Q t q I q S R W c 2 G u a B x R 1 P J t R b D 1 L p x a 0 S 0 q / J i X C G 9 L k V K c J 2 9 3 Q 5 S v

Results

We now present results of numerical experimentations which address the questions stated in the introduction. For most of them a choice of parameter values had to be made. Unless stated otherwise, we chose 𝑝 " = 0.03, 𝑝 # = 0.02, 𝑝 $ = 0.7. Although this appears arbitrary, we have carried out enough experimentations to be convinced that the phenomena we emphasize hold for other "typical'' values of the parameters. Moreover, the available software enables the reader to check by himself for any set of values.

Reliability of the inference: confidence intervals

In what follows we will compare different predictors on the basis of their ROC curves and in particular the area under these curves (AUC). These ROC curves are obtained after hours of intensive computation but are still submitted to the randomness of the samples drawn to construct them. In other words, the AUC is affected by a random noise which has two distinct origins: first the sampling of families on which the predictors are tested, second the computation of the predictor itself when using simulation to evaluate the conditional probability. Note that this second source of randomness disappears when using a brute force method.

Consequently, it would be safe to obtain precise confidence intervals for all the AUC that we will provide but that seems too ambitious regarding the computation time needed. In particular, the noise on the AUC does not seem to be close to Gaussian with the current sample sizes and real confidence intervals are much wider than Gaussian ones. Consequently we have chosen to evaluate confidence intervals in the specific case of the Reg(3,2) tree which allows rather quick computation both by brute force or simulation method. The size of the confidence intervals in this case will allows us to calibrate the software in terms of numbers of simulations, taking also into account the limits in computing time and it will guide us on the conclusions that can be drawn from the other ROC curves.

Below are the graphs produced by the demonstration software in its original form. It gives the performance of the mutation predictor for a regular family tree with 3 generations and 2 children per couple. It shows the ROC curve i.e. resulting from the sensitivity and specificity of the predictor to assign individual 0 to its correct group (i.e. mutated or not). Since this family contains 10 individuals, there are 2 $# genotypes (including the inadmissible ones) so a brute force method for the computation of the conditional probability is possible. Also, as we have already mentioned, conditioning the family by a certain number of occurrences of K makes sense and we will tell more on this subject in the next subsection. For the moment we consider both the unconditioned case and the case of at least one occurrence of K. Hence the figure below shows four graphs corresponding to the different cases with brute force method represented by continuous lines and simulation method represented by crosses. First, we see that the brute force and the simulation methods give almost exactly the same results. This is a reassuring fact about the accuracy of the simulation method recalling that, for most of the cases we will study, it will be the only method available. We performed 50 000 simulations of the family tree and we imposed, for the simulation method based on Formula (2), that both evaluations of the denominator and the numerator in (2) get a least 5000 contributions. These are the numbers we will use by default for the other ROC curves in this article. For large families, it is close to the limit imposed by the computation time.

In the unconditioned case represented above, the empirical confidence interval for the AUC based on 100 attempts is [0.949, 0.959] for a confidence level of 90%. The conditioned case shows a confidence interval of similar length. This means that ranking predictors can safely be made when AUC show differences typically larger than 0.02. This will often be the case in the following experiments, but not always. For tighter differences on AUCs, we may simply conclude that the corresponding predictors have quite similar performances which is certainly sufficient in practice but if a ranking is nevertheless desired, we will raise the number of simulations or multiply the number of parallel computations.

Influence of conditioning

The performance of a predictor evolves when we condition the total number of ill people in the family. For instance, for a family Reg(4, 2) with 22 individuals, that allows very good performances (as we will see later), we "force" the minimum number of patients successively to 0, 1, 2 or 3 and we obtain the following ROC curves: The first two curves overlap and hardly stand out. We notice that beyond 1 case of disease, conditioning decreases the effectiveness of the predictor, which in the end is quite logical because it becomes more and more difficult for the model to distinguish the mutated families from the others as they all systematically present cases of illness.

Influence of the generation number

We consider regular pedigrees Reg(ng, 1) for a number of ng generations varying from 2 to 6 and a fixed number of children ne per couple equal to 1. The ROC curves for the associated predictors are as follows: Obviously, the data from only 2 generations, i.e. the individual and his parents are not sufficient to obtain a successful prognosis. Rising to 3 generations allows a significant gain and the addition of the 4 th generation still provides a small improvement. Beyond, any addition of information from previous generations is useless or even harmful: performance seems to regress slightly with the 6 th generation. This is confirmed by performing 20 independent attempts: the AUC for 6 generations is less than the AUC for 5 generation (yet infinitesimally) for all attempts. This phenomenon also holds when the families are conditioned to at least one occurrence of the illness and essentially the same thing happens when we set the number of children per couple at 2. Our interpretation is that too distant individuals bring more randomness than information.

Number of children

Investigating the performance according to the number of children in pedigrees is more theoretical than the previous issue. Indeed, in the previous case, the doctor / statistician wishing to make a prognosis can choose the number of generations that he provides. For children, their existence, a fact, cannot be the object of a choice, even if the practitioner can choose to include only one child: the one who is in the ascending line of the individual studied (ex. the index case1 ). The ROC curves obtained for the families Reg(4, ne) where ne varies between 1 and 5 are as follows: We observe that the larger the pedigree, the better the prediction with a significant gain added between one child and two children. Beyond two children the gain is lower. This means that adding information from individuals directly connected to the lineage is beneficial. However, since the AUC increases slightly less as the number of children increases, there may be a limit at which the performance of the predictors peaks.

Height against width

We have just seen that the predictors performances increase with the number of generations up to 4 (height) and also with the number of children per couple (width). Hence the question: for the same size (same number of individuals), is it better to have a large pedigree or a tall one? If we refer to the following ROC curves relating to pedigrees "Grt Grd Par." and "Wide 3 gen." of 15 individuals each, it would seem that a large pedigree does slightly better than a tall one. However, conditioning the presence of at least one disease case per pedigree is sufficient to produce a markedly opposite result: Given the presence of cancer cases, either sporadic or familial, in pedigrees usually observed in oncogenetics, we can estimate that the second case is more relevant. Therefore, with an identical number of members, the height of a pedigree seems more informative than its width, even if, as seen previously, beyond 4 generations compiling information does not bring anymore benefit.

Utility of cousins

We have previously shown the value of the information corresponding to individuals directly connected to the ancestral line of the individual studied. It makes sense now to study the case of somewhat more distant relatives, starting with cousins. We thus compare the predictors associated with the family "uncles" and "cousins" among themselves but also with the pedigree Reg(3, 1) without uncle or cousin: This figure is interesting because it shows that indirect links can also be very informative -this is echoing the conclusion as to the number of children per couple. If the addition of uncles does not improve ROC AUC, the addition of cousins when you have uncles, brings clearly contributing elements. This could be explained by the fact that if a cancer is found in cousins, usually rather young, these cancers are very informative about a possible familial mutation.

Asymmetric pedigrees

We study the extreme case where all phenotypic information on either the father or the mother and all their ancestry is missing. Is it then illusory to hope for a mutation prognosis? We compare the associated predictors of the following three pedigrees:

• The "Asym" family made up of the father and the paternal grand and greatgrandparents • The "Grt grd par." with grand and great-grandparents on both sides.

• The "Grd par." with grandparents but no great-grandparents on both sides. See Figures 3 and6 Depending on the conditioning, the loss of symmetry (i.e. a loss of one familial branch) may induce or not a significant loss of information. When we do not force the presence of a case of disease (0-conditionning), the loss of a branch can remove the pivotal genetic information. This is the case of figure 15-A: the blue ROC curve is much less sensitive than the orange one although the number of members is similar. On the contrary, when the information has to appear in the pedigree, thus in the remaining branch, the other one becomes less necessary, and the pedigree works as well as a pedigree containing its two branches (figure 15-B). For the oncogeneticist, knowing only one branch of a pedigree may be sufficient as long as this branch includes at least one cancer case.

Ascendants against descendants

We previously suggested that the prediction using the descendants should be worse than the one including the ancestors since recent generations are less exposed to the disease risk. We tested this hypothesis and here is the result (A) in the standard context then (B) when we condition to at least one disease case per pedigree:

A B

Figure 16: Comparison of predictors whether pedigrees contain rather ascending or descending generations and whether we condition (B) or not (A) the occurrence of at least one disease by pedigree.

These results contradict our initial guess. The prediction using children or grandchildren is better at identical number of members per pedigree. Perhaps this is happening because a pedigree containing children and/or young adults imposes a high age requirement on individual 0. Also, conditioning one cancer at least per family implies that they occur more likely at younger ages and therefore must be caused by a deleterious mutation.

3.9 Influence of parameters 𝑝 " , 𝑝 # and 𝑝 $ For Reg(4, 2)-type pedigree containing 22 individuals, the five ROC curves below correspond to increasing mutation frequencies but respectively weaker and weaker areas under the curve (AUC). Thus, we see that the predictor is all the better, both in sensitivity and in specificity, as the prevalence of the mutation is low. Another parameter of the model is the penetrance of disease K in mutation carriers. We tested a penetrance 𝑝 $ equal to 20%, 40% and 70% (this latter similar to that of BRCA mutations). For its part, the 𝑝 # risk of occurrence of the disease in non-mutated individuals (wild-type) during its lifetime can vary greatly and we have tested values ranging from 0.1% to 10%: It appears that a high penetrance in case of mutation improves the quality of the predictor while conversely, a higher incidence of the disease in the general population decreases this quality. Clearly, the greater the gap between penetrance and incidence of disease at any age, the better the predictor becomes.

Overall, these effects are in line with what we would expect: the more the disease K is confined to mutated individuals, the easier the prediction. However, these parameters do not all have the same influence. The incidence of the disease in wild-type individuals is not very influential, while penetrance in mutation carriers has a strong impact on the performance of the prediction. The frequency of mutation is also very influential and its direction of variation seems less intuitive: the more frequent the mutation, the more difficult its detection.

Conclusion on the cost / effectiveness trade-off

Previous experiments with family pedigrees of growing size show that there is an optimal size for the prediction. This notion of optimum is reinforced by the increase in the computational cost with the pedigree size, which might even shift the practical optimum towards smaller sizes. A difference should also be made regarding the cost of collecting phenotypes, depending on whether the data is preexisting in a database or collected specifically for mutation diagnosis. Typically, we

A B

will seek to get information regarding parents, grandparents and greatgrandparents. This ancestral line can be completed, albeit with a relative gain in performance. One should probably avoid going too far, as illustrated by the three ROC curves below: We can observe in this figure that by adding the last 20 members no gain is obtained, which suggests that a very exhaustive search for collateral relatives may cost more than it helps. The three pedigrees concerned are represented nested in the following diagram: 4 Generalization and extensions

Model with two mutations

In this model, two cross-effect mutations are implemented. We could probably call them "polymorphisms" because separately they have minimal penetrance and no impact on the age of onset, but when they are present together in the genome, the risk of cancer increases considerably as well as its precocity, with a penetrance similar to previous single mutation models. The phenotype is therefore regulated in a similar way, but it is not the presence of a mutation that generates the increased risk but the simultaneous presence of both.

The difference between the single and double models lies in the transmission of mutations which, of course, obeys Mendel's laws, but results in different characteristics. Thus, in the double model we consider that the genome of an individual takes its values from a set of 4 elements:

{ wild-type, mutated 1, mutated 2, mutated 1 and 2 } = { (0, 0), (1, 0), (0, 1), (1, 1) }

For each mutation 1 or 2, the probability of transmission from parents to their children is still governed by Table 1. In the single mutation model, a child of two parents carrying the same heterozygous mutation has the probability 3/4 to be a carrier. In the double mutation model, a child whose two parents are doubly mutated (always heterozygous), now has only a probability (3/4)² » 0.56 of being in turn doubly mutated. Note that we continue to consider only heterozygous mutations by ignoring homozygous cases, previously thought to be lethal. But in the case of polymorphisms that are not frankly deleterious, the possibilities of homozygosity could legitimately also be considered. This would complicate the models2 so much without likely providing any additional information that we have ruled them out.

The parameters of this double-mutation model are therefore slightly more numerous, namely:

• The frequency of mutation 1 in the general population, fmut1

• The frequency of mutation 2 in the general population, fmut2

• The penetrance of K for non-mutated subjects, p0

• The penetrance of K for carriers of mutation 1 only, p1

• The penetrance of K for the carriers of mutation 2 only, p2

• The penetrance of K for doubly mutated subjects, p1,2

In the cases studied below, the values of p0, p1 and p2 are close and low while the value of p1,2 is high. Besides, the law of the diagnosis age was assumed identical, for carriers of a double mutation, to the case of carriers of a single mutation in the previous sections. The problem therefore comes down to the prediction of the doubly mutated state which is the only one to constitute a significant risk of disease. The double mutation probability formula is an obvious adaptation of the single mutation model (1).

( In this two-mutation model, the number of genotypes to be scanned in the bruteforce method is 4 N where N denotes the size of the pedigree. So, except for extremely small trees, the simulation method should be preferred and can be written as an obvious variant of (2).

The performances obtained in this model with two mutations are overall a little lower than in the one-mutation model. But the phenomena described above persist. In particular, performance reaches a practical optimum around fifteen individuals as shown by the little gain in efficiency when going to 22 members, after what a loss is observed when going to 42 members:

Other predictors

Calculating the probability of the genotype conditional on the phenotype is the most reliable estimate of an individual's risk of mutation. However, other predictors naturally come to mind. Given the popularity of likelihood maximization methods, one predictor might be the genome value that maximizes the likelihood of a phenotype. For small families, this maximization calculation can be done by bruteforce, but for more realistic sizes, it is necessary to move towards a simulated annealing algorithm. Let us consider the case of the single mutation model and, as always, denote by F [resp. G] the variable which compiles the phenotypes [resp. genotypes] of individuals in the family. We still rate 0 the individual for whom the prognosis is to be made. Let us denote by 𝐺 [ 0 and 𝐺 [ 1 respectively the genotypes which maximize the likelihood of the phenotype with the individual 0 non-mutated or mutated respectively, i.e.:

𝐿(𝐹|𝐺 = 𝐺 [ X ) = max{ 𝐿(𝐹|𝐺 = 𝑔); 𝑔(0) = 𝑖} 𝑓𝑜𝑟 𝑖 ∈ {0,1}
A natural predictor of mutation is:

𝑝 $ = 1 { cd𝐹e𝐺 [ $ f g cd𝐹e𝐺 [ # f }
That is, the mutational state of individual 0 would be the one found in the genome that makes the phenotype most likely. Unfortunately, the performance of this predictor is catastrophic. It can be modified to obtain a usable variant by introducing:

𝑝 h = 1 { cd𝐹e𝐺 [ $ f g h cd𝐹e𝐺 [ # f }
where s is a threshold to adjust. By varying s, we thus obtain a family of predictors based on the likelihood ratio. In the case of small family pedigrees, these predictors yield results very close to the thresholding prediction of the conditional probability that we have used so far. Indeed, on the graph below, the ROC curves obtained by the brute-force method for each of the two types of predictor merge. The same is true for the curve corresponding to the predictor based on thresholding the conditional probability, but this time calculated using simulations. In the case of maximizing the likelihood of the genome by simulated annealing, the result is not as good and we must certainly conclude that this simulated annealing should be improved either in the choice of the kernel enabling to explore the space of the genomes, or in the setting of the temperature scheme (number and duration of temperature steps).

If we condition at least one case of disease, we obtain, as reported above, a lower performance but with an identical hierarchy:

Additions considered

We have limited the statement of our results to cases that we suppose the most interesting. The computing time is also a limiting factor because to get regular andhopefully -reliable ROC curves as above, the computing time is counted in hours or even days. However, various extensions are possible:

• To carry out for the double mutation model all the experiments done for the single mutation model • To evaluate a double mutation model taking into account homozygous mutations • To improve the annealing method based on the likelihood-ratio predictor so that it equals the brute-force method (when the latter is feasible).

• To examine the sensitivity of predictors to parameter estimation errors or even to model change, i.e. simulate families according to the one-mutation model and analyze them with the model for two mutations, and vice versa.

Limitations

General approaches like the one used in the article face several limitations. The first one concerns the type of disease that was implemented in our models: this disease was uniform and occurred according to two laws of age, one for mutation carriers and another for non-mutated individuals. When we consider the cancer predispositions encountered in the oncogenetic consultation, for example those caused by BRCA mutations, several cancer locations are favored by the mutations but with an incidence varying a lot by location: breast and ovaries are the organs that are the most often affected, but case of cancer striking other organs can also occur. Besides, breast cancers do not constitute a single and simple category since there are different histological types and a classification according to hormone receptors. All these subclasses influence both the age at disease detection and the prognosis. Therefore, the oncogeneticist needs not only to consider in a pedigree the number of cancers, but also their type in detail.

Conclusion

This study, although based on rather simple algorithms, yields practical conclusions about the daily work of the oncogeneticist. In particular, it suggests that there is an optimal size for familial pedigrees: bigger is not always better. It also quantifies the prediction quality for other shapes and sizes of pedigrees. Our computations have focused on a single-mutation model for computational feasibility, but we have also considered a double-mutation model with interaction, with similar results. This second type of model fits present oncogenetics issues since nowadays, many of the hereditary cancer risks seem to originate from the interaction between several nonpathogenic variants.

Figure 1 :

 1 Figure 1: probability laws used to model the onset of cancer according to age.

Figure 2 :

 2 Figure 2: Diagram of a regular pedigree of 4 generations and 3 children per couple. Women are colored in pink, men in blue and regardless of gender in grey. Individual quoted 0 is the one considered for the calculation of the mutation risk.

Figure 3 :

 3 Figure 3: Regular pedigree with 4 generations and 1 child per couple

Figure 4 :

 4 Figure 4: Design of the pedigree « wide 3 gen. » with 15 members.

Figure 5 :

 5 Figure 5 : design of the pedigree "Cousins" (19) after the addition of 8 cousins (in green) to the pedigree "Uncles" (11)

Figure 6 :

 6 Figure 6: Asymmetry created between two pedigrees Reg(3, 1) by moving individual 0 down of one generation and omitting his father.

Figure 7 :

 7 Figure 7: Two kinds of pedigree with one or two descending generations

  n b N m / e Y + d p z 3 b i N b E e w 2 J N X B J 7 L 9 0 9 5 X / q 7 O 9 G O j D B 9 e D o p 4 K x 9 j u h H e p 3 K 3 Y k + N v X R l y K I i z u E d 5 T V g 4 5 f 0 9 o 9 O U r n d 7 t 9 z l v 7 l K y 9 p Y + N o K 7 u w p 6 Y G j h 8 / 5 G J y 2 m t G 7 Z u v 4 f X 3 v o 3 / q O V i H 1 9 C g 9 9 y G P T i A I 2 i T 9 1 2 w G K w G a + H X 8 H v 4 I / w 5 L Q 0 D r 1 m B P 8 Z M + A u e F L Y s < / l a t e x i t >

Figure 8 :

 8 Figure 8: ROC curves for the mutation predictor; context: Reg(3,2) family tree conditioned by at least one case of K or unconditioned.

Figure 9 :

 9 Figure 9: Performance of predictor depending on the conditioning of a minimal number of disease case per family (context: Reg(4, 2) pedigree, one deleterious mutation)

Figure 10 :

 10 Figure 10: Influence of the generation number on the model performance (context: one mutation, conditioning = 0 and 1 children per couple).

Figure 11 :

 11 Figure 11: Influence of children number per couple on model performance; context: one mutation and 4 generations.

Figure 12 :

 12 Figure 12: Performance comparison of two pedigrees of 15 members, one large and the other one high

Figure 13 :

 13 Figure 13: Performance comparison of two pedigrees of 15 members, one large and the other one high but constraining the presence of at least one disease case per pedigree

Figure 14 :

 14 Figure 14: Influence on predictors of the presence of uncles and cousins in pedigrees (context: 1 mutation, 3 generations, 1 child per couple, except when adding uncles and cousins, and condit. of 1 disease case)

Figure 15 :

 15 Figure 15: Influence of asymmetry on predictors performance (context: 1 mutation and one child per couple; A -no condit. and B -condit. of 1 case)

Figure 17 :

 17 Figure 17: Predictors performance according to the mutation prevalence (context: 4 generation pedigree and 2 children per couple, thus 22 individuals, condit. of 1 case of disease))

Figure 18 :

 18 Figure 18: Influence on the prediction of the deleterious mutation penetrance (A) or its incidence in non-mutated individuals (B).

Figure 19 :

 19 Figure 19: Predictor performance when pedigree size increases from 15 members to 22, then to 42.

Figure 20 :

 20 Figure 20: The 3 embedded pedigrees corresponding to the 3 previous ROC curves: 15 members in pink for the basal pedigree, + 7 members in blue for the second (N = 22) and finally + 20 members in grey for the last one (N = 42)

  = 𝑔) 𝑃(𝐺 = 𝑔) B∶B(X)Z($,$)

Figure 21 :

 21 Figure 21: Comparison of the predictor performance according to the number of members included in the pedigrees; context: 2 mutations interacting and conditioning of at least one disease case per family

Figure 22 :

 22 Figure 22 : Comparison of performance of two predictors (thresholding the conditional probability or likelihood ratio) calculated either by brute-force, or by simulated annealing algorithm ("ann." for simulated annealing); context: families Reg(3, 1) of 7 individuals, one deleterious mutation and no conditioning of disease case

Figure 23 :

 23 Figure 23 : Comparison of the performance of both predictors (thresholding the conditional probability or the likelihood ratio) obtained either by brute-force, or simulated annealing (context: Reg(3, 1) families of 7 individuals thus and one deleterious mutation and conditioning of at least one disease case per pedigree)

Table 1 :

 1 Probability for an individual to inherit a deleterious mutation from his parents depending on whether or not they are carriers

	¯ mother / father ®	wild-type mutated
	wild-type	0	0.50
	mutated	0.50	0.75

A case in a family is called "index" if he is the first family member who comes to consultation. Most often, he has reported a cancer that corresponds to the familial syndrome, although it is not mandatory.

In particular, it would be necessary to rule on the lethal nature or not of a double homozygosity, on a difference in penetrance depending on whether one of the two mutations is homozygous or not, etc.