A fundamental issue in multiscale materials modeling and design is the consideration of traction-separation behavior at the interface, which generally influences the failure properties. This paper develops a physics-based machine learning model based on the deep material network (DMN) enriched by cohesive layers, which enables the accurate and efficient prediction of multiscale responses for heterogeneous materials with interfacial effect. New fitting parameters are invoked in the cohesive building block and have physical meanings related to the length scale and orientation of the cohesive layer. It is shown that the enriched material network can be effectively optimized via a multi-stage training strategy, with training data generated from linear elastic direct numerical simulation (DNS). The extrapolation capability of the method to unknown spaces is demonstrated through the debonding analysis of a unidirectional fiber-reinforced composite, where the interface behavior is governed by an irreversible softening mixed-mode cohesive law. Its predictive accuracy is validated against the nonlinear DNS results, and the reduction in computational time is particularly significant.