Details of data source and the study population were described in our previous study (9). In brief, we used nationally representative data from the Korean national health insurance service (NHIS) (14). The NHIS is the single insurer managed by the Korean government and covers all Koreans. All people who enrolled in the NHIS are recommended to receive standardized medical examinations every 2 years. Data on sociodemographic characteristics, a medical treatment history, a health examination results and death information as de-identified and anonymous format are available in the NHIS database.

Inclusion criteria for this study was those who had undergone the health examinations provided by the South Korean NHIS between January 1, 2013, and December 31, 2014 (index year). Of 17,664,058 people with age ≥ 20 years with health examination data in the index year, 8,522,376 underwent three or more health examinations during the past 5 years (including 2013 and 2014). Exclusion criteria was those who had a history of BD (International Classification of Disease, 10th Revision [ICD-10] codes: F30 and F31) before the index year (n = 53,857). Individuals with missing data for one or more variables were also excluded in the analyses. A total of 8,326,953 people were finally included in the analyses (Figure 1). This study was approved by the Institutional Review Board (IRB) of Samsung Medical Center (IRB No. 2020–10-005). We did not have an informed consent because all data was de-identified before they became publicly available.

The detailed general medical examination process and definitions of each medical condition were described elsewhere (15). Each hospital where the health examinations were took place had certification by the NHIS for the quality control. In brief, it includes medical history taking, vital sign check-up and blood sampling. A trained clinician measured blood pressure after at least 5 min’ rest with the arm at an appropriate position. Blood samples were collected after an overnight fast. Measurements included in the blood analysis were total cholesterol, triglyceride, high-density lipoprotein (HDL) and low-density lipoprotein (LDL) cholesterol, hemoglobin, serum creatinine, aspartate aminotransferase (AST), alanine transaminase (ALT), gamma glutamyl transferase (GTP) (14).

The diagnosis of diabetes mellitus (DM) was defined according to the presence of ≥1 claim per year under ICD-10 codes (I11–14) and ≥ 1 claim per year for the prescription of antidiabetic medication or fasting glucose levels ≥126 mg/dL. The diagnosis of hypertension was defined according to the presence of ≥1 claim per year under ICD-10 codes (I10–13 or I15) and ≥ 1 claim per year for the prescription of antihypertensive agents or systolic/diastolic BP of ≥140 mmHg/90 mmHg.

We defined the diagnosis of dyslipidemia as the presence of ≥1 claim per year under ICD-10 code E78 and ≥ 1 claim per year for the prescription of a lipid-lowering agent or a total cholesterol (TC) level ≥ 240 mg/dL. Obesity was defined as a body mass index (BMI) ≥ 25 kg/m² (16) according to the Korean classification. We used survey data in evaluating smoking and alcohol consumption history. Heavy alcohol consumption was defined as ≥30 g/day. Regular exercise was defined as performing >30 min of moderate physical activity ≥5 times per week or > 20 min of strenuous physical activity ≥3 times per week. Low household income level was defined as the income level below 20%.

We used the National Cholesterol Education Program Adult Treatment Panel guideline (NCEP-ATP III) definition (17) to define metabolic syndrome.

We used the same definitions for metabolic variability with those used in our previous study (9). In brief, metabolic variability was defined as intra-individual variability of metabolic variables recorded in health examinations. We used variability independent of the mean as in previous studies (15). Variability independent of the mean was calculated by dividing the standard deviation (SD) by the mean to the power x. Power x was modeled as SD = k × mean^x. It was derived from fitting curves by nonlinear regression analysis. For the analyses, we used the PROC NLIN procedure of the SAS package (18).

Highly variable metabolic parameters were defined as the highest quartile (Q4) of the variability. Low variability was defined as the lower three quartiles (Q1– Q3) of variability. Participants were classified further based on the number of high-variability metabolic parameters (waist circumference, blood pressure, triglyceride, HDL, and fasting blood sugar levels) using a score ranging from 0 to 5, meaning a score of 0 indicated no high-variability parameter, 1 to 5 indicated the number of high-variability parameters out of a total of five parameters.

The primary outcome measure was the new diagnosis of BD during the follow-up period. New onset BD was defined by a primary diagnosis of BD (F30 and F31) based on the International Statistical Classification of Disease and Related Health Problems 10th revision (ICD-10). The study population was followed from baseline to the date of new onset BD, death, or until December 31, 2018, whichever came first.

