The author gives a systematic study of the Hardy spaces of functions with values in the noncommutative $L^p$-spaces associated with a semifinite von Neumann algebra $\mathcal{M}.$ This is motivated by matrix valued Harmonic Analysis (operator weighted norm inequalities, operator Hilbert transform), as well as by the recent development of noncommutative martingale inequalities. In this paper noncommutative Hardy spaces are defined by noncommutative Lusin integral function, and it is proved that they are equivalent to those defined by noncommutative Littlewood-Paley G-functions. The main results of this paper include:

(i) The analogue in the author's setting of the classical Fefferman duality theorem between $\mathcal{H}^1$ and $\mathrm{BMO}$.

(ii) The atomic decomposition of the author's noncommutative $\mathcal{H}^1.$

(iii) The equivalence between the norms of the noncommutative Hardy spaces and of the noncommutative $L^p$-spaces $(1 < p < \infty )$.

(iv) The noncommutative Hardy-Littlewood maximal inequality.

(v) A description of $\mathrm{BMO}$ as an intersection of two dyadic $\mathrm{BMO}$.

(vi) The interpolation results on these Hardy spaces.