Local interaction patterns play an important role in the structural and functional organization of complex networks. Here we propose a method for obtaining parsimonious decompositions of networks into higher order interactions which can take the form of arbitrary motifs. The method is based on a class of analytically solvable generative models which in combination with non-parametric priors allow us to infer higher order interactions from dyadic graph data without any prior knowledge on the types or frequencies of such interactions. We test the presented approach on simulated data for which we recover the set of underlying higher order interactions to a high degree of accuracy. For empirical networks the method identifies concise sets of atomic subgraphs from within thousands of candidates that cover a large fraction of edges and include higher order interactions of known structural and functional significance. Being based on statistical inference the method also produces a fit of the network to analytically tractable higher order models opening new avenues for the systematic study of higher order interactions.