Because we used the diagnostic code based on the ICD-10 when defining new onset BD, there could be patients previously diagnosed with depression (F32 and F33) who later had a diagnosis of BD (F30 and F31). Considering this as a diagnostic conversion from unipolar depression to BD within mood spectrum disorder, we conducted subgroup analysis by dividing study participants into two groups based on history of depression prior to examination. A previous history of depression was defined as a primary diagnosis of depressive disorder (F32 and F33).

Continuous variables were presented as mean ± standard deviation (SD), while categorical variables were presented as n (%). Participants were classified into six groups based on the number of high-variability metabolic variables. Characteristics of these six groups were compared. The incidence rate of the primary outcome was calculated by dividing the number of incident cases by the total follow-up duration (person-years). The cumulative incidence of new onset BD according to the number of parameters with high variability was presented using unadjusted Kaplan–Meier curves. We conducted the log-rank test to analyze differences between groups. Hazard ratio (HR) and 95% confidence interval (CI) for new onset BD were analyzed using the Cox proportional-hazards model. A multivariable-adjusted proportional-hazards model was applied. Model 1 was not adjusted. Model 2 was adjusted for age and sex. Model 3 was adjusted for age, sex, smoking, alcohol consumption, regular exercise, and income status. Model 4 was adjusted for age, sex, smoking, alcohol consumption, regular exercise, low-income status, baseline diabetes, hypertension, and dyslipidemia. We additionally explored whether the number of metabolic parameters with high variability showed different effects on the development of BD depending on gender and age group. All statistical analyses were performed using SAS version 9.4 (SAS Institute Inc). A p-value of less than 0.05 indicated significance.

Table 1 shows characteristics of participants. Participants with greater numbers of highly variable metabolic parameters (metabolic parameters in the 4th quartile of variability) were older and more likely to be females. Those with greater numbers of highly variable metabolic parameters had diagnoses of hypertension, DM, and dyslipidemia more frequently than those with fewer highly variable metabolic parameters. They also smoked less and drink less. They performed less regular exercise. They had low household income more frequently. Those with greater numbers of highly variable metabolic parameters more frequently had metabolic syndrome (≥ 3 metabolic parameters met the criteria for metabolic syndrome (MS)), while those with fewer highly variable metabolic parameters had MS less frequently (p < 0.001 for all variables displayed in Table 1).

Table 2 displays hazard ratio of new onset BD in the presence of metabolic syndrome. There were a total of 11,027 cases of new onset BD during a mean follow-up of 4.22 years (standard deviation (SD): 0.5 years) in the entire cohort. The presence of metabolic syndrome was associated with the development of BD even after adjusting for aforementioned covariates (HR: 1.16, 95% CI: 1.11–1.21). The presence of the condition that met the criteria for each criterion of metabolic syndrome was also independently associated with the development of depression (HR: 1.05–1.26).

Those with highly variable metabolic parameters (i.e., 4th quartile) had a greater risk of having new onset BD compared to those with less variable metabolic parameters (i.e., 1st quartile) in all models (from model 1 to model 4) (Table 3). All metabolic parameters showed similar results. The number of highly variable metabolic parameters also showed a significant linear association with new onset BD (Figure 1 and Table 4). The association was significant in model 4 (for those with 5 highly variable metabolic parameters: HR = 2.5, 95% CI: 1.9–3.3).

As previously discussed, the experience of depression before the diagnosis of BD can be seen as a diagnostic switch from unipolar to BD. To adjust for the effect of pre-existing depression (F32 and F33), we conducted subgroup analyses depending on the history of pre-existing depression. In both groups, the risk of having diagnosis of BD increased as the number of highly variable metabolic parameters increased. No interaction effect was detected (Figure 2).

Figure 1 shows the effect of sex on the association between newly developed BD and the number of highly variable metabolic parameters. Overall patterns and hazard ratios for metabolic parameters with high variability were similar between males and females (HR = 2.4, 95% CI: 1.6–3.5 for males with five highly variable metabolic parameters; HR = 2.5, 95% CI: 1.7–3.6 for females with five highly variable metabolic parameters).

We also divided participants into two groups based on their age with a cut-off of 50 years old. Overall patterns of hazard ratios for metabolic parameters with high variability were similar between those aged <50 years and those aged ≥50 years, although the hazard ratio was greater in those aged ≥50 years (HR: 1.67, 95% CI: 0.90–3.12 for those aged <50 years; HR: 2.65, 95% CI: 1.97–3.57 for those aged ≥50 years). No interaction effect was detected.