In Magnetic Resonance Imaging (MRI), the sequential acquisition of raw complex-valued image data in Fourier space, also known as k-space, results in extended examination times. To speed up the MRI scans, k-space data are usually undersampled and processed using numerical techniques such as compressed sensing (CS). While the majority of CS-MRI algorithms primarily focus on magnitude images due to their significant diagnostic value, the phase components of complex-valued MRI images also hold substantial importance for clinical diagnosis, including neurodegenerative diseases. In this work, complex-valued MRI reconstruction is studied with a focus on the simultaneous reconstruction of both magnitude and phase images. The proposed algorithm is based on the nonsubsampled contourlet transform (NSCT) technique, which offers shift invariance in images. Instead of directly transforming the complex-valued image into the NSCT domain, we introduce a wavelet transform within the NSCT domain, reducing the size of the sparsity of coefficients. This two-level hierarchical constraint (HC) enforces sparse representation of complex-valued images for CS-MRI implementation. The proposed HC is seamlessly integrated into a proximal algorithm simultaneously. Additionally, to effectively minimize the artifacts caused by sub-sampling, thresholds related to different sub-bands in the HC are applied through an alternating optimization process. Experimental results show that the novel method outperforms existing CS-MRI techniques in phase-regularized complex-valued image reconstructions.