Objective. To determine the optimal approach for identifying and mitigating batch effects in PET/CT radiomics features, and further improve the prognosis of patients with head and neck cancer (HNC), this study investigated the performance of three batch harmonization methods. Approach. Unsupervised harmonization identified the batch labels by K-means clustering. Supervised harmonization regarding the image acquisition factors (center, manufacturer, scanner, filter kernel) as known/given batch labels, and Combat harmonization was then implemented separately and sequentially based on the batch labels, i.e. harmonizing features among batches determined by each factor individually or harmonizing features among batches determined by multiple factors successively. Extensive experiments were conducted to predict overall survival (OS) on public PET/CT datasets that contain 800 patients from 9 centers. Main results. In the external validation cohort, results show that compared to original models without harmonization, Combat harmonization would be beneficial in OS prediction with C-index of 0.687–0.740 versus 0.684–0.767. Supervised harmonization slightly outperformed unsupervised harmonization in all models (C-index: 0.692–0.767 versus 0.684–0.750). Separate harmonization outperformed sequential harmonization in CT_m+clinic and CT_cm+clinic models with C-index of 0.752 and 0.722, respectively, while sequential harmonization involved clinical features in PET_rs+clinic model further improving the performance and achieving the highest C-index of 0.767. Significance. Optimal batch determination especially sequential harmonization for Combat holds the potential to improve the prognostic power of radiomics model in multi-center HNC dataset with PET/CT imaging.