D -Kondo N

Moreno-Barbosa

Celis De

Alonso

Ramos-Méndez

DNA damage modeled with Geant4-DNA: effects of plasmid DNA conformation and experimental conditions

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Monte Carlo Track Structure (MCTS) codes have allowed us to deepen our understanding of radiation damage in biological tissue. MCTS codes provide detailed information of each and every interaction of an energetic ionizing particle with matter along the particle track. A list of MCTS codes and past works on DNA damage modeling was published along the last decades, see e.g., [START_REF] Nikjoo | Track structure in radiation biology: theory and applications[END_REF][START_REF] Berger | ETRAN, Monte Carlo code system for electron and photon transport through extended media[END_REF], [START_REF] Booz | A microdosimetric understanding of lowdose radiation effects[END_REF], Chatterjee et. al.(1985), [START_REF] Friedland | Simulation of DNA damage after proton and low let irradiation[END_REF], [START_REF] Leenhouts | Radiation Energy Deposition In Water Calculation of DNA Damage and Its Association with RBE[END_REF], [START_REF] Nikjoo | Track structure in radiation biology: theory and applications[END_REF]. [START_REF] Toburen | Development of Monte Carlo Track Structure Codes[END_REF], [START_REF] Turner | Physical and chemical development of electron tracks in liquid water[END_REF], [START_REF] Wright | Calculations of Physical And Chemical Reactions Produced in Irradiated Water Containing DNA[END_REF], [START_REF] Zaider | The applications of track calculations to radiobiology I. Monte Carlo simulation of proton tracks[END_REF]. Radiation damage to biological tissue has been associated with the breakage of cellular DNA caused by either direct or indirect action. In the direct action, the damage is caused by energy transfer to DNA molecules, whereas in indirect action, the creation of reactive chemical species by radiolysis of the medium in close proximity to the DNA induces the damage by reacting with the DNA molecule [START_REF] Hall | Radiobiology for the radiologist: Seventh edition[END_REF]. Simulation with MCTS has historically allowed for detailed studies of the interaction of ionizing radiation with biological structures such as DNA damage through the direct and indirect pathways.

The Geant4-DNA project [START_REF] Bernal | Track structure modeling in liquid water: A review of the Geant4-DNA very low energy extension of the Geant4 Monte Carlo simulation toolkit[END_REF]Incerti et al., 2010;Incerti et al., 2010;[START_REF] Incerti | Geant4-DNA example applications for track structure simulations in liquid water: A report from the Geant4-DNA Project[END_REF] extended the Geant4 toolkit to include track structure capabilities in 2010. More recently, Geant4-DNA added radiation chemistry models to simulate all stages of water radiolysis including the physical, pre-chemical and chemical stages. For the physical and pre-chemical stages, a comprehensive description of the processes and models available in Geant4-DNA has been given elsewhere [START_REF] Karamitros | Modeling Radiation Chemistry in the Geant4 Toolkit[END_REF][START_REF] Karamitros | Diffusion-controlled reactions modeling in Geant4-DNA[END_REF], [START_REF] Incerti | Geant4-DNA example applications for track structure simulations in liquid water: A report from the Geant4-DNA Project[END_REF], [START_REF] Shin | Geant4-DNA simulation of the pre-chemical stage of water radiolysis and its impact on initial radiochemical yields[END_REF]. For the chemical stage, two methods were provided; one based on the Step-By-Step approach [START_REF] Karamitros | Modeling Radiation Chemistry in the Geant4 Toolkit[END_REF][START_REF] Karamitros | Diffusion-controlled reactions modeling in Geant4-DNA[END_REF], the other based on the Independent Reaction Times (IRT) approach ( Clifford et al., 1981;[START_REF] Clifford | Stochastic models of multi-species kinetics in radiation-induced spurs[END_REF]Greenet. al1990;[START_REF] Pimblott | Stochastic modeling of partially diffusioncontrolled reactions in spur kinetics[END_REF]Ramos-Méndez et al., 2020) Geant4-DNA has been used widely for DNA damage simulations. Applications included the simulation of a fibroblast cell nucleus irradiated with protons [START_REF] Meylan | Simulation of early DNA damage after the irradiation of a fibroblast cell nucleus using Geant4-DNA[END_REF], chromatin fibers irradiated with protons and alpha particles (de la [START_REF] De La Fuente | Accounting for radiation-induced indirect damage on DNA with the Geant 4-DNA code[END_REF], a bacterial cell irradiated with electrons and protons (Lampe et al., 2018), a fractal cell nucleus irradiated with protons [START_REF] Sakata | Evaluation of early radiation DNA damage in a fractal cell nucleus model using Geant4-DNA[END_REF] and straight DNA segments irradiated with electrons (Lampe et al., 2018). In these works, direct and indirect contributions to the DNA damage in complex structures were simulated. In all cases, Geant4-DNA reproduced the trend of the experimental results as a function of the particle linear energy transfer (LET). However, the use of such complex geometries had the disadvantage of demanding longer times (hours to days) to prepare the geometry for irradiation together with the particle source. In addition, experimental procedures to measure DNA breakage in cells has proven to be a complicated process with multiple techniques available [START_REF] Figueroa-González | Strategies for the evaluation of DNA damage and repair mechanisms in cancer[END_REF]) not providing the sufficient spatial resolution to directly compare calculated to measured DNA breakage. For instance, in most cases, the presence of repair proteins is used as a site marker providing an indirect method to estimate DNA damage.

However, MCTS accuracy for DNA damage estimation can be questioned, with SSB efficiencies ranging from 25% to 65% (Lampe et. al., 2018;[START_REF] Meylan | Simulation of early DNA damage after the irradiation of a fibroblast cell nucleus using Geant4-DNA[END_REF][START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF]Friedland et. al., 2011) compared with the experimental results of 12% [START_REF] Milligan | Variation of single-strand break yield with scavenger concentration for plasmid DNA irradiated in aqueous solution[END_REF]. This could lead to an over estimation not only of SSBs, but also of more complex break types that depend on SBs distributions. With the recent improvements to the physical and chemical models of Geant4-DNA/TOPAS-nBio [START_REF] Shin | Geant4-DNA simulation of the pre-chemical stage of water radiolysis and its impact on initial radiochemical yields[END_REF][START_REF] Shin | Geant4-DNA simulation of the pre-chemical stage of water radiolysis and its impact on initial radiochemical yields[END_REF], which have been observed to have a great impact on •OH and eaq -, with differences around 35% and 22% were observed at 1µs (Peukert et. al., 2019). Therefore, a decrease in DNA SBs is expected, which is of particular interest to the radiation chemistry area moving forward into temperature and oxygen depletion effects on DNA SBs.

In this work, the IRT approach implemented in Geant4-DNA was benchmarked for its use in predicting plasmid DNA damage response in experimental conditions. This was achieved with a less complex geometrical approach than a cell nucleus consisting of plasmids pUC19 (2686 base pairs) irradiated by low LET radiation. Due to their relatively low number of base pairs compared with DNA contained in a cell nucleus [START_REF] Meylan | Geant4-DNA simulations using complex DNA geometries generated by the DnaFabric tool[END_REF](3)(4)(5)(6), plasmids were easier to manage computationally than cellular DNA. Experimentally, their trademark "supercoiling" property allowed scientists to identify strand breaks experimentally and classify them as a SSB or a DSB [START_REF] Casali | E.coli Plasmid Vectors[END_REF][START_REF] Hempel | Determination of g-values for single and double Strand break induction in plasmid DNA using agarose gel electrophoresis and a curve-fitting procedure[END_REF][START_REF] Higgins | Topological Behavior of Plasmid DNA[END_REF]. The plasmid supercoiling property was measured by the super helix density. A super helix density different than zero meant that the DNA was under stress creating higher order helix geometries [START_REF] Vologodskii | Conformational And Thermodynamic Properties of Supercoiled DNA[END_REF]. The use of computational plasmid models allowed us to efficiently compare different scenarios of DNA damage response. The investigated scenarios included the comparison with published measured data of SSB yields as a function of DNA concentration, scavenger concentration and plasmid superhelix density (σ). In addition, comparison of calculated SSB and DSB yields as a function of the absorbed dose with published data was performed.

Methodology

A Geant4-DNA application using the IRT method was developed to model radiation damage in supercoiled plasmid DNA. The application used the low energy physical models available in the "G4EmDNAPhysics_option2" physics list and the chemical processes modeled by the IRT provided by Geant4-DNA version 10.7.p01 (Ramos-Méndez et al., 2020). The application facilitated the control of the radiation source, loading of plasmid data files, the number of plasmids geometries and radiation chemistry conditions.

Chemical stage simulation

The IRT implementation of Geant4-DNA includes more than 70 reactions to help model radiation chemistry processes (Ramos-Méndez et al., 2020). The heterogenous chemistry process that follow water radiolysis, has been successfully simulated using a reduced number of chemical reactions [START_REF] Pimblott | Stochastic simulation of the electron radiolysis of water and aqueous solutions[END_REF][START_REF] Ramos-Méndez | Monte Carlo simulation of chemistry following radiolysis with TOPAS-nBio[END_REF][START_REF] Shin | Evaluation of the influence of physical and chemical parameters on water radiolysis simulations under MeV electron irradiation using Geant4-DNA[END_REF]. Reactions used with the IRT in this work are presented in Table 1. The use of a reduced number of chemical reactions was done in order to decrease simulation times as much as possible while maintaining precision.

