Pure gravity in three dimensions is conjecturally dual to an ensemble of 2d CFT data. This ensemble is highly constrained by modular invariance and crossing symmetry, which allows us to formulate a precise version of the generalized eigenstate thermalisation hypothesis (ETH) for 2d CFTs. In the dual gravitational description, higher statistical moments in this ensemble correspond to multiboundary Euclidean wormholes. I will show how to compute the exact partition functions of these wormholes using the recently proposed Virasoro TQFT. This leads to a precise match between non-Gaussian statistics of OPE coefficients and on-shell wormhole solutions in pure gravity. Based on WIP with Diego Liška and Jan de Boer (UvA).