We develop a new 3‐D shear velocity model for the lithosphere and sublithospheric mantle under Saudi Arabia by jointly inverting P wave receiver functions, S wave receiver functions, and fundamental‐mode Rayleigh wave group velocities. P and S receiver functions are calculated from earthquakes recorded between 2012 and 2015 at 156 Saudi National Seismic Network stations operated by the Saudi Geological Survey. Rayleigh wave dispersion data are extracted from independent tomographic studies. Our model reveals significant lateral variations in crustal and upper‐mantle S velocity below Saudi Arabia. Particularly, a low‐velocity zone, with minimum S velocity of ~4.0 km/s in the depth range of 70‐190 km, is observed under the Arabian Shield coinciding with Cenozoic surface volcanism. The low‐velocity zone is found consistent with the presence of partial melts in the mantle and interpreted as a potential deep magma source for the volcanism in western Arabia. We propose that magmas responsible for the Arabian volcanism may be derived from multiple sources. Both lithospheric thinning and the Afar plume trigger magma production beneath southwestern Arabia, while decompression melting caused by the lithospheric thinning may be the main factor in the central and northern portions. Furthermore, we find evidence for localized crustal low shear velocity anomalies (2–4% reductions) that appear spatially correlated to the volcanism; these may be due to fractures caused by magma ascent and small amounts of partial melt in the crust. This spatial distribution of S velocity reductions may indicate the plumbing system for magma ascent underneath western Arabia.