Experimental irradiation conditions of DNA plasmids included a cryopreservative which acted as a scavenger for • OH, thus protecting the DNA from damage after irradiation, e.g., dimethyl sulfoxide (DMSO), tris(hydroxymethyl)aminomethane (TRIS) or ethylenediaminetetraacetic acid (EDTA). The presence of cryopreservatives was modeled by first-order reactions considering their scavenging capacity. The reactions involving such scavengers also occurred at longer times (around 10 µs for 10 -4 M DMSO) than the heterogeneous chemistry stage (1µs). Thus, classes in Geant4-DNA dealing with the IRT method were modified to allow reactions between chemical species and buffered solutions until 20 µs. Reactions with dissolved O2 were also considered to recreate experimental aerobic conditions [START_REF] Milligan | The effect of superhelical density on the yield of single-strand breaks in γ-irradiated plasmid DNA[END_REF][START_REF] Milligan | Variation of single-strand break yield with scavenger concentration for plasmid DNA irradiated in aqueous solution[END_REF][START_REF] Milligan | Characterization of the Reaction Rate Coefficient of DNA with the Hydroxyl Radical[END_REF][START_REF] Tomita | Strand Break Formation in Plasmid DNA Irradiated in Aqueous Solution: Effect of Medium Temperature and Hydroxyl Radical Scavenger Concentration[END_REF][START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF] and these were obtained from [START_REF] Buxton | Critical Review of rate constants for reactions of hydrated electrons, hydrogen atoms and hydroxyl radicals (⋅OH/⋅O-in Aqueous Solution[END_REF], see Table 1. In reactions 11, 12, and 13, DNA referred to a deoxyribose-phosphate pair center. This was done to reduce the number of chemical species in the IRT, reducing computing time.

Reactions for simulation of radiolysis in pure liquid water

DNA damage reactions

No.

Reaction kobs (M

-1 s -1 ) No. Reaction kobs (M -1 s -1 ) 1 eaq -1 + • OH → OH - 2.95×10 10 11 • OH + DNA → OHDNA Equation 1 2 eaq -1 + H2O2 →OH -+ • OH 1.10×10 10 12 H • + DNA → HDNA 3.00×10 7 3 • OH + H • → H2O 1.55×10 10 13 eaq -1 + DNA 1.00×10 7 4 H • + H2O2 → • OH + H2O 9.00×10 7 5 • OH + • OH → H2O2 5.50×10 9 Oxygen Reactions 6 eaq -1 + eaq -1 → H2 + OH - 6.36×10 9 No. Reaction kobs (M -1 s -1 )
The reaction rate coefficient of reaction 11, k11 (M -1 s -1 ), was measured by [START_REF] Milligan | Characterization of the Reaction Rate Coefficient of DNA with the Hydroxyl Radical[END_REF], and is shown in equation 1 as a function of • OH scavenging capacity, s:

𝑘 !! = 1.32 × 10 " 𝑠 #.%& (1)
2.2. DNA Model: pUC19 Plasmids

7 H3O + + OH -→ 2H2O 1.13×10 11 14 eaq -1 + O2 → O2 - 1.90×10 10 8 eaq -1 + H3O + →H • + H2O 2.11×10 10 15 H • + O2 → HO2 2.10×10 10 9 eaq -1 + H • → OH -+ H2 2.50×10 10 16 • OH + HO2 → O2 + H2O 7.90×10 9 10 H • + H • →OH -+ H2 5.03×10 10 DMSO TRIS and EDTA No. Reaction kobs (M -1 s -1 )
No.

Reaction kobs (M -1 s -1 ) 17 In this work, plasmid "skeletons" for pUC19 were generated using a Brownian Dynamic (BD) approach [START_REF] Ermak | Brownian dynamics with hydrodynamic interactions[END_REF][START_REF] Huang | Dynamics of site juxtaposition in supercoiled DNA[END_REF].

Starting from an open circle composed of n evenly spaced vertexes and slowly deforming the configuration until a total simulation time of 5 µs was achieved. The main parameters for the plasmid generation algorithm were temperature, super helix density, and number of vertexes (Vologodskii & Cozzarelli, 1994). Plasmid super helix density is a measure of plasmid supercoiling. It is the ratio between the difference in the number of turns of the DNA helix in a relaxed state and the actual number of turns divided by the number of turns in a relaxed state (Vologodskii & Cozzarelli, 1994). This difference causes an energy change that plasmid generation algorithms use to perform small plasmid vertex displacements. The temperature and super helix density were used to compute the total energy of the configuration, which affects the grade of supercoiling of the plasmid. While the number of vertexes corresponds to the length of the plasmid and as so, the type of plasmid (pUC19, pB322, pVC40, etc). The initial parameters were 20 °C and a superhelix density (σ) ranging from -0.06 to -0.03 and 2686 vertexes. In total, a set of 40 pUC19 plasmids were generated, 10 for each σ value: -0.03, -0.04, -0.05 and -0.06 (Figure 1a). The mean writhe of these plasmids were: -11.83 ± 0.37, -9.94 ± 0.46, -7.85 ± 0.45, and -5.92 ± 0.51, respectively, calculated with the Gaussian integration method reported by [START_REF] Klenin | Irradiation of plasmid and phage DNA in water-alcohol mixtures: Strand breaks and lethal damage as a function of scavenger concentration[END_REF].

A basic DNA model consisting of semi-spheres was used to wrap the plasmid skeleton as shown in Figure 1b. The 2686 base pairs of pUC19 were accommodated in a deformed polygon of 179 straight sides of the same length. Each side contained 14 to 15 base pairs accommodated to avoid overlapping. The number of straights sides and base pairs could be changed by the user. However, we found that between 10 and 20 base pairs (3.4 to 6.8 nm) produced a reasonably well behaved figure.

The spatial and angular separations between two base pairs were 0.34 nm and 36 o , respectively. Each base pair was modeled using six cut spheres, two for the nucleobases, two for the deoxyriboses and two for the phosphates. The radii were 2.45 Å, 2.93 Å, and 2.7 Å, respectively. These values were obtained to achieve the same geometric volume as that provided by the DNAfabric software previously reported in [START_REF] Meylan | Geant4-DNA simulations using complex DNA geometries generated by the DnaFabric tool[END_REF].

DNA damage efficiency.

Due to different reaction pathways, not all hydrogen atom abstractions involved in the reaction between • OH and H • with DNA produced an indirect Strand Break (SB) [START_REF] Dizdaroglu | Mechanisms of free radical-induced damage to DNA[END_REF]. Thus, for modeling purposes an indirect SB is produced with an efficiency given by a probability λ, after one reaction between • OH or H • with DNA occurred. Several efficiencies values has been used ranging from 25% to 65% (Lampe et al., 2018;[START_REF] Meylan | Simulation of early DNA damage after the irradiation of a fibroblast cell nucleus using Geant4-DNA[END_REF][START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF]Friedland et. al., 2011). In this work, we calibrated the indirect SB efficiency by using a Nelder-Mead (NM) optimization algorithm [START_REF] Nelder | A Simplex Method for Function Minimization[END_REF]. We compared the simulation results with experimental data and minimized their difference using a twoparameter method shown by the following function:

𝐹(𝜆 '( , 𝜆 ( ) = 8 9𝑦 )*++ (𝑥) -(𝜆 '( 𝑂𝐻𝐷𝑁𝐴(𝑥) + 𝜆 ( 𝐻𝐷𝑁𝐴(𝑥))9𝑑𝑥 , ! , " (2)
where yreff was the measured reference data, OHDNA and HDNA the yields of the products from reactions 11 and 12, respectively, while the integration limits xi and xf are dependent on the specific simulation setup. The fitted parameters λOH and λH, were the SB efficiencies for • OH and H • , respectively. Initial values of λOH and λH were randomly taken from 0 to 1 and 0 to 0.1, respectively, with neither value allowed to be below 0 or above 1 during the NM procedure. The NM algorithm was conducted until the difference between the best and worst value of F(λOH, λH) converged within 0.1%. The NM procedure was repeated 100 times with different random number seeds and the results were used to obtain the mean and standard deviation of the SB efficiencies. This was done to mitigate the effect of the algorithm finding a local minimum instead of the global minimum. This calibration procedure was applied to scenario S1 (Table 2) and the obtained SB efficiencies were used for the remaining scenarios.

Strand Break Scoring

A direct SB was registered if at least 17.5 eV was accumulated in the sugar and phosphate volumes on one side of a single base pair [START_REF] Charlton | Calculation of Initial Yields of Single-Strand and Double-Strand Breaks in Cell-Nuclei From Electrons, Protons and Alpha-Particles[END_REF]Lampe et. al., 2018). For indirect SB, we performed the following.

1. Each individual base pair (nucleotide, phosphate and deoxyribose) was taken as one whole chemical species or molecule and inputted into the IRT, with its position being centered between the sugar-phosphate volumes. This was done based on the use of observed reaction rate constants between chemical species and DNA [START_REF] Perry | DNA condensation with a boron-containing cationic peptide for modeling boron neutron capture therapy[END_REF][START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF]). 2. OHDNA and HDNA yields were converted to SB yields by using the relationship 𝑆𝐵 = 𝑂𝐻𝐷𝑁𝐴 • 𝜆 '( + 𝐻𝐷𝑁𝐴 • 𝜆 ( . 3. SB were accumulated and classified as SSBs and DSBs in order to be compared with experimental data. In this regard, we defined one DSB as two SSBs located on opposite strands with a distance of 10 base pairs or less.

Cell survival fraction is attributed to complex DNA strand breaks like the ones studied in [START_REF] Charlton | Calculation of Initial Yields of Single-Strand and Double-Strand Breaks in Cell-Nuclei From Electrons, Protons and Alpha-Particles[END_REF][START_REF] Georgakilas | Induction and repair of clustered DNA lesions: What do we know so far?[END_REF][START_REF] Watanabe | Spectrum of Radiation-Induced Clustered Non-DSB Damage -A Monte Carlo Track Structure Modeling and Calculations[END_REF]. However, due to the lack of experimental resolution, these high complexity breaks are hard to measure. Most plasmid experimental results on the other hand are presented using only SSBs and DSBs by making use of the loss of super coiling and DNA linearization fractions [START_REF] Edel | Computer simulation of strand break yields in plasmid pBR322: DNA damage following 125I decay[END_REF][START_REF] Vyšín | Proton-induced direct and indirect damage of plasmid DNA[END_REF] and do not include DNA repair processes. In this validation work, Geant4-DNA simulations were compared against different experimental results by using the SSBs and DSBs yields, without considering the repair kinetics. However, Geant4-DNA has been successfully used alongside DNA repair models [START_REF] Sakata | Fully integrated Monte Carlo simulation for evaluating radiation induced DNA damage and subsequent repair using Geant4-DNA[END_REF].

Given the presence of a geometrical model like bound DNA volumes (filled with water for all simulations), we considered two scenarios involving the origin of primary chemical species. In the first scenario (A1) only chemical species created outside of the DNA volumes were inputted into the IRT. Species created inside the DNA were discarded without contributing to the chemistry. For the second scenario (A2), the simulation included both, species created inside and outside of the DNA volumes into the IRT. A2 was conducted only to serve as a comparison with A1 as inside the DNA volumes there should not be any free water molecules to be dissociated.

Simulation Setup

Four different setups were used to compare simulation results with experimental data. All setups shared the same geometry consisting of a 1 µm diameter water sphere (ρ=1 g/cm 3 ) containing up to 20 plasmids positioned uniformly (Figure 1b), as well as the same oxygen concentration ([O2]) of 0.27 mM. One pUC19 plasmid in the sphere was equivalent to a concentration of 5.62 µg/mL DNA. Simulation scenarios were configured to recreate the experimental conditions from their respective references (Table 2). Scenario S1 was conducted using 30 Gy of dose with a [DNA] of 50.6 µg/mL, [DMSO] ranged from 10 -4 to 1 M and the plasmids super helix density was -0.03. Scenario S2 configuration used 30 Gy of dose with [DNA] ranging from 5.62 -100.12 µg/mL, [DMSO] of 10 [START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF] Table 2: Parameters of the simulation setups; The scenarios were set to determine the dependence of SSB on DMSO concentration (S1), DNA concentration (S2), pUC19 superhelix density (S3), and for SSB and DSB on absorbed dose (S4). Each scenario was configured to recreate the experimental conditions of their respective references.

Results from scenarios S1, S2, and S3 were compared with data from [START_REF] Milligan | Variation of single-strand break yield with scavenger concentration for plasmid DNA irradiated in aqueous solution[END_REF][START_REF] Milligan | The effect of superhelical density on the yield of single-strand breaks in γ-irradiated plasmid DNA[END_REF][START_REF] Milligan | Characterization of the Reaction Rate Coefficient of DNA with the Hydroxyl Radical[END_REF]. Results from scenario S4 were compared with data from Tomita et al., [START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF][START_REF] Tomita | Strand Break Formation in Plasmid DNA Irradiated in Aqueous Solution: Effect of Medium Temperature and Hydroxyl Radical Scavenger Concentration[END_REF]. In Tomita et al., the authors used pBR322 (4362 bps) plasmids at a concentration of 30 μg/mL and the • OH scavenging capacity of the environment was 1.62×10 6 s -1 (0.1 Tris EDTA buffer). Thus, the concentration of DMSO for scenario S4 was calculated to match the • OH scavenging capacity used by Tomita in order to allow the comparison. The reaction kinetics of the products from reactions with scavengers was not considered in this work. Simulations were repeated using different random seeds; 100 times for scenarios S1, S2 and S4, and 24 times for S3. From the results, the mean yields and standard deviations were calculated. In all the scenarios, the energy deposited by particle tracks was accumulated in the water sphere until an absorbed dose was achieved. We considered that each individual track was independent from any the other.

Due to the low interaction probability of γ-rays in such a small volume (1 μm), the spectrum of secondary electrons produced by 60 Co or 137 Cs γ-rays interacting in water (neglecting pair production) was obtained to increase computational efficiency (Figure 2). For that, each spectrum was calculated in a separated simulation using condensed history Monte Carlo with TOPAS [START_REF] Perl | TOPAS: An innovative proton Monte Carlo platform for research and clinical applications[END_REF]. A γ-ray point source was positioned at the center of a spherical water phantom of 5 cm radius. The source was isotropic, emitting γ-rays of 0.662 MeV for 137 Cs, or two γ-rays of 1.33 MeV and 1.17 MeV for 60 Co. The spectrum of initial kinetic energy (at the point of energy transfer) of the secondary electrons produced by γ-rays was scored at 5-5.2 cm depth in water to account for photon scatter. For the SSB and DSB simulations (scenarios S1 to S4), the electron source was isotropic and uniformly distributed in the spherical water phantom of 1 µm diameter. The initial electron energy was sampled from the appropriate spectrum ( 137 Cs or 60 Co). The calculated yield of SSB, G(SSB), as a function of the DMSO concentration (scenario S1) is shown in Figure 3 for setup A1, see section (2.4). The values λOH of 13.77% ± 0.64% and λH of 0.74% ± 0.07% were obtained with the minimization algorithm. The results obtained agreed with the experimental measurements within one standard deviation. For setup A2 which included radiolitic species created in all regions of the simulation setup (see section 2.3), the resulted efficiency is 23.77%. Subsequent results in this work use the 13.77% efficiency value. The calculated dependence of SSB yield on DNA concentration (scenario S2) is shown in Figure 4. Results show reasonable agreement within the statistical uncertainties. The mean statistical uncertainty of the simulation results was 2.4%. Simulation agrees with measured data with a slope difference of 4.47 ± 2.3%. The G(SSB) response to superhelix density obtained in scenario S3 is presented in Figure 5. Simulation results are shown alongside experimental data.A reasonable agreement was achieved between simulation and experimental measurements within one standard deviation. Simulation mean statistical uncertainties were of 4.6%. Simulation results were obtained by considering correction factors from [START_REF] Milligan | The effect of superhelical density on the yield of single-strand breaks in γ-irradiated plasmid DNA[END_REF]. These correction factors account for the ethidium binding to supercoiled DNA, not simulated, which goes from 3.3 at 0 σ to 1.1 at -0.08 σ. Hence, a half-point value for ethidium correction of 2.2 was used resulting in satisfactory agreement with the measured data. Results for DNA damage as a function of absorbed dose (scenario S4) are shown in figure 6. Both SSB and DSB agreed with experimental data within statistical uncertainties with differences of 5.84% and 14% for SSB and DSB respectively. The curves of plasmid fraction are shown in Figure 7. Calculated results were fitted using the Cowan model [START_REF] Cowan | Breakage of double-stranded DNA due to single-stranded nicking[END_REF]. Fitting parameters for the Cowan model resulted in µ = 2.92×10 -2 SSB Plasmid -1 Gy -1 and Φ = 4.45×10 -4 DSB Plasmid - 1 Gy -1. For the fitting procedure, it was assumed that no strand break yields were produced at zero doses (µ0 = 0 and Φ0 = 0).

Results

Direct DNA damage results for different DNA concentrations are shown in Table 3. These results were obtained with 50,000 source electrons interacting in the 1 micrometer diameter sphere.

[DNA] (μg mL -1 ) 60 Co SSB (μmol J -1 )

137 Cs SSB (μmol J -1 ) 5.62 3.58×10 -8 ± 1.0×10 -8 3.11×10 -8 ± 7.7×10 -9 11.24 7.20×10 -8 ± 1.45×10 -8 6.21×10 -8 ± 1.1×10 -8 16.86 9.89×10 -8 ± 1.74×10 -8 1.06×10 -7 ± 1.45×10 -8 22.48

1.47×10 -7 ± 2.1×10 -8 1.28×10 -7 ± 1.6×10 -8 28.1

1.59×10 -7 ± 1.8×10 -8 1.74×10 -7 ± 2.4×10 -8 33.72

2.18×10 -7 ± 2.6×10 -8 1.92×10 -7 ± 1.9×10 -8

Discussion

In this work we extended and benchmarked the Geant4-DNA IRT implementation for its use as a DNA SSB and DSB prediction tool. An overall good agreement was achieved for all simulated scenarios, recreating the trend of published experimental results within one standard deviation. The results were applicable to low LET regimes while only considering radiolitic species created outside of the DNA volumes (setup A1).

For the simulation scenario S1, DNA damage efficiencies exceed experimental values of 12% for • OH (versus 13.77 ± 0.64% from our simulations) as reported by [START_REF] Milligan | Variation of single-strand break yield with scavenger concentration for plasmid DNA irradiated in aqueous solution[END_REF] and 0.8% for H • (versus 0.74 ± 0.07%) obtained by [START_REF] Aydogan | Monte Carlo simulations of site-specific radical attack to DNA bases[END_REF]. For • OH yields calculated in liquid water with the IRT implemented in Geant4-DNA (Ramos-Méndez et al., 2020, 2018), an exceess in the G value • OH radicals of about 10% exists at around 1 ns, around the time corresponding to the scavenging capacities used in this work. By reducing the • OH yields discrepancies, we expect an increase of up to 10% in SB efficiency, which will give us 23.77%, well in excess of 12%. The DNA geometry model might assist in resolving that issue having an impact on SB yields. For example, a 10% difference in SB yields was found when comparing A2 with A1. This difference is based purely on DNA geometry volumes. Hence, by using more detailed DNA models, it is possible to compensate the increase in SB efficiencies due to Both experimental and simulated results for scenario S3 showed a trend of monotonic increasing SSB yield with increasing superhelix density. This could be explained by a "protective" effect caused by the DNA volume [START_REF] Takata | Chromatin Compaction Protects Genomic DNA from Radiation Damage[END_REF].

In other words, due to the larger region of condensed DNA where chemical species cannot be created, the separation distance between inner DNA and chemical species increased, reducing the likelihood of reaction. The increase from -0.06 to -0.03 accounts for a difference of 4.5%, which suggest that the more packed the DNA is, the less likely it is to break, as mentioned by [START_REF] Takata | Chromatin Compaction Protects Genomic DNA from Radiation Damage[END_REF]. However, these effects in plasmids are not as accentuated as with chromatin fibers in which the SB yield can go up to 50 times between condensed and decondensed chromatin.

For scenario S4, simulated SSB yields had a linear response and eventually reach curvature towards 120 Gy. Although it was reported that SSB should behave linearly at low doses [START_REF] Klenin | Irradiation of plasmid and phage DNA in water-alcohol mixtures: Strand breaks and lethal damage as a function of scavenger concentration[END_REF][START_REF] Tomita | Monte Carlo simulation of DNA strand-break induction in supercoiled plasmid pBR322 DNA from indirect effects[END_REF], non-linear behavior of SSB at high doses had also been observed experimentally [START_REF] Kassis | Strand breaks in plasmid DNA after positional changes of auger electron-emitting iodine-125: Direct compared to indirect effects[END_REF](Kassis et. al. , 1999a)). As shown by Kassis et. al., SSBs exhibited a threshold of dose (around 30 Gy in absence of DMSO) after which it reached a plateau and remained constant. However, DSBs continued to increase nonlinearly. This suggested that as SSBs continued being produced with the dose increase, these also were being created closer leading to more DSBs. Such behavior was accelerated as DMSO decreased (as shown in Figure 8). That behavior was consistent with data from Kassis et al. These effects are expected to be seen in the scenario of high doses in cellular scavenging conditions (equivalent to around 0.1 to 1 M [DMSO]). Thus Geant4-DNA was capable of reproducing that behavior satisfactorily. By using the results from scenario S4, we provided a means to obtain curves of plasmid fraction response to dose. The obtained results behave consistently with published data. For instance, [START_REF] Vyšín | Proton-induced direct and indirect damage of plasmid DNA[END_REF] reported plasmid fraction response curves for 20 MeV protons. They obtained parameters for SSB and DSB per plasmid and Gy of 3.27×10 -1 and 1.26×10 -2 whereas we obtained 2.92×10 -2 and 4.45×10 -4 . The observed difference may be attributed to the difference in LET from ≈ 5.0 keV/µm to 0.3 keV/µm for 20 MeV protons, and 60 Co γ-ray, respectively.

