Suicide mortality among youth in southern Brazil: a spatiotemporal evaluation of socioeconomic vulnerability.

OBJECTIVE: To conduct a geospatial analysis of suicide deaths among young people in the state of Paraná, southern Brazil, and evaluate their association with socioeconomic and spatial determinants.
METHODS: Data were obtained from the Mortality Information System and the Brazilian Institute of Geography and Statistics. Data on suicide mortality rates (SMR) were extracted for three age groups (15-19, 20-24, and 25-29 years) from two 5-year periods (1998-2002 and 2008-2012). Geospatial data were analyzed through exploratory spatial data analysis. We applied Bayesian networks algorithms to explore the network structure of the socioeconomic predictors of SMR.
RESULTS: We observed spatial dependency in SMR in both periods, revealing geospatial clusters of high SMR. Our results show that socioeconomic deprivation at the municipality level was an important determinant of suicide in the youth population in Paraná, and significantly influenced the formation of high-risk SMR clusters.
CONCLUSION: While youth suicide is multifactorial, there are predictable geospatial and sociodemographic factors associated with high SMR among municipalities in Paraná. Suicide among youth aged 15-29 occurs in geographic clusters which are associated with socioeconomic deprivation. Rural settings with poor infrastructure and development also correlate with increased SMR clusters.

Introduction

Global suicide rates have increased 60% over the last 45 years,1 making suicide a public health emergency as a major cause of death worldwide.2 It is estimated that approximately 800,000 people die by suicide each year; in other words, a suicide is completed every 40 seconds.3 In Brazil, 9.7 deaths from suicide per 100,000 population were recorded in 2016, a relatively low rate compared to other countries in the Americas (exceeding only Venezuela, Peru, and Mexico) and worldwide.4 However, suicide is increasing in all age groups, with young people aged 15-29 years at the most risk.3 In this age group, suicide is the second leading cause of death worldwide,3 with an alarmingly high rate of suicide attempts (5-10%).5

Most studies on youth suicide focus on identifying individual risk factors.6,7 However, these individual factors may lose relevance when considered within a geographic and context-specific perspective. Thus, understanding the geospatial patterns and changes over time of suicide mortality clusters can provide context-specific, macro-political evidence to inform public policy and public health. However, no global studies have conducted geospatial analyses of the distribution of suicide fatalities in this age group over time, nor of how these deaths are associated with geopolitical and socioeconomic development.

Recent research has reported geospatial associations between suicide rates and indicators of socioeconomic deprivation, such as unemployment and low household income,8 rural residence,9 material deprivation,10 and educational attainment.11 However, these studies used cross-sectional designs, and did not consider the spatiotemporal evolution of suicide. Moreover, these results are from general age groups, not specifically young people, who may have their own unique, specific predictors. In fact, studies focusing on adolescent and young-adult populations have found that, unlike other age groups, young people have an increased risk of suicide due to depression,12 frequent experiences of bullying,13 belonging to marginalized groups,14 and substance use or mood disorders.15 Given the known spatial association of suicide rates with socioeconomic deprivation, and the fact that youth have unique risk factors for suicide which also are impacted by socioeconomic status, evaluating and understanding the spatiotemporal variance of suicide clusters can support public policies.

Spatial analyses are crucial in the epidemiologic assessment of the impact of social processes and structures in health events,16 especially in low- and middle income countries.17 Understanding the spatial structure and dynamics of the population is a priority for the characterization of health conditions.18 This is especially true where there are socioeconomic disparities that might explain preventable outcomes, such as suicide rates; geographic information systems and spatial analysis allow an evaluation of the context in which the event takes place and describe the surrounding environmental setting.

Thus, the aim of this study was to conduct a spatial analysis of suicide deaths among young people aged 15-29 years in the state of Paraná, Brazil, during two distinct time periods (1998-2002 and 2008-2012), and evaluate their association with socioeconomic determinants.

