We characterize Gaussian estimates for transition probability of a discrete time Markov chain in terms of geometric properties of the underlying state space. In particular, we show that the following are equivalent: (1) Two sided Gaussian bounds on heat kernel (2) A scale invariant Parabolic Harnack inequality (3) Volume doubling property and a scale invariant Poincaré inequality. The underlying state space is a metric measure space, a setting that includes both manifolds and graphs as special cases. An important feature of our work is that our techniques are robust to small perturbations of the underlying space and the Markov kernel. In particular, we show the stability of the above properties under quasi-isometries. We discuss various applications and examples.