  We propose a method to incorporate the coupling between shape and pairing
collective degrees of freedom in the framework of the interacting boson model
(IBM), based on the nuclear density functional theory. To account for pairing
vibrations, a boson-number non-conserving IBM Hamiltonian is introduced. The
Hamiltonian is constructed by using solutions of self-consistent mean-field
calculations based on a universal energy density functional and pairing force,
with constraints on the axially-symmetric quadrupole and pairing intrinsic
deformations. By mapping the resulting quadrupole-pairing potential energy
surface onto the expectation value of the bosonic Hamiltonian in the boson
condensate state, the strength parameters of the boson Hamiltonian are
determined. An illustrative calculation is performed for $^{122}$Xe, and the
method is further explored in a more systematic study of rare-earth $N=92$
isotones. The inclusion of the dynamical pairing degree of freedom
significantly lowers the energies of bands based on excited $0^+$ states. The
results are in quantitative agreement with spectroscopic data, and are
consistent with those obtained using the collective Hamiltonian approach.
