K-means Clustering Reveals Long-term Thyrotropin receptor antibody Patterns in Graves’ Disease: Insights from a 10-Year Study with Implications for Graves’ Orbitopathy

1. Introduction

Graves’ disease (GD) is an autoimmune disorder primarily characterized by excessive thyroid hormone secretion. Thyroid-stimulating hormone receptor antibodies (TRAb) play a key etiological role in this process by stimulating thyroid-stimulating hormone (TSH) receptors, thereby leading to hyperthyroidism [1]. The TRAb assay is a vital tool for diagnosing GD owing to its high sensitivity and specificity. [2] In addition, TRAb levels typically decline as GD is managed with antithyroid medications or thyroidectomy, although radioactive iodine (RAI) treatment can transiently elevate TRAb levels during the first year [3,4,5]. During the clinical course of GD treatment, the specific cutoff values of TRAb serve as good predictors of relapse [6]. By evaluating whether TRAb remains elevated or decreases to normal ranges, clinicians can make more informed decisions regarding the continuation or cessation of antithyroid therapy [7,8].

TRAb also plays an important role in the pathogenesis of Graves’ orbitopathy (GO) by increasing hyaluronic acid synthesis and promoting adipogenesis in orbital fibroblasts [9]. Patients with a severe course of GO exhibit significantly higher TRAb levels at diagnosis and during follow-up [10]. High TRAb is also considered a risk factor for GO progression following RAI treatment and guides clinicians in determining whether steroid prophylaxis should be prescribed to prevent GO exacerbation after RAI treatment [11].

As TRAb is an important factor in both the pathogenesis and treatment guidelines for GD and GO, investigating its long-term changes in GD is worthwhile. While several studies have explored the overall average changes in TRAb levels during different GD treatments, there is a notable research gap in classifying the distinct patterns of TRAb changes over time during GD treatment and their relationship with GO comorbidity [3,5,12]. Identifying distinct patterns of TRAb changes could provide deeper insights into disease prognosis and help optimize individualized treatment strategies.

Predicting TRAb trends can also aid in clinical decision-making by guiding the duration of antithyroid medication and monitoring for extrathyroidal manifestations, including GO [13]. Therefore, the primary objective of this study was to characterize the long-term patterns of TRAb changes during GD treatment and to evaluate their association with the presence of GO. We also aimed to identify key factors related to TRAb normalization. By investigating these aspects, we hope to enhance our understanding of TRAb dynamics and their implications for the management and potential progression of GO.

2. Materials and Methods

2.1. Patient Recruitment

A retrospective analysis was conducted on the medical records of adult patients (aged 19 years and above) diagnosed with GD at Tertiary National University Hospital between the years 2010 and 2021. The inclusion criteria for the study required patients to have a TRAb level exceeding 1.5 IU/L and a minimum of 3-year follow-up on TRAb levels. The exclusion criteria included the absence of laboratory test results, prior diagnosis of GO, or missing treatment records. The Institutional Review Board (IRB) of Pusan National University Hospital approved this study, as indicated by the assigned IRB number 2306-030-128.

2.2. Clinical and Laboratory Assessment

TRAb levels were measured using the TRAK human radio receptor assay (Thermo Fisher Scientific, Germany), which has a reference range of 0–1.5 IU/L. Yearly TRAb values were collected for each patient, and missing measurements were estimated by linear interpolation [14]. Linear interpolation is a commonly used method for time-series data in which a straight line is drawn between two known data points, and the missing value is assumed to lie along this line.

We applied k-means clustering to identify distinct longitudinal patterns of TRAb changes (“Patterns A, B, C, D”) based on baseline TRAb levels, changes over time, and whether normalization occurred. In addition, we performed a cox-regression analysis using a median-based threshold of 6.14 IU/L to divide patients into “lower-TRAb” (<6.14 IU/L) and “higher-TRAb” (≥6.14 IU/L) groups. Although no formal cutoff analysis (e.g., ROC curve) was performed, 6.14 IU/L was chosen because it represented the median TRAb value in our study population.

