Computational prediction of the pKa of ionizable groups remains a central challenge in biomolecular modeling. Although all-atom fixed-charge force fields could be accurate to describe the interaction network within the biomolecules, proper sampling techniques are required to obtain the thermodynamic information in the (de)protonation event. Sufficient sampling requires an ensemble of structures from simulations with proper treatments of the acid-base equilibria, and the grand canonical simulation technique could be used to model the growth/annihilation of hydrogen atoms by merging Hamiltonians of different protonation states into one simulation ensemble. The electrostatic feature of nucleotide systems is especially difficult to model, and the situation becomes more challenging when the ionizable site is highly perturbed. Although there are many successful predictions obtained from the grand canonical constant pH simulations, few reports focus on highly perturbed nucleotide systems with unconventional base-pair features. In this work, with the discrete constant pH method, we investigate the titration thermodynamics of an adenine in the catalytic triad in a 35-nucleotide single-stranded RNA hairpin, featuring an unconventional GA mismatch and a substantially shifted pKa value. Validation tests are performed with two system setups, both of which provide pKa predictions in good agreement with the experimental value. A single-configuration-based technique is used to calculate the pKa for comparison. The current success indicates the predictive power of the current nucleotide modeling framework.