The contribution of SSB produced by direct action under low LET radiation conditions (< 0.91 keV/um) [START_REF] Meesungnoen | Monte Carlo Calculation of the Primary Radical and Molecular Yields of Liquid Water Radiolysis in the Linear Energy Transfer Range 0.3-6.5 keV/?m: Application to 137 Cs Gamma Rays1[END_REF]) was expected to be low [START_REF] Klenin | Irradiation of plasmid and phage DNA in water-alcohol mixtures: Strand breaks and lethal damage as a function of scavenger concentration[END_REF][START_REF] Vyšín | Proton-induced direct and indirect damage of plasmid DNA[END_REF]. The calculated total direct contribution was less than 0.01%, whereas the experimental value measured by [START_REF] Tomita | Strand Break Formation in Plasmid DNA Irradiated in Aqueous Solution: Effect of Medium Temperature and Hydroxyl Radical Scavenger Concentration[END_REF] for indirect DNA damage contribution was 98.5 ± 2.5% (1.5 ± 2.5% direct DNA damage contribution). The differences in contribution values fell within the experimental uncertainty. The result also depends on the size of the spherical components of the DNA plasmid and the contribution of very low energy electrons. Yingxia et. al. (2021), showed that low energy electrons irradiating iced plasmids produced considerable DNA damage (4.86 -160 SSB µmol J -1 ). This data might provide an improved benchmark which is subject of a future work.

