We consider identification of linear dynamical systems comprising of high-dimensional signals, where the output noise components exhibit strong serial, and cross-sectional correlations. Although such settings occur in many modern applications, such dependency structure has not been fully incorporated in existing approaches in the literature. In this paper, we explicitly incorporate the dependency structure present in the output noise through lagged values of the observed multivariate signals. We formulate a constrained optimization problem to identify the space spanned by the latent states, and the transition matrices of the lagged values simultaneously, wherein the constraints reflect the low rank nature of the state information, and the sparsity of the transition matrices. We establish theoretical properties of the estimators, and introduce an easy-to-implement computational procedure for empirical applications. The performance of the proposed approach, and the implementation procedure is evaluated on synthetic data, and compared with competing approaches, and further illustrated on a data set involving weekly stock returns of 75 US large financial institutions for the 2001-2017 period.