We present a data-driven approximation framework of a phase field model for brittle fracture. The proposed approach aims to alleviate the computational cost associated with full-scale simulations in studying brittle fracture in heterogeneous materials; where the fracture toughness of the material is spatially varying. Our strategy relies on the combination of dimensionality reduction and learning between function spaces. We first explore optimal strategies to encode and decode smooth and non-smooth physical fields. Neural network surrogates are subsequently introduced and trained to learn the solution map between ad hoc lower-dimensional latent spaces. We also explored a probabilistic framework - manifold learning, introduced to explore the geometrical structure exhibited by the latent variables and achieve further dimensionality reduction for reduced latent data. New sets of latent variables are sampled from manifold and subsequently lifted to corresponding full scale values. The accuracy of the methods proposed is finally demonstrated by comparing the predicted results and with those obtained from high-fidelity computations.