Dr François Mazuel

Dr Samuel Mathieu

Prof. J.-C Riccardo Di Corato

Jean-Claude Bacri

Thierry Meylheuc

Teresa Pellegrino

Myriam Reffay
email: myriam.reffay@univ-paris-diderot.fr

Claire Wilhelm
email: claire.wilhelm@univ-paris-diderot.fr

Dr M Bacri

Prof C Reffay

Wilhelm

Dr R Di Corato

Prof T Pellegrino

Forced and Self Rotation of Magnetic Nanorods Assembly at the Cell Membrane: a Bio-Magnetic Torsion Pendulum

Keywords: magnetic nanorods, assembly, membranes, rotating magnetic field magnetic nanorods, assembling, membranes, rotating magnetic field

In order to give insights into how anisotropic nano-objects interact with living cell membranes, and possibly self-assemble, we designed magnetic nanorods with average size around 100 nm x 1µm by assembling iron oxide nanocubes within a polymeric matrix under a magnetic field.

We then explored the nano-bio interface at the cell membrane under the influence of a rotating magnetic field. We observed a complex structuration of the nanorods intertwined with the membranes. Unexpectedly, after a magnetic rotating stimulation, the resulting macrorods were able to rotate freely for multiple rotations, revealing the creation of a biomagnetic torsion pendulum.

Table of contents entry

Assembling nanorods at the cell membrane by a remote spinning magnetic field leads to a complex bio-magnetic macrorod constituted of both nanorods and intertwined cell membrane. This new object behaves as a torsion pendulum able to store elastic energy: it is able to freely rotate over multiple turns in response to the spinning stimulation.

In addition to their recognized potential for advancing diagnostics and therapies, nanomagnetic materials were recently described as tools to probe or act on the cell membrane.

In magnetic manipulation of cell membranes, an external rotating magnetic field is used to act at a distance, with precisely controlled intensity, direction and localization of the applied magnetic torque. Magnetically driven rotating or oscillating magnetic nanoparticles could thus be used to probe cells' mechanical properties, [1][2][3] to deliver drugs, [4][5][6][7][8] or to kill cells by physical membrane rupture. [9][10][11][12][13][14][15][16] In this last application, anisotropic nanoparticles in disk-or rod-like shapes, stimulated with low-frequency magnetic fields, could compromise the cell membrane and thereby trigger apoptotic or necrotic programmed cell death. Another promising and recent field of research takes advantage of field-mediated magnetic nanoparticle bioassembly to activate biochemical signaling mechanisms. [17][18][19][20][21][22][23] When localized to the cell membrane, [17][18][19][20][21] the resulting clustering or orientation of targeted receptors can be seen as a nanomagnetic switch to trigger cell responses.

In parallel, chemical synthesis has yielded magnetic nanomaterials with anisotropic geometries, often rod-like shapes. [15,24] Synthesis or self-assembly of magnetic particles under a magnetic field is then frequently used to form super-organized anisotropic structures.

