We have investigated the waveform distortion of energetic particle driven off-axis fishbone mode (OFM) in tokamak plasmas with kinetic magnetohydrodynamic (MHD) hybrid simulations. We extended our previous simulations (Li et al 2022 Nucl. Fusion62 026013) by considering higher-n harmonics in the MHD fluid, where n is toroidal mode number. The waveform distortion is successfully reproduced in the simulation for both magnetic fluctuations and temperature fluctuations. It is clarified that the waveform distortion arises from the superposition of the n = 2 harmonics on the fundamental n = 1 harmonics of OFM, where the n = 2 harmonics are generated by the MHD nonlinearity from the n = 1 OFM. Two types of waveform distortion can occur depending on the phase relationship between the n = 1 and n = 2 harmonics and the relative amplitude of the n = 2 harmonics to the n = 1 harmonics. Lissajous curve analyses indicate that the wave couplings between the n = 1 and n = 2 harmonics with phase-lock and lead to 'rising distortion' and 'falling distortion', respectively. The two types of waveform distortion can be attributed to the strong shearing profile of radial MHD velocity with n = 2 around the q = 2 magnetic flux surface. The dependence of waveform distortion on viscosity is investigated. It is found that the viscosity which is needed to reproduce the waveform distortion is larger than that in the experiment.