Methods

Study design and sample site

We conducted an ecological observational study of secondary data obtained from the DATASUS Mortality Information System and the Brazilian Institute of Geography and Statistics (Instituto Brasileiro de Geografia e Estatística [IBGE]), applying geospatial analysis techniques to suicide mortality rate (SMR) data. Data were collected from two 5-year periods (1998-2002 and 2008-2012) for the state of Paraná, southern Brazil (Figure S1, available as online-only supplementary material). According to the latest census (2010), the state of Paraná had a population of 10,444,526 (representing 5.5% of the population of Brazil) distributed across an area of 199,880 km2. It is bordered by the states of São Paulo and Mato Grosso do Sul to the North, Santa Catarina to the South, the Atlantic Ocean to the East, and the countries of Paraguay and Argentina to the West (Figure S1). The state is socially, economically, and racially diverse, and has experienced considerable economic growth. In 2010, it ranked sixth place among Brazilian states in Human Development Index (0.749). By 2012, the state’s gross domestic product (GDP) had increased 170.5% compared to 2002.19

Data sources

Data was obtained from open-access secondary sources at the municipality level (Table S1, available as online-only supplementary material). We obtained and analyzed the annual average number of suicides and the municipal populations in the two five-year periods (1998-2002 and 2008-2012). The two time periods were chosen on the basis of availability of information from the Brazilian Census, which is conducted every 10 years.

Suicide deaths rate

Data on suicide death rates in the 399 municipalities of Paraná were obtained from the Mortality Information System of the Brazilian Ministry of Health.20 Suicide death rates were defined according to the ICD-10 categories of intentional self-harm (X60-X84). Specific mortality rates per 100,000 population were obtained for each of the 399 municipalities. Data was extracted for the three age groups defined a priori (15-19, 20-24, and 25-29 years). We used an empirical Bayes spatial estimator to smooth variance in mortality rates by municipalities.

Socioeconomic indicators

Population sizes for the 15-19, 20-24, and 25-29 age groups were obtained from the IBGE. Socioeconomic data were extracted from the 2000 and 2010 population censuses, available online at IBGE,21 and the Paraná Institute for Economic and Social Development (Instituto Paranaense de Desenvolvimento Econômico e Social [IPARDES]).22 Six socioeconomic and demographic indicators for each city were analyzed: education, income, unemployment, informality, human development index, and GDP.

The data sources for analysis are shown in Table S2, available as online-only supplementary material.

Data analysis

Geospatial data were analyzed by geographic locations (areas) to evaluate the geospatial distribution of SMRs with higher densities of occurrences (clusters). We performed exploratory spatial data analysis in the open-source software environments GeoDa™ version 0.9.5i, to determine measures of global and local spatial autocorrelation,23,24 and QGIS version 2.14.

To evaluate the existence of spatial autocorrelation, we used a Queen-type matrix that allows measurement of nonrandom associations between the value of a variable observed in a given geographical unit and the value of variables observed in neighboring units.24 Using the Global Moran Index,23,24 we calculated global and local spatial autocorrelations evaluating SMR for each municipality in the state of Paraná, stratified by age groups (15-19, 20-24, and 25-29 years). We calculated univariate associations for SMR and bivariate associations with Census-based socioeconomic indicators.21 The Global Moran Index identifies whether the value of the SMR tends to be clustered (positive Moran I) or dispersed (negative Moran I) among geographical areas.23-25

To graphically depict spatial autocorrelation, we applied the local indicators of spatial association (LISA) clustering method. LISA choropleth maps were plotted to identify significant spatial clusters throughout the state of Paraná with high or low association values for the SMR.24,25 Clustered areas are categorized according to the pattern of characteristics in adjacent districts. High/high areas are those districts with high suicide death rates surrounded by other districts with high suicide deaths on univariate analysis. Likewise, low/low districts are those with low rates surrounded by districts which also have low values for the variables of interest. When districts with low rates of suicide deaths are surrounded by districts with high suicidality, LISA maps categorize them as low/high; when the opposite pattern is found, the districts are categorized as high/low.

