Dwarf spheroidal galaxies are among the most promising targets for indirect dark matter (DM) searches in $\gamma$-rays. The $\gamma$-ray flux from DM annihilation in a dwarf spheroidal galaxy is proportional to the $J$-factor of the source. The $J$-factor of a dwarf spheroidal galaxy is the line-of-sight integral of the DM mass density squared times $\langle \sigma_{\rm ann} v_{\rm rel} \rangle/(\sigma_{\rm ann} v_{\rm rel})_0$, where $\sigma_{\rm ann} v_{\rm rel}$ is the DM annihilation cross-section times relative velocity $v_{\rm rel}=|{\bf v}_{\rm rel}|$, angle brackets denote average over ${\bf v}_{\rm rel}$, and $(\sigma_{\rm ann} v_{\rm rel})_0$ is the $v_{\rm rel}$-independent part of $\sigma_{\rm ann} v_{\rm rel}$. If $\sigma_{\rm ann} v_{\rm rel}$ is constant in $v_{\rm rel}$, $J$-factors only depend on the DM space distribution in the source. However, if $\sigma_{\rm ann} v_{\rm rel}$ varies with $v_{\rm rel}$, as in the presence of DM self-interactions, $J$-factors also depend on the DM velocity distribution, and on the strength and range of the DM self-interaction. Models for self-interacting DM are increasingly important in the study of the small scale clustering of DM, and are compatible with current cosmological observations. Here we derive the $J$-factor of 20 dwarf spheroidal galaxies from stellar kinematic data under the assumption of Yukawa DM self-interactions. $J$-factors are derived through a profile Likelihood approach, assuming either NFW or cored DM profiles. We also compare our results with $J$-factors derived assuming the same velocity for all DM particles in the target galaxy. We find that this common approximation overestimates the $J$-factors by up to one order of magnitude. $J$-factors for a sample of DM particle masses, self-interaction coupling constants and density profiles are provided electronically, ready to be used in other projects.