Three-dimensional (3D) genome organization, which determines how the DNA is packaged inside the nucleus, has emerged as a key regulatory mechanism of cellular processes. High-throughput chromosomal conformation capture (Hi-C) technologies have enabled the study of 3D genome organization by experimentally measuring interactions among genomic regions in 3D space. Changes in 3D genome organization have been associated with gene expression changes during cell fate-specification as well as with diseases such as cancer. Therefore, a key challenge in regulatory genomics is to systematically detect higher-order structural changes across Hi-C datasets from multiple conditions. Existing computational methods either do not model higher-order structural units or do not take into account complex relationships among the conditions of interest. We address these limitations with Tree-Guided Integrated Factorization (TGIF), a new multi-task Non-negative Matrix Factorization (NMF) approach. TGIF models arbitrary relationships among multiple Hi-C datasets as a tree such that closely related Hi-C datasets have similar lower-dimensional representation. TGIF provides a statistically significant set of differential TAD boundaries with higher precision than existing approaches. When applied to cardiomyocyte differentiation time-course data, TGIF identifies known and novel boundary elements specific to each stage of development. We also find that cardiovascular-disease-associated SNPs are depleted in TGIF boundaries, and differentially expressed genes are enriched near differential boundaries identified by TGIF. Finally, application to tissue-specific endothelial cell Hi-C datasets further demonstrates TGIF's ability to identify cell type-specific regulatory elements of interest. Taken together, TGIF is a tool to interrogate the dynamics of 3D genome organization and its impact on normal and disease phenotypes.