Bayesian networks

Based on geographic information system results, we modeled the causal structure of those data with geospatial significance. To discover the relationships between variables we applied Bayesian network (BN) structure learning algorithms, which allowed us to learn the structure of networks using the data set. BNs are directed acyclic graphs in which each node (circles in the graph) represents a study variable (socioeconomic indicator or SMR), and each edge (lines connecting nodes) represents the conditional dependent relationship between two nodes.26

The quality of BNs was assessed by a group of investigators specializing in suicide.

To build the BN model, we first created a combination of 32 data sets (Table S2), variously including or excluding the Child Labor (CLB) and Human Development Index (HDI) variables and using either the age variables or their average.

We used the bnpa package in R27 to learn the BN structure of each of the 32 data sets created. The bnpa then used four (“gs,” “iamb,” “fast.iamb,” “inter.iamb”) constrained-based algorithms, combined with four conditional independence tests (“cor,” “zf,” “mi-g,” “mi-g-sh”), to learn the structure of the BN, as well as two (“hc,” “tabu”) score-based algorithms combined with three (“aic-g,” “bic-g,” “bge”) network scores. All of these conditional independence tests and network scores are suitable for working with continuous data, and were selected according to discovery of type of variable for data sets by a function of bnpa. Detailed information about these algorithms, the conditional independence tests, and network scores is available elsewhere.28 As our proposal was to first explore the data, this first step in learning the BN structure was a purely data-driven process. Then, two investigators specializing in suicide analyzed the learned BN structures to incorporate prior knowledge based on the literature and to remove spurious relations, using whitelists and blacklists, both implemented in bnlearn. Briefly, arcs in the whitelist are always included in the network, while arcs in the blacklist are never included in the network. The bnpa package also includes a parametrized process to obtain unbiased estimates of model goodness of fit, using k-fold cross-validation. For this study, we used 100 replications and 10 splits, generating 1,000 BN structures. After the execution of k-fold cross-validation, bnpa also executed the model averaging process,29 whereby all significant arcs from each BN structure are selected. To be considered significant, the arc strength must exceed a suitable threshold and its direction must be greater than 0.5. These values are calculated by bnlearn. All BN arcs that failed this criterion were dropped. Finally, to choose the best networks, the two suicide experts independently analyzed each BN structure. When interpretations differed, these two investigators discussed until a consensus was achieved. If no consensus could be reached, a third investigator analyzed the results. BN analysis data are available online.30

Results

Due to a decrease in overall suicide rates over time, there was an increasing proportion of municipalities with lower suicide rates and a decreasing proportion of municipalities with higher suicide rates (Figure 1). This was seen in all but the 25-to-29-year age groups, where an increase in municipalities with high mortality rates was seen (Figure 1). Among all subregions that showed a high SMR (≥ 15/100,000), the South Central subregion seemed to be most affected in both evaluated periods across all ages (Figure 1). As seen in Figure 1, in the youngest age group (15-20 years), high-SMR clusters shifted from municipalities in the western subregions to the center of the state. The 20-to-24-year age group had a change in geospatial pattern of high suicide rate, moving from municipalities in the center of the state to the southern (Southeast and South Central) and western subregions, with one cluster in the East center subregion. Similarly, the 25-to-29-year age group pattern condensed to high mortality rates in the center and southeast subregions of the state. In all three age groups, we also observed a concentration of high SMR in the western part of the state during the first period of analysis, which subsequently decreased over time.

Figure 1

Spatial distribution of suicide mortality rate (SMR) across municipalities, stratified by age group, with ranges of standard deviation from the average for delimitation of class intervals; Paraná, 1998-2002 and 2008-2012.

