Results on $\phi$ meson production in inelastic p+p collisions at CERN SPS energies are presented. They are derived from data collected by the NA61/SHINE fixed target experiment, by means of invariant mass spectra fits in the $\phi \to K^+K^-$ decay channel. They include the first ever measured double differential spectra of $\phi$ mesons as a function of rapidity $y$ and transverse momentum $p_T$ for proton beam momenta of 80 GeV/c and 158 GeV/c, as well as single differential spectra of $y$ or $p_T$ for beam momentum of 40 GeV/c. The corresponding total $\phi$ yields per inelastic p+p event are obtained. These results are compared with existing data on $\phi$ meson production in p+p collisions. The comparison shows consistency but superior accuracy of the present measurements. The emission of $\phi$ mesons in p+p reactions is confronted with that occurring in Pb+Pb collisions, and the experimental results are compared with model predictions. It appears that none of the considered models can properly describe all the experimental observables.