Study design and population

This study uses a nested matched case-control design. Cases are defined as adults (ages 18 and over) living in Ontario who died by suicide between April 1, 2007 and December 31, 2015. Controls are randomly selected age- and sex-matched adults living in Ontario who were alive at the time of the matched control’s death by suicide. Four controls were matched for every suicide case.

Data sources

Data include administrative databases held at ICES (formerly known as the Institute for Clinical Evaluative Sciences) in Toronto, Ontario. ICES is an independent, non-profit research institute whose legal status under Ontario’s health information privacy law allows it to collect and analyze health care and demographic data, without consent, for health system evaluation and improvement. The databases used for this analysis include the Registered Persons Database (RPDB), which includes data on individual’s age, sex and eligibility for public health care insurance and used to capture cases and controls; the Canadian Institute for Health Information National Ambulatory Care Reporting System (CIHI-NACRS), which includes information on emergency department (ED) visits, including prior suicide attempts; the 2016 Postal Code Conversion File (PCCF) for geocoding of hospital postal codes; 2018 road data accessed as open data through Land Information Ontario portal to calculate drive times; Census data, which provide information on neighbourhood-level income expressed as quintiles and the Ontario Marginalization Index (ON-MARG), which includes data on area-level residential instability and dependency; the Immigration, Refugees and Citizenship Canada (IRCC)’s Permanent Resident Database, which provides information on migrant status; and Vital Statistics Death data from the Office of the Registrar General (ORG-D), which provide information on date and cause of death, including suicide-related deaths. The Vital Statistics cause of death field has been found to have over 95% sensitivity when compared to coroner-confirmed suicides, [9] and the algorithm for suicide attempts has also been validated [10]. These datasets were linked using unique encoded identifiers and analyzed at ICES, in accordance with Ontario privacy legislation.

Outcome

The outcome of interest is death by suicide, which is ascertained using International Classification of Diseases (ICD), Ninth and Tenth Revision (ICD-9 and ICD-10), codes for suicide (ICD-9 E950-E959; ICD-10 × 60-X84) or possible suicide (ICD-10 Y10-Y19, Y28) and validated using coroner’s data [9]. This definition is based on the standardization document for ICES at the time of study commencement, taking into account the sensitivity and specificity of these codes in the Ontario population.

Aim 1 exposure: travel time to care

Travel time to general and psychiatric hospitals are the primary exposures. Drive travel time is calculated from the resident’s home to the nearest hospital based on the postal codes of both locations. This analysis uses ArcGIS v. 10.2 by ESRI [11]. No traffic congestion or other travel impedances were assumed. This measure considers road networks, rather than measuring direct point-to-point distance measures from the individual’s postal code to the hospital. There are 5 psychiatric hospitals in Ontario, and 171 general hospitals. Travel time is measured in minutes and then also converted into hours.

Aim 2 exposure: rurality

The exposure is rurality, which is measured using Metropolitan Influence Zones (MIZ) [12]. MIZs are a measure of rurality that considers population size and access to urban centre labour markets. A MIZ where the population is larger than 10,000 people is considered urban [12]. Levels of rurality are determined by those commuting to work in a census metropolitan agglomeration (CMA). A metropolitan influence zone (MIZ) is considered strong when 30% or more individuals work in the CMA, moderate (5–30%), weak (1–5%), and remote ( = < 40 people) [12].

Aim 2 mediation variable: travel time to care

For the second aim, travel time to care is the mediator. It is defined based on road networks between the individual’s residence and the nearest hospital, as described above.

Covariates

The Johns Hopkins Adjusted Clinical Groups (ACGs) ® system version 10 is used to categorize chronic conditions into 32 Aggregate Diagnosis Groups (ADGs) [13]. These 32 disease classifications have been validated for use in predicting mortality among a population-based cohort of adults with schizophrenia in Ontario, Canada [13]. Examples of ADGs include dermatologic conditions, chronic stable and unstable conditions, and time-limited conditions. We also included migrant status, the Ontario Marginalization Indices of instability and dependency, [14] neighbourhood income quintile, and any prior suicide attempt. Prior suicide attempts were defined based on emergency department presentation with self-harm codes X60-X84, Y10-Y19, Y28 [10].

Statistical analyses

We first compared patient characteristics by cases and controls using standardized differences, because we have a large sample size and standardized differences measures the effect size independent of sample size [15]. We then calculated the median travel times and interquartile ranges (IQRs) to general hospitals and psychiatric hospitals among cases and controls of both sexes, and compared them using the Wilcoxon test. We then used multiple logistic regression to model the relationship between travel time to hospitals and death by suicide. In addition, we used multiple logistic regression to analyze the relationship between travel time to psychiatric hospitals and death by suicide. All covariates are included in the model to allow for comparability between effect estimates.

For the mediation analysis, the exposure is rurality, the outcome is death by suicide, and the mediating variable is travel time to general or psychiatric hospital (Fig. 1). To examine travel time as a mediator of the relationship between rurality and suicide, we used multiple logistic regression and the CAUSALMED procedure in SAS to derive the natural direct effect (NDI), the natural indirect effect (NDE), and the percent mediation [16]. This procedure is based on VanderWeele’s mediation approach based on the counterfactual framework. The natural direct effect is the conditional association between the exposure and the outcome, whereas the natural indirect effect is the combination of exposure’s effect on the mediator, and the mediator’s effect on the outcome [17, 18]. We used the case control option in SAS to fit the mediator model by only using observations for those in the control group [16]. We generated a 95% bootstrap Wald confidence interval for the effect estimates, using 1,000 bootstrap samples with a seed value of 740,404 (based on a random number generator).

Fig. 1

Finally, we also examined the possibility of effect modification, by including the interaction term between rurality and travel time in the model.

All analyses were completed using SAS Enterprise version 7.1.

Role of the funding source

The analyses, conclusions, opinions, and statements expressed herein are solely those of the authors and do not reflect those of the funding or data sources.