  We present an algorithm for the stochastic simulation of gene expression and
heterogeneous population dynamics. The algorithm combines an exact method to
simulate molecular-level fluctuations in single cells and a constant-number
Monte Carlo method to simulate time-dependent statistical characteristics of
growing cell populations. To benchmark performance, we compare simulation
results with steadystate and time-dependent analytical solutions for several
scenarios, including steadystate and time-dependent gene expression, and the
effects on population heterogeneity of cell growth, division, and DNA
replication. This comparison demonstrates that the algorithm provides an
efficient and accurate approach to simulate how complex biological features
influence gene expression. We also use the algorithm to model gene expression
dynamics within "bet-hedging" cell populations during their adaption to
environmental stress. These simulations indicate that the algorithm provides a
framework suitable for simulating and analyzing realistic models of
heterogeneous population dynamics combining molecular-level stochastic reaction
kinetics, relevant physiological details and phenotypic variability.