SMRs suggest positive spatial autocorrelations in the two 5-year periods in all age groups (p = 0.0010). In the first 5-year period, the age groups 15-19, 20-24, and 25-29 had municipalities with high SMR surrounded by other municipalities with high SMR (Global Moran Index: I = 0.551923, I = 0.457190, I = 0.542864, respectively, per age group), demonstrating that the respective Moran indices showed a positive spatial autocorrelation, that is, for all age groups in the two periods analyzed, municipalities with high suicide rates were surrounded by other municipalities with high rates. For the second 5-year period, the same clustering pattern was observed for each age group (Global Moran Index: I = 0.495064, I = 0.525674, I = 0.386640) (Figure 2).

Figure 2

Local indicators of spatial association, univariate analysis: cluster formation according to suicide mortality rate (SMR) according to age group (types of cluster: high-high; low-low; low-high, high-low).

Bivariate analysis demonstrated that all socioeconomic and demographic indicators used for our analysis were significantly associated with the SMR (p < 0.05). During each 5-year period, SMR correlated negatively with higher socioeconomic status: specifically, educational level (≥ 8 years), income, HDI, and GDP. Unemployment and informal employment had a more complex relationship. Higher unemployment (more common in urban areas) correlated with lower SMR clusters, and informal employment (more common in rural areas), with higher SMR clusters (Table 1). The relationship between informal employment and SMR clusters appeared stronger during the second time period, suggesting that increasing informal jobs (and decreasing formal jobs) led to higher rates of suicide mortality (Table 1).

Table 1

Global Moran Index (I) bivariate coefficient of suicide rates according to municipality of residence and socioeconomic and demographic indicators in the 1998-2002 and 2008-2012 time periods

	1998-2002		2008-2012
Variable	Moran’s I	p-value	Moran’s I	p-value
≥ 8 years of schooling	-0.16	< 0.001	-0.19	0.001
Income	-0.18	0.001	-0.18	0.001
Unemployment	-0.15	0.001	-0.11	0.001
Informality	0.01	0.033	0.18	0.001
Human development index	-0.17	0.001	-0.16	0.001
Gross domestic product	-0.06	0.009	-0.08	0.002

By combining cross-validation results (Figure 3), model averaging, and their expertise, the researchers selected the BN structure learned with the hc algorithm for both study periods. All scores for this algorithm presented the same loss rate.

Figure 3

Loss rate of cross-validation procedure.

The BN of the first studied period showed the positive impact of educational attainment on income and unemployment, i.e., additional years of study result in improvement of incomes and lower unemployment, which leads to a decline in mortality rates. On geospatial analysis, these variables showed a positive association with suicide rates. However, regarding educational level, the spatial analysis has also identified that, during this period, SMR correlated negatively with higher educational attainment (Figure 4).

Figure 4

Relationship between possible suicide predictors as of the year 2000, using hill climbing and aic-g score network test. EDU = education; INF = informality; INC = income; UNP = unemployed; SCD = suicide.

During the second studied period, educational attainment again showed a positive relationship with income and unemployment; higher educational attainment, higher income, and lower unemployment were associated with consequent reductions in mortality rates. Likewise, the geospatial analysis also identified during a negative correlation of suicide rates with higher educational level during this period. Informality had a marked positive association with suicide rates, as well as a negative relationship with income, i.e., the higher the income, the lower the informality (Figure 5).

Figure 5

Relationship between possible suicide predictors as of the year 2010. c = regression coefficient; EDU = education; INF = informality; INC = income; UNP = unemployed; SCD = suicide.

Discussion

This is the first study to conduct a spatiotemporal analysis of suicide mortality among young people and its association with socioeconomic indicators. While overall we found that SMRs are decreasing, we found a significant clustering of suicide mortality in both periods of analysis (1998-2002 and 2008-2012), with positive associations between socioeconomic deprivation and mortality rates.