[ [25][26][27][28][29] However, how magnetic anisotropic nano-objects interact and assemble with cell membranes under the action of a remote magnetic field has rarely been investigated. Here, we addressed this field by producing highly magnetic nanorods, by assembling them at the cell membrane into larger macrorods through remote spinning by a rotating magnetic field, and by exploring the macrorod bio-structuration and dynamic response.

Magnetic macrorod formation at the cell membrane

Nanorods were prepared by self-aggregation of iron oxide nanocubes (20-nm edge) embedded within a polymeric matrix (loaded with FITC dye) in the presence of a 0.05 T magnetic field.

The resulting fluorescent magnetic nanorods (Figure 1A1) are (1.1±0.1) µm long and (110±15) nm wide (Figure 1A1) and exhibit strong magnetization (80 emu/g). The nanorods spontaneously attach to the cell membrane (Figure 1A2) and are then internalized when incubated overnight in absence of a magnetic field (Figure 1A3). However, if incubation takes place in a magnetic field (Figures 1B and Figure 1C), internalization is prevented. The nanorods are then forced to assemble at the membrane into larger anisotropic structures aligned in the direction of the magnetic field, which we call macrorods. The formation of such macrorods is a two-step process. First, dipolar magnetic interactions (created by a 0.2-T static vertical magnetic field) forces the assembly of multiple nanorods into thin micron-long elongated clusters (Figure 1B) during their first stage of attachment to the cell membrane. In the second step, cells are subjected to a spinning 0.28-T magnetic field rotating around a 45°C cone, for 30 min at 1 Hz. This conical rotating field is created by two permanent magnets that generate a horizontal 0.2-T magnetic field in-between the magnets and spin around the sample, itself within the vertical 0.2-T magnetic field (see experimental. The thin nanorod clusters then merge, forming larger ellipsoid "macrorods" (average length: 22±11μm and diameter: 4±1μm). SEM and confocal images of typical macrorods are shown in Figures 1C1-3.

First step was to evaluate the impact of the magnetic stimulations on cell viability (Figure 1D).

It demonstrated that neither the static magnetic field nor the rotating one result in any cell damage, as quantified by measuring the cells metabolic activity (relative to untreated control cells) the day after the magnetic stimulations. The nanorods / cell membrane interactions thus do not induce any significant alteration of cells viability and proliferation capacity, even in the rotating setting, contrary to some other studies reporting membrane physical rupture triggered by rotating magnetic nanoparticles, generally targeted to a specific receptor. [9][10][11][12][13][14][15][16] Here it shows that cells can also adapt to a rotating stress and avoid massive harm.

Free rotation

Right after the rotating stimulation, we removed the cells from the magnetic fields (rotating and static) and observed them under the microscope. As expected from the cell viability measurements, cells morphologies were similar to control, with no cell rounding or detachment. However, and surprisingly, some of the macrorods continued rotating, with no external energy source (magnetic or other) to explain this motion (Figure 1D, see also supplementary Movie S1accelerated 5 times). This "spontaneous" rotation reveals in fact the remarkable adaptability of cell membrane materials. Indeed, the macrorods rotated in the direction opposite to the magnetic field rotation to which the rods were exposed. This phenomenon is impressive, as it implies that the cell membranes confer a very efficient elastic energy storage ability on the macrorods. This energy accumulated during the forced magnetic excitation is release by free rotations in the opposite direction when the stimulation is stopped (hundreds of free rotations were sometimes recorded).

Macrorods: intertwined nanorods and membrane filaments

The explanation for this elastic response comes first from the microscopic structure of the macrorods. Macrorods are composed not only of nanorods but are mingled with membrane structures all along their axis. This is clearly illustrated by the typical confocal image of a macrorod in Figure 2A, with the cell membranes labeled with the red Pkh26 membrane marker, and nanorod fluorescence collected in the green channel (FITC). Note that the membrane labels are present all along the macrorod. Looking at the same samples with scanning electron microscopy (SEM, Figure 2B) an intertwining of membrane filaments and nanorods are observed (Figure 2B2). Membrane fragments are also detected on transmission electron microscopy (TEM) images (Figure 2C) of macrorods detached from cells.

The use of the rotating magnetic conical field is decisive to structure the anchor point of the nanorods within the membrane. It is responsible for membrane-rods entanglement and for radial organization of membrane filaments around the macrorod. Indeed, these hybrid membrane/nanorod structures are absent when rods are assembled without rotation (Figure 2B1).

Biomagnetic macrorods at the cell membrane: a torsion pendulum?

To model and quantify the observed free rotation, the rotating stimulation process was explored. A programmable electromagnetic rotating field was used to impose the same rotating cone under the microscope as that generating macrorod formation. Both the magnetic rotating stimulation and the free relaxation could then be video-monitored (supplementary Movie S2accelerated 5 times). Because the field intensity was now lower (30 mT) and thus not sufficient to form the macrorods, they were first formed as previously described with rotating permanent magnets (130 mT). The electromagnetic field was programmed to spin the macrorods counterclockwise, for a given number of rotations (generally 20), at 2 Hz. When the forcing was stopped, and thus in absence of a magnetic field, the macrorods rotated clockwise, some of them up to 20 times, while others rotated only a few degrees (Figure 3).

We propose a torsion pendulum model to describe this free rotation.

The membrane structure at the anchor point during the forced rotation generates an elastic torque C and the rod is subjected to a rotational viscous friction  when it moves. As inertia can be neglected (low Reynolds number ≈ 10 -2 ), the rod's angle  obeys the equation

  R C dt d       
, with R the total of relaxation rotation (in degrees), and can thus be

simply written as                    t = t R exp 1 with C =  
. The fitting in Figure 3A demonstrates the correct matching of this torsion pendulum model. The parameters  and R reflect the peculiar coiling of cellular membranes around and within the macrorod. Figure 3B shows the plot of as a function of R for different rods. Two groups emerge: the first one corresponds to rods relaxing only few tens of degrees (black symbols), while the second corresponds to rods relaxing through more than a thousand degrees (> 3 rotations). For this second group, the rod's anchor point is organized such that it enables extremely efficient energy storage during forced rotations as the membrane filament tension increases. In contrast, for the first group, either the organization is inappropriate or the coiling structure may be damaged during the forced rotation, resulting in an abrupt energy release that would prevent free rotation.

Interestingly, the typical relaxation curve shown in Figure 3A presents a succession of angular jumps that follow the exponential trend described above. This could be explained by different coiling modes of the membrane at the anchor point. Jumps would be associated with twisted membrane filaments, whereas the overall relaxation would correspond to membrane filaments coiled around the rod. The supplementary movie S2 illustrate these two modes.

Within the framework of the torsion pendulum model, successive excitations of a rod is expected to damage the anchor point organization as tension increases and exceeds the membrane filament resistance. This is illustrated by Figure 3C. In this case, rods were subjected to successive periods of magnetic field rotation and immobilization. The intensity of the magnetic field in the immobilized phase is half the intensity of the rotating field, so that a rod with a proper coiling at the anchor point can relax if the elastic restoring force is stronger than the magnetic torque. In the particular example shown in Figure 3C, this is the case for the two first excitations: the rod is able to relax for one rotation despite the magnetic field (indicated with black stars). The increase in tension against the magnetic forcing at the anchor point is materialized in Figure 3C with events (indicated with red crosses) when the rod slows down and stops before being caught up by the magnetic field after an additional rotation. The occurrence of such events increases with successive excitations until the structure at the anchor is damaged (red arrow). After that the rod follows perfectly the magnetic field: the damaged structure cannot store anymore elastic energy efficiently.

Taken together, all those observations confirm the simple torsion pendulum model based on a specific membrane coiling that enable energy storage and restitution.

In summary, by remote assembly of magnetic nanorods during their first stage of attachment to cell membranes through the application of a rotating magnetic field, we managed to create a super-assembly (macrorod) trapping within and wrapping around membrane filaments. As a result, this biomagnetic structure exhibited an elastic behavior which provided an impressive numbers of free rotations, when the rotating field is released and no other stimulation applied.

This movement, evidences a remarkable pool of membranes available and the possibility to arrange them in a biomagnetic torsion pendulum.

Experimental section

Magnetic nanorods preparation

The preparation method is based on iron oxide nanocubes assembly into a polymeric matrix, in presence of a static magnetic field, as adapted by a procedure previously described for spherical magnetic nanoparticles [29]. Briefly, the polymeric anisotropic construct is obtained by mixing nanocubes with the poly(maleic anhydride alt-1 octadecene) polymer, in chloroform. The nanoparticles dispersion was placed in an ultrasound bath under the influence of two opposite permanent magnets placed against each other on the vial wall, generating a 0.05 T magnetic field. The slow and controlled addition of acetonitrile induces a change in the solubility of both the nanocubes and polymer and promotes their aggregation. The application of a permanent magnetic field promoted the formation of elongated superstructure. The magnetic nanorods were finally magnetically sorted and resuspended in water.

The iron oxide nanocubes, about 20 nm in cube-edge, were synthesized by thermal decomposition technique detailed as detailed in [30]. Briefly, they were prepared by mixing 1 mmol of iron(III) acetylacetonate and 4 mmol of decanoic acid in 25 mL of dibenzyl ether, heating the solution to 200°C (5°C/min) for 2.5 h. The temperature is then increased further to reflux temperature (300°C at a rate of 10°C/min) for 1 h. At the end of the process, nanocubes are dispersed in chloroform.

Magnetic field stimulations

Different magnetic devices were designed and fabricated. The simplest one dedicated to macrorods formation (Figure 4A) consists in a permanent magnet (neodymium, 50x40x20 mm, Supermagnet) placed below the cell sample. It creates a magnetic field of 130 mT in the sample region. The second one (Figure 4B) is adapted to this first one by adding a set of two permanent magnets facing one to the other and fixed to a motor controlling a gear wheel (frequency up to 2 Hz). This set of magnets creates a horizontal field of 130 mT, so that the resulting field (about 180 mT) makes an angle of 45° with the vertical axis. As a result, when the motorized magnets are rotated, the magnetic field spins, describing a 45° cone. Both devices contain a cylindrical 20 mm cradle to welcome the cell sample, and are thermostated by water circulation. The third device (Figure 4C) was designed to be adapted to a microscope and to allow switching on and off the magnetic field at will. It is composed of four coils which cores are made of soft iron, connected by pairs and supplied by an alternative current. The space between each core is about 1cm and it creates in between a 30 mT magnetic field. To generate a rotating 45°C rotating cone, this magnetic device is placed 2 cm beneath the cell sample. The magnetic field on the cells is then reduced to 10 mT. To generate the rotating field in the plane of the magnetic coils, the two pairs of coils are supplied with sinusoidal currents (2A amplitude) displaying the same frequency (up to 5 Hz) but 90° out of phase. Finally, the whole set-up was mounted on the vertical arm of a Leica DMIRB with control of the z position, and the microscope was thermostated at 37°C by cube&box (Life Imaging Services).

Nanorods incubation with cells and macrorods formation

PC3 cancer cells were cultured in DMEM supplemented with 10% fetal bovine serum, Lglutamine and penicillin/streptomycin. Prior to experiments, cells were transferred into 20 mm cylindrical home-made cuves with glass bottom which fit within the magnetic fields devices. When cells reached 80-90 % of confluency, nanorods were incubated at an iron concentration of 1 mM, corresponding roughly to 10 6 nanorods per ml, or equivalently to 10 3 nanorods per cell. This concentration was adjusted by testing a range of concentration, observing the cells with confocal microscopy, and selecting concentration where nanorods were numerous on the cells membrane, but still quasi-individual.

After the 2-hour incubation period, the cells were immediately transferred into the thermostated magnetic field devices, first the permanent static field for 30 min, then the spinning magnetic field for 30 min.

Cell viability was assessed by Alamar Blue metabolic assay (Thermo Fisher Scientific).

Fluorescence (appearing post-metabolization of the active ingredient, Resaruzin), was quantified with a microplate reader (excitation 550 nm, detection 590 nm). In brief, 100 000 cells were first seeded in the cuves (3 per conditions). 24 hours after, nanorods incubation was performed (2-hours at 1mMFe, except for the control conditions) and the cuves were submitted either to the static magnetic field only (30 min, condition B in Figure 1), or to both the static and rotating magnetic field (30 min + 30 min, condition C in Figure 1). 24 hours after, all cuves were incubated (800 µl total) with 10 % Alamar Blue in DMEM for 2 hours, and the reagent was transferred to 96-well plate for analysis (200 µl per well). All values are expressed relative to control (normalized at 100% viability).

Imaging

Nanorods and macrorods were imaged by scanning (SEM) and transmission (TEM) electron microscopy, by confocal microscopy, and by conventional transmission microscopy. TEM was used to observe nanorods in aqueous suspension, or macrorods extracted from the cells.

For this second case, the cells were first fixed with paraformaldehyde (2%) for 30 min right after the whole stimulation protocol. Culture medium was then washed and replaced with ultra-pure water to avoid salts contaminations, and cells were scratched in order to detach cells and macrorods. For both nanorods and macrorods TEM observation, a small drop of suspension (5µl) was pipetted onto a copper grid, and observed after total evaporation of the liquid with a Phillips Tecnai 12.

SEM was used to observe whole cells. Cells were fixed with glutaraldehyde (2% in 0.1 M cacodylate buffer), and dehydrated by soaking in a graded series of ethanol before criticalpoint drying under CO2. Samples were mounted on aluminum stubs with conductive silver paint and sputter coated with gold palladium for 200 s at 10 mA. Samples were then imaged with a Hitachi S4500 instrument. The corresponding movie can be seen as supplementary movie S1. D. Cell viability measured by quantifying the metabolic activity (Alamar blue) of cells incubated with the magnetic nanorods in presence of the static magnetic field (condition "static", similar to part B), and under the influence of the rotating magnetic field (condition "rotating", similar to part C), and compared to control cells (seeded at the same exact number of cells, see Experimental Section).

For

  confocal microscopy, cells were labeled either by Pkh26 label (20 min incubation, according to manufacturer's instruction, Sigma) which binds to the plasma membrane, or with phalloidin to see actin filaments, or with DAPI to image the nuclei. Cells were observed by means of an Olympus JX81/ BX61 Device/Yokogawa CSU Device spinning disk microscope (Andor Technology plc, Belfast, Northern Ireland), equipped with a 60x Plan-ApoN oil objective lens. Conventional transmission microscopy was carried out with a DMIRB Leica microscope, equipped with a Cube & Box device to maintain the cells at 37°C during experiment.Macrorods imposed and free rotations were captured with an ultra-fast camera. Time sampling was used in between 100 to 10 images per second. Image J home-made plugins (succession of thresholding, selecting and measuring) were used to analyze the macrorods rotational movements.

Figure 1 (

 1 Figure 1 (2 columns / 16 cm): Assembly process of magnetic nanorods into a macrorod at

Figure 2 ( 1 .

 21 Figure 2 (1.5 column / 12 cm): A. Confocal microscopy image of a macrorod in the red

Figure 3 ( 1 ..

 31 Figure 3 (1.5 column / 12 cm): A. Typical angular evolution of two macrorods after the same

Figure 4 (

 4 Figure 4 (one column / 8 cm): Magnetic set-ups. A. For macrorods formation: a home-made

Acknowledgements

This work was supported by the European Union (ERC-2014-CoG project MaTissE 648779).

We acknowledge the ImagoSeine core facility of the Institut Jacques Monod, member of IBiSA and France-BioImaging (ANR-10-INBS-04) infrastructures.

Received: ((will be filled in by the editorial staff))

Revised: ((will be filled in by the editorial staff)) Published online: ((will be filled in by the editorial staff))


Antoine Avignon
email: avignon.antoine@gmail.com

Pierre Fesler

Claire Daien

David Costa

Marie-Christine Picot

François Roubille

Ariane Sultan

Christine Viarouge-Reunier

Vincent Attalin

Mélanie Badin

Catherine Boegner

Pascal Demoly

Michel Dauzat

Michel David

Béatrice Lognos

Jacques Morel

Jean-Luc Pasquié

Jean Ribstein

Sophie Granier

Bernard Combe

Jacques Mercier

Rodolphe Bourret

Jean Bousquet

Living Lab MACVIA Maladies chroniques

à la diffusion de documents scientifiques 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.

Àl 'origine de 85 % de la mortalité et de plus de 70 % des dépenses de santé en France, les maladies chroniques non transmissibles (MCNT) ont des conséquences médico-économiques considérables. Intimement liées à l'environnement, au mode de vie et au vieillissement, elles s'associent fréquemment chez une même personne conduisant à des situations de multiet comorbidité. Il est impératif que nos systèmes de santé s'adaptent au plus vite pour assurer un vieillissement actif en bonne santé. Le Living Lab MACVIA Maladies chroniques, intégré au réseau des Living Lab MACVIA-LR s'inscrit dans le Partenariat européen d'innovation pour un vieillissement actif et en bonne santé (European Innovation Partnership on Active and Healthy Ageing ; EIP on AHA, DG Sanco et DG CNECT) et a pour objectif de favoriser le dépistage et la prise en charge précoce et intégrative des MCNT et de leurs comorbidités au niveau des soins primaires et de définir un parcours de soins intégré « Maladies Chroniques & Comorbidités ». Il s'appuie fortement sur les TIC-santé, outils connectés et capteurs qui représentent un formidable levier pour approcher la complexité inhérente aux MCNT. Il intègre une web application (plateforme Aviitam ® ), à l'interface du médecin et du patient qui, outre ses fonctions de coordination et d'optimisation des soins, vise à mesurer le niveau de risque des patients de façon à les orienter vers le niveau de soins le plus adapté (1 er , 2 e ou 3 e recours). La plateforme est coordonnée avec les cliniques des comorbidités (Centre ambulatoire des maladies chroniques, CAMAC) qui permettent d'évaluer les patients ayant les index de gravité les plus élevés sur une journée pour proposer une approche thérapeutique coordonnée au niveau des soins primaires. Le Living Lab MACVIA Maladies chroniques s'intègre dans une dynamique régionale portée par la région LR, le CHRU de Montpellier et le CHU de Nîmes, le réseau de soins de santé primaire, un tissu scientifique, les acteurs institutionnels, ainsi que le tissus associatif et industriel. Il propose une approche innovante et créative pour favoriser l'émergence de pistes d'actions possibles permettant d'appréhender la complexité des MCNT.

Le Living Lab dans le contexte international, national et régional

Santé publique

La seconde moitié du XX e siècle a vu s'opérer une mutation profonde dans la santé avec la transition des maladies aiguës (MA) aux MCNT, actuellement responsables dans notre pays de 85 % de la mortalité et de plus de 70 % des dépenses de santé [START_REF] Has | Maladies chroniques : une vision rénovée du parcours de soins[END_REF]. Les conséquences médico-économiques sont considérables avec un impact financier substantiel pour les entreprises [2,[START_REF] Bousquet | System medicine approaches for the definition of complex phenotypes in chronic diseases and ageing. From concept to implementation and policies[END_REF]. Les MCNT sont intimement liées à l'environnement, au mode de vie et au vieillissement, leur nombre allant croissant à mesure que la population vieillit. Avec l'allongement de la durée de vie, l'association de plusieurs MCNT chez une même personne est ainsi de plus en plus fréquente et les situations de multi-et comorbidité [START_REF] Barnett | Epidemiology of multimorbidity and implications for health care, research, and medical education: a cross-sectional study[END_REF] façonnent désormais le quotidien des consultations de médecine générale [START_REF] Bousquet | System medicine approaches for the definition of complex phenotypes in chronic diseases and ageing. From concept to implementation and policies[END_REF][START_REF] Oms | 2008-2013 Action plan for the global strategy for the prevention and control of non communicable diseases. Prevent and control cardiovascular diseases, cancers, chronic respiratory diseases, diabetes. Geneva[END_REF]. La prévention, le dépistage, la prise en charge précoce et intégrée des MCNT sont au coeur des stratégies visant à assurer un vieillissement actif en bonne santé. Sous peine d'effondrement, il est impératif que nos systèmes de santé s'adaptent au plus vite à ces nouveaux challenges pour lesquels nous n'avons pas de solution ni même de stratégie claire et bien définie. Alors que la prise en charge des MA repose sur un modèle mécaniste biomédical, les MCNT nécessitent une approche holistique prenant en compte le mode de vie (alimentation, activités physiques, sommeil, relations familiales et sociales, etc.). Le spécialiste en médecine générale apparaît comme le mieux placé pour saisir de façon unitaire les conditions de multi-et comorbidité, les autres spécialistes ayant par définition une vision segmentaire. Il est également le mieux placé pour appréhender les patients dans leur globalité et intégrer leur mode de vie dans une prise en charge holistique. Il existe de multiples freins à l'adaptation des soins primaires à ce changement de paradigme, les principaux étant : la contrainte médico-économique qui ne permet pas de consacrer le temps nécessaire à l'intégration des données du mode de vie dans la prise en charge, à la mise en place des modifications thérapeutiques du mode de vie Le projet est né des difficultés rencontrées par des praticiens pour la prise en charge ambulatoire des patients ayant un diabète de type 2 et/ou une obésité. Des focus groupes animés par un étudiant en marketing et impliquant des spécialistes en médecine générale ont été organisés de façon à recueillir les attentes de ces derniers dans le cadre de la prise en charge des MCNT. Un premier site prototype a été construit et la preuve de concept validée cliniquement par une étude préliminaire conduite sur 2 ans auprès de 500 patients. Nous avons montré une amélioration très significative de la qualité de vie des patients et une bonne appropriation de l'outil. Le retour des utilisateurs patients a également été pris en compte et une plateforme SaaS (Software as a Service) plus aboutie a alors été développée, regroupant plusieurs outils (encadré 1) au sein de quatre modules :

préparation des consultations (synthèse des données du mode de vie) ; aide à la décision pour le spécialiste en médecine générale (thérapeutiques non médicamenteuses, gradation des soins, orientation vers le niveau de recours adapté) ; implication du patient dans la prise en charge de sa maladie ( ETP, adhésion thérapeutique) ; formation du médecin (développement professionnel continu). Pour optimiser le temps médical, le patient alimente lui-même la plateforme, seul ou à l'aide d'objets connectés. Il est également en lien avec la communauté des patients alors que le médecin est en lien avec la communauté des médecins permettant de générer un savoir communautaire patients et médecins. L'ensemble des activités figure sur une « Time Line » accessible par les professionnels de santé, leurs permettant de trois habitants sur dix (29 %) en France métropolitaine. Certains territoires comme le département de la Lozère sont particulièrement touchés par ce phénomène de vieillissement qui se retrouve couplé à la problématique des déserts médicaux [START_REF]Projection de population à l'horizon 2030 en Languedoc-Roussillon[END_REF]

Objectifs

Le Living Lab MACVIA Maladies chroniques a pour objectif de favoriser le dépistage et la prise en charge précoce et intégrative des MCNT et de leurs comorbidités au niveau des soins primaires et de définir un parcours de soins intégré « Maladies Chroniques & Comorbidités ». Il s'appuie fortement sur les TIC-santé, outils connectés et capteurs qui représentent un formidable levier pour approcher la complexité inhérente aux MCNT. L'innovation portée par ce Living Lab s'attache particulièrement au respect et au renfort de la place centrale des spécialistes en médecine générale dans la coordination des soins et à leur implication dans la recherche, à la création de nouveaux modes de développement professionnel continu plus intégrés à leur pratique quotidienne, à développer des solutions innovantes pour la mise en place des recommandations de santé ainsi qu'à fournir des outils pour une éducation thérapeutique de proximité. Il s'efforce également de promouvoir des approches plus globales en prenant en compte le mode de vie, en promouvant les thérapeutiques non médicamenteuses et en rendant les patients acteurs de leur propre prise en charge. Conformément au concept des Living Lab, une participation accrue et précoce de l'utilisateur profanele patientou professionnelle spécialiste en médecine généraleest mise à profit pour relever ces défis.

Le Living Lab MACVIA Maladies chroniques

Clinique

Plateforme web d'aide à la prise en charge holistique des maladies chroniques pour les médecins et les patients Les consultations en médecine générale sont complexes et doivent traiter dans un court laps de temps de problèmes multiples ayant trait à un large éventail de domaines thérapeutiques [START_REF] Salisbury | The content of general practice consultations: cross-sectional study based on video recordings[END_REF]. L'intégration des comorbidités et des données du de suivre leurs patients. Une « Map » permet au patient de repérer tous les services disponibles dans son environnement (professionnels de santé, salles de sport, groupes de soutien, etc.) pour l'aider dans sa prise en charge. Toutes les données collectées sont stockées dans une base de données sur serveur sécurisé et agréé pour données de santé. L'analyse des données aussi bien médecins que patients (Data Mining) permettra de générer des algorithmes évolutifs adaptés à chaque nouveau patient. À terme, la plateforme web Aviitam ® vise à être un véritable « GPS Santé » (Gestion personnalisée de la santé) adapté aux principales maladies chroniques (figure 1).

Clinique des comorbidités au CHRU de Montpellier

Outre ses fonctions de coordination et d'optimisation des soins, la plateforme web vise à mesurer le niveau de risque des patients ayant une MCNT [START_REF] Fuhrman | Surveillance épidémiologique de la multimorbidité[END_REF], de façon à les orienter vers le niveau de soins le plus adapté (1 er , 2 e ou 3 e recours). Les patients ayant les index de gravité les plus élevés seront orientés vers la cliniques des comorbidités au sein du CHU de Montpellier (Centre ambulatoire des maladies chroniques, CAMAC). À partir de la maladie index (maladie cardiovasculaire, diabète, rhumatisme inflammatoire, maladie respiratoire chronique, etc.), le patient bénéficie d'une évaluation des comorbidités à travers des questionnaires communs aux différentes MCNT ou spécifiques de chacune d'entre elles, ainsi que d'examens complémentaires spécifiques de chaque maladie. La clinique CAMAC est informatisée avec le système du CHU de Montpellier [START_REF] Bourret | An integrated approach to telemonitoring noncommunicable diseases: best practice from the European innovation partnership on active and healthy ageing[END_REF]. À travers ce bilan complet, la clinique des comorbidités permet d'évaluer un patient atteint de MCNT en une journée et de proposer une approche thérapeutique qui sera

Encadr e 1

Principaux outils proposés par la plateforme Aviitam ® Questionnaires d'évaluation psycho-socio-comportementale. Semainiers originaux alimentaires, activités physiques et sommeil. Outil de prescription automatisé des thérapeutiques non médicamenteuses (Vidal ® des « non médicaments »). Recommandations de bonne pratique clinique avec mise à disposition automatique au cours de la consultation. Modules de formation du médecin PENDANT le temps de consultation avec un outil de mise en pratique des formations (réseau de médecins, partage, etc.). Moteur d'aide à la décision en temps réel, à la manière des GPS routiers, capables d'adapter l'itinéraire à la volée en fonction du trafic et des incidents. Module de création d'ordonnances automatisée de thérapeutiques non médicamenteuses en 1 clic (Module One-Click Prescription).

Figure 1

Résumé de la prise en charge à l'aide de la plateforme web Aviitam ® L'augmentation de fréquence de ces comorbidités s'explique par l'utilisation de certains traitements en particulier les corticoïdes, par la présence de facteurs de risque tels que le tabac mais surtout par la présence d'une inflammation systémique chronique. L'accroissement du risque cardiovasculaire et ostéoporotique est également constaté dans d'autres rhumatismes chroniques comme les spondylo-arthrites et notamment le rhumatisme psoriasique [START_REF] Jamnitski | Cardiovascular comorbidities in patients with psoriatic arthritis: a systematic review[END_REF]. Ceci a conduit la Ligue européenne contre les rhumatismes (EULAR) à recommander entre autres une évaluation annuelle des risques cardiovasculaires et une vaccination contre le pneumocoque tous les cinq ans chez les patients atteints de PR [START_REF] Peters | EULAR evidencebased recommendations for cardiovascular risk management in patients with rheumatoid arthritis and other forms of inflammatory arthritis[END_REF]. Un programme de dépistage des comorbidités dans les maladies chroniques a été mis en place dans le Languedoc-Roussillon (MACVIA-LR), dans le cadre du partenariat européen d'innovation pour un vieillissement actif et en bonne santé (EIP-AHA) [START_REF] Bousquet | System medicine approaches for the definition of complex phenotypes in chronic diseases and ageing. From concept to implementation and policies[END_REF][START_REF] Bousquet | The French reference sites of the European Innovation Partnership on active and healthy ageing[END_REF][START_REF] Bousquet | MACVIA-LR Reference Site of the European Innovation Partnership on Active and Healthy Ageing[END_REF]. Il semblait donc important d'y intégrer les rhumatismes inflammatoires chroniques. Ce dépistage vient de débuter en hospitalisation de jour du département de rhumatologie du CHRU de Montpellier. Il est centré sur l'évaluation du rhumatisme, le dépistage des facteurs de risque et atteintes cardiovasculaires, de l'ostéoporose, des anomalies respiratoires, la vérification du statut vaccinal et des dépistages néoplasiques recommandés (figure 3). Les anomalies dépistées concernant le rhumatisme inflammatoire et les troubles osseux sont prises en charge au moment de l'hospitalisation, les anomalies cardiovasculaires et autres sont indiquées au médecin traitant à qui les conclusions du bilan sont envoyées. Les premiers résultats montrent un taux particulièrement élevé de comorbidités non détectées et, en particulier, de risque d'évènement cardiovasculaire létal avec un SCORE (Systematic COronary Risk Evaluation) modifié selon les recommandations de l'EULAR supérieur à 5 % à 10 ans [START_REF] Peters | EULAR evidencebased recommendations for cardiovascular risk management in patients with rheumatoid arthritis and other forms of inflammatory arthritis[END_REF][START_REF] Conroy | Estimation of ten-year risk of fatal cardiovascular disease in Europe: the SCORE project[END_REF].

Coordination avec le milieu associatif et activité physique La sédentarité, ou l'absence d'exercice physique, est un facteur de risque pour les maladies chroniques [START_REF] Booth | Lack of exercise is a major cause of chronic diseases[END_REF] ; l'activité physique représente un élément essentiel de la prise en charge des patients atteints de MCNT et participe au vieillissement actif en bonne santé [START_REF] Pahor | Effect of structured physical activity on prevention of major mobility disability in older adults: the LIFE Study Randomized Clinical Trial[END_REF]. Le succès dans la mise en place des activités physiques nécessite des approches combinées associant des interventions cliniques avec des activités à domicile ou de proximité et éventuellement sur le lieu de travail. Les fédérations et associations sportives (fédérations EPMM [Entraînement physique dans le monde moderne], EPGV [Éducation physique et de gymnastique volontaire], etc.) assurent une offre de pratique d'activité physique adaptée complémentaire sur l'ensemble du territoire régional, au profit d'un public qui n'est pas toujours ciblé par d'autres disciplines, comme les personnes ayant des MCNT. Un travail en collaboration avec le milieu associatif permettra la mise en place d'activité physique par la suite coordonnée au niveau des soins de santé primaire. Une attention particulière est portée à l'accessibilité pour les médecins de ville et aux délais de rendez-vous, ainsi qu'à la production de compte rendus automatisés dans des délais extrêmement courts (figure 2) [START_REF] Bousquet | System medicine approaches for the definition of complex phenotypes in chronic diseases and ageing. From concept to implementation and policies[END_REF]. L'objectif du CAMAC est d'établir la preuve de concept de la modélisation du dépistage des comorbidités chez les personnes ayant une MCNT index et de l'orientation de la prise en charge en favorisant les MTMV. En cas de démonstration de la faisabilité et de l'efficience de ce mode de prise en charge, l'expérience sera répétée au sein des territoires dans le LR, et notamment au niveau des maisons de santé. Le couplage entre la plateforme web et CAMAC permettra le partage des informations et le suivi du patient au plus proche de son domicile pour la mise en place des MTMV et l'adhésion thérapeutique. Cette clinique CAMAC servira de projet pilote à d'autres cliniques dans des services hospitaliers et en dehors.

DeProPass

Une clinique CAMAC mobile a été initiée (DeProPass : dépistage des maladies chroniques de proximité) par le CHU de Nîmes. Les examens seront effectués dans les maisons médicales pluriprofessionnelles des Cévennes.

Clinique des comorbidités en rhumatologie

Les pathologies musculo-squelettiques affectent directement la mobilité, l'autonomie et la qualité de vie de plus de 100 millions d'européens [START_REF] Vos | Years lived with disability (YLDs) for 1160 sequelae of 289 diseases and injuries 1990-2010: a systematic analysis for the Global Burden of Disease Study 2010[END_REF]. Elles augmentent notamment le risque de maladies cardiovasculaires, de diabète et diminuent l'espérance de vie. Parmi les affections articulaires, la polyarthrite rhumatoïde (PR) est le rhumatisme inflammatoire le plus fréquent et le plus sévère. Les patients atteints de PR ont un risque accru d'évènements cardiovasculaires [START_REF] Daien | Rheumatoid arthritis: a cardiovascular disease?[END_REF][START_REF] Peters | EULAR evidencebased recommendations for cardiovascular risk management in patients with rheumatoid arthritis and other forms of inflammatory arthritis[END_REF], d'infections, de pneumopathie interstitielle, de cancers ou encore de fractures ostéoporotiques.

Figure 2 de proximité en complément de l'accompagnement proposé par la plateforme web.

Acteurs impliqués

Associations de patients L'appropriation des outils proposés est un des facteurs clés de la réussite du projet. L'implication des patients est à ce titre une étape essentielle dès la phase initiale de la mise en place du Living Lab et de tous les outils associés. Les principales associations de patients atteints de maladies chroniques du Languedoc-Roussillon seront donc impliquées (Association française des diabétiques, associations de personnes obèses, association BPCO, Association nationale de défense contre l'arthrite rhumatoïde, association française des malades et opérés cardiovasculaires Languedoc-Roussillon, etc.). Liens avec les acteurs institutionnels L'expérience initiale a été présentée à de nombreux acteurs du monde de la santé : spécialistes en médecine générale au cours de « focus group » ; Maisons de santé pluriprofessionnelles, Union régionale des professionnels de santé Languedoc-Roussillon, département de médecine générale de l'UFR de médecine de Montpellier -Nîmes, Agence régionale de la santé Languedoc-Roussillon, Caisse primaire d'assurance maladie de l'Hérault, CHU de Montpellier, coordination des centres spécialisés obésité (coordonnateur national), École centrale. L'accueil a été très favorable, soulignant notamment :

Enseignement

Références

  Les Universités de Montpellier 1, 2, 3 et l'Université de Nîmes sont des centres d'enseignement et de formation reconnus en Europe. De nombreux diplômes d'université sont dispensés dans le cadre des principales MCNT liées à l'environnement (diabète, obésité, etc.). Outre les UFR de médecine et de pharmacie, naturellement impliquées dans la formation autour des maladies chroniques, l'UFR STAPS propose une filière Réhabilitation par les activités physiques adaptées. L'université Montpellier 1 propose un Master 2 Santé Humanité Autonomie avec une spécialité Maladies chroniques et société. Cette formation permet de qualifier des professionnels de la recherche alliant un savoir de type biomédical ou soignant à des perspectives innovantes en sciences humaines et sociales dans la considération des pathologies chroniques. Les universités Montpellier 1 et Montpellier 2 proposent des formations supérieures dans le domaine des Technologie de l'Information et de la Communication appliquée au domaine de la santé. Ces formations permettent de fournir des compétences concernant le suivi à distance de patients souffrant de maladies chroniques et l'assistance des personnes âgées dépendantes.

Figure 3

 3 Figure 3 Évaluation des rhumatismes inflammatoires et de ses comorbidités en rhumatologie BPCO : bronchopneumopathie chronique obstructive évaluée par IPAG (International Primary care Airways Guidelines) ; SAOS : syndrome d'apnée obstructive du sommeil évalué par le score de BERLIN ; risque d'alcoolo-dépendance évalué par le questionnaire CAGE (Cut-down Annoyed Guilty Eye-opener).

  dans la prise en charge est ainsi rendue extrêmement difficile. La plateforme Aviitam ® est une web application à l'interface du médecin et du patient, conçue et développée dans l'objectif de mettre à disposition des spécialistes en médecine générale un outil permettant d'appréhender les MCNT et leurs comorbidités en intégrant les données du mode de vie. Plus précisément, il s'agit d'améliorer la prise en charge des patients ayant une ou plusieurs MCNT associées liées au mode de vie (obésité, anomalie de la glycémie à jeun/diabète, HTA, dyslipidémies, troubles du sommeil, maladie cardiovasculaire, etc.) en :optimisant la relation patient-médecin et la consultation de soins primaire ; formant le médecin en temps réel à la prise en charge des MCNT liées au mode de vie ;

	(MTMV) et à l'éduca-
	tion thérapeutique du patient (ETP), problématique qui devrait
	s'exacerber dans les années à venir avec la diminution de la
	démographie médicale [2] ;
	la formation médicale (initiale et continue) qui souffre d'iner-
	tie et reste basée sur la résolution linéaire de problèmes
	simples -ou compliqués réduits en une somme de problèmes
	simples -dont relèvent les pathologies aiguës alors que les
	MCNT relèvent des problèmes complexes, non linéaires et
	divergents (données et dimensions multiples, démarche
	à inventer, nombreuses solutions possibles) ;
	le manque d'outils efficaces et de moyens de coordination
	entre les professionnels.
	Intégration dans l'EIP on AHA
	Le Living Lab MACVIA Maladies chroniques fait partie du réseau
	des Living Lab MACVIA-LR et s'inscrit parfaitement dans le
	Partenariat européen d'innovation pour un vieillissement actif
	et en bonne santé (European Innovation Partnership on Active
	and Healthy Ageing ; EIP on AHA, DG Sanco et DG CNECT) [6] qui
	vise à mobiliser les acteurs et les ressources en Europe pour
	identifier et éliminer les obstacles persistants à l'innovation en
	faveur d'un vieillissement actif et en bonne santé.
	Le projet est centré sur les technologies de l'information et de la
	communication en santé (TIC-santé) pour préparer les consul-
	tations médicales, former les médecins et mettre en oeuvre
	l'éducation thérapeutique. Il s'intègre donc prioritairement dans
	le plan d'action B3 sur la mise en place et le suivi de modes de
	prise en charge intégrés des maladies chroniques, y compris les
	systèmes de surveillance à distance au niveau régional (Imple-
	menting and tutoring integrated care for chronic diseases,
	including remote monitoring at regional level). En favorisant
	la mise en place des thérapeutiques non médicamenteuses
	(TNM) ainsi que leur adhésion par les patients, il relève égale-
	ment du plan d'action A1 sur l'adhésion aux prescriptions et aux
	programmes médicaux (adherence to prescription and medical
	plans).
	Problématique régionale

L'allongement de l'espérance de vie, associé à l'arrivée aux âges élevés des générations nombreuses du baby-boom en Languedoc-Roussillon (LR), devrait entraîner une augmentation du nombre de personnes âgées de plus de 60 ans de 602 000 en 2005 à plus de 1 million en 2030. Celles-ci représenteraient ainsi près du tiers (31,5 %) de la population totale du LR en 2030, contre près mode de vie favorisant la prescription des TNM et l'adhésion thérapeutique ; participant à l'ETP du patient et en le rendant acteur de sa prise en charge ; créant, des chemins cliniques individualisés, véritables « feuilles de route SANTÉ » à partir des données récupérées et analysées à l'aide d'un outil statistique (Data Mining).

  . La région LR est un acteur d'excellence dans la lutte contre les maladies chroniques grâce aux services cliniques des CHRU de Montpellier et du CHU de Nîmes, au réseau de soins de santé primaire, et au tissu scientifique reconnu au plan national (Inserm, CNRS, LABEX, EQUIPEX, CR2I, SAS) et international. Collaborativement, le CHRU de Montpellier, le CHU de Nîmes et l'université Montpellier 1 ont investi un ensemble de forces hospitalo-universitaires de recherche et de soins afin de développer une stratégie commune de prévention et de lutte contre les maladies chroniques. Le Living Lab MACVIA Maladies chroniques s'intègre dans cette dynamique régionale et propose une approche innovante et créative pour favoriser l'émergence de pistes d'actions possibles permettant d'appréhender la complexité des MCNT.

  Le développement du Living Lab MACVIA Maladies chroniques se fait en partenariat avec des acteurs de l'industrie implantés dans le Languedoc-Roussillon. On peut notamment noter dans les partenaires IBM et BioRad à travers l'Unité mixte de recherche SysDiag. La plateforme web Aviitam ® est développée à travers une startup qui bénéficie du soutien d'incubateurs régionaux tels que le Business & Innovation Centre (BIC) de Montpellier et Languedoc-Roussillon Incubation (LRI). Aviitam ® est également accompagné à travers un programme de maturation de la société d'accélération du transfert de technologies (SATT) du Languedoc-Roussillon, AxLR. Antoine Avignon et Vincent Attalin sont actionnaire de Aviitam SAS ; Vincent Attalin en est le dirigeant. mode d'appropriation des nouvelles technologies dans le champ de la santé de la part des patients et des médecins ; analyses statistiques et Data Mining pour la définition de chemins cliniques ; aspects juridiques. Ces voies de recherche permettent d'envisager une approche transdisciplinaire de la maladie chronique avec l'intervention de discipline aussi variées que la biologie et la médecine (Inserm, Université Montpellier 1 & 2), le marketing (Montpellier Research in Management), les sciences humaines et sociales (Universités Montpellier 1 & 3), l'électronique (Institut d'électronique, Unité mixte de recherche CNRS-UM2) ou encore le droit et la santé (Centre d'étude et de recherche droit et santé, Université Montpellier 1), les statistiques et l'informatique (CNRS, BioRad, IBM).

	le respect et le renfort de la place centrale du MG dans le
	système de soin ;
	le caractère formateur vis-à-vis des MG en cours de
	consultation ;
	la capacité d'impliquer les MG dans la recherche et la promo-
	tion de la constitution de laboratoires territoriaux « vie quo-
	tidienne & santé » ;
	la pertinence économique appuyée sur le Développement
	professionnel continu lui donnant un caractère pérenne ;
	une solution innovante dans la mise en place des recomman-
	dations de santé.
	Liens avec l'industrie

Déclaration d'intérêts :


Thilo Meyer-Brandis

Peter Tankov
email: tankov@math.jussieu.fr

Multi-factor jump-diffusion models of electricity prices *

Keywords: electricity prices, multi-factor models, Lévy-Driven Ornstein-Uhlenbeck type processes, statistical estimation, nonlinear filtering

The recent deregulation of electricity markets has led to the creation of energy exchanges, where the electricity is traded like any other commodity. In this paper, we study the most salient statistical features of electricity prices with a particular attention to the European energy exchanges. These features can be adequately reproduced by the sum-OU model: a model representing the price as a sum of Lévy-driven Ornstein-Uhlenbeck (OU) processes. We present a new method for filtering out the different OU components and develop a statistical procedure for estimating the sum-OU model from data.

Introduction

Similarly to how the end of the Bretton Woods system in 1973 triggered the appearance of currency exchanges, the recent deregulation of electricity markets has led to the creation, in Europe, United States and other countries of a network of energy exchanges, where the electricity is quoted almost as any other commodity. The power prices, which were fixed by governments before competition was introduced, are now determined via an equilibrium of supply and demand, and present a much higher volatility than equity prices. A precise statistical model of electricity spot price behavior is therefore necessary for energy risk management, pricing of electricity-related options and evaluation of production assets.

This paper focuses on statistical modeling of spot electricity prices with a particular attention to the (more or less recently established) European energy exchanges. The electricity prices are usually quoted on hourly basis, however, they cannot be seen as one long hourly series, because in most markets, the delivery prices for all 24 hours of a given day are determined by the system operator simultaneously on the previous day (day-ahead prices). Depending on the market, this may happen between 11AM and 1PM of the previous day. After this, the day-ahead market closes and only small changes can be made to delivery contracts using the adjustment market in case of emergency such as plant outage. There is thus no causality relationship between different hourly prices on the same day, and the hourly electricity data should rather be seen as 24 dependent daily series than one long hourly series. Given the structure of the day-ahead electricity market, a possible model for hourly prices would be

X h t = Y t f (t, h) + ε h t ,
where • X h t , h ∈ {1, . . . , 24} is the price for day t and hour h;

• Y t is the common factor (average daily price);

• f (t, h) is a slowly varying daily pattern;

• ε h t is a white noise.

A principal component analysis on the Powernext hourly data shows that the first component corresponding to a constant daily pattern f (t, h) = f (h) explains 70% of variance. Therefore, most of the variability of electricity prices, as well as all interesting statistical features like mean reversion are contained in the average daily price series, and the daily pattern is mostly related to seasonality. In this paper, we therefore concentrate on modeling the daily electricity price series. The rest of this paper is structured as follows. In Section 2 we elucidate the qualitative features of electricity prices common for different data sets that we study. Section 3 discusses different models for electricity spot prices found in the literature and how well they describe the stylized features of the data. At the end of this section, we introduce our model of choice: the sum-OU model, where the spot price is described as sum of several non-Gaussian Ornstein-Uhlenbeck processes. This type of model was introduced in [START_REF] Benth | A non-Gaussian Ornstein-Uhlenbeck process for electricity spot price modeling and derivatives pricing[END_REF], and it both captures well the stylized facts of electricity spot prices and is analytically tractable for derivatives pricing and risk management. However, the important question of how to estimate and implement the model has not yet been addressed in the literature. The rest of this paper is therefore dedicated to the estimation of this model from data: Section 4 shows how different components can be separated, and Section 5 presents a detailed case study where the sum-OU model is estimated completely from the EEX spot price series.

While the main objective of this work is to provide an estimation procedure for the multi factor non-Gaussian OU model, the original contribution in this paper is threefold:

• We discuss the stylized features of electricity prices some of which are being reported here for the first time. In addition, we focus on quotes from European energy exchanges which are much less studied than the American ones (except maybe the Nord pool market) and present qualitatively different characteristics.

• We propose a novel method for estimating the sum-OU model from price data, based on recent findings in nonparametric statistics.

• We give a complete data-driven specification of the model and present a detailed case study where all parameters of this specification are estimated from EEX time series.

The principal application of our model lays in simulating market evolution scenarios for computing risk measures, testing hedging strategies and evaluating investment policies.

Stylized features of electricity prices

This section discusses qualitative features of electricity prices common for different data sets that we study. Weekend effect As shown in Table 1, the distribution of electricity index returns from Friday to the following Monday is significantly different from the distribution of intraweek daily returns for all six series under study. First of all, the variances of weekend returns are 1.2 to 2 times larger than those of intraweek returns. This is explained by the fact that unlike stocks, electricity is traded on the weekends, although in much smaller volumes (see Fig. 1). Secondly, for all six series under study, the expected Friday to Monday return is reliably positive whereas the intraweek returns are negative or indistinguishable from zero. This is behavior is opposite to what has been observed for equity returns [START_REF] French | Stock returns and the weekend effect[END_REF][START_REF] Schwert | Anomalies and market efficiency[END_REF]. While weekends introduce a lot of seasonality into the price series they don't account for interesting statistical features (e.g., no spikes during weekends). In order to make deseasoning easier and to concentrate on statistical aspects, we consider time series without weekends.

Description of data sets

Series Affine trend R Yearly seasonality Many studies (see for instance [START_REF] Geman | Understanding the fine structure of electricity prices[END_REF]) have established the presence of annual seasonality in electricity spot prices. Following [START_REF] Geman | Understanding the fine structure of electricity prices[END_REF], we remove the seasonality by representing it as a linear combination of sinusoids with periods 12 and 6 months. We estimate the seasonality by fitting a 6parameter trend function

f (t) = a + bt + c 1 sin(2πt) + c 2 cos(2πt) + d 1 sin(4πt) + d 2 cos(4πt) (1)
to the log-price series by ordinary least squares. A linear trend and both harmonics are statistically significant at a 1% level in all 6 six series under study. However as seen from the Table 2, the seasonality contributes very little to the overall variability of the series. The proportion of variance explained by seasonality is highest for the Nord Pool series (whose huge share of hydro energy highly depends on weather cycles) and lowest for the COB series. These two series are represented in Figure 2. From now on, we concentrate on the deseasonalized series, that is, on the residuals of the least squares regression.

Stationarity and mean reversion Contrary to stock prices, the electricity prices tend to exhibit stationary behavior. For all six series under study, the unit root hypothesis was confidently rejected by the Dickey-Fuller test, both on the deseasonalized data and on the non-deseasonalized one (in the latter case, a variant of the Dickey-Fuller test with a linear trend was used). The test statistics are reproduced in Table 3.

The electricity prices, similarly to other commodities, are therefore mean reverting to a trend which may exhibit slow stochastic variations as in Figure 3 Spikes A fundamental property of electricity spot prices, already observed by many authors, and present in all 6 time series studied in this paper, is the presence of spikes, that is, rapid upward price moves followed by a quick return to about the same level. Typical spike profile for Nord Pool price series is shown in Figure 4. During the peak period, the price process has different properties. In particular, the rate of mean reversion is much higher than during normal evolution. The presence of spikes is a fundamental feature of electricity prices due to the non-storable nature of this commodity and any relevant spot price model must take this feature into account.

The presence of long-term stochastic variation in the mean level of electricity prices makes it difficult to establish a range of prices, for which the price process is in peak mode: what is considered as peak level now, may become normal in 3 years. The rate of mean reversion is therefore determined not only by the current price level but also by the previous evolution of the price, which suggests that the behavior of electricity spot prices may be non-markovian. As observed in [START_REF] Geman | Understanding the fine structure of electricity prices[END_REF], the intensity of spikes can exhibit yearly seasonal patterns. We will come back to this issue in Section 4 devoted to spike detection.

Multiscale autocorrelation

The examination of autocorrelation structure of the six electricity price series (Figure 5) reveals the presence of two very distinct types of behavior: the autocorrelations of COB, MEAD and NP series decay very slowly, while the autocorrelation function of APX, EEX, and UKPX series present a two-scale behavior with a quickly decaying component and a slow one. This striking difference may be due to differences in the organization of electricity markets [START_REF] Geman | Commodities and Commodity Derivatives[END_REF]: the APX, EEX and UKPX markets are organized as power exchanges where the reported price is simply the price of the last transaction. On the other hand, the price index we use for the Nord Pool market is the so-called system price resulting from the bids of all market participants, and the COB and MEAD price series are indices compiled by Platts as average of all transactions reported by market participants [START_REF]Platts methodology and specifications guide[END_REF]. This suggests that APX, EEX and UKPX prices will have a higher volatility and faster mean reversion than COB, MEAD and NP, which is exactly what is observed in Figure 5 and Table 4. Let us now focus on the study of the multiscale autocorrelation structure in APX, EEX and UKPX time series. The observed structure is described quite precisely with a sum of two exponentials:

ρ(h) = w 1 e -h/λ1 + w 2 e -h/λ2 .
Figure 6 and Table 5 show the quality of the fit and the fitted parameters. This type of correlation structure arises in a model where the price is a sum of two independent mean-reverting components with two different rates of mean reversion. We conjecture that in this context, the fast rate of mean reversion, corresponding to correlation length between 3 and 5 days corresponds to the "spike mode" while the second slower rate of mean reversion, with correlation length of about 100 days, corresponds to normal price evolution. This idea will be further developed in the subsequent sections.

Non-gaussian distribution of returns

The descriptive statistics of return distributions (skewness and kurtosis) shown in Table 6 highlight the non-gaussian character of electricity spot returns. This is explained by the presence of spikes in the data since these low-probability large-amplitude moves clearly cannot arise in a Gaussian framework. First of all, all six distributions are strongly leptokurtic, with excess kurtosis values ranging from 9 to 29. For comparison, the excess kurtosis of daily returns of the S&P 500 index is about 3. Secondly all six series are slightly positively skewed. This is related to the fact that the downward price move at the end of a spike is not quite as fast as the upward move in the beginning and can last several days. Figure 7 plotting the logarithm of the empirical distribution function for the right and left tail of EEX return distribution shows that these tails are fatter than exponential. In case of exponential tail decay, the empirical distribution function

F n (x) = 1 n n i=1 1 Xi≤x satisfies F n (x) ∼ e -λ1|x| as x → -∞ and 1 -F n (x) ∼ e -λ2x as x → ∞.
Therefore log F n (x) must behave like a linear function when x → -∞ and similarly for the right tail. However, this is not what is observed in Figure 7: the graph is convex and we conclude that both tails are fatter than exponential. The same conclusion holds for all 6 series that we study.

To detect the possible power-law behavior in the tails of return distribution we use the Hill estimator. To fix the ideas, let us concentrate on the upper tail of the distribution and suppose that the distribution function satisfies

1 -F (x) ∼ x -α L(x), x → ∞, α > 0,
where L is a function satisfying lim t→∞ L(tx)/L(t) = 1 for all x (slowly varying function). Then the extreme value index γ = α -1 can be estimated with the Hill estimator defined by

H k,n := 1 k k i=1 log X (i) X (k+1) ,
where (X (i) ) are order statistics. This estimator is consistent under a broad set of assumptions (see [START_REF] De Haan | How to make a Hill plot[END_REF] for a review). In practice, if H k,n remains stable over a large range of k, this is an indication of power-like tail behavior.

Looking at the data (see Figure 8), we see that for all six series under study, in the range between k = 30 and k = 150, the Hill plot stabilizes around the value between 0.3 and 0.5 (for both tails). This suggests that the tail behavior of electricity prices may be well described by a power law with decay parameter 2 < α < 3. Electricity returns have thus finite second moment which justifies the use of autocorrelation functions above in this section.

Literature review

In this section we will give a short overview of the main types of continuous time reduced form spot price models from the existing literature. We mention that there also exist a variety of econometric time series models (see e.g. [START_REF] Escribano | Modelling electricity prices: International evidence[END_REF][START_REF] Knittel | An empirical examination of deregulated electricity prices[END_REF][START_REF] Leon | Testing for weekly seasonal unit roots in the spanish power pool[END_REF][START_REF] Misiorek | Point and interval forecasting of spot electricity prices: Linear vs. non-linear time series models[END_REF] and references therein). However, we don't go further into this line of research here.

Structural models

In structural or equilibrium models [START_REF] Barlow | A diffusion model for electricity prices[END_REF][START_REF] Kanamura | A structural model for electricity prices with spikes: Measurement of jump risk and optimal policies for hydropower plant operation[END_REF], one tries to mimic the price formation in electricity market as a balance of supply and demand. The (very inelastic) demand for electricity is described by a stochastic process:

D t = D t + X t , dX t = (µ -λX t )dt + σdW t ,
where D t describes the seasonal component and X t corresponds to the stationary stochastic part. The price is obtained by matching the demand level with a deterministic supply function which must be nonlinear to account for the presence of price spikes. Barlow [START_REF] Barlow | A diffusion model for electricity prices[END_REF] proposes

P t = a 0 -D t b 0 1/α
for some α > 0 while Kanamura and Ohashi [START_REF] Kanamura | A structural model for electricity prices with spikes: Measurement of jump risk and optimal policies for hydropower plant operation[END_REF] suggest a "hockey stick" profile

P t = (a 1 + b 1 D t )1 Dt≤D0 + (a 1 + b 1 D t )1 Dt>D0 .
A typical spot price trajectory in the latter model is shown in Figure 9.

While structural models reproduce price spikes on a qualitative level, the assumption of deterministic supply function is probably too restrictive. It implies in particular that spikes can only be caused by surges in demand, while in electricity markets spikes can also be due to sudden changes in supply, such as plant outages. In addition, prices in such models do not have the multiscale autocorrelation property observed in the data, and, with a deterministic supply function, do not allow for a stochastic base level as in Figure 3.

Markov models To make estimation easier, Geman and Roncoroni [START_REF] Geman | Understanding the fine structure of electricity prices[END_REF] model the electricity log-price as a one-factor Markov jump diffusion.

dP t = θ(µ t -P t )dt + σdW t + h(t)dJ t
The spikes are introduced by making the jump direction and intensity leveldependent: if the price is high, the jump intensity is high and downward jumps are more likely, whereas if the price is low, jumps are rare and upward-directed. This approach produces realistic trajectories and reproduces the seasonal intensity patterns observed in American price series, but once again, the process reverts to a deterministic mean level rather than the stochastic pre-spike value.

Regime-switching models In the one-factor Markov specification of Geman and Roncoroni [START_REF] Geman | Understanding the fine structure of electricity prices[END_REF], the 'spike regime' is distinguished from the 'base regime' by a deterministic threshold on the price process: if the price is higher than a given value, the process is in the 'spike regime' otherwise it is in the 'base regime'. This threshold value may be difficult to calibrate and it is not very realistic to suppose that it is determined in advance. The regime-switching models of Deng [START_REF] Deng | Stochastic models of energy commodity prices and their applications: mean-reversion with jumps and spikes[END_REF] and Weron [START_REF] Weron | Heavy tails and electricity prices[END_REF] alleviate this problem by introducing a two state unobservable markov chain which determines the transition from "base regime" to "spike regime" with greater volatility and faster mean reversion:

dP t = θ 1 (µ t -P t ) + σ 1 dW t (base regime) dP t = θ 2 (µ t -P t ) + σ 2 dW t (spike regime).
This nonlinear model is more difficult to estimate than a one-factor Markov specification (filtering must be used) and once again, the stochastic base level is not taken into account since in the spike regime, the process quickly reverts to the seasonal mean and not to the pre-spike base level.

Multifactor models To allow mean reversion to a stochastic base level, several authors have suggested to model the price as sum of several factors. In the simplest case of a two-factor model, the first factor corresponds to the base signal with a slow mean reversion, and the second factor represents the spikes and has a high rate of mean reversion. The first models of this type (see [START_REF] Barlow | Calibration of multifactor models in electricity markets[END_REF] and references therein) reprensented the electricity price (or log-price) as sum of Gaussian Ornstein-Uhlenbeck processes. However, while the first factor can in principle be a Gaussian process, the second factor is close to zero most of the time (when there is no spike) and takes very high values during a spike. This behavior is hard to describe with a Gaussian process, and several authors [START_REF] Benth | A non-Gaussian Ornstein-Uhlenbeck process for electricity spot price modeling and derivatives pricing[END_REF][START_REF] Deng | Lévy process-driven mean-reverting electricity price model: the marginal distribution analysis[END_REF][START_REF] Villaplana | Pricing power derivatives: A two-factor jump-diffusion approach[END_REF] have therefore explored a jump-diffusion specification.

In Benth, Kallsen and Meyer-Brandis [START_REF] Benth | A non-Gaussian Ornstein-Uhlenbeck process for electricity spot price modeling and derivatives pricing[END_REF], the deseasonalized spot price is a sum of independent Lévy-driven Ornstein-Uhlenbeck components:

X(t) = n i=1 Y i (t) (2)
dY i (t) = -λ -1 i Y i (t)dt + dL i (t) , Y i (0) = y i , (3)
and processes L i (t) are independent, possibly time inhomogeneous Lévy processes with E[L 2 i (1)] < ∞. In this model, the autocorrelation function ρ t (h) defined as the correlation between X(t) and X(t -h) satisfies

ρ t (h) := n i=1 w i (t)e -h/λi n i=1 w i (t)
In the stationary case (i.e. the underlying processes are time homogeneous Lévy processes), we have simply

ρ(h) := n i=1 w i e -h/λi n i=1 w i with w i = λ i Var L i (1).
This modeling approach thus allows to explain the multiscale autocorrelation phenomenon observed in European spot price series in Section 2

Figure 10 shows a typical price trajectory in this model, with two factors, one of which is driven by a Brownian motion and another one by a compound Poisson process (with Gaussian jump sizes).

This model, coupled with a good description of price seasonality, provides a precise characterization of electricity spot price behavior, especially well suited to European energy markets. In addition, due to its arithmetic structure, it is analytically tractable when it comes to futures and other derivatives pricing. In the remaining part of this paper we will discuss the question of fitting this model to data.

Detecting and filtering the spikes

For the remaining parts of the paper we will deal with the question of how to estimate a multi-factor model as presented in (2)-(3). In particular, we will focus on the following two factor specification for the deseasonalized spot price (eventually log-spot price) X(t):

X(t) = Y 1 (t) + Y 2 (t) ( 4
)
dY i (t) = -λ -1 i Y i (t)dt + dL i (t) , i = 1, 2 . (5)
Here Y 2 (t), the Ornstein-Uhlenbeck component responsible for spikes, is driven by a Lévy processes L 2 (t) with possibly time varying jump intensity. The Ornstein-Uhlenbeck component Y 1 (t), responsible for the remaining price variations, is assumed to be Gaussian, i.e. driven by a Brownian motion L 1 (t). We recall that the Ornstein-Uhlenbeck components (5) are explicitly given by

Y i (t) = e -(t-s)/λi Y i (s) + t s e -(t-r)/λi dL i (r) , s < t , i = 1, 2 . (6)
In this section we will discuss the methods to identify and filter out the spikes in the time series of spot prices. First of all, since our two factor framework is non-Markovian, this exercise is necessary to estimate the parameters for the two components of the model. Especially when spike behavior (frequency and/or size) appears in a non-stationary seasonal way (as for example seems to be the case for the Nord pool series), effective estimation of the model requires the separation of the non-stationary spikes path from the stationary path of the remaining daily variations. Further, besides serving estimation purposes, the identification of the individual factors is essential to apply the model to futures and other derivatives pricing (see Benth et al. [3]).

In the literature, most of the existing spot price models driven by processes including jumps are of jump diffusion type. When estimating a model of this type, the difficulty is to determine which price variations are caused by jumps and which ones are caused by the continuous part of the process. The easiest and most common way to deal with this problem is to fix a threshold (according to some criteria) above which price variations are considered to be caused by jumps.

For our purposes, these threshold methods are not appropriate. The reason is that they are designed to filter out the path of jumps whereas we need to identify the spike component Y 2 (t), i.e., we need to filter out the path of complete spikes (jumps followed by rapid mean reversion).

One potential method to filter out the spikes is stochastic filtering. The aim of this procedure is to estimate the spike component at time t given the price observations X(s), s ≤ t, until time t via the optimal filter

E[f (Y 2 (t))|F X t ]
, where F X t denotes the filtration generated by the observations X(t). In the case where both factors are Gaussian, the solution is explicitly given by the well known Kalman filter, which was used for calibrating Gaussian multifactor models of electricity prices in [START_REF] Barlow | Calibration of multifactor models in electricity markets[END_REF]. In our framework, since the spike component obviously has to be non-Gaussian, there is no explicit solution and Monte Carlo methods (particle filters) have to be employed in order to approximate the filter. However, when implementing the particle filter in our situation one encounters several difficulties, in particular:

• the filter is not easy to design when parameters are unknown;

• rare events such as spikes lead to sample impoverishment;

• sequential filtering makes less sense when the complete series is available for estimation.

Because of the above objections, we don't pursue the particle filter approach any further. The approach that we want to focus on in more detail in this section is to filter out spikes by adopting methods from nonparametric statistics. In this setting, the general idea is to consider the spike path as deterministic data given through a function f (t) which is disturbed by Gaussian (autoregressive) noise given through the base component Y 1 (t):

X(t) = Y 1 (t) + f (t) .
We will assume the spikes path is of the form as proposed by our model, i.e. f (t) ∈ OU where OU is the space of all Ornstein-Uhlenbeck functions

f (t) = M i=1 α i 1 t≥τi e -(t-τi)/λ2 .
Here, M ∈ N denotes the number of spikes, α i ∈ R are the spike sizes, and τ i ∈ R are the starting times of the spikes in the spike path f (t).

Let N be the sample size of our price series and X(t j ), j = 1, ..., N , the deseasonalized (log) price observations, equally spaced in time (one day). For notational convenience we introduce for a function g(t) the notation g(j) := g(t j ). Then, according to our model,

(X(j) -f (j)) = Y 1 (j), j = 1, ..., N ,
is a sequence of correlated Gaussian random variables. In order to pass to a sequence of iid Gaussian random variables suited for estimation purposes we can consider the transformation (∆X(j) -∆f (j)), j = 2, ..., N ,

where we define for a function g(t)

∆g(j) = g(j) -e -1/λ1 g(j -1).

Then

∆X(j) = Y 1 (j) + f (j) -e -1/λ1 (Y 1 (j -1) + f (j -1)) = ∆Y 1 (j) + ∆f (j) ,
and from the explicit solution of an Ornstein-Uhlenbeck process we get that ∆Y 1 (j) is a sequence of iid Gaussian random variables. Taking [START_REF] Escribano | Modelling electricity prices: International evidence[END_REF] as a starting point, we propose in the following two methods to filter out the spikes path f (t) from our data, i.e. to determine the spike number M and spike sizes and times (α i , τ i ): hard thresholding and the adapted Potts filter.

Hard thresholding Given the number M of spikes and assuming the two mean reversion coefficients λ 1 and λ 2 to be known (estimated from the autocorrelation function), we know that the maximum likelihood estimator (α * i , τ * i ) of (α i , τ i ), i=1,...,M, is given by least squares estimation, i.e.

(α * i , τ * i ) = arg inf (αi,τi) N j=1 (∆X(j) -∆f (j)) 2 .
However, the complexity of this optimization problem is of order N M which makes it unfeasible for implementation. Instead, the hard thresholding method suggests to approximate the problem by placing the spikes one by one. To this end suppose we want to optimally place one spike given through

f (t) = α1 t≥τ e -(t-τ )/λ2 .
In this case, the maximum likelihood estimator

(α * , τ * ) = arg inf (α,τ ) N j=1 (∆X(j) -∆f (j)) 2
can be computed as

α * = ∆X(j)∆f (j) ∆f (j) 2 , τ * = arg sup ( ∆X(j)∆f (j)) 2 ∆f (j) 2 .
It follows that the algorithm in order to place M spikes optimally one by one can be formulated as follows: In contrast to the global maximum likelihood estimation, the hard thresholding algorithm is of order M N and thus well suited for implementation. The procedure stops, that is the optimal number M of spikes is found, when the spike size becomes small according to some criteria (see the following section for an example of such a criteria).

1. Set X (1) = X and X (m+1) = X (m) -f (m) for m = 1, ...M . 2. For every m = 1, ...M , the m'th spike f (m) (t) = α m 1 t≥τ m e -(t-τ m )/λ2 is determined by α m = ∆X (m) (j)∆f (j) ∆f (j) 2 , τ m = arg sup ∆X (m) (j)∆f (j) 2 ∆f (j) 2 .
In Figure 11 the hard thresholding method is tested on some simulated data, and in Figure 12 we applied hard thresholding to Powernext data. The procedure appears satisfactory, however, two remarks have to be made. Typically, the performance of the filter improves when a spike correlation length λ 2 is chosen that is shorter than the one obtained from the autocorrelation function. For Powernext data the autocorrelation function yields λ 2 ∼ 2.5 -4, whereas λ 2 ∼ 1 -2 is a better choice for the filter. This indicates that λ 2 captures some additional mean reversion not belonging to spike behavior when fitting the autocorrelation function. As for the first correlation length λ 1 , the filter behaves in a rather robust way in the range λ 1 = 10 -100. Secondly, we see that the hard thresholding algorithm works better for non-interacting spikes. When there is a cluster of spikes closely following each other, hard thresholding encounters problems in separating the spikes. In this situations, the adapted Potts filter presented in the following might be a better choice. [START_REF] Winkler | Smoothers for discontinuous signals[END_REF] the authors introduce a Potts filter which estimates a piecewise constant signal from noisy data by penalized least squares. The idea we propose in this paragraph is to adapt the Potts filter from [START_REF] Winkler | Smoothers for discontinuous signals[END_REF] in order to detect spikes rather than piecewise constant signals disturbed by some noise. In the situation of spike detection the penalized least squares optimization problem takes the following form:

Adapted Potts filter In Winkler and Liebscher

H(f * ) = arg inf f ∈OU   γ|{t : f t -= f t }| + N j=1 (∆X(j) -∆f (j)) 2   . (8)
Here f t -denotes the left limit of f in t (where we set f 0 -= 0) and |{t : f t -= f t }| is the number of spikes in the path of f . We thus want to find the spike path f * that is closest to our data in the least squares sense (after accounting for autoregressive noise) where, however, every spike is penalized by a weight γ.

In order to solve problem [START_REF] French | Stock returns and the weekend effect[END_REF] we first consider the case of uncorrelated noise which corresponds to λ 1 = 0. In this situation problem [START_REF] French | Stock returns and the weekend effect[END_REF] takes the form

H(f * ) = arg inf f ∈OU   γ|{t : f t -= f t }| + N j=1 (X(j) -f (j)) 2   , (9)
which can be explicitly solved by employing a dynamic programming argument.

For this purpose we denote by

B(r) := arg inf f ∈OU   γ|{t ≤ n : f t -= f t }| + t≤r (X(j) -f (j)) 2  
the solution of the problem on the interval [1, r] (such that B(N ) = H(f * )). Further, we denote the cost of placing a spike on an interval [r, s], 1 ≤ r ≤ s ≤ N , as

H [r,s] := H [r,s] (µ) := γ + s j=r X(j) -µe -(j-r)/λ2 2 ,
and respectively we denote as

H [0,s] = s j=1 X(j) 2
the cost if now spike has been placed on [1, s]. Note here that µ does not denote the spike size (which above was denoted by α) but the level f (r) of the spike path at time r (i.e. spike size α plus f (r -)). Then, with the same arguments as in ( [START_REF] Winkler | Smoothers for discontinuous signals[END_REF]), we can represent B(n), 1 ≤ n ≤ N , recursively by employing dynamic programming:

B(n) = min 0≤r≤n B(r -1) + min µ∈R H [r,n] (µ) . (10)
Here we set B(-1) := 0 and B(0) := 0. The algorithm to solve problem (9) can then be formulated as follows:

1. For all 1 ≤ r ≤ s ≤ N determine µ * In the general case of autoregressive noise, i.e. λ 1 ∈ (0, ∞), problem (8) can no longer be solved explicitly with the above algorithm. In this situation the cost of placing a spike on an interval [r, s], 1 ≤ r ≤ s ≤ N , takes the form

H [r,s] (µ) = γ + s j=r (∆X(j) -∆f (j)) 2 = ∆X(r) -µ -e -1/λ1 f (r -1) 2 + s j=r+1 ∆X(j) -µ∆e -(j-r)/λ2 2 . ( 11
)
We see that in order to solve for an optimal level µ on the interval [r, s] we need to know the spike path level f (r -1) at time r -1. In other words we can no longer place spikes optimally using local considerations only on the interval [r, s], and therefore the dynamic programming argument is not applicable. Instead we propose to approximate the problem by substituting the optimal cost function H

[r,s] (µ * [r,s]
) in 1. of the above algorithm with ). Right: spike detection using the hard thresholding procedure. The procedure is applied to the deseasonalized price series directly (rather than its logarithm), and only part of the data series is shown.

K [r,s] = γ + s j=r+1 ∆X(j) -X(r)∆e -(j-r)/λ2 2 . ( 12
Step 1 (Deseasonalizing the data): We fit the seasonality function (1) to logarithmic spot prices (multiplicative seasonality). This is carried out as described in Section 2. The deseasonalized spot price series is

X des t = X t e f (t) .
The seasonal component of the EEX series is plotted in Figure 14, left graph.

Step 2 (Separation of components): The spike component in the deseasonalized spot price series can be identified by using one of the two methods described in the previous section (hard thresholding or Potts filter). Each of these methods has 3 tuning parameters: the correlation lengths λ 1 and λ 2 , and a third parameter which determines how many spikes we want to remove.

Estimation from the autocorrelation function described in Section 2 gives λ 1 = 99.4 and λ 2 = 2.9 for the EEX series (this is different from the numbers given in Table 5 because the latter were computed for the logarithm of the spot price). As explained in Section 4, in practice, the hard thresholding algorithm is insensitive to the value of λ 1 in the range of 10 -100 and relatively sensitive to the value of λ 2 , so that best results (the smallest residual volatility for the same number of spikes removed) are obtained with λ 2 between 1 and 2.

The number of spikes to remove can be determined by fixing some "target noise level" and then removing spikes with the hard thresholding method, or decreasing the penalty parameter γ in the Potts filter until the standard deviation of returns of the residual signal (the signal from which spikes have been removed) reaches this target noise level. The target noise level represents our a priori knowledge of the noise level of the base signal. One possible way to do assess this quantity is to remove some percentage ǫ of returns with highest absolute value, and compute the standard deviation of the remaining returns. The idea is that ǫ biggest returns correspond to spikes. This procedure is not sensitive to outliers (occasional very big returns).

For the actual estimation, we used the hard thresholding method taking λ 1 = 100, λ 2 = 1. To compute the target noise level, the percentage of returns to be removed was set to ǫ = 5%. This gives the target standard deviation σ * = 0.154 (the standard deviation computed using all available returns is 0.38). With this criterion, the adapted Potts filter algorithm detected 115 spikes in a series of 6.5 years with the penalty parameter γ = 0.13. The hard thresholding algorithm on the other hand needed to place only 62 jumps to achieve the prescribed standard deviation. Despite this, the residual signal for the Potts filter had a higher excess kurtosis than that for the hard thresholding algorithm, suggesting that in the former case some spikes/jump persist after the filtering procedure. We therefore concentrate on the hard thresholding algorithm for the rest of this study. Figure 14, right graph, illustrates the spike detection using the hard thresholding procedure.

Step 3 (Distribution of the base signal): In the case of EEX series, the increments of the base signal have a skewness of -0.008 (down from 2.62 for the initial deseasonalized series including spikes) and an excess kurtosis of 1.05 (down from 113). We conclude that the base signal is sufficiently close to Gaussian to be modeled by an AR(1) process with Gaussian increments. The parameters of this process can be estimated by maximum likelihood, which gives e -1/λ1 = 0.85 ⇒ λ 1 = 6.2 and σ = 0.148. The difference between this value and the value of λ 1 estimated from the autocorrelation function may be explained by the presence of not two but three exponential components in the autocorrelation function. Indeed, fitting the autocorrelation function with three exponentials, we obtain λ 2 = 1.69, λ 1 = 14.4 plus a third factor with a very long decay rate. However, for the sake of simplicity, we concentrate on the two-factor model.

Step 4 (Distribution of the spike signal): We need to estimate the spike intensity and the distribution of the spike amplitude. Since we don't observe any seasonal behavior in spike occurrence on EEX data we assume a constant spike intensity denoted by λ. Then λ is determined as λ = number of spike detected total number of data points , which gives λ = 0.037 in the case of the EEX series. To model the spike size distribution, we choose the Pareto law with distribution function

P (Z > z) = z z 0 -α .
There are several reasons for this choice:

• Since the spike amplitudes are from bounded below by a threshold, we need a distribution supported on the interval [z 0 , ∞).

• The number of data points for estimation is limited (62 spikes in the case of the EEX series) therefore we cannot estimate a distribution with many parameters.

• Figure 15 shows that the empirical CDF of the spike amplitude distribution is very close to a straight line on a log-log scale, which is an indication of a power-law decay.

For the EEX series, the threshold z 0 was taken equal to 0.45 (smallest spike amplitude), and the decay rate α is equal to 1.44 (this was found by fitting a straight line to the empirical CDF graph in the log-log scale).

In summary, we have the following model for the EEX price: α. This model specification can be used to simulate price trajectories. Figure 16 confirms that indeed the model mimics the path behavior of the spot price series very well.

X(t) = e f (t) (Y 1 (t) + Y 2 (t)) dY 1 (t) = -λ -1 1 Y 1 (t) + σdW (t) dY 2 (t) = -λ -

Figure 2 :

 2 Figure 2: Left: Nord Pool series (strongest seasonal effect). Right: COB series (weakest seasonality).

Figure

  Figure 3: UKPX time series (after the seasonality has been removed) and its 6-month moving average.

Figure 5 :

 5 Figure 5: Autocorrelation functions of electricity price series.

Figure 6 :

 6 Figure 6: Fitting the autocorrelation function with a sum of two exponentials. From left to right: APX, EEX, UKPX.

Figure 7 :

 7 Figure 7: Tails of the empirical distribution function of EEX returns. The left graph shows log F n (x) and the right graph log(1 -F n (x))

Figure 8 :Figure 9 :

 89 Figure 8: Hill plots for MEAD (left) and UKPX (right) series.

Figure 10 :

 10 Figure 10: Sample spot price trajectory in a two-factor non-Gaussian Ornstein-Uhlenbeck model.

Figure 11 :

 11 Figure 11: Performance of hard thresholding on simulated data.

Figure 12 :

 12 Figure 12: Performance on Powernext data. Left: in the log-price series both positive and negative spikes are present. Right: in the price series mostly positive spikes are present.

  [r,s] = arg min µ H [r,s] (µ) and H [r,s] (µ * [r,s] ).Determine further H [0,s] and set µ * [0,s] = 0 for all 1 ≤ s ≤ N .2. Determine B(n) for all 1 ≤ n ≤ N , keeping track of (at least one) r * n giving the minimum in[START_REF] Geman | Understanding the fine structure of electricity prices[END_REF].

3 .

 3 Determine recursively through N, r * N , r * r * N , ... the spike number and spike times of the optimal path f * . The optimal spike sizes are then implied by the corresponding µ * [r * N ,N ] , µ * [r * r * N ,r * N ] , ....

Figure 13 :Figure 14 :

 1314 Figure 13: Comparison of the adapted Potts filter and the hard thresholding algorithm for spike detection.

1 2 Y 2 - 1 Figure 15 :

 22115 Figure 15: Empirical CDF of spike size in log-scale.

Figure 16 :

 16 Figure 16: Comparison of the real EEX series and the simulated series with estimated parameters. Which is which? (The market series is the lower one).

Table 1 :

 1 Expectations and variances of intraweek and Friday to Monday returns. The standard deviations of the estimators, computed by bootstrap, are given in brackets.

	FM var.

Our study uses the following six sets of daily data (weekends removed), provided by Datastream. The last date in all series is 21/11/2006.

• Dow-Jones California-Oregon border firm on-peak electricity price index (COB) starting 19/05/1997. • Dow-Jones Mead/Marketplace firm on-peak electricity price index (MEAD) starting 08/08/2000. • European Power Exchange Phelix Base electricity price index (EEX) starting 16/06/2000. • Amsterdam Power Exchange daily average (APX) starting 26/05/1999. • United Kingdom Power Exchange daily average (UKPX) starting 27/03/2001. • Nord Pool system price daily average (NP) starting 1/1/1994. For some illustrations we also use the Powernext (French energy exchange) daily data downloaded from Powernext website. Series Intraweek exp. ret. FM exp. ret. Intraweek var. Figure 1: Daily trade volume, in MWh, on the EEX market, from 24/09/2006 to 10/11/2006. Bright-colored bars correspond to weekends. Data downloaded from EEX website at www.eex.de.

Table 2 :

 2 R 2 (proportion of initial variance explained by the trend function) with seasonality and without it (in this case, the trend function is affine).

			2 Seasonal trend R 2
	COB	0.0990	0.1031
	MEAD	0.0274	0.0385
	EEX	0.5160	0.5373
	APX	0.1378	0.1715
	UKPX	0.5615	0.5900
	NP	0.2504	0.3106

Table 3 :

 3 Dickey-Fuller unit root test statistics for the six price series under study.

	3: UKPX time series (after the seasonality has been removed) and its
	6-month moving average.
	Series	APX COB EEX MEAD NP UKPX
	Standard DF test on desea-	485	59.5	505	43.4	71.9	304
	sonalized series
	DF test with linear trend on	466	59.4	483	43.2	66.4	284
	raw series
	For all six series the non-stationarity hypothesis is rejected at 1% confidence
	level. For the standard test, the non-stationarity hypothesis is rejected at 1%
	level if the statistics is greater than 13.8. For the test with linear trend, it is
	rejected if the statistics is greater than 29.5.

Table 4 :

 4 Annualized volatilities of electricity prices.

	0.8
	0.6
	0.4
	0.2
	0.0
	-0.2
	-0.4
	-0.6
	-0.8
	Dec 1999	Jan 2000	Feb 2000
	APX COB EEX MEAD	NP	UKPX
	513% 231% 363%	197%	133%	282%

Figure 4: Spikes in the Nord Pool price series (after logarithmic transformation and seasonal trend removal).

Table 5 :

 5 Fitting the autocorrelation function with a sum of two exponentials.

		APX EEX UKPX
	w 1 0.300 0.161	0.468
	λ 1	81.5	94.3	112.6
	w 2 0.688 0.761	0.483
	λ 2	2.9	4.3	3.7

Table 6 :

 6 Descriptive statistics of electricity return distributions.

		Series				APX COB EEX MEAD	NP	UKPX
		Skewness			0.11 0.163 0.641	0.025		0.431	0.801
		Excess kurtosis 14.8	15.4	12.7	9.5		29.0	15.0
	-1									0
		Log of empirical CDF									Log of empirical CDF
	-2	LS fit								-2			LS fit
	-3
										-4
	-4
	-5									-6
	-6
										-8
	-7
										-10
	-8
	-9									-12
	-1.6	-1.4	-1.2	-1.0	-0.8	-0.6	-0.4	-0.2	0.0	0.0	0.5	1.0	1.5	2.0	2.5

* This research project was supported by the Europlace Institute of Finance

This amounts to assuming that the size of the optimally placed spike on the interval [r, s] is equal to ∆X(r) (all the variation is explained by the spike) and then spike function at time r is equal to X(r). Obviously, the cost of placing spike of this type can then be computed locally. Also we set

for the cost of not placing a spike. As before we then introduce the recurrence relation

Note that in this case B(N ) is not the solution of problem ( 8) but of the corresponding problem where instead of optimizing over paths f ∈ OU we optimize over paths resulting from the type of spikes in [START_REF] Knittel | An empirical examination of deregulated electricity prices[END_REF]. We use this optimization procedure in order to determine the spike number and times (but not spike sizes) by the following algorithm: Given the spike number M and spike times τ 1 , ..., τ M obtained from the above algorithm, the optimal spike sizes α * 1 , ...α * M are then determined by global least squares optimization:

In Figure 13 the performance of the autoregressive Potts filter is compared to the hard thresholding method on Powernext data. We see that both filters deliver comparable and acceptable results.

Case study: Modeling EEX spot prices

We conclude this paper by estimating the two-factor Ornstein-Uhlenbeck model defined by ( 4)-( 5) on the EEX time series. We choose to estimate our model on the spot price rather than on its logarithm to avoid introducing artificial downward jumps by the logarithmic transformation. We discuss both model specification and parameter estimation in four steps.