Tissue mechanics have been widely studied for the past two decades because they are significantly associated with disease states. The variations of mechanical properties in soft tissues are considerable indicators used for clinical diagnosis and disease monitoring. Optical coherence elastography (OCE) has been extensively developed to characterize the mechanical properties of various tissues. However, these methods are generally based on time-domain data and measure the time-of-flight of the localized shear wave propagations. Here, we propose a method to evaluate phase velocity in four-dimensional (4D) space (x, y, z, f), called 4D-OCE phase velocity, and applied it to a heterogeneous phantom with a 4 mm diameter inclusion. A 7.5 MHz highly focused transducer transmitting a toneburst of 500 µs was used to provide an acoustic radiation force (ARF) for transient elastography. A customized 4D acquisition pattern was developed under Microsoft Visual C++ 2019 development environment. The 3D-OCE group velocity was reconstructed to compare with 4D-OCE phase velocity as cross-verification. The experimental result from the heterogeneous phantom showed that 4D-OCE phase velocity is capable of accurately reconstructing a 4D-OCE phase velocity map and provides good contrast between two gelatin materials with different concentrations, 7% versus 14%. The proposed 4D-OCE phase velocity can provide further information to better understand the delicate change of mechanical properties of various biological tissues in each voxel over a range of temporal frequencies.