Maximally localized Wannier functions (MLWFs) are widely used to construct first-principles tight-binding models that accurately reproduce the electronic structure of materials. Recently, robust and automated approaches to generate these MLWFs have emerged, leading to natural sets of atomic-like orbitals that describe both the occupied states and the lowest lying unoccupied ones (when the latter can be meaningfully described by bonding/anti-bonding combinations of localized orbitals). For many applications, it is important to instead have MLWFs that describe only certain target manifolds separated in energy between them — the occupied states, the empty states, or certain groups of bands. Here, we start from the full set of MLWFs describing simultaneously all the target manifolds, and then mix them using a combination of parallel transport and maximal localization to construct orthogonal sets of MLWFs that fully and only span the desired target submanifolds. The algorithm is simple and robust, and it is applied to some paradigmatic but non-trivial cases (the valence and conduction bands of silicon, the top valence band of MoS₂, the 3d and t2g/eg bands of SrVO₃ and to a mid-throughput study of 77 insulators.