Comparison of both periods shows a decline in SMR over time; this runs counter to the global trend of increasing suicide rates, specifically among the young population.31 We attribute this phenomenon to significant investments in social programs designed to combat poverty in Brazil.32 This is evident from the increase in HDI and GDP between the first and second time periods.20 Similarly, these social programs and economic development can explain the displacement of high suicide mortality to the less developed central regions of the state of Paraná.

Our findings that suicide mortality clusters in regions of socioeconomic deprivation and limited population development mirrors the international literature.33 As an example, the South Central region of Paraná, in which we identified high-SMR clusters, showed low population density, high levels of poverty, low HDI, high unemployment, and substantial social inequality.23 All of these conditions have been established as risk factors for suicide.8 One possible explanation for the association between low socioeconomic development and suicide is limited availability of essential services, such as healthcare, education, culture, employment and housing, which can lead to hopelessness, persistent stress, and anxiety.34

National and international studies have found associations of high SMR with income and unemployment.33,34 Our finding of higher youth SMR in areas with these phenomena correlates with international data that high poverty is a risk factor for suicide attempts.35 This association between income and suicide may vary according to the socioeconomic development of the geographic area analyzed. A study in Brazil found the highest suicide rates in the poorer areas of São Paulo state, and richer areas of Brazil and São Paulo city.36 However, a spatial analysis study showed that suicide mortality in Brazil has a weak spatial correlation and low or no spatial relationship with socioeconomic factors.37

Informal employment (more common in rural areas) was associated with clusters of higher SMR. The rural area has particular relevance in the context of suicide, given the significant, growing relation between suicide and pesticides. An increasing number of suicides have been observed due to ingestion of these chemicals, especially in Asian countries, but also in Central and South America, including Brazil.38 The overall number of suicides is probably higher than estimated; pesticide ingestion is the most common method of suicide worldwide,38 and should be the subject of further studies.

Our findings mirror previous research demonstrating an association between suicide and low educational level.39 Individuals with higher educational attainment are more receptive to health promotion and prevention interventions, and likely have higher rates of employment and, therefore, income or earning potential.40 In our findings, the South Central region of the state of Paraná had high uneducated proportions of the population and low rates of formal employment, highlighting the limited socioeconomic development of the region and explaining higher-SMR clusters. Similarly, our data supports findings suggesting that urbanization is protective against suicide, as no high SMR clusters occurred in large metropolitan areas.41 In Brazil, urban regions also have higher per capita income, less economic inequality, and higher educational levels, providing more opportunities.34 For young people, urban areas provide infrastructure and possibilities for educational, professional, and cultural development, as well as leisure opportunities, all of which can be protective factors against suicide.31

While suicide in young people is multifactorial, there are predictable geospatial and sociodemographic factors associated with high rates of suicide mortality among municipalities in Brazil.42 Targeted interventions in these geospatial regions or for those with known risk factors can further reduce the rates of suicide among younger age groups. We hope these findings can serve as inputs to increase social programs and health infrastructure support to those regions most in need of reducing youth suicide.

Some limitations of our study must be borne in mind when interpreting our results. First, we performed an analysis of secondary data, which is thus subject to limitations of the data source. For instance, underreporting of suicide is common, as is incomplete data collection. These limitations are likely more pronounced in rural settings, which could lead to underrepresentation of the magnitude of our findings. Again, performing a secondary data analysis restricted the depth of the socioeconomic variables available; further data on violence, access to and quality of health care, mental health, and social inequality could enhance our understanding of associations. Similarly, many of our variables are collinear (education, income, employment), which might make relationships more opaque.

While youth suicide is multifactorial, there are predictable geospatial and sociodemographic factors associated with high SMR among municipalities in Paraná, Brazil. Suicide among youth aged 15-29 occurs in geographic clusters which are associated with socioeconomic deprivation. Rural settings with poor infrastructure and development also correlate with increased SMR clusters.

Disclosure

The authors report no conflicts of interest.