GD treatment status during the follow-up period was collected, including the duration of methimazole prescription, statin use, RAI treatment history, and thyroidectomy. GD duration and GO comorbidities were also reviewed. GO diagnosis was confirmed by a specialized ophthalmologist at the same hospital whose primary responsibility was to manage patients with GO. The clinical parameters were compared between the different TRAb patterns and investigated to determine whether they were related to TRAb normalization.

2.3. Statistical Analysis

Clustering analysis: To identify homogeneous patient trajectories, we used the KmL algorithm (Christophe and Bruno, 2009) for functional cluster analysis. KmL, a non-parametric extension of k-means for longitudinal data, clusters based on the overall shape and temporal changes of trajectories, rather than individual points. For handling missing values, we utilized the 'linearInterpol.local' imputation method, which involves linear interpolation for intermittent missing values and considers the global slope at the end [14,15]. The number of clusters was determined based on clinical judgment and interpretation. Clinical parameters between the different patterns were compared using Pearson's chi-square test, Kruskal–Wallis rank sum test.

Survival analysis: We estimated the median time to normalization of TRAb values using Kaplan-Meier curves. Survival data from the original dataset were utilized, with "survival time" defined as the duration from the index date to either normalization time or the last follow-up observation. The "status" variable was coded as 1 for normalization and 0 for cases without normalization (censored). Univariate and multivariate Cox regression analyses were conducted to identify factors affecting the time to TRAb value normalization. For the Cox models, we chose to present results as “Survival ratio” rather than the conventional hazard ratio, aiming to emphasize the proportion of patients who achieve TRAb normalization over time. Specifically, the “Survival ratio” reflects the relative rate of normalization events in one group compared to another at a given time point.

In all analyses, statistical analyses were performed using the R statistical language (version 4.0.5; R Core Team, 2021) and additional packages (kml, longitudinalData). A p-value of <0.05 was considered statistically significant.

3. Results

A total of 403 patients (111 men) were included in this study. The mean age at enrollment was 47.65 ± 15.58 years, and the mean follow-up duration was 4.88 ± 2.33 years. At diagnosis, the mean disease duration of GD was 0.49 ± 0.86 years, and 63 patients (16%) presented with comorbid GO. Antithyroid medication (methimazole) was prescribed to 359 patients (89%) for an average of 44 ± 33.1 months (1,320 ± 993 days). Statin therapy was used in 40 patients (9.9%), radioactive iodine (RAI) therapy was performed in 47 (12%), and 5 (1.2%) underwent thyroidectomy during the follow-up (Table 1).

The mean baseline TRAb level was 8.98 ± 8.00 IU/L. Over a 10-year period, the mean of all measured TRAb values (“10-year average”) was 4.76 ± 6.38 IU/L, with a within-patient standard deviation (“10-year SD”) of 4.10 ± 3.74 IU/L. The annual rate of TRAb change was -0.29 ± 1.26 IU/L/year, indicating a general decline in TRAb levels across most patients, although some exhibited stable or even increasing trends. In total, 333 patients (83%) achieved TRAb normalization (<1.5 IU/L), and the mean time to normalization was 2.60 ± 1.87 years.

3.1. Patterns of TRAb Over Time and Overall Clinical Information

Out of 4,433 total TRAb measurements, 567 interpolated data points were generated after excluding measurements obtained following TRAb normalization (with no subsequent follow-up). Four distinct patterns of TRAb were identified over time and were designated as A, B, C, and D based on initial TRAb levels, annual rates of change, and final normalization outcomes (Figure 1, Table 1).

Pattern A (n = 291, 72.2%) had the lowest baseline TRAb (6.8 ± 5.42 IU/L) and showed a mean TRAb change of -0.41 ± 0.41 IU/L/year. This group attained a 96% normalization rate, with a median survival time of 2 years (p < 0.001).

Pattern B (n = 51, 12.6%) had the highest baseline TRAb (20.39 ± 10.21 IU/L) and exhibited a steep negative trend (-1.98 ± 0.77 IU/L/year). TRAb normalization occurred in 80% of these patients, and the median survival time was 8 years (p < 0.001).

