Subspace-based methods have been effectively used to estimate multi-input/multi-output, discrete-time, linear-time invariant systems from uniformly spaced spectrum samples. A critical step in these methods is the splitting of causal and noncausal invariant subspaces of a Hankel matrix built from spectrum measurements via singular-value decomposition in order to determine model order. Quite often, in particular when signal-to-noise ratio is low, unmodelled dynamics is present, and the number of measurements is small, this step is inconclusive since the assumed mirror image symmetry with respect to the unit circle between the eigenvalues of the invariant spaces is lost. In this paper, we propose a new model order selection and invariant space splitting scheme based on the regularized nuclear norm optimization in combination with a particular subspace method that is insensitive to noise and undermodelling. By a detailed numerical study, efficacy of the proposed scheme is shown for a broad range of signal-to-noise ratio and short data records.