Print
facebook
twitter
Linkedin
Introduction
Thyroid cancer is the most common malignancy of the endocrine system, and its incidence continues to rise worldwide [1]. Differentiated thyroid carcinoma (DTC) represents more than 90% of thyroid cancers and generally carries a favorable prognosis. However, a subset of DTCs progresses to locally advanced disease with increased aggressiveness and poor clinical outcomes [2-4]. This progression follows a genetic multistep process: early-stage DTCs harbor initiating driver mutations, while the subsequent accumulation of secondary aggressive mutations promotes invasion and dedifferentiation [5-7]. Therefore, accurate risk stratification requires assessing both the initial driver mutations and these secondary alterations in order to optimize individualized management.
Most thyroid cancers originate from mutually exclusive initiating mutations in BRAF, RAS, or RET rearrangements [8,9]. The BRAFV600E mutation has long been regarded as a key biomarker for classification and therapeutic guidance [6,10,11]. However, in Asian populations, BRAF mutations are highly prevalent (50%-90%), even in micro-papillary thyroid carcinomas (PTCs), yet only a minority of these cases develop adverse outcomes [12,13]. This discrepancy indicates that single initiating mutations (e.g., BRAF alone) are insufficient to predict aggressiveness; rather, the interaction between initiating drivers and secondary mutations may determine clinical behavior.
Growing evidence shows that additional alterations, such as mutations in the TERT promoter, TP53, EIF1AX, and genes in the phosphoinositide 3-kinase (PI3K)/AKT/mammalian target of rapamycin (mTOR) pathway, can enhance invasiveness and reduce differentiation, leading to poorer prognosis [14-18]. In particular, BRAFV600E combined with TERT promoter mutations is linked to iodine resistance, higher recurrence, and increased mortality [19,20]. Pappa et al. [21] further demonstrated that BRAF mutations in combination with PIK3/AKT/mTOR pathway alterations are associated with more aggressive clinicopathologic features and greater progression risk compared to BRAF mutation alone. These findings highlight the importance of identifying multigenic co-mutations (initiating+secondary) for accurate prognostication and treatment planning.
Nevertheless, comprehensive next generation sequencing (NGS) panels to detect co-mutations remain inaccessible in non-tertiary centers (e.g., community hospitals) and cost-prohibitive in resource-limited regions. Given the low incidence of multigenic co-mutations, universal NGS testing for all DTC patients not only wastes medical resources but also exposes low-risk individuals to unnecessary anxiety and additional follow-up procedures [14,21]. Therefore, a cost-effective pre-screening method is needed to identify high-risk individuals who would benefit most from NGS, thereby optimizing resource allocation and personalized care. Contrast-enhanced ultrasonography (CEUS) offers a promising non-invasive alternative: it enables dynamic, real-time evaluation of tumor microvascular features, which cannot be assessed with grayscale ultrasonography (US) alone. This capability supports more precise lesion characterization and has demonstrated value in diagnosing thyroid cancer and predicting prognosis [22-25]. Although US characteristics have been applied to predict single-gene mutations (e.g., BRAF, TERT promoter) [26-28], the CEUS features associated with multigenic co-mutations—especially those involving key secondary mutations such as TERT promoter and the PIK3/AKT/mTOR pathway—remain unclear.
This study aimed to construct a non-invasive predictive nomogram for multigenic co-mutations in DTC using CEUS imaging, thereby exploring the role of CEUS as a pre-screening tool to refine patient selection for NGS testing and to facilitate targeted genetic evaluation and optimized clinical management.
Materials and Methods
Compliance with Ethical Standards
This retrospective observational study was approved by the Institutional Review Board (grant number SYSKY-2024-289-01). The requirement for informed consent was waived because no personal health information was disclosed.
Patients
From September 2021 to December 2023, patients were included if they met the following criteria: (1) NGS was performed for suspected thyroid carcinoma, and (2) preoperative grayscale US and CEUS were available. Exclusion criteria were as follows: (1) failed NGS due to inadequate samples, (2) pathology confirmed as non-DTC, and (3) non-primary or non–treatment-naïve disease. The sample size was determined by available cases. Ultimately, 116 consecutive patients were enrolled (Fig. 1).
US Examination
All grayscale US and CEUS examinations were performed by experienced US physicians using various instruments, including Mindray Resona, Siemens Sequoia, and Samsung systems. For patients with multiple lesions, each suspicious nodule was examined, but only the lesion that underwent fine-needle aspiration and NGS was analyzed in this study. Patients were positioned supine with maximum neck extension. After routine grayscale US and color Doppler evaluation, a bolus of 1.2-2.4 mL of sulfur hexafluoride microbubble suspension (SonoVue, Bracco, Milano, Italy) was injected intravenously, followed by 5.0 mL of 0.9% saline. Lesions were observed until contrast wash-out, during which patients were instructed to breathe steadily and avoid swallowing or speaking. All grayscale US and CEUS images were digitally stored on a cloud platform for subsequent analysis.
US Imaging Analysis
All grayscale US and CEUS images were independently reviewed by two certified sonographers blinded to participants’ clinicopathologic features and NGS results. In cases of disagreement, consensus was reached through discussion. Extracted imaging features included tumor size, shape, margin, microcalcifications, extrathyroidal extension, contrast agent arrival time, enhancement direction, enhancement pattern, peak intensity, and rim enhancement. Definitions of these features are provided in Supplementary Table 1. Contrast agent arrival time was assessed relative to surrounding normal thyroid parenchyma: "earlier" indicated enhancement preceding normal tissue, while "synchronous/later" denoted simultaneous or delayed enhancement. Enhancement homogeneity was classified as "homogeneous" for uniform uptake and "heterogeneous" for patchy, non-uniform uptake.
Additional Mutations
NGS of surgical and biopsy specimens was performed using the Illumina NovaSeq 6000 or NextSeq 550 platforms with a thyroid cancer–specific 88-gene panel. For this study, additional mutations were defined as secondary alterations with established associations with aggressive phenotypes and relatively high incidence in DTC (based on prior literature [15,21]) that were included in the panel. These comprised: (1) TERT promoter mutations (c.-124C>T and c. -146C>T) and EIF1AX A113 splice mutation; (2) key mutations in the PI3K/AKT/mTOR pathway (AKT1 p.E17K, PIK3CA p.E542K, p.E545K, p.H1047R/L); and (3) loss-of-function mutations in the tumor suppressor TP53. Because of the low frequency of individual secondary mutations, a multigenic co-mutation was defined as the presence of a primary driver mutation together with any of these additional alterations.
Study Design
The analytical workflow consisted of the following steps. First, clinical and US features were compared between groups to identify potential predictive variables. Second, collinearity analysis was performed to exclude highly correlated variables and minimize multicollinearity. Third, a Bayesian logistic regression model with 6,000 iterations (including 3,000 burn-in iterations) was applied to identify independent predictors and construct the final prediction model. Given the absence of prior evidence linking clinical and CEUS features to multigenic co-mutations in DTC, non-informative normal priors (mean, 0; standard deviation, 10) were used to avoid subjective bias. Independent variable selection was based on posterior inclusion probabilities (PIP >0.9). Model convergence was assessed using Markov-chain Monte Carlo (MCMC) diagnostics, and posterior predictive checks were conducted to evaluate the model’s ability to replicate observed data. Internal validation was then performed using bootstrap resampling (1,000 replicates) to test model stability and generalizability. Finally, the optimal cutoff probability was determined using the Youden index from the receiver operating characteristic (ROC) curve, and confusion matrix analysis at this cutoff was used to assess clinical screening performance.
Statistical Analysis
Continuous variables were compared using the Student t-test for normally distributed data or the Mann-Whitney U test for non-normally distributed data. Categorical variables were analyzed using the chi-square test or Fisher exact test. Multicollinearity among candidate predictors was assessed using correlation coefficients and variance inflation factors (VIFs). Independent predictors of multigenic co-mutations were identified using Bayesian logistic regression, with convergence assessed by MCMC diagnostics, including trace plot stationarity and potential scale reduction factor (PSRF <1.01 across four independent chains). Model performance was evaluated using (1) posterior predictive checking with coverage (≥0.85, proportion of observed values within 95% predictive intervals), (2) ROC curves, and (3) confusion matrix analysis to obtain sensitivity, specificity, and accuracy. For Bayesian inference, regression coefficients were considered meaningful if the 95% confidence interval (CI) excluded zero. For frequentist analyses, two-tailed tests were applied, with P<0.05 considered significant. Statistical analyses were performed with SPSS Statistics version 25.0 (IBM Corp., Armonk, NY, USA) and R version 4.3.0 (R Foundation for Statistical Computing, Vienna, Austria) using the brms package (for Bayesian modeling) and pROC package (for ROC analysis).
Results
Baseline Characteristics
As shown in Table 1, a total of 116 patients (83 women and 33 men) were included, with a mean age of 39.84 years and a median tumor size of 11.00 mm. The majority were diagnosed with PTC (111 cases, 96.5%), while five patients (4.3%) had follicular thyroid carcinoma (FTC). Lymph node metastasis was identified in 64 patients (55.2%), with only four cases at advanced TNM stage (Ⅲ -Ⅳ). With respect to genetic alterations, BRAF mutation was the most frequent driver mutation (90 cases, 77.6%), followed by RAS mutation (6 cases, 5.2%) and RET/PTC rearrangement (4 cases, 3.5%). Additional mutations included four cases each of TP53 and TERT promoter mutations, three cases of PIK3CA mutation, and one case of AKT1 mutation. Overall, 12 patients (10.3%) were identified as positive for multigenic co-mutations. The distribution of genetic alterations across all patients is presented in Supplementary Fig. 1.
Clinical Characteristics of DTC with Multigenic Co-mutations
The distribution of clinical and pathological characteristics between patients with and without multigenic co-mutations is summarized in Table 2. No significant differences in demographic variables were observed between the two groups. However, a higher proportion of aggressive histological subtypes, including the tall cell and columnar cell variants of PTC and the widely invasive variant of FTC, was found in the multigenic co-mutation group (25.0% vs. 1.9%, P=0.008). Although the incidence of capsular invasion did not differ significantly between groups, lymph node metastasis occurred more frequently in patients with multigenic co-mutations (83.3% vs. 51.9%, P=0.039). Furthermore, this group exhibited a higher proportion of advanced TNM stages (16.7% vs. 1.9%, P=0.053), though the difference did not reach statistical significance.
Bayesian Predictive Model for Multigenic Co-mutations
The grayscale US and CEUS features of DTC with multigenic co-mutations are presented in Table 3. Differential analyses identified four features with significant differences between groups (P<0.05): echogenicity, tumor size, enhancement homogeneity, and contrast agent arrival time. As shown in Supplementary Fig. 2A and B, multicollinearity was minimal, with all pairwise correlation coefficients <0.8 and VIF values well below 10.
Bayesian regression identified three key predictors: tumor size, enhancement homogeneity, and contrast agent arrival time, each with PIP >0.9. These predictors demonstrated stable posterior distributions (Supplementary Fig. 2C). MCMC trace plots indicated good convergence across chains, with all PSRF values <1.01 (Supplementary Fig. 2D). The model incorporating these three predictors achieved an area under the ROC curve (AUC) of 0.873 (Fig. 2A). Posterior predictive checks confirmed strong concordance with observed data, with the 95% predictive interval covering more than 85% of observed values (Fig. 2B). Internal validation using bootstrap resampling demonstrated consistent performance, with a mean AUC of 0.880 (95% CI, 0.745 to 0.978) (Supplementary Fig. 3).
The optimal cutoff probability was determined as 0.174 using the Youden index. Confusion matrix analysis at this cutoff produced the following outcomes: 9 true positives, 94 true negatives, 10 false positives, and 3 false negatives (Fig. 2C). The model achieved an accuracy of 0.888, sensitivity (recall) of 0.750, specificity of 0.904, precision of 0.473, and an F1 score of 0.581.
Representative CEUS images of DTC with and without multigenic co-mutations are presented in Fig. 3.
Discussion
In this study, multigenic co-mutations were identified in 10.3% of patients with DTC and were associated with more aggressive clinicopathological features. Using three key US features (tumor size, enhancement homogeneity, and contrast agent arrival time) this study developed a Bayesian predictive model that achieved an AUC of 0.873. Internal validation by bootstrap resampling demonstrated a mean AUC of 0.880 (95% CI, 0.745 to 0.978), with accuracy of 0.888, sensitivity of 0.750, and specificity of 0.904. These findings suggest that US-based features can serve as effective non-invasive predictors of multigenic co-mutations in DTC, supporting personalized guidance for NGS panel selection and clinical management.
In this study, BRAF mutations were present in 90 patients (77.6%), while 12 patients (10.3%) harbored multigenic co-mutations, defined as primary driver mutations accompanied by secondary alterations. These included TERT promoter mutations (4 cases, 3.5%), TP53 mutations (4 cases, 3.5%), PIK3CA mutations (3 cases, 2.6%), and AKT1 mutation (1 case, 0.9%). Such genetic patterns may be influenced by ethnicity [8,29,30]. Specifically, BRAF mutation rates in Asian thyroid cancer populations range from 50% to 90%; this is consistent with the present study’s findings and significantly higher than in Western populations [8,13]. The characteristics of this cohort therefore align with the genetic profile commonly observed in Asian populations.
Beyond genetic prevalence, the clinical implications of these co-mutations are consistent with established mechanisms of thyroid cancer progression. Previous studies show that well-differentiated thyroid cancers may acquire aggressive phenotypes through key events such as PI3K/AKT pathway activation (e.g., PIK3CA), TP53 loss-of-function, and TERT reactivation, all of which were observed in the present co-mutation cohort [14]. Consistent with this, patients with multigenic co-mutations demonstrated a higher prevalence of aggressive histological subtypes (25.0% vs. 1.9%, P=0.008) and lymph node metastasis (83.3% vs. 51.9%, P=0.039), echoing earlier reports [14,15]. Although the co-mutation group also had a greater proportion of advanced TNM stages (16.7% vs. 1.9%), the borderline P-value (P=0.053) may be attributable to the limited number of events (only 2 advanced cases in the co-mutation group and 4 overall), which likely reduced statistical power. Larger cohorts, particularly with more advanced-stage cases, may help to confirm this trend.
Tumor size, contrast agent arrival time, and enhancement homogeneity were identified as independent predictors of multigenic co-mutations in DTC using Bayesian regression, with each predictor showing high inclusion probabilities (>0.9) and stable posterior distributions. Prior research has validated US as a valuable tool for preoperative prediction of thyroid cancer molecular alterations. For instance, BRAF mutations have been associated with features such as microcalcifications, irregular margins, and hypo-enhancement [26,31,32], while TERT promoter mutations have been linked to age, non-parallel orientation, microlobulated margins, multifocality, and capsule involvement [28,33]. Unlike these studies, which largely focused on single-gene alterations, the present research extends this approach to multigenic co-mutations, demonstrating the feasibility of preoperative US-based prediction for this aggressive genetic subset.
Because the number of patients with multigenic co-mutations was small (n=12), Bayesian regression was adopted for model development. Compared with conventional logistic regression, Bayesian methods are more robust in limited-sample contexts and allow probabilistic interpretation of predictor importance. This methodological rigor increases confidence in the findings despite modest sample size. The Bayesian model demonstrated strong performance, with an AUC of 0.873 in the primary analysis and consistent results in bootstrap resampling (AUC, 0.880; 95% CI, 0.745 to 0.978). Reliability was further supported by MCMC trace plots and posterior predictive checks. While the sensitivity of 0.750 reflects reasonable ability to detect true positives, the moderate precision (0.473) indicates challenges with false-positives—likely influenced by the low prevalence of multigenic co-mutations (10.3%) and resulting class imbalance. Adjusting the optimal threshold (0.174 by the Youden index) to increase sensitivity reduced specificity, reflecting the inherent trade-off in rare-event prediction. This limitation may be addressed in future studies using larger, better-balanced cohorts.
With the growing emphasis on comprehensive genetic evaluation in DTC, this Bayesian model demonstrates potential clinical value in streamlining preoperative risk stratification for NGS testing, particularly as a screening tool to exclude low-risk patients and avoid unnecessary referrals. By identifying individuals at high risk of multigenic co-mutations (based on the 0.174 threshold), it helps prioritize NGS testing for those most likely to benefit from comprehensive genetic assessment, while effectively sparing patients with very low predicted risk from unnecessary procedures. In this way, the model supports more precise risk stratification, enables personalized management of DTC, and promotes more efficient allocation of medical resources.
This study has several limitations. First, its retrospective, single-center design with a relatively small sample size may have introduced bias, and the absence of external validation limits the credibility and generalizability of the findings. Future studies should therefore expand sample sizes and involve multi-center collaborations with more diverse populations to mitigate these constraints. Second, the limited sample size prevented us from incorporating rarer co-mutations (such as the EIF1AX A113 splice mutation) into the analysis, underscoring the need for larger cohorts to achieve a more comprehensive genetic assessment. Finally, the CEUS parameters included in the model are operator-dependent, which may introduce subjective variability. To address this, future work should assess interobserver and intraobserver variability in larger patient cohorts and further explore the role of quantitative parameters derived from multi-modal US (including grayscale US, elastography, and CEUS) in predicting multigenic co-mutations in DTC.
In conclusion, DTC patients with multigenic co-mutations exhibited greater aggressiveness compared with those without. The proposed preoperative CEUS-based Bayesian model demonstrated preliminary potential for identifying these co-mutations; however, external validation is essential before its clinical applicability can be confirmed.
Download Citation
PDF Links
PubReader
ePub Link
XML Download
Full text via DOI
Full text via PMC