Conclusions

In this work the IRT approach implelemented in Geant4-DNA was succesfully modified to simulate DNA strand breaks and water radiolisis simulations at times greater than 1 μs. The algorithm was benchmarked for its use in estimating SSB and DSB in pUC19 plasmids using low LET radiation using both 60 Co and 137 Cs equivalent sources. We achieved agreement between simulated results and experimental data for all simulated scenarios within one standard deviation of SSB and DSB yields. The obtanied SSB efficiency values were close to experimental results, albeit outside of experimental uncertainty for • OH. This Geant4-DNA application will soon be incorporated as part of the Geant4 examples and can then be used as a reliable tool to calculate SSB and DSB in plasmids. This application may also be used in a future work to test higher LET regimes by using heavy ions.

Figure 1 :

 1 Figure 1: Plasmid DNA geometries used in Geant4-DNA. a) Examples of four pUC19 plasmids generated with a σ ranging from -0.03 to -0.06 used to calculate the DNA geometry for both direct and indirect DNA damage. b) Basic simulation geometry setup illustrating the use of plasmid skeletons to wrap the simple DNA model used in this work. The coordinate positions of DNA components for radiation chemistry simulation were located at the center of the semi-spheres.

Figure 2 :

 2 Figure 2: Secondary electron spectrum for 60 Co and 137 Cs. Energy spectrums were calculated using condensed history MC using TOPAS.

Figure 3 :

 3 Figure 3: S1) Response of G(SSB) (μmol J -1 ) to [DMSO]. Error bars represent statistical uncertainties, one standard deviation. Geant4-DNA IRT results (gray ◇), reference experimental data from[START_REF] Milligan | Variation of single-strand break yield with scavenger concentration for plasmid DNA irradiated in aqueous solution[END_REF] corresponds to pBR322 (☆), pUC18 (△), and pEC (○) plasmids.

Figure 4 :

 4 Figure 4: [DNA] (μg/mL) response of G(SSB) (μmol J -1 ). Error bars represent statistical errors, one standard deviation. Geant4-DNA IRT results are shown (◇), experimental data are from Milligan, 1993 (○, □, △).

Figure 5 :

 5 Figure 5: Response of G(SSB) (μmol J -1 ) to σ. Error bars represent statistical errors, one standard deviation. Geant4-DNA IRT results (□), experimental data from Milligan et. al. (1992) is shown (○) with a correspondent fit. Both, the experimental (solid line) data and simulation (dashed line) were fitted using a double logarithmic function to better represent their behaviors

Figure 6 :

 6 Figure 6: Dose response of SSB (left), and DSB (right). Error bars represent one standard deviation. Simulation results for SSB and DSB (○) are shown alongside Tomita's (1998) reference data (□). Simulation data 2 nd order polynomial fit is shown with the dotted line.

Figure 7 :

 7 Figure 7: Fraction of pUC19. Error bars represent statistical uncertainty, one standard deviation. Results separated in supercoiled (▽), relaxed (△) and linearized fractions (◇) fitted to the Cowan (1987) model and.

Figure 8 :

 8 Figure 8: DMSO effect on SSB and DSB as dose increases. Error bars represent one standard deviation. 1×10 -6 M [DMSO] is shown (□) against 2.28×10 -4 M [DMSO] (○), from MCTS simulations with pUC19 plasmids.

Table 1 : List of chemical reactions used in this work with Geant4-DNA IRT to perform DNA damage simulations

 1

	• OH + DMSO	7.10×10 9

. Reactions involving O2 were treated as pseudo-first-order. Products from reactions with DMSO were not included. 1 M = 1 mole dm -3

Table 3 : Direct DNA damage against DNA concentration for 60 Co and 137 Cs γ sources respectively.

 3 Results are shown alongside their statistical uncertainties, one standard deviation. and normalized per μmol J -1 .

  • OH yields once the chemical models are updated, maintaining a good agreement with the experimental value of 12%. Thus, reducing the DNA geometry dimensions might compensate the excess in efficiency. This can be explained as follows. By accepting chemical species produced in the DNA volumes, the chance of reactions at the contact with DNA molecules increased and hence increased the number of SB. The study corresponding to the DNA geometry dimension is a subject of a future work. Scenario S2 recreated the behavior of the measured SSB yields as a function of the DNA concentration within one standard deviation. Slope values differ slightly with 4.47 ± 2.3%. Simulated results tend to diverge from the fitted experimental data results as concentration of DNA increases; with an initial difference of 4.07% at 5.62 µg/mL and 8.04% at 101.19 µg/mL. However, since this fit is only taking into account the reported values without considering statistical errors from those, we can assume a reasonable agreement for concentration of DNA response of SSB yields.

Acknowledgements

This work is part of D-Kondo N. doctoral studies from the "Programa de Doctorado en Ciencias Fisica Aplicada, Benemerita Universidad Autonoma de Puebla" Mexico. The first author received financial support for his studies from fellowship 2019-000002-01NACF-05024 from CONACYT, México. Ramos-Méndez J.,Schuemann J. and Faddegon B. were supported by the NIH/NCI R01 CA187003 (TOPAS-nBio) grant. To Wook-Geun Shin from Department of Radiation Oncology, Seoul National University Hospital, Seoul 03080, Republic of Korea for his contribution in the development of Geant4-DNA IRT. All the simulation work was made with the support of the computer resources provided by the Laboratorio Nacional de Supercómputo del Sureste de México trough grant number 202001012C.


Farida Zehraoui
email: farida.zehraoui@univ-evry.fr

Naziha Sendi
email: sendi.neziha@gmail.com

Nadia Abchiche-Mimouni
email: nadia.abchichemimouni@univ-evry.fr

MS-LSTMEA: Predicting clinical events for Hypertension using Multi-Sources LSTM Explainable Approach

Keywords: Deep Learning, Explainable Artificial Intelligence, Attention mechanism, Hypertension, Electronic Health Records

Hypertension is a major risk factor for cardiovascular disorders and diseases.

Since it is expected to increase dramatically, effective Hypertension management becomes more and more critical. Early detection of patients with uncontrolled Hypertension would allow to employ personalized medicine for anticipating and providing the adequate drug class. In this paper, we propose MS-LSTMEA, a Multi-Source and Explainable Hypertension prediction Approach based on Long Short-Term Memory algorithm, for predicting both the drug class for patients and the date of their next medical appointment. MS-LSTMEA can successfully combine different sources of medical information about patients, represented by tabular data, while processing them separately and differently, and taking into account the temporal aspect of Electronic Health Records. In addition, it integrates an Attention mechanism which allows to improve the model results and to explain its outcomes.

Introduction

Hypertension is a major factor behind one of the most important causes of death in the world due to Cardiovascular diseases [START_REF] Roth | Global burden of cardiovascular diseases and risk factors, 1990-2019: Update from the gbd 2019 study[END_REF]. Improving upon the detection and the control of Hypertension, and working towards its early detection is of great importance in reducing the impact of cardiovascular diseases, their progression, and the associated risks. In order to confirm the diagnostic of Hypertension, repeated medical measurements are collected from patients.

With the availability of home kits for self-measurement allowing patients to be involved in the management of their disease, adherence to classic management rules could be advantageously replaced by personalized management. This enables to deeply interpret the characteristics of the patients and a better understanding of disease mechanisms, leading to better treatment choice and prevention. The collected measurements added to the available data sources are recorded in Electronic Health Records (EHRs).

