The Forward-backward semiclassical treatment of ensemble averaged quantities is combined with a discretized path integral description of the Boltzmann operator describing the initial density. We present a practical Monte Carlo methodology for calculating time-dependent expectation values and time correlation functions, applicable to polyatomic systems.