Accurate classification of skin diseases using deep learning models is crucial in medical diagnostics, significantly impacting patient care. This study addresses the challenges of class imbalance and uneven data distribution in the HAM10000 dataset by proposing a novel approach that integrates advanced data augmentation, model customization, and statistical validation techniques. We employed comprehensive data preprocessing and augmentation strategies to balance underrepresented skin disease classes, creating a more equitable training set. The EfficientNetB0 architecture was customized specifically for skin disease classification, with key layers fine-tuned through Bayesian Optimization, which optimized hyperparameters such as the number of trainable layers and the handling of class weights, boosting classification accuracy. To ensure clinically relevant predictions, we evaluated the model using top-k accuracy metrics (top-1, top-2, and top-3), achieving over 99% top-3 accuracy, thus providing valuable shortlists of potential diagnoses to assist healthcare professionals. Furthermore, we introduced a stratified data splitting method based on clustering PCA-reduced features to ensure similar distributions across training and validation sets, reducing the risk of overfitting. Finally, the Kolmogorov-Smirnov test was used to statistically validate the similarity between training and validation distributions, ensuring the model’s robustness and minimizing biases introduced by uneven data splits. These results demonstrate the model’s reliability and clinical utility, offering precise and actionable predictions in real-world medical contexts.