The model reduction problem for networks of interconnected dynamical systems is studied in this paper. In particular, networks of identical passive subsystems, which are coupled according to a tree topology, are considered. For such networked systems, reduction is performed by clustering subsystems that show similar behavior and subsequently aggregating their states, leading to a reduced-order networked system that allows for an insightful physical interpretation. The clusters are chosen on the basis of the analysis of controllability and observability properties of associated edge systems, representing the importance of the couplings and providing a measure of the similarity of the behavior of neighboring subsystems. This reduction procedure is shown to preserve synchronization properties (i.e., the convergence of the subsystem trajectories to each other) and allows for the a priori computation of a bound on the reduction error with respect to external inputs and outputs. The method is illustrated by means of an example of a thermal model of a building.