The theoretical approach we outlined here can be widely applicable to any mixed quantum-classical methods to simulate polariton dynamics (population, coherence, different types of correlation functions). The scaling reduction in the single excitation subspace is the most relevant part of the linear spectra simulation in this paper, whereas the scaling reduction in the double excitation subspace will be relevant for 2DES simulations in Paper II.