By contrast, Patterns C (n = 38, 9.4%) and D (n = 23, 5.7%) showed positive mean TRAb changes (0.88 ± 1.28 IU/L/year and 2.95 ± 1.65 IU/L/year, respectively) and lower normalization rates (29% and 13%). Fewer than half of the patients in Patterns C and D achieved normalization within the study period (p < 0.001). Patterns B and D demonstrated steeper annual TRAb changes than A and C (P < 0.001), and Patterns A and B had higher final normalization rates (96% and 80%) than C and D (29% and 13%) (P < 0.001).

Among these four patterns, Pattern B had the longest mean disease duration (1.02 ± 1.48 years, p = 0.018), whereas Pattern D showed the highest thyroidectomy rate (8.7%, p = 0.006). Comorbidity with GO did not significantly differ across patterns (p = 0.876). MMI was prescribed for the shortest duration in Pattern A (1,118 ± 910 days, p < 0.001), while RAI therapy was more frequently used in Patterns B (35%), C (16%), and D (30%) than in Pattern A (5.5%, p < 0.001). Statin prescription rates did not differ significantly among the four groups (p = 0.759). Detailed demographic and clinical variables for each pattern are provided in Table 1.

3.2. Two Main Patterns Based on Normalization Rate

When Pattern A was compared with the combined Patterns B, C, and D (hereafter “BCD group”), Pattern A showed a significantly higher normalization rate. Specifically, 278 of 291 patients (96%) in Pattern A achieved TRAb normalization, whereas only 55 of 112 patients (49%) in the BCD group did so (p < 0.001). Consistent with these findings, the baseline TRAb level in Pattern A was markedly lower (6.80 ± 5.42 IU/L) than in the BCD group (14.65 ± 10.50 IU/L; p < 0.001; Table 2).

Disease duration was also shorter in Pattern A (0.39 ± 0.69 years) compared with the BCD group (0.73 ± 1.16 years; p = 0.004). Methimazole (MMI) therapy was prescribed for an average of 1,118 ± 910 days in Pattern A, significantly less than the 1,830 ± 1,013 days in the BCD group (p < 0.001). In addition, the BCD group had higher rates of radioactive iodine (RAI) therapy (28% vs. 5.5%, p < 0.001) and thyroidectomy (3.6% vs. 0.3%, p = 0.022). No significant differences in MMI usage frequency (p = 0.427) or statin prescription (p = 0.246) were noted between the two groups.

The mean TRAb change rate was -0.41 ± 0.41 IU/L/year in Pattern A, compared with 0.00 ± 2.28 IU/L/year in the BCD group (p < 0.001). Likewise, the time to TRAb normalization was shorter in Pattern A (2.39 ± 1.64 years) than in the BCD group (4.27 ± 2.64 years; p < 0.001). Detailed demographic and clinical variables for both Pattern A and the BCD group are summarized in Table 2.

3.3. Survival Analysis of TRAb

The median time to TRAb normalization in the overall cohort was 3 years (Figure 2A). When patients were classified into the four predefined patterns (A, B, C, D), Pattern A had a median survival time of 2 years, whereas Pattern B showed a notably longer median survival time of 8 years. In Patterns C and D, fewer than half of the patients achieved normalization within the 10-year observation period.

An additional survival analysis compared Pattern A with the combined Patterns B, C, and D (hereafter “BCD group”). In this comparison, the BCD group demonstrated a median survival time of 9 years (Figure 2B,C), indicating a slower normalization rate than that observed in Pattern A.

3.4. Prognostic Factors Associated with Patterns

We performed survival analyses using a Cox regression model but present the results as “Survival ratio” (SR) rather than the conventional hazard ratio. These analyses demonstrated that Pattern A was significantly associated with a higher probability of TRAb normalization compared with Patterns B, C, and D (p < 0.001; Table 3). In the multivariable model, the SR for TRAb normalization in Pattern A was 8.391 (95% CI: 3.087–22.81), whereas Patterns B and C did not show significant differences relative to Pattern D.

