The inner workings of a biological cell or a chemical reactor can be rationalized by the network of reactions, whose structure reveals the most important functional mechanisms. For complex systems, these reaction networks are not known a priori and cannot be efficiently computed with ab initio methods; therefore, an important goal is to estimate effective reaction networks from observations, such as time series of the main species. Reaction networks estimated with standard machine learning techniques such as least-squares regression may fit the observations but will typically contain spurious reactions. Here we extend the sparse identification of nonlinear dynamics (SINDy) method to vector-valued ansatz functions, each describing a particular reaction process. The resulting sparse tensor regression method "reactive SINDy" is able to estimate a parsimonious reaction network. We illustrate that a gene regulation network can be correctly estimated from observed time series.