We propose a local, non -intrusive model order reduction technique to accurately approximate the solution of coupled multi -component parametrized systems governed by partial differential equations. Our approach is based on the approximation of the boundary response maps arising from a non -overlapping domain decomposition method. To construct the surrogate, we combine dimensionality reduction techniques with interpolation or regression approaches, such as kernel interpolation methods and artificial neural networks. Two alternative training strategies, making use of the full coupled problem or an artificial parametrization of the boundary conditions, are proposed and discussed. We show the potential of our approach in a series of test cases, ranging from linear diffusion -like models to nonlinear multi -physics problems. High levels of accuracy and computational efficiency are achieved in all cases.