Among the patient characteristics examined, female sex increased the likelihood of TRAb normalization by approximately 1.3-fold (p = 0.044), and having an initial TRAb level below 6.14 IU/L was associated with a 1.7-fold higher probability of normalization (p < 0.001). Age, presence of GO, methimazole (MMI) treatment, statin use, and radioactive iodine (RAI) therapy were not significant predictors of TRAb normalization in the multivariable analysis. Although RAI therapy showed a reduced likelihood of TRAb normalization in the univariate model (SR, 0.415; 95% CI: 0.272–0.632; p < 0.001), this effect did not persist in the multivariable model (SR, 0.715; 95% CI: 0.459–1.113; p = 0.138). Thyroidectomy also did not exhibit a statistically significant association with TRAb normalization in the final analysis.

In a separate comparison of Pattern A (Group A) versus the combined BCD group (Patterns B, C, and D), similar results were obtained: female sex and lower baseline TRAb (< 6.14 IU/L) were significantly linked to faster TRAb normalization (data not shown). Detailed SRs and confidence intervals for each factor are provided in Table 3.

4. Discussion

Our findings demonstrate that the median time to TRAb normalization in this cohort was 3 years, underscoring the heterogeneous clinical courses of GD. Based on baseline TRAb levels, rate of changes, and normalization rates, we identified four distinct patterns. These could be further grouped into two main categories depending on whether the median time to normalization was within 2 years. Pattern A rapidly achieved normalization (median 2 years), whereas the remaining patterns generally required longer durations or exhibited normalization rates below 50%. Although these discrepancies were linked to factors such as disease duration, treatment modalities, and baseline TRAb concentrations, GO comorbidity rates did not differ significantly across the patterns. Notably, female sex and lower initial TRAb (<6.14 IU/L) were associated with greater probabilities of achieving normalization.

Among these patterns, Pattern A—observed in the majority of our cohort—appears to reflect the most favorable disease trajectory.Patients in Pattern A had the lowest baseline TRAb, the shortest time to normalization, less frequent use of RAI, and shorter durations of MMI therapy. This aligns with earlier reports indicating that an early and sustained response to a 12–18-month course of MMI often leads to remission in GD [7,16]. Patterns B, C, and D demonstrated relatively higher initial TRAb levels, longer methimazole prescription durations, and higher RAI and thyroidectomy prescription rates than pattern A. Notably, a high baseline TRAb level is a recognized risk factor for relapse [17,18], and discontinuation of antithyroid medication has been associated with relapse rates of 40–50% in Europe and Japan and approximately 37% in South Korea [19,20]. Consistent with these observations, the 2016 American Thyroid Association guidelines recommend either extending antithyroid therapy or opting for definitive treatments (RAI or surgery) in patients who demonstrate persistently elevated TRAb levels or experience clinical relapse [7].

Patterns B, C, and D likely represent clinical trajectories in which patients failed to achieve remission with an initial MMI course or experienced relapse after medication withdrawal.These patterns showed higher baseline TRAb levels, more frequent second-line or ablative interventions, and extended MMI prescriptions, consistent with earlier reports that prolonged or repeated antithyroid therapy—combined with elevated TRAb—often necessitates RAI or thyroidectomy [21,22]. Notably, Pattern B achieved a higher rate of TRAb normalization (80%) than Patterns C (29%) and D (13%), despite exhibiting the highest baseline TRAb level. Although Pattern B required a longer median time (8 years) to normalization than Pattern A, its eventual success highlights the role of additional or extended therapies for patients who fail to achieve remission after the initial MMI course. Indeed, previous studies have noted that a substantial proportion of GD patients who relapse following first-line MMI can still attain euthyroidism with prolonged medical therapy or RAI [23]. Still, the median survival time of 8 years and a mean TRAb change rate of -1.98 ± 0.77 IU/L/year in Pattern B highlight the impact of baseline TRAb levels on long-term outcomes [24]. Even with comprehensive management—including pharmacological and surgical interventions—achieving TRAb normalization within a widely acceptable timeframe (e.g., 2–3 years) can remain elusive for many patients [24].