An EHR can be viewed as a real-time patient health record (or a digital version of a patient's health information) with access to evidence-based decision support tools that can be used to aid physicians in their medical care task. An EHR contains the complete medical history of a patient. Such history can contain diagnoses, drug class plans, medications, comorbidities, allergies, laboratory results, radiology images, etc.

In this study, we aim to improve appointment times and drug class quality to each patient. To achieve this goal, we use a real dataset provided by a medical prescription platform "CEGEDIM" (used by 23,000 doctors in France), that has been collected from 3,000 physicians. This dataset contains 429,087 patients. Each patient is associated with a sequence of visits. A visit contains several patient's characteristics (age, gender, height, diastolic pressure, systolic pressure, pulse rate, etc.) represented as tabular data. At each visit, a drug class is associated to each patient.

The purpose of the work presented in this paper is to propose a Machine Learning (ML) tool for predicting drug class of Hypertension and the next medical appointment for a patient, by using EHRs represented as tabular data that comes from different sources. Among ML approaches, Deep Learning (DL) showed impressive results in various areas like image recognition, translation, etc. DL is actually widely used to improve the quality of healthcare.

Recently, recurrent neural networks (RNNs) [START_REF] Sherstinsky | Fundamentals of recurrent neural network RNN and long short-term memory LSTM network[END_REF] have encountered a great success in solving several ML problems due to their ability to learn hierarchical features in temporal domains [START_REF] Chimmula | Time series forecasting of covid-19 transmission in canada using lstm networks[END_REF]. It is why RNNs have been applied in different tasks such as representation learning, predictive learning, and classification. Nevertheless, simple RNNs are limited in their own since, for example, it is not easy to train them on long time sequences. This is due to the problems of vanishing and exploding gradients that are caused by back-propagation of errors over many time steps. The Long Short Term Memory models (LSTM) [START_REF] Sherstinsky | Fundamentals of recurrent neural network RNN and long short-term memory LSTM network[END_REF] [START_REF] Van Houdt | A review on the long short-term memory model[END_REF] deal with such problem by integrating memory units that allow to learn long temporal dynamics. LSTMs showed their effectiveness in various areas, especially in healthcare [START_REF] Iqbal | Covid-19 patient count prediction using lstm[END_REF] [START_REF] Chimmula | Time series forecasting of covid-19 transmission in canada using lstm networks[END_REF]. Despite their performance, and like most other ML algorithms, LSTMs are limited in their ability to explain results. They suffer from what is called the "black-box" problem.

In almost all approaches, the tabular feature sources that characterize the patients are fused before using the ML models (with or without preprocessing step). The problem is that these features generally come from different sources and their difference must be captured by the ML algorithm. This can be achieved by learning the tabular source representation when building the ML model.

To deal with the problems cited above, we propose a novel deep architecture through a Multi-Sources LSTM Explainable Approach (MS-LSTMEA) applied to Hypertension management. For that purpose, the approach successfully treats separately and differently several sources of tabular medical information, and takes into account long temporal trajectories of events embedded in EHRs. In addition, MS-LSTMEA uses an Attention mechanism [START_REF] Niu | A review on the attention mechanism of deep learning[END_REF] which consists in highlighting the information that has been most crucial in the predictive results. Such a mechanism can improve the model's explanation; a very important property in critical areas such as healthcare.

MS-LSTMEA is able to predict the drug class and the next medical appointment for a hypertensive patient. Indeed, the choice of the relevant and personalized drug class for each patient is very important since a bad choice can lead to disastrous effects that can conduct to the death of the patient. The next medical appointment time is also crucial for physicians when making decisions, because it provides information on the gravity of the patient's condition. The scientific contributions and challenges of this work are summarized below:

1. Novel multi-source LSTM architecture for dealing with the temporal trajectory of multi-sources events embedded in EHRs (temporal and heterogeneous aspects). 4. Advising support system with the potential to improve patient healthcare by predicting the date of his/her next visit (self-management of Hypertension).

In the next section, we will discuss related works in terms of temporal aspect, multi-sources of EHRs, and use of Attention mechanism that improves explainability of DL approaches. Then, Section 3 describes the proposed approach (MS-LSTMEA). Section 4 presents the experiment design, the results and a discussion allowing to highlight the main contributions of this work. Section 5 concludes and draws some improvements of our work and its perspectives.

Related works

Researchers have recently applied DL methods to EHRs in order to predict medical events. We can distinguish the proposed approaches using several criteria [START_REF] Shickel | Deep ehr: a survey of recent advances in deep learning techniques for electronic health record (ehr) analysis[END_REF][8]: the kind of features they consider, the feature representation they used, the type of DL architectures they used, the tasks they performed, their methods attempting to explain the obtained results, etc.

We focus on related works that consider the temporal aspects and the multisources nature of EHRs and those that use Attention mechanism for model explanation.

Temporal aspect of EHRs

There are different neural network architectures in the literature that deal with the temporal aspects of EHRs. These architectures are mainly composed of RNNs or Convolutional Neural Networks (CNNs).

RNNs are specially designed to deal with sequences and have been very successful in the NLP field. Since EHRs data are similar to text documents from many perspectives, most of the RNNs that represent EHRs were inspired from Natural Language Processing (NLP) applications.

Pham et al. [START_REF] Pham | DeepCare: A deep dynamic memory model for predictive medicine[END_REF] proposed Deepcare, that is based on LSTMs. It used temporal parameterization that takes into account the dynamics of the diseases.

DeepCare deals with the variable size of the input sequences, their irregular timing and the confounding interactions between the intervention and disease progression. It was used to predict unscheduled readmission of patients with mental disorders and diabetes in an Australian hospital. Lipton et al.

[10] used LSTM with target replication to represent the sequential aspect of data to perform diagnoses prediction. They evaluated their approach for predicting different diseases on EHRs data from patients in intensive care units.

Doctor AI [START_REF] Choi | Doctor AI: Predicting clinical events via recurrent neural networks[END_REF] was proposed to predict the diseases diagnoses, the treatments, and the date of the next visit. It projected the high dimensional inputs to a lower dimensional space and then used a Gated Recurrent Unit (GRU) network that is based on the same principle as LSTM with less training parameters. Doctor AI used a joint loss function combining the cross entropy for diagnoses and drug class prediction with squared loss for the next visit prediction. In [START_REF] Ashfaq | Readmission prediction using deep learning on electronic health records[END_REF], authors addressed the class imbalance problem using LSTMs for the prediction of unplanned readmission risk within 30-day. The proposed approach is based on contextual embedding and medical features of clinical concepts to generate visits representations. The authors succeeded to obtain good results on real data from a Swedish Hospital. Dasilva et al. [START_REF] Da Silva | Deepsigns: A predictive model based on deep learning for the early detection of patient health deterioration[END_REF] developed a model based on LSTM to assess the severity of patient's condition in order to anticipate its appropriate treatment. The authors showed good results when predicting vital signs of patients and therefore predict the prognostic indices to treat patients before the deterioration of their health condition.

In [START_REF] Sha | Interpretable predictions of clinical outcomes with an attention-based recurrent neural network[END_REF] and [START_REF] Zhang | Patient2vec: A personalized interpretable deep representation of the longitudinal electronic health record[END_REF], authors proposed the use of bidirectional RNNs (BiGRU or BiLSTM) that can be considered as variants of RNNs. These models consist first in applying the RNNs from the first to the last visit then, from the last to the first visit. The resulting representations are concatenated before performing predictions.

CNNs have also been adapted to deal with EHRs sequential data. EHRs are represented as matrices, with the horizontal dimensions representing the time and the vertical dimensions representing medical events [START_REF] Cheng | Risk prediction with electronic health records: A deep learning approach[END_REF][17] [START_REF] Li | Intelligent diagnosis with chinese electronic medical records based on convolutional neural networks[END_REF].

Other works have chosen to combine different DL algorithms for the purpose of improving the prediction results. For example, a CNN was combined with LSTM in [START_REF] Lauritsen | Early detection of sepsis utilizing deep learning on electronic health record event sequences[END_REF] and with BiLSTM in [START_REF]An attention-based cnn-bilstm hybrid neural network enhanced with features of discrete wavelet transformation for fetal acidosis classification[END_REF]. More recently, several works were published about the prediction of COVID-19 comparing the different architectures that deal with temporal data [21][22]. In [START_REF] Verma | Temporal deep learning architecture for prediction of covid-19 cases in india[END_REF], authors compared different architectures (CNN, BiLSTM, LSTM, and combined CNN-LSTM) for the prediction of COVID-19 and find that LSTM and combined CNN-LSTM give better results.

Multi-Sources DL and EHR

Medicine practices rely on information provided by different sources: clinical data, laboratory data, imaging data, narrative data, etc. In recent DL applications in the medical field, a fusion of complementary data sources was proposed aiming to improve the performance results [START_REF] Huang | Fusion of medical imaging and electronic health records using deep learning: a systematic review and implementation guidelines[END_REF]. Most recent works combined the tabular EHR data with medical imaging or textual data [START_REF] Huang | Fusion of medical imaging and electronic health records using deep learning: a systematic review and implementation guidelines[END_REF].

In [START_REF] Spasov | A multimodal convolutional neural network framework for the prediction of alzheimer's disease[END_REF][25], authors used CNNs to represent the images, combined the extracted high level features with clinical characteristics and then used a Feed forward neural network. These models did not deal with the temporal aspect of EHRs. In [START_REF] Reda | Deep learning role in early diagnosis of prostate cancer[END_REF] [START_REF] Qiu | Deep learning for automated extraction of primary sites from cancer pathology reports[END_REF], authors aggregated the prediction results using different approaches. Reda et al. [START_REF] Reda | Deep learning role in early diagnosis of prostate cancer[END_REF] built a new classifier that takes the prediction probabilities from the different sources as input for early diagnosis of prostate cancer. Qiu et al. [START_REF] Qiu | Deep learning for automated extraction of primary sites from cancer pathology reports[END_REF] used different ways to aggregate the predic-tions of the different sources: majority voting, max and min. The results of the three methods were merged and then fused with clinical models. Qiao et al. [START_REF] Qiao | MNN: Multimodal attentional neural networks for diagnosis prediction[END_REF] proposed a multimodal approach based on Attention that combines CNN and Bi-GRU to represent the textual and clinical data for diagnosis prediction. The authors captured the correlations between the two sources before using Joint fusion, where features interactions were performed before using bidirectional RNNs with Attention.

Bagheri et al. [START_REF] Bagheri | Multimodal learning for cardiovascular risk prediction using ehr data[END_REF] combined clinical temporal data with textual data for cardiovascular risk prediction using Joint fusion approach. They used a BiL-STM model that provided word embedding. These embedding were then concatenated with clinical data, resulting in a larger embedding which was then fed to dense layers before performing the prediction.

The multi-sources models proposed in the literature deal only with heterogeneous data (tabular EHR with images or/and text). Even if the data in EHR comes from different sources, these sources are merged and processed in the same way.

Our approach learns the tabular source representation during the model construction by simultaneously taking into account their temporal aspect. To our knowledge, it is the first approach that deals with this problem.

As a critical application area, healthcare requires explanation of the results provided by DL tools to be thoroughly explained. The interpretation of the models and/or the predictions is crucial since it can highlight the features that are involved in the predictions and therefore allows models to be validated or trusted and/or the discovery of unknown relationships between the features and the predictions. Therefore, different approaches that aim to ex-plain DL models used in healthcare were proposed in the literature [START_REF] Shen | An interpretable deep hierarchical semantic convolutional neural network for lung nodule malignancy classification[END_REF] [START_REF] Zeiser | Deepbatch: A hybrid deep learning model for interpretable diagnosis of breast cancer in whole-slide images[END_REF].

Recently, Attention mechanism has received a lot of interest in several application areas like Natural Language Processing (NLP), where it has been used to improve model results and to explain prediction. Most of the models applied to EHRs are inspired by those used in the field of NLP. We focus on the related works that use Attention mechanism for model or predictions interpretation.

Attention based interpretation

Attention mechanism in DL allows to highlight elements/features that are involved in the prediction. Self-attention provides a representation of the same sequence by linking its different positions.

Attention mechanisms are mainly combined with different RNNs models at different levels in the architectures. Shickel et al. [START_REF] Shickel | Deepsofa: a continuous acuity score for critically ill patients using clinically interpretable deep learning[END_REF] proposed DeepSOFA for the disease severity prediction. It used a GRU to handle the temporal aspect of EHRs and Self-attention to point up the most relevant time points when predicting hospital mortality. After embedding EHRs, the model assigns weights to time points that precede the current time point. The Attention scores allow to identify and highlight the time periods that influences the predictions. In [START_REF] Ma | Dipole: Diagnosis prediction in healthcare via attention-based bidirectional recurrent neural networks[END_REF] proposed in [START_REF] Sha | Interpretable predictions of clinical outcomes with an attention-based recurrent neural network[END_REF], based on a Bi-GRU and hierarchical Attention, performs mortality prediction. It is able to compute two levels of Attention: one level for medical codes and another one for the visits of a patient. Another model that introduced hierarchical Attention mechanisms, called Patient2Vec, was proposed by Zhang et al. [START_REF] Zhang | Patient2vec: A personalized interpretable deep representation of the longitudinal electronic health record[END_REF]. Patient2Vec used a BiGRU with a two-level hierarchical Self-attention in order to predict the future risk of hospitalization.

Authors proposed to group visits into times periods called sub-sequences and proposed hierarchical self-attention mechanisms at two levels: medical codes (variables) and sub-sequences. Kaji et al. [START_REF] Kaji | An attention based deep learning model of clinical events in the intensive care unit[END_REF] applied an Attention mechanism at the level of input variables, then used LSTM in order to predict several medical outcomes (daily sepsis, vancomycin antibiotic, myocardial infarction, and administration) over two-week-long patient Intensive Care Unit (ICU) courses using the MIMIC-III dataset. They provided, from a simple Attention mechanism, an interpretation at the level of variables and at the level of visits. More recently, Men et al. [START_REF] Men | Multi-disease prediction using lstm recurrent neural networks[END_REF] proposed to predict multiple diseases based on EHR data, from a hospital in Southeast China, using LSTM extended with time-aware (for handling temporal irregularities in EHRs) and Attention mechanism (for prediction interpretation). The latter allows to identify the relevance of each visit for the disease prediction.

The authors showed that their approach outperforms classical ML and DL methods for the prediction of the diagnoses of future disease.

In our work, we propose a novel approach, which is based on a joint fusion of data sources, that takes into account the temporal aspect of EHRs using LSTMs and the different sources of EHRs tabular data even though they can be merged directly at the level of inputs (early fusion) as done in related works. We showed that this representation improves the results.

Moreover, unlike existing approaches, the construction of the data sources embedding is not separated from the overall architecture, it is learned during the training step of the whole NN. Our approach deals simultaneously with multiple sources over time in the RNN since at each visit, the learned highlevel representation of the sources is transmitted to the LSTM units. We also use a Self-attention mechanism that helps to improve the results and to highlight the most important variables, sources, and patient visits. Unlike existing approaches that embed the Attention mechanism inside RNN units, or that use complex multi-level Attention, our approach implements a simple attention mechanism at the input level. This mechanism allows a direct interpretation of the model predictions.

Method

As already said, we present a novel explainable Multi-Source Approach that is based on LSTM architecture, called MS-LSTMEA. This approach is dedicated to predicting clinical events for each patient, such as drug class and the date of the next medical appointment.

We chose LSTM architecture due to its ability to remember previous entries, which makes it well suited for time series prediction. MS-LSTMEA successfully integrates multiple tabular data sources containing medical patients information and takes into account the temporal aspect of EHRs. The architecture of MS-LSTMEA is showed in Figure 1. Every data source is passed to fully connected layers that build individual embedding for the source. These embeddings constitute abstract high level representations of the data sources that are concatenated before being processed by an LSTM unit. The goal is to predict the drug class at each visit and the date of the next medical appointment for each patient.

In addition, we introduced a Self-attention mechanism at the level of the input. This improves both the results and the explanation of the predictions and the model. Indeed, we are able to determine the most relevant features that are used for each prediction thanks to the attention mechanism that is located at the input level. The importance of the sources and the visits can then be deduced from the scores of the important features.

We represent each patient p ∈ P by a succession of n p medical events (visits)

v t (1 ≤ t ≤ n p ): p = [v 1 , . . . , v np ].
The goal is to predict the date of the next medical appointment (see Sec- Guided by the EHR data that we used, we defined the three following sources 1 :

• s(1): refers to the personal features (gender, age, etc.)

• s(2): corresponds to clinical features (clinical codes related to the diagnoses, the procedures of medications, etc.);

• s(3): contains comorbidity features that correspond to the presence of other diseases in addition of hypertension (disease of interest).

Each source of information is transmitted to fully connected layers after learning attention scores. 1 Our architecture is enough general since we can add easily other data sources that contain additional information about patients.

Attention mechanism

A visit v t at time t of a patient p is represented by a vector that contains the values of the sources of information (see Figure 1):

v t = [s t (1), s t (2), s t (3))]
, where s t (1), s t (2) and s t (3) are the three sources that describe the patients p at the visit v t . The sources s(k) t (k ∈ {1, 2, 3}) are represented by numerical vectors of dimension d k . We introduce a set of weights W k (k ∈ {1, 2, 3}) that allows to calculate the feature-level selfattention scores α t (k) associated to the source s t (k) (k ∈ {1, 2, 3}). These weights are learned during the training step and the attention scores α(k) t (k ∈ {1, 2, 3}) are computed as follows:

α t (k) = sof tmax(W k • s t (k)) (1)
The Attention scores α t (k) associated to the input patient visit v t and then transmitted to the fully connected layers in order to build the high level of representation as follows:

ŝt (k) = α t (k) ⊙ s t (k), k ∈ {1, 2, 3} (2)
where ⊙ denotes element-wise product.

Note that, unlike most of existing approaches that introduce Attention scores inside the LSTM unit (mostly with the aim to improve the prediction results), we add an Attention layer at the level of input sources to preserve interpretability. The computed score is based on learned weights W k (k ∈ {1, 2, 3}), that provide variables importance like in [START_REF] Kaji | An attention based deep learning model of clinical events in the intensive care unit[END_REF].

After this step, we obtain, for each patient p, high level representations s(k), of the sources s(k) (k ∈ {1, 2, 3}) using, a specific number l k of connected hidden layers.

at time t, these embeddings are concatenated in order to form one multidimensional vector ṽt defined as follows 2 :

ṽt = [s t (1); st (2); st (3)] (3)
For the rest of the paper, the components of the high level representation of the visit at time t, ṽt , are denoted as ṽtj , i.e., ṽt = [ṽ t1 , . . . , ṽtd ], where d is the dimension of ṽt (4)

Visits representation

The embeddings built for the visit of each patient are fed sequentially to the LSTM units that take into account the temporal aspect of these visits.

We use a one-to-one LSTM in order to provide a prediction at each time step. Figure 2 illustrates the temporal aspect of our model.

A visit v t ∈ V , at time t, of patient p ∈ P is embedded into a representation space where it is transformed to a vector ṽt . The resulting vector is fed to an LSTM unit that produces the prediction of the drug class or the next medical appointment. An LSTM unit (see Figure 3) is composed of a memory cell c, that keeps the observed information from the input to the current time step, and three gates that control the state of the cell: forget gate f that controls the retention or deletion of old information; input gate i, which selects the information that will be added to the cell state c; output gate o that selects the information that will be returned by the cell state.

+ tanh σ + + tanh + + + σ σ x x ht-1 ct-1 Ct ht ht 0 3 2 1 i t c' t o t f t Figure 3: LSTM Unit
The LSTM unit receives three inputs: c t-1 , h t-1 and ṽt . c t-1 represents the memory cell and h t-1 corresponds to the state of the LSTM unit at time step (t -1). The value of the initial LSTM state h 0 is equal to 0.

At time step t, three gates f t , i t , and o t are defined as follows:

f t = sigmoid(W f [ṽ t ; h t-1 ] + b f ) i t = sigmoid(W i [ṽ t ; h t-1 ] + b i ) o t = sigmoid(W o [ṽ t ; h t-1 ] + b o ) (5)
Where

(W f , b f ), (W i , b i ),
ct = tanh(W c [ṽ t ; h t-1 ] + b c ) c t = f t • c t-1 + i t • ct h t = o t • tanh(c t ) (6)
A each time step t, the state h t of the LSTM unit is updated before being propagated to the following output layer in order to predict the drug class or the date of the next medical appointment. For the classification task, we used a softmax activation function in the output layer in order to compute the probabilities π t associated to the predictions as following:

π t = sof tmax(W.h t + b) (7)
Where W represents the learned weights and b is the bias.

The prediction ỹt at time t of our model is defined by: ỹt = argmax φ i=1 π ti . A categorical cross entropy function E is used to train the model3 , it is defined a follows:

E = - 1 #P p∈P 1 n p np t=1 φ i=1 [y ti log(π ti )] (8)
Where #P represents the patients number in the training set, π t is the model output probabilities and y t is the true target (class)4 .

Experimentation

In order to validate our approach, we performed different experiments whose results are presented in this section. For these experiments, we have considered Hypertension, which is known as one of the major factors of heart disease and stroke. Gathering regular medical measurements helps to consolidate the diagnosis of Hypertension. These measurements are exploited and combined into EHRs with different sources in order to improve the prediction of the right drug class for a patient and the date of patient's next medical appointment. This allows to help patients to self-manage their disease.

The Dataset

The dataset considered in the experiments comes from a real data collected by 3,000 doctors, provided by CEGEDIM (medical prescription platform used by 23,000 doctors in France). This dataset contains the characteristics of 429,087 patients which are under follow-up for Hypertension.

Each patient is represented by a sequence of visits and received an average of twenty visits. We have retained, for each visit, 15 characteristics (Table to proceed to the predictions by using LSTM units. The model has been trained for 1000 epochs 5 . Dropout has been used between LSTM layers to avoid overfitting. All model hyperparameters are fixed using grid search.

Note that we tested different variations of LSTM: LSTM with stacked layers, Bi-LSTM and we did not see any improvement in results compared to vanilla LSTM.

Experimental results

Our algorithm has been validated by comparing the results given by MS-LSTMEA with those of other proposed methods that are based on different DL architectures. We considered: Doctor AI [START_REF] Choi | Doctor AI: Predicting clinical events via recurrent neural networks[END_REF] (based on (GRU), Med2vect [START_REF] Choi | Multi-layer representation learning for medical concepts[END_REF] (based on Auto-encoders), Deepr [START_REF] Nguyen | Deepr: A convolutional net for medical records[END_REF] (based on CNN), DeepCare [START_REF] Pham | DeepCare: A deep dynamic memory model for predictive medicine[END_REF] (based on LSTM), DeepPatient [START_REF] Miotto | Deep patient: An unsupervised representation to predict the future of patients from the electronic health records[END_REF] (based on Multilayer Perceptrons). The following aspects have been evaluated:

• Prediction results.

• Influence of source type on the prediction results.

• Influence of the high level representation of sources on the prediction results.

• Determining relevant clinical factors.

• Determining individual relevant clinical factors.

Prediction results

The and the date of the next medical appointment, are reported in Tables 3 and4. Table 3 shows that our approach greatly outperformed the existing ones while predicting drug class. Indeed, MS-LSTMEA obtained an accuracy by up to 80.4% while Doctor AI, which resulted in the second best score, obtained barely 68.3%. This can be attributed to the use of LSTM instead of GRU in DoctorAI since LSTMs remember longer sequences than GRUs.

While this is the main reason why LSTMs outperform other algorithms in tasks requiring modeling long sequences, we will see hereafter that the introduction of Attention mechanism further improves the results.

The prediction results are in the favor of our approach, compared to other approaches. This confirms the efficiency of the proposed visit representation, the use of an Attention mechanism and fully connected layers.

The prediction of the next medical appointment should be modeled as a regression problem. The time between visits (time period) could be biased because it depends not only on the availability of doctors and patients, but also on patient's health condition. In our approach, we decided to discretize the time period. Thus, the regression problem was transformed into a clas-sification problem.

After having explored several discretization values, and having had several discussions with physicians, we opted for a time period covering one month.

In the used EHRs, the smallest time between two visits is equal to 10 days and the largest one is up to 200. So, seven classes C k have been fixed, depending on the time period as follows:

Let k be the number of months since the previous visit, so seven classes are fixed, depending on the time period as follows:

• For k = 1 to 6:

Class C k = k, if (k -1) months ≤ time period < k months • Class C 7 = 7 if time period ≥ 6 months
In practice, physicians6 plan the next visit of their patient according to the patient's situation. When a patient's Systolic Blood Pressure is very high, physicians usually set their next visit very close (less than 25 days) for therapeutic adjustments. The visit is planned within three months if the Hypertension is severe, i.e., with associated risk factors (e.g., a patient who had a previous stroke) or with comorbidities (e.g., diabetes). Renewal of the drug class requires close monitoring. In the case of a renewal of a drug class for a person with few health problems and for young patients whose blood pressure has been stabilized for a long while, the delay goes beyond six months. The use of classes incremented by month applies well to patients in self-measurement at home as it creates an opportunity to detect any prob- lematic situations beyond the above framework.

Table 4 also illustrates the superiority of our MS-LSTMEA based approach, which gives an accuracy of up to 79.5%.

The influence of source type on the prediction results

For studying the influence of combining multiple sources of information on the performance, we have tested the different combinations 7 . Indeed, determining the significant sources of features plays an important role in predicting Hypertension. Choosing the relevant combination of features plays also a crucial role in improving the performance of the prediction model. when they are included to predict drug class. However, using only Personal features gives a low accuracy (score of 24.5%). This can be interpreted by the fact that the age, gender, weight etc. without Blood pressure measurements do not really make sense for drug class prediction. What is important to note is that the results show that we can produce more accurate predictions (accuracy=80.4% for the prediction of the drug class) when merging the different information sources, which provide measurements that can be different and complementary in their nature.

Including solely Clinical features for predicting the time before the next medical appointment gives 51.8% of accuracy that is higher than the accuracy obtained with personal data only (19.4%) or with only Comorbidity features (13.7%).

Influence of high level representation of sources on the prediction results

We evaluated the influence of high level representation of sources (construction of embeddings). For that purpose, we tested our model with and without construction of embeddings. In Table 6, it clearly appears that construction of embeddings from EHRs provides new data representations and improves prediction results whether when predicting the drug class or the date of a patient's next medical appointment. This is due to the ability of our approach to build relevant high level representation from the initial data sources using the fully connected layers.

Determining relevant clinical factors

The Attention scores over visits were used to construct heatmaps for the models trained to predict the drug class or the date of a patient's next medical appointment. The softmax scores for each feature were averaged over test patients who were predicted to belong to the same class in order to obtain patients-averaged Attention maps. This is helpful for highlighting the individual characteristics that had the most influence on the prediction of a class of interest.

Figure 4(A) shows a heatmap while averaging the Attention scores of patients when predicting the drug class (BB) at each visit. This shows, for each visit, the most relevant features that lead to the prediction. Figure 4(B) shows the averaged Attention scores obtained in Figure 4(A) on the features for each patients visit. This shows the visits which contributed the most to the prediction. We can see that Visits 4 and 10 are the most efficient visits for the prediction. Finally, the Attention mechanism has been also applied to each data source to observe the importance of sources. By averaging the scores of the features that constitute each data source, we measured their importance scores. The second heatmap of Figure 5 shows the patients-averaged Attention scores when predicting the date of the next medical appointment (class C 1 ).

The results shows that we can extract the features with the highest Attention scores for each visit from Figure 5 features with an Attention score up to 0.7.

Similarly, the heatmap of the

Discussion

This study presents a deep learning model based on LSTM and an Attention mechanism for performing two prediction tasks that are important in Hypertension management: the drug class for patients and the date of their next medical appointment. The model gives high accuracy of up to 80.4%

while predicting the drug class and an accuracy of up to 72.3% while predicting the date of the next medical appointment. The results showed that high level representation of Personal, Clinical, and Comorbidity features realized separately and differently, from EHRs using a fully connected layers, with a simple Attention mechanism, is a strategy which outperforms models that ignore any of these aspects. Moreover, we showed that combining different features results gives better performance than using either of them alone. In addition, our model proved itself to be able to capture the significant features, sources or visits that contribute the most to the predictions thanks to the use of the Attention mechanism. Thus, our model can easily provide an explainable prediction.

There are a many works regarding the prediction of high Blood pressure and cardiovascular risk, but few address the prediction of the drug class for Hypertension. Close to this is the work of Xiangyang and his colleagues [START_REF] Farhan | A predictive model for medical events based on contextual embedding of temporal sequences[END_REF] on the prediction of the optimal therapeutic pathway for Hypertension. They used RNNs including LSTM to predict the probability of reaching a target BP (Blood Pressure) control with different drug class pathways. As in our case, LSTMs reveal a high degree of accuracy in prediction. We think that the strategies of predicting the drug class and the optimal pathway can be considered as two complementary ways that can help patients and doctors collaborate and improve Hypertension management. Another related work concerns the personalization of hypertensive drugs [START_REF] Chunyu | Characterizing the critical features when personalizing antihypertensive drugs using spectrum analysis and machine learning methods[END_REF]. The authors considered five antihypertensive drug prescriptions for a 14,581 inpatient population where each patient is described by 15 features and belongs to the successful group if its DBP (Diastolic BP) <90 and SDB (Systolic BP) <140, to the unsuccessful group otherwise. To reveal the difference in characteristics of the features between the two groups, the authors applied statistical and machine learning methods to assign a weight to each feature. Each weight represents the importance of a feature for a patient regarding the effectiveness of a drug. The results of this study can help doctors in choosing the most adequate drugs for patients. Unfortunately, unlike in our study, the considered dataset does not take into account the temporal aspects related to Hypertension drug class. Moreover, the study could have gained in reliability if the considered dataset had been larger.

A limitation in the construction of the predictive model of current research is related to the use of only demographic features [START_REF] Basu Roy | Dynamic hierarchical classification for patient risk-of-readmission[END_REF] and clinical features [START_REF] Farhan | A predictive model for medical events based on contextual embedding of temporal sequences[END_REF].

The former discards a huge proportion of information in each patient records, whereas the latter ignores the knowledge and guidelines coming from human intelligence. Our experiments showed that using information regarding comorbidities improves the performance of deep learning algorithms. Such expert knowledge is considered as a domain-specific knowledge which is also a big challenge since it impacts the progression of the disease and the state of the patient.

When dealing with deep learning models for predicting clinical events with EHRs, it is important to ensure that the model is flexible. Such flexibility ensures the possibility of adding new data sources. In fact, the combination of multiple data inputs (biological, pathological, medical images, etc.) raises expectations and hopes in terms of understanding the causes and mechanisms of diseases as well as for the personalization of medical monitoring. Some of existing works did not take into account the integration of new data sources such as the work in [START_REF] Nguyen | Deepr: A convolutional net for medical records[END_REF]. Our MS-LSTMEA approach is integrative since it can easily handle new data sources, either being discrete or categorical, regardless of heterogeneity. This is because we build for each data source differently and separately an embedding to represent a visit.

Conclusion

In this work, we have proposed MS-LSTMEA, which is an LSTM-based model that can learn efficient patient representation from a large amount of electronic patient records and predict future events for patients, such as the drug class or the date of their next medical appointment. We have tested MS-LSTMEA on a large dataset of real EHRs, which achieved 82.90% accuracy and significantly outperformed many existing approaches. We have also shown that high level representation of Personal, Clinical, and Comorbidities features from EHRs, with an LSTM network highly improves the performance.

In the conducted experimentation, we have used a real set of data collected from 3,000 doctors, provided by CEGEDIM (medical prescription platform used by 23,000 doctors in France). The provided dataset describes the characteristics of 429,087 patients, each of whom is represented by a series of visits. Our method is generic enough to integrate several kinds of knowledge sources (Ex. genomics, images...). The approach and the results can be used to design a Decision Support System (DSS) that could be used by hospitals and physicians for Telemonitoring (follow-up) of Hypertensive pa-tients in terms of predicting the drug class and the date of their next medical appointment. At the last, but not the least, the proposed approach can help patients self-manage their Hypertension while providing explanations for the decisions of the model, contributing to the challenges of explainable Artificial Intelligence (AI), a subject of high interest in the current AI landscape [43][44].

This work raises different research perspectives. Indeed, we used the LSTM to model the temporal trajectory of patient visits. However, other models like the transformer [START_REF] Vaswani | Attention is all you need[END_REF] have recently proved their efficiency in the processing of sequential information and especially the EHRs [START_REF] Rasmy | Med-BERT: pretrained contextualized embeddings on large-scale structured electronic health records for disease prediction[END_REF]. Future work will focus on these models. In order to improve the explainability and transparency of our approach, we plan to combine the neural networks with symbolic AI approaches. We will focus on applying a multi-agent argumentation protocol that will be able to deal with prior knowledge injection from the official guidelines behind Hypertension. This has already been successfully done with an ensemble method based on multilayer perceptrons in [START_REF] Sendi | Towards a transparent deep ensemble method based on multiagent argumentation[END_REF].

CRediT authorship contribution statement Farida Zehraoui (FZ) is the main contributor and the supervisor of the work. She designed and formalized the method. She supervised the experiments. She also wrote the state of the art and highly contributed to the writing of the paper. Naziha Sendi realized the implementation and performed the experiments. She also wrote the first draft of results analysis.

Nadia Abchiche-Mimouni (NAM) participated to the conception of the method and to the writing and proof-reading of the paper. The Funding

2 .

 2 Specific processing of data sources of patient characteristics (personal, clinical and comorbidity), that are usually merged in the literature, leading to a high level patient representation learned during the model construction at each time step.

3 .

 3 Attention mechanism allowing to: (a) Highlight the information that has contributed most to the provided results (decision explainability). (b) Adapt the drug class by identifying important features that are specific for each patient (drug class personalization).

  , Ma et al. proposed a model that is based on bidirectional RNN, called Dipole. It is used to predict the future health information of a patient. To improve and interpret the prediction results, authors proposed three variants of their algorithm by introducing three Attention mechanisms that allow to compute different context vectors. These mechanisms differ by the way they capture the relationships between the hidden states of the RNN. Edward Choi and al. [34] developed a DL model, called RETAIN (REverse Time AttentIoN) that predicts heart failure in the outpatient setting. It is able to compute the importance of a variable or a visit through two RNNs and a two-level Attention mechanism which can identify important visits and the relevant clinical features in the highlighted visits. The tool GRNN-HA

Figure 1 :

 1 Figure 1: MS-LSTMEA architecture

  tion 3.3) or the drug class among a set of drug classes: Y = {y 1 , . . . , y m },where m is the number of drug classes. The different knowledge sources are usually concatenated (early fusion) and not separated in the literature since they share a tabular representation. In our approach, we used joint fusion were the high-level representations (embeddings) of sources are combined.

Figure 2 :

 2 Figure 2: Unfolded MS-LSTMEA Architecture

Figure 4 (

 4 C) shows the averaged Attention scores on the visits to highlight the most significant features that are used for the prediction. As we can see, the most significant features are 'Pulse' and 'Diastolic pressure' with an Attention score equals respectively to 0.63, and 0.62. The feature 'Systolic pressure' has an Attention score up to 0.58. In summary, we can confirm, once more, that the most significant features are the clinical ones.

  Figure 4(D) shows that the most important source type for prediction of the drug class is the Clinical features with an Attention score of up to 0.53.

Figure 4 :

 4 Figure 4: Attention heatmaps averaged over all test patients with visit Attention, features, and sources with the highest mean activation per visit while predicting the drug class BB.

  Figures6 and 7show the same information as in the previous section for each patient. Heatmaps are built for each patient individually in (A). Then, we plotted the Attention scores related to visits (B), features (C), and data sources (D). Figure6(B) shows that visit 5 is the most relevant one for pre-

Figure 6 (

 6 B) shows that visit 5 is the most relevant one for predicting the drug class for this patient.

Figure 6 (

 6 C) shows that 'Diastolic pressure', 'Pulse', and 'Prescription blood sugar' are the most important features with respectively 0.78, 0.74, and 0.74 as Attention scores. Figure6(D)confirms, once more, that the most important features belong to the clinical

Figure 5 :

 5 Figure 5: Attention heatmaps averaged over all test patients with visit Attention, features, and sources with the highest mean activation per visit while predicting the date of the next medical appointment (class C 1 ).

Figure 7 (

 7 B) shows that for the same patient, the visit 4 is the most determinant one, the features 'Duration' and 'Diastolic pressure' are the most determinant features regarding the prediction of the next medical appointment with respectively 0.71 and 0.60 as Attention scores (Figure 7(C)). This figure illustrates the sources scores and definitely confirms that clinical features represent the most determinant ones with an Attention score of up to 0.47.

Figure 6 :

 6 Figure 6: Personalized Attention heatmap while predicting the drug class BB for one patient.

Figure 7 :

 7 Figure 7: Personalized Attention map for the prediction of the next medical appointment for one patient (class 1 ).

  and (W o , b o ) are the (weights, biases) couples associated to the three gates.The unit state h t at the current time step t depends on the input ṽt , the states of the three gates and the previous LSTM state h t-1 . It is calculated by:

Table 1 :

 1 Most important features of a visit.

	1 describes the most important ones). A class of drug is associated to each
	layer and a mean squared error loss function.

Table 2 :

 2 Datasets description.

3.2. Experiment setup

An LSTM model has been trained for predicting the drug class and the date of their next medical appointment for a patient. The proposed model receives a multidimensional vector v t representing the features related to the visit at time t for a patient p. As explained previously, a visit v t is composed of Personal features (s(1)), Clinical features (s(2)), and Comorbidities features (s(3)), as shown in

Table 1 .

 1 First, the data is split up into training, validation, and test sets. 60% of patients data were used for training the models, 20% as the validation set and 20% as the test set. An Attention mechanism is applied to each data source. The output of the Attention mechanism is used, for embeddings'

	construction, as input to a two fully connected layers with sigmoid activation
	function. Then, each visit of each patient is represented by concatenating
	the high level representations sources. It is this representation which is used

  results of the comparison of our MS-LSTMEA based Attention algorithm with different algorithms, regarding the prediction of the drug class

	Approach	Accuracy	Precision	Recall	F1-Measure
	MS-LSTMEA	80.4 ±0.03	78.4 ±0.02	83.9 ±0.2	81.0 ±0.04
	Doctor AI	68.3 ±0.07	69.5 ±0.01	65.2 ±0.02	67.3 ±0.03
	M2vect	65.4 ±0.02	59.5 ±0.07	69.2 ±0.11	64.0 ±0.05
	DeepCare	55.5 ±0.11	54.6 ±0.03	64.9 ±0.01	59.3 ±0.09
	Deepr	51.9 ±0.04	52.2 ±0.02	46.6 ±0.03	49.2 ±0.02
	DeepPatients	58.8 ±0.01	57.2 ±0.03	69.3 ±0.01	62.7 ±0.13

Table 3 :

 3 Comparative results while predicting the best drug class.

Table 4 :

 4 Comparative results while predicting the time before the next medical appointment.

	Approach	Accuracy	Precision	Recall	F1-Measure
	MS-LSTMEA	79.5 ±0.04	78.6 ±0.01	80.9 ±0.01	79.7 ±0.03
	Doctor AI	68.2 ±0.02	66.7 ±0.11	72.2 ±0.04	69.4 ±0.03
	M2vect	57.6 ±0.01	56.2 ±0.01	66.8 ±0.1	61.1 ±0.03
	DeepCare	74.7 ±0.02	73.4 ±0.01	77.1 ±0.05	75.2 ±0.01
	Deepr	59.4 ±0.03	58.2 ±0.05	65.6 ±0.07	61.7 ±0.01
	DeepPatients	62.2 ±0.02	60.3 ±0.03	70.9 ±0.02	65.2 ±0.01

Table 5

 5 shows that Clinical features outperform the Personal and the Comorbidity features when predicting the drug class or the date of the next medical appointment. The prediction accuracy is improved by up to 44.7%

	(s1)	(s2)	(s3)	Accuracy (T)	Accuracy (V)
	0	0	1	22.2 ±0.03	13.7 ±0.04
	0	1	0	49.3 ±0.01	51.8 ±0.03
	0	1	1	64.3 ±0.02	59.1 ±0.05
	1	0	0	24.5 ±0.04	19.4 ±0.03
	1	0	1	36.3 ±0.06	46.2 ±0.02
	1	1	0	67.7 ±0.07	57.1 ±0.01
	1	1	1	80.4 ±0.03	72.3 ±0.14

Table 5 :

 5 Influence of source type (1: included, 0: not included) on the prediction Accuracy.

Table 6 :

 6 Influence of embeddings construction on the prediction results.

Note that the temporal aspect of the visits is not taken into account in the fully connected layers. It is handled by the LSTM.

In the case of regression tasks, we can use a linear activation function in the output

Note that we represent the classes using One-Hot Encoding.

One epoch corresponds to one iteration over the entire training data.

Source: interview with Dr François Teboul, Emergency Doctor, Paris, May 2021

The fact of having three distinct input sources, necessitated to set up seven experiments.

acquisition and the Project administration have been ensured by NAM and FZ. All authors approved the manuscript.

Conflict of interest

The authors declare that they have no conflict of interest.