In this study, comparative analysis of calculated (GIAO method, DFT level) and experimental P-31 NMR shifts for a wide range of model palladium complexes showed that, on the whole, the theory reproduces the experimental data well. The exceptions are the complexes with the P=O phosphorus, for which there is a systematic underestimation of shielding, the value of which depends on the flexibility of the basis sets, especially at the geometry optimization stage. The use of triple-zeta quality basis sets and additional polarization functions at this stage reduces the underestimation of shielding for such phosphorus atoms. To summarize, in practice, for the rapid assessment of P-31 NMR shifts, with the exception of the P=O type, a simple PBE0/{6-311G(2d,2p); Pd(SDD)}//PBE0/{6-31+G(d); Pd(SDD)} approximation is quite acceptable (RMSE = 8.9 ppm). Optimal, from the point of view of "price-quality" ratio, is the PBE0/{6-311G(2d,2p); Pd(SDD)}//PBE0/{6-311+G(2d); Pd(SDD)} (RMSE = 8.0 ppm) and the PBE0/{def2-TZVP; Pd(SDD)}//PBE0/{6-311+G(2d); Pd(SDD)} (RMSE = 6.9 ppm) approaches. In all cases, a linear scaling procedure is necessary to minimize systematic errors.