In recent years there has been a considerable drive towards data-driven analysis, discovery and control of dynamical systems. To this end, operator theoretic methods, namely, Koopman operator methods have gained a lot of interest. In general, the Koopman operator is obtained as a solution to a least-squares problem, and as such, the Koopman operator can be expressed as a closed-form solution that involves the computation of Moore-Penrose inverse of a matrix. For high dimensional systems and also if the size of the obtained data-set is large, the computation of the Moore-Penrose inverse becomes computationally challenging. In this paper, we provide an algorithm for computing the Koopman operator for high dimensional systems in a time-efficient manner. We further demonstrate the efficacy of the proposed approach on two different systems, namely a network of coupled oscillators (with state-space dimension up to 2500) and IEEE 68 bus system (with state-space dimension 204 and up to 24,000 time-points).