Detecting citrus diseases at an early stage is very important for ensuring fruit quality and minimizing production losses, as well as for raising awareness about sustainable agriculture. As a solution to this problem, the authors of this paper propose a hybrid feature-based citrus disease classification system that integrates deep learning representations, handcrafted descriptors, feature selection, and ensemble learning, all of which are tuned via Bayesian optimization. We perform tests on two real-world citrus disease datasets that differ greatly in nature: a four-class lemon dataset and a two-class orange dataset. Both datasets were collected under quite different environmental conditions, so they show diverse disease symptoms and varying background complexity. Deep semantic features were obtained by running a pretrained ResNet50 network. In addition to those, other complementary handcrafted features, such as color, texture, spatial, and statistical features, were extracted from the segmented infected areas. Together, the hybrid feature vector of 2082 dimensions was subjected to an optimization method known as Neighborhood Component Analysis (NCA). This technique selects 300 features that are most effective for discrimination while at the same time ensuring the preservation of class separability and the minimization of redundancy. For the classification task, two classifiers, namely Random Forest (RF) and Bayesian-Optimized Random Forest (BORF), were employed. The latter is based on Bayesian optimization to locate the model hyperparameters. To measure the model's performance in an unbiased manner, five-fold cross-validation was performed. Based on the experimental results, BORF can improve classification accuracy on the lemon dataset from 90.42% to 93.75% and on the orange dataset from 95.42% to 95.92% compared to the baseline RF classifier. Cross-validation mean accuracies of the proposed system were 93.75 ± 0.82% and 95.92 ± 0.47% for the lemon and orange datasets, respectively. Receiver Operating Characteristic (ROC) analysis provided class-specific area under the curve (AUC) values of 0.975 and 0.972 for the orange dataset, with a macro-averaged AUC of about 0.94 for the multi-class lemon dataset. The combination of hybrid feature fusion, NCA-based feature optimization, and Bayesian-optimized ensemble classification leads to enhanced discriminative power, greater robustness, and better generalization performance for the system in citrus disease identification in a real-world agricultural setting, as demonstrated by the results.