Patterns C and D further illustrate the complexity of TRAb dynamics by offering distinct pathways to persistent elevation.Pattern C began with a relatively elevated baseline TRAb but exhibited only a modest positive slope, resulting in a 29% normalization rate. In contrast, Pattern D, despite starting from a moderate baseline TRAb, showed both fluctuations and an overall upward trend, even after multiple treatments, including RAI and, in 8.7% of cases, thyroidectomy. Some individuals in Pattern D maintained persistently elevated TRAb levels following surgery, potentially reflecting residual thyroid tissue or ongoing autoimmune activity [25]. This resistance to conventional treatments underscores the need for more aggressive or alternative therapeutic approaches and suggests possible immunologic or genetic underpinnings [5,26,27,28].

Regarding GO, our data showed no significant difference in GO comorbidity across these TRAb patterns, despite notable variations in RAI use. Specifically, 28% of patients in the combined BCD group underwent RAI compared with 14% in Pattern A (p < 0.001), yet GO incidence was 14% vs. 16%, respectively. This may be partly due to the concurrent and prolonged use of MMI in nearly 90% of the BCD group. In line with these findings, Azizi et al. [28] demonstrated that concurrent MMI administration reduces the risk of GO exacerbation in patients receiving RAI. Nevertheless, prospective studies are warranted to clarify whether extended antithyroid therapy truly mitigates GO risk in individuals presenting with high baseline TRAb.

The observed GO comorbidity rates of 11–17% are lower than the 27% reported by Tanda et al. [29] potentially reflecting differences in study design and patient selection. We required at least three years of TRAb follow-up, which may have excluded some patients with severe GO who lacked the requisite longitudinal data or were referred early for specialized care. Although no significant differences in GO comorbidities emerged among the four patterns, other investigations have found that higher TRAb levels correlate with more severe GO [10]. Future work examining GO severity rather than mere presence or absence may provide deeper insights.

Our results also corroborate previous evidence that female sex and lower initial TRAb levels predict a higher likelihood of TRAb normalization. Women have been reported to achieve higher remission rates than men [30], and an elevated baseline TRAb is a well-established marker for relapse. Cappelli et al. [24] similarly demonstrated that TRAb titers at GD diagnosis can help forecast long-term outcomes, where a high titer predicts persistent disease.

Despite identifying four distinct patterns of TRAb changes, clinicians should note that significant variability can exist within each cluster. Although these patterns offer a valuable framework for understanding general trajectories, individual patient courses may deviate substantially from the cluster average. Such heterogeneity reflects the multifactorial nature of GD, encompassing genetic predispositions, immune-regulatory differences, and variable treatment responses [31]. Accordingly, these patterns should be treated as broad guides rather than definitive predictions of each patient’s course, reinforcing the need for personalized approaches, possibly incorporating additional biomarkers and clinical variables.

Our study has some limitations. First, the retrospective design and single-center setting may introduce biases and underestimate GO comorbidities compared with previous reports, as not all newly diagnosed GD patients were referred for ophthalmological evaluation. Second, we did not classify disease states (relapse vs. remission) in our TRAb analysis, potentially overlooking important clinical nuances. Lastly, while linear interpolation is practical for handling missing data in a long-term dataset, it does not account for abrupt fluctuations that individual patients may experience. Nonetheless, this study provides valuable insights into the longitudinal behavior of TRAb and its clinical relevance in GD management.

In conclusion, we identified four distinct TRAb patterns in GD, which could be consolidated into two broader groups based on median survival times for normalization. Each pattern reflects varying TRAb trajectories and therapeutic intensities, offering clinicians a potential tool to anticipate disease courses and customize treatment strategies. Prospective, multicenter studies that incorporate a broader range of clinical, immunologic, and imaging parameters are necessary to validate these findings and further elucidate the interplay between TRAb patterns, GO severity, and long-term remission rates.