Robin Recours
email: robin.recours@univ-montp1.fr

François Nick Trujillo

A E Franc

Aussaguel Ae

Metal Music and Mental Health in France

Keywords: Music, Rock metal music, Mental health, Anxiety, Depression

teaching and research institutions in France or abroad, or from public or private research centers.

Introduction

Performed by bands with evocative names such as Megadeth, Black Sabbath and Suicidal Tendencies, ''metal music'' is characterized by references to death, suicide and Satanism. Critics argue that this type of music has negative impacts on the health and well-being of people who listen to it, especially on young fans. Defenders counter that metal music is a form of entertainment similar to horror films and does not have harmful consequences. This study examines the potential relationship of metal music to the mental health of fans.

A Brief History of Metal Music

Metal music began in the 1960s and 1970s, primarily in the United Kingdom and the United States (see [START_REF] Binder | Constructing Racial Rhetoric: Media Depictions of Harm in Heavy Metal and Rap Music[END_REF][START_REF] Campion | the face of death[END_REF][START_REF] Scheel | Heavy Metal Music and Adolescent Suicidality: An Empirical Investigation[END_REF][START_REF] Wass | Adolescents and Destructive Themes in Rock Music: A Follow-up[END_REF][START_REF] Weinstein | Heavy Metal: A Cultural Sociology[END_REF]. Although the roots of this music were blues and the first hard rock bands, Cream followed by Led Zeppelin, early metal bands like Black Sabbath, Judas Priest, Iron Maiden and Motorhead radicalized blues and rock 'n roll by incorporating distortion, hard beats and explosive sounds into the music and dark images into the lyrics. Subsequent metal bands fused elements of punk, hip hop and even classical music into their songs with an even harder edge, as in the works of Megadeth, Iron Maiden, The Scorpions and Evanescence. This music gained in popularity with the introduction of new bands with new sounds, as well as with popular films such as Sam Dunn's Metal: A Headbanger's Journey.

Many variations of metal music have been introduced over the last few decades, including ''thrash, '' ''death'' and ''black'' metal (see Christe 2007;[START_REF] Hein | Hard Rock, Heavy Metal, Metal. Histoire[END_REF]. Developed in the early 1980s, thrash metal music is known for its fast and aggressive beat and critical lyrics about social issues such as drugs, suicide, globalization and consumption, as illustrated in the work of Anthrax, Metallica and Slayer. Developed in the mid to late 1980s, death metal music is characterized by melodic distortion, harsh guttural vocals and morbid lyrics, as performed by bands such as Death, Morbid Angels and Possessed. Developed in the late 1980s and early 1990s, black metal music is known for its shrieking vocals, experimentation and anti-Christian and Satanic symbolism, as in the music of Barzum, Mayhem and Venom. Metal music continues to evolve into various fragmented forms.

Critics and Defenders of Metal Music

Critics of metal music argue that it has negative impacts on society, pointing to media accounts about band members and metal music fans who, individually and/or in cults, have desecrated graves and churches, engaged in satanic rituals and committed suicide and murder (see [START_REF] Balducci | Adorateurs du diable et rock satanique[END_REF]Bobineau 2004;[START_REF] Domergue | Culture barock et gothic flamboyant, un e ´cho surgi des abı ˆmes: la musique extre ˆme[END_REF][START_REF] Regimbald | Le Rock n'Roll, viol de la conscience par les messages subliminaux[END_REF]. For example, Marilyn Manson, a superstar of metal music, chose his name because of the reference to serial killer Charles Manson, and he has used images of burning crosses and the number 666 in his music and videos. The lead singer from the band Korn was reported to have mental health problems. In the 1990s, the lead singer of the Norwegian band Mayhem, known for throwing animal heads in public, committed suicide, and the guitarist of the band was brutally murdered. The media blamed the ''diabolical'' nature of the band and of the heavy metal counter culture as being responsible for these deaths [START_REF] Campion | the face of death[END_REF].

Parent groups and government authorities are especially concerned because children and adolescents have easy access to metal music, lyrics, images and cults through Internet Web sites, file sharing, and other means that are not under the control of parents or officials. Not surprisingly, then, the American Academy of Child and Adolescent Psychiatry has recommended that parents pay close attention to children who listen to metal music because of the potential negative effects (see [START_REF] Alessi | The Influence of Music and Rock Videos[END_REF].

In France, the situation is somewhat unique because, unlike the United States, cults are outlawed by the Constitution. (For example, the Church of Scientology is prohibited in France.) No country in the world has a religious diversity as extensive as that found in the United States. This religious diversity stems in significant part from the First Amendment to the United States Constitution, which guarantees that all Americans have the freedom to practice the religion(s) of their choice. Recent editions of J. Gordon Melton's Encyclopedia of American Religions lists some 2,300 religious groups active in the United States; most of them have at least 2,000 members. In France, 172 religious groups are officially registered. French officials are particularly concerned about Satanic cults related to metal music. In fact, some black metal music as well as earlier forms have been described as ''Catholic rock'' because of its images of sin, ghosts and, of course, the devil. As French journalist Robert Belleret (2007) put it with respect to metal music, ''It is on the Internet that today the Devil spins its web'' (p. 3). Thus, the French government's Interministral Mission of Vigilance and Struggle against the Influence of Cults (La Mivilude 2004) has warned parents to limit their children's exposure to metal music and their access to Internet sites associated with this music. Some Catholic groups have also protested during metal music concerts and festivals, such as Slayer's 2006 concert in Paris and the Hellfest Festival.

In opposition to these critics are various members of metal music fan clubs and musicians who argue that listeners of metal music can distinguish between fiction and reality and that the rituals associated with metal music concerts are fantasies that rarely have serious consequences to the health of the participants. In addition, these groups say that the impact of metal music may be weaker in France because the majority of lyrics are in English, even for French bands. Thus, fans and musicians of metal music believe that the criticism of reporters, politicians and other opponents is unfounded.

Given the controversy surrounding metal music, we were interested in determining if there is a link between the immersion in metal music culture and the mental health of the individuals who participate. Is an individual who listens to metal music more likely to have serious mental health problems? These are empirical questions, and they form the basis of this study.

Variables and Research Questions

Generally speaking, we agree with [START_REF] Satcher | Natural and Manmade Disasters and Mental Health[END_REF], who define mental health as ''the successful performance of mental function, resulting in productive activities, fulfilling relationships, the ability to adapt to change, and successful coping with adversity'' (p. 2540). Although there are several dimensions of mental health, we chose to focus on the anxiety and depression of metal music fans because of their exposure to images of death, suicide and Satanism in this music and culture.

Dependent Variables

Anxiety and Depression

The two dependent variables in the study were anxiety and depression. Conceptually, anxiety is defined as a feeling of apprehension and fear characterized by physical symptoms such as sweating and palpitations and/or feelings of stress, restlessness and tension. Depression is defined as a more serious illness that affects the physical and emotional well-being of a person and may be characterized by symptoms such as loss of appetite and sleep, feelings of despair and possible suicidal tendencies, negative thoughts and the notion of the uselessness of life. We wondered if fans of metal music subjected to its dark symbolism would exhibit greater levels of anxiety and depression.

Independent Variables

Metal music fans differ according to their preferences for musical styles, the type of clothing they wear, their age, and many other factors [START_REF] Mombelet | Les concerts de musique metal: approche anthropologique[END_REF]. Accordingly, we took into account variables that we believed might mediate the anxiety and depression of fans of metal music.

Age, Number of Years as a Fan and Musical Preference

Metal music has already had three generations of musicians and listeners, and we do not know if the first generation of metal music fans shares a common culture with more recent generations. Thus, we identified as independent variables the age of the participant, the number of years the participant has identified with metal music culture and the type of metal music preferred (thrash, death, black).

Gender

In general, males tend to score higher on scales that measure sensation-seeking [START_REF] Zuckerman | Dimensions of Sensation Seeking[END_REF]. Inasmuch as metal music is a leisure activity that is associated with high levels of stimulation (e.g., loud music, very physical dancing), we identified gender as an independent variable.

Employment Status, Education and Field of Study

Some studies have shown that employment status and education are strong markers of human behavior [START_REF] Baudelot | L'Ecole capitaliste en France[END_REF][START_REF] Berthelot | Le pie `ge scolaire[END_REF][START_REF] Boudon | Pierre 1984 Distinction: A Social Critique of the Judgement of Taste[END_REF]. We broke employment status into three categories based on Bourdieu (1984): ''high'' (factory owners, managers, professionals), ''medium'' (teachers, business employees, office workers) and ''low'' (factory workers, welfare recipients). In addition, an individual's field of study (e.g., literary and arts, scientific, economic, professional, technical) has been shown to influence one's perspective as well (Dubet 1991). We believed that each of these factors might contribute to the mental health of metal music fans, and thus, we identified them as independent variables.

Body Modification and Concert Attendance and Participation

Metal music fans and performers create a unique subculture through various symbolic images and performances [START_REF] Canetti | Masse et puissance[END_REF][START_REF] Walser | Running with the Devil. Power, Gender and Madness in Heavy Metal Music[END_REF]. For example, many fans of metal music engage in body modification practices, ranging from discreet tattooing to severe scarring. In addition, fans differ in terms of the number of concerts they attend each year and the type of performances they engage in at concerts. For example, the most common forms of participation include ''pogoing'' (pushing and jumping into each other), ''slamming'' (passing a person overhead), ''stage diving'' (leaping from the stage into the crowd), ''moshing'' (hitting or smashing into each other with high impact) and ''Bravehearting'' (when fans separate into two lines and rush at each other, as in the film Braveheart). For this reason, we identified as independent variables types of body modification and concert attendance and participation.

Motivation for Participation

Metal music fans attend concerts and participate in metal music rituals for different reasons, such as to listen to their favorite bands, to enjoy the atmosphere of concerts and to drink and/or use drugs, among others. Thus, we included motivation as an independent variable.

Research Questions

R1 Do metal music fans in France exhibit great levels of anxiety and depression? R2 What variables mediate the levels of anxiety and depression for metal music fans?

Methods

Instrument

We assessed a group of metal fans' mental health with the Hospital Anxiety and Depression Scale (HADS). This scale was developed by [START_REF] Zigmond | The Hospital Anxiety and Depression Scale[END_REF] to detect anxiety and depressive states, the most common form of psychological disturbance in hospitals and in the general population. The HADS instrument has been used widely by researchers to study various populations (e.g., patients with cancer and heart problems, medical students, the elderly).

The HADS instrument consists of 14 items, each using a 4-point scale from 0 to 3, and it is simple to use. In addition, because of the subtle nature of the questions, the subject taking the instrument is usually unaware that it is measuring anxiety or depression. This was a distinct advantage in this study because we believe that fans of metal music are sensitive to the prejudice against them and, therefore, might be nervous or defensive in taking a questionnaire about their mental health.

The HADS instrument is considered to be effective in terms of factor structure, intercorrelation, homogeneity, and reliability (see [START_REF] Herrmann | International Experiences with the Hospital Anxiety and Depression Scale: A Review of Validation Data and Clinical Results[END_REF][START_REF] Mykletun | Hospital Anxiety and Depression (HAD) Scale: Factor Structure, Item Analyses and Internal Consistency in a Large Population[END_REF]. In most countries where this scale has been administered, factor analysis has yielded a model of two factors, although a third factor has emerged in some studies of populations with good health [START_REF] Brandberg | Anxiety and Depressive Symptoms at Different Stages of Malignant Melanoma[END_REF]). After orthogonal (varimax) rotation, factor correlations with the noncorresponding subscales are low (\0.40) indicating that the subscales can be conceptualized as two theoretically independent dimensions. Most studies have shown satisfactory item-total correlations within the two subscales. Internal consistency has been acceptable, with Cronbach's a between 0.80 and 0.93 for anxiety and between 0.81 and 0.90 for depression. Retest reliability has shown a high correlation (r [ 0.80) after up to 2 weeks, although this correlation decreases after longer time intervals. Generally, the scores of the HADS instrument are better distributed (in terms of normality) than the State-Trait Anxiety Inventory (STAI) and the Beck Depression Inventory (BDI). Due to its simplicity and neutrality, most studies have had response levels between 95 and 100 percent, even in sensitive populations (see [START_REF] Herrmann | International Experiences with the Hospital Anxiety and Depression Scale: A Review of Validation Data and Clinical Results[END_REF].

Although it is not used to make definitive diagnoses, the HADS instrument can reveal pathological levels of anxiety and depression. To identify subjects likely to be anxious or depressed, the ordinal scale score must be dichotomized at a particular cutoff point. There is no single, generally accepted cutoff score for the HADS. In their original study, Zigmond and Snaith recommended two cutoff scores for both subscales: 7/8 for possible and 10/11 for probable anxiety or depression. In subsequent studies, the same authors proposed a cutoff score of 14 for ''severe'' mental disorders [START_REF] Snaith | HADS: Hospital Anxiety and Depression Scale[END_REF].

In this study, we followed research on the quality of life in France [START_REF] Boyer | Les aidants familiaux des personnes myopathes: liaisons entre charge de soins, qualite ´de vie et sante ´mentale[END_REF]) and chose the conservative cutoff score of 11 for each dimension (anxiety score C11, probable anxiety; depression score C11, probable depression). A study using the same scale and the same cutoffs [START_REF] Herrmann | HADS-D: Hospital Anxiety and Depression Scale-Deutsche Version; Ein Fragebogen zur Erfassung von Angst und Depressivita ¨t in der somatischen Medizin[END_REF] showed that an abnormal HADS score was observed in 10 percent of subjects with good mental health (n = 152) and 71 percent with psychiatric patients (n = 69). In another study with people suffering from depression (n = 2669), 91 percent were found to have abnormal HADS scores [START_REF] Friedman | Three-Dimensional Structure of the Hospital Anxiety and Depression Scale in a Large French Primary Care Population Suffering from Major Depression[END_REF].

We used all 14 items of the HADS instrument, translated and validated in French (Razavi et al. 1989). The scale was divided into two subscales with alternate items: anxiety (items 1, 3, 5, 7, 9, 11, 13) and depression (items 2, 4, 6, 8, 10, 12, 14). In addition, we included items to measure each of the independent variables, including age, gender, employment status, education and field of study, body modification, motivation and number of and type of participation at metal concerts.

The survey was composed of three parts: the first dealt with the demographic characteristics of the subjects, the second dealt with their motivation and participation in metal music and the third was composed of 14 items from the HADS instrument.

Pretest of Instrument

We expected that adolescents and young adults would comprise a majority of our sample of metal music fans. Therefore, we tested the instrument on a group of 20 adolescents in school to make sure that all the questions were easily understood and accepted by a young population. After administration of the questionnaire, subjects were interviewed in order to isolate and modify any problematic items or questions. This new version was then administered to the sample.

Sample

We needed a sample of metal music fans to complete the HADS instrument. Our first inclination was to administer the instrument at places where fans gather: at concerts or at practice areas. However, we believed that group dynamics or concert euphoria might bias the results. We also considered administering the survey to fans of metal music that we could identify on the streets (through clothing such as teeshirts with morbid themes and studded necklaces) but realized that we would miss metal music fans who do not/did not wear such attire.

Thus, we administered the instrument via the Internet by posting it on the largest French Internet forums dedicated to metal music. This decision had two distinct advantages. First, it permitted us to contact metal music fans individually, as participation in Internet forums is usually conducted by individuals in front of a computer screen at home or in a cybercafe ´. As such, the results are not biased by group dynamics. Second, using the Internet permitted us to reach fans with high levels of interest in metal music, because those who are most interested are more likely to visit and participate in Internet Web sites and forums devoted to metal music. We believed that a sample of metal music fans who are most devoted to this type of music was the best sample to assess the potential effects of this music.

We posted the instrument on the 10 Internet sites for a two-week period. Questionnaires with incoherent or incomplete responses were removed from the analyses. This could eliminate some of the ill, we concede, as the most depressed would be least likely to fill out a questionnaire.

Statistical Analyses

We conducted a factor analysis to verify the validity of the scale. We used multiple linear regression analysis and ANOVA to asses the effects of the independent variables. Unless otherwise indicated, significance levels for all analyses was the p \ 0.05 level.

Results

Among the 333 metal enthusiasts who responded to the survey, 9 gave partial or incoherent responses. We eliminated these individuals from the study before conducting the analyses, which resulted in a 97.3 percent response rate.

Descriptive Results

The sample consisted of 282 males (87.85 percent) and 39 females (12.15 percent). The mean age was 22.67 years old (SD = 5.29). Approximately one-half (49.53 percent) were students and 41.74 percent were employed (8.73 percent were both employed and students). The sample was diversified with respect to other variables (see Table 1).

The average number of years in metal music culture was 9.22 years (SD = 5.77, maximum = 30 years). The average number of concerts each year was 15.65. However, there was a large disparity between the participants: SD = 18.60, maximum = 100 concerts. A little more than one-third of the respondents (35.20 percent) had permanent body modifications (e.g., tatoos, piercings, scarifications) whereas 5.30 percent had a combination of two types of body modifications (e.g., tatoo and piercing).

The ritual most widely performed during concerts (Table 2) was ''pogo'' (M = 1.86, SD = 0.69), whereas the least widely performed was ''stage diving'' (M = 1.21, SD = 0.47). Almost half of the sample (44.55 percent) belonged to a metal rock band or club. Respondents expressed varying preferences for different types of metal music, including ''death metal'' (37.69 percent), ''black metal'' (22.74 percent) and ''thrash metal'' (18.07 percent).

The motivations for metal concert participation were also diverse (see Fig. 1). In the survey, respondents were given the option to select multiple motivations for attending concerts. The majority of the reasons for attending concerts included to listen to music they enjoy (95.93 percent), for the ambiance (84.62 percent), to party (66.97 percent) and/or to unwind (55.66 percent); a smaller percentage of reasons included to drink a lot (33.94 percent), while a very small percentage included to take drugs (4.98 percent), to sell drugs (2.72 percent) and/or to fight (0.91 percent).

Factor Analysis

Factor analysis yielded a three-dimensional structure that explained 47.60 percent of the variance. A three-dimensional structure was not unexpected since we had no evidence a priori that our sample would have high levels of anxiety or depression. Following [START_REF] Friedman | Three-Dimensional Structure of the Hospital Anxiety and Depression Scale in a Large French Primary Care Population Suffering from Major Depression[END_REF], we performed an orthogonal rotation (see Table 3) to test the two-factor model (depression and anxiety). These two factors derived from the HADS factor analysis explained 38.71 percent of the variance (anxiety explaining 27.55 percent of the variance, eigenvalue = 3.86; depression capturing 11.17 percent of the variance, eigenvalue = 1.56). The reliability factors using Cronbach's a were 0.70 for anxiety and 0.67 for depression.

In general, the results indicated that the respondents had low levels of anxiety and depression. The average scores (anxiety, M = 7.26, SD = 3.44; depression, M = 3.76, SD = 2.91) were far below the cutoff levels used to determine mental disorders (score C11). However, as in all populations, there were some individuals who scored above this cutoff; specifically 15.58 percent of the subjects had a high probability of anxiety (score C11) and 3.43 percent of them had a high probability of depression (score C11). So-called ''severe'' subjects (as measured by Zigmond and Snaith) constituted \5 percent for anxiety (n = 16) and\1 percent for depression (n = 2).

Multiple Regression Analyses

The multiple regression analyses revealed that gender, age, number of concerts per year, preferred musical style and body modification practices of tattooing and piercing had no link to mental health, in terms of either anxiety or depression. On the other hand, these same analyses revealed that grade level and employment status had links to mental health. The results indicated that as the grade level of the respondents goes up, the anxiety level of respondents goes down. The least anxious respondents were those who had a master's degree, whereas the most anxious respondents were those with only a middle school or a high school degree (Newman-Keuls, p \ 0.05). Students coming from literary and arts backgrounds were found to be more anxious than those of scientific backgrounds [t(91) = 2.36, p \ 0.05]. The results also showed that the lower the employment status, the higher the level of depression and of anxiety. Lower classes of workers were more anxious and depressed than members of the middle class (Newman-Keuls, p \ 0.05), whereas middle-class individuals were more anxious and depressed than upper-class individuals (p \ 0.05).

The analyses also revealed that writing lyrics for a band, consuming high levels of alcohol during concerts and the specific body modification practice of scarification have links to mental health. Specifically, fans of metal music who write song lyrics scored significantly higher than others on the subscales of anxiety [t(320) = 2.75, p \ 0.05] and depression [t(320) = 2.46, p \ 0.05]. Metal music enthusiasts who drink a lot of alcohol at concerts scored significantly higher than the others on the subscale of anxiety [t(210) = 3.24, p \ 0.05] but not on the subscale of depression. In addition, metal music fans who had scarifications scored significantly higher than the others on the depression subscale [t(320) = 2.38, p \ 0.05] but not on the anxiety subscale. However, the higher scores on the subscales of depression or anxiety do not necessarily mean that these metal music fans are abnormally anxious or depressed because the scores were below the danger cutoff (score C11).

A multiple linear regression analysis (see Table 4) for concert rituals showed that ''stage diving,'' ''moshing'' and ''Bravehearting'' were not significantly linked to mental health. This same analysis showed, however, that individuals who practiced ''pogoing'' were significantly less anxious than others and individuals who practiced ''slamming'' were significantly more anxious than the others. Nevertheless, the R 2 remained very low, meaning that if there is a significant relation between these kinds of rituals and mental health, the relation is very weak.

R = 0.25, R 2 = 0.06 R = 0.

Discussion

Our sample consisted of mostly male subjects, ages 13-44 that came from a wide range of education and employment statuses. In this regard, although the sample may appear to be biased toward males, it is very representative of the culture of metal music as described by the scientific literature in France [START_REF] Mombelet | Les concerts de musique metal: approche anthropologique[END_REF] and in other Western countries [START_REF] Arnett | Metalheads: Heavy Metal Music and Adolescent Alienation[END_REF][START_REF] Walser | Running with the Devil. Power, Gender and Madness in Heavy Metal Music[END_REF]. In addition, because this sample consisted of individuals who were more likely to be serious fans of metal music, they were more likely to be influenced by this music and culture. The use of the HADS instrument proved to be effective for this sample in that 97.3 percent of respondents filled it out completely. However, the factor analysis of the instrument yielded three factors instead of two (anxiety and depression), although three factors also emerged in previous studies of populations in good health [START_REF] Brandberg | Anxiety and Depressive Symptoms at Different Stages of Malignant Melanoma[END_REF].

The results indicate that fans of metal music are in good health with respect to anxiety and depression. The percentage of anxious and depressed individuals in the sample was low, and very few were judged to be severely pathological. These results indicate that, contrary to critics who suggest that images of death and destruction in metal music have harmful consequences, the mental health of fans of this type of music is generally good.

However, this study provided evidence that certain variables have effects on the mental health of metal music fans. For example, with respect to variables external to metal music culture, the study showed some links between the educational and employment status of metal music fans and their mental health. We believe, however, that individuals who are less educated and less established in their career will feel more anxiety and depression in general, regardless of the type of music they listen to (see [START_REF] Santelmann | Formation professionnelle et emplois non qualifie ´s[END_REF]. These individuals have less qualified/skilled jobs and may see themselves as valued less in social life than those in more qualified/skilled professions (see [START_REF] Baudrillard | La socie ´te ´de consommation, ses mythes, ses structures[END_REF]Bourdieu 1984).

With respect to the influence of a literary and arts background on anxiety scores, some researchers suggest that literary and arts schools and universities are sometimes perceived as less professional than scientific and engineering schools (Dubet 1991), so any associated levels of anxiety could come from the difficulty (real or presumed) experienced by of individuals with literary or arts backgrounds who must enter the working world. It is necessary to conduct further research to determine if this result would be found in a larger and more heterogenous sample.

With respect to the influence of gender, the small sample size of females prevented a meaningful analysis. It is unclear if female fans of metal music are more like other fans (both female and male) or like other women in general. Additional research on this question should be undertaken.

Concerning variables internal to metal music culture, the results show that several factors had little or no influence on mental health, including preferred style of music, number of concerts attended each year and whether fans have tatoos or piercings. However, other factors did seem to influence the anxiety and depression of metal music fans, including whether they wrote metal music lyrics, consumed lots of alcohol during concerts and engaged in scarification. With respect to writing songs, it would be interesting to determine whether writing metal music songs increases anxiety or if being anxious motivates one to write metal music songs. After all, history is filled with writers who spent time in psychiatric hospitals, prompting Kaufmann (1990:194) to suggest a ''circular'' relationship between writing and sanity. The connection between alcohol use at concerts and anxiety is fraught with similar ambiguity: it may be that people who are anxious will drink more, regardless of the type of music they prefer. However, previous studies have found links between the use of alcohol and depression [START_REF] Deykin | Adolescent Depression, Alcohol and Drug Abuse[END_REF], and our results did not confirm this connection with our sample of metal music fans.

Finally, the results show a link between depression and scarification, the most extreme form of body modification. Because the lesions are self-inflicted, scarification is usually considered to be an act of self-mutilation. Some studies have found that, of all the middle and high school students needing to see the school nurse, 11.30 percent of girls and 6.60 percent of boys admit to having voluntarily injured themselves (e.g., cuts and burns) in the past 12 months [START_REF] Choquet | Les e ´le `ves a `l'infirmerie scolaire: identification et orientation des jeunes a `haut risque suicidaire. Enque ˆte re ´alise ´e aupre `s de 21 e ´tablissements scolaires du de ´partement de la Gironde[END_REF]). Among the students that admitted to having tried to commit suicide, almost three-quarters of them (72.60 percent) reported cutting or burning themselves, whereas only 15.90 percent of those who have not tried to commit suicide reported similar cuts and burns. Clearly, if scarification is an indicator of suicide risk, it also may be an indicator of anxiety and depression [START_REF] Conwell | Relationships of Age and Axis I Diagnoses in Victims of Completed Suicide: A Psychological Autopsy Study[END_REF]. The results of the present study confirm this finding, as we found a link between scarifications and depression. However, our results also indicated that, although fans of metal music who engage in scarification may have a higher depression score than those who do not engage in this body modification practice, this score remains well below the danger cutoff.

Concerning the dance rituals, the results revealed links between practicing ''pogoing'' or ''slamming'' and mental health. It appears that metal music enthusiasts who show signs of more anxiety participate in the more collective ''slamming,'' whereas metal music enthusiasts who show less anxiety participate in the more individual ''pogoing.'' It may be that the individual nature of ''pogoing'' is enacted by those with a stronger sense of individuality and, thus, with less anxiety about themselves. In contrast, the more collective ''slamming,'' in which participants are passed above the crowd, may be a way for those with more anxiety to learn to trust others and decrease their anxiety. More research could be done to examine these dynamics.

Conclusion

Considering band names such as Suicidal Tendencies, Megadeth and Black Sabbath, and references to death, suicide and self-destruction, it is not surprising that critics of metal music portray it as sick, dangerous, psychopathic and even Satanic, and argue that it has harmful consequences to fans (see [START_REF] Balducci | Adorateurs du diable et rock satanique[END_REF][START_REF] Domergue | Culture barock et gothic flamboyant, un e ´cho surgi des abı ˆmes: la musique extre ˆme[END_REF][START_REF] Regimbald | Le Rock n'Roll, viol de la conscience par les messages subliminaux[END_REF]). This study, however, showed that fans of metal music in France exhibit overall low levels of anxiety and depression. We believe it is quite possible that metal music fans' fascination with images of death may have positive consequences. Many scholars concur that images and open discussions of death have decreased over the past several decades in most Western societies (see [START_REF] Baudry | Le corps extre ˆme. Approche sociologique des conduites a `risque[END_REF][START_REF] Campbell | Puissance du mythe[END_REF]Le Breton 1991;[START_REF] Leroi-Gourhan | Dictionnaire de la pre ´histoire[END_REF]. The mass media in these cultures emphasize youth and fitness rather than aging and death, especially in entertainment and advertising images. Previous generations of people were far more likely to die at home in front of family and close friends, but now the end of life usually occurs behind the closed doors of hospitals and nursing homes, where the sights, sounds and smells of death are camouflaged. Most of us are not taught how to die or how to care for someone who is dying. Many people now do not even talk of death but, instead, use euphemisms such as ''passed away'' and ''lost.''

Obviously, the mass media present many images of death to various audiences. Horror and slasher films offer ghastly and grotesque images of death and destruction, but these images are presented as a form of entertainment and escape rather than reality. In addition, television news and crime shows are often graphic in their presentation of dead and decomposing bodies, some of which may have been tortured or mutilated, though these cases are considered to be deviant, illegal and outside the norm.

In contrast, metal music may be a form of entertainment, but it presents images of death and destruction as more typical occurrences that are not outside the norm. Indeed, the ''death'' and ''black'' forms of metal music, in particular, portray murder and suicide as everyday events that happen to everyday people, even if such portrayals may be thought to be macabre by some. As the band Metallica concluded in their 1984 hit ''Fade to Black,'' ''Death greets me warm, now I will just say goodbye.'' Kubler-Ross and other scholars who have studied death and dying argue that people who confront death more openly and directly are more likely to deal with it in a more healthy manner and to be more comfortable with their own mortality (see [START_REF] Dickinson | Death, Dying and Bereavement[END_REF][START_REF] Kubler-Ross | On Death and Dying[END_REF][START_REF] Berg | Cancer and Death: A Love Story in Two Voices[END_REF]. Perhaps the constant exposure to symbols of death helps fans of metal music to feel less anxiety and depression and, in some ways, to have better mental health than members of the general population, who do not confront death and mortality so openly. Additional research could shed more light on this interesting possibility.

Fig. 1

 1 Fig. 1 Motivations for metal concert attendance/participation

  11, R 2 = 0.13 F(5,293) = 4.04** F(5,293) = 0.75* Note: *p \ 0.05; **p \ 0.01

Table 1

 1 Descriptive results

	a The English term ''piercing''
	in France refers not to
	traditional ear piercings of
	many women and men, but to
	piercings in nontraditional
	areas such as the nose, tongue,
	eyelid, and nipple

Pushing or smashing into each other (when music stops, dancers separate into two opposite camps; when music resumes, dancers in the two camps rush at each other)

Table 3

 3

	HADS factor analysis (after two factors for		Factor 1 (anxiety)	Factor 2 (depression)
	orthogonal rotation were predetermined)	Mean	7.26	3.76
		SD	3.44	2.91
		Variance explained	27.55%	11.17%
		Eigenvalue	3.86	1.56
		a	0.70	0.67

Table 4

 4 Multiple linear regression analysis for rituals in concerts

		Anxiety			Depression
		b	B	t(293)	b	B	t(293)
	Intersect		6.564	9.327		3.482	5.612
	Pogoing	-0.205**	-1.023	-3.036	-0.072	-0.308	-1.035
	Slamming	0.208**	1.341	2.612	0.067	0.372	0.821
	Stage diving	0.067	0.505	0.860	0.067	0.436	0.841
	Moshing	0.034	0.191	0.529	0.002	0.011	0.034
	Bravehearting	-0.001	-0.005	-0.013	-0.029	-0.133	-0.405

Acknowledgments The authors would like to thank Professor Marika Moisseeff (Psychiatry) and Professor Ve ´ronique Campion-Vincent (Folklorism) for their many useful discussions and incisive comments on the first and subsequent drafts of the manuscript. They also would like to thank Fre ´de ´ric Rocheron and Jean Griffet for their support, Thomas McQuaid and Pamela Ruttum for their translation and the two anonymous reviewers for their significant insights and suggestions for improvement.


Augustin Touron

Consistency of the maximum likelihood estimator in seasonal hidden Markov models

In this paper, we introduce a variant of hidden Markov models in which the transition probabilities between the states, as well as the emission distributions, are not constant in time but vary in a periodic manner. This class of models, that we will call seasonal hidden Markov models (SHMM) is particularly useful in practice, as many applications involve a seasonal behaviour. However, up to now, there is no theoretical result regarding this kind of model. We show that under mild assumptions, SHMM are identifiable: we can identify the transition matrices and the emission distributions from the joint distribution of the observations on a period, up to state labelling. We also give sufficient conditions for the strong consistency of the maximum likelihood estimator (MLE). These results are applied to simulated data, using the EM algorithm to compute the MLE. Finally, we show how SHMM can be used in real world applications by applying our model to precipitation data, with mixtures of exponential distributions as emission distributions.

Introduction

Hidden Markov models (HMM) have had various application fields during the last decades: finance [START_REF] Rogemar | Hidden markov models in finance[END_REF], ecology [START_REF] Patterson | Statistical modelling of individual animal movement: an overview of key methods and a discussion of practical challenges[END_REF], climate modelling [START_REF] Wilks | Multisite generalization of a daily stochastic precipitation generation model[END_REF], speech recognition [START_REF] Gales | The application of hidden markov models in speech recognition[END_REF], genomics [START_REF] Yoon | Hidden markov models and their applications in biological sequence analysis[END_REF]) and many more. Let X be a finite set and (Y, Y) a measurable space. A hidden Markov model (HMM) with state space X is a X × Y-valued stochastic process (X t , Y t ) t≥1 where (X t ) t≥1 is a Markov chain and (Y t ) t≥1 are Y-valued random variables that are independent conditonnally on (X t ) t≥1 and such that for all j ≥ 1, the conditionnal distribution of Y j given (X t ) t≥1 only depends on X j . The law of the Markov chain (X t ) t≥1 is determined by its initial distribution π and its transition matrix Q. For all k ∈ X, the distribution of Y 1 given X 1 = k is called the emission distribution in state k. The Markov chain (X t ) t≥1 is called the hidden Markov chain because it is not accessible to observation. The process (Y t ) t≥1 only is observed. See [START_REF] Rabiner | An introduction to hidden Markov models[END_REF] for an introduction to HMM and [START_REF] Cappé | Inference in hidden Markov models[END_REF] for a more general formulation. One very common approach to fit such models is to give a parametric form to the emission distributions and to infer the parameters by maximizing the likelihood function. The asymptotic properties of such an estimator have been widely studied. In [START_REF] Leonard | Statistical inference for probabilistic functions of finite state markov chains[END_REF], the authors proved the consistency and the asymptotic normality of the maximum likelihood estimator (MLE) when the emission distributions have a finite support. Since then, these results have been extended to more general hidden Markov models: see e.g. [START_REF] Douc | Asymptotic properties of the maximum likelihood estimator in autoregressive models with markov regime[END_REF] and references therein.

Motivation In many applications, especially in climate modeling, simple stationary HMM are not adapted because the data exhibit non-stationarities such as trends and seasonal behaviours. Obviously, temperature has a seasonal component, but this is also the case of precipitations or wind speed for example. It is sometimes possible to preprocess the data in order to obtain a stationary residual. However, finding the right form for the non-stationarity can be very tricky, as well as testing the stationarity of the residual. In the case where the non-stationarity is caused by the presence of a seasonality, a popular solution to avoid this pitfall is to split the period into several sub-periods, and to assume stationarity over each sub-period. For example, if there is an annual cycle, one may consider each month separately and fit twelve different models. However, this is not entirely satisfactory, for several reasons.

• A choice has to be made for the length of the time blocks.

• The stationarity assumption over each sub-period may not be satisfied.

• We have to fit independently several sub-models, which requires a lot of data.

• The time series used to fit each of the sub-models is obtained by concatenation of data that do not belong to the same year. For exemple, if the time periods are months, the 31st of January of year n will be followed by the first of January of year n + 1. This is a problem if we use a Markovian model, which exhibits time dependence.

• If the purpose is simulation, it is preferable to be able to simulate a full year using only one model.

Therefore, we prefer to use a extension of HMM that allows seasonality, both in the transition probabilities and in the emission distributions. This will be refered to as a Seasonal Hidden Markov Model (SHMM).

Our contribution In this paper, we will first detail the mathematical framework of a general SHMM. Although this generalization of hidden Markov models is very useful in practice, as far as we know there exists no theoritical result concerning this class of models. The first question that arises is the identifiability of such models. The identifiability of stationary hidden Markov models is not obvious and has been solved only recently.

In [START_REF] Gassiat | Finite state space non parametric hidden markov models are in general identifiable[END_REF], it is proved that the transition matrix and the emission distributions of a hidden Markov model are identifiable from the joint distribution of three consecutive observations, provided that the transition matrix is non-singular and that the emission distributions are linearly independent. [START_REF] Alexandrovich | Nonparametric identification and maximum likelihood estimation for hidden markov models[END_REF] proved that the identifiability can be obtained with the weaker assumption that the emission distributions are distinct. In this paper, we extend the result of [START_REF] Gassiat | Finite state space non parametric hidden markov models are in general identifiable[END_REF] by proving that the SHMM are identifiable up to state labelling, under similar assumptions. To achieve this, we use a spectral method as described in [START_REF] Hsu | A spectral algorithm for learning hidden Markov models[END_REF]. Once we have proved that SHMM are identifiable, it is natural to seek a consistent estimator for their parameters. Regarding HMM, it has been proved by [START_REF] Douc | Consistency of the maximum likelihood estimator for general hidden markov models[END_REF] that under weak assumptions, the maximum likelihood estimator in the framework of HMM is strongly consistent. In this paper, we generalize this result to SHMM. This is done by applying the result of Douc et al. to a well chosen hidden Markov model. We then give several examples of specific models for which our consistency result can be applied. The practical computation of the maximum likelihood estimator for HMM, and a fortiori for SHMM, is not straightforward.

We address this problem by adapting the classical EM algorithm to our framework. Then we run this algorithm on simulated data to illustrate the convergence of the MLE to the true parameters. Finally, we successfully fit a SHMM to precipitation data and show that our model is able to reproduce the statistical behaviour of precipitation.

Outline In Section 2, we give the general formulation of a SHMM in a parametric framework. Using a spectral method [START_REF] Hsu | A spectral algorithm for learning hidden Markov models[END_REF] and a general identifiability result for HMM [START_REF] Gassiat | Finite state space non parametric hidden markov models are in general identifiable[END_REF], we show that under weak assumptions, SHMM are identifiable up to state labelling. Then, we generalize the existing results on the convergence of the MLE in HMM to prove our main result, that is the strong consistency of the MLE in SHMM.

We also give two examples of models (using mixtures of exponential distributions and mixtures of Gaussian distributions) for which this convergence theorem is applicable. In Section 3, we describe the EM algorithm used for the numerical computation of the maximum likelihood estimator and we illustrate the convergence of the MLE for a simple SHMM using simulated data. Finally, in Section 4, we fit a SHMM to precipitation data and we show that such a model can be used to simulate realistic times series of weather variables.

2 Consistency result

Model description

Let K a positive integer and (X t ) t≥1 a non-homogeneous Markov chain with state space X = {1, . . . , K} and initial distribution π, that is π k = P(X 1 = k). For t ≥ 1 and 1 ≤ i, j ≤ K, let

Q ij (t) := P(X t+1 = j | X t = i). Q(t) is the transition matrix from X t to X t+1 . Let us assume that Q(•) is a T -periodic function, so that there exists some integer T ≥ 1 such that for all t ≥ 1, Q(t + T ) = Q(t).
Let Y a Polish space, Y its Borel σ-algebra, and (Y t ) t≥1 a Y-valued stochastic process. Assume that conditionally to (X t ) t≥1 , the (Y t ) t≥1 are independent and that the distribution of Y s conditionally to the (X t ) t≥1 only depends on X s and s. We shall denote by ν k,t the distribution of Y t given X t = k. We also assume that for all t ≥ 1 and for all k ∈ X, ν k,t+T = ν k,t .

Then the process (X t , Y t ) t≥1 is called a seasonal hidden Markov model (SHMM) and the ν k,t 's are its emission distributions.

The law of the process (X t , Y t ) t≥1 is determined by the distribution π of X 1 , the transition matrices Q(1), . . . , Q(T ) and the emission distributions ν k,1 , . . . , ν k,T for 1 ≤ k ≤ K.

Choosing T = 1, we retrieve the classical hidden Markov model (HMM). A general SHMM includes periodicity in both the transition probabilitites and the emission distributions. However, for some applications, it may be enough to consider periodic transitions and constant emission distributions, or vice versa (see e.g. Section 4).

The following remark will be the key to the proof of our main result. Given (X t , Y t ) t≥1 a seasonal hidden Markov model, two (classical) hidden Markov models naturally arise.

• For any t ∈ {1, . . . , T } the process (X jT +t , Y jT +t ) j≥0 ∈ (X × Y)

N is a hidden Markov model with transition matrix Q(t)Q(t + 1) . . . Q(T )Q(1) . . . Q(t -1) and emission distributions (ν k,t ) k∈X .

• For j ≥ 0, let U j := (X jT +1 , . . . , X jT +T ) and

W j := (Y jT +1 , . . . , Y jT +T ). Then (U j , W j ) j≥0 ∈ X T × Y T N is a hidden Markov model. For u = (u 1 , . . . , u T ) ∈ X T , the conditional distribution of W 0 given U 0 = u is T t=1 ν ut,t .
The transition matrix of the homogeneous Markov chain (U j ) j≥0 is given by Quv :

= P(U 1 = v | U 0 = u) = Q u T v1 (T )Q v1v2 (1) . . . Q v T -1 v T (T -1).
Parametric framework Assume that there exists a measure µ defined on Y such that all the emission distributions are absolutely continuous with respect to µ and let f k,t := dν k,t dµ be the emission densities. We consider that all the emission distributions and the transition matrices depend on a parameter θ ∈ Θ, where Θ is a compact subset of some finite-dimensional vector space, e.g. R q . Let us precise the structure of Θ.

• The function t → Q(t) belongs to a (known) parametric family of T -periodic functions indexed by a parameter β that we wish to estimate.

Example

Q ij (t) ∝ exp d l=0 a ijl cos 2πlt T + b ijl sin 2πlt T (1)
In this example, β = (a ijl , b ijl ) 1≤i,j≤K,0≤l≤d and all the transition matrices are entirely determined by β.

• For any t ≥ 1, the emission densities f k,t belong to a (known) parametric family (which does not depend on t) indexed by a parameter θ Y (t). In addition, we assume that the T -periodic function t → θ Y (t) itself belongs to a parametric family indexed by a parameter δ.

Example Denoting by E(α) the exponential distribution with parameter α,

ν k,t = E δ k 1 + cos 2πt T In such a case, θ Y (t) = δ k 1 + cos 2πt
T k∈X and δ = (δ 1 , . . . , δ K ) ∈ R K . • Hence we can define θ = (β, δ).

We shall denote by P π,θ the law of the process (Y t ) t≥1 when the parameter is θ and the distribution of X 1 is π, and by E π,θ (•) the corresponding expected value. If the initial distribution π is the stationary distribution associated with the transition matrix Q(1) • • • Q(T ), then the two HMM described above are stationary. In such a case, we will simply write P θ and E θ (•) for the law of (Y t ) t≥1 and the corresponding expected value. Our purpose is to infer θ from a vector of observations (Y 1 , . . . , Y n ) ∈ Y n . We assume that the model is well specified, which means that there exists a true initial distribution π * and a true parameter θ * = (β * , δ * ) in the interior of Θ such that the observed vector (Y 1 , . . . , Y n ) is generated by the SHMM defined by π * and θ * . We denote by Q * (t) and ν * k,t the corresponding transition matrices and emission distributions. Note that we consider the number K of hidden states to be known.

Identifiability

In [START_REF] Gassiat | Finite state space non parametric hidden markov models are in general identifiable[END_REF], the authors prove that the transition matrix and the emission distributions of a stationary HMM are identifiable (up to state labelling) from the law of three consecutive observations, provided that the transition matrix has full rank and that the emission distributions are linearly independent. Still in the context of HMM, the authors of [START_REF] Alexandrovich | Nonparametric identification and maximum likelihood estimation for hidden markov models[END_REF] use the weaker assumption that the emission distributions are all distinct to obtain identifiability up to state labelling. However, it requires to consider the law of more than three consecutive observations. In this paragraph, we show that under similar assumptions, the SHMM described above is identifiable: we can retrieve the transition matrices and the emission distributions from the law of the process (Y t ) t≥1 up to permutations of the states. We will use the following assumptions: (A1). For 1 ≤ t ≤ T , the transition matrix Q * (t) is invertible and irreducible.

(A2). The matrix Q * (1) • • • Q * (T ) is ergodic and its unique stationary distribution π * is the distribution of X 1 .

(A3). For 1 ≤ t ≤ T , the K emission distributions (ν * k,t ) k∈X are linearly independent.

Remark These assumptions only involve the data generating process and its corresponding parameter θ * , they are not about the whole parameter space Θ.

Theorem 1. Assume that the set of true parameters (π * , Q * (t), ν * k,t ) 1≤t≤T,k∈X satisfies Assumptions (A1)-(A3). Let (π, Q(t), νk,t ) 1≤t≤T,k∈X another set of parameters such that the joint distribution of (Y 1 , . . . , Y T ) is the same under both sets of parameters. Then there exist σ 1 , . . . , σ T permutations of X such that for all k ∈ X, πk = π σ1(k) and for all t ∈ {1, . . . , T }, k, l ∈ X, νk,t = ν σt(k),t and Q(t) kl = Q(t) σt(k),σt+1(l) , with σ T +1 = σ 1 .

Proof. We shall follow the spectral method as presented in [START_REF] Hsu | A spectral algorithm for learning hidden Markov models[END_REF] (see also [START_REF] De Castro | Minimax adaptive estimation of nonparametric hidden Markov models[END_REF][START_REF] De Castro | Consistent estimation of the filtering and marginal smoothing distributions in nonparametric hidden Markov models[END_REF]). Before going through the spectral algorithm, let us first introduce some notations. Let (φ n ) n∈N be a sequence of measurable functions such that for any probability measures ν 1 , ν 2 on (Y, Y),

Y φ n dν 1 n∈N = Y φ n dν 2 n∈N =⇒ ν 1 = ν 2 .
Such a sequence exists because Y is a Polish space. For t ≥ 2 and N ≥ 1, we shall consider the following matrices:

• Let O t = O (N ) t ∈ R N ×K be the matrix defined by (O t ) ak = E [φ a (Y t ) | X t = k] = Y φ a dν * k,t .
• Let L(t) ∈ R N be the vector such that, for 1

≤ a ≤ N , L a (t) = E [φ a (Y t )].
• Let N (t) ∈ R N ×N be the matrix defined by

N ab (t) = E [φ a (Y t )φ b (Y t+1 )].
• Let P (t) ∈ R N ×N be the matrix defined by

P ac (t) = E [φ a (Y t-1 )φ c (Y t+1 )]. • For b ∈ {1, . . . , N }, let M t (•, b, •) be the matrix defined by M t (a, b, c) = E[φ a (Y t-1 )φ b (Y t )φ c (Y t+1 )].
Notice that all these quantities can be computed from the law of (Y t ) 1≤t≤T , except O t which requires the emission distributions. Although they all depend on N , we do not indicate it, for the sake of readability. Using Assumption (A3), we see that there exists an integer N 0 > K such that for all N ≥ N 0 , the matrices O (N ) t have full rank. From now on, we will consider that N ≥ N 0 . We denote by π * (t) the (unconditional) distribution of X t . Elementary calculations show that the following equalities hold:

L(t) = O t π * (t)
(2)

N (t) = O t diag(π * (t))Q * (t)O T t+1
(3)

P (t) = O t-1 diag(π * (t -1))Q * (t -1)Q * (t)O T t+1 (4) M t (•, b, •) = O t-1 diag(π * (t -1))Q * (t -1)diag[O t (b, •)]Q * (t)O T t+1 , 1 ≤ b ≤ N (5)
where diag(v) is the diagonal matrix whose diagonal entries are those of the vector v. Thanks to Assumptions (A1)-(A2), all the entries of π * (t) are positive, so that diag(π * (t)) is invertible. In addition, Assumption (A1) and equations ( 3) and (4) show that the matrices P (t) and N (t) also have rank K.

Let P (t) = U ΣV T be a singular value decomposition (SVD) of P (t): U and V are matrices of size N × K whose columns are orthonormal families being the left (resp. right) singular vectors of P (t) associated with its K non-zero singular values, and Σ = U T P (t)V is an invertible diagonal matrix of size K containing these singular values. Note that such a decomposition is not unique, as we may choose arbitrarily the order of the diagonal entries of Σ, which is equivalent to swapping the columns of U and V , using the same permutation of {1, . . . , K}. Let us define, for 1 ≤ b ≤ N :

B(b) := (U T P (t)V ) -1 U T M t (•, b, •)V.
Hence we have:

B(b) = U T O t-1 diag(π * (t -1))Q * (t -1)Q * (t)O T t+1 V -1 U T O t-1 diag(π * (t -1))Q * (t -1)diag[O t (b, •)]Q * (t)O T t+1 V = Q * (t)O T t+1 V -1 diag[O t (b, •)] Q * (t)O T t+1 V ,
so that there exists an invertible

K × K matrix R := Q * (t)O T t+1 V -1 such that for all b ∈ {1, . . . , N }, diag[O t (b, •)] = R -1 B(b)R.
Besides, as O t has rank K, there exists (α 1 , . . . , α N ) ∈ R N such that the eigenvalues of B := N b=1 α b B(b) are distinct. Hence the eigenvalue decomposition of B is unique up to permutation and scaling. As R diagonalizes B, we obtain R (up to permutation and scaling) by diagonalizing B. Then we can deduce O t (b, •) for all b ∈ {1, . . . , N }, up to a common permutation that we denote by σ t . It follows that for all t ≥ 2, the matrix O t is computable from M t , N (t) and P (t), up to permutation of its columns, corresponding to the states. Then, since O t has full rank, we obtain π * (t) from O t and L(t) thanks to equation (3). Again, π * (t) is only determined up to permutation of its entries: if R is replaced by RP σt where P σt is the matrix of the permutation σ t , we get P T σt π * (t) instead of π * (t). We finally obtain the transition matrix:

Ũ T O t diag(π * (t)) -1 Ũ T N (t)V O T t+1 V -1 = Q * (t),
where Ũ is the matrix whose columns are the left singular vectors of N (t

P T σt Q * (t)P σt+1 instead of Q * (t)
, which means that Q * (t) is only determined up to permutations of its lines and columns, those permutations being possibly different.

Therefore, we have proved that if the law of (Y t ) t≥1 is the same under both sets of parameters (π * , Q * (t), ν * k,t ) 1≤t≤T,k∈X and (π, Q(t), νk,t ) 1≤t≤T,k∈X with the first one satisfying Assumptions (A1)-(A3), then there exists σ

(N ) 1 , . . . , σ (N ) T , σ (N ) T +1 , permutations of {1, . . . , K} such that for all k ∈ X, πk = π * σ (N ) 1 (k)
and for all t ∈ {1, . . . , T },

Q(t) kl = Q * (t) σ (N ) t (k),σ (N ) t+1 (l) , Õ(N) t (•, k) = O (N ) t •, σ (N ) t (k) , k, l ∈ X
where Õt is the analog of O t with respect to the distributions νk,t . The relationship between Q(t) and Q * (t) shows that for all N ≥ N 0 and for all t ∈ {1, . . . , T }, σ

(N ) t = σ (N0) t
. Hence we denote by σ t this permutation. Thus, for all N ≥ N 0 , Õ(N)

t = O (N ) t
P σt . This implies that for all t ∈ {1, . . . , T } and for all k ∈ X, νk,t = ν * σt(k),t and the theorem is proved.

Remarks

• In Theorem 1, we need not assume that the second set of parameters (π, Q(t), νk,t ) 1≤t≤T,k∈X satisfies Assumptions (A1)-(A3) because it has to be the case. Indeed, as the two sets of parameters induce the same distribution of (Y 1 , . . . , Y T ), we have, using the notations of the above proof,

N (t) = O t diag(π(t))Q(t)O T t+1 = Õt diag(π(t)) Q(t) ÕT t+1
Therefore, if the second set of parameters does not satisfy Assumptions (A1)-(A3), there exists some t such that Õt diag(π(t)) Q(t) ÕT t+1 has not full rank. This is a contradiction since

O t diag(π(t))Q(t)O T t+1 has full rank.
• We proved that we could identify the emission distributions and the transition matrices up to permutations, these permutations depending on the time step. However it is not possible to prove that there exist a single permutation of the states that is common to all the times steps. Indeed, if we choose another labelling at time t, the matrix O t is replaced by O t P σ and the matrix Q(t) is replaced by Q(t)P σ , with P σ a permutation matrix. Then we see, using (5), that the matrices M (•, b, •) remain unchanged, which means that the permutation σ cannot be identified from the distribution of (Y t-1 , Y t , Y t+1 ).

• The spectral method presented above provides a non-parametric moment estimator for the parameters of a HMM. The properties of such an estimator are studied in De Castro et al. ( 2017).

• The identifiability of the parameters β and δ depends on the parametric form chosen for t → Q(t) and t → θ Y (t). Hence this should be studied separately for each particular version of the SHMM.

Consistency

In this paragraph, we prove our main result, that is the strong consistency of the maximum likelihood estimator for SHMM. Assume that we have observed (Y 1 , . . . , Y n ) (recall that X 1 , . . . , X n are not observed). For a probability distribution π on X and θ ∈ Θ, let L n,π [θ; (Y 1 , . . . , Y n )] be the likelihood function when the parameter is θ and the distribution of X 1 is π. We define the maximum likelihood estimator by θn,π := arg max

θ∈Θ L n,π [θ; (Y 1 , . . . , Y n )] .
We will need the following assumptions:

(A4). The parameter β can be identified from the transition matrices Q(1), . . . , Q(T ) and the parameter δ can be identified from the emission distributions ν k,t .

(A5).

α := inf θ∈Θ inf 1≤t≤T inf i,j∈X Q ij (t) > 0 (A6).
The transition probabilities (resp. the emission densities) are continuous functions of β (resp. δ).

(A7). For all y ∈ Y, k ∈ X and t ∈ {1, . . . , T },

inf θ∈Θ f θ k,t (y) > 0, sup θ∈Θ f θ k,t (y) < ∞ (A8). For t ∈ {1, . . . , T } and y ∈ Y, we define c t (y) := inf θ∈Θ k∈X f θ k,t (y), d t (y) := sup θ∈Θ k∈X f θ k,t (y),
and we assume that

E ν θ * k,t [-log c t (Y )] < ∞, E ν θ * k,t [log d t (Y )] < ∞ Let S K be the set of permutations of {1, . . . , K}. For σ = (σ 1 , . . . , σ T ) ∈ (S K )
T and θ ∈ Θ, let us denote by σ(θ) the parameter obtained from θ by swapping the states according to the permutations σ 1 , . . . , σ T . More precisely, we compute the transition matrices and the emission distributions corresponding to the parameter θ, we swap them using the permutations σ 1 , . . . , σ T , and using (A4), we identify the parameter corresponding to the swapped matrices and emission distributions. This parameter is denoted by σ(θ). It follows from Theorem 1 that under Assumptions (A1) to (A4),

Θ * := {θ ∈ Θ : P θ = P θ * } ⊂ {σ(θ * ) : σ ∈ (S K ) T }.
Note that due to the parametric form of the transition matrices and the emission distributions, the set {σ ∈ (S K )

T : P θ * = P σ(θ * ) } may actually be much smaller than (S K )

T . However, it contains at least

{σ ∈ (S K ) T : σ 1 = • • • = σ T }.
Theorem 2. Under Assumptions (A1) to (A8), for any initial distribution π and P θ * -a.s., there exists σ

(n) n∈N a (S K ) T -valued sequence such that lim n→∞ σ (n) θn,π = θ * .
Proof. Let us consider (U j , W j ) j≥0 the HMM defined in paragraph 2.1. Recall that its transition matrix is given by

Quv = Q u T v1 (T )Q v1v2 (1) . . . Q v T -1 v T (T -1)
and its emission densities are

g θ (w | u) = T t=1 f θ Y ut,t (w t ).
Thus the law of the process (U j , W j ) j≥0 is entirely determined by θ and π and it is stationary under Assumption (A2). Denoting by Q π,θ the law of the process (W j ) j≥0 when the parameter is θ and the distribution of X 1 is π, we notice that for any θ 1 , θ 2 ∈ Θ,

Q π,θ1 = Q π,θ2 =⇒ P π,θ1 = P π,θ2 .
Therefore, using Theorem 1 and Assumption (A4), we have

Q θ = Q θ * =⇒ ∃ σ ∈ (S K ) T , θ = σ(θ * ). (6)
We notice that for all θ ∈ Θ, J ≥ 0 and initial distribution π, we have :

LJ,π [θ; (W 0 , . . . , W J )] = L (J+1)T,π [θ; (Y 1 , . . . , Y (J+1)T )],
where LJ,π is the likelihood function corresponding to the model (U j , W j ) j≥0 when the distribution of X 1 is π. Let θJ,π be a maximizer of LJ,π . If we are able to prove the strong consistency of θJ,π , by the same arguments, for all s ∈ {0, . . . , T -1}, we can prove that the estimator θs J,π := arg max θ∈Θ L (J+1)T +s,π [θ; Y 1 , . . . , Y (T +1)J+s ] is strongly consistent. From there we easily deduce that θn,π is strongly consistent. Therefore it is sufficient to prove the strong consistency of the maximum likelihood estimator for the HMM (U j , W j ) j≥0 . To this end, we shall use the consistency result stated in Theorem 13.14 in [START_REF] Douc | Nonlinear time series: theory, methods and applications with R examples[END_REF]. The following properties must hold in order to apply this theorem to the HMM (U j , W j ) j≥0 :

1. α := inf θ∈Θ inf u,v∈X T Qθ uv > 0 2. For any u, v ∈ X T and w ∈ Y T , the functions θ → Qθ (u, v) and θ → g θ (w | u) are continuous. 3. For all w ∈ Y T , b -(w) := inf θ∈Θ u∈X T g θ (w | u) > 0, b + (w) := sup θ∈Θ u∈X T g θ (w | u) < ∞ 4. E θ * [| log b + (W 0 )|] < ∞, E θ * [| log b -(W 0 )|] < ∞
The three first properties are straightforward consequences of Assumptions (A5) to (A7). Let us prove that

E θ * [| log b + (W 0 )|] < ∞. The proof that E θ * [| log b -(W 0 )|] < ∞
follows the same lines. We have:

E θ * [| log b + (W 0 )|] = u∈X T πθ * (u) g θ * (w | u)| log b + (w)|µ ⊗T (dw),
where πθ * is the stationary distribution associated with Qθ * . Hence it is enough to prove that for all u ∈ X T ,

g θ * (w | u)| log b + (w)|µ ⊗T (dw) < ∞
We have:

g θ * (w | u)| log b + (w)|µ ⊗T (dw) = b+>1 g θ * (w | u) log b + (w)µ ⊗T (dw) + b+<1 g θ * (w | u)(-log b + (w))µ ⊗T (dw)
Using Assumption (A7), we get inf w∈Y T b + (w) > 0. Therefore, as g θ * is a probability density function, the second term is finite. In order to show that the first one is finite, it is enough to find a function C such that for all w ∈ {b

+ > 1}, C(w) ≥ b + (w) and b+>1 g θ * (w | u) log C(w)µ ⊗T (dw) < ∞. Let C(w) = T t=1 d t (w t ).
For all w ∈ Y T , we get:

C(w) = T t=1 sup θ∈Θ K ut=1 f θ Y ut,t (w t ) ≥ sup θ∈Θ T t=1 K ut=1 f θ Y ut,t (w t ) = sup θ∈Θ u∈X T T t=1 f θ Y ut,t (w t ) = sup θ∈Θ u∈X T g θ (w | u) = b + (w)
In addition:

b+>1 g θ * (w | u) log C(w)µ ⊗T (dw) ≤ Y T g θ * (w | u) T t=1 log d t (w t ) µ ⊗T (dw) = T t=1 Y T T s=1 f θ * Y us,s (w s ) log d t (w t )µ ⊗T (dw) (Fubini) → = T t=1 Y f θ * Y ut,t (w t ) log d t (w t )µ(dw t ) = T t=1 E ν θ * u t ,t [log d t (Y )] (A8) < +∞. Hence E θ * [| log b + (W 0 )|] < ∞.
Thus we can apply Theorem 13.14 in [START_REF] Douc | Nonlinear time series: theory, methods and applications with R examples[END_REF] to get that for any initial distribution π, P θ * -a.s., lim

n→∞ d θn,π , θ ∈ Θ : Q θ = Q θ * = 0,
where d is a distance on Θ. Combining this theorem with (6), we obtain the strong consistency of θJ,π and then Theorem 2 is proved.

Remarks

• The strong consistency of the MLE does not depend on the choice of the initial distribution in the computation of the likelihood function. This is due to the forgetting property of HMM under Assumptions (A5) and (A7).

• Here we just remind the outline of the proof of Theorem 13.14 in [START_REF] Douc | Nonlinear time series: theory, methods and applications with R examples[END_REF]. Note that in the original proof, the authors use a stronger assumption than us. They assume that sup w,u sup θ g θ (w | u) < ∞. However the proof is still valid with our slightly weaker assumption. Let p θ,π (W t | W t-1 , . . . W t-m ) the conditional probability density function of W t given (W t-1 , . . . , W t-m ) when the parameter is θ and the initial distribution is π. The first step of the proof consists in showing that almost surely, for any θ ∈ Θ, lim

m→∞ p θ,π (W t | W t-1 , . . . W t-m ) =: ∆ t,∞ (θ)
is well defined and does not depend on the initial distribution π. This result relies on the geometric forgetting rate of the initial distribution by the HMM. This is where (A5) is crucial. Then we prove that (∆ t,∞ (θ)) t is an ergodic process and that for any θ ∈ Θ and any initial distribution π,

1 J LJ,π [θ; (W 0 , . . . , W J )] -→ J→∞ E θ * [∆ 0,∞ (θ)] =: (θ).
In the last step of the proof, we show that the function is continuous, that for any θ ∈ Θ, (θ) ≤ (θ * ), and

(θ) = (θ * ) =⇒ Q θ = Q θ * .

Applications

In this section, we introduce two examples of SHMM and we show that under weak assumptions, Theorem 2 can be applied. We have seen in section 2.2 that one of the conditions to obtain identifiability in SHMM is the linear independence of the emission distributions. To obtain this, we will use the following lemma.

Lemma 1. 1. Let λ 1 , . . . , λ n be pairwise distinct positive numbers and let us denote by E(λ) the exponential distribution with parameter λ. Then the distributions E(λ 1 ), . . . , E(λ n ) are linearly independent.

2. Let m 1 , . . . , m n be real numbers and σ 2 1 , . . . , σ 2 n be pairwise distinct positive numbers. For k ∈ {1, . . . , n}, let us denote by µ k the Gaussian distribution with mean m k and variance σ 2 k . Then the distributions µ 1 , . . . , µ n are linearly independent.

Proof. Without loss of generality, we can assume that

λ 1 < λ 2 < • • • < λ n . Let (a 1 , . . . , a n ) ∈ R n such that a 1 E(λ 1 ) + • • • + a n E(λ n ) = 0.
For x > 0, we apply this inequality to the Borel set (-∞, x] and we take the derivative with respect to x. We get that for any x > 0,

a 1 λ 1 e -λ1x + • • • + a n λ n e -λnx = 0. ( 7
)
This implies that for any x > 0,

a 1 λ 1 + a 2 λ 2 e -(λ2-λ1)x + • • • + a n λ n e -(λn-λ1)x = 0
Hence, letting x go to infinity, we obtain that a 1 = 0, and equation ( 7) reduces to

a 2 λ 2 e -λ2x + • • • + a n λ n e -λnx = 0.
Thus, step by step, we show that a 1 = • • • = a n = 0, which ends the proof of the first statement. The second one can be proved using the same arguments. Now let us prove a second lemma that we will need to prove the invertibility of transition matrices.

Lemma 2. For n ≥ 1, let us denote by λ n the Lebesgue measure on R n . Let f : R n -→ R be an analytic function and Z(f

) := {x ∈ R n : f (x) = 0} its zero-set. If λ n (Z(f )) > 0 then f ≡ 0.
Proof. We proceed by induction. Let f : R -→ R a real analytic function such that λ 1 (Z(f )) > 0. Then Z(f ) is uncountable, hence it has an accumulation point. As f is analytic, this implies that f ≡ 0. Now let n ≥ 2 and assume that the result holds for analytic functions on R n-1 . Let f : R n -→ R be an analytic function such that λ n (Z(f )) > 0. By Fubini's theorem, we have:

0 < λ n (Z(f )) = R n 1 Z(f ) (x)dx = R dx n R n-1 1 Z(f ) (x 1 , . . . , x n-1 , x n )dx 1 . . . dx n-1 .
This implies that there exists A ⊂ R with λ 1 (A) > 0 such that:

∀x n ∈ A, R n-1 1 Z(f ) (x 1 , . . . , x n-1 , x n )dx 1 . . . dx n-1 > 0, that is ∀x n ∈ A, λ n-1 Z(f ) ∩ (R n-1 × {x n }) > 0.
Thus, for all x n ∈ A, the function f (•, x n ) : (x 1 , . . . , x n-1 ) → f (x 1 , . . . , x n-1 , x n ) is analytic on R n-1 and vanishes on a set with positive Lebesgue measure. Therefore, for all

x n ∈ A, f (•, x n ) ≡ 0. Let y ∈ R n-1 . The function f (y, •) : x n → f (y, x n
) is real analytic and vanishes on A with λ 1 (A) > 0. Hence f (y, •) ≡ 0. As this holds for any y, it shows that f ≡ 0.

Mixtures of exponential distributions

Let (X t , Y t ) t≥1 be a SHMM whose transitions are given by equation ( 1), and whose emission distributions are mixtures of exponential distributions. Denoting by E the exponential distribution, the emission distribution in state k and time t is

ν k,t = M m=1 p km E λ km 1 + σ k (t)
,

where (p k1 , . . . , p kM ) is a vector of probability, the λ km are positive and σ k (t) is a trigonometric polynomial whose constant term is zero and whose degree d is known. This polynomial can be interpreted as a periodic scaling factor. Let δ k the vector of coefficients of σ k . Note that the emission densities are given by

f k,t (y) = M m=1 p km λ km 1 + σ k (t) exp - λ km 1 + σ k (t) y , y > 0. ( 8
)
The vector of parameters of this model is θ = (β, p, Λ, δ) where : 2d+1) is the parameter of the transition probabilities, given by:

• β = (β ijl ) i,j,l ∈ R K×(K-1)×(
Q ij (t) ∝ exp β ij1 + d l=1 β ij,2l cos 2π T lt + β ij,2l+1 sin 2π T lt • p = (p km ) k,m ∈ [0, 1] K×(M -1)
is the set of weights of the mixtures.

• Λ = (λ km ) k,m ∈ (0, +∞) K×M is the set of parameters of the exponential distributions.

• δ = (δ kl ) k,l ∈ R 2d+1 are the coefficients of the trigonometric polynomials (σ k ) k∈X .

We shall make the following assumptions about the true parameter θ * . They ensure that the mixtures have exactly M components.

(A9). For all k ∈ X, λ * k1 < • • • < λ * kM .
(A10). For all k ∈ X and m ∈ {1, . . . , M }, p * km > 0.

Let us first show that the transition matrices and the emission distributions are identifiable. Following the result of paragraph 2.2, it suffices to show that Assumptions (A1)-(A3) are satisfied.

• Clearly, for any time t, the transition matrix Q * (t) is irreducible as all its entries are positive.

• Let us show that for almost every β ∈ R K×(K-1)×(2d+1) , Q * (t) is invertible for all t ∈ {1, . . . , T }. Recall Leibniz's formula for the determinant:

det Q(t) = σ∈S K (σ) K i=1 Q σ(i),i (t),
where S K is the set of permutations of {1, . . . , K} and (σ) is the signature of the permutation σ. Looking at the definition of Q(t), we see that det Q(t) = 0 if and only if

σ∈S K ε(σ) exp K-1 i=1 Z(t) • β σ(i),i = 0,
where • is the standard inner product in R 2d+1 ,

Z(t) = 1, cos 2π T t , sin 2π T t , • • • , cos 2dπ T t , sin 2dπ T t ∈ R 2d+1 ,
and

β σ(i),i = a σ(i),i,0 ; a σ(i),i,1 ; b σ(i),i,1 ; • • • , a σ(i),i,d ; b σ(i),i,d ∈ R 2d+1 .
For all t ∈ {1, . . . , T }, the function

φ t : β → σ∈S K ε(σ) exp K-1 i=1 Z(t) • β σ(i),i
is analytic. Therefore, by Lemma 2, either φ t (β) = 0 for all β, either the Lebesgue measure of its zero-set is zero. Let us define β by βijl = 1 i=j 1 l=1 log K. The corresponding transition matrix for any time t is

1 2K -1      K 1 • • • 1 1 K • • • 1 . . . . . . . . . 1 • • • 1 K     
, which is invertible as it is a diagonally dominant matrix. Hence, for all t ∈ {1, . . . , T }, φ t ( β) = 0. Thus, for any time t, the zero-set of φ t is negligible. As there is a finite number of such functions, we get the desired result and Assumption (A1) is satisfied.

• The following lemma shows that Assumption (A3) is easily satisfied.

Lemma 3. Assume that for all k ∈ X, λ * k1 < λ * k2 < • • • < λ * kM .
For t ∈ {1, . . . , T }, let us define the set

E t = λ * km 1 + σ k (t) : 1 ≤ k ≤ K, 1 ≤ m ≤ M .
If is linearly independent, using Lemma 1. Hence, the family (ν k,t ) k∈X is linearly independent if and only if the rank of the matrix B(t) is K (this requires that p(t) ≥ K). This holds true except if the family (p * km ) belongs to a set of roots of a polynomial. Moreover, the entries of B(t) are among the p * km . This implies that the range of the map t → B(t) is finite. Thus we obtain linear independence of the (ν k,t ) k∈X for all t, except if the p * km belong to a finite union of roots of polynomials.

Then, we can identify the parameter themselves from the emission distributions and the transitions matrices.

• Let us denote by V the variance operator. For t ∈ {1, . . . , T } and k ∈ X, let

s(t) := 1 + σ k (t) 1 + σ k (1) = V(Y t | {X t = k}) V(Y 1 | {X 1 = k})
If the emission distributions are known, s(t) can be computed for any time step t. Let c be the constant coefficient of the trigonometric polynomial s. It follows from the above that c = 1 1+σ k (1) . Hence we get σ k (t) = s(t) c -1. From this we can obtain δ * k .

• Using Assumptions (A9)-(A10) and equation ( 8), we have

lim y→∞ log f * k,1 (y) y = - λ * k1 1 + σ k (1) ,
from which we find λ * k1 . Then we can determine p * k1 using the following formula:

p * k1 = exp lim y→∞ log f * k,1 (y) + λ * k1 y 1 + σ k (1) -log λ * k1 1 + σ k (1) .
Step by step, following the same method, we identify the rest of the parameters of the emission distributions.

• In order to retrieve β from the transition matrices (Q(t)) 1≤t≤T , notice that for any t ∈ {1, . . . , T } and for any i ∈ X,

K-1 j=1 Q ij (t) = K-1 j=1 exp d l=0 a ijl cos 2πlt T + b ijl sin 2πlt T 1 + K-1 j=1 exp d l=0 a ijl cos 2πlt T + b ijl sin 2πlt T , so that K-1 j=1 exp d l=0 a ijl cos 2πlt T + b ijl sin 2πlt T = K-1 j=1 Q ij (t) 1 - K-1 j=1 Q ij (t)
.

Then, for any j ∈ {1, . . . , K -1},

exp d l=0 a ijl cos 2πlt T + b ijl sin 2πlt T = Q ij (t) K-1 j =1 Q ij (t) 1 - K-1 j =1 Q ij (t) .
As a trigonometric polynomial of degree d has at most 2d zeros over a period, this implies that we can retrieve β, as long as T > 2d. Hence Assumption (A4) is satisfied.

In order to prove the strong consistency of the maximum likelihood estimator, it remains to check that Assumptions (A5) to (A8) are satisfied.

• Assuming that there exists β min and β max such that for all θ ∈ Θ and for all i, j, l,

β ijl ∈ [β min , β max ],
Assumption (A5) is satisfied.

• Assumption (A6) is clearly satisfied.

• Assumption (A7) is satisfied provided that:

inf θ∈Θ inf k∈X,m∈{1,...,M } p km > 0 there exists positive numbers λ min and λ max such that for all θ ∈ Θ and for all k, m, λ km ∈ [λ min , λ max ]

there exists positive numbers σ min and σ max such that for all θ ∈ Θ, for all k ∈ X and for all t ∈ {1, . . . , T }, σ k (t) ∈ [σ min , σ max ]. Implicitly, this is a boundedness condition on the parameter δ.

• Under the same boundedness assumptions, we obtain Assumption (A8).

Thus, under weak conditions on the parameters, we can apply our identifiability and convergence results to this particular model.

Mixtures of Gaussian distributions

We choose the same transition matrices as in the previous example, and the emission distribution in state k at time t writes:

ν k,t = M m=1 p km N m k (t), σ 2 km
where m k is a trigonometric polynomial with (known) degree d, (p k1 , . . . , p kM ) is a vector of probability and N (m, σ 2 ) refers to the Gaussian distribution with mean m and variance σ 2 . The following lemma ensures that Assumption (A3) is easily satisfied.

Lemma 4. Assume that for all k ∈ X, σ 2 k1 , . . . , σ 2 kM are pairwise distinct. Let

E = {σ 2 km , 1 ≤ m ≤ M, 1 ≤ k ≤ K}.
If E has at least K elements, then, for all t ∈ {1, . . . , T }, the distributions ν 1,t , . . . , ν K,t are linearly independent for almost every choice of p.

Proof. Using the second statement of Lemma 1, the proof is the same as in Lemma 3.

Hence we can identify the transition matrices and the emission distributions, up to state labelling. Thus, using the fact that finite Gaussian mixtures are identifiable, we can identify, for each state k and each time step t, the vector (p k1 , . . . , p kM ), the mean m k (t) and the variances σ 2 km , up to permutation of the components of the mixture. Finally, for each k, we can identify the coefficients of the trigonometric polynomial m k (•) from its values (m k (1), . . . , m k (T )), so that Assumption (A4) is satisfied. Then, under boundedness conditions on the parameters that are very similar to those of the previous example, we see that Assumptions (A5) to (A8) are satisfied. Hence the strong consistency of the maximum likelihood estimator.

3 Simulation study

Computation of the maximum likelihood estimator

We have already shown the strong consistency of the maximum likelihood estimator. This paragraph deals with its practical computation. Assume that we have observed a trajectory of the process (Y t ) t≥1 with length n. Let X := (X 1 , . . . , X n ) and Y := (Y 1 , ..., Y n ) and recall that X is not observed. The likelihood function with initial distribution π is then

L n,π [θ; Y ] = x π x1 f θ Y x1,1 (Y 1 ) n t=2 Q xt-1xt (t -1)f θ Y xt,t (Y t ),
where x = (x 1 , . . . , x n ). As X is not observed, we use the Expectation Maximization (EM) algorithm to find a local maximum of the log-likelihood function. The EM algorithm is a classical algorithm to perform maximum likelihood inference with incomplete data. See [START_REF] Arthur P Dempster | Maximum likelihood from incomplete data via the EM algorithm[END_REF] for a general formulation of the EM algorithm and [START_REF] Baum | A maximization technique occurring in the statistical analysis of probabilistic functions of Markov chains[END_REF] for its application to HMM. For any initial distribution π, we define the complete log-likelihood by:

log L n,π [θ; (X, Y )] := log π X1 + n-1 t=1 log Q XtXt+1 (t) + n t=1 log f θ Y Xt,t (Y t ) .
This would be the log-likelihood function if X were observed. The algorithm starts from an initial vector of parameters (θ (0) , π (0) ) and alternates between two steps to construct a sequence of parameters θ (q) , π (q) q≥0 .

• The E step is the computation of the intermediate quantity defined by:

Q (θ, π) , θ (q) , π (q) := E π (q) ,θ (q) [log L n,π (θ; (X, Y )) | Y ] .
This requires to compute the smoothing probabilties, that are the a posteriori distributions of X given Y . More precisely, the following quantities are to be computed:

π (q) t|n (k) := P π (q) ,θ (q) (X t = k | Y )
for all k ∈ X and 1 ≤ t ≤ n, and

π (q) t,t+1|n (k, l) := P π (q) ,θ (q) (X t = k, X t+1 = l | Y )
for k, l ∈ X and 1 ≤ t ≤ n-1. The computation of the smoothing probabilities can be done efficiently using the forward-backward algorithm. See [START_REF] Rabiner | An introduction to hidden Markov models[END_REF] or [START_REF] Cappé | Inference in hidden Markov models[END_REF] for a description of this algorithm in the framework of HMM. The adaptation of the forward-backward algorithm for SHMM is straightforward. The intermediate quantity writes:

Q (θ, π) , θ (q) , π (q) = E π (q) ,θ (q) [log L n,π (θ; (X, Y )) | Y ] = K k=1 π (q) 1|n (k) log π k + n-1 t=1 K k=1 K l=1 π (q) t,t+1|n (k, l) log Q kl (t) + n t=1 K k=1 π (q) t|n (k) log f θ Y k,t (Y t ).
• The M step consists in finding θ (q+1) , π (q+1) maximizing the function (θ, π) → Q (θ, π) , θ (q) , π (q) , or at least increasing it. Depending on the specific models chosen for t → Q(t) and t → θ Y (t) it is sometimes possible to find an analytic formula for the solution of this maximization problem. However, in most cases, a numerical optimization algorithm is required.

It can be shown that the sequence of likelihoods L n θ (q) ; Y is increasing and that under regularity conditions, it converges to a local maximum of the likelihood function [START_REF] Cf | On the convergence properties of the EM algorithm[END_REF]. We alternate the two steps of the EM algorithm until we reach a stopping criterion. For example, we can stop the algorithm when the relative difference

L n,π (q) (θ (q+1) ;Y )-L n,π (q) (θ (q) ;Y ) L n,π (q) (θ (q) ;Y )
drops below some threshold ε. The last computed term of the sequence θ (q) q≥0 is then an approximation of the maximum likelihood estimator. However, if the EM algorithm does converge, it only guarantees that the limit is a local maximum of the likelihood function, which may not be global. Therefore it is a common practice to run the algorithm a large number of times, starting from different (e.g. randomly chosen) initial points and select the parameter with the largest likelihood. In [START_REF] Biernacki | Choosing starting values for the em algorithm for getting the highest likelihood in multivariate gaussian mixture models[END_REF], the authors compare several procedures to initialize the EM algorithm, using variants such as SEM [START_REF] Broniatowski | Reconnaissance de mélanges de densités par un algorithme d'apprentissage probabiliste[END_REF]. Introducing randomness in the EM algorithm provides a way to escape from local maxima.

Example

Model Let us consider the following SHMM (X t , Y t ) t≥1 :

• Two states: K = 2 • T = 365
• Q(t), the transition matrix at time t is determined by its first column, given by:

Q i1 (t) = exp β i1 + β i2 cos 2πt T + β i3 sin 2πt T 1 + exp β i1 + β i2 cos 2πt T + β i3 sin 2πt
T Thus the law of the Markov chain (X t ) t≥1 is determined by π (the distribution of X 1 ) and the β il for 1 ≤ i ≤ K and 1 ≤ l ≤ 3.

• The emission distributions are Gaussian:

Y t | {X t = k} ∼ N m k (t), σ 2 k
where the mean m k (t) is given by

m k (t) = µ k + δ k1 cos 2πt T + δ k2 sin 2πt T
Thus the parameters of the emission distributions are the µ k , δ k1 , δ k2 and σ 2 k , for 1 ≤ k ≤ K. Note that this is a special case of the model introduced in paragraph 2.4.2, with M = d = 1. θ is the vector containing both the parameters of the transitions and the parameters of the emission distributions. Then for any choice of θ and π, it is easy to simulate a realization of (X t , Y t ) 1≤t≤nmax . First we simulate the Markov chain (X t ), then we simulate (Y t ) conditionally to (X t ). We chose n max = 200000.

True parameters

• The transition probabilities are given by π * = 0.5 0.5 , β *

1 := (β * 1l ) 1≤l≤3 = 1 0.7 0.5 , β * 2 := (β * 2l ) 1≤l≤3 = -1 -0.6 0.7
The graphs of the functions t → Q * 11 (t) (black) and t → Q * 22 (t) (red) for 1 ≤ t ≤ 365 are depicted in Figure 1. Estimation In order to compute the MLE in this model, we use the EM algorithm described in paragraph 3.1. The estimation procedure is described below:

1. Select a random initial parameter and run the EM algorithm using Y 1 , . . . , Y 500 , with 50 iterations starting from this initial point.

2. Repeat the first step 30 times.

3. Among the 30 initial points candidates, select the one that led to the largest log-likelihood after the EM.

4. Run a long EM using the selected initial point and Y 1 , . . . , Y n . Stop when the relative difference in log-likelihoods drops below 10 -7 . Then the result of the last M-step is θn .

Hence we can compute a sequence of MLE θnp p≥1

. We chose n p = 100p, for 1 ≤ p ≤ 150.

Results

The following graphs show the estimated parameters for the emission distributions, i.e. the µ k , δ k and σ 2 k . The dashed lines represent the true parameters. We can also check that the periodic means m k (t) have been well estimated. In Figure 6 are drawn the graphs of the true means (solid line) and their estimated counterparts (dashed line), for each state. The coefficients β kl for the transition matrices, as well as the transition matrices themselves, are well estimated, as shown in Figures 7 and8.

Remarks

• We had to swap the two states before producing these graphs. This illustrate the fact that the emission distributions and transitions matrices are identifiable only up to label swapping.

• It is not necessary to have such a long times series (200000 observations) to produce good estimates. The figures above show that we can obtain good results with only 20000 observations.

Application to precipitation data

This work was motivated by the design of a stochastic weather generator. A stochastic weather generator [START_REF] Katz | Use of conditional stochastic models to generate climate change scenarios[END_REF]) is a statistical model used whenever we need to quickly produce synthetic time series of weather variables. These series can then be used as input for physical models (e.g. electricity consumption models, hydrological models...), to study climate change, to investigate on extreme values [START_REF] Yaoming | Stochastic modeling of daily precipitation in china[END_REF]... A good weather generator produces times series that can be considered realistic. By realistic we mean that they mimic the behaviour of the variables they are supposed to simulate, according to various criteria. For example, a temperature generator may need to reproduce daily mean temperatures, the seasonality of the variability of the temperature, its global distribution, the distribution of the extreme values, its temporal dependence structure... and so on. The criteria that we wish to consider largely depend on applications. In this section, we introduce a univariate stochastic weather generator that focuses on rainfall.

Data We use data from the European Climate Assessment and Dataset (ECA&D project: http://www.ecad. eu). It consists in daily rainfalls measurements (in millimeters) at the weather station of Bremen, Germany, from 1/1/1950 to 12/31/2015. We remove the 16 February 29 so that every year of the period of observation has 365 days. Thus there are 24090 data points left. Missing data are replaced by drawing at random a value among those corresponding to the same day of the year. The distribution of daily precipitation amounts naturally appears as a mixture of a mass at 0 corresponding to dry days, and a continuous distribution with support on R + corresponding to the intensity of precipitations on rainy days. Also, the data exhibits a seasonal behaviour with an annual cycle. Rainfalls tend to be less frequent and heavier in summer than in winter. Hence, a simple HMM cannot be used to model this data, but using a SHMM seems appropriate.

Model We use a model that is very similar to the one introduced in paragraph 2.4.1. To account for dry days, for each state, we replace the first component of the mixture of exponential distributions by a Dirac mass at 0, so that the emission distribution in state k is

ν k = p k1 δ 0 + M m=2 p km E (λ km ) .
Notice that the emission densities do not depend on t, as introducing periodic emission distributions is not necessary to generate realalistic time series of precipitations. It is enough to consider periodic transitions. Here the dominating measure is µ = δ 0 + λ where λ is the Lebesgue measure over (0, +∞). Hence the emission densities are given by

f k (y) = p k1 1 y=0 + M m=2
p km λ km e -λ km y 1 y>0 .

Recall that (p k1 , . . . , p kM ) is a vector of probability and the transition probabilities between the hidden states are given by equation (1).

Results

We estimate the parameters of the model by maximum likelihood inference, using the EM algorithm described in paragraph 3.1. We choose K = 4 (four states), M = 3 (mixtures of two exponential distributions and a Dirac mass at 0) and d = 2 (complexity of the seasonal components) as these parameters give good validation results. The physical interpretation of the four states is straightforward. State 1 is mostly a dry state: in this state, it rarely rains and when it does rain, the rainfalls amounts are small. On the opposite, state 4 is a rainy state, with heavy rainfalls. Between these two extremes, state 2 and state 3 are intermediate states with moderate precipitation amounts. However, they differ by their precipitations frequency, as state 2 is dry most of the time whereas state 3 is almost always rainy. Figure 10 shows the transition probabilites between the four states as functions of time. Transitions from state 4 ... In particular, we see that in summer, dry states 1 and 2 are less frequent whereas state 4 is the most visited state. It means that in summer, we observe either dry days or heavy rain. It is the opposite in winter, where rainfalls are more frequent but also lighter compared to the summer. These interpretations of the states are consistent with climatology.

Validation As we wish to generate realistic simulations of daily rainfall amounts, i.e. simulations whose statistical properties mimic those of the real data, we evaluate the model by comparing the simulations produced by the model using the estimated parameters to the observed time series. To be specific, 1000 independent simulations are produced, each of them having the same length as the observed series. To perform a simulation, we first simulate a Markov chain X sim t t≥1 with transition matrices Q(t). Then we simulate the observation process Y sim t t≥1 using the estimated densities f θY X sim t ,t . Several criteria can be considered to carry out the comparison: daily statistics (moments, quantiles, maxima, rainfall occurrence), overall distribution of precipitations, distribution of annual maximum, interannual variability, distribution of the length of dry and wet spells... The choice of the criteria mostly depends on the specific application of the model. Each of these statistics is computed from the simulations, which provides an approximation of the distribution of the quantity of interest under the law of the generator (in other words, we use parametric bootstrap), hence a 95% prediction interval. Then this distribution is compared to the value of the same statistic computed using the data. Let us first compare the overall distributions of the real precipitation amounts and the simulated ones by looking at the quantile-quantile plot (see Figure 12). The match is correct, except in the upper tail of the distribution. The last point corresponds to the maximum of the simulated values, which is much larger than the maximum observed value. This should not be considered as a problem: a good weather generator should be able to (sometimes) generate values that are larger than those observed. We then focus on daily distributions. Figure 13 shows the results obtained for the first four daily moments and for the daily frequency of rainfall. It shows that these statistics are well reproduced by the model. Even though we did not introduce seasonal coefficients in emission densities, seasonalities appear both in the frequency of rainfall and the amounts. This is only due to the seasonality of the transition probabilities between the states. The distribution of the duration of dry and wet spells is another quantity of interest when one studies precipitation. A wet (resp. dry) spell is a set of consecutive rainy (resp. dry) days. This statistic provides a way to measure the time dependence of the occurrence process. The results are presented in Figure 14. The dry spells are well modelled, whereas there is a slight underestimation of the frequency of 2-day wet spells while the single day events frequency is slightly overestimated.

Conclusion

We introduced a variant of hidden Markov models called SHMM, adapted to data with a seasonal behaviour. In these models, the transition probabilities between the states are periodic, as well as the emission distributions. We gave sufficient conditions of identifiabiity for SHMM and we proved that under reasonable assumptions on the parameter space, the maximum likelihood estimator is strongly consistent, thus generalizing previous results on HMM. Two specific models for which those conditions are satisfied were given as examples. In order to compute the maximum likelihood estimator, we described the EM algorithm adapted to the framework of SHMM and we used it with simulated data to illustrates our consistency result. In the last section, we applied a SHMM with zero-inflated mixtures of exponential distributions as emission laws to precipitation data. We showed that such a model provides a good example of a stochastic weather generator, as the statistical properties of time series generated by the model are close to those of the observations. In this paper, we considered the number of states of the hidden process as a known parameter. However, in most real world applications, it is unknown. This model selection problem has yet to be adressed. In many applications, the data exhibit trends (e.g. climate change) in addition to seasonalities. However, the techniques presented in this paper cannot be directly adapted to deal with trends, so that this case requires further investigation.

Figure 1 :Figure 2

 12 Figure 1: Transition probabilities Q * 11 (t) (black) and Q * 22 (t) (red)

Figure 2 :

 2 Figure 2: Realization of Y 1 , . . . , Y 1000

Figure 3 :

 3 Figure 3: Estimators μ1,n (left) and μ2,n (right)

Figure 6 :

 6 Figure 6: True means m k (t) and their estimators mk (t) (here computed with n = n max )

Figure 8 :

 8 Figure 7: Estimators βkl

Figure 9 :

 9 Figure 9: Estimated emission densities

Figure 10 :Figure 11 :

 1011 Figure 10: Estimated transition probabilities

  Figure 12: Quantile-quantile plot

Figure 13 :

 13 Figure13: Daily moments and frequency of precipitations. The black line relates to observations, the red solid line is the mean over all simulations, and the dashed lines depict an estimated 95% prediction interval under the model.

Figure 14 :

 14 Figure 14: Distribution of the lengths of dry (left plot) and wet (right plot) spells: observed (bars) versus simulated (error bars). The dots represent the means of the simulations

  ). Replacing O t by O t P σt , π * (t) by P T σt π * (t) and O t+1 by O t+1 P σt+1 in the last equation, we obtain

  For all k ∈ {1, . . . K}, ν k,t is a linear combination of exponential distributions whose parameters belong to E t . Hence,

	p(t)
	ν k,t =	B kj (t)E( λ * j ),
	j=1
	where the B kj (t) are among the p * km . The λ * j being pairwise distinct, the family E( λ * j )	1≤j≤p(t)

E t has cardinality at least K, Assumption (A3) is generically satisfied. Proof. For t ∈ {1, . . . , T }, let us denote by p(t) the cardinality of E t . Thus K ≤ p(t) ≤ KM and we can write E t = { λ * 1 , . . . , λ * p(t) }.

Acknowledgements

The author would like to thank Yohann De Castro, Élisabeth Gassiat, Sylvain Le Corff and Luc Lehéricy from Université Paris-Sud for fruitful discussions and valuable suggestions. This work is supported by EDF. We are grateful to Thi-Thu-Huong Hoang and Sylvie Parey from EDF R&D for providing this